Mass Exchange Dynamics of Surface and Subsurface Oil in Shallow-Water Transport

TL;DR

This paper develops a large-scale model for oil spill transport in shallow water, focusing on the mass exchange between surface and subsurface oil driven by droplet dynamics and buoyancy, addressing microscale challenges.

Contribution

It introduces an upscaled, environmentally relevant model for surface-subsurface oil exchange that maintains mass conservation and is applicable at large spatio-temporal scales.

Findings

Model captures oil exchange mechanisms due to droplet dynamics and buoyancy.

Addresses microscale challenges through upscaling for large-scale predictions.

Ensures mass conservation in shallow-water oil spill simulations.

Abstract

We formulate a model for the mass exchange between oil at and below the sea surface. This is a particularly important aspect of modeling oil spills. Surface and subsurface oil have different chemical and transport characteristics and lumping them together would compromise the accuracy of the resulting model. Without observational or computational constraints, it is thus not possible to quantitatively predict oil spills based upon partial field observations of surface and/or sub-surface oil. The primary challenge in capturing the mass exchange is that the principal mechanisms are on the microscale. This is a serious barrier to developing practical models for oil spills that are capable of addressing questions regarding the fate of oil at the large spatio-temporal scales, as demanded by environmental questions. We use upscaling to propose an environmental-scale model which incorporates the mass exchange between surface and subsurface oil due to oil droplet dynamics, buoyancy effects, and sea surface and subsurface mechanics. While the mass exchange mechanism detailed here is generally applicable to oil transport models, it addresses the modeling needs of a particular to an oil spill model [1]. This transport model is designed to capture oil spills at very large spatio-temporal scales. It accomplishes this goal by specializing to shallow-water environments, in which depth averaging is a perfectly good approximation for the flow, while at the same time retaining mass conservation of oil over the whole oceanic domain.

Equations

$$S=\sum _{i=1}^N{s}_i,$$

$$M_s(T)={\int }_{{\Omega }_T}\sum _{i=1}^N{s}_i(\mathbf{x},T){\rho }_{i}d\mathbf{x},$$

$$\frac{\partial s_i}{\partial T}+{\nabla }_{\perp }\cdot \left(C_{xs}\mathbf{V}{s}_i+\frac{\mathbf{\tau }}{2{\mu }_i}{s}_i^2\right)-{\nabla }_{\perp }\cdot [\Psi {\nabla }_{\perp }{s}_i]=R_i+E_i^s(s_i)+P_i^s+G_i^s.$$

$$M_c={\int }_{{\Omega }_T}{\int }_{-h}^{\eta }\sum _{i=1}^N{c}_i(\mathbf{x},z,T)dzd\mathbf{x}.$$

$$C_i(\mathbf{x},T)=\frac{1}{H(\mathbf{x})}{\int }_{-h}^{\eta (T)}{c}_i(\mathbf{x},z,T)dz$$

$$\frac{\partial H{C}_i}{\partial T}+{\nabla }_{\perp }\cdot (H\mathbf{V}{C}_i)-{\nabla }_{\perp }\cdot [H\Psi {\nabla }_{\perp }{C}_i]=-{\rho }_i{R}_i+{\rho }_i\left[E_i^C(C_i)+P_i^C+G_i^C\right].$$

$$M(\mathbf{x},T)$$

$$={\int }_{-h}^{\eta (T)}c(\mathbf{x},z,T)dz+\rho s(\mathbf{x},T)$$

$$={\int }_{-h}^{\eta (T)}{\int }_0^{\infty }\rho V(r)n(r,\mathbf{x},z,T)drdz+\rho s(\mathbf{x},T)$$

$$={\int }_0^{\infty }\rho V(r)\bar{n}(r,\mathbf{x},T)dr+\rho s(\mathbf{x},T)$$

$$\frac{\partial HC}{\partial T}+\mathbf{V}\cdot \nabla (HC)+\rho R=-\rho \left[\frac{\partial s}{\partial T}+\mathbf{V}\cdot \nabla s-R\right]$$

$$R=-{\int }_0^{\infty }V(r)w_b(r)n(r,\mathbf{x},\eta ,T)dr+{\int }_{-h}^{\eta }{\int }_0^{\infty }V(r)S_{b}drdz.$$

$$b(T)=\sum _{i=1}^{N_b(T)}{b}_i{\delta }_{T_i}(T)$$

$$N(r):=\frac{3(4-s)}{4\pi (r_m^{4-s}-r_1^{4-s})}{r}^{-s},r\in [r_1,r_m],$$

$${\int }_0^{\infty }{g}_b(z)dz=1.$$

$$S_b(T_i,z,r)=\frac{1}{dT}N(r)g_b(z)b_{i}s(T_i).$$

$$R(\mathbf{x},T)=-{\int }_0^{\infty }V(r)w_b(r)n(r,\mathbf{x},\eta ,T)dr+\sum _{i=1}^{N_b(T)}s(T_i)b_i{\delta }_{T_i}(T).$$

$$\frac{\partial n}{\partial T}=D_z(r)\frac{{\partial }^{2}n}{\partial z^2}-(W(\mathbf{x},z,T)+w_b(r))\frac{\partial n}{\partial z}+S_b+\mathcal{K}[n],$$

$${\int }_0^{\infty }\mathcal{K}[n](r)dr=0.$$

$$\frac{\partial n}{\partial z}(t,\eta ,r)=0,t>0,r\in [0,\infty ).$$

$$\frac{\partial n}{\partial T}$$

$$\frac{\partial s}{\partial T}$$

$$n(r,z,T)$$

$$n(0,z,r)=\frac{1}{\rho }HC(0)N(r)g_b(z)$$

$$G(r,z,z^{\prime },T)=$$

$$\times \exp \left[-\frac{w_b(2(z^{\prime }-z)+T{w}_b)}{4D_z}\right]$$

$$-\frac{1}{2D_z}{w}_b\exp \left[\frac{(\eta -z^{\prime })w_b}{D_z}\right]\text{erfc}\left[\frac{2\eta -(z+z^{\prime })+T{w}_b}{2\sqrt{T{D}_z}}\right].$$

$$s(T)$$

$$F_0(T)$$

$$\bar{n}(r,T)={\int }_{-h}^{\eta (T)}n(r,z,T)dz,$$

$$\bar{n}(r,T)$$

$$\bar{G}(r,z^{\prime },T)={\int }_{-h}^{\eta (T)}G(r,z^{\prime },z,T)dz.$$

