LARGE EDDY SIMLIATION WAVE BREAKING OVER MUDDY SEABED*

HU Yi, NIU Xiao-jing, YU Xi-ping

State Key Laboratory of Hydro-Science and Engineering, Department of Hydraulic Engineering, Tsinghua University, Beijing 100084, China, E-mail: y-hu07@mails.tsinghua.edu.cn

LARGE EDDY SIMLIATION WAVE BREAKING OVER MUDDY SEABED*

HU Yi, NIU Xiao-jing, YU Xi-ping

State Key Laboratory of Hydro-Science and Engineering, Department of Hydraulic Engineering, Tsinghua University, Beijing 100084, China, E-mail: y-hu07@mails.tsinghua.edu.cn

The Large Eddy Simulation (LES) of the wave breaking over a muddy seabed is carried out with a Coupled Level Set and Volume Of Fluid (CLSVOF) method to capture the interfaces. The effects of the mud on the wave breaking are studied. The existence of a mud layer beneath an otherwise rigid bottom is found to have a similar effect as an increase of the water depth. As compared with the case of a simple rigid bottom, the inception of the wave breaking is evidently delayed and the breaking intensity is much reduced. The dissipation of the wave energy is shown to have very different rates before, during and after the breaking. Before and after the breaking, the mud plays an important role. During the breaking, however, the turbulence as well as the entrainment of the air also dissipate a large amount of energy.

wave breaking, muddy seabed, Large Eddy Simulation (LES), a Coupled Level Set and Volume Of Fluid (CLSVOF) method

Introduction

The wave dissipation over a muddy seabed is usually significant. Sheremet and Stone[1]compared the field measurements of the wave height over a muddy seabed with those over a sandy one, and found that the waves dissipate much more energy over the muddy seabed. Elgar and Raubenheimer[2]found that more than 70% of the wave energy flux was dissipated across the Louisiana muddy shelf according to their measurements. In addition, the fluid mud can be easily transported under the nearshore hydrodynamic environment. The related issues are often directly relevant to the maintenance of navigation channels or the assessment of a coastal ecosystems. However, both the wave dissipation over a muddy seabed and the wave induced mud transport are closely related with the interaction between the wave and the muddy seabed.

Analytical studies were extensive in the past decades to investigate the interaction between waves and muddy seabed. Dalrymple and Liu[3]developed a linear theory of viscous water propagating over a viscous mud layer. Based on the boundary layer assumption, Ng[4]obtained an explicit solution of the small amplitude wave propagation over a thin mud layer. The Non-Newtonian properties of the mud were also considered by many authors, such as in the pseudoplastic model of Zhang et al.[5], and in the viscoelastic model of Liu and Chan[6]. In addition to the analytical studies, there were also a great number of numerical simulations for wave propagations over a muddy bed. In the numerical models, the effects of the muddy bed were usually represented by an extra energy dissipation term in the existing wave models, such as in the energy balance equation and in the mild slope equation, etc.[7,8]. The analytical studies mentioned above were frequently considered in formulating the energy dissipation term.

Most of the existing studies on the interaction between waves and muddy seabed focused on the non-breaking cases. Investigations on the influences of the mud on the wave breaking are few and far between in spite of its practical importance. Among the very few related studies, Azam and Mokhatar[9]experimentally studied the surf zone dynamics over the muddy seabed. The surface roller, the velocity, and the turbulent energy were measured. Soltanpour et al.[10]proposed a numerical model for the wave propagation over a muddy seabed, in which the energy dissipation in the surf zone was considered as the sum-mation of those due to the mud and due to the wave breaking, but the interaction between the mud and the breaking wave was not taken into account.

With an ever-increasing capacity of computers, the Large Eddy Simulation (LES) for the wave breaking over a rigid bottom becomes possible and some results have been reported in the literature[11,12]. In the present article, we first perform the LES for the wave breaking over a muddy seabed based on the filtered Navier-Stokes equations. As a first step, the muddy seabed is modeled as a viscous fluid with much greater density and viscosity than water. A Coupled Level Set and Volume Of Fluid (CLSVOF) method[13,14]is used to capture the water surface and the water-mud interface.

