Vertical mixing of oil droplets by breaking waves

Vertical mixing of oil droplets by breaking waves

Marine Pollution Bulletin 44 (2002) 1219–1229 www.elsevier.com/locate/marpolbul Vertical mixing of oil droplets by breaking waves Pavlo Tkalich *, En...

242KB Sizes 0 Downloads 40 Views

Marine Pollution Bulletin 44 (2002) 1219–1229 www.elsevier.com/locate/marpolbul

Vertical mixing of oil droplets by breaking waves Pavlo Tkalich *, Eng Soon Chan Tropical Marine Science Institute, The National University of Singapore, 14 Kent Ridge Road, Singapore 119223

Abstract Oil spilled on a sea surface can be dispersed by a variety of natural processes, of which the influence of breaking waves is dominant. Breaking waves are able to split the slick into small droplets, facilitating oil mixing in the water column. Vertical dynamics of the droplets plays a major role in the oil mass exchange between the slick and the water column. In this paper a mathematical model of oil droplet mixing by breaking waves is developed. The model uses a kinetic approach to describe the vertical exchange of the droplets at the interface between the slick and the water column. The majority of the coefficients and parameters are conveniently combined into a single ‘‘mixing factor’’. The model is verified using sensitivity analysis and empirical formulae of other authors. The model permits a rapid estimation of the amount of dispersed oil under the breaking waves. The ultimate goal of the research is to parameterise influence of breaking waves on vertical mixing of oil droplets to be used in a general 3-D oil spill model. Ó 2002 Elsevier Science Ltd. All rights reserved. Keywords: Oil spills; Breaking waves; Droplet kinetics; Entrainment; Mixing processes; Dispersion

1. Introduction Marine oil pollution has dramatic consequences for the aquatic and onshore environment. In spite of global efforts, current methods of oil storage, transportation and treatment do not allow the complete elimination of technogenic spills. Once spilled on the sea surface, oil forms a thin layer, known as the slick. Composition of the oil changes from the initial time of the spill. Light (low-molecular weight) fractions evaporate (20–40%), water-soluble components dissolve in the water column (1–5%), and immiscible components become emulsified and dispersed in the water column as small droplets. Wind, waves and currents speed up the droplet generation, the smallest of which may propel deep into the water column to be dispersed eventually by the currents. Certain chemicals, known as dispersants, can reduce the oil–water interfacial tension significantly, thus facilitating much smaller droplets to shear from the slick. A larger number of droplets become small enough to experience the excessive mixing force of turbulence over buoyancy, increasing the dispersed oil mass. Although the above description is commonly accepted, majority of

*

Corresponding author. E-mail address: [email protected] (P. Tkalich).

models do not account accurately for the droplet kinetics at the interface between the slick and the water column. Reviews (Huang, 1983; Spaulding, 1988; ASCE, 1996; Reed et al., 1999) have provided a summary of data and models, available for the evaluation of oil vertical mixing by breaking wind-waves (the phenomenon is often referred to as ‘‘dispersion’’ in the oil spill literature). The earliest models, such as Blaikley et al. (1977), employed the tabulated dispersion rate depending on oil type, sea state, and time after the spill. Audunson (1979), using the square of the wind speed as a scale for the breaking wave energy, suggested a first order kinetics for the oil dispersion. Even though the approach looked promising, it does not differentiate types of oil, making it difficult to apply to real spills. In a later approach of Spaulding et al. (1982), a decay function was included in the model to account for oil weathering and ‘‘mousse’’ formation. However, according to Spaulding (1988), the ‘‘technique is not sufficiently developed for use in practical spill models’’. Johansen (1982), using the same first order kinetics, scaled the oil entrainment rate with the ‘‘whitecap coverage’’ parameter, W. The parameter was estimated as being a linear function of the wind speed, U, i.e. W  U / , where / ¼ 1. If one takes into account that the exponent / may vary within the interval 0 to 3.3,

0025-326X/02/$ - see front matter Ó 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 2 5 - 3 2 6 X ( 0 2 ) 0 0 1 7 8 - 9

1220

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

depending on the wind speed range (Monahan, 1971), then Johansen’s approach looks similar to the Audunson’s model in terms of wind speed dependence, and exclusion of major properties of the spilled oil. Much more sophisticated is the semi-empirical model of Mackay et al. (1980) that considers a slick consisting of two parts: thick and thin slicks, both having dispersion rates proportional to U 2 . In addition, the dispersion rate is inversely proportional to the oil–water interfacial surface tension coefficient, r, and oil viscosity. For a decade after the publication, Mackay’s model was often regarded as the best available model. One of the first successful attempts to look at the oil spill as the three-dimensional dynamic phenomenon was done by Elliott (1986). The random walk technique was used to follow the motion of individual oil droplets in the water column. The action of a vertical shear on the droplets was accounted utilizing vertical velocity profiles due to tide and wind forcing, and the Stokes drift for linear waves. The latter phenomenon is relatively important for non-breaking waves. With the formation of breaking waves, upper layer of the water column (together with oil contents) becomes well mixed to the depth of approximately significant wave height. Mathematical description of this phenomenon along with the empirical formula for oil droplet entrainment (Delvigne and Sweeney, 1988), being combined with the Elliott’s (1986) particle approach, laid basics for the Reed et al.’s (1994) model. One of the key issues in oil vertical mixing is a simple and yet correct quantitative description of the phenomena of droplet formation, dynamics and size distribution. The exact mechanism of droplet formation is not fully understood yet. Theoretical approaches to the problem have been proposed by Hinze (1955), Raj (1977), and Li and Garrett (1998), using different assumptions on the dominant forces during the droplet break-up. The authors presented relationships for a maximal possible radius of the droplets as a function of various oil parameters and available mixing energy. The droplet size spectra and the droplet vertical distribution were studied by empirical (Forrester, 1971; Delvigne and Sweeney, 1988; Delvigne and Hulsen, 1994) and theoretical (Aravamudan et al., 1980; Bouwmeester and Wallace, 1986) means. Based on extensive laboratory experiments, Delvigne and Sweeney (1988), and Delvigne and Hulsen (1994) derived a relatively simple oil droplet entrainment model, elements of which are now widely used in commercial and research models, such as ADIOS (ADIOS, 1994), OILMAP (Spaulding et al., 1992), OSCAR (Reed et al., 1995), and SINTEF (Daling et al., 1997). The success of the model may be attributed to the careful account of major physical parameters during the oil entrainment. Until now the experiments of Delvigne and co-workers are regarded as the most extensive, and the