$$n(r,z^{\prime },T^{\prime }+d{T}^{\prime })-n(r,z^{\prime },T^{\prime })\approx \mathcal{K}[n(\cdot ,z^{\prime },T^{\prime })]d{T}^{\prime }.$$

$${\int }_{-h}^{\eta (T^{\prime })}G(r,z^{\prime },z,T-T^{\prime })\mathcal{K}[n(\cdot ,z^{\prime },T^{\prime })](r)d{T}^{\prime }d{z}^{\prime }.$$

Full text

Oil Spill; Ocean pollution; Oil transport; Oil slick; Ocean oil fate model

I Introduction

A critical component of a mass-conserving oil fate model for oceanic oil spills is a model of the mechanics of oil exchange between surface oil and sub-surface oil. Our primary focus in this work is on capturing the interplay between buoyancy, shearing and advection, as they affect the overall balance of oil in the surface and in the subsurface. In order to do so one has to include many of the details of the microscopic dynamics, for instance the changes in surface tension and in viscosity over time as a result of chemical reactions. In particular, it is very important to track the droplet size distribution, since the size governs how micro-scale processes like merging, diffusion/dispersion, advection, shearing and buoyancy affect a dense population of droplets.

Our focus is on oil transport and spills at environmental scales (e.g., at the small end these would be kilometers and hours, at the larger scales it would be estuarine/basin scales and months/years). Our long-term goal is to develop a model for oil spill dynamics that could be implemented as a simulation platform. This is a daunting challenge for the following reasons:

The environmental scales are several orders of magnitude larger than the microscopic scales characterizing of the droplet dynamics alluded to previously. A model that resolves the dynamics on the microscopic scale will end up with an enormous number of degrees of freedom. 2. 2.

Oil is a composite, consisting of hundreds or thousands of different chemical species. This further compounds the “curse of dimensionality”, if one seeks to accurately capture the multitude of chemical reactions within the oil. An alternative is to lump the individual species into a few chemical complexes. This approach now requires modeling the important chemical processes by empirical rules rather than from first principles. 3. 3.

Finally, there is the matter of coupling the oil transport to models for the ocean, the atmosphere, and land.

In order to reach environmental scales we have embarked on the development of an oil transport model that was first introduced by Restrepo, et al Restrepo et al. (2015) (hereon, RRV). In developing this model, we have placed an emphasis respecting the underlying physics and in reducing its overall dimensionality, possibly at the expense of its quantitative accuracy. The rationale for our methodology is the following –

In order for a model to be useful over long time scales, we need to respect “universal" physical laws. For instance, our model should be matter conserving and should not violate the laws of thermodynamics. Even tiny violations of these laws on short time scales can be amplified over time and lead to large errors or unphysical results over longer time scales. (Mass conservation is a non-trivial and basic requirement in oil reservoir simulations, particularly in chemically-reacting flows Chen (2007)). This is a particular concern for models that are entirely built from empirical rules obtained from measurements on short time scales. 2. 2.

As we argued above, we need to reduce the number of degrees of freedom to develop a model that can form the basis of a numerical simulation platform for oil spills. 3. 3.

While items (1) and (2) from above are certainly desirable, one has to pay a price to achieve them. This price is in terms of a tradeoff between computational feasibility and quantitative accuracy. In particular, empirical models with higher number of degrees of freedom might potentially run quicker and be more accurate for short times, than the mass conserving, reduced dimension model we develop in this paper. However, there are many situations where one needs to study the dynamics of oil spills over long periods of time, and it it then crucial that the model not give unphysical results.

In its current incarnation, the oil transport model is a shallow-water model, i.e., the horizontal spatio-temporal scales of the oil spill are much larger than the vertical scale. This does not imply that there are no significant vertical variations in distribution of the oil within the water column. Indeed, there are significant vertical variations due to the various processes, buoyancy, merging and splitting of oil droplets, turbulent dispersion, that act on the micro-scale. Our key contributions in this paper are that (i) we solve the relevant balance equations for the vertical distribution of oil in the water column; and (ii) we use the solutions to parameterize the vertical distribution in terms of quantities (depth averages) that are slowly varying functions of the horizontal variables. These form the basis of an upscaling method that leads to effectively 2 dimensional modes for the oil transport with significantly fewer degrees of freedom.

In order to couple this oil transport model with models of the ocean, we also demand that the oceanic flow be well captured by 2 dimensional models given by depth-averaging. This is not a significant restriction for practical application of our model. Indeed, familiar depth-averaged circulation models are successfully used to capture oceanic flows in estuaries, the continental shelf and shallow basins. These are regions of the ocean that abut fragile shore environments as well as large cities and population centers. These are also where large scale oil exploration takes place, as does complex shipping traffic.

As we argue above, we view as a fundamental requirement that any practical oil spill model with environmental relevance should conserve mass and describe the geographical distribution of oil over time. In essence should answer the most basic question of where and how much of the oil is to be found at any given time. It is impractical for existing observer network to capture how much oil is found because measurements sample the ocean interior and the surface very sparsely. Consequently, even in a data assimilation framework, it is crucial that the model conserve mass. The crucial ingredient in a mass conserving model is a proper accounting of the vertical exchange of oil between the surface and subsurface domains; this is the subject of this paper.

Although we develop our parameterization for the vertical distribution in the context of our shallow-water oil spill model, the parametrization itself or the ideas we use in the proposal could be ported to other oil models, even Lagrangian models.

This paper is organized as follows: We begin by briefly summarizing the shallow-water oil model and identify where the mass exchange term appears in the model. In RRV we focused on transverse dispersive and advective effects on oil spill dynamics. The mass exchange model, which is mentioned but not described in RRV is proposed next. The mass exchange model is formulated in terms of droplet size distributions in the bulk and it accounts for many distinct processes including buoyancy, vertical mixing due to wave breaking and boundary turbulence and merging/splitting of droplets. The size of a droplet plays a key role in its dynamics and thereby, in the oil model. Our focus is on large scale dynamics, meaning over spatio-temporal expanses that cover hundreds of kilometers and thousands of seconds. At first glance micro-mechanics in the oil could be perhaps ignored. However, surface tension effects, viscous effects, chemistry, and the background hydrodynamics leads to the breakup of oil into droplets. Because these droplets have a dynamic that is size dependent, it is not possible to ignore microscopic physics; for example, droplets that are very small take a very long time to buoyantly reach the sea surface. Whether droplets advect with the flow or are affected by shearing and compression also depend on the size of the droplets is accounted for in our model, as is merging and splitting in the subsurface oil. The model of the process of droplet merging/splitting will describe how size distribution and droplet mechanics can be accounted for at the larger scales of the oil fate model. In the Results section we revisit a nearshore phenomenon that we originally proposed by Restrepo et al. Restrepo et al. (2014) (hereon, RVD) and have since used as a benchmark (see RRV) to test our refinements of the oil transport model. The nearshore sticky-waters phenomenon concerns the dynamics of oil approaching a nearshore environment and the conditions whereupon the oil slows down and even parks itself just beyond the break zone. The original model for nearshore sticky-waters Restrepo et al. (2014) had a crude mass conserving rule for the exchanges of oil mass exchanges between the slick and the sub-surface. The calculation shown here suggests that the original conceptual model gives results that are consistent with the refined mass exchange model that we propose here. In the Summary section we the mass exchange parametrization.