1. Numerical approach

1.1 Governing equations

A LES is to simulate the large scale eddies of a fluid flow and represent the small eddies by a sub-grid scale stress model[15]. The governing equations are the filtered Navier-Stokes equations and the continuity equation:

where u is the velocity vector, Π is the total pressure,ρ is the density, μt=μm+μewith μmbeing the molecular viscosity of the water, μethe eddy viscosity, and μtthe total viscosity, D is the deformation tensor with Dij=(?ui/?xj+?uj/?xi)/ 2,g is the gravity acceleration vector, γ is the surface tension coefficient, κ is the curvature of the interface, δ is the Dirac delta function, and φ is the level-set function defined according to the normal distance to the interface. We define φ in such a way: negative in the upper fluid and positive in the lower fluid in this study.

The eddy viscosity is calculated following the renormalization group (RNG) sub-grid scale model proposed by Yakhot and Orszag[16], i.e.,

where H is the Heaviside step function, or,H(x )= 0if x≤0 and H(x)=x if x＞0,is the Smagorinsky constant,is the filter scale, C is a constant (C=75in the present study). By this model,μtis reduced to the molecular viscosity of the fluidmμ in a relatively low strain field andsμ in a high strain field.

Equations (1) and (2) are solved with the projection method. Discretization is performed on a staggered grid with the finite difference schemes. The second-order Runge-Kutta scheme is used for time integration. The viscous terms are discretized by the second-order central difference scheme, while the convection terms by the second-order Essentially Non-Oscillatory (ENO) upwind scheme.

To avoid the unphysical tearing of the free surface due to the large density difference of air and water, a density-weighted velocity smoother with the projection method in the air regime is applied[17]. The smoother is set to work every 20 steps following Fu et al.[17]

1.2 CLSVOF method

The level-set function φ used to represent the interface is defined by the signed distance to the interface. φ＜0 in the upper fluid and φ＞0 in the lower fluid. The governing equation of the level-set function is

The curvature of the interface κ may be derived from the level-set function as

The level-set function can also be used to represent the properties such as the density and the viscosity of different fluids in a unified form. To avoid the sharp jump of fluid properties across the interface, a smoothed Heaviside function Hε(φ,ε) may be introduced to allow a continuous variation of the relevant parameters in a small region near the interface

where the subscriptsu and ldenote the up per and lower parts of the fluid, respectively, f denotes the property of the fluid, ε is the half thickness of the smoothed region, Hε(φ,ε) is defi- ned as

It is convenient to use Eq.(4) to represent the interface, but the conservation of mass separated by φ=0 is not guaranteed. To solve this problem, a CLSVOF method was proposed[13,14], where one also considers the VOF fraction F, which is defined such that a cell is filled with only the upper fluid if F=0, with only the lower fluid if F=1, and with both if 0＜F＜1. The governing equationof F is

A second order operator-split scheme is used to solve φ and F as functions of the computational time from Eqs.(4) and (8). To ensure the mass conservation, the discrete fluxes of F need be determined from the reconstructed geometry of the interface. By a piecewise function of F, asφ is used to determine the slope, the interfaceis reconstructed. After obtaining the updatedφ and F, the geometry of the interface can also be updated. The re-evaluation of the level-set function is thus possible through computing the signed distance to the reconstructed interface.

1.3 Validation

To verifythe numerical model, the breaking of a deep-water Stokes wave is computed. The computa-computational domain is set to be λ×λ (with λ being the wavelength) with a 512×512grid. The com-to capture the water surface while the CLSVOF method is employed in the present study. The two results are shown to be in a good agreement.

2. Wave breaking over muddy seabed

2.1 Computational conditions