obtained data and the empirical model are the most complete. The entrainment rate initially was found to be inversely proportional to the oil kinematic viscosity, mo , as 1/mo . More recent data of Delvigne and Hulsen (1994) did not confirm the relationship for low-viscous oil. For high viscosities, the dispersion decreased at a much higher rate. Since the flume experiments were performed for a limited temperature range and for a few types of oil, the derived relationships may not be accurate beyond the experimental conditions. There is no differentiation of oil entrainment on the interfacial surface tension or oil density, which make difficult the model application for the cases of dispersant spraying. To use the model for a greater variety of oil parameters and environmental conditions, it has to be generalized utilising some additional theoretical considerations, such as energy/mass conservation, scale analysis, etc. In this paper a mathematical model of oil droplets vertical kinetics due to the breaking waves is developed. To satisfy basic physical laws, and to ensure accuracy of the derived relationships between the sparse data, general theoretical considerations are used for the model development, and then the equations are verified using the experimental formulae. Since the research was motivated by the need for a more accurate and physically relevant parameterisation of wave breaking phenomena in general oil spill models, such as ADIOS (ADIOS, 1994), OILMAP (Spaulding et al., 1992), particle model of Reed et al.’s (1994), and MOSM (Tkalich et al., 1999), the kinetic model takes care only of vertical exchange of the droplets between the slick and the mixing layer. Rest of numerous phenomena, taking place during or after the spill, are out of concern of the study, even though some of them are used in the paper for illustrative purposes. The structure of the paper is as follows: in the first two sections, energy of waves is budgeted for the most typical breaking scenarios, and the dimensionless droplet size spectrum is introduced. Then mechanisms of the droplet entrainment and resurfacing are considered, leading to the kinetic model of oil droplet vertical mixing. Last two sections of the paper concern with the model verification using other empirical formulae and sensitivity analysis.

2. Energy of breaking waves Breaking waves can introduce significant energy into the upper ocean layer. The rate of the wave energy dissipation, dE=dt, is different for different wave and environmental conditions. However, most cases are expected to follow the decay law (Hasselmann, 1974) as dE ¼ cxE: dt

ð1Þ

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

Here, E is the wave energy, x is the wave frequency, c is the dimensionless damping coefficient, t is the time. The right-hand-side term of the Eq. (1) is linear for the pressure-induced decay, and quasi-linear for wave–wave and wave–whitecap interaction (Hasselmann, 1974). Further details on different forms of the damping coefficient can be found in Komen et al. (1984, 1994), SWAMP (1985), and Young (1999). In the case of the equilibrium energy transfer from the wind to waves and from the waves to whitecaps, the wave energy E will not deviate significantly from the mean wave energy, Ew , as given by (Komen et al., 1994) Ew ¼ gqH 2 =16;

ð2Þ

where H is the significant wave height, q is the water density, g is the gravity acceleration. The wave height depends on the wind speed, U, and local hydrological and meteorological conditions. The general relationship H  U 12 is expected at the different growing stages of the wave spectra (Shore Protection Manual, 1984; Komen et al., 1994), namely,  U for limited fetch or duration; H ð3Þ U 2 for \fully developed sea": In the case of the ‘‘fully developed sea’’ approach for the Pierson–Moskowitz wave spectrum, the wave height is H  0:2U 2 =g (Neumann and Pierson, 1966). Taking Eq. (2) into account, the relationship (1) can be rewritten as dE  cgqxH 2 =16; ð4Þ dt where one can chose the damping coefficient c according to the white-capping (Hasselmann, 1974), swell decay (SWAMP, 1985) or any other suitable model, i.e.  5 0:25 10 xEw ; for \white-capping"; c¼ ð5Þ 1:8  107 x3 ; for swell decay: Influence of oil slicks on wave spectra can be neglected, provided the waves are formed outside the affected field, and dimensions of the slick are significantly smaller than the fetch length. There are neither reliable data nor a model on how an oil film affects breaking waves. The slick impact is presumably maximal for the small-scale capillary waves, and is small for the larger waves. 3. Oil droplet population under breaking waves It was noted from observations that breaking waves develop a well-defined mixing layer in the upper part of the water column (Fig. 1). It has been reported (Donelan, 1978) that the immediate mixing of near-surface waters due to turbulence, produced by whitecaps, penetrates to a depth of the order of the wave height. According to Li and Garrett (1998), the region of strong dissipation extends from a depth aH to a few meters

1221

Fig. 1. Sketch of the vertical structure of the mixing layer.

(where a ¼ 1:2–1:6 is a dimensionless scaling factor depending on the sea state). Laboratory measurements of oil entrainment from the slick to the water column (Delvigne and Sweeney, 1988) led to a similar conclusion, i.e. a ¼ 1:50  0:35. Generalizing the empirical data, one can assume uniform mixing of the droplets within the mixing layer zm ¼ aH ;

ð6Þ

where zm is the layer thickness (m). Below the mixing layer, the vertical distribution of the droplets is governed by the advection and turbulent diffusion phenomena. According to Delvigne and Sweeney (1988), the shape of the steady-state oil droplet size spectra can be expressed as the power relationship N^ ðrÞ  rs , where N^ ðrÞ is the number of the droplets (per unit water volume), r is the droplet radius and s is an empirical parameter. The exponent is found to be s ¼ 2:3  0:06 in the experiments, though it may vary with the oil density and viscosity, weathering state, temperature, etc. The droplet size distribution was observed within the range ðrmin ; rmax Þ  ð1; 1000Þ lm, where rmin and rmax are the minimal and maximal radius, respectively. The droplet spectra can be normalized according to N^ ðrÞ ¼ N 0 =rs ;