II The Shallow-Water Oil Fate Model

The oil fate model itself is described in RRV. The general transport model is formulated for a surface slick layer $S$, in units of length, and a subsurface component of oil, described by a concentration field $C$ in units of density times volume over volume. While there may be thousands of chemicals in a particular sample of oil and each oil spill may have different ratios of these chemicals, we describe the slick oil in terms of

[TABLE]

where $S({\bf x},T)$, so that the total oil slick mass is given by

[TABLE]

${\bf x}$ is the transverse coordinate, $\Omega_{T}$ is the (possibly time dependent) transverse domain. Time is denoted by $T$. The densities of each chemical complex are denoted by the $\rho_{i}$.

For the $i^{th}$ chemical complex, the transport equation is

[TABLE]

Here, $\nabla_{\perp}$ is the transverse gradient operator, ${\bf V}$ are the ${\bf x}$-components of the transport velocity, $C_{xs}$ is a slip parameter, ${\bm{\tau}}$ is the wind stress, $\mu_{i}$ the viscosity of the $i^{th}$ oil species, $\Psi$ is the dispersion. The last term on the left hand side is the eddy/dispersion term. $R_{i}$ is the rate associated with mass exchanges between the oil slick and the interior of the ocean and is the focus of this study. The reaction term $E^{s}_{i}(b_{i},s_{i},s_{j})$ encompasses chemical reactions among chemical complexes within the slick as well as other processes modeled by rates, such as evaporation. $P^{s}_{i}$ is the rate of biodegradation, emulsification, photodegradation and $G^{s}_{i}$ is the source rate term. Unlike $R_{i}$ the terms with $s$ superscript are particular to the surface oil.

The total mass of the oil in the sub-surface is defined as

[TABLE]

where $c_{i}$ is the concentration of the $i$-th complex in the sub-surface layer, in units of mass per unit volume. Define the depth average concentration as

[TABLE]

$H$ is the total water column depth.

The evolution of the $i^{th}$ chemical complex in the subsurface is given by

[TABLE]

On the right hand side, the first term is the mass exchange term $\rho_{i}R_{i}$, the reaction/chemistry term is $E_{i}^{C}$, $P_{i}^{C}$ represents the sedimentation and biodegradation mechanisms. The subsurface source term is $G_{i}^{C}$.

This paper concerns itself with developing a model for $R_{i}$, which we connote as the mass exchange term and has units of length per unit time. To do so we do not need to carry around reference to chemical complexes, so we will just assume that there is a single chemical complex (we can then drop the subscript $i$ from $s$ and $C$).

III The Mass Exchange Model

The mass exchange $R$ represents oil lost to the slick by turbulent mixing and wave breaking and oil gained from the subsurface by buoyancy. Because buoyancy forces depend on the radius of the droplets, it is necessary to keep track of the droplet distribution in the subsurface, which in turn, requires that droplet merging and splitting be incorporated in the dynamics. The consequence of this is that micromechanics need to be resolved or they must be upscaled if the scales of interest are far larger than the droplet coherence scale. The oil model Restrepo et al. (2015) is being designed to capture the dynamics of oil spills at environmental spatio-temporal scales, which are large. Hence, we will develop the rationale for an upscaled version of the mass exchange process, applicable to environmental scales.

In summary, the mass exchange model takes into account mass exchanges between the slick oil and the subsurface oil as well as the dynamics of droplet merging and splitting in the subsurface. We will be making frequent reference to the study by Tkalich and Chan Tkalich and Chan (2002) and thus we will refer to it as TC, hereon.

III.1 Slick/Subsurface Oil Exchanges

Per unit area, at a fixed horizontal location ${\bf x}$, the total mass of oil is

[TABLE]

where $V(r)=\frac{4}{3}\pi r^{3}$ and $n(r,{\bf x},z,t)\,dr$ is the total number of droplets of radius in $(r,r+dr)$ per unit volume of water at depth $z$. $\bar{n}$ is the total number of droplets of radius in $(r,r+dr)$ per unit area of water surface water, in the subsurface. If we ignore the dispersion, reaction/chemistry, sedimentation, biodegradation and subsurface sources, the total derivative of $M$ must be zero, and therefore

[TABLE]

where

[TABLE]

The first term of $R$ is the rate of oil gained by the slick (in units of length/time), coming from the subsurface by means of buoyancy. The second term is the rate of oil lost (primarily) due to wave breaking (see TC for an alternative derivation). In the first term $w_{b}\geq 0$ is the buoyancy of oil, which is known to depend on the radius of the droplets. The droplet distribution is classified by $V(r)$ the droplet volume, which in turn depends on the radius $r$. The largest resolved radii is $r_{0}$ and the smallest is $r_{m}$. The second term (again, in units of length/time), is associated with the oil lost at the surface to the subsurface, with $S_{b}$ being the wave breaking source term.

III.2 The wave breaking source term $S_{b}$

Wave breaking events occur at the location of interest at a rate $\lambda_{b}$ and have a ‘strength’ that may be characterized, for example, in terms of the slope $hk$ of the breaking wave. For the time scales of interest, a wave breaking event occurring at time $T_{i}$ can be regarded as an instantaneous transport event, in which a fraction $b_{i}$ of the volume $s(T_{i})$ is distributed into the water column. Let $\delta_{T_{i}}(T)$ denote the delta function at $T_{i}$, $\delta_{T_{i}}(T)=\frac{1}{dT}$ for $T\in[T_{i},T_{i}+dT)$, and zero everywhere else. The breaking inputs can then be written in terms of a compound Poisson $N_{b}$ process of rate $\lambda_{b}$,