We consider the breakingof the sur face wave on a layer of water with depth h over a layer of mud with depth d (see Fig.2(a), which will be referred to as the water-mud case hereinafter). For a comparison, we also consider the water waves with the same period breaking over a rigid bottom (see Fig.2(b)). Two different water depths, i.e., h+d and h, which will be referred to as the pure-water Cases1 and 2, respectively, are assumed for detailed discussions on the effects of the mud. It is worthwhile to mention that, for the water-mud case, the two level-set functions are necessary to capture both the air-water interface and the water-mud interface. In all cases, a wave of the laboratory scale with period T=1s is considered. In the water-mud case, the depths ofthe air, the water and the mud, when normalized by the wavelengthλ, are set to be 0.21, 0.21 and 0.063, respectively.The density and viscosity ratios of air to water are ρa(bǔ)/ρw=0.0012 and μa/μw=0.0185, respectively, while the corresponding ratios of mud to water are chosen as ρm/ρw=1.21 and μm/μw= 4840. Governed by the dispersion relation, the wavelengthis λ＇=1.01λ in the pure-water case 1 and λ＇=0.96λ in the pure-water Case 2. In all computations, the domain is set to be 1 wavelength in the horizontal direction. The grid size is 1/384 of the wavelength. Periodic condition is assumed at the lateral boundaries while the non-slip condition is assumed at the bottom, and the free-slip condition at the top.

Fig.1 Computational results of a deep-water wave plunging

For the water-mud case, the initial velocity and the interface profiles are prescribed based on thelinear theory of Dalrymple and Liu[3]. For the purewater case, the initial conditions are given according to the small amplitude wave theory. According to Lubin et al.[12], Chen et al.[18], Iafrati[19], the initial wave steepness is set to be larger than usual (being 0.36 in the present study) so that the wave will break.

Fig.2 Schematic diagrams

2.2 Wave profile and flo w velocity

Figure 3 shows the computedwater surface profileat different time instants. The breaking process can be summarized as follows. Starting from a large initial slope, the wave becomes steeper and steeper (Fig.3(a)). Then, a jet is formed near the wave crest (Fig.3(b)). The jet shoots forward and plunges onto the water surface (Fig.3(c)).

Compared with the water-mud case and the purewater Case 1, the onset of the breaking is evidently earlier in the pure-water Case 2. This is not unexpected because in the pure-water Case 2, the total fluid depth is smaller. Since the mud modeled here is soft, the onset of the breaking is close in the water-mud case and the pure-water Case 1. But still, the breaking in the water-mud case is a little delayed as compared with the pure-water Case 1, which might be attributed to the energy transfer from water to mud.

The computed streamwise velocities attimet= 0.78T for the water-mud case and the pure-water Case 1are shown in Fig.4. Near the wave crest, the velocity is very large and the maximum value reaches 1.32, or, 32% larger than the initial wave celerity in the water-mud case. Within the mud, the streamwise velocity under the crest of the water-mud interface is much larger than that under the trough. For the purewater Case 1, a generally similar distribution of the streamwise velocity can be observed, except in the region occupied by mud in the other case.

Fig.3 Computed surface profile of a breaking wave. Note that a vertical shift of 0.05 is made to the surface profile of the pure-water Case I in (c) to make the comparison clear

To study the vertical variation of the streamwise velocity, the vertical distribution of the streamwise velocity at x=0.9 are shown in Fig.5 as an example. For the water-mud case, the streamwise velocity decreases as the depth increases, but it increases slightly at the water-mud interface. It decreases again to the bottom within the mud layer. For the pure-water case, the streamwise velocity decreases downward from the water surface monotonically. Near the bottom, it drops to zero.

2.3 Energy dissipation

The total energy of the water-mud system can be expressed aswhere the potential energy is defined with the still water level as the datum. The computed time variations of the total energy for the water-mud case, the pure-water Case 1, and the pure-water Case 2 are shown in Fig.6. In all cases, the total energy is shown to have very different decaying rates at three different stages, i.e., before, during and after the breaking, which are in a general agreement with the results of the existing studies, including Lubin et al.[12], Chen et al.[18], Iafrati[19]. Before the breaking, the energy decreases slowly. Viscous dissipation is responsible for the decay at this stage. The energy decreases faster in the water-mud case as compared with the purewater cases, because the mud has a much larger viscosity. During the breaking, the energy decays very quickly. At this stage, in addition to the viscous and turbulent dissipations, the work done against the buoyancy during the downward transport of the entrapped air also dissipates a large amount of energy[20,21]. Though the viscous dissipation is larger in the watermud case as compared with the pure-water cases, the total energy decay rates during the breaking in the water-mud case and the pure-water Case 1 are nearly the same, indicating the breaking intensity in the water-mud case is weakened. The energy decay rate during the breaking is the largest in the pure-water Case 2, because the breaking is more violent when the water depth is smaller. After the breaking, since most of the wave energy is dissipated, the energy decay process slows down. The energy dissipation is relatively faster in the water-mud case at this stage because the viscous and turbulent dissipations dominate the phenomenon again.