ð7Þ 0

where the normalization function, N , can be found from the expression Z rmax 3s 3s N^ ðrÞ dv ¼ 4pN 0 ðrmax  rmin Þ=ð3  sÞ ðs 6¼ 3Þ: V ¼ rmin

ð8Þ Here, V is the volume of oil in all droplets within unit volume of the mixing layer, and v ¼ 4pr3 =3 is the volume of a single droplet. Combination of Eqs. (7) and (8) leads to the dimensionless droplet size spectra N^ ðrÞ ¼

ð3  sÞV 1 : 3s 3s 4pðrmax  rmin Þ rs

ð9Þ

The formulation (9) allows the spectra parameters s, rmin , and rmax to be different for each particular case. One has to note that Eq. (9) is valid only for oil-in-water emulsions, having the droplet size distribution in the

1222

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

(PB0 , PB3 , and PB10 will refer to unweathered, 3-day, and 10-day evaporation) and ‘‘Ekofisk’’ oil. One can easily include the effect of weathering and emulsification into the model by adopting one of the available empirical formulae qo ¼ qo ðtÞ, mo ¼ mo ðtÞ, and r ¼ rðtÞ (see reviews of Huang (1983), Spaulding (1988), ASCE (1996) and Reed et al. (1999)), leading to the unsteady-state versions of Eq. (10). Since the present paper focuses only on vertical kinetics of the droplets, the above phenomena should be considered separately, if needed. Fig. 2. The normalized oil droplet size spectrum, and volume of droplets, versus the normalized droplet radius.

form of Eq. (7). Fig. 2 shows the normalized spectra, N^ ðrÞ=N ðrmax Þ, and the total volume of similar droplets, N^ ðrÞvðrÞ=N ðrmax Þvðrmax Þ, versus the dimensionless droplet radius, r=rmax , for s ¼ 2:3, V ¼ 106 m3 . According to Fig. 2, the most amount of oil is contained in a relatively small number of large droplets. It is important to have an accurate estimation of the radii rmin and rmax in order to calculate the volume of the entrained oil correctly. For large droplets, break-up may be attributed to the dynamic pressure force of turbulent flows (Hinze, 1955, hereafter HZ). For small droplets, viscous shear plays a dominant role in the break-up (Li and Garrett, 1998, hereafter LG). Based on the scale analyses of Hinze, Li and Garrett, and experiments of Delvigne and Sweeney (1988, hereafter DS), the maximum size is given by 8 ðHZÞ 0:363ðr=qÞ3=5 e2=5 ; > > > 9 > 3 1 > 5 0 0 3=8 > > 2:6  10 ðr=mqÞ e ðq m Þ ; > > > > > > < = for rmax > g=2 rmax ¼ ; ðLGÞ > 1:9  102 ðr=mqÞm1=2 e1=2 ðq0 m0 Þ1=8 ; > > > > > > > > ; > > for rmax < g=2 > > : : ðDSÞ 1818e0:500:1 m0:340:05 o ð10Þ Here, m0 ¼ mo =m, where mo and m are the kinematic viscosity of oil and water; q0 ¼ qo =q, where qo is the oil density; e is the mean energy dissipation rate; r is the oil–water interfacial surface tension coefficient; g ¼ ðm3 =eÞ1=4 is the Kolmogorov micro-length scale. Eq. (10, DS) was fine-tuned by Spaulding et al. (1992). Apart from the discrepancy in the exponent, the major difference above is that the formulae of Hinze (Eq. (10, HZ)) and Li and Garrett (Eq. (10, LG)) include the surface tension coefficient r and oil density qo , but the formula (10, DS) does not. One explanation (Delvigne and Sweeney, 1988) is that the measurements were performed with a narrow range of oil types, and it was not possible to determine the dependence of rmax on r or qo . The oil used in the experiments included ‘‘Prudhoe Bay’’

4. Oil droplet kinetics To define the droplet kinetics in stormy conditions it is necessary to account for major phenomena taking place during wave breaking. 4.1. Oil droplet entrainment The rate of near-surface oil mixing in open seas depends strongly on the energy of the breaking waves. It is estimated (Lamarre and Melville, 1991) that during wave breaking up to 50% of the dissipated wave energy is expended for air bubbles to penetrate in the water column against buoyancy. Similarly, it can be assumed that during the wave breaking a certain part of the dissipated wave energy, kb dE=dt, is expended to entrain the oil droplets from the slick into the water column. The coefficient kb may be evaluated from experiments. Currently, it is reasonable to assume that kb  0:3–0:5 (Lamarre and Melville, 1991). Using the relationships (4) and (6), the rate of oil entrainment from slick to the mixing layer, can be scaled as kow ¼ 

kb dE kb xcH ¼ : 16aLow qgzm Low dt

ð11Þ

Here, kow is the entrainment rate (s1 ); Low is the vertical length-scale parameter (m), depending on type of the breaking wave. The central term of the relationship (11) can be interpreted as the ratio of the wave energy flux, ðkb =qÞdE=dt, and the gravity force, gzm . The negative sign is taken to depict the oil mass gain in the mixing layer during the wave energy loss. Utilizing Eq. (3), Eq. (11) can be rewritten with a wind speed dependence. In this case, for a certain wave environment, the entrainment rate acquires a dependence kow  U 2 as in the model of Audunson (1979). The oil droplet mass transfer from the slick into the mixing layer is prescribed by the kinetic equation dMe ¼ kow Ms : dt

ð12Þ

Here, Ms and Me are the mass of oil in the slick and the mixing layer, respectively (both per unit surface area A).

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

1223

Total mass of oil in the slick and the mixing layer (per area A) is Ma ¼ Ms þ Me . At the instant of spill, t0 , Ma is equal to the mass of the spilled oil, M0 .

dMe ¼ kwo ðMe  Mes Þ: dt

4.2. Oil droplet resurfacing