[TABLE]

We will refer to $b_{i}$ as the strength of the $i$th wave breaking event.

According to Delvigne and Sweeney Delvigne and Sweeny (1988), if a unit of oil per unit area is available to be distributed into droplets, the number of droplets with radius in $(r,r+\,dr)$ will be $N(r)\,dr$ with

[TABLE]

zero otherwise, which is such that $\int_{0}^{\infty}V(r)N(r)\,dr=1$.

$S_{b}(T,r,z)\,dz$ is the rate, per unit time, at which droplets of radius in $(r,r+\,dr)$ become available in the volume $(z,z+\,dz)$ at time $T$. This term is nonzero only at breaking times $\{T_{i}:i=1,2,\dots\}$. At $T_{i}$, the unit of ocean surface area we are considering suffers a breaking event displacing a volume $b_{i}s(T_{i})$ of oil to be distributed along the depth of water. Assume that such distribution is independent of everything else and is given by $g_{b}(z)$ with

[TABLE]

The volume of oil available in the volume $(z,z+\,dz)$ of water is $b_{i}s(T_{i})g_{b}(z)\,dz$, and by (11), the resulting rate of available droplets of radius in $(r,r+\,dr)$ is

[TABLE]

By (12) and (11) we can write the mass exchange term as

[TABLE]

III.3 Subsurface Droplet Distribution Dynamics

We must now specify a dynamical model for the droplets in the subsurface. In accordance with the oil model Restrepo et al. (2015), droplets are assumed to undergo a linear vertical advection-diffusion process. Furthermore, the time scales at which the vertical distribution of oil sorts itself out, given fast disturbances at the sea surface, are assumed to be very short compared to the horizontal advection scales (the shallow water wave approximation).

We denote the turbulent diffusivity by $D(r)\sim 10^{-2}-10^{-3}$\text{,}\mathrm{m}\text{,}{\mathrm{s}}^{-2} which can be parametrized in terms of the turbulent description of the flow, while $w(r)\sim 10^{-6}-10^{-3}$\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1} is the limiting velocity of a oil droplet moving through water as described by TC (Equation 13).

In what follows we focus on vertical balances, an adiabatic approximation, in which for given radius $r$, depth $z$ and time $T$, two flux sources affect $n(r,{\bf x},z,t)$: the droplets that are entrapped from the surface into depth $z$ after a wave breaking event, and changes on the distribution of radii due to merging and splitting of droplets. Putting all together, gives

[TABLE]

Furthermore, $W\ll w_{b}(r)$, hence we will ignore $W$. (We showed in RRV that depth integration effectively dispenses with $W$). The operator $\mathcal{K}$ encodes the effect of merging and splitting of droplets. This operator is discussed in detail in section IV.1, but for now, it suffices to remark that it should be conservative, namely:

[TABLE]

At the surface $z=\eta(T)$ we assume that there is no diffusive flux as the slick is saturated, while all incoming droplets from the upward effect of buoyancy are added to the volume $s$ and removed from the concentration $n$. The boundary condition at the subsurface/slick interface is simply:

[TABLE]

At $z=-h$ we allow for sedimentation, which we can ignore for now.

In what follows we fix ${\bf x}$ (which we omit in the subsequent notation), and solve for the amount of oil in the slick and subsurface droplets occupying the one-dimensional moving water column $(-h,\eta(T))$, ignoring as well horizontal advection/dispersion. Specifically, we solve

[TABLE]

The solution to (18) is:

[TABLE]

where we have made the convenient choice for the initial condition

[TABLE]

and set $T_{0}=0,\,b_{0}=\frac{1}{\rho}\frac{HC(0)}{s(0)}$. The function $G$ in (20) is the Green’s function associated with the vertical terms of the parabolic equation (15) of $n$ which may be approximated with the Greens’ function corresponding to the operator $w_{b}\partial_{z}+D_{z}\partial^{2}_{z}$ on $(-\infty,\eta)$, with boundary condition (17). That is to say,

[TABLE]

In the simplest of possible situations, we illustrate the dynamics of the mass exchange and its impact on the slick and subsurface oil. For simplicity we assume $\mathcal{K}\equiv 0$ and propose a simple backwards difference scheme for he approximate solution of $s$, via (19):

[TABLE]

Here we use the parametrization of $w(r)$ and $D(r)$, suggested by TC, $w(r)=\frac{2g(1-\rho/\rho_{w})}{9\nu_{w}}r^{2}$, $D(r)=D_{0}\left(\frac{D_{0}}{D_{m}}\right)^{-\frac{r}{r_{m}}}$, with $\rho=$800\text{,}\mathrm{kg}\text{,}{\mathrm{m}}^{-3}, $D_{0}=10^{-2}$\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}\text{\,}{\mathrm{}}^{2}, $D_{m}=10^{-3}$\text{,}\mathrm{m}\text{,}{\mathrm{s}}^{-1}\text{,}{\mathrm{}}^{2}, $r_{1}=10^{-6}$\text{\,}\mathrm{m}, $r_{m}=10^{-4}$\text{,}\mathrm{m}. We picked the rate of breaking as $\lambda_{b}=10^{-3}$\text{\,}{\mathrm{s}}^{-1} and the fraction of entrapped oil randomly distributed as $b_{i}\sim\text{Uniform}(0,10^{-1}).$ After a breaking event, the entrapped volume is distributed exponentially over the depth so 95% of it is above $z_{m}=$1\text{,}\mathrm{m}$$, $g_{b}(z)=\kappa e^{-\kappa z}$ with $\kappa=$3\text{,}{\mathrm{m}}^{-1}$.$

For this illustrative example, the evolution of the fraction of slick and subsurface oil, as well as the distribution of the radii in the subsurface is shown in Figure 1. There is an increase in the droplet number for small radii.

IV Depth-Averaged Dynamics

The total amount of submerged particles of radius in $(r,r+dr)$ is

[TABLE]

and is given by

[TABLE]

Here,

[TABLE]

We would like to derive an expression for the dynamics of merging and splitting of droplets that is independent of $z$ and hence arrive to an approximate version of (24) in terms of solely $\bar{n}$. We may assume that the merging and breaking processes occur at time scales much shorter than those of $T$. In particular, during an interval $dT^{\prime}$, the change in the distribution of radii is

[TABLE]

