A numerical study on the evolution and structure of a stress-driven free-surface turbulent shear flow

30  Download (0)

全文

(1)

 2005 Cambridge University Press doi:10.1017/S0022112005007044 Printed in the United Kingdom

A numerical study on the evolution and

structure of a stress-driven free-surface

turbulent shear flow

By W U -T I N G T S A I1,2, S H I -M I N G C H E N2

A N D C H I N -H O H M O E N G3

1Institute of Hydrological Sciences, National Central University, Taoyuan 32001, Taiwan 2Department of Civil Engineering, National Chiao Tung University, Hsinchu 30050, Taiwan

3National Center for Atmospheric Research, Boulder, Colorado 80307, USA

(Received7 October 2004 and in revised form 20 June 2005)

Turbulent shear flow beneath a flat free surface driven by a surface stress is simulated numerically to gain a better understanding of the hydrodynamic processes governing the scalar transfer across the air–water interface. The simulation is posed to mimic the subsequent development of a wind-driven shear layer as in a previous experiment except that the initiation of the surface waves is inhibited, thus focusing on the boundary effect of the stress-imposed surface on the underlying turbulent flow and vice versa. Despite the idealizations inherent in conducting the simulation, the computed flow exhibits the major surface features, qualitatively similar to those that appear in the laboratory and field experiments. Two distinct surface signatures, namely elongated high-speed streaks and localized low-speed spots, are observed in the simulated flow. Including temperature as a passive tracer and describing an upward heat flux at the surface, we obtain high-speed streaks that are colder and low-speed spots that are warmer than the surrounding regions. The high-low-speed streaks, arranged with somewhat equal cross-spacing of centimetres scale, are formed by an array of streamwise jets within the viscous sublayer immediately next to the surface. Beneath the streaks, counter-rotating streamwise vortex pairs are observed among other prevailing elongated vortices. However, they are significantly shorter in length and more irregular than their corresponding high-speed streaks at the surface. Accompanying the more organized high-speed streaks, localized regions of low streamwise velocity emerge randomly on the surface. These low-speed spots are attributed to strong upwelling flows which disrupt the viscous sublayer and also bring up the submerged fluids of low streamwise velocity. The occasional interruptions of the streamwise high-speed jets by the upwelling flows account for bifurcation or dislocation of the surface streaks. Statistics of the turbulence are presented and their implications for the formation of the flow structures are discussed.

1. Introduction

The layer of water adjoining and including the air–sea interface controls the exchange of slightly soluble gases, such as carbon dioxide, between the atmosphere and the ocean. Such a liquid-phase controlled process is determined by the near-surface turbulence in the ocean side, of which the dominant forcing is the wind stress on the interface. The energetic, coherent vortical motions within the aqueous turbulent boundary layer, or known as the surface renewal eddies, therefore account for the

(2)

effectiveness of transport across the air–water interface. Through interacting with the boundary, these coherent eddies also induced structural patterns on the water surface. A popular approach in recent years to estimate the air–water transfer rate is to infer the renewal events from remote-sensing images of water surface (e.g. Haussecker et al. 2002; Garbe, Schimpf & J¨ahne 2004; Schimpf, Garbe & J¨ahne 2004). The accuracy of the technique thus relies on relating quantitatively the small-scale surface signatures to the underlying turbulence processes.

The observation by Woodcock (1941) is among the pioneering experiments to study the small-scale motions on the water surface. He observed that when winds blow, high-speed surface water movements occur that subsequently form surface streaming in lines and streaks roughly parallel to the wind direction. The water in the streams moves down-wind much faster than that between the streams. The streaks are approximately a centimetre across and with a distance of the order of centimetres between the streaks. The subsequent visualization experiments of Gemmrich & Hasse (1992) and Kenny (1993) all observed similar small-scale surface streaks on natural water surfaces. These high-speed streaks, which are the most frequently appearing patterns of surface water movement at wind speeds of 2 or 3 cm s−1 or higher, were taken by Gemmrich & Hasse (1992) as indication of surface renewals that might enhance air–sea gas exchange. They interpreted the appearance of these surface streaks with a cross-spacing of the order of centimetres either as a result of wave breaking or as an indication of the existence of counter-rotating helical vortex rolls similar to those in Langmuir circulations, which have a typical spacing of metres to kilometres (Langmuir 1938; or see the reviews by Leibovich 1983 and Thorpe 2004). The observation of Kenney (1993), nevertheless, indicated the contrary. No evidence was found of fixed circulations beneath the streaky surface. In addition, no significant surface waves were observed and the streaks could appear on a calm water surface at low wind condition, indicating the absence of wave–current interaction for generating Langmuir circulations.

In a laboratory experiment by Melville, Shear & Veron (1998) (and also Veron & Melville 2001), predominant streaks aligned in the wind direction were also observed at a wind-driven water surface. The streaky surface undergoes several stages in its evolution which begins with a quiescent water surface at the onset of the wind. The appearance of the surface streaks was shown to follow closely the initial growth of the surface waves and the inception of the Langmuir circulations. (In another laboratory experiment by Caulliez, Ricci & Dupont 1998, transition of the initial laminar flow to a turbulent field was observed, which then leads to an explosive growth of the initial wind-generated ripples and the formation of surface streaks.) In addition to the surface streaks, local upwellings and divergence of the velocity field at the surface were also observed by Melville et al. (1998), which caused bifurcation or dislocation of the streaks.

The acting of wind stress on the water surface induces both surface waves and turbulent shear flow. It is possible that the coherent surface signatures observed in various experiments might result from several different mechanisms, including the interaction among surface waves, the interaction between waves and shear current, and also the development of the turbulent boundary layer itself. None of these mechanisms has been demonstrated to the exclusion of the others. The formation of a wind-driven shear current always accompanies the inception of the surface waves. Accordingly, it is difficult, if not impossible, to distinguish various mechanisms in the experiments, in particular those pertinent to surface waves. Attempts to separate the generation of surface waves from the wind-driven shear current by supplied

(3)

surfactants to inhibit the wave initiation (Faller & Perini 1984; Veron & Melville 2001), however, would also modify the sub-surface flow structures and, consequently, affect the induced surface signatures (McKenna & McGillis 2004; Tsai & Yue 1995; Tsai 1996).

In this study, the difficulty encountered in the experiments is resolved by numerical simulation of a canonical problem. Specifically, turbulent shear layer driven by a shear stress acting at a water-surface boundary is simulated numerically without the formation of the surface waves. Such a model problem has previously been studied for stress-driven oceanic mixed layers by large-eddy simulations (e.g. Skyllingstad & Denbo 1995; McWilliams, Sullivan & Moeng 1997; Sander et al. 2000), which require parameterization of the near-surface flow using the law-of-the-wall theory. The approach adopted here is a direct numerical simulation of the wind-driven turbulent shear flow using neither any turbulence closure to model the subgrid flow processes nor a ‘vortex force’ to initiate the mechanism of wave–current interaction as in these large-eddy simulations. Therefore, with direct numerical simulation we can resolve the flow scales down to that of the viscous sublayer and examine the detailed flow field immediately next to the surface. The present study thus attempts to deduce the possible correlation between the observed centimetre-scale surface signatures and the coherent flow processes induced by shear turbulence. The simulation also focuses on providing further insight into the detailed hydrodynamic mechanism, particularly the coherent flow structures that govern the process of scalar exchange across the air–sea interface.