kwo ¼ wðrl Þ=Lwo ;

Due to the buoyancy, the droplets may rise back to the water surface with the terminal velocity (Raj, 1977)

where Lwo is the vertical length-scale parameter (m), depending primarily on oil buoyancy, and vertical components of current velocity and diffusion.

wðrÞ ¼ kw r p :

Here, p ¼ 2 and kw ¼ 2gð1  q Þ=9m for Re < 50; p ¼ 0:5 and kw ¼ ð16gð1  q0 Þ=3Þ1=2 otherwise; Re ¼ 2rw=m is the droplet Reynolds number. This division of the rise velocity according the Reynolds number is due to application of a Stokes law for small droplets and a Reynolds law for large ones. Eq. (13) shows that larger droplets would rise more rapidly, whereas smaller droplets may remain in the water column for a longer time, increasing the chance to be dispersed by currents. It is assumed that the larger droplets from the droplet spectra (within the range rc < r < rmax ) contribute to the oil resurfacing, and droplets smaller than some threshold radius, rc , will permanently remain in the water column. Generally, the value rc is a function of water column state and oil characteristics; however, a typical empirical value of 50 lm is used sometimes for the simplification (Delvigne and Sweeney, 1988; Forrester, 1971; Proctor et al., 1994). Calculating a mean terminal velocity of the rising p droplets using relationship wðrl Þ ¼ 0:5kw ðrmax þ rcp Þ, the 1=p respective ‘‘mean’’ radius follows from rl ¼ ½wðrl Þ=kw . Alternatively, one can use the simplified estimation rl ¼ 0:5rmax (Spaulding et al., 1992). The characteristic radius of the small droplets is taken as rs ¼ rmin . The characteristic radii, rchar ¼ ðrs ; rc ; rl Þ, depend on the physical–chemical properties of the oil and water, and the mixing energy. Any commonly accepted functions rchar ¼ rchar ðe; r; m0 ; q0 Þ, can be substituted into the model. The volume of oil that may resurface, Vl , and the volume that remains submerged, Vs , (both per unit volume of the mixing layer) can be derived using Eq. (8) as 3s Vl ¼ 4pN 0 ðrmax  rc3s Þ=ð3  sÞ ¼ Bl V ; 3s Vs ¼ 4pN 0 ðrc3s  rmin Þ=ð3  sÞ ¼ Bs V :

ð14Þ

Here Bl and Bs ð¼ 1  Bl Þ are the fractions of the resurfaced and submerged oil, respectively, given by 3s rmax  rc3s ; 3s 3s rmax  rmin

Here, kwo is the resurface rate (s1 ), defined by ð17Þ

ð13Þ 0

Bl ¼

ð16Þ

Bs ¼

3s rc3s  rmin : 3s 3s rmax  rmin

ð15Þ

If one defines the highest mass of the droplets in the mixing layer (per unit area A) during the oil entrainment (storm is on) as Me ¼ qo Vzm A, then when the mixing force of the breaking waves is discontinued (storm is over), and the large droplets are resurfaced, the remaining oil mass is Mes ¼ Bs Me . The kinetics of the droplet resurfacing can be described using the first-order equation

4.3. Governing equations One can combine the two competing processes: the downward directed mixing, Eq. (12), and the upward oil buoyancy, Eq. (16), to yield dMe ¼ kow ðMa  Me Þ  kwo ðMe  Bs Me Þ: ð18Þ dt According to a definition, Me is the highest mass of droplets in the mixing layer for the specified wave parameters and oil properties. This is also a quasi-equilibrium value, because the downward mixing flux and the upward buoyancy flux are equalized. The value Me can be found by prescribing the extremum condition, dMe =dt ¼ 0, to the Eq. (18), yielding Me ¼

K Ma ; 1þK

ð19Þ

where the dimensionless ‘‘mixing factor’’, K, is defined as K¼

kow kb Lwo cxH ¼ : kwo Bl 16a Low wðrl ÞBl

ð20Þ

Using Eq. (19), and introducing the oil mass exchange rate coefficient as K¼

ðkow þ kwo Þkwo Bl ; kow þ kwo Bl

ð21Þ

the Eq. (18) is simplified to dMe ¼ KðKMs  Me Þ: ð22Þ dt Eq. (22) shows that the mixing factor K plays a role of oil partition coefficient between the slick and the mixing layer. To maintain a mass balance in the system of ‘‘slick–mixing layer’’, oil mass in the slick has to follow the kinetics dMs ¼ KðKMs  Me Þ: ð23Þ dt Derived Eqs. (22) and (23) are the main result of the paper. They offer a simple and consistent way to compute vertical exchange of oil droplets mass, if all other losses are neglected. An example of inclusion of losses for evaporation and vertical diffusion is given in the next section. To be compatible with a general oil spill model, the kinetic Eqs. (22) and (23) can be rewritten in terms of

1224

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

the oil slick thickness, h, and the oil droplet concentration, Ce , using the relationships h ¼ Ms =qo A and Ce ¼ Me =zm A. The combined model is expected to be an efficient tool for simulation of the spilled oil fate. It can be used to assess the slick combating techniques in order to select the most efficient countermeasure. 4.4. Oil dispersion and evaporation The total oil mass Ma may change with time due to losses for hydrolysis, evaporation, photolysis, oxidation, biodegradation, etc. These phenomena have to be taken into account in a general oil spill model. In this paper, simplified expressions for evaporation and vertical diffusion are utilised for the illustrative purpose only. Other losses can be easily included in a similar fashion, if required. The turbulent mixing force might propel the smaller oil droplets out of the mixing layer deeper into the water column, reducing their opportunity to resurface back. Denoting the concentration of the droplets at the depth Ld as Ced (see Fig. 1), kinetics of diffusion can be approximated with dMe D ¼  2 ðMe  Med Þ; Ld dt

ð24Þ