Out of these added droplets, the following amount is then transported to depth $z$

[TABLE]

We approximate the second integral in (24) by conmuting these mechanism, and writing

[TABLE]

which essentially disregards the depth-wise structure of $n$ when it comes to the merging and splitting mechanism. The equation for $\bar{n}$ is now

[TABLE]

When the oil spill is young, or the sea state is rough, droplet distributions change quickly due to merging and splitting. Further, it is well known that droplets in the subsurface smaller than a critical radius, will remain in the subsurface (in TC the critical radius for oil is quoted as $r_{c}=50$ $\mu$m). For computational reasons it is convenient to bin droplet sizes for each chemical complex and formulate the subsurface oil transport in terms of a finite, discrete set of droplet classes. (Tracking the oil slick droplet sizes is unnecessary). At the same time, we will incorporate depth-averaging and pose the dynamics of subsurface droplets in terms of fluxes into and out of the depth-averaged horizontal computational cells due to advection/dispersion, interactions with the slick via buoyancy and wave mixing, and finally via changes to droplet class distributions within each cell.

Let $[r_{j},r_{j+1})$, $j=1,\dots,N_{r}:=m-1$ denote different radius ranges to be considered and denote:

[TABLE]

Ignoring advection and dispersion, an equation for the evolution of $\bar{n}_{j}$ can be obtained via equation (15):

[TABLE]

We now make the following approximations: we use $w_{b}(r)\approx w_{b}(r_{j})$ for all $r\in[r_{j},r_{j+1})$. Lastly, we assume that at the time scale of interest, we can write

[TABLE]

The operator $\mathcal{K}$ can be simplified, by ignoring depth dependence and we can thus propose a dynamic between size classes that takes into account droplet merging/splitting and interactions with the eddies in the flow, within the subsurface compartment. We then arrive to

[TABLE]

where $N_{j}:=\int_{r_{j}}^{r_{j+1}}N(r)\,dr$, and $j=1,2..,m$. And similarly for the slick:

[TABLE]

The terms $B_{j}$ and $D_{j}$ are birth and death rates for droplets in class size $j$. These two terms will be developed fully in Section IV.1.

IV.1 Subsurface Oil Droplet Merging and Coalescence

We will focus now on the merging/splitting mechanism and ignore all other effects. Let us consider the situation for a subsurface cell at position ${\bf x}$. The changes to the total number of droplets of size class in $(V(r),V(r+dr))$ occupying a fixed volume of water, due to breaking and merging obey the balance equation

[TABLE]

where $D_{V},B_{V}$ are death and birth rates, collectively defined as $\mathcal{K}[\bar{n}]$, that approximate the Smoluchowsky operator $\mathcal{K}[n]$. The death $D_{V}$ and birth $B_{V}$ terms, in turn, are defined as

[TABLE]

where $V:=V(r)$, $V^{\prime}:=V(r^{\prime})$. In the above expressions $\lambda(V,V^{\prime})$ is the coalescence efficiency for collisions between droplets of classes $V$ and $V^{\prime}$, $h(V,V^{\prime})$ is the frequency of droplets of class $V$ colliding with droplets of class $V^{\prime}$, $g(V^{\prime})$ is the breakage frequency of droplets of class $V^{\prime}$, $\nu(V^{\prime})$ represents the number of droplets formed by breakage of droplets of class $V^{\prime}$ and $\beta(V,V^{\prime})$ is the probability of forming droplets of class $V$ from breakage of droplets of class $V^{\prime}$. To first order, it can be assumed that only binary collisions occur, and thus $\nu(V^{\prime})=2$.

Once the droplet distribution is binned (assuming the total number of bins $m$ is even), and the distribution of droplet sizes $r_{1}$ and $r_{m}$, is labelled by $j=1,..,m$. The $HD_{j}(\cdot,T)\approx D_{V(r_{j})}$ and $HB_{j}(\cdot,T)\approx B_{V(r_{j})}$ are:

[TABLE]

IV.1.1 Parametrizing Breakup Mechanics

Droplet breakup is caused when an external force in the surrounding fluid exceeds the surface and internal forces in the droplet (Liao and Lucas, 2009). Specific breakup mechanisms are droplet elongation in sheared flow, turbulent pressure fluctuations, relative velocity fluctuations, and droplet-eddy collisions. The breakup frequency, $f(d)$, assuming droplet-eddy collisions in isotropic turbulence conditions, droplet sizes in the inertial sub-range, and droplet breakup due solely by eddy-droplet collisions with smaller or the comparably-sized eddies is defined as:

[TABLE]

$K_{b}$ is a system dependent parameter determined as calibration factor (Zhao et al., 2014).

The eddy-droplet collision frequency, $h$, is defined as:

[TABLE]

where $S_{ed}=\pi(r_{e}+r)^{2}$ is the collision cross-section area, $r_{e}$ is the eddy size which collide with a droplet of radius $r$, $dn_{e}$ is number of eddies of size between $r_{e}$ and $r_{e}+dr_{e}$. Here, $u_{e}^{2}=3.8(\epsilon/k)^{\frac{2}{3}}$ is an estimate of the velocity-squared of the eddy, $k=\pi/r_{e}$ is wavenumber, and $u_{r}^{2}=1.7(\epsilon r)^{\frac{2}{3}}$ is the droplet velocity. The estimated number of eddies per unit volume of the fluid is found from $dn_{e,m}(k)=3\times 10^{-3}k^{2}dk$ (see (Azbel, 1981)).

The breakup efficiency $B_{e}$ is the probability of breakup if a drop-eddy collision takes place:

[TABLE]

The minimum energy for droplet breakup is related to the amount of required energy to keep the surface of a droplet intact. For binary breakage, this energy is assumed equal to the mean energy required to break the drop into two equally-sized droplets (maximum energy) and into a one small and one big droplets (minimum energy) (Tsouris and Tavlarides, 1994). Hence,

[TABLE]

where $e=\frac{1}{2}mu_{e}^{2}=\frac{2}{3}\pi\rho_{w}r_{e}^{3}u_{e}^{2}$, $r_{\min}=r_{1}$, $r_{\max}=r_{m}$, and $\sigma_{d}$ is the droplet surface tension (Zhao et al., 2014).

Oil droplets are subjected to surface tension effects as well as internal viscous forces. Hence we need to take into account the balance of these (see (Calabrese et al., 1986)). The energy associated with this internal force is

[TABLE]