Fig. 4 Contours of the streamwise velocity at the timet= 0.78T. The water-mud interface is marked bythe dash line

Fig.5 Vert ical profile of streamwise velocity atx=0.9for the water-mud case and the pure-water Case 1 in Fig.4

Fig.6 Time variation of the total energy for the water-mud case, the pure-water Cases 1 and 2. The results are normalized by the initial energy. The dots correspond tothe time instants in Fig.3

Fig.7 Variation of viscous and turbulent dissipations.The dissipations are normalized by σE(0) with σbeing the wave frequency

In order to investigat the role played by the viscous and turbulent dissipations, we consider the integrals of viscous and turbulent dissipations in the waterand in the mud, separately

where ε=μt(?ui/?xj+?uj/?xi)?ui/?xj. Figure 7 shows the time variation of Dw(t) and Dm(t) for the water-mud case. It is seenthat Dw(t)increases sharply as the jet touches downontothe water surface where a strong rotational flow is generated. Then, it fluctuates with the splash-ups during the breaking. The decay slows down after the breaking as most of the energy is dissipated. A similar variation of Dw(t) was also observed by Iafrati[19,21]. From Fig.7(b), it is seen that Dm(t) decreases rather continuously with time. However, as compared with Dw(t), the decrease of Dm(t) is an order larger due to the large viscosity of the mud.

Fig.8 Variation of viscous and turbulent dissipationsin the water Dw(t)in the water-mud case and in the purewater Case 1. The dissipations are normalized by their corresponding σE(0) with σ being the wave frequency

In Fig.8,Dw(t)in the water-mud case is compared with that in the pure-water Case 1. The variationsof Dw(t) in the two cases are shown to be rather similar, except Dw(t) has smaller fluctuation in the water-mud case, which again indicates that the breaking intensity is weakened when the mud is present.

3. Conclusions

In the present study, the LES of the wave breaking over a viscous muddy seabed is carried out for thefirst time, based on the filtered Navier-Stokes equations and a coupled level-set and VOF method for capturing the interfaces. Simulations of the wave breaking over a rigid bottom are also conducted for comparison.

It is found that the muddy bed has significant influences on the wave breaking. The existence of a mud layer beneath an otherwise rigid bottom is shown to have a similar effect as the increase of the water depth. When compared with that on a rigid bottom, the inception of the wave breaking is evidently delayed and the breaking intensity is weakened. The results are also compared with the wave breaking over a rigid bottom with the water depth equal to the mud depth plus the water depth over the muddy bed. Still, the onset of breaking is slightly delayed and the breaking intensity is weakened due to the dissipation within the mud. The variation of the total energy has different decay rates before, during and after the breaking. Before and after the breaking, the energy decays more quickly in the water-mud case due to the large viscosity of the mud. During the breaking, in addition to the fluid turbulence and the viscous effect of the mud, the air entrainment also dissipates a large amount of energy. The energy dissipation in the water increases sharply with the breaking, while that in the mud varies continuously. But, the dissipations in the mud always are in a dominant amount.

Acknowledgement

The first author thanks Professor Shen L. and Professor Dalrymple R. A. for their kind help during her stay in Johns Hopkins University.

[1]SHEREMET A., STONE G. W. Observation of nearshore wave dissipation over muddy sea beds[J]. Journal of Geophysical Research, 2003, 108(C11): 3357- 3368.

[2]ELGAR S., RAUBENHEIMER B. Wave dissipation by muddy seafloors[J]. Geophysical Research Letters, 2008, 35: L07611.

[3]DALRYMPLE R. A., LIU P. L. F. Waves over soft muds: A two-layer fluid model[J]. Journal of Physical Oceanography, 1978, 8(6): 1121-1131.