where Med ¼ Azm Ced ; D is the vertical diffusion coefficient (m2 /s); Ld is the vertical length-scale (m). The typical values for the near-surface diffusion coefficient are about D  0:001–0.01 m2 /s (Thorpe, 1984; Proctor et al., 1994). The value Ld depends on wave and hydrodynamics conditions and oil parameters, and generally it can be estimated from the measurements (Forrester, 1971; Thorpe, 1984) as being around 10–20 m. Dispersants may increase the length-scale to the higher values. The similar magnitude is expected for the parameter Lwo . In the majority cases one can assume zero concentration of the droplets at the depth below Ld , except for the shallow waters or pycnocline, requiring different boundary conditions. Oil slicks are subjected to rapid evaporation. In this paper a simplified first order kinetics model is employed dMs ¼ bMs ; ð25Þ dt where b is the evaporation rate (s1 ). One can easily substitute any other kinetic terms, using a comprehensive review of the available oil evaporation models (Fingas, 1995). 5. Model verification

a measure of the entrained oil mass. The mixing factor can be related to the empirical oil entrainment expression of Delvigne and co-workers (1988, 1994). The analysis of the Eq. (20) is complicated by a multiple choice of expressions for the terminal velocity, w, and the maximal droplet size, rmax . After generalization of 0 0 n0 the Eqs. (10) as rmax  rl mm o qo , the Eq. (13) assumes the form wðrl Þ  rl mmo qno :

Here, the exponents are ðl; m; nÞ ¼ pðl0 ; m0 ; n0 Þ, where the values of l0 , m0 , and n0 depend on the chosen droplet formation hypothesis, and the parameter p depends on the droplet Reynolds number. The exponents l, m, n and p, being chosen to match conditions of the Delvigne and Sweeney (1988) experiments, are summarized in the Table 1. The Eq. (10, DS) for maximal droplet size leads to m ¼ 0:68 for Re 6 50, and m ¼ 0:17 for Re > 50. Parameters l and n can be estimated using Eqs. (10, HZ and LG). For small Re, of the two formulae (10, LG) the one for rmax P g=2 is chosen to match the wave breaking conditions of the Delvigne and co-workers’ experiments, giving l ¼ 6 and n ¼ 0:75. For large Re, Hinze’s definition leads to l ¼ 0:3 and n ¼ 0. Substituting the Eq. (26), and Bl ¼ Vl =V (derived from the Eq. (14)) into Eq. (20), one gets 8 cHV > > < r6 m0:680:1 q0:75 V for Re < 50; l o o ð27Þ K cHV > > : 0:3 0:170:025 for Re > 50: r mo Vl In this form, the mixing factor can be compared with the similar empirical quantity introduced by Delvigne and 0:7 Sweeney (1988) as QDS ðro Þ ¼ Ew0:570:06 m1 o ðr0 drÞ. Here, QDS (kg=m2 ) was defined as ‘‘the entrained mass of the oil droplets, having radius in an interval dr around r0 , per unit surface area per breaking event’’. The empirical relationship can be converted into the dimensionless form: Z A Ew0:570:06 rmax 0:7 H 1:140:12 V KDS ¼ r dr  ; ð28Þ M0 mo mo rmin where the Eqs. (2) and (8) were used in the derivation. Spaulding et al. (1992, hereafter OM), using the same Delvigne and Sweeney (1988) data fitted somewhat difTable 1 Estimation of parameters using different hypotheses on the maximal droplet size LG

5.1. Model comparison with empirical formulae The mixing factor K, Eq. (20), being a ratio of the downward wave mixing and upward buoyancy forces, is

ð26Þ

l m n

DS

HZ

(small Re)

Re < 50

Re > 50

(large Re)

6.00 0.75 0.75

– 0.68  0.1 –

– 0.17  0.025 –

0.3 – –

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

ferent power relationship QOM  mo0:5 , leading to KOM  H 1:140:12 V m0:5 . o The theoretically derived mixing factor, given by Eq. (20) or (27), considers more oil properties, than the respective empirical relationship (28). Among newly included parameters are: the interfacial surface tension coefficient r and oil density qo . Rest parameters, such as wave height H and oil viscosity mo are in a good agreement with the empirical Eq. (28). Note that the combination of the relationships K  w1 (Eq. (20)) and w  mm o (Eq. (26)) lead to the observed correlation K  mm o . The derived exponent value m ¼ 0:68 for small Re is similar to that in the Delvigne and Sweeney (1988) experiments. More recent data of Delvigne and Hulsen (1994, hereafter DH) did not confirm the correlation for low-viscous types of oil (m0 < 100), however for the higher viscosities the entrainment decreased much faster with increasing viscosity. One explanation for the data discrepancy follows from the presented theoretical model: smaller droplets originating from the lower viscosity oil may not resurface fast enough to be sampled, thus masking the dependence of the mixing factor K on the terminal velocity w (and mo , respectively). Higher viscosity oil, in contrast, generate much larger droplets which may resurface at much shorter time-scale, leading to the strong correlation of K and mo . Alternative explanation is that some of the new experiments may be associated with large droplet Reynolds number, indicating a weak dependence K and mo (m ¼ 0:17, in Eq. (27, second term)); or droplets formation may followed the Hinze regime, showing no dependence of the mixing factor and the oil viscosity. Eq. (27) predict that the oil droplet mixing intensity is inversely proportional to the surface tension coefficient r. This phenomenon is routinely used in the practice of slick combating to facilitate the slick breaking-up and oil dispersion by means of the dispersant spraying. Similar correlation was obtained from the measurements of Mackay et al. (1984) and postulated in the empirical model of Mackay et al. (1980). Delvigne and co-workers (1988, 1994) did not find any reliable correlation of the maximal droplet size rmax with the oil density qo , as it is done theoretically (Li and Garrett, 1998). To assess an influence of qo on the terminal velocity wðrmax Þ, the Eq. (13) is computed for Re 6 50 using the values rmax given by Li and Garrett (1998), Hinze (1955) and Delvigne and Sweeney (1988). All the approaches reveal almost identical dependence of w versus qo (Fig. 3). Eq. (5) show that, depending on the wave origin and type of breaking, the damping coefficient could be c  H 0 for the swells decay, or c  H 0:5 for the white-capping. Generalization of the last two expressions for an arbitrary wave type yields the relationship K  cH  H 1:250:25 . Similar correlations were observed by Delvigne and co-workers (1988, 1994) in the large-scale