where $\rho_{w}$ is the fluid density, $\rho$ is the density of the particular chemical complex oil density, $T_{tot}$ is interaction simulation time (one time step of fluid-oil interaction), and $\mu_{d}$ is the oil species viscosity. The constant $a=1.1$ is a calibration constant (see Zhao et al. Zhao et al. (2014), Equation 15).

In Tsouris and Tavlarides (Tsouris and Tavlarides, 1994) it is argued that, assuming binary breakage (recall it was assumed that $\nu(V^{\prime})$ is 2), the amount of energy required for the generation of two like-sized droplets is larger than what is required to produce one small and one big droplets. They suggested the breakup PDF of a drop of class size $V$, from a parent drop of class size $V^{\prime}$ is given by:

[TABLE]

where:

[TABLE]

IV.1.2 Parametrizing Coalescence Mechanics

The merging process in which two (or more) colliding droplets form a bigger droplet is called coalescence. The coalescence of droplets normally evolves in three stages: (1) collision; (2) interfacial drainage of liquid at touching interfaces; and (3) the rupture of the droplet surface. In this study we assume that background turbulence is the main reason for the collisions. The volume absorbed by a moving droplet when colliding with other droplets is the basis for calculating the collision frequency. The swept volume rate $h(V,V^{\prime})$ is modeled following the kinetic gas theory. Assuming both droplets moving with turbulent eddies with relative velocity of $c$, a collision tube with cross section of $\pi(r+r^{\prime})^{2}$ is swept. Therefore, the $h(V,V^{\prime})$ could be defined as:

[TABLE]

Assuming isotropic turbulence Kolmogorov (1941), in Kuboi et al. Kuboi et al. (1972) the relative colliding velocity in stirred and pipe flows was estimated as:

[TABLE]

Hence the swept volume rate for droplets $V$ and $V^{\prime}$ reads:

[TABLE]

The droplets need to stay in contact for a certain time to complete the coalescence process. Turbulence promotes collisions, however, it can also move the droplets away from each other before the coalescence process takes place. Therefore the collision efficiency is defined as the ratio between the average drainage time $\bar{t}$ and average contact time $\bar{\tau}$:

[TABLE]

(See Prince and Blanch (Prince and Blanch, 1990)). The drainage time of the thin interface between two droplets depends on the viscosities of each droplet and the main phase (could be different for air bubbles in comparison with oil droplets). As a result we need to be more careful about implementation of this parameter. The average drainage time for fully mobile interfaces is given by Prince and Blanch (1990), Equation 18 :

[TABLE]

where $\rho_{l}$ is the density of the continuous phase, $\sigma_{d}$ is the surface tension and $h_{0}$ is the initial film thickness. The critical film thickness $h_{c}$, as estimated in (Chesters, 1991), is

[TABLE]

where $A$ is the Hamaker constant and $R_{eq}$ is the equivalent droplet radius, given by:

[TABLE]

where $r$ and $r^{\prime}$ are the radii of the colliding droplets. The average contact time $\bar{\tau}$ is given by

[TABLE]

(see Levich Levich (1962)).

IV.2 The Complete Depth-Averaged Model

We now summarize the physical/computational shallow-water oil fate model. We use a split-step strategy, with time step $\Delta T$.

• •

For each chemical complex, at some time $T+\Delta T/2$, the slick is updated and denoted by $s^{*}$. This is done using $s(\cdot,T)$ and (33) (the sum runs for $j=1,\ldots,m$). For the subsurface, we produce the update $\bar{n}^{*}_{j}$ via (32) knowing $\bar{n}_{j}(\cdot,T)$.

• •

To step to $T+\Delta T$, we proceed as follows: Using $s^{*}$ we update $s(\cdot,T+\Delta T)$ using (1) with no $R$ (it has already been taken into account in the half step). We also update the subsurface concentration: we replace $HC(\cdot,T+\Delta T)=\sum_{j=1}^{m}\rho V_{j}\bar{n}_{j}(\cdot,T+\Delta T)\Delta r$ in (3), with $R=0$, and use $C^{*}:=\rho\sum_{j=1}^{m}V_{j}\bar{n}^{*}_{j}\Delta r$ with $\Delta r=R_{r}/m$, where $R_{r}$ is the range of radii of the oil droplets, and $m$ is the number of bins used in partitioning the distribution of oil droplets, classified by their radii.

In total, for each chemical complex there will be 1 evolution equation for the slick, $m$ for the subsurface. All told $N\times(m+1)$ degrees of freedom, where $N$ is the number of chemical complexes.

V Results

We first benchmark the mechanics of merging and coalescence in our model by comparing its predictions to the results shown in Tsouris and Tavlarides (Tsouris and Tavlarides, 1994), since this work largely inspired the droplet merging and breaking dynamics of our oil model. We then explore the dynamics within a single stack. A stack is composed of a surface compartment, and a subsurface compartment, directly below. The surface compartment corresponds to the slick. The subsurface compartment has oil and water. In the stack computations that follow, we will feature the birth/death dynamics as well as the mechanics between the slick and the subsurface compartments. Since it is a single stack, no advection and dispersion takes place.

Finally, we explore a simple nearshore ocean/oil dynamics problem related to nearshore sticky waters. In RVD we proposed a conceptual model for the slowing down and possible parking of buoyant material, such as oil, as it approaches the nearshore. The model suggests that oil being advected by ocean currents will be largely subjected to dispersive effects that gain prominence near the shore. We specifically claim that the dynamics of nearshore oil cannot be properly captured if dispersion is ignored. In fact, it is essential in explaining how oil reaches the shore even when the depth-averaged velocity normal to the shore goes to zero, if the oil cannot establish a conveyor dynamic near the shore. We have been replacing simple conceptual components of the original model proposed in RVD to show that, (1) as the mechanics of the model are refined, the proposed dynamics for sticky waters still holds; (2), how the oil model proposed in RRV captures the phenomenon of nearshore sticky waters. In the calculations that will follow the stack benchmark we will highlight how dispersion, advection, and now droplet dynamics play out in the nearshore context; (3) we also want to suggest that the nearshore sticky water problem, which is in itself an important dynamics problem for pollution transport close to the shore, could be used as a field problem to test many aspects of oil transport dynamics and models.

V.1 Benchmarking the Merging and Coalescence Mechanics in our Model