The numerical aspects of the simulation, including the model formulation, the numerical method and the simulation implementation are described in § § 2, 3 and 4, respectively. Statistics of the turbulence are first discussed in § 5. The numerical results of the first-order mean profiles are compared with those derived from the surface renewal theory (Liu & Businger 1975) and the experimental measurements (Wu 1984). Turbulent energetics budgets of the present shear turbulence bounded by a stress-driven boundary, compared to those next to either a shear-free boundary or a no-slip wall, are discussed. Two-dimensional surface signatures of the velocity and scalar fields arising from the development of the turbulent boundary layer are then presented in § 6, and compared with those observed in the experiment of Melville et al. (1998). The underlying three-dimensional flow structures, which induce the surface signatures, are finally examined and discussed in§ 7. Some of the preliminary results from this investigation were reported in letters (Tsai 2001; Tsai et al. 2003).

2. Mathematical formulation

We consider a turbulent shear flow beneath a flat free surface at z = 0. The coordinate axes x = (x, y, z) are in the streamwise, spanwise and upwards vertical direction, respectively, with the corresponding velocities v = (u, v, w). The turbulent shear flow is driven by a shear stress τs acting at the boundary surface. The fluid is assumed to be incompressible and Newtonian, such that the velocity v and the dynamic pressure p are governed by the solenoidal condition and the Navier–Stokes equations: ∇ · v = 0, (2.1) v ∂t +v · ∇v = − 1 ρ∇p + ν∇ 2v, (2.2)

(4)

With the focus of this study being on the boundary effect of the stress-imposed surface on the underlying shear turbulence and vice versa, an idealized boundary of a flat free surface is considered (e.g. Walker, Leighton & Garza-Rios 1996; Calmet & Magnaudet 2003) in this study. By a flat free surface, what is meant is that the exerted stress is balanced by local surface tangential stress and that the surface is free to flow horizontally, but not vertically. A previous study (Tsai 2001) using more realistic free-surface boundary conditions (Borue, Orszag & Staroselsky 1995; Tsai 1998) indicates only a minor impact of the small-amplitude surface motions on both the surface and subsurface turbulent structures. At the surface, the tangential surface stress is balanced by the mean exerting shear stress, whereas the normal surface stress vanishes. Accordingly, the corresponding kinematic and dynamic conditions on the mean surface, z = 0, are

w= 0, ∂u ∂z = τs ρν, ∂v ∂z = 0, ∂p ∂z = 0. (2.3)

Note that in contrast to the Dirichlet pressure condition on a free surface with deformation (e.g. Tsai 1998), the normal-stress condition on a flat free surface results in a homogeneous Neumann pressure condition. With w = 0, the surface gravity waves are excluded from our simulations. Although the idealized free surface differs from that of the experiment of Melville et al. (1998) and also from the real ocean surface, the exclusion of surface waves allows us to isolate the effect of shear turbulence on small-scale signatures at the air–sea interface without the complicating effects of the mean-flow/surface-wave interaction. Our simulations also differ from that of a wall-bounded shear turbulent boundary layer (e.g. Mansour, Kim & Moin 1988; Moin & Mahesh 1998). In contrast to shear turbulence over a no-slip wall, our simulated flows allow for tangential motions on the surface in response to the driven surface stress and the underlying turbulence.

To study the effects of the near-surface flow processes on the transport of passive scalars, advection–diffusion equations governing the evolutions of the temperature field θ (x, t) and the gas concentration c(x, t) are also integrated in time with the flow simulation. The boundary condition at the surface, however, depends on the transport process considered. A given uniform heat flux Q gives rise to a Neumann condition for the temperature field at the surface as

∂θ(x, y, 0, t)

∂z =−

Q

ρcpνθ

, (2.4)

where cp and νθ are specific heat and thermal diffusivity, respectively. In our simulations, the temperature field is treated as a passive tracer; the buoyancy effect due to temperature fluctuations and the imposed surface heat flux does not modify the vertical momentum equation in (2.2).

In the case of transfer of sparingly soluble gas, the resistance is dominated by the subsurface aqueous flow and the atmosphere can be regarded as an infinite reservoir with a constant concentration. This results in a Dirichlet condition for the concentration field at the surface as

c(x, y, 0, t) = cs. (2.5)

3. Numerical method

The numerical method used in the present simulation resembles that in Tsai (1998) but with modifications to improve the near-surface resolution and to impose the

(5)

solvability condition for the pressure equation. A mixed scheme of the pseudo-spectral method in the two horizontal directions and second-order finite differencing in the vertical is adopted for approximation of the spatial-differential operators in the equations and boundary conditions. A stretched grid system is employed such that the discretization grids cluster when approaching the surface so as to allow proper resolution of the viscous sublayer adjacent to the upper surface. The vertical domain of physical coordinate −h 6 z 6 0 is discretized according to the transformation function zk=−h  1− tanh aζk tanh a  , 06 ζk 6 1, (3.1)

where the grids of the numerical coordinate ζ are uniformly distributed. The degree of cluster is determined by the parameter a. A typical value used in the present computations is a = 1.8417 (Gavrilakis 1992).

A low-storage Runge–Kutta scheme (Spalart, Moser & Rogers 1991) is adopted for temporal integration of the Navier–Stokes equation (2.2) for the velocity field and the transport equations for the scalar fields. The pressure field is first solved at each interval step in the Runge–Kutta time integration before updating towards the next time step for the new velocities. The Poisson equation for the pressure field is obtained by taking the divergence of the temporal-discretized Navier–Stokes equation and applying the solenoidal constraint to the velocity field at the new time interval. Since the pressure conditions on both the flat free surface and the free-slip bottom are of the Neumann type, an additional solvability condition is required to make the solution unique when solving the pressure Poisson equation. In the present numerical implementation, a convenient value for pressure is assigned to a chosen grid point. For spectral representation of the pressure field, this corresponds to explicitly specifying the constant mode. Such a solvability condition, however, would not affect the solution as only the pressure gradients are required in integrating the velocity equations.

4. Posing the simulation

We pose the simulations to mimic the subsequent development of the flow of water which has been driven by a wind field accelerating constantly from rest to a reference speed of ua in a time interval of t0 as in the experiment of Melville et al. (1998). Accordingly, numerical simulation of the flow is initiated with a mean streamwise velocity distribution (in the unit of cm s−1):

U(z) = At0  (1 + 2η2)erfc(η)− 2 π1/2ηexp(−η 2)  , (4.1)

where η = − z/2(νt)1/2. Equation (4.1) corresponds to the velocity profile of a plane laminar shear layer driven by a surface stress, which induces a linear increase in the surface velocity with a constant acceleration A for a period of t0. Assuming continuity in the evolution of the air flow, the shear stress acting at the surface,

τs = ρν  dU dz  z=0 = 2 π1/2ρA(νt) 1/2, (4.2)

therefore increases with time in the rate proportional to t1/2 from the start of the simulation. The surface velocity develops freely and the maximum turbulence production occurs near the water surface. The constant acceleration A depends on

(6)

the reference wind speed ua. For the three reference wind speeds considered in the experiments of Melville et al. (1998), ua= 3, 4 and 5 m s−1, the estimated acceleration constants are A = 0.27, 0.5 and 1 cm s−2 with the corresponding acceleration periods