1225

Fig. 3. Calculated terminal velocity for the maximal droplet sizes using different droplet break-up hypotheses.

wave flume experiments (KDS  H 1:14 ) and in the smallscale experiments (KDH  H 1:4 ). Both experimentally fitted exponents fall into the theoretically predicted interval ð1; 1:5Þ. 5.2. Model sensitivity to the parameters variation To explore the model sensitivity to the variation of oil density and viscosity, the expression (19) for quasiequilibrium droplet mass is rewritten in the dimensionless form as Me =M0 ¼ K=ð1 þ KÞ. The following parameters are used in the test: q ¼ 1000 kg=m3 ; m ¼ 106 m2 =s; rmin ¼ 1 lm; 2

A¼1m ;

qo ¼ 800–1000 kg=m3 ; mo ¼ 106 –103 m2 =s;

rc ¼ 50 lm; Low ¼ 1 m;

s ¼ 2:3;

h0 ¼ 0:001 m;

Ld ; Lwo ¼ 20 m;

3

e ¼ 1000 J=m s; x ¼ 1:91 l=s; H ¼ 2 m; c ¼ 0:0013; kb ¼ 0:3; a ¼ 1:5: According to the scenario, oil is instantly spilled on the water surface, forming a slick with the initial thickness h0 . The actual area of the slick and the total mass of the oil are not important for the study, because results are analysed using normalized quantities. To compare the computations against the Delvigne and Sweeney (1988) data, the ‘‘Prudhoe Bay’’ and ‘‘Ekofisk’’ oils are used along with the authors’ definition of rmax Eq. (10, DS). The computed mass of droplets in the mixing layer is shown in the Figs. 4 and 5 as functions of the oil/water density ratio, q0 , and the viscosity ratio, m0 , respectively. The calculations clearly demonstrated that the droplet mixing rate is proportional to the oil density, and inversely proportional to the oil viscosity. For the case of neutral buoyancy, q0 ¼ 1, the surface oil rapidly mixes with water, as expected (Fig. 4). Fig. 5 exhibits decrease of slope of the droplet entrainment curve from m0:68 to o m0:17 , when the droplet Reynolds number exceeds the o

1226

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

Fig. 4. Computed oil droplet mass in the mixing layer versus the oil/ water density ratio for the different oil/water viscosity ratios.

Fig. 6. Computed oil droplet mass in the mixing layer versus the oil/ water viscosity ratio using different definitions for the maximal droplet size.

expected to be dominant at the final stage or after chemical dispersants application. Speculating that the experiments of Delvigne and Sweeney (1988) might have depicted the both features, the droplet mass is computed for the selected oil types using the suggested formulae for rmax Eqs. (10). The results, plotted in Fig. 6, indeed suggest that the empirical Eq. (10, DS) exhibits features of the both hypotheses.

6. Efficiency of chemical dispersants in combating of oil slicks

Fig. 5. Computed oil droplet mass in the mixing layer versus the oil/ water viscosity ratio for the different oil/water density ratios.

critical value Re ¼ 50. Similar reduction in the slope of the oil entrainment curve was shown by Delvigne and Sweeney (1988, p. 22), in the transition from lower viscosity oil to higher viscosity oil. If one requires a uniform exponent for the arbitrary droplet Reynolds number, then the fit QOM  m0:5 is within the theoretio cally predicted range (Fig. 5). On the whole, Figs. 4 and 5 qualitatively agree well with the empirical data of Delvigne and Sweeney (1988, pp. 21–22). More detailed quantitative comparison is not possible at this stage, since the authors did not provide the reference scale for the empirical plots. The two hypotheses on the dominant mechanism of oil droplet formation, namely the pressure force and the viscous shear yield different estimates for rmax . Li and Garrett (1998) suggested that the transition between the regimes is not expected to be abrupt. The Hinze mechanism must be more important at the early stage of the wave breaking, and the Li and Garrett phenomenon is

One of the major techniques for the slick combating is the application of dispersants, which are able to reduce the oil–water interfacial tension coefficient from the usual 0.02–0.03 N/m to possibly less than 0.001 N/m (Mackay et al., 1984), thus facilitating the droplets to shear from the slick. As the result of complex physical– chemical processes between the oil, water and dispersant, the fractured droplets become so small, that their rising velocity cannot overcome the mixing force of turbulence near the surface. To be efficient, the dispersant has to be spread over the slick prior to a storm that could provide energy for the treated slick to break and for the droplets to disperse. It has to be noted, that dispersants are not equally effective for all kind of oil, and they are not useful for the high-viscosity or weathered oil. To demonstrate the phenomenon using the model, Eqs. (22) and (23), an instant oil spill on a water surface is assumed. Oil losses due to the dispersion into a deep water and evaporation are considered using Eqs. (24) and (25). A storm, started 5 h after the spillage, supposedly lasts for one whole day (24 h). The three cases of a dispersant application are considered: (1) 2 h prior to the storm; (2) half a day (12 h) after the storm arriving;

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

and (3) 3 h after the storm is over. The following parameters have been used for the simulation:

1227

Actual percentage of evaporated and residual fractions depends strongly on type of oil. In the test case it is assumed that only 30% of the oil mass is available for the evaporation, and 10% is going to form the ‘‘chocolate mousse’’, persistently floating near the surface. Weathering phenomenon is simulated by linear increase of oil viscosity from the initial value 105 m2 /s to the intermediate value 102 m2 /s within first five days after the spill. A dispersant presumably changes the surface tension coefficient from the value r ¼ 0:01 N/m down to the 0.001 N/m. The maximal droplet size is calculated using Li and Garrett (1998) formulae, Eq. (10, LG).