We benchmarked the birth/death mechanism in our model by comparing our results to those of in Tsouris and Tavlarides (Tsouris and Tavlarides, 1994). Rather than replicating here what Tsouris and Tavlarides did in the comparison of their model to data, we suggest the reader refer to Tsouris and Tavlarides (Tsouris and Tavlarides, 1994) for the details. In order to accomplish a similar comparison between our model and data we omit the slick dynamics, and suppressed all of the mechanics in the subsurface dynamics, except for the terms associated with $B$ and $D$. Hence, we are using only (34) in this benchmark comparison. In the calculations we adopted the same parameter values as those used to generate Figures 3 and 4 in Tsouris and Tavlarides (Tsouris and Tavlarides, 1994): For continuous phase (water): the viscosity is $\mu_{c}=1.787\times 10^{-3}$ Ps, water density $\rho_{w}$ = 1025.0 kg/m3, and the kinematic viscosity $\nu_{c}$ = $\mu_{c}/\rho_{w}$ m2/s. For the dispersed phase (Toluene), we choose a surface tension $\sigma_{d}$ = 28.4 $\times 10^{-3}$ N/m, a dispersed phase molecular viscosity $\mu_{d}$ = 0.59 $\times 10^{-3}$ Ps , and a density $\rho$ = 865 kg/m3. The droplets ranged in size from $r_{1}$ = 0.0001 m to $r_{m}$ = 0.00038, and we used $m-1=25$ size classes. In Figure 2 we superimpose the output of our model (lines) with digitized disconnected dots representing the experimental data used in Tsouris and Tavlarides (Tsouris and Tavlarides, 1994) and originally obtained by Bae and Tavlarides Bae and Tavlarides (1989).

Our results show similar consistency between model and data as was reported by Tsouris and Tavlarides, when they compared their model to the same data, especially when we compare the evolution of the smaller size classes. We note that our model is considerably simpler than theirs, and apparently the simplification degrades the performance of our model for the larger size classes. This degradation of fidelity is acceptable and deliberate: for the large spatio-temporal scales demanded by the applications of our oil fate model, gains in computational efficiency at the expense of small degradations in fidelity is an acceptable tradeoff.

V.2 Mass Exchanges in a Single Stack

We consider next the dynamics of a stack, under equations (32) and (33). For this case the size classes were chosen using $r_{1}$ = 0.0001 m and $r_{m}$ = 0.00038 m and a volume ratio of $V_{i}/V_{i-1}$ = 1.5 resulting in 9 different size classes. We simulated one hour of random wave breaking given by (10) with $\lambda_{b}=\frac{1}{6T_{w}}$ ($T_{w}$ being the wave period, see Restrepo et al. Restrepo et al. (2011)). The material properties of the oil and the water are the same as was used in the preceding benchmark calculation. We also initialized oil at the surface and inside the mixing layer with equal volume to follow the dynamical interaction between surface and sub-surface oil wave breaking strength $b_{i}$ of varying ranges.

Figure 3 shows the total oil in the slick and the subsurface as a function of time for wave breaking strength $b_{i}\sim\text{Unif}(0,10^{-2})$. The wave breaking strength takes high enough values to force a substantial portion of the oil in the slick to go to the subsurface. The fluxes associated with wave breaking nearly balance, in this case, with the fluxes associated with buoyancy. Larger $b_{i}$ lead to wave breaking fluxes dominating.

The situation changes for a sufficiently small wave breaking strength. See Figure 4. In this figure we plot the mass balance as a function of time for different wave breaking strengths $b_{i}$ of maximum possible value 0.001, 0.01, and 0.1. For small $b_{i}$ the buoyancy fluxes dominate and the mass changes smoothly after a very short transient. As the breaking strength increases, the breaking fluxes dominate, and for very large breaking strengths the mass reaches a quasi-equilibrium configuration.

In Figure 5 we plot the equivalent thickness of each size class as it evolves driven by the wave breaking regimes of Figure 4. The main conclusion drawn from this calculation is that the oil in the large size classes changes the overall balance of oil in the slick and the subsurface, while the small size classes vary in time more dramatically for any breaking strength.

When transversely coupled stacks are considered, we can introduce advection and dispersion mechanics, which along with transversal heterogeneity of the wave forcing, would alter considerably the dynamics of coupled stacks. We will highlight next a special environmentally-inspired case in which we can explore the stack dynamics coupled with advection and dispersion.

V.3 Mass Exchanges In The Nearshore Sticky Waters Problem

In RVD we proposed a conceptual model for the slowing down and possible parking of buoyant pollutants, such as oil, approaching a shore due to ocean circulation. When oil is brought toward the shore by ocean currents, it was found that the advective effects and the dispersive effects can balance each other in this nearshore so as to slow down or park a buoyant tracer. We coined the phrase nearshore sticky waters to refer to the phenomenon of slowing or parking of pollutants. The basic mechanics are proposed conceptually in RVD. In RRV the nearshore sticky waters problem was revisited in order to examine how the shallow-water oil transport model that is being developed here captures the phenomena.

Here we incorporate the mass exchange mechanism proposed in earlier sections of this paper into the shallow water oil fate model. The scenario is presented in Section 6.1 of RRV, and we are only going to point out in what ways the present calculation differs from the ones in RRV.

The changes to the original model, given by Equations (39), (40), (43) in RRV, and the one captured by (1), (3) and (14), in the algorithmic process described in Section IV.2, are as follows:

• •

For simplicity, we will use a constant mixed layer depth $P$ in the calculations;

• •

the subsurface oil will be captured by 4 size classes, chosen in the range $r_{1}=0.0001$m and $r_{m}=0.00038$m. The only size class that can buoyantly reach the slick is the oil contained in the largest subsurface size class. In this regard we are following TC, who report that only oil droplets larger than about $0.0005$m are affected significantly by buoyancy forces.

• •

The mass exchange term in the dynamics of the slick and the subsurface oil, used in the older paper (which reads $\frac{(1-\gamma)s-\gamma PS}{\varsigma(x)}$), are replaced by the mass exchange terms derived in this paper; further,

• •