t0= 45, 30 and 20 s. The characteristic depths Ls= 2(νt0)1/2 and velocities us= At0 of the shear layers at the beginning of the simulation are 1.4, 1.1 and 0.9 cm, and 12.15, 15 and 20 cm s−1, respectively. These result in the corresponding Reynolds numbers of the initial driven shear flows, Re = usLs/ν= 1630, 1650 and 1800, for the three reference wind speeds ua= 3, 4 and 5 m s−1, respectively. From (4.2), the initial surface friction velocities are determined as u= (τs/ρ)1/2= 0.45, 0.56 and 0.71 cm s−1, which result in the wall Reynolds numbers Re= uLs/ν= 60.4, 61.9 and 63.6 for ua= 3, 4 and 5 m s−1, respectively.

In the experiment of Melville et al. (1998), predominant streamwise-elongated high-speed streaks (see also the present simulation results in§ 6) were observed. The streaks reveal a nearly equal spacing, which varies with the wind speeds. Accordingly, the length and width of the computational domain are chosen to be at least six times that of the observed characteristic steak spacing; the domain size in both x and

y is therefore 45, 33.6 and 23 cm, respectively, for the three reference wind speeds

ua= 3, 4 and 5 m s−1. All the simulations are carried out with 1283 grid points. This means that the lateral region between two streaks is resolved by at least twenty grids. The corresponding non-dimensional streamwise and spanwise grid sizes in wall unit,

x+ and y+, are 15.89, 14.54 and 12.87 for u

a= 3, 4 and 5 m s−1, respectively. In the vertical direction near the upper surface, the stretched grid distribution allows sufficient resolution of the viscous layer. For the three simulations considered, there are at least ten grids in the near-surface region (z+6 10).

In the depth direction, the computational domain is closed by imposing a free-slip bottom condition at a depth of a quarter of the length (and the width) of the domain to emulate the infinite depth. Thus, the depth of the turbulent layer will grow in time through entrainment in proportion to the applied wind stress. The thickness of the shear layer develops to about a quarter of the depth of the computational domain at the end of the simulation. For all the simulations conducted in the present study, the numerical depth is sufficiently large not to affect the near-surface flow.

In order to generate a realistic turbulent velocity field, a fluctuating but solenoidal vector field is first generated by taking the curl of another random vector field. The simulation is then conducted using the fluctuating solenoidal vector field as the initial velocity condition without the mean streamwise velocity. As the integration of the momentum equations proceeds, the simulated flow field adjusts itself from an initial artificial random vector field to a realistic turbulence. The ‘spin-up’ simulation is performed until the simulated flow reaches an equilibrium state in which the inertial and dissipating ranges of the energy spectra converge. During this stage, the turbulence field near the surface become effectively anisotropic, with a rapid attenuation of the vertical fluctuating velocities and an accompanying increase in the horizontal turbulence intensities. The turbulent velocity fluctuations, which are homogeneous in the horizontal plane ((x, y)-plane), is then added to the mean flow (4.1) to form the initial velocity field for the simulation.

The abrupt superposition of the turbulent velocity field onto the mean streamwise velocity would require additional spin-up time before meaningful flow statistics and structures can be obtained from the simulations. Evolutions of three averaged flow properties, the surface-averaged mean velocity us=uz= 0 and turbulent kinetic en-ergyu2+v2z= 0, and the volume-averaged turbulent kinetic energyu2+v2+w2V,

(7)

tUs/Ls 0 100 2000 0 100 200 t (s) f ug z=0 /us f KE gz=0 /us 2 f KE g /us 2 0 200 300 0 10 0.4 0.6 0.8 1.0 f ug z=0 /us 0.4 0.6 0.8 1.0 f ug z=0 /us 0.4 0.6 0.8 1.0 0 0.05 0.10 f KE gz=0 /us 2 f KE g /us 2 0 0.05 0.10 f KE gz=0 /us 2 f KE g /us 2 0.05 0.10 (a) (b) (c) 5 0 5 10 0 5 10 15 20 100

Figure 1. Evolutions of the surface-averaged streamwise velocity uz= 0/us (solid curves),

the surface-averaged turbulent kinetic energy u2+ v2z= 0/u2s (dashed curves), and the

volume-averaged turbulent kinetic energy u2+ v2V/u2

s (dash-dotted curves) for the wind

speeds ua= (a) 5 m s−1, (b) 4 m s−1and (c) 3 m s−1.

are shown in figure 1, where the turbulent fluctuating quantity q= q− q, and · and ·V represent the averaging over the horizontal plane and the entire volume, respectively. Immediately after the start of the simulation, both the surface mean velocity and fluctuation intensity decrease drastically showing a rapid adjustment of the imposed initial flow to turbulence. The subsequent evolution of the flow after the surface mean velocity and turbulent intensity reaching the minima can be divided into two stages. During the first stage (non-dimensional time tus/Ls≈ 30 to 120), the mean surface velocity rises in a higher rate while the turbulent intensity evolves gradually. After the developing stage (tus/Ls & 120), the flow can be considered as that of fully

(8)

developed turbulence in which the mean surface velocity increases in a reduced rate while the surface turbulent intensity approaches a constant level. To determine if the simulated flow has reached an equilibrium state and has become representative of turbulence, the evolutions of the one-dimensional turbulent energy spectra at various depths are examined. (The details are the same as those in Tsai 1998.) Once the surface turbulent intensities reach the fully-developed stage, the distributions of the inertial and dissipation ranges of the energy spectra also reach the steady states at both the shallow and submerged depths, indicative of the equilibrium of turbulence. The turbulence statistics, which are discussed in§ 5, are space and time averages (over the whole horizontal plane and the large-scale time tus/Ls= 120∼ 200) of the flow fields in this fully developed stage.

5. Flow statistics

5.1. Mean profiles

For the temporal mean current beneath a wind-blown air–water interface, a two-layer velocity profile, analogous to that of the turbulent wall layer, has been observed in the laboratory (Wu 1975, 1984; McLeish & Putland 1975) and field (Churchill & Csanady 1983) experiments. The mean current near, but not immediately below, the interface is found to be characterized by a logarithmic velocity profile:

us− u(z) u = 1 κ ln z ++ ψ+= 1 κln z z0 , (5.1)

where κ is the von K ´arm ´an constant, the friction velocity u= (τs/ρ)1/2, the wall coordinate z+= zu

/ν, and the constant of integration ψ+ which is related to surface roughness length z0 by z0= (ν/u∗) exp(−κψ+). The greater the constant ψ+, the smaller the surface roughness length z0. The universal function of (5.1) is also appli-cable to the profiles of the mean temperatureθ(z) and dissolved gas concentration c(z) as θ(z) − θs θ = 1 κθ ln z++ ψ+ θ (5.2) and c(z) − cs c = 1 κc ln z++ ψ+ c, (5.3)

where θs=θ(0), cs=c(0), θ∗= − Q/(ρcpu), c∗= − J/u, and Q and J are the mean heat and gas fluxes, respectively.

Immediately below the interface, the vertical component of the turbulence is largely suppressed so that a thin viscous sublayer exists where molecular viscosity dominates the momentum transport. This viscous sublayer was found to be thinner than the analogous sublayer adjacent to a solid no-slip wall (Csanady 1984; Wu 1984). Adopting the surface renewal model that the viscous sublayer undergoes cyclic growth and disruption, Liu & Businger (1975) proposed an exponential profile for the mean velocity distribution within the sublayer as

us− u(z) us− ub = 1− exp  −z δu  , (5.4)

where ub denotes the bulk mean velocity. δu is a scaling depth equivalent to the thickness of the viscous sublayer, and relates to the average residence duration of

(9)

renewal fluid elements tr∗ through δu= (νtr∗) 1/2= νus− ub u2 ∗ . (5.5)