[4]NG C. O. Water waves over a muddy bed: A two layer Stokes’ boundary layer model[J]. Coastal Engineering, 2000, 40(3): 221-242.

[5]ZHANG Qing-he, WAI O. W. H. and LEE J. H. W. Wave attenuation over mud bed: A pseudo-plastic model[J]. Journal of Hydrodynamics, 2003, 15(6): 32- 38.

[6]LIU P. L. F., CHAN I.-C. A note on the effects of a thin visco-elastic mud layer on small amplitude water-wave propagation[J]. Coastal Engineering, 2007, 54(3): 233- 247.

[7]KAIHATU J. M., SHEREMET A. and HOLLAND K. T. A model for the propagation of nonlinear surface waves over viscous muds[J]. Coastal Engineering, 2007, 54(10): 752-764.

[8]NIU Xiao-jing, YU Xi-ping. A practical model for the Decay of random waves on muddy beaches[J]. Journal of Hydrodynamics, 2008, 20(3): 288-292.

[9]AZAM M. H., MOKHTAR N. An experimental investigation of broken wave over mud bed[J]. Journal of Coastal Research, 2000, 16(4): 965-975.

[10]SOLTANPOUR M., SHIBAYAMA T. and NOMA T. Cross-shore mud transport and beach deformation model[J]. Coastal Engineering Journal, 2003, 45(3): 363-386.

[11]WATANABE Y., SAEKI H. and HOSKING R. J. Three-dimensional vortex structures under breaking waves[J].Journal of Fluid Mechanics, 2005, 545: 291-328.

[12]LUBIN P., VINCENT S. and ABADIE S. et al. Threedimensional large eddy simulation of air entrainment under plunging breaking waves[J]. Coastal Engineering, 2006, 53(8): 631-655.

[13]SUSSMAN M., PUCKETT E. G. A coupled level set and volume-of-fluid method for computing 3D and Axisymmetric incompressible two-phase flow[J]. Journal of Computational Physics, 2000, 162(2): 301-337.

[14]COYAJEDD E., BOERSMA B. J. Numerical simulation of drop impact on a liquid-liquid interface with a multiple marker front-capturing method[J]. Journal of Computational Physics, 2009, 228(12): 4444-4467.

[15]HUAI Wen-xin, WU Zhen-lei and QIAN Zhong-dong et al. Large eddy simulation of open channel flows with non-submerged vegetation[J]. Journal of Hydrodynamics, 2011, 23(2): 258-264.

[16]YAKHOT V., ORSZAG S. A. Renormalization group analysis of turbulence. (I. Basic theory)[J]. Journal of Scientific Computing, 1986, 1(1): 3-51.

[17]FU T. C., RATCLIFFE T. and O’Shea T. T. et al. A comparison of experimental measurements and computational predictions of a deep-V planninghull[C]. Proceeding of the 28th Symposium on Naval Hydrodynamics. Pasadena, CA, USA, 2010.

[18]CHEN G., KHARIF C. and ZALESKI S. et al. Twodimensional Navier-Stokes simulation of breaking waves[J]. Physics of Fluids, 1999, 11(1): 121-133.

[19]IAFRATI A. Numerical study of the effects of the breaking intensity on wave breaking[J]. Journal of Fluid Mechanics, 2009, 622: 371-411.

[20]BLENKINSOPP C. E., CHAPLIN J. R. Void fraction measurements in breaking waves[J]. Proceedings of the Royal Society A, 2007, 463(2088): 3151-3170.

[21]IAFRATI A. Energy dissipation mechanisms in wave breaking process: Spilling and highly aerated plunging breaking events[J]. Journal of Geophysical Research, 2011, 116: C07024.

January 28, 2012, Revised February 18, 2012)

* Project supported by the National Natural Science Foundation of China (Grant No. 51109119), the State key Laboratory of Hydro-science and Engineering, Tsinghua University (Grant No. 2011-KY-1).

Biography: HU Yi (1986-), Female, Ph. D. Candidate

YU Xi-ping,

E-mail: yuxiping@tsinghua.edu.cn