- Email: [email protected]
Wave deformation and vortex generation in water waves propagating over a submerged dike
Coastal Engineering 37 Ž1999. 123–148 www.elsevier.comrlocatercoastaleng
Wave deformation and vortex generation in water waves propagating over a submerged dike Ching-Jer Huang ) , Chih-Ming Dong Department of Hydraulics and Ocean Engineering, National Cheng Kung UniÕersity, 1 Dah Shyue Road, Tainan 70101, Taiwan Received 12 August 1998; received in revised form 25 January 1999; accepted 16 February 1999
Abstract The Navier–Stokes equations and the exact free surface boundary conditions are solved to simulate wave deformation and vortex generation in water waves propagating over a submerged dike. Incident waves are generated by a piston-type wavemaker set up in the computational domain. Numerical results are compared with experimental data in order to confirm the validity of the numerical model. The fast Fourier transform and a wave resolution technique are applied to decompose the transformed waves and the higher harmonics. Effects of different parameters on wave transformation and vortex generation are studied systematically. These parameters include the Ursell number, the Keulegan–Carpenter number, the water depth ratio, the Reynolds number, the length aspect ratio of the dike, and the type of dike. q 1999 Elsevier Science B.V. All rights reserved. Keywords: Submerged dike; Wave deformation; Higher harmonics; Flow separation; Vortex generation
1. Introduction There are two main physical phenomena associated with the propagation of water waves over a submerged dike. The first is the generation of higher harmonics and the second is flow separation and vortex generation on the weather side and on the lee side of the obstacle. The former phenomenon has long been known from direct field observations and laboratory measurements Že.g., Johnson et al., 1951; Jolas, 1960.. Earlier researchers have tried to predict the generation of higher harmonics using )
Corresponding author. Tel.: q886-6-2757575 ext. 63252; Fax: q886-6-2741463; E-mail: [email protected] 0378-3839r99r$ - see front matter q 1999 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 8 - 3 8 3 9 Ž 9 9 . 0 0 0 1 7 - 4
124
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
nonlinear shallow-water wave theories, such as the Boussinesq and the KdV equations. Since the derivation of the Boussinesq equations is based on the assumptions of both weak nonlinearity and weak dispersivity of waves, these equations may not be valid for the prediction on the lee side of the dike, where higher harmonics may arise in the form of deep-water waves. In order to overcome this defect, improvements on the Boussinesq equations have been developed. Peregrine Ž1967. developed equations of motion for long waves in water of varying depth. These equations can be viewed as the Boussinesq equations for water of varying depth. Madsen et al. Ž1991. and Madsen and Sorensen Ž1992. improved the Boussinesq equation by adding a term to improve its dispersion characteristics. Battjes and Beji Ž1991. derived a new set of Boussinesq equations by combining the method used by Madsen et al. Ž1991. and the equation proposed by Peregrine Ž1967.. Beji and Battjes Ž1994. also employed these equations to study the deformation of water waves passing over a dike. The numerical results of wave height were in good agreement with experimental data. The boundary element method has also been applied to solve the Laplace equation and the nonlinear free surface boundary conditions in order to study the decomposition phenomena of waves passing over a rectangular step ŽKittitanasuan et al., 1993. and a rectangular submerged dike ŽOhyama and Nadaoka, 1992, 1994.. In Ohyama and Nadaoka, the incident waves were generated by a numerical wave tank model developed by Ohyama and Nadaoka Ž1991.. This is similar to the present method of setting up a piston-type wavemaker in the computational domain to generate the incident waves. The viscosity effect has been neglected in the aforementioned theories. Most of the studies concentrated on the transformation of waves passing over submerged dikes. As far as the wave transformation is concerned, the results are satisfactory. However, the prediction of the flow fields without taking the viscosity effect into account may in some situations fail to provide useful knowledge about the real flow fields. For example, the vortices exiting in front of and behind dikes in the nearshore region may have an important effect on the safety and stability of the dikes, but they cannot be predicted by earlier theories. Recently, Ting and Kim Ž1994. conducted laboratory experiments to investigate flow separation and vortex generation induced by waves propagating over a submerged rectangular dike. Velocity measurements were performed using a two-component laser-doppler anemometer ŽLDA.. Significant results were obtained, which showed that the Keulegan–Carpenter ŽK–C. number is the most important parameter in determining the formation and development of vortices. Based on the dimensional analysis, it was also expected that the zone of flow separation at the corners of the dikes would be larger at a higher Reynolds number. Ting and Kim’s measurements showed a clear picture of the flow fields around a submerged rectangular dike. Due to the complex nature of the flow around the dike, Ting and Kim suggested that it would be very difficult to determine flow separation effects unambiguously without solving the viscous flow equations in the near field. The Navier–Stokes equations and the exact free surface boundary conditions must be solved in order to predict the realistic flow fields around the submerged dikes. In the present study, the numerical wave tank model developed by Huang et al. Ž1998. was used to simulate wave deformation and the flow fields as waves propagated over submerged trapezoidal and rectangular dikes. In this model the unsteady two-dimen-
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
125
sional Navier–Stokes equations were discretized by means of a finite-analytical scheme developed by Chen and Chen Ž1982.. The SUMMAC method ŽChan and Street, 1970. was used to determine the location of the free surface. In order to verify the accuracy of the numerical calculations, numerical results were compared with the experimental data of Beji and Battjes Ž1994. for a trapezoidal dike and with those of Driscoll et al. Ž1992. for a rectangular dike. The fast Fourier transform ŽFFT. and Grue’s method ŽGrue, 1992. were also applied to decompose the transformed waves and the higher harmonics. Effects of different parameters on wave transformation, flow separation, and vortex generation were studied and compared systematically. These parameters include the Ursell number, the K–C number, the water depth ratio of the dike, the Reynolds number, the length aspect ratio of the dike, and the type of dike.
2. Governing equations and boundary conditions The wave transformation and flow fields in water waves passing over a submerged dike are studied in the present paper. A schematic diagram of a trapezoidal dike is shown in Fig. 1. For a rectangular dike, only the rectangular part Žindicated by dotted lines. of the trapezoid is kept. A piston-type wavemaker with stroke S 0 is located at x s 0 and generates the incident waves. The still water depth is h 0 . The trapezoidal dike consists of an upslope of S 1 and a downslope of S 2 . The shallow water depth above the dike has a depth of qh 0 Ž0 - q - 1.. For an incompressible, viscous fluid, the dimensionless continuity equation in the Cartesian coordinate system is Eu
EÕ q
Ex
Ey
s0
Ž 1.
and the Navier–Stokes equations are given by Eu Et
Eu qu
Ex
Eu qÕ
Ep sy
1 q
Ey
Ex
Re
EÕ
Ep
1
ž
E2 u Ex2
E2 u q
E y2
/
Ž 2.
/
Ž 3.
and EÕ Et
EÕ qu
Ex
qÕ
sy Ey
q Ey
Re
ž
E2 Õ Ex2
E2 Õ q
E y2
where u and Õ are the horizontal and vertical velocity components, t is the time, p is the hydrodynamic pressure, and R e is the Reynolds number, defined as Re s
u0 h0
y
Ž 4.
where y is the kinematic viscosity of the fluid, and u 0 is the velocity–amplitude of the wavemaker, which moves back and forth with a velocity of u s u 0 sin v t. In the present case, u 0 and h 0 are used to non-dimensionalize the velocity and length, while t 0 s h 0ru 0 is chosen to non-dimensionalize the time.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
126
Fig. 1. Schematic diagram of submerged trapezoidal and rectangular dikes.
The boundary conditions which must be satisfied to solve the problem shown in Fig. 1 are: Ž1. the kinematic and the dynamic free surface boundary conditions at the free surface, Ž2. the no-slip boundary condition at the channel bottom, Ž3. the upstream boundary condition on the wavemaker plate, and Ž4. the radiation condition on the outgoing boundary. The initial conditions of the velocities, hydrodynamic pressure, and surface displacements are set at zero at time t s 0. Conditions Ž1., Ž3., and Ž4. mentioned above are explained in more detail below. 2.1. Free surface boundary conditions The kinematic condition states that fluid particles at a free surface remain on the free surface, and can be expressed as Eh Et
Eh qu
Ex
Ž 5.
sÕ
where h s h Ž x,t . is the location of the free surface. The dynamic condition requires that, along the free surface boundary, the normal stress is equal to the atmospheric pressure and the tangential stress is zero. These conditions are expressed by Huang et al. Ž1998. as the following Eh
2 1q
h p0 s
Fr
2
q
2
ž / ž / Ex
Eh
Re 1 y
EÕ 2
Ey
Ž 6.
Ex
and Eu
EÕ sy
Ey
Ex
EÕ Eh
4 q Eh
ž / Ex
2
Ey Ex y1
Ž 7.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
127
where p 0 s pŽ x,h . is the hydrodynamic pressure at the free surface and Fr is the Froude number, defined as u0 Fr s Ž 8. gh 0
(
In numerical computations, Eq. Ž6. is used to determine the pressure at the free surface and Eq. Ž7. is used to extrapolate the horizontal velocity at the free surface from the flow domain. The vertical velocity component Õ is then calculated from the continuity equation using the known velocity component u, obtained from Eq. Ž7.. 2.2. The upstream boundary condition The no-slip boundary condition at the wavemaker is released at two nodal points beneath the free surface. The slip velocities Õ along the wavemaker are extrapolated from the flow inside the domain. 2.3. The downstream boundary condition The downstream boundary condition requires that, at a large distance from the wavemaker, the wave is outgoing without reflection. According to the wave equation and continuity equation, we can set Ep Et
qc
Ep Ex
Eu s 0,
Et
Eu qc
Ex
Eu s 0,
EÕ q
Ex
Ey
s0
Ž 9.
where c is the phase speed of the wave and is determined from the dispersion relation.
3. Numerical method 3.1. Discretization of goÕerning equations In the present study, the finite-analytic ŽFA. method is applied to the discretization of the unsteady two-dimensional Navier–Stokes equations. In this method, local analytical solutions obtained from the linearized Navier–Stokes equations for the discretized computational elements are incorporated into the numerical method. The value of a variable at a particular node is expressed in terms of the values of neighboring nodes. The simplified four-point version of the finite-analytic method developed by Chen and Patel Ž1987. was used here to discretize the Navier–Stokes equations. Please refer to Chen and Patel for more information about the derivation of the FA coefficients. The Navier–Stokes Eqs. Ž2. and Ž3. can be written in linearized form as E2 f Ex
2
E2 f q Ey
2
s2 A
Ef Ex
q2 B
Ef Ey
q Re
Ef Et
q Sf
Ž 10 .
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
128
where f represents the velocity components u or Õ, and As
R eU 2
Su s R e
,
Ep Ex
Bs ,
R eV 2
SÕ s R e
Ep
Ž 11 .
Ey
where U, V are mean values of u and Õ in a given computational element. Eq. Ž10. can be rewritten as E2 f Ey
2
s2 B
Ef Ey
q Gf
Ž 12 .
where Gf is the nonhomogeneous term Gf s 2 A
Ef Ex
q Re
Ef Et
q Sf y
E2 f Ex2
Ž 13 .
Eq. Ž12. is a nonhomogeneous second-order ordinary differential equation and can be easily solved to obtain the analytic solution, namely f s a Ž e 2 B y y 1 . q by q c
Ž 14 .
By using the values of the variables at a node P and at the two neighboring nodal points N and S in the y-direction ŽFig. 2., the constants a, b, and c can be determined. An equation similar to Eq. Ž12. can be obtained by substituting Eq. Ž14. into Eq. Ž12., but with variable y replaced by x. Following procedures similar to those used to solve the differential equation in Eq. Ž12., one obtains f p s Ž Cs fs q Cn f n q C w f w q Ce fe q Ct f pny 1 . q Cp S f
Ž 15 .
where the coefficient C with different subscripts are the FA coefficients, which relate to the time step, the grid size, and the parameters appearing in Eq. Ž10., such as A, B, and
Fig. 2. Symbols used in the present numerical method.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
129
the Reynolds number. In obtaining Eq. Ž15., the time derivative is approximated by an implicit, backward difference scheme. From Eq. Ž15. we see that the value of f at node P is related to the value of f at the neighboring four nodes, to the value of f in the previous time step, and to the pressure gradient. The variable f can be determined by applying Eq. Ž15. to all the grid nodes in the flow domain and using the given initial and boundary conditions. The SIMPLER algorithm developed by Patankar Ž1980. was used to calculate the coupled velocity and pressure fields. A staggered numerical grid was used. The velocity components u and Õ are defined as located at the boundaries of a control volume, while the pressure p is located at the center. 3.2. Treatment of the boundary The coordinate transformation technique was not applied in the computation. In order to make a rectangular grid system for a trapezoidal dike, the smoothly varying bed topography is divided into a series of small rectangular steps with a horizontalrvertical ratio of 20:1 or 10:1 ŽFig. 3.. The Marker and Cell method ŽMAC. and its modified version SUMMAC are used in combination to calculate the free surface boundary. The major concept underlying the MAC method is the use of marker particles to identify the location of the free surface. By tracking the positions of these marker particles, the transient location of the free surface can be determined. Pressure at the free surface is calculated by the normal dynamic free surface boundary condition ŽEq. Ž6.., using the known velocities u and Õ. The velocity component u at a free surface cell is extrapolated from the velocity of the main body of the fluid by means of the tangential dynamic boundary condition ŽEq. Ž7... The value of Õ at a free surface cell can be extrapolated from the fluid domain using an extrapolation formula. However, the velocity components obtained in this way may result in a violation of the continuity equation at local cells. To avoid this risk, we use the continuity equation to determine the velocity component Õ. A method which satisfies the discretized Poisson equation for the pressure and is called ‘irregular star’ is applied to determine the pressure at the cells near the free surface.
Fig. 3. Numerical grids in the computational domain of a trapezoidal dike.
130
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
In the numerical model, the Eulerian grid system is used, which means the form and size of the computational grids are fixed. However, the computational domain near the wavemaker changes as the wavemaker moves back and forth. In the computation the grids in the neighborhood of the wavemaker are arranged such that the boundary of the computational domain matches the location of the wavemaker. The numerical results shown in the following were calculated with dimensionless D x and D y Žnormalized with respect to h 0 . equal to 0.1, except for the region near the free surface where the finest grid size D y was 0.025. One example of the grid system is shown in Fig. 3. The increment of dimensionless time, normalized with respect to t 0 , was chosen to be 0.001. 4. Water waves passing over a trapezoidal dike The numerical methods described in Section 3 are used to simulate the deformation of waves propagating over submerged trapezoidal and rectangular dikes. The incident waves are generated by a piston-type wavemaker imposed in the computation domain. To generate the small-amplitude waves, the strokes of the numerical wavemaker are determined by linear wavemaker theory ŽDean and Dalrymple, 1984.. To generate relatively long waves of finite amplitude and to make sure that the generated waves are of permanent form, Madsen’s second-order wavemaker theory ŽMadsen, 1971. is applied to determine the wavemaker motion. The accuracy of the numerical wave tank model applied in the present study has been verified and discussed in detail in Huang et al. Ž1998.. In an experimental study of water waves propagating over a rectangular submerged dike, Ting and Kim Ž1994. indicated by dimensional analysis that the normalized water surface elevation is a function of water depth ratio, q, the length aspect ratio of the dike, lrh 0 , the Ursell number, Hi L2rh 30 , the K–C number, Hi Lrh20 , and the Reynolds number, which was defined to be gl 3 ry , in which h 0 is the stillwater depth, and Hi and L are the incident wave height and wavelength. The Ursell number represents the nonlinearity of the incoming waves. The K–C number is interpreted as the ratio of maximum displacement of a fluid particle from its neutral position to the body length, and is an important parameter in determining the formation and development of vortices in the flow. In accordance with the above ideas, we have chosen three different incident waves passing over a trapezoidal and a rectangular dike in order to study the main characteristics of waves passing over submerged dikes. Detailed information about the incident waves and the dikes are summarized in Table 1 for case 1 through case 10. Cases 1 to 3 in Table 1 refer to a trapezoidal dike and are discussed here, while the other cases are for a rectangular dike and are left to Section 5. It should be noted that we have confined ourselves to nonbreaking waves in all the numerical simulations. In the following sections, we simulate the transformation of wave passing over a trapezoidal dike and compare the results with the experimental data of Beji and Battjes Ž1993, 1994.. The FFT and Grue’s method are then applied to decompose the transformed waves and the higher harmonics. Effect of the Ursell number on the wave transformation is also discussed.
(
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
131
Table 1 Numerical conditions Variable Case 1
Case 2
Case 3
Case 4
Case 5
Case 6
Case 7
Case 8
Case 9
Case 10
h 0 Žm. Hi Žcm. L Žm. T Žs. Ur K–C Re
0.4 2.0 3.695 2.0 4.27 0.46 3.4=10 5
0.4 4.0 3.695 2.0 8.53 0.92 6.8=10 5
0.5 2.5 3.327 1.7 2.21 0.33 3.3=10 5
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
0.4 2.0 3.695 2.0 4.27 0.46 3.4=10 5
0.4 4.0 3.695 2.0 8.53 0.92 6.8=10 5
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
0.6096 5.3 9.529 4.0 21.24 1.36 2.0=10 6
0.25 2.0 1r20 1r10
0.25 2.0 1r20 1r10
0.24 0.79
0.25 2.0
0.25 2.0
0.25 2.0
0.5 2.0
0.25 3.8
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
Dike type A q l Žm. S1 S2
0.25 2.0 1r20 1r10
B 0.5 0.6096
A: trapezoidal; B: rectangular. Ur s Hi L2 r h 30 ; K–C s Hi Lr h20 ; R e s Hi L2 ry h 0 T.
4.1. Comparison with experiments In the experiment of Beji and Battjes Ž1993, 1994. the stillwater depth was set at 0.4 m. A piston-type wavemaker was installed at one end of the flume. The model of the submerged trapezoidal dike consisted of an upslope of 1:20 and a downslope of 1:10. The water depth was 0.1 m in the shallowest region above the horizontal part Ži.e., q s 0.25.. The wave profiles were measured at 7 different locations as shown in Fig. 4. The incident wave height Hi and period T were 2.0 cm and 2.0 s, respectively. The incident wavelength L was 3.695 m based on linear theory. The Ursell number Hi L2rh 30 was 4.27. Fig. 5 shows the numerical and the experimental results of the wave profiles
Fig. 4. Definition sketch of wave tank and locations of wave gages Žfrom Beji and Battjes, 1994..
132
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 5. Wave profiles at several stations for monochromatic incident waves with Hin s 2.0 cm, T s 2.0 s, . Results of present numerical model; q s 0.25. Žv . Experimental data ŽBeji and Battjes, 1994.; Ž Ž- - -. BEM ŽOhyama et al., 1994..
obtained at stations 1 to 6. The numerical results obtained by Ohyama et al. Ž1994. using the BEM model are also presented in Fig. 5 for comparison. From Fig. 5 we see that both our numerical results and those obtained by Ohyama et al. are in good agreement with the experimental data. The reason that our results are not as good as those of Ohyama may be due to the simplified treatment of the boundary. The wave profiles shown in stations 1 to 6 indicated that waves undergo rapid variation as they propagated over the dike. The main reason for this is the generation of higher harmonics and energy transfer from the bound waves to the free waves of the higher harmonics. This will be discussed in more detail in Section 4.2. 4.2. Basic mechanism of waÕe decomposition The numerical results for free surface elevations along the wave channel at different times are shown in Fig. 6. It can be seen from Fig. 6 that as waves propagate over the dike, the primary wave crests become steeper and a dispersive tail gradually develops. A small wave, indicated by an arrow in Fig. 6, appears at the trailing edge of a primary
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
133
Fig. 6. Surface elevations along the wave channel at several times for waves with Hin s 2.0 cm, T s 2.0 s, and q s 0.25 passing over a submerged trapezoidal dike.
wave. Because the small wave propagates more slowly than the main crest, it gradually detaches from the main crest and is overtaken by the next wave. To understand the basic mechanism of wave evolution and decomposition, the method of Grue Ž1992. is used to determine the amplitude of high frequency components. This method is based on the resolution technique of two-wave gage systems and is described briefly as follows.
134
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
We assume that on the lee side of the dike the waves are given by Ž1. Ž1. hq Ž x , t . s aq cos Ž Kx y v t q d q . q Ý aŽlqn. cos n Ž Kx y v t q d qŽ 1 . .
n)1 Ž n. q Ý aq cos Ž K n x y n v t q dqŽ n. .
Ž 16 .
n)1 n. where K represents the wavenumber of the first-harmonic wave, aŽlq , n s 2,3, . . . , Ž n. denote the amplitudes of the nth-harmonic bound waves, aq and K n with n s 2, 3, . . . denote the amplitudes and the wavenumber of the transmitted nth-harmonic free waves,
Fig. 7. Comparison of wave number–frequency spectra for monochromatic incident waves with Hin s 2.0 cm, T s 2.0 s, q s 0.25.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
135
and dqŽ n., n s 1, 2, . . . , are the phase angles. From Eq. Ž16. we see that the transformed waves on the lee side of the dike can be decomposed into first harmonic waves and nth-harmonic bound and free waves. The bound waves have the same phase velocity as the first harmonic wave and are bound to it. Hence, they cannot be distinguished simply from the wave profiles. The free waves of higher harmonics satisfy the dispersion relation and propagate at different velocity from the first harmonic and will be detached from it as the wave propagates. If we define the Fourier transform of the surface elevation h Ž x,t . as
hˆ Ž n. Ž x . s
v
2p r v
H 2p 0
h Ž x , t . exp Ž yin v t . d t ,
n s 1, 2, . . . ,
Ž 17 .
then the amplitudes of the free and bound waves can be determined as follows by giving hˆ Ž n. at two different locations, namely x 1 and x 1 q D x: Ž n. aq s
1 sin Ž 1r2 Ž K n y nK . D x .
hˆ Ž n . Ž x 1 . y hˆ Ž n. Ž x 1 q D x . exp Ž inKD x . ,
n s 2, 3, . . . n. aŽlq s
1 sin Ž 1r2 Ž K n y nK . D x .
n s 2, 3, . . .
hˆ Ž n . Ž x 1 . y hˆ Ž n. Ž x 1 q D x . exp Ž iK n D x . ,
Ž 18 .
The wavenumber–frequency spectra of waves above the submerged dike and on the lee side of the dike are shown in Fig. 7Ža. and Žb., respectively. The results shown in Fig. 7Ža. are obtained using the surface elevations simulated at stations 3 and 4, while those shown in Fig. 7Žb. use the data of stations 7 and 8, which is at a distance of 0.4 m shoreward from station 7. In Fig. 7, a n and a1 are the nth-harmonic and first harmonic
Fig. 8. Spatial evolution of harmonic amplitudes for a monochromatic incident wave with Hin s 2.9 cm, T s1.25 s, Ur s1.84, q s 0.25. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
136
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 9. Spatial evolution of harmonic amplitudes for a monochromatic incident wave with Hin s 2.0 cm, T s 2.0 s, Ur s 4.27, q s 0.25. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
wave amplitudes, f 0 is the frequency of the incident wave, K is the wave number, and h is the local water depth. In each frf 0 y Kh plane, a curved line ŽF. represents the dispersion relation, and a straight line ŽB. is drawn from the origin to a point at which the ratio of wave amplitude for the first harmonic Ž frf 0 s 1. is maximum. The wave components lying on line ŽF. are free waves, and we can determine their phase velocities from the dispersion relation Ž4p 2 f 2 s gk tanh kh.. The phase velocity of the wave components on line ŽB. is the same as that of the incident wave. The free waÕes haÕe shorter periods than incident waÕes; thus, from the dispersion relation, their phase Õelocities are slower than those of incident waÕes. Hence, the free wave is detached from the main crest and is overtaken by the next wave as shown in Fig. 6. From Fig. 7 we noted also that above the submerged dike, the wave amplitude of the free wave and bound wave are about the same. However, on the lee side of the dike, the amplitude of the bound wave almost disappears. The main reason for this is that as the
Fig. 10. Spatial evolution of harmonic amplitudes for a monochromatic incident wave with Hin s 4.0 cm, T s 2.0 s, Ur s8.53, q s 0.25. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 11. Surface elevation and velocity field for case 1; Hin s 2.9 cm, T s1.25 s.
137
138
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
waves propagate into the deep water region where wave nonlinearity is so weak that the bound wave is reduced substantially, a large amount of energy is transferred from the bound wave to the free wave. 4.3. Spatial amplitude modulations of waÕes and the flow field Figs. 8–10 show the spatial evolution of the first, second and third harmonic wave amplitudes Ž n s 1, 2, 3. for cases 1, 2, and 3. The Ursell numbers in cases 1 through 3 are 1.91, 4.27, and 8.53, and the K–C numbers are 0.37, 0.46, and 0.92, respectively. We noted from Fig. 8 through 10 that on the lee side of the dike the higher harmonic wave amplitudes become larger as the Ursell number increases. In case 3, the second harmonic wave amplitude is even 55% larger than that of the first harmonic wave. The instantaneous velocity fields for case 1 are shown in Fig. 11Ža,b. for relative time trT s 0 and 1r2. In this figure, there is no flow separation or vortex generation around the dike. This may be due to the gentle slopes located in front of and behind the dike. The flow fields for cases 2 and 3 are very similar to that of case 1, and no flow separation or vortex generation occurs. 5. Waves passing over rectangular dikes Before we proceed to study the effect of different parameters on the main characteristics of water waves passing over a rectangular dike, the accuracy of the numerical results are verified. Fig. 12 shows a comparison of numerical and experimental results at several stations and times for monochromatic incident waves passing over a rectangular dike. Information about the incident wave and the dike are specified in case 4 of Table 1. Experimental results were obtained by Driscoll et al. Ž1992.. Results obtained by Driscoll et al. Ž1992. using the BEM model are also presented in the same figure for comparison. We see from Fig. 12 that our numerical results are very close to the experimental data and are better than those obtained by the BEM model. When water waÕes propagate oÕer a rectangular dike, flow separates at the corner of the dike and Õortices are generated. These phenomena cannot be predicted by a Laplace equation solÕer and this may be the reason that causes the results from BEM to be not as good as the present numerical results. Fig. 13Ža–e. shows the spatial evolution of the first, second, and third harmonic wave amplitudes for cases 5 to 9 in Table 1. In cases 5 to 7, the incident waves are exactly the same as those in cases 1 to 3, except that the dikes are rectangular rather than trapezoidal. The instantaneous velocity fields around the dike at the relative time trT s 0 and 1r2 for cases 5 to 7 are shown in Figs. 14–16. The Ursell numbers of the incident waves in cases 5 to 7 are 1.91, 4.27 and 8.53, and the K–C numbers are 0.37, 0.46, and 0.92, respectively. Comparing Fig. 13Ža–c. with Figs. 8–10, it may be noted that the surface elevations in front of a rectangular dike behave like a partial standing wave, while the reflected waves in front of a trapezoidal dike seem to be smaller, except in case 1, where due to the relatively small period and small wavelength of the incident wave, the reflection is also apparent. We noted also under the same incident wave conditions and the same water depth ratio, the wave amplitudes of higher harmonics on the lee side of a rectangular dike are larger than those for a
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
139
Fig. 12. Wave profiles at several stations and times for monochromatic incident waves with Hin s 2.5 cm, T s1.7 s, q s 0.24 passing over a rectangular dike Žcase 4.. Žv . Experimental data ŽDriscoll et al., 1992.; Ž . Results of present numerical model; Ž- - -. BEM ŽDriscoll et al., 1992..
140
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 13. Spatial evolution of the first, second, and third harmonic wave amplitudes for waves passing over a rectangular dike as specified in cases 5 to 9 of Table 1. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
141
Fig. 14. Surface elevation and velocity field for case 5; Hin s 2.9 cm, T s1.25 s.
trapezoidal dike. From Figs. 14–16, we found that flow separation occurs at the up-wave edge Žleft corner. and down-wave edge Žright corner. of the dike. The separations are followed by a reattachment due to the inertial effect of the flow. In Fig. 15Ža., a vortex formed on the weather side of the dike, while in Fig. 15Žb. a vortex
Fig. 15. Surface elevation and velocity field for case 6; Hin s 2.0 cm, T s 2.0 s.
142
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 16. Surface elevation and velocity field for case 7; Hin s 4.0 cm, T s 2.0 s.
formed on the dike’s lee side. As the K–C number increases from 0.46 in case 6 to 0.92 in case 7, the complete vortices are shed from the edges of the dike, as can be seen from Fig. 16. We see from these figures that the vortex becomes larger as the K–C number increases. When the height of the dike decreases to q s 0.5 Žsee case 8 of Table 1., while the values of other parameters are the same as in case 5, the spatial evolution of wave amplitudes is shown in Fig. 13Žd. and the velocity fields around the dike are shown in Fig. 17. We see from these two figures that at the trailing side of the dike, the higher harmonics almost disappear and the region of flow separation is smaller than that in Fig. 14. On the other hand, if the length of the dike is increased from 2.0 m to 3.8 m Žsee case 9 of Table 1., the corresponding results are shown in Fig. 13Že. and Fig. 18. The behavior of the spatial evolution of wave amplitudes in Fig. 13Že. is similar to that of Fig. 13Ža., except that in Fig. 13Že. the ‘recurrence length’ or the ‘beat length’ of the second and third harmonics appears. The beat length is due to the interplay between the free and bound higher harmonics and can be determined as follows ŽDingemans, 1997.: Assuming weak nonlinearity of the wave fields, the beat length of the second harmonic is determined by 2prŽ k 2 y 2 k 1 ., where k 2 and 2 k 1 are the wavenumbers of the free and bound waves in the second harmonic and both k 1 and k 2 satisfy the linear dispersion relation. The beat length of the second harmonic is 3.3 m from the figure and 3.8 m from calculations. This difference arises because wave nonlinearity is so strong over the dike that the assumption of weak nonlinearity in the theory is no longer adequate. A similar difference has also been reported by Ohyama and Nadaoka Ž1994.. The separation zone in Fig. 18 is about the same as in Fig. 14. As mentioned in the beginning of Section 4, Ting and Kim defined the Reynolds number for this flow to be gl 3 ry and expect that the separation zone would be larger at higher Reynolds
(
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
143
Fig. 17. Surface elevation and velocity field for case 8; Hin s 2.9 cm, T s1.25 s.
numbers. Accordingly, the Reynolds number in case 9 is about 2.6 times larger than in case 5. However, the separation zones in Figs. 14 and 18 are about the same. This indicates that the separation zone seems not to be influenced by the Reynolds number defined by Ting and Kim. The flow separation is the result of the interaction between the inertial force, the viscous force, and the pressure force. The size of separation zone does depend on the Reynolds number and other parameters which represents the effect of the pressure gradient, such as the water depth ratio, q. To reflect the inertia effect of the incident wave, we suggest the Reynolds number for this problem be defined as R e s Ž Hi L2 .rŽy h 0 T ., where h 0 is the stillwater depth, Hi , L, and T are the incident wave height, wavelength, and period, respectively. This Reynolds number is not the same as that defined in Eq. Ž4. for non-dimensionalizing the Navier–Stokes equations. The reason why the Reynolds number is so defined is explained in Appendix A. The Reynolds numbers specified in Table 1 are calculated from this definition. From Table 1 we noted that the Reynolds numbers in cases 5 and 9 are identical. The present numerical results confirm that the magnitude of flow separation does depend on the Reynolds number and the water depth ratio. In order to show that the vortex grows with the K–C number, we performed another numerical experiment with relatively larger values of the K–C number, the Ursell number, and the Reynolds number. The values of these parameters are 1.36, 21.24, and 2.0 = 10 6 , and detailed information is specified in case 10 of Table 1. The water surface elevation and the velocity fields at different times are shown in Fig. 19. The velocity fields shown in Fig. 19 are very similar to the Test 2 experimental results in Ting and Kim. The Ursell number and the K–C number in Ting and Kim’s case are respectively 23.45 and 1.5, while in our case they are 21.24 and 1.36. The variation of the flow fields
144 C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 18. Surface elevation and velocity field for case 9; Hin s 2.9 cm, T s1.25 s.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 19. Surface elevations and velocity fields for case 10; Hin s 5.3 cm, T s 4.0 s.
145
146
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
with respect to different wave phases has been discussed in detail by Ting and Kim. For completeness we will also briefly discuss them here. We have noted from Fig. 19 that as a wave crest propagates over the leading edge of the dike Ž trT s 0., flow separation with reattachment occurs at the left corner of the dike due to the large horizontal velocity near the crest, and it develops into a separation bubble Ž trT s 1r8.. The effects of the separation seem not to reach the free surface. As the wave crest propagates downstream, the separation bubble diminishes in the decelerated flow and eventually disappears at trT s 2r8. Meanwhile, the horizontal velocities of the flow particle on the lee side of the dike accelerate. Due to the strong velocity gradient, a clockwise vortex is generated in this region. The clockwise vortex continues to grow until the horizontal velocity changes the direction. As the wave trough propagates over the dike the flow direction is reversed Ž trT s 4r8.. The flow separates at the down-wave edge of the dike and develops into a separation bubble at trT s 5r8. A counter-clockwise vortex is generated on the weather side of the dike at trT s 6r8. As the wave trough propagates downwards, at trT s 7r8, the separation bubble at the down-wave edge diminishes and flow separates at the up-wave edge of the dike, which develops into flow separation with reattachment as the wave crest passes over this point. 6. Conclusions A numerical scheme has been developed to solve the Navier–Stokes equations and the exact free surface boundary conditions in order to simulate the wave deformation and flow fields in water waves propagating over a submerged dike. The numerical results have been compared with experimental data to confirm the validity of the numerical model. Effects of different parameters on wave deformation, flow separation at the corners of the dike, and vortex generation have been studied systematically. These parameters include the water depth ratio, the length aspect ratio of the dike, the Ursell number, the K–C number, the Reynolds number, and the type of dike. Based on the present results, we may conclude the following. Ž1. The present numerical scheme is able to accurately predict wave deformation and the detailed flow fields, including flow separation and vortex generation, as waves propagate over a dike. This will provide useful information for real applications. Ž2. The first harmonic wave amplitude on the lee side of the dike becomes smaller as the Ursell number of the incident wave increases and the water depth ratio of the dike decreases. However, the shape of the dike plays also an important role. Due to the gentler depth change, the amplitudes of higher harmonics on the lee side of a trapezoidal dike are smaller than those for a rectangular dike. On the other hand, the length aspect ratio of the dike seems not to affect the result very much, except that as the dike becomes longer, a ‘recurrence length’ of higher harmonics becomes apparent. Ž3. Flow separation occurs at the up-wave edge Žleft corner. of a rectangular dike as a wave crest passes over. Conversely, it also occurs at the down-wave edge Žright corner. of the dike as the wave trough propagates over this point. The region of flow separation expands as the K–C number or the Reynolds number of the incident wave increases. In the extreme cases of our study, the separation zone develops into a separation bubble. Even so, the effect of the flow separation does not reach the free surface.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
147
Ž4. An increase in the dike length seems not to affect the magnitude of the flow separation. For this reason, we suggest that the Reynolds number for a propagating wave field be defined as Ž Hi L2 .rŽy h 0 T .. Ž5. The strength of the vorticity increases as the K–C number increases. It should be noted that formation of vorticity was also influenced by the water depth ratio and the shape of the dike. If the height of the dike is small, then the vorticity is weak. Under the same K–C number and the same water depth ratio, voticity was generated around a rectangular dike, but this did not appear in the case of a trapezoidal dike with gentle slopes. Ž6. The generation of vorticity with respect to the variation of free surface has also been explained in this paper. Acknowledgements This research was supported by the National Science Council, Taiwan, R.O.C. under Contract No. NSC 86-2611-E-006-017. Appendix A. Definition of the Reynolds number for viscous wave fields When the viscosity of the fluid is taken into account in wave fields, we often need to use the Reynolds number, which is defined as the ratio of the inertial force to the viscous force acting on the fluid particles. From the order of magnitude analysis, we can estimate that the inertial force term uŽEurE x . is on the order of u˜ 2rL, where u˜ is the maximum horizontal velocity of the fluid and L is the wavelength. Similarly, the order of magnitude of the viscous force term y ŽE 2 urE x 2 . is y urL ˜ 2 . In the above estimation, we assume the wave is long wave, hence the dominant velocity gradient is in the horizontal direction. The Reynolds number is then uL ˜ Re s Ž A1. y The horizontal velocity of fluid particles under the progressive wave is: Hi
Re s
v
cosh k Ž z q h .
cos Ž kx y v t . Ž A2. 2 sinh kh The maximum value of u is: Hi v Hi c u˜ s ; Ž A3. 2 kh h In the above equation, we have assumed that the wave is a shallow water wave Žsinh kh f kh.. Substituting Eq. ŽA3. into Eq. ŽA1. we have: us
LHi c
Hi L2
s Ž A4. yh y hT This is the Reynolds number that we suggested for use with viscous wave fields.
148
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
References Battjes, J.A., Beji, S., 1991. Spectral evolution in waves traveling over a shoal. In: Peregrine, D.H. ŽEd.., Proc. of the Nonlinear Water Waves Workshop. Univ. of Bristol, Bristol, pp. 11–19. Beji, S., Battjes, J.A., 1993. Experimental investigation of the wave propagation over a bar. Coastal Eng. 19, 151–162. Beji, S., Battjes, J.A., 1994. Numerical simulation of nonlinear wave propagation over a bar. Coastal Eng. 23, 1–16. Chan, R.K.C., Street, R.L., 1970. A computer study of finite-amplitude water waves. J. Comput. Phys. 6, 68–94. Chen, C.J., Chen, H.C., 1982. The finite-analytic method. IIHR Report 232-IV, Iowa Institute of Hydraulic Research, The Univ. of Iowa. Chen, H.C., Patel, V.C., 1987. Laminar flow at the trailing edge of a flat plate. AIAA J. 25, 920–928. Dean, R.G., Dalrymple, R.A., 1984. Water wave mechanics for engineers and scientists. Prentice-Hall, Englewood Cliffs, NJ. Dingemans, M.W., 1997. Water wave propagation over uneven bottoms: Part II. World Scientific Publishing, Singapore. Driscoll, A.M., Dalrymple, R.A., Grill, S.T., 1992. Harmonic generation and transmission past a submerged rectangular obstacle. Proc. 23rd Int. Coastal Eng. Conf., Venice, ASCE, pp. 1142–1152. Grue, J., 1992. Nonlinear water waves at a submerged obstacle or bottom topography. J. Fluid Mech. 244, 455–476. Huang, C.J., Zhang, E.C., Lee, J.F., 1998. Numerical simulation of nonlinear viscous wavefields generated by a piston-type wavemaker. J. Eng. Mech. 124, 1110–1120. Johnson, J.W., Fuchs, R.A., Morison, J.R., 1951. The damping action of submerged breakwaters. Trans. Am. Geophys. Union 32, 704–718. Jolas, P., 1960. Passage de la houle sur un seuil. Houille Blanche 2, 148–152. Kittitanasuan, W., Goda, Y., Shiobara, T., 1993. Deformation of nonlinear waves on a rectangular step. Coastal Eng. in Japan 36, 133–153. Madsen, O.S., 1971. On the generation of long waves. J. Geophys. Res. 76, 8672–8683. Madsen, P.A., Murray, R., Sorensen, O.R., 1991. A new form of the Boussinesq equations with improved linear dispersion characteristics, Part 1. Coastal Eng. 15, 371–388. Madsen, P.A., Sorensen, O.R., 1992. A new form of the Boussinesq equations with improved linear dispersion characteristics, Part 2. Coastal Eng. 18, 183–204. Ohyama, T., Nadaoka, K., 1991. Development of a numerical wave tank for analysis of nonlinear and irregular wave field. Fluid Dynamics Research 8, 231–251. Ohyama, T., Nadaoka, K., 1992. Modeling the transformation of nonlinear waves passing over a submerged dike. Proc. 23rd Int. Coastal Eng. Conf., Venice, ASCE, pp. 526–539. Ohyama, T., Nadaoka, K., 1994. Transformation of a nonlinear wave train passing over a submerged shelf without breaking. Coastal Eng. 24, 1–22. Ohyama, T., Beji, S., Nadaoka, K., Battjes, J.A., 1994. Experimental verification of numerical model for nonlinear wave evolution. J. Waterway, Port, Coastal, and Ocean Eng. 120, 637–644. Peregrine, D.H., 1967. Long waves on a beach. J. Fluid Mech. 27, 815–827. Patankar, S.V., 1980. Numerical heat transfer and fluid flow. McGraw-Hill, New York. Ting, F.C.K., Kim, Y.K., 1994. Vortex generation in water waves propagating over a submerged obstacle. Coastal Eng. 24, 23–49.
Wave deformation and vortex generation in water waves propagating over a submerged dike Ching-Jer Huang ) , Chih-Ming Dong Department of Hydraulics and Ocean Engineering, National Cheng Kung UniÕersity, 1 Dah Shyue Road, Tainan 70101, Taiwan Received 12 August 1998; received in revised form 25 January 1999; accepted 16 February 1999
Abstract The Navier–Stokes equations and the exact free surface boundary conditions are solved to simulate wave deformation and vortex generation in water waves propagating over a submerged dike. Incident waves are generated by a piston-type wavemaker set up in the computational domain. Numerical results are compared with experimental data in order to confirm the validity of the numerical model. The fast Fourier transform and a wave resolution technique are applied to decompose the transformed waves and the higher harmonics. Effects of different parameters on wave transformation and vortex generation are studied systematically. These parameters include the Ursell number, the Keulegan–Carpenter number, the water depth ratio, the Reynolds number, the length aspect ratio of the dike, and the type of dike. q 1999 Elsevier Science B.V. All rights reserved. Keywords: Submerged dike; Wave deformation; Higher harmonics; Flow separation; Vortex generation
1. Introduction There are two main physical phenomena associated with the propagation of water waves over a submerged dike. The first is the generation of higher harmonics and the second is flow separation and vortex generation on the weather side and on the lee side of the obstacle. The former phenomenon has long been known from direct field observations and laboratory measurements Že.g., Johnson et al., 1951; Jolas, 1960.. Earlier researchers have tried to predict the generation of higher harmonics using )
Corresponding author. Tel.: q886-6-2757575 ext. 63252; Fax: q886-6-2741463; E-mail: [email protected] 0378-3839r99r$ - see front matter q 1999 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 8 - 3 8 3 9 Ž 9 9 . 0 0 0 1 7 - 4
124
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
nonlinear shallow-water wave theories, such as the Boussinesq and the KdV equations. Since the derivation of the Boussinesq equations is based on the assumptions of both weak nonlinearity and weak dispersivity of waves, these equations may not be valid for the prediction on the lee side of the dike, where higher harmonics may arise in the form of deep-water waves. In order to overcome this defect, improvements on the Boussinesq equations have been developed. Peregrine Ž1967. developed equations of motion for long waves in water of varying depth. These equations can be viewed as the Boussinesq equations for water of varying depth. Madsen et al. Ž1991. and Madsen and Sorensen Ž1992. improved the Boussinesq equation by adding a term to improve its dispersion characteristics. Battjes and Beji Ž1991. derived a new set of Boussinesq equations by combining the method used by Madsen et al. Ž1991. and the equation proposed by Peregrine Ž1967.. Beji and Battjes Ž1994. also employed these equations to study the deformation of water waves passing over a dike. The numerical results of wave height were in good agreement with experimental data. The boundary element method has also been applied to solve the Laplace equation and the nonlinear free surface boundary conditions in order to study the decomposition phenomena of waves passing over a rectangular step ŽKittitanasuan et al., 1993. and a rectangular submerged dike ŽOhyama and Nadaoka, 1992, 1994.. In Ohyama and Nadaoka, the incident waves were generated by a numerical wave tank model developed by Ohyama and Nadaoka Ž1991.. This is similar to the present method of setting up a piston-type wavemaker in the computational domain to generate the incident waves. The viscosity effect has been neglected in the aforementioned theories. Most of the studies concentrated on the transformation of waves passing over submerged dikes. As far as the wave transformation is concerned, the results are satisfactory. However, the prediction of the flow fields without taking the viscosity effect into account may in some situations fail to provide useful knowledge about the real flow fields. For example, the vortices exiting in front of and behind dikes in the nearshore region may have an important effect on the safety and stability of the dikes, but they cannot be predicted by earlier theories. Recently, Ting and Kim Ž1994. conducted laboratory experiments to investigate flow separation and vortex generation induced by waves propagating over a submerged rectangular dike. Velocity measurements were performed using a two-component laser-doppler anemometer ŽLDA.. Significant results were obtained, which showed that the Keulegan–Carpenter ŽK–C. number is the most important parameter in determining the formation and development of vortices. Based on the dimensional analysis, it was also expected that the zone of flow separation at the corners of the dikes would be larger at a higher Reynolds number. Ting and Kim’s measurements showed a clear picture of the flow fields around a submerged rectangular dike. Due to the complex nature of the flow around the dike, Ting and Kim suggested that it would be very difficult to determine flow separation effects unambiguously without solving the viscous flow equations in the near field. The Navier–Stokes equations and the exact free surface boundary conditions must be solved in order to predict the realistic flow fields around the submerged dikes. In the present study, the numerical wave tank model developed by Huang et al. Ž1998. was used to simulate wave deformation and the flow fields as waves propagated over submerged trapezoidal and rectangular dikes. In this model the unsteady two-dimen-
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
125
sional Navier–Stokes equations were discretized by means of a finite-analytical scheme developed by Chen and Chen Ž1982.. The SUMMAC method ŽChan and Street, 1970. was used to determine the location of the free surface. In order to verify the accuracy of the numerical calculations, numerical results were compared with the experimental data of Beji and Battjes Ž1994. for a trapezoidal dike and with those of Driscoll et al. Ž1992. for a rectangular dike. The fast Fourier transform ŽFFT. and Grue’s method ŽGrue, 1992. were also applied to decompose the transformed waves and the higher harmonics. Effects of different parameters on wave transformation, flow separation, and vortex generation were studied and compared systematically. These parameters include the Ursell number, the K–C number, the water depth ratio of the dike, the Reynolds number, the length aspect ratio of the dike, and the type of dike.
2. Governing equations and boundary conditions The wave transformation and flow fields in water waves passing over a submerged dike are studied in the present paper. A schematic diagram of a trapezoidal dike is shown in Fig. 1. For a rectangular dike, only the rectangular part Žindicated by dotted lines. of the trapezoid is kept. A piston-type wavemaker with stroke S 0 is located at x s 0 and generates the incident waves. The still water depth is h 0 . The trapezoidal dike consists of an upslope of S 1 and a downslope of S 2 . The shallow water depth above the dike has a depth of qh 0 Ž0 - q - 1.. For an incompressible, viscous fluid, the dimensionless continuity equation in the Cartesian coordinate system is Eu
EÕ q
Ex
Ey
s0
Ž 1.
and the Navier–Stokes equations are given by Eu Et
Eu qu
Ex
Eu qÕ
Ep sy
1 q
Ey
Ex
Re
EÕ
Ep
1
ž
E2 u Ex2
E2 u q
E y2
/
Ž 2.
/
Ž 3.
and EÕ Et
EÕ qu
Ex
qÕ
sy Ey
q Ey
Re
ž
E2 Õ Ex2
E2 Õ q
E y2
where u and Õ are the horizontal and vertical velocity components, t is the time, p is the hydrodynamic pressure, and R e is the Reynolds number, defined as Re s
u0 h0
y
Ž 4.
where y is the kinematic viscosity of the fluid, and u 0 is the velocity–amplitude of the wavemaker, which moves back and forth with a velocity of u s u 0 sin v t. In the present case, u 0 and h 0 are used to non-dimensionalize the velocity and length, while t 0 s h 0ru 0 is chosen to non-dimensionalize the time.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
126
Fig. 1. Schematic diagram of submerged trapezoidal and rectangular dikes.
The boundary conditions which must be satisfied to solve the problem shown in Fig. 1 are: Ž1. the kinematic and the dynamic free surface boundary conditions at the free surface, Ž2. the no-slip boundary condition at the channel bottom, Ž3. the upstream boundary condition on the wavemaker plate, and Ž4. the radiation condition on the outgoing boundary. The initial conditions of the velocities, hydrodynamic pressure, and surface displacements are set at zero at time t s 0. Conditions Ž1., Ž3., and Ž4. mentioned above are explained in more detail below. 2.1. Free surface boundary conditions The kinematic condition states that fluid particles at a free surface remain on the free surface, and can be expressed as Eh Et
Eh qu
Ex
Ž 5.
sÕ
where h s h Ž x,t . is the location of the free surface. The dynamic condition requires that, along the free surface boundary, the normal stress is equal to the atmospheric pressure and the tangential stress is zero. These conditions are expressed by Huang et al. Ž1998. as the following Eh
2 1q
h p0 s
Fr
2
q
2
ž / ž / Ex
Eh
Re 1 y
EÕ 2
Ey
Ž 6.
Ex
and Eu
EÕ sy
Ey
Ex
EÕ Eh
4 q Eh
ž / Ex
2
Ey Ex y1
Ž 7.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
127
where p 0 s pŽ x,h . is the hydrodynamic pressure at the free surface and Fr is the Froude number, defined as u0 Fr s Ž 8. gh 0
(
In numerical computations, Eq. Ž6. is used to determine the pressure at the free surface and Eq. Ž7. is used to extrapolate the horizontal velocity at the free surface from the flow domain. The vertical velocity component Õ is then calculated from the continuity equation using the known velocity component u, obtained from Eq. Ž7.. 2.2. The upstream boundary condition The no-slip boundary condition at the wavemaker is released at two nodal points beneath the free surface. The slip velocities Õ along the wavemaker are extrapolated from the flow inside the domain. 2.3. The downstream boundary condition The downstream boundary condition requires that, at a large distance from the wavemaker, the wave is outgoing without reflection. According to the wave equation and continuity equation, we can set Ep Et
qc
Ep Ex
Eu s 0,
Et
Eu qc
Ex
Eu s 0,
EÕ q
Ex
Ey
s0
Ž 9.
where c is the phase speed of the wave and is determined from the dispersion relation.
3. Numerical method 3.1. Discretization of goÕerning equations In the present study, the finite-analytic ŽFA. method is applied to the discretization of the unsteady two-dimensional Navier–Stokes equations. In this method, local analytical solutions obtained from the linearized Navier–Stokes equations for the discretized computational elements are incorporated into the numerical method. The value of a variable at a particular node is expressed in terms of the values of neighboring nodes. The simplified four-point version of the finite-analytic method developed by Chen and Patel Ž1987. was used here to discretize the Navier–Stokes equations. Please refer to Chen and Patel for more information about the derivation of the FA coefficients. The Navier–Stokes Eqs. Ž2. and Ž3. can be written in linearized form as E2 f Ex
2
E2 f q Ey
2
s2 A
Ef Ex
q2 B
Ef Ey
q Re
Ef Et
q Sf
Ž 10 .
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
128
where f represents the velocity components u or Õ, and As
R eU 2
Su s R e
,
Ep Ex
Bs ,
R eV 2
SÕ s R e
Ep
Ž 11 .
Ey
where U, V are mean values of u and Õ in a given computational element. Eq. Ž10. can be rewritten as E2 f Ey
2
s2 B
Ef Ey
q Gf
Ž 12 .
where Gf is the nonhomogeneous term Gf s 2 A
Ef Ex
q Re
Ef Et
q Sf y
E2 f Ex2
Ž 13 .
Eq. Ž12. is a nonhomogeneous second-order ordinary differential equation and can be easily solved to obtain the analytic solution, namely f s a Ž e 2 B y y 1 . q by q c
Ž 14 .
By using the values of the variables at a node P and at the two neighboring nodal points N and S in the y-direction ŽFig. 2., the constants a, b, and c can be determined. An equation similar to Eq. Ž12. can be obtained by substituting Eq. Ž14. into Eq. Ž12., but with variable y replaced by x. Following procedures similar to those used to solve the differential equation in Eq. Ž12., one obtains f p s Ž Cs fs q Cn f n q C w f w q Ce fe q Ct f pny 1 . q Cp S f
Ž 15 .
where the coefficient C with different subscripts are the FA coefficients, which relate to the time step, the grid size, and the parameters appearing in Eq. Ž10., such as A, B, and
Fig. 2. Symbols used in the present numerical method.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
129
the Reynolds number. In obtaining Eq. Ž15., the time derivative is approximated by an implicit, backward difference scheme. From Eq. Ž15. we see that the value of f at node P is related to the value of f at the neighboring four nodes, to the value of f in the previous time step, and to the pressure gradient. The variable f can be determined by applying Eq. Ž15. to all the grid nodes in the flow domain and using the given initial and boundary conditions. The SIMPLER algorithm developed by Patankar Ž1980. was used to calculate the coupled velocity and pressure fields. A staggered numerical grid was used. The velocity components u and Õ are defined as located at the boundaries of a control volume, while the pressure p is located at the center. 3.2. Treatment of the boundary The coordinate transformation technique was not applied in the computation. In order to make a rectangular grid system for a trapezoidal dike, the smoothly varying bed topography is divided into a series of small rectangular steps with a horizontalrvertical ratio of 20:1 or 10:1 ŽFig. 3.. The Marker and Cell method ŽMAC. and its modified version SUMMAC are used in combination to calculate the free surface boundary. The major concept underlying the MAC method is the use of marker particles to identify the location of the free surface. By tracking the positions of these marker particles, the transient location of the free surface can be determined. Pressure at the free surface is calculated by the normal dynamic free surface boundary condition ŽEq. Ž6.., using the known velocities u and Õ. The velocity component u at a free surface cell is extrapolated from the velocity of the main body of the fluid by means of the tangential dynamic boundary condition ŽEq. Ž7... The value of Õ at a free surface cell can be extrapolated from the fluid domain using an extrapolation formula. However, the velocity components obtained in this way may result in a violation of the continuity equation at local cells. To avoid this risk, we use the continuity equation to determine the velocity component Õ. A method which satisfies the discretized Poisson equation for the pressure and is called ‘irregular star’ is applied to determine the pressure at the cells near the free surface.
Fig. 3. Numerical grids in the computational domain of a trapezoidal dike.
130
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
In the numerical model, the Eulerian grid system is used, which means the form and size of the computational grids are fixed. However, the computational domain near the wavemaker changes as the wavemaker moves back and forth. In the computation the grids in the neighborhood of the wavemaker are arranged such that the boundary of the computational domain matches the location of the wavemaker. The numerical results shown in the following were calculated with dimensionless D x and D y Žnormalized with respect to h 0 . equal to 0.1, except for the region near the free surface where the finest grid size D y was 0.025. One example of the grid system is shown in Fig. 3. The increment of dimensionless time, normalized with respect to t 0 , was chosen to be 0.001. 4. Water waves passing over a trapezoidal dike The numerical methods described in Section 3 are used to simulate the deformation of waves propagating over submerged trapezoidal and rectangular dikes. The incident waves are generated by a piston-type wavemaker imposed in the computation domain. To generate the small-amplitude waves, the strokes of the numerical wavemaker are determined by linear wavemaker theory ŽDean and Dalrymple, 1984.. To generate relatively long waves of finite amplitude and to make sure that the generated waves are of permanent form, Madsen’s second-order wavemaker theory ŽMadsen, 1971. is applied to determine the wavemaker motion. The accuracy of the numerical wave tank model applied in the present study has been verified and discussed in detail in Huang et al. Ž1998.. In an experimental study of water waves propagating over a rectangular submerged dike, Ting and Kim Ž1994. indicated by dimensional analysis that the normalized water surface elevation is a function of water depth ratio, q, the length aspect ratio of the dike, lrh 0 , the Ursell number, Hi L2rh 30 , the K–C number, Hi Lrh20 , and the Reynolds number, which was defined to be gl 3 ry , in which h 0 is the stillwater depth, and Hi and L are the incident wave height and wavelength. The Ursell number represents the nonlinearity of the incoming waves. The K–C number is interpreted as the ratio of maximum displacement of a fluid particle from its neutral position to the body length, and is an important parameter in determining the formation and development of vortices in the flow. In accordance with the above ideas, we have chosen three different incident waves passing over a trapezoidal and a rectangular dike in order to study the main characteristics of waves passing over submerged dikes. Detailed information about the incident waves and the dikes are summarized in Table 1 for case 1 through case 10. Cases 1 to 3 in Table 1 refer to a trapezoidal dike and are discussed here, while the other cases are for a rectangular dike and are left to Section 5. It should be noted that we have confined ourselves to nonbreaking waves in all the numerical simulations. In the following sections, we simulate the transformation of wave passing over a trapezoidal dike and compare the results with the experimental data of Beji and Battjes Ž1993, 1994.. The FFT and Grue’s method are then applied to decompose the transformed waves and the higher harmonics. Effect of the Ursell number on the wave transformation is also discussed.
(
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
131
Table 1 Numerical conditions Variable Case 1
Case 2
Case 3
Case 4
Case 5
Case 6
Case 7
Case 8
Case 9
Case 10
h 0 Žm. Hi Žcm. L Žm. T Žs. Ur K–C Re
0.4 2.0 3.695 2.0 4.27 0.46 3.4=10 5
0.4 4.0 3.695 2.0 8.53 0.92 6.8=10 5
0.5 2.5 3.327 1.7 2.21 0.33 3.3=10 5
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
0.4 2.0 3.695 2.0 4.27 0.46 3.4=10 5
0.4 4.0 3.695 2.0 8.53 0.92 6.8=10 5
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
0.6096 5.3 9.529 4.0 21.24 1.36 2.0=10 6
0.25 2.0 1r20 1r10
0.25 2.0 1r20 1r10
0.24 0.79
0.25 2.0
0.25 2.0
0.25 2.0
0.5 2.0
0.25 3.8
0.4 2.9 2.052 1.25 1.91 0.37 2.4=10 5
Dike type A q l Žm. S1 S2
0.25 2.0 1r20 1r10
B 0.5 0.6096
A: trapezoidal; B: rectangular. Ur s Hi L2 r h 30 ; K–C s Hi Lr h20 ; R e s Hi L2 ry h 0 T.
4.1. Comparison with experiments In the experiment of Beji and Battjes Ž1993, 1994. the stillwater depth was set at 0.4 m. A piston-type wavemaker was installed at one end of the flume. The model of the submerged trapezoidal dike consisted of an upslope of 1:20 and a downslope of 1:10. The water depth was 0.1 m in the shallowest region above the horizontal part Ži.e., q s 0.25.. The wave profiles were measured at 7 different locations as shown in Fig. 4. The incident wave height Hi and period T were 2.0 cm and 2.0 s, respectively. The incident wavelength L was 3.695 m based on linear theory. The Ursell number Hi L2rh 30 was 4.27. Fig. 5 shows the numerical and the experimental results of the wave profiles
Fig. 4. Definition sketch of wave tank and locations of wave gages Žfrom Beji and Battjes, 1994..
132
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 5. Wave profiles at several stations for monochromatic incident waves with Hin s 2.0 cm, T s 2.0 s, . Results of present numerical model; q s 0.25. Žv . Experimental data ŽBeji and Battjes, 1994.; Ž Ž- - -. BEM ŽOhyama et al., 1994..
obtained at stations 1 to 6. The numerical results obtained by Ohyama et al. Ž1994. using the BEM model are also presented in Fig. 5 for comparison. From Fig. 5 we see that both our numerical results and those obtained by Ohyama et al. are in good agreement with the experimental data. The reason that our results are not as good as those of Ohyama may be due to the simplified treatment of the boundary. The wave profiles shown in stations 1 to 6 indicated that waves undergo rapid variation as they propagated over the dike. The main reason for this is the generation of higher harmonics and energy transfer from the bound waves to the free waves of the higher harmonics. This will be discussed in more detail in Section 4.2. 4.2. Basic mechanism of waÕe decomposition The numerical results for free surface elevations along the wave channel at different times are shown in Fig. 6. It can be seen from Fig. 6 that as waves propagate over the dike, the primary wave crests become steeper and a dispersive tail gradually develops. A small wave, indicated by an arrow in Fig. 6, appears at the trailing edge of a primary
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
133
Fig. 6. Surface elevations along the wave channel at several times for waves with Hin s 2.0 cm, T s 2.0 s, and q s 0.25 passing over a submerged trapezoidal dike.
wave. Because the small wave propagates more slowly than the main crest, it gradually detaches from the main crest and is overtaken by the next wave. To understand the basic mechanism of wave evolution and decomposition, the method of Grue Ž1992. is used to determine the amplitude of high frequency components. This method is based on the resolution technique of two-wave gage systems and is described briefly as follows.
134
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
We assume that on the lee side of the dike the waves are given by Ž1. Ž1. hq Ž x , t . s aq cos Ž Kx y v t q d q . q Ý aŽlqn. cos n Ž Kx y v t q d qŽ 1 . .
n)1 Ž n. q Ý aq cos Ž K n x y n v t q dqŽ n. .
Ž 16 .
n)1 n. where K represents the wavenumber of the first-harmonic wave, aŽlq , n s 2,3, . . . , Ž n. denote the amplitudes of the nth-harmonic bound waves, aq and K n with n s 2, 3, . . . denote the amplitudes and the wavenumber of the transmitted nth-harmonic free waves,
Fig. 7. Comparison of wave number–frequency spectra for monochromatic incident waves with Hin s 2.0 cm, T s 2.0 s, q s 0.25.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
135
and dqŽ n., n s 1, 2, . . . , are the phase angles. From Eq. Ž16. we see that the transformed waves on the lee side of the dike can be decomposed into first harmonic waves and nth-harmonic bound and free waves. The bound waves have the same phase velocity as the first harmonic wave and are bound to it. Hence, they cannot be distinguished simply from the wave profiles. The free waves of higher harmonics satisfy the dispersion relation and propagate at different velocity from the first harmonic and will be detached from it as the wave propagates. If we define the Fourier transform of the surface elevation h Ž x,t . as
hˆ Ž n. Ž x . s
v
2p r v
H 2p 0
h Ž x , t . exp Ž yin v t . d t ,
n s 1, 2, . . . ,
Ž 17 .
then the amplitudes of the free and bound waves can be determined as follows by giving hˆ Ž n. at two different locations, namely x 1 and x 1 q D x: Ž n. aq s
1 sin Ž 1r2 Ž K n y nK . D x .
hˆ Ž n . Ž x 1 . y hˆ Ž n. Ž x 1 q D x . exp Ž inKD x . ,
n s 2, 3, . . . n. aŽlq s
1 sin Ž 1r2 Ž K n y nK . D x .
n s 2, 3, . . .
hˆ Ž n . Ž x 1 . y hˆ Ž n. Ž x 1 q D x . exp Ž iK n D x . ,
Ž 18 .
The wavenumber–frequency spectra of waves above the submerged dike and on the lee side of the dike are shown in Fig. 7Ža. and Žb., respectively. The results shown in Fig. 7Ža. are obtained using the surface elevations simulated at stations 3 and 4, while those shown in Fig. 7Žb. use the data of stations 7 and 8, which is at a distance of 0.4 m shoreward from station 7. In Fig. 7, a n and a1 are the nth-harmonic and first harmonic
Fig. 8. Spatial evolution of harmonic amplitudes for a monochromatic incident wave with Hin s 2.9 cm, T s1.25 s, Ur s1.84, q s 0.25. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
136
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 9. Spatial evolution of harmonic amplitudes for a monochromatic incident wave with Hin s 2.0 cm, T s 2.0 s, Ur s 4.27, q s 0.25. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
wave amplitudes, f 0 is the frequency of the incident wave, K is the wave number, and h is the local water depth. In each frf 0 y Kh plane, a curved line ŽF. represents the dispersion relation, and a straight line ŽB. is drawn from the origin to a point at which the ratio of wave amplitude for the first harmonic Ž frf 0 s 1. is maximum. The wave components lying on line ŽF. are free waves, and we can determine their phase velocities from the dispersion relation Ž4p 2 f 2 s gk tanh kh.. The phase velocity of the wave components on line ŽB. is the same as that of the incident wave. The free waÕes haÕe shorter periods than incident waÕes; thus, from the dispersion relation, their phase Õelocities are slower than those of incident waÕes. Hence, the free wave is detached from the main crest and is overtaken by the next wave as shown in Fig. 6. From Fig. 7 we noted also that above the submerged dike, the wave amplitude of the free wave and bound wave are about the same. However, on the lee side of the dike, the amplitude of the bound wave almost disappears. The main reason for this is that as the
Fig. 10. Spatial evolution of harmonic amplitudes for a monochromatic incident wave with Hin s 4.0 cm, T s 2.0 s, Ur s8.53, q s 0.25. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 11. Surface elevation and velocity field for case 1; Hin s 2.9 cm, T s1.25 s.
137
138
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
waves propagate into the deep water region where wave nonlinearity is so weak that the bound wave is reduced substantially, a large amount of energy is transferred from the bound wave to the free wave. 4.3. Spatial amplitude modulations of waÕes and the flow field Figs. 8–10 show the spatial evolution of the first, second and third harmonic wave amplitudes Ž n s 1, 2, 3. for cases 1, 2, and 3. The Ursell numbers in cases 1 through 3 are 1.91, 4.27, and 8.53, and the K–C numbers are 0.37, 0.46, and 0.92, respectively. We noted from Fig. 8 through 10 that on the lee side of the dike the higher harmonic wave amplitudes become larger as the Ursell number increases. In case 3, the second harmonic wave amplitude is even 55% larger than that of the first harmonic wave. The instantaneous velocity fields for case 1 are shown in Fig. 11Ža,b. for relative time trT s 0 and 1r2. In this figure, there is no flow separation or vortex generation around the dike. This may be due to the gentle slopes located in front of and behind the dike. The flow fields for cases 2 and 3 are very similar to that of case 1, and no flow separation or vortex generation occurs. 5. Waves passing over rectangular dikes Before we proceed to study the effect of different parameters on the main characteristics of water waves passing over a rectangular dike, the accuracy of the numerical results are verified. Fig. 12 shows a comparison of numerical and experimental results at several stations and times for monochromatic incident waves passing over a rectangular dike. Information about the incident wave and the dike are specified in case 4 of Table 1. Experimental results were obtained by Driscoll et al. Ž1992.. Results obtained by Driscoll et al. Ž1992. using the BEM model are also presented in the same figure for comparison. We see from Fig. 12 that our numerical results are very close to the experimental data and are better than those obtained by the BEM model. When water waÕes propagate oÕer a rectangular dike, flow separates at the corner of the dike and Õortices are generated. These phenomena cannot be predicted by a Laplace equation solÕer and this may be the reason that causes the results from BEM to be not as good as the present numerical results. Fig. 13Ža–e. shows the spatial evolution of the first, second, and third harmonic wave amplitudes for cases 5 to 9 in Table 1. In cases 5 to 7, the incident waves are exactly the same as those in cases 1 to 3, except that the dikes are rectangular rather than trapezoidal. The instantaneous velocity fields around the dike at the relative time trT s 0 and 1r2 for cases 5 to 7 are shown in Figs. 14–16. The Ursell numbers of the incident waves in cases 5 to 7 are 1.91, 4.27 and 8.53, and the K–C numbers are 0.37, 0.46, and 0.92, respectively. Comparing Fig. 13Ža–c. with Figs. 8–10, it may be noted that the surface elevations in front of a rectangular dike behave like a partial standing wave, while the reflected waves in front of a trapezoidal dike seem to be smaller, except in case 1, where due to the relatively small period and small wavelength of the incident wave, the reflection is also apparent. We noted also under the same incident wave conditions and the same water depth ratio, the wave amplitudes of higher harmonics on the lee side of a rectangular dike are larger than those for a
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
139
Fig. 12. Wave profiles at several stations and times for monochromatic incident waves with Hin s 2.5 cm, T s1.7 s, q s 0.24 passing over a rectangular dike Žcase 4.. Žv . Experimental data ŽDriscoll et al., 1992.; Ž . Results of present numerical model; Ž- - -. BEM ŽDriscoll et al., 1992..
140
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 13. Spatial evolution of the first, second, and third harmonic wave amplitudes for waves passing over a rectangular dike as specified in cases 5 to 9 of Table 1. Žv . ns1; Ž'. ns 2; ŽB. ns 3.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
141
Fig. 14. Surface elevation and velocity field for case 5; Hin s 2.9 cm, T s1.25 s.
trapezoidal dike. From Figs. 14–16, we found that flow separation occurs at the up-wave edge Žleft corner. and down-wave edge Žright corner. of the dike. The separations are followed by a reattachment due to the inertial effect of the flow. In Fig. 15Ža., a vortex formed on the weather side of the dike, while in Fig. 15Žb. a vortex
Fig. 15. Surface elevation and velocity field for case 6; Hin s 2.0 cm, T s 2.0 s.
142
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 16. Surface elevation and velocity field for case 7; Hin s 4.0 cm, T s 2.0 s.
formed on the dike’s lee side. As the K–C number increases from 0.46 in case 6 to 0.92 in case 7, the complete vortices are shed from the edges of the dike, as can be seen from Fig. 16. We see from these figures that the vortex becomes larger as the K–C number increases. When the height of the dike decreases to q s 0.5 Žsee case 8 of Table 1., while the values of other parameters are the same as in case 5, the spatial evolution of wave amplitudes is shown in Fig. 13Žd. and the velocity fields around the dike are shown in Fig. 17. We see from these two figures that at the trailing side of the dike, the higher harmonics almost disappear and the region of flow separation is smaller than that in Fig. 14. On the other hand, if the length of the dike is increased from 2.0 m to 3.8 m Žsee case 9 of Table 1., the corresponding results are shown in Fig. 13Že. and Fig. 18. The behavior of the spatial evolution of wave amplitudes in Fig. 13Že. is similar to that of Fig. 13Ža., except that in Fig. 13Že. the ‘recurrence length’ or the ‘beat length’ of the second and third harmonics appears. The beat length is due to the interplay between the free and bound higher harmonics and can be determined as follows ŽDingemans, 1997.: Assuming weak nonlinearity of the wave fields, the beat length of the second harmonic is determined by 2prŽ k 2 y 2 k 1 ., where k 2 and 2 k 1 are the wavenumbers of the free and bound waves in the second harmonic and both k 1 and k 2 satisfy the linear dispersion relation. The beat length of the second harmonic is 3.3 m from the figure and 3.8 m from calculations. This difference arises because wave nonlinearity is so strong over the dike that the assumption of weak nonlinearity in the theory is no longer adequate. A similar difference has also been reported by Ohyama and Nadaoka Ž1994.. The separation zone in Fig. 18 is about the same as in Fig. 14. As mentioned in the beginning of Section 4, Ting and Kim defined the Reynolds number for this flow to be gl 3 ry and expect that the separation zone would be larger at higher Reynolds
(
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
143
Fig. 17. Surface elevation and velocity field for case 8; Hin s 2.9 cm, T s1.25 s.
numbers. Accordingly, the Reynolds number in case 9 is about 2.6 times larger than in case 5. However, the separation zones in Figs. 14 and 18 are about the same. This indicates that the separation zone seems not to be influenced by the Reynolds number defined by Ting and Kim. The flow separation is the result of the interaction between the inertial force, the viscous force, and the pressure force. The size of separation zone does depend on the Reynolds number and other parameters which represents the effect of the pressure gradient, such as the water depth ratio, q. To reflect the inertia effect of the incident wave, we suggest the Reynolds number for this problem be defined as R e s Ž Hi L2 .rŽy h 0 T ., where h 0 is the stillwater depth, Hi , L, and T are the incident wave height, wavelength, and period, respectively. This Reynolds number is not the same as that defined in Eq. Ž4. for non-dimensionalizing the Navier–Stokes equations. The reason why the Reynolds number is so defined is explained in Appendix A. The Reynolds numbers specified in Table 1 are calculated from this definition. From Table 1 we noted that the Reynolds numbers in cases 5 and 9 are identical. The present numerical results confirm that the magnitude of flow separation does depend on the Reynolds number and the water depth ratio. In order to show that the vortex grows with the K–C number, we performed another numerical experiment with relatively larger values of the K–C number, the Ursell number, and the Reynolds number. The values of these parameters are 1.36, 21.24, and 2.0 = 10 6 , and detailed information is specified in case 10 of Table 1. The water surface elevation and the velocity fields at different times are shown in Fig. 19. The velocity fields shown in Fig. 19 are very similar to the Test 2 experimental results in Ting and Kim. The Ursell number and the K–C number in Ting and Kim’s case are respectively 23.45 and 1.5, while in our case they are 21.24 and 1.36. The variation of the flow fields
144 C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 18. Surface elevation and velocity field for case 9; Hin s 2.9 cm, T s1.25 s.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
Fig. 19. Surface elevations and velocity fields for case 10; Hin s 5.3 cm, T s 4.0 s.
145
146
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
with respect to different wave phases has been discussed in detail by Ting and Kim. For completeness we will also briefly discuss them here. We have noted from Fig. 19 that as a wave crest propagates over the leading edge of the dike Ž trT s 0., flow separation with reattachment occurs at the left corner of the dike due to the large horizontal velocity near the crest, and it develops into a separation bubble Ž trT s 1r8.. The effects of the separation seem not to reach the free surface. As the wave crest propagates downstream, the separation bubble diminishes in the decelerated flow and eventually disappears at trT s 2r8. Meanwhile, the horizontal velocities of the flow particle on the lee side of the dike accelerate. Due to the strong velocity gradient, a clockwise vortex is generated in this region. The clockwise vortex continues to grow until the horizontal velocity changes the direction. As the wave trough propagates over the dike the flow direction is reversed Ž trT s 4r8.. The flow separates at the down-wave edge of the dike and develops into a separation bubble at trT s 5r8. A counter-clockwise vortex is generated on the weather side of the dike at trT s 6r8. As the wave trough propagates downwards, at trT s 7r8, the separation bubble at the down-wave edge diminishes and flow separates at the up-wave edge of the dike, which develops into flow separation with reattachment as the wave crest passes over this point. 6. Conclusions A numerical scheme has been developed to solve the Navier–Stokes equations and the exact free surface boundary conditions in order to simulate the wave deformation and flow fields in water waves propagating over a submerged dike. The numerical results have been compared with experimental data to confirm the validity of the numerical model. Effects of different parameters on wave deformation, flow separation at the corners of the dike, and vortex generation have been studied systematically. These parameters include the water depth ratio, the length aspect ratio of the dike, the Ursell number, the K–C number, the Reynolds number, and the type of dike. Based on the present results, we may conclude the following. Ž1. The present numerical scheme is able to accurately predict wave deformation and the detailed flow fields, including flow separation and vortex generation, as waves propagate over a dike. This will provide useful information for real applications. Ž2. The first harmonic wave amplitude on the lee side of the dike becomes smaller as the Ursell number of the incident wave increases and the water depth ratio of the dike decreases. However, the shape of the dike plays also an important role. Due to the gentler depth change, the amplitudes of higher harmonics on the lee side of a trapezoidal dike are smaller than those for a rectangular dike. On the other hand, the length aspect ratio of the dike seems not to affect the result very much, except that as the dike becomes longer, a ‘recurrence length’ of higher harmonics becomes apparent. Ž3. Flow separation occurs at the up-wave edge Žleft corner. of a rectangular dike as a wave crest passes over. Conversely, it also occurs at the down-wave edge Žright corner. of the dike as the wave trough propagates over this point. The region of flow separation expands as the K–C number or the Reynolds number of the incident wave increases. In the extreme cases of our study, the separation zone develops into a separation bubble. Even so, the effect of the flow separation does not reach the free surface.
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
147
Ž4. An increase in the dike length seems not to affect the magnitude of the flow separation. For this reason, we suggest that the Reynolds number for a propagating wave field be defined as Ž Hi L2 .rŽy h 0 T .. Ž5. The strength of the vorticity increases as the K–C number increases. It should be noted that formation of vorticity was also influenced by the water depth ratio and the shape of the dike. If the height of the dike is small, then the vorticity is weak. Under the same K–C number and the same water depth ratio, voticity was generated around a rectangular dike, but this did not appear in the case of a trapezoidal dike with gentle slopes. Ž6. The generation of vorticity with respect to the variation of free surface has also been explained in this paper. Acknowledgements This research was supported by the National Science Council, Taiwan, R.O.C. under Contract No. NSC 86-2611-E-006-017. Appendix A. Definition of the Reynolds number for viscous wave fields When the viscosity of the fluid is taken into account in wave fields, we often need to use the Reynolds number, which is defined as the ratio of the inertial force to the viscous force acting on the fluid particles. From the order of magnitude analysis, we can estimate that the inertial force term uŽEurE x . is on the order of u˜ 2rL, where u˜ is the maximum horizontal velocity of the fluid and L is the wavelength. Similarly, the order of magnitude of the viscous force term y ŽE 2 urE x 2 . is y urL ˜ 2 . In the above estimation, we assume the wave is long wave, hence the dominant velocity gradient is in the horizontal direction. The Reynolds number is then uL ˜ Re s Ž A1. y The horizontal velocity of fluid particles under the progressive wave is: Hi
Re s
v
cosh k Ž z q h .
cos Ž kx y v t . Ž A2. 2 sinh kh The maximum value of u is: Hi v Hi c u˜ s ; Ž A3. 2 kh h In the above equation, we have assumed that the wave is a shallow water wave Žsinh kh f kh.. Substituting Eq. ŽA3. into Eq. ŽA1. we have: us
LHi c
Hi L2
s Ž A4. yh y hT This is the Reynolds number that we suggested for use with viscous wave fields.
148
C.-J. Huang, C.-M. Dong r Coastal Engineering 37 (1999) 123–148
References Battjes, J.A., Beji, S., 1991. Spectral evolution in waves traveling over a shoal. In: Peregrine, D.H. ŽEd.., Proc. of the Nonlinear Water Waves Workshop. Univ. of Bristol, Bristol, pp. 11–19. Beji, S., Battjes, J.A., 1993. Experimental investigation of the wave propagation over a bar. Coastal Eng. 19, 151–162. Beji, S., Battjes, J.A., 1994. Numerical simulation of nonlinear wave propagation over a bar. Coastal Eng. 23, 1–16. Chan, R.K.C., Street, R.L., 1970. A computer study of finite-amplitude water waves. J. Comput. Phys. 6, 68–94. Chen, C.J., Chen, H.C., 1982. The finite-analytic method. IIHR Report 232-IV, Iowa Institute of Hydraulic Research, The Univ. of Iowa. Chen, H.C., Patel, V.C., 1987. Laminar flow at the trailing edge of a flat plate. AIAA J. 25, 920–928. Dean, R.G., Dalrymple, R.A., 1984. Water wave mechanics for engineers and scientists. Prentice-Hall, Englewood Cliffs, NJ. Dingemans, M.W., 1997. Water wave propagation over uneven bottoms: Part II. World Scientific Publishing, Singapore. Driscoll, A.M., Dalrymple, R.A., Grill, S.T., 1992. Harmonic generation and transmission past a submerged rectangular obstacle. Proc. 23rd Int. Coastal Eng. Conf., Venice, ASCE, pp. 1142–1152. Grue, J., 1992. Nonlinear water waves at a submerged obstacle or bottom topography. J. Fluid Mech. 244, 455–476. Huang, C.J., Zhang, E.C., Lee, J.F., 1998. Numerical simulation of nonlinear viscous wavefields generated by a piston-type wavemaker. J. Eng. Mech. 124, 1110–1120. Johnson, J.W., Fuchs, R.A., Morison, J.R., 1951. The damping action of submerged breakwaters. Trans. Am. Geophys. Union 32, 704–718. Jolas, P., 1960. Passage de la houle sur un seuil. Houille Blanche 2, 148–152. Kittitanasuan, W., Goda, Y., Shiobara, T., 1993. Deformation of nonlinear waves on a rectangular step. Coastal Eng. in Japan 36, 133–153. Madsen, O.S., 1971. On the generation of long waves. J. Geophys. Res. 76, 8672–8683. Madsen, P.A., Murray, R., Sorensen, O.R., 1991. A new form of the Boussinesq equations with improved linear dispersion characteristics, Part 1. Coastal Eng. 15, 371–388. Madsen, P.A., Sorensen, O.R., 1992. A new form of the Boussinesq equations with improved linear dispersion characteristics, Part 2. Coastal Eng. 18, 183–204. Ohyama, T., Nadaoka, K., 1991. Development of a numerical wave tank for analysis of nonlinear and irregular wave field. Fluid Dynamics Research 8, 231–251. Ohyama, T., Nadaoka, K., 1992. Modeling the transformation of nonlinear waves passing over a submerged dike. Proc. 23rd Int. Coastal Eng. Conf., Venice, ASCE, pp. 526–539. Ohyama, T., Nadaoka, K., 1994. Transformation of a nonlinear wave train passing over a submerged shelf without breaking. Coastal Eng. 24, 1–22. Ohyama, T., Beji, S., Nadaoka, K., Battjes, J.A., 1994. Experimental verification of numerical model for nonlinear wave evolution. J. Waterway, Port, Coastal, and Ocean Eng. 120, 637–644. Peregrine, D.H., 1967. Long waves on a beach. J. Fluid Mech. 27, 815–827. Patankar, S.V., 1980. Numerical heat transfer and fluid flow. McGraw-Hill, New York. Ting, F.C.K., Kim, Y.K., 1994. Vortex generation in water waves propagating over a submerged obstacle. Coastal Eng. 24, 23–49.