In developing (5.4), Liu & Businger (1975) assumed that the probability of renewal of any fluid element is independent of the duration for which it has resided in contact with the surface (Danckwerts 1951), and the probability distribution of the residence time tr is described by an exponential function, exp(−tr/tr)/tr∗. Equation (5.4) can be re-written in terms of the wall coordinate as

us− u(z) u = ξu  1− exp  −z+ ξu  , (5.6)

where ξu= δuu/νis the sublayer parameter. The asymptotic limit of (5.6) as (z+/ξu)→ 0 becomes the familiar linear profile (us− u)/u= z+.

For heat and gas transport across the conduction and diffusion sublayers, similar exponential profiles can be obtained for temperature and concentration distributions within the sublayers as

θ(z) − θs θ = ξθ  1− exp  −P r z+ ξθ  , (5.7) c(z) − cs c = ξc  1− exp  −Scz+ ξc  , (5.8)

where the Prandtl number P r = ν/νθ and the Schmidt number Sc = ν/νc, νθ and

νc are the diffusivities of heat and dissolved gas, respectively, the non-dimensional parameters ξθ= δθu/νθand ξc= δcu/νc, the thicknesses of the thermal and molecular sublayers δθ= (νθtθ)1/2 and δc= (νctc)1/2.

Comparisons between the theoretical and simulated vertical mean distributions of the streamwise velocity, the temperature and the gas concentration are shown in figure 2 for the three wind speeds ua= 5, 4 and 3 m s−1. As shown in figure 2, the simulated mean velocity, temperature and concentration distributions are all well represented by the two-layer profiles. Our simulated flows reveal an exponential viscous sublayer above z+ ≈ 10 and a logarithmic inertial layer between z+ ≈ 20 and 100. For the velocity comparison, the available laboratory measurements by Wu (1984) at the wind speeds of 4.8, 4.0 and 3.0 m s−1 are also plotted. The experimental measurements, the numerical simulations and the parameterized profiles are comparative for the whole sublayer range extending to the buffer layer.

In determining the non-dimensional constants κ and ψ+ in (5.1), least-squares fitting to the simulated mean streamwise velocities for the depths z+= 20 to 100 is employed. The exponential and logarithmic profiles are then matched by equating the right-hand sides and their first derivatives of (5.1) and (5.4) at the matching depth z+= ζ+ (Liu, Katsaros & Businger 1979; Kraus & Businger 1994). This yields two nonlinear equations for the sublayer parameter ξu and the matching coordinate

ζ+. The same procedure is employed to reach the matched temperature and gas concentration profiles.

The corresponding parameters of the matched mean velocity profiles in figure 2 for the three wind speeds are listed in table 1. For logarithmic fitting of the averaged velocity distributions, the typical computed values of κ are between 0.35 and 0.40, which are close to the universal K ´arm ´an constant, 0.4. The typical values of the computed ψ+, however, are much less than 5.0 of a no-slip smooth wall, indicating

(10)

0 10 15 100 101 102 (a) (d) (g) (c) ( f ) (i) (b) (e) (h) z+ 5 0 10 15 100 101 102 5 10 15 100 101 102 5 0 10 15 100 101 102 z+ 5 0 10 15 100 101 102 5 10 15 100 101 102 5 0 10 15 100 101 102 z+ 5 0 10 15 100 101 102 5 10 15 100 101 102 5

Figure 2.(a–c) Mean profiles of the streamwise velocity (us−u(z))/u, (d–f ) the temperature

(θ(z)−θs)/θ, (g–i) the gas concentration (c(z)−cs)/cfor various wind speeds: ua= 5 m s−1

(a, d, g), 4 m s−1 (b, e, h), and 3 m s−1 (c, f, i). The computed results are represented by solid

curves. The dash-dotted curves are the exponential profiles (5.6) of the viscous sublayer developed based on the surface renewal model, and the dash-dot-dotted curves are the linear

approximation. The dashed curves denote the matched logarithmic profiles (5.1).䊊, laboratory

measurements by Wu (1984) at wind speeds of 4.8, 4.0 and 3.0 m s−1.

ua us uz0 δu tr∗ (m s−1) (cm s−1) (cm s−1) κ ψ+ (mm) ζ+ ξ u (cm) (s) 5.0 13.0 0.74 0.39 1.1 0.089 20.9 9.8 0.13 1.8 4.0 10.7 0.61 0.34 1.5 0.092 23.3 10.5 0.17 2.9 3.0 8.2 0.47 0.36 1.9 0.098 26.2 11.2 0.24 5.6

Table 1. Computed parameters in the mean velocity profiles (5.1) and (5.4) for the three

wind speeds.

the increase in equivalent surface roughness z0. For the flow beneath a wind-blown sea surface, this is caused by direct energy supply to the near-surface turbulence through influences such as breaking wavelets (e.g. Csanady 1984). In the present case of stress-driven shear flow, the increase infers the enhancement of near-surface

(11)

horizontal turbulence intensity by the presence of the upper surface. (This will be discussed in§ 5.2.) The increase in surface roughness z0has also been noted in a recent direct numerical simulation of an oceanic boundary layer driven by wave breaking (Sullivan, McWilliams & Melville 2004). For the limiting case with no breaking surface forcing, their simulation also indicates a higher effective roughness z0 than that of wall-bounded turbulence. Their simulations with breaking surface forcing further reveal that wave breaking would effectively increase the surface roughness.

The computed non-dimensional viscous sublayer parameters ξu are around 10, smaller than that immediately next to a smooth no-slip wall (≈16), and increase slightly (within the range of O(10−1)) with decreasing wind speeds. This is quanti-tatively similar to that observed in the experiments (Csanady 1984; Wu 1984). The increase in ξu with decreasing wind speed also implies the thickening of the viscous sublayer δu and the growth of resistance against the vertical momentum transport.

5.2. Turbulence intensities and production

The vertical distributions of the turbulent velocity variance normalized by the surface friction velocity, ui

2(z)/u2

∗, are shown in figure 3 for the three reference wind speeds

ua= 5, 4 and 3 m s−1. Away from the viscous layer, the magnitudes of the turbulent velocity variances are close to those next to a no-slip solid wall (e.g. Kim, Moin & Moser 1987). In particular, the upper bound of w2/u2

∗ is about 1; it does not exhibit a large increase like that in the large-eddy simulations of oceanic mixed layers (Skyllingstad & Denbo 1995; McWilliams et al. 1997), in which the Craik– Leibovich ‘vortex force’ attributed to the wave–current interaction mechanism (Craik & Leibovich 1976) is incorporated into the simulations.

The turbulence field near the interface is effectively anisotropic and two-dimensional, with a rapid attenuation of the vertical velocity fluctuation and an accompanying increase in the streamwise component on approaching the interface. In comparison with the streamwise and vertical components, the spanwise turbulence intensity shows relatively mild variation. This is contrary to the observations of the turbulent shear flows bounded by shear-free surfaces (such as the free-surface jet flow in Walker, Chen, & Willmarth 1995; the free-surface wake flow in Tsai 1998; and the open-channel flow in Calmet & Magnaudet 2003), in which the transfer of the turbulence energy to the streamwise component is much more significant than that to the spanwise one. Immediately beneath the upper surface and within the viscous sublayer (z+. 10), both streamwise and spanwise turbulence intensities exhibit salient features, in particular for the case of ua= 3 m s−1(figure 3c). (The depths of the viscous sublayers, z+= ξ+≈ 10, are marked in the figures.)