During the storm the wave height increases from the flat value H ¼ 0:1 m up to the 3 m, and the energy dissipation rate (according to recommendations of Li and Garrett, 1998) from e ¼ 0:1 to 10 m2 /s3 , respectively. The actual area of the slick and the total mass of the spill are not important for the study, because the oil kinetics is analysed in terms of fraction of the total spill. Fig. 7 presents the first case, where the dispersant is applied 2 h prior to the storm. The computation exhibits the known phenomenon, that application of the chemical alone, without wave mixing, does not increase rate of oil dispersion. In contrast, after the storm started, a large fraction of surface oil (up to 59%) is effectively dispersed. In the second case, when the dispersant is applied during the storm (Fig. 8), about 40% of initial oil mass would be permanently dispersed in the deeper water. If the dispersant is applied after the storm, as in the third case (Fig. 9), only a few percent of oil can be dispersed. Figs. 7–9 depict rapid increase of the oil entrainment rate kow during the storm. The same time, the resurface rate kwo falls sharply when the dispersant is added, and rebounds back after the storm is over. Gradual decrease of the mixing factor K is caused by the oil viscosity increase due to weathering. It is clear from the results, that the chemical dispersants and the mixing

Fig. 7. Computation of oil droplet dispersion by a storm. The dispersant added 2 h prior to the storm.

Fig. 8. Computation of oil droplet dispersion by a storm. The dispersant added 12 h after the storm arriving.

q ¼ 1000 kg=m3 ; m ¼ 106 m2 =s;

qo ¼ 800 kg=m3 ; mo ¼ 105 –102 m2 =s;

r ¼ 0:001–0:01 N=m; s ¼ 2:3;

h0 ¼ 0:001 m;

Low ¼ 1 m;

a ¼ 1:5;

D ¼ 0:002–0:02 m2 =s; c ¼ 0:004;

rmin ¼ 1 lm;

rc ¼ 50 lm;

Ld ; Lwo ¼ 20 m; e ¼ 0:1–10 m2 =s3 ; x ¼ 1:7 s1 ;

kb ¼ 0:3; b ¼ 2:0  10

5

H ¼ 0:1–3 m; 1

s :

1228

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229

Fig. 9. Computation of oil droplet dispersion by a storm. The dispersant added 3 h after the storm is over.

due to breaking waves are equally important factors for efficient combating of oil slicks. The dispersants, being applied prior to the storm, significantly increase production of small droplets, thus enhancing the rate of oil dispersion. The higher the energy of the waves, the more intensive is the breaking of the slick and the droplets dispersion. Increase of oil viscosity due to weathering causes larger droplets (or even blobs) to shear from the slick, thus making the dispersant application and wave mixing inefficient. The simulation results depict the typical observed behaviour of oil spills. In calm waters, oil practically does not mix with water, sustaining in the slick for a long time. In rough seas, a high percentage of the oil can still be located within the slick. The mixing is drastically improved when the chemicals are applied before the storm, leading to the intensive dispersion of the spilled oil.

7. Conclusions In this paper, a new theoretical model of oil droplet vertical mixing by breaking waves is developed. Considering dominant forces affecting the droplet formation and vertical distribution, the model describes the droplet

mass kinetics in the upper part of the water column. The most important characteristics of oil, waves and water column are conveniently combined into the single ‘‘mixing factor’’ K, (Eq. (20) or (27)) that is shown to be consistent with the previously empirically obtained ‘‘droplet entrainment mass’’ Q of Delvigne and Sweeney (1988). Having similarity with the experimental formula on oil viscosity and wave height dependence, the mixing factor includes additionally the interfacial surface tension coefficient and oil density. The quantity K (as well as Q) may represent a special case of oil quasi-equilibrium distribution between the slick and the water column, leading to a simple formula (19) for computing the droplet mass in the mixing layer. However, majority of cases would require usage of the entire kinetic model, Eqs. (22) and (23), to allow for oil redistribution between the slick and the water column according to the changes in physical–chemical and hydro-meteorological characteristics. These new features extend the model applicability to a wider range of oil parameters and environmental conditions. A series of test computations have been performed to ensure adequate sensitivity of the model to parameter variation, and to verify the model using available empirical formulae. The model can be used for a quick assessment of oil mixing rate under breaking waves. However, real significance of the research has become evident when the authors worked on the development of Multiphase Oil Spill Model (Tkalich et al., 1999; Tkalich and Chan, 2001), which required an accurate and physically relevant quantification of oil droplet mixing in the upper part of the water column. In the context of a general oil spill model, the droplet mixing kinetics provides a convenient parameterisation of wave breaking phenomena at the oil slick–water column interface.

References ADIOS, 1994. Automated Data Inquiry for Oil Spills, User’s Manual. NOAA/Hazardous Materials Response and Assessment Division, Seattle Washington. Aravamudan, K., Raj, P., Newman, E., Tucker, W., 1980. Breakup of oil on Rough Seas––simplified models and step-by-step calculations. Report no. CG-D-69-79 for US Department of Transportation, US Coast Guard Office of Research and Development, Washington DC 20590. ASCE Task Committee, 1996. State-of-the-art review of modeling transport and fate of oil spills. Journal of Hydraulic Engineering 122, 594–609. Audunson, T., 1979. Fate of oil spills on the Norwegian continental shelf. In: Proceedings of the 1979 Oil Spill Conference, Los Angeles, California. API, Washington DC, p. 675. Blaikley, D.R., Dietzel, F.F.L., Glass, A.W., van Kleef, P.J., 1977. Sliktrak––a computer simulation of offshore oil spills, Cleanup, Effects and Associated Costs. In: Proceedings of the 1977 Oil Spill Conference. API, Washington DC, pp. 45–53. Bouwmeester, R.J.B., Wallace, R.B., 1986. Oil entrainment by breaking waves. In: Proceedings of the Ninth Arctic Marine Oil