to the mass exchange component of the subsurface (32) associated with oil entrainment due to wave breaking $\sum_{i=1}^{N_{b}(T)}N_{j}s(T_{i})b_{i}\delta_{T_{i}}(T)$, we add $s(T)\frac{2\epsilon}{\rho_{w}gA^{2}}$. Here, $\epsilon$ is the wave action damping (see Equation (33) of RRV, $\rho_{w}$ is the ocean density, $g$ is gravity, $A$ is the wave amplitude. In the calculations that follow the wave amplitude in the deeper end of the domain is 1 m, and the wave period is 9.1 s. The new term accounts for the intense turbulence associated with the shoaling process of the waves. In general it tends to dominate the mass exchanges due to wave breaking in the shoaling region;

• •

similarly, to the mass exchange component of the slick (33) associated with oil fluxes to the subsurface, $\sum_{i=1}^{N_{b}(T)}s(T_{i})b_{i}\delta_{T_{i}}(T)$, we add $s(T)\frac{2\epsilon}{\rho_{w}gA^{2}}$ to account for vertical oil fluxes due to wave shoaling;

Each computational grid location in the transverse coordinate $x$ has a stack composed of a multi-size subsurface compartment, and an upper slick compartment. Variations along the direction $y$ of the shore, are ignored. The calculation is performed as follows: the model performs mass exchange calculations using an adiabatic approximation in the transverse distribution of oil among the many coupled stacks. The mass exchanges within each stack is run for about 1.5 hrs. Then the advection and dispersion across the coupled stacks is performed. This process is repeated. The total time lapsed is in the order of 1000 hrs. The parameters used in the calculation are the same as those that were used to generate Figure 7 of RRV. Figure 6a shows the initial oil distribution. The total amount of oil is $v_{tot}$= 2.334 $\times 10^{-3}$ m3/m/m, and all of the oil is initially the slick. Figures 6b-c show the final oil distribution, roughly after 1000 hrs. The oil is reported on a logarithmic scale as the volume ratio with respect to $v_{tot}$ in a size class or in the slick.

The cases shown in Figure 6b-c differ only in the thickness of the mixed layer, with (b) corresponding to $P=1$m, and (c) to $P=6$m.

As was shown in RVD and RRV the $P=1$m case shows little evidence of slowing or parking. In the earlier papers we reported that the oil that initially started in the slick stayed in the slick even very close to the shore. This is no longer the case with the mass exchange term proposed in this paper. The slick still contains most of the oil away from the shore, but there is a redistribution of the oil between the slick and the subsurface due to shoaling. The $P=6$m shows that there is a parking of the oil away from the shore: the center of mass is roughly 200m away from the shore. The oil has also redistributed between the slick and the subsurface both close as well as far away from the coast. In both of these cases the bulk of the subsurface oil is found in the smaller droplet size class.

VI Summary

In this paper and in RRV we propose an Eulerian, mass-conserving, semi-empirical, oil fate model that is specialized to shallow water environments. The model, whose development is far from complete at this stage, is being designed to capture oil spills at very large spatio-temporal scales required for environmental applications, in the shallow reaches of the ocean. Achieving the large scales is largely due to our ability to parametrize the vertical distribution in terms of depth-averaged quantitites. In RRV we focused on a proposal for the horizontal transport through advection and the dispersion mechanics. In this paper we focused on the mass exchanges between the slick, which is formally the oil present on the surface, and the subsurface oil. The key to modeling vertical exchanges between surface and subsurface oil is that microscale mechanics are not ignored. However, in order to handle large scale dynamics we cannot resolve micromechanics. In the vertical direction the waves and the turbulence in the oceanic water column, particularly close to the surface, can send oil to the deeper reaches of the water column. At the same time, most oil is buoyant and this buoyancy counteracts the downward flux due to waves and turbulence. Surface tension, viscous effects, and shearing easily breaks up oil into tiny droplets. Stokes balances on droplets tend to keep smaller droplets from resurfacing. Furthermore, oil droplets in high concentrations will merge/coalesce and breakup due to fluid stretching and shearing, affecting the overall distribution of oil droplets. This means that advection, dispersion, droplet dynamics, buoyancy and near surface eddy and shearing forces are critical components of an oceanic oil fate model. In this and in RRV we have proposed models for these phenomena. Other effects, such as reaction, sedimentation, and biodegradation will be added to this model in the future. These bring into the model the inherent challenges associated with handling the dynamics of a buoyant material made up of thousands of simpler chemical components. Within the model development program proposed in RRV we will next proceed to incorporate weathering or aging. Oil is a chemically-reacting substance and it is composed of a large number of chemicals. Clearly, tracking many chemicals increases substantially the computational expense of a general oil transport model. The main thrust of tackling oil aging is thus to derive equations for a few lumped chemical complexes that mimic the behavior of the thousands of chemicals present in real oil.

Neither the mass exchange model we propose in this paper and the subsequent weathering model will be specific to our transport model. Hence, in principle, the outcome of these modeling efforts will have a bearing and perhaps find application in other oil-fate models.

The mass exchange model proposed in this study upscales the microscale mechanics and captures merging, splitting, vertical transport due to buoyancy and turbulence. The exchanges preserve mass conservation, a feature we consider essential to an oil-fate model. The upscaling is achieved by exploiting the disparate time scales of (slower) advection/dispersion and (faster) microscale mechanics. This assumption has obvious computational implications: the faster microscale mechanics can be upscaled to the slower scales of transverse oil motion. Properly exploited, this leads to a model that can handle large transverse spatio-temporal scales without resorting to small scale stationarity that could allow for the use of homogeneization of the small scale dynamics.

The merging and coalescence mechanics in the model were benchmarked against data and shown to have reasonable consistency. When compared with our intuition, the vertical fluxes due to wave breaking and buoyancy were shown to be reasonably well parametrized with minimal free parameters and with enough phenomenology to capture general oceanic situations. Specialized experiments and field campaigns would be necessary to test and refine the model further.

In this study we also applied the mass exchange model to a conceptual model that explains how oil and other buoyant pollutants behave dynamically in the near shore zone. Nearshore sticky waters still persist, as proposed in RVD: the added micromechanics and mass exchanges do not change the basic premise behind the suggested explanation for the slowing down or parking of oil and other buoyant pollutants that can reach the shores via oceanic circulation. The mass exchange model suggests that the mixed layer depth has a significant effect on the overall distribution of oil: when the mixed layer is thin the distribution remains fairly constant in time, far from the shore. However, this is not the case when the mixed layer is large. Our model suggests that oil droplet distributions play a critical role in the dynamics of nearshore oil dynamics, and a more careful and extensive analysis of the model and of field data merits attention.

Acknowledgements

This work was supported by GoMRI and by NSF DMS grants 1524241, 1434198 and 1109856. JMR also acknowledges the hospitality of Stockholm University, where some of this work was done, and the support provided by a  Stockholm University Rossby Fellowship.