For turbulence bounded by a shear-free surface, the normal turbulent fluctuations are forced to decay towards the boundary. This is a manifestation of the blocking effect brought about by the kinematic condition of vanished normal velocity at the surface. The shear-free boundary conditions, nevertheless, give rise to much weaker constraints on the tangential velocities and lead to the dissipation rate close to the surface actually becoming smaller than in the bulk of the flow (Teixeira & Belcher 2000). Accordingly, the horizontal velocity fluctuations increase upon approaching a shear-free surface (e.g. Perot & Moin 1995; Walker et al. 1996; Tsai 1998; Calmet & Magnaudet 2003). On the present stress-imposed surface, the constraints of the horizontal fluctuating velocities, ∂u/∂z= ∂v/∂z= 0, are identical to those of a shear-free surface although the mean velocity field is subjected to a constant vertical gradient, ∂u/∂z = τs/µ. As figure 3 shows, the intensities of the horizontal turbulence are also amplified as they approach a stress-imposed surface. This is also because of the decrease in the

(12)

0 5 10 0 50 100 (a) (b) (c) z+ 0 5 10 0 50 100 z+ 0 5 10 0 50 100 z+ — — — — — — — — — — — —

Figure 3.Vertical distributions of the normalized turbulent velocity variances, u2(z)/u2

(solid curves),v2(z)/u2

∗ (dash-dotted curves), w2(z)/u2∗ (dashed curves) and the averaged

turbulent kinetic energy,u2+ v2+ w2(z)/u2

∗ (thick solid curves) for various reference wind

speeds (a) ua= 5 m s−1, (b) 4 m s−1 and (c) 3 m s−1. The depths of the viscous sublayers,

z+= ξ+≈ 10, are marked by dashes in the figures.

near-surface dissipation profiles, which will be discussed later. Such an enhancement of the horizontal turbulence intensities has also been observed in the recent simulation by Sullivan et al. (2004).

The drastic changes within the viscous sublayer are also discernible in vertical variations of the vorticity variances,2

x(z), ω2y(z) and ωz2(z), as well as the vertical strain rate (horizontal divergence), (∂w/∂z)2 = (∂u/∂x + ∂v/∂y)2, as shown in figure 4. Both the vertical and horizontal vorticity variances, in particular the spanwise component, reach the maxima within the viscous sublayer. At the interface, the streamwise and spanwise vorticities decrease to negligible levels, whereas the vertical component remains finite. Away from the interface (z+& 40; the characteristic depth of the shear layer is approximately at z+≈ 70), the horizontal vortical field becomes essentially isotropic and the vorticity variances continue reducing. No significant

(13)

0 0.5 1.0 0 50 100 (a) (d ) (b) (e) (c) ( f ) z+ — — — — — — — —— — — — 0.1 0.2 0 50 100 0 0.5 1.0 0 50 100 z+ 0.1 0.2 0 50 100 0 0.5 1.0 0 50 100 z+ 0.1 0.2 0 50 100

Figure 4. Vertical distributions of the vorticity variances (a–c): streamwise vorticity,ωx2(z)

(solid curve), fluctuating spanwise vorticity,2y(z) (dash-dotted curve), and vertical vorticity,

2

z(z) (dashed curve); and the vertical velocity–strain variance, (∂w/∂z)2(z) (d–f ) for the

wind speeds of ua= 5 m s−1(a, d), 4 m s−1(b, e) and 3 m s−1(c, f ). All the terms are normalized

by us/Ls. The depths of the viscous sublayers, z+= ξ+≈ 10, are marked by dashes.

change in the streamwise vorticity variance is observed in these submerged depths for all the wind speeds considered throughout the simulations, implying the absence of predominant streamwise vortices.

To further elucidate the impact of the sheared surface on the transport of near-surface turbulence energy, the energetics budgets are examined next. Assuming equilibrium, the mean budget equations for turbulent flow with two-dimensional mean shear ∂u(z)/∂z are

∂u2 ∂t =−2u w∂u ∂z    P11 −∂u2w ∂z    T11 +2 ρ p∂u  ∂x    Π11 +ν∂ 2u2 ∂z2    D11 −2ν ∂u ∂xk ∂u ∂xk    11 , (5.9)

(14)

–0.001 0 0.001 0 20 40 (a) (b) (c) z+ –0.0005 0 0.0005 –0.0002 0 0.0002 0 20 40 z+ 0 20 40 z+

Figure 5.Vertical distributions of the various production/consumption terms of (a)

stream-wise, (b) spanwise and (c) vertical turbulent intensities in the Reynolds-stress equations (5.9),

(5.10) and (5.11) for the reference wind speed ua= 5 m s−1: P11, thick solid curves; Tii,

dash-dot-dotted curves; Πii, solid curves; Dii, dash-dash-dot-dotted curves; ii, dashed curves. All the terms are

normalized by u3 s/Ls. ∂v2 ∂t =− ∂v2w ∂z     T22 +2 ρ p∂v  ∂y    Π22 +ν∂ 2v2 ∂y2    D22 −2ν ∂v ∂xk ∂v ∂xk    22 , (5.10) ∂w2 ∂t =− ∂w3 ∂z    T33 −2 ρ ∂pw ∂z + 2 ρ p∂w  ∂z    Π33 +ν∂ 2w2 ∂z2    D33 −2ν ∂w ∂xk ∂w ∂xk    33 , (5.11)

where P11 is the shear production rate, Tii the turbulent advection rate, Πii the velocity–pressure correlation term, Diithe viscous diffusion rate, and iithe dissipation rate. Figure 5 shows the vertical variations of the various terms in (5.9), (5.10) and (5.11) for the reference wind speed ua= 5 m s−1. The distributions for other wind speeds resemble those in figure 5.

(15)

The distributions of the energy budgets shown in figure 5 reveal some general features known for shear turbulence bounded by a no-slip wall (e.g. Mansour et al. 1988): the shear production rate, the P11term, is the major turbulence production for the streamwise budget, which reaches its maximum level within the viscous sublayer before diminishing at the surface. The pressure–strain term, Π11, converts the u2 component to v2 and w2 and hence is a consumption for u2 and a major source for the other two components, except near the viscous sublayer for w2. The dissipation terms, ii, account for the major loss for all three energy components. The contributions of turbulent advection, Tii terms, and viscous diffusion, Dii terms, are in general negligibly small except within the viscous sublayer. These general features of the turbulent kinetic energy budget in the region of submerged water have also been observed in the simulation of a stress-driven oceanic boundary layer by Sullivan et al. (2004).

The main differences in the energetics budgets between the present stress-driven turbulence and the wall-bounded turbulence are within the viscous sublayer attributed to the different boundary conditions. For a free surface (sheared or shear-free), the vertical fluctuating strain rate, ∂w/∂z= − (∂u/∂x+ ∂v/∂y), increases rapidly approaching the boundary (figure 4). Concurring with the amplification of near-surface streamwise turbulent velocities (figure 3), the convective transport of streamwise turbulence, T11, increases from consumption in the submerged water to production within the sublayer, and reaches its maximum at the upper surface. The near-surface streamwise turbulent transport is of a magnitude comparable to that of shear production P11 in the submerged region. Such a significant increase in near-surface T11 is balanced by pressure transport (Π11) and an enhancement in viscous diffusion (D11). The present near-surface turbulent energy budget is somewhat different from that in a similar simulation by Sullivan et al. (2004), but with a higher Reynolds number, in which negligible contributions from turbulent transport and viscous diffusion have been observed.