P. Tkalich, E.S. Chan / Marine Pollution Bulletin 44 (2002) 1219–1229 Spill Program Technical Seminar. Environment Canada, pp. 39– 49. Daling, P.S., Aamo, O.M., Lewis, A., Strøm-Kristiansen, T., 1997. SINTEF/IKU oil-weathering model: predicting oil properties at sea. In: Proceedings of the 1997 Oil Spill Conference. API, Washington DC, pp. 297–307. Delvigne, G.A.L., Hulsen, L.J.M., 1994. Simplified laboratory measurements of oil dispersion coefficient––application in computations of natural oil dispersion. In: Proceedings of the 17th Arctic and Marine Oil Spill Program Technical Seminar. Environment Canada, pp. 173–187. Delvigne, G.A.L., Sweeney, C.E., 1988. Natural dispersion of oil. Oil and Chemical Pollution 4, 281–310. Donelan, M.A., 1978. Whitecaps and momentum transfer. In: Proceedings of the NATO Conference on turbulent fluxes through the sea surface, wave dynamics, and prediction. Plenum, NY, pp. 273–287. Elliott, A.J., 1986. Shear diffusion and the spread of oil in the surface layers of the North Sea. Deutsche Hydrographische Zeitschrift 39, 113–137. Fingas, M.F., 1995. A literature review of the physics and predictive modelling of oil spill evaporation. Journal of Hazardous Materials 42, 157–175. Forrester, W.D., 1971. Distribution of suspended oil particles following the grounding of the tanker Arrow. Journal of Marine Research 29, 151–170. Hasselmann, K., 1974. On spectral dissipation of ocean waves due to whitecapping. Boundary-Layer Meteorology 6, 107–127. Hinze, J.O., 1955. Fundamentals of the hydrodynamic mechanism of splitting in dispersion processes. AICHE Journal 1, 289–295. Huang, J.C., 1983. A review of the state-of-the-art of oil spill fate/ behaviour models. In: Proceedings of the 1983 Oil Spill Conference. API, Washington DC, pp. 313–322. Johansen, Ø., 1982. Dispersion of oil from drifting slicks. Spill Technology Newsletter, November–December, pp. 134–149. Komen, G.J., Hasselmann, S., Hasselmann, K., 1984. On the existence of a fully developed wind–sea spectrum. Journal of Physical Oceanography 14, 1271–1285. Komen, G.J., Cavaleri, L., Donelan, M., Hasselmann, K., Hasselmann, S., Janssen, P.A.E.M., 1994. Dynamics and Modelling of Ocean Waves. Cambridge University Press, New York, p. 532. Lamarre, E., Melville, W.K., 1991. Air entrainment and dissipation in breaking waves. Nature 351, 469–472. Li, M., Garrett, C., 1998. The relationship between oil droplet size and upper ocean turbulence. Marine Pollution Bulletin 36, 961–970. Mackay, D., Paterson, S., Trudel, K., 1980. A mathematical model of oil spill behavior. Report EE7, Environment Canada, Ottawa. Mackay, D., Chau, A., Hossain, K., Bobra, M., 1984. Measurements and prediction of the effectiveness of oil spill chemical dispersants. In: Allen, T.E. (Ed.), Oil Spill Chemical Dispersants: Research,

1229

Experience, and Recommendations, STP 840. American Society for Testing and Materials, Philadelphia, pp. 38–54. Monahan, E.C., 1971. Oceanic whitecaps. Journal of Physical Oceanography 1, 139–144. Neumann, G., Pierson, W.J., 1966. In: Principles of Physical Oceanography. Prentice-Hall, Englewood Cliffs, NJ, p. 545. Proctor, R., Elliott, A.J., Flather, R.A., 1994. Forecast and hindcast simulations of the Braer oil spill. Marine Pollution Bulletin 28, 219–229. Raj, P.P.K., 1977. Theoretical study to determine the sea state limit for the survival of oil slicks on the ocean. Report CG-D-90-77, US Coast Guard, Washington DC. Reed, M., Turner, C., Odulo, A., 1994. The role of wind and emulsification in modelling oil spill and surface drifter trajectories. Spill Science and Technology Bulletin 12, 143–157. Reed, M., French, D., Rines, H., Rye, H., 1995. A three-dimensional oil and chemical spill model for environmental impact assessment. In: Proceedings of the 1995 Oil Spill Conference, Long Beach, California. API, Washington DC, p. 61. Reed, M., Johansen, Ø., Brandvik, P.J., Daling, P., Lewis, A., Fiocco, R., Mackay, D., Prentki, R., 1999. Oil spill modelling toward the close of the 20th century: overview of the state of the art. Spill Science and Technology Bulletin 5, 3–16. Shore Protection Manual, 1984. Prepared for Department of the Army, US Army Corps of Engineers. Dept. of the Army, Coastal Engineering Research Center, Washington DC, 1. Spaulding, M.L., Saila, S.B., Reed, M., Lorda, E., Walker, H., Isaji, T., Swanson, J., Anderson, E., 1982. Assessing the impact of oil spills on a commercial fishery. Final report to Minerals Management Service. Contract no. AA851-CTO-75; NTIS no. PB83149104. Spaulding, M.L., 1988. A state-of-the-art review of oil spill trajectory and fate modelling. Oil and Chemical Pollution 4, 39–55. Spaulding, M.L., Howlett, E., Anderson, E., Jayko, K., 1992. OILMAP: a global approach to spill modelling. In: Proceedings of the 15th Arctic and Marine Oil Spill Program Technical Seminar. Environment Canada, pp. 15–21. The SWAMP Group, 1985. Ocean Wave Modeling. Plenum Press, NewYork, p. 256. Thorpe, S.A., 1984. On the determination of Kv in the near-surface ocean from acoustic measurements of bubbles. Journal of Physical Oceanography 14, 855–863. Tkalich, P., Huda, K., Gin, K., 1999. A multiphase model of oil spill dynamics. In: Proceedings of XXVIII IAHR Congress. Graz, Austria. Tkalich, P., Chan, E.S., 2001. A High Resolution Oil Spill Model. In: Proceedings of International Conference on Port and Maritime R&D and Technology. Singapore, pp. 647–652. Young, I.R., 1999. Wind Generated Ocean Waves. Elsevier, Amsterdam, New York, p. 288.