The variations in the dissipation rates within the sublayer next to a stress-imposed surface are completely opposite to those next to a wall boundary. The dissipation rates of both the streamwise and spanwise turbulent energies, 11 and 22, decrease drastically when approaching the upper surface. This results in a corresponding increase in the horizontal fluctuating velocities as explained by the theoretical model of Teixeira & Belcher (2000). The dissipation rate of the vertical turbulent energy,

33, increases rapidly and reaches its maximum level at the surface, indicative of an abrupt reduction in vertical velocity towards the surface.

For shear turbulence near a no-slip wall, viscous diffusion dominates the productive transport of u2 and v2, and is negligibly small for vertical component w2. However, for turbulence bounded by a sheared surface, both the streamwise and spanwise diffusion rates, D11 and D22, become the primary consumption of the horizontal components within the viscous sublayer and reach their negative minima at the surface. On the other hand, viscous diffusion, D33, is the major term for vertical transport immediately beneath the surface. This confirms that the vertical exchange across the interface is controlled by viscous diffusion within the viscous sublayer.

6. Surface characteristics

6.1. Streaky appearance

The first clue of the existence of organized structures within the stress-driven surface turbulent boundary layer is the emergence of elongated filaments of high-speed streaks

(16)

with somewhat equal transverse spacing on the water surface as have been observed in both field and laboratory experiments. In the experiment of Melville et al. (1998) the typical spanwise spacing between the streaks ranges from∼3 cm to ∼7 cm as the reference wind speed decreases from 5 m s−1 to 3 m s−1. These streamwise elongated streaks with transverse spacing of O(1 cm) to O(10 cm) have also been observed in the ocean (Gemmrich & Hasse 1992) and lake (Kenney 1993; Woodcock 1941). Such a surface signature arises from the wind-driven shear layer and is geometrically similar to the low-speed streaks near a stationary wall, which develop within the viscous sublayer and extend into the logarithmic region (e.g. Kline et al. 1967; Smith & Metzler 1983).

In the present simulation of surface-driven turbulent shear flow, high-speed surface streaming, similar to that observed in the laboratory and field experiments, emerges on the surface boundary. Instantaneous distributions of the surface streamwise velocity and the temperature contours for the three wind speeds are shown in figure 6. The distribution patterns clearly demonstrate the existence of elongated surface streaks carrying cooler and faster-moving fluids. The major streaks are arranged with somewhat equal transverse spacing, and the spacing increases with decreasing wind speed.

To further visualize the structures and also the inception and evolution of the surface streaming, 1282 uniformly distributed floating Lagrangian particles are released on the water surface and their trajectories are tracked. Figure 7 depicts the evolution of the Lagrangian particles for the wind speed of ua= 5 m s−1. Immediately after the particles are released, they aggregate and form a ‘fish-skin’ pattern as shown in figure 7(a). The merged particles further accumulate into several elongated groups and align along the fast-moving velocity streaks (figures 7b and 7c). The rapidity of the formation of this surface streaky pattern is the same as the observation of Melville et al. (1998) and as that described by Woodcock (1941): ‘When the powder is dusted over the water it moves into lines so quickly that, even in moderate breeze, the eye can hardly detect the transition from random scattering to the linear pattern.’ It takes less than two seconds for the particles to evolve from the initial uniform distribution to the streaky pattern in figure 7(c). The aligned particles continue travelling along the high-speed streaky regions although some particles of the individual streets may diverge and join other ones.

As time proceeds, some segments of the particle streets break, and the particles further cluster in a streamwise fashion and form shorter streets, as shown in figures 7(e) and 7(f ). These shorter streets, however, still travel along the same velocity streaks. To confirm that the elongated velocity streaks persist though the streets of Lagrangian particles seem to have dislocated, a new set of uniformly distributed particles are released on the surface immediately after the time of figure 7(d) (5.364 s from the release of the first set of particles). The distributions of the newly released particles at the corresponding escaping times in figure 7 are shown in figure 8. Similar evolution patterns are observed for the two sets of surface floating particles released at different times.

Accompanying the more organized fast-moving cooler streaks are slowly-moving warmer localized spots appearing intermittently at the surface, as shown in figure 6. The appearance of these random surface spots can be educed from the initial formation of the ‘fish-skin’ pattern in figure 7(a). At the early stage after the particles are released, the floating particles bypass the localized regions of lower streamwise velocities and begin to accumulate around the upstream edges of the slowly moving spots. The process is focused on by showing in figure 9(a) the distribution of the

(17)

10 (a) 5 0 –5 –10 –10 –5 0 5 10 10 5 0 –5 –10 –10 –5 0 5 10 y (cm) 10 (b) 5 0 –5 –10 –10 –5 0 5 10 5 0 –5 –10 –10 –5 0 5 y (cm) 10 (c) (d) (e) ( f ) 5 0 –5 –10 –10 –5 0 5 10 5 0 –5 –10 –10 –5 0 5 y (cm) x (cm) x (cm)

Figure 6. Representative contour distributions of the streamwise velocity (a–c) and

tempera-ture (d–f ) at the free surface for the three reference wind speeds: (a, b) 5 m s−1; (c, d) 4 m s−1;

(e, f ) 3 m s−1. The flow travels from left to right. The warm/cold-colour areas represent regions

(18)

(a) (b) 10 5 0 –5 –5 0 5 10 –10 10 5 0 –5 –10 –10 (c) 10 5 0 –5 –5 0 5 10 –10 –10 (d) 10 5 0 –5 –5 0 5 10 –10 –10 ( f ) 10 5 0 –5 –5 0 5 10 –10 –10 (e) 10 5 0 –5 –5 0 5 10 –10 –10 –5 0 5 10 –10

Figure 7.Instantaneous distribution of the surface floating Lagrangian particles at times

(a) 0.224 s, (b) 0.671 s, (c) 0.894 s, (d) 1.788 s (e) 3.129 s and (f ) 5.364 s after the release of the

particles. 1282uniformly distributed particles are released.

Lagrangian particles superposed on the contours of the streamwise velocity shortly after the release of the particles. A high-speed streak bifurcates into branches when it encounters a localized low-speed spot. This can be seen from the streamwise velocity contours (figure 6) as well as the distribution of Lagrangian particles (figure 7c, d). A close-up superposed distribution of streamwise velocity and particles revealing the bifurcation events is shown in figure 9(b). Such a feature is identical to the streak dislocation or bifurcation observed by Melville et al. (1998). The merging and bifurcation of Lagrangian particles also look similar to those shown in McWilliams et al. (1997, figure 20), which they claimed are due to Langmuir cells.

(19)

(a) (b) 10 5 0 –5 –5 0 5 10 –10 10 5 0 –5 –10 –10 (c) 10 5 0 –5 –5 0 5 10 –10 –10 (d) 10 5 0 –5 –5 0 5 10 –10 –10 ( f ) 10 5 0 –5 –5 0 5 10 –10 –10 (e) 10 5 0 –5 –5 0 5 10 –10 –10 –5 0 5 10 –10

Figure 8. Instantaneous distribution of the surface floating Lagrangian particles at times

(a) 0.224 s, (b) 0.671 s, (c) 0.894 s, (d) 1.788 s (e) 3.129 s and (f ) 5.364 s after the release of the second set of uniformly distributed particles. The particles are released immediately after the time of figure 7(f ).

6.2. Spanwise distribution of streaks

For wall-bounded turbulent flows, low-speed streaks similar in geometry to the present high-speed counterparts on a stress-imposed surface, have been identified within the viscous sublayer in close proximity to the no-slip wall. The non-dimensional mean spanwise spacing, λ+= λu

, between these low-speed streaks next to a no-slip wall is generally accepted to be around the value of 100 (e.g. Kline et al. 1967), where

λ is the dimensional mean spanwise spacing. Further study by Smith & Metzler

(1983) indicated that such a non-dimensional mean value is essentially invariant with Reynolds number, and the streak spacing exhibits a log-normal distribution as suggested previously by Nakagawa & Nezu (1981).

(20)

(a) (b) 0 –5 –5 0 –10 –5 0

Figure 9.Close-up snapshots of the Lagrangian particles superposed on the streamwise

velocity contours. Immediately after the particles are released (a), the floating particles bypass the localized regions of lower streamwise velocities (blue areas). The particles eventually accumulate and travel along the fast-moving streaks (red areas) as shown in (b), where the fast-moving streak bifurcates into branches when it encounters a localized slowly moving spot.

-y+ 10 20 30 40 50 60 70 80 90 100 0 500 1000 1500 Frame number

Figure 10.Sequence of the streak marks for the reference wind speed ua= 5 m s−1. Each

vertical line corresponds to one image frame for streak identification. The time interval between two subsequent frames is about 0.0224 s. The dashes indicate the transverse locations of the identified streaks.

Following the study of Smith & Metzler (1983), we determine the streak-spacing values from visual identification using the image sequence of the contours of streamwise velocity and the distributions of Lagrangian particles, as presented in figures 6 and 7. The images of particle distributions are superposed on the velocity contours to allow a cross-examination when the streak positions are determined. For each frame of the images, a reference spanwise axis at a fixed streamwise coordinate is chosen where the streak positions are to be marked. A streak intersecting the reference spanwise axis is identified when the Lagrangian particles aggregate along an elongated well-defined region of high streamwise velocity with a streamwise extent of at least x+= 300. The locations of the streak intersections along the transverse axis are then marked using a graphic digitizer software, and the transverse coordinates are stored. An example of the sequence of streak marks for the reference wind speed

(21)

P(λ+) P(λ+) P(λ+) 200 400 600 0 0.002 0.004 0.006 0.008 (a) (b) (c) λ+ = 183 σ+ λ = 90 200 400 600 0 0.002 0.004 0.006 0.008 λ+ 0 200 400 600 0.002 0.004 0.006 λ+ = 161 σ+ λ = 63 λ+ = 201 σ+ λ = 77

Figure 11. Probability density histograms of the spanwise streak spacing and the log–normal

distributions with the corresponding mean value λ+ and standard deviation σ+

λ for the wind

speed (a) ua= 5 m s−1, (b) 4 m s−1and (c) 3 m s−1.

frame for streak identification, and the dashes indicate the transverse locations of the identified streaks. The transverse distances between the neighbouring streaks, λ, are then measured from the sequence of streak marks.

Distribution histograms of the spanwise spacing between two streamwise streaks at the water surface determined from the described procedure of streak identification and counting are shown in figure 11 for the three reference wind speeds considered. The streak-spacing values are presented in non-dimensional form λ+= λu

. The three distributions all skew toward values lower than the means, λ+, similar to that of low-speed streaks within a wall boundary layer. The determined non-dimensional values of the mean spacing λ+ are 183, 161 and 201, with the corresponding standard deviations σ+

λ = 90, 63 and 77, respectively, for the three wind speeds ua= 5, 4 and 3 m s−1, as plotted in figure 12(a). These results are about twice the values of λ+≈ 100 and σ+

λ ≈ 40 for the mean spacing of low-speed streaks near a no-slip wall (Smith & Metzler 1983).

(22)

ua (cm s–1) λ (cm) 2 4 6 8 10 3 4 5 3 4 5 0 100 200 300 (a) (b) λ+

Figure 12.(a) Variation of the non-dimensional mean streak spacing, λ+, determined by the

visual procedure with the wind speed. (b) Comparison of the dimensional streak spacings

determined by both the visual procedure () and the spectral analysis () in the present

simulations with those from the experiments of Melville et al. (1998) () for the three wind

speeds. The corresponding spacing of the second peaks for ua= 5 m s−1 from the spectral

analysis is beyond the upper range and is not shown in the figure.

As suggested by Nakagawa & Nezu (1981) and also Smith & Metzler (1983), the skewed distribution of the streak spacing can be characterized statistically using a log–normal probability density function determined from the values of λ+ and σ+

λ in the form: P(λ+) = √ 1 2πλ+ψ 0 exp −1 2  1 ψ0 lnλ + λ+0 2 , (6.1) where λ+0 = λ+(1 + ψ2

λ)−1/2 is the median value of λ+, ψ0=

ln(1 + ψ2 λ)

1/2 is the coefficient of variation of ln λ+, and ψ

λ= σλ++. The fitted log–normal probability density functions are superposed over the corresponding histograms in figure 11. The comparison between the histograms and the log–normal distribution functions appear to be quite good. These results also suggest the possible similarity between the small-scale high-speed streaks on a wind-blown water surface and the low-speed streaks next to a wall boundary.

The characteristics of the surface streak spacing can also be estimated by examining the spanwise spectrum of the streamwise-averaged surface fluctuating streamwise velocity, us(y) =



u(x, y, 0) dx. The distributions of the normalized average spectra, us(λ) /u2, for the three reference wind speeds are shown in figure 13. The power spectra are obtained from the ensemble average of the streamwise-averaged velocity

us(y) at various times within the interval when the surface turbulence reaches the stationary state. Despite the initial broadband spectra without any predominant span-wise wavelengths, the streamspan-wise fluctuating velocities are dominated by components

(23)

200 1000 1800 0.01 0.03 0.05 (a) (b) (c) 200 1000 1800 0.01 0.03 0.05 200 1000 1800 0.01 0.03 λ+

Figure 13. Distributions of the normalized average spectra, us(λ) /u2, for the three

reference wind speeds of (a) 5 m s−1, (b) 4 m s−1 and (c) 3 m s−1. The power spectra are

obtained from ensemble averages of the streamwise-averaged velocity us(y) at various times

within the intervals when the surface turbulence reaches the stationary state.

with preferential cross-spacing wavelengths when surface streaks appear, as shown in figure 13. The spectra all display two characteristic peaks for the three wind speeds considered. The non-dimensional average streak spacing λ+, which can be considered as the shorter wavelength of the first spectrum peak, are 219, 249 and 270, respectively, for the wind speeds ua= 5, 4 and 3 ms−1. Consistent with the finding for low-speed streaks next to a no-slip wall (Kline et al. 1967), the spacing values determined in this more objective manner are greater than those determined by the visual counting scheme. The longer wavelength implies possible spanwise modulation of the streaks as shown in figure 6.

While the present numerical simulations are not completely compatible with the conditions in the experiments of Melville et al. (1998), the simulated surface signatures, however, quantitatively resemble the observations. It is therefore worth comparing the characteristic cross-spacing of the surface streaks arising from the simulated boundary layer with that from a wind-driven surface shear layer. Summary of the surface streak spacing determined by both the visual procedure and the spectral analysis in the present study, and from the experiment of Melville et al. (1998) is shown in figure 12(b). The spacing values observed in the experiment

(24)

of Melville et al. (1998) are in general slightly larger than those determined in the present numerical simulations. For the wind speed of ua= 5 m s−1, the numerical predictions are close to the measurement of Melville et al. (1998). Both numerical and experimental results clearly show the same trend: the streak spacing decreases as the wind speed increases.

7. Coherent flow structures

The surface features discussed in§ 6 are strongly linked to the coherent structures in the subsurface flow. To reveal the underlying structures, three-dimensional cut-away views of instantaneous streamwise velocity u and temperature θ distributions are shown in figures 14 and 15 for the wind speed ua= 5 m s−1. The corresponding distributions of vertical velocity w and gas concentration c on the facing cross-stream planes are also depicted in figures 14 and 15, respectively.

It is discernible that the longitudinal streaks are the surface signatures of an array of high-speed jets strongly contiguous to the surface. The vertical depths of the jets are limited within the viscous sublayer. Looking into the distribution of the vertical velocity, downdrafts occur beneath the longitudinal streaks (marked with downward arrows in the figure). The cold water collects along the high-speed jets which are characterized by cold streaks at the surface. The well-known cool-skin effect (e.g. Katsaros 1980), therefore, is attributed to these coherent thermal structures within the surface sublayer. Despite the different surface boundary conditions between the temperature and the gas concentration fields, the two scalar distributions are virtually identical except in the region immediately next to the water surface.

Beneath the sublayer, distinct ‘tongues’ of warmer water with lower streamwise velocities stick vertically out of the well-mixed submerged region. In contrast to the more organized cold high-speed jets, these warm low-speed tongues appear randomly with no spatial pattern. Comparing the distributions of the streamwise and vertical velocities on the cross-stream vertical plane, it can be seen that the coherent tongues are formed by upward advection of the warm and low-speed water underneath. When the upwelling is strong enough to penetrate the viscous sublayer (one particular event is highlighted by a circle in figure 14), it results in a localized low-speed spot at the surface. The bifurcation of the elongated surface streak occurs (figure 9) when a high-speed jet encounters such a localized low-speed updraft.

An explanation of this local upwelling flow is the emergence of Ω-shaped horseshoe vortices in the turbulent shear layer by turning and stretching of the spanwise vortices (Tsai 1998). These coherent vortices move upwards and impinge on the surface. The upper heads of the Ω-shaped vortices thereby induce the local upwelling flows and bring up the submerged slowly moving fluids. Such an impinging process can also be deduced from the enhancement of horizontal divergence of the velocity field immediately beneath the surface (figure 4).

To assess the presence of circulatory flow that forms the surface streaks, we plot in figure 16 the instantaneous iso-surfaces of streamwise vorticity ωx, the contours of

ωx on various vertical cross-sections (y, z-plane) along with the streamwise velocity contour on the surface. In sharp contrast to the more organized high-speed streaks at the surface, the streamwise vortical field display is irregularly spaced and is much shorter in length. The spanwise spacing between the vortex tubes is also narrower than that between the streaks. Furthermore, the predominant streamwise vortices are strongly confined to near the surface, highly consistent with the vertical distribution of streamwise vorticity (figure 4). While the pairing of the counter-rotating vortices

(25)

–15 –10 –5 0 5 10 15 15 10 y/L x/L 5 0 –5 –10 –15 0 20 40 z+ 60

Figure 14. Three-dimensional cut-away view of streamwise velocity u for the reference wind

speed ua= 5 m s−1. The corresponding contours of vertical velocity w on the facing cross-stream

plane are also depicted on an additional vertical section. Note that the non-dimensional friction

coordinate z+ is used for the vertical axis to demonstrate better the vertical scale of the flow

structures. The downward arrow indicates the occurrence of downdraft beneath the high-speed streaks. The circle highlights a particular region where the upwelling tongue penetrates the viscous sublayer, and the upward advection carries the submerged slowly-moving water and induces a low-speed spot on the free surface.

–15 –10 –5 0 5 10 15 15 10 y/L x/L 5 0 –5 –10 –15 0 20 40 z+ 60

Figure 15. Three-dimensional cut-away view of temperature θ for the reference wind speed

ua= 5 m s−1. The corresponding contours of gas concentration c on the facing cross-stream

plane are also depicted on an additional vertical section. Note that the non-dimensional friction

coordinate z+ is used for the vertical axis to demonstrate better the vertical scale of the flow

structures. The circle highlights the region where the upwelling penetrates the conductive sublayer and induces a warm spot on the free surface.

數據

Figure 1. Evolutions of the surface-averaged streamwise velocity uz = 0 /us (solid curves),

Figure 1.

Evolutions of the surface-averaged streamwise velocity uz = 0 /us (solid curves), p.7
Figure 2. (a–c) Mean profiles of the streamwise velocity (us −u(z))/u ∗ , (d–f ) the temperature

Figure 2.

(a–c) Mean profiles of the streamwise velocity (us −u(z))/u ∗ , (d–f ) the temperature p.10
Table 1. Computed parameters in the mean velocity profiles (5.1) and (5.4) for the three

Table 1.

Computed parameters in the mean velocity profiles (5.1) and (5.4) for the three p.10
Figure 3. Vertical distributions of the normalized turbulent velocity variances, u 2 (z)/u 2 ∗

Figure 3.

Vertical distributions of the normalized turbulent velocity variances, u 2 (z)/u 2 ∗ p.12
Figure 4. Vertical distributions of the vorticity variances (a–c): streamwise vorticity, ω x 2 (z)

Figure 4.

Vertical distributions of the vorticity variances (a–c): streamwise vorticity, ω x 2 (z) p.13
Figure 5. Vertical distributions of the various production/consumption terms of (a) stream-

Figure 5.

Vertical distributions of the various production/consumption terms of (a) stream- p.14
Figure 6. Representative contour distributions of the streamwise velocity (a–c) and tempera-

Figure 6.

Representative contour distributions of the streamwise velocity (a–c) and tempera- p.17
Figure 7. Instantaneous distribution of the surface floating Lagrangian particles at times

Figure 7.

Instantaneous distribution of the surface floating Lagrangian particles at times p.18
Figure 8. Instantaneous distribution of the surface floating Lagrangian particles at times

Figure 8.

Instantaneous distribution of the surface floating Lagrangian particles at times p.19
Figure 9. Close-up snapshots of the Lagrangian particles superposed on the streamwise

Figure 9.

Close-up snapshots of the Lagrangian particles superposed on the streamwise p.20
Figure 10. Sequence of the streak marks for the reference wind speed ua = 5 m s −1 . Each

Figure 10.

Sequence of the streak marks for the reference wind speed ua = 5 m s −1 . Each p.20
Figure 11. Probability density histograms of the spanwise streak spacing and the log–normal

Figure 11.

Probability density histograms of the spanwise streak spacing and the log–normal p.21
Figure 12. (a) Variation of the non-dimensional mean streak spacing, λ + , determined by the

Figure 12.

(a) Variation of the non-dimensional mean streak spacing, λ + , determined by the p.22
Figure 13. Distributions of the normalized average spectra,   u  s (λ)/u 2 , for the three

Figure 13.

Distributions of the normalized average spectra,  u  s (λ) /u 2 , for the three p.23
Figure 15. Three-dimensional cut-away view of temperature θ for the reference wind speed

Figure 15.

Three-dimensional cut-away view of temperature θ for the reference wind speed p.25
Figure 14. Three-dimensional cut-away view of streamwise velocity u for the reference wind

Figure 14.

Three-dimensional cut-away view of streamwise velocity u for the reference wind p.25
Figure 16. Representative iso-surfaces of streamwise vorticity ωx beneath a selective area

Figure 16.

Representative iso-surfaces of streamwise vorticity ωx beneath a selective area p.26

參考文獻