THREE-DIMENSIONAL NUMERICAL MODELLING OF THE WAVE-INDUCED RIP CURRENTS UNDER IRREGULAR BATHYMETRY*

XIE Ming-xiao

School of Civil Engineering, Tianjin University, Tianjin 300072, China

Tianjin Research Institute for Water Transport Engineering, Ministry of Transport, Tianjin 300456, China, E-mail: crabsaver@163.com

(Received March 30, 2012, Revised May 20, 2012)

THREE-DIMENSIONAL NUMERICAL MODELLING OF THE WAVE-INDUCED RIP CURRENTS UNDER IRREGULAR BATHYMETRY*

XIE Ming-xiao

School of Civil Engineering, Tianjin University, Tianjin 300072, China

Tianjin Research Institute for Water Transport Engineering, Ministry of Transport, Tianjin 300456, China, E-mail: crabsaver@163.com

(Received March 30, 2012, Revised May 20, 2012)

A process-based 3-D hydrodynamic model is established to simulate the rip current structures under irregular bathymetry. The depth-varying wave-induced residual momentum, the surface rollers, the turbulent mixing and the wave-current interactions are considered. Experimental datasets are used to validate the model, and it is shown that the model can effectively describe the 3-D structures of the rip currents in both normal and oblique wave incident cases. The flow patterns of the rip currents see various characteristics for different incident wave directions. In the normal incident case, pairs of counter-rotating primary circulation cells are formed, and an offshore rip flow occurs in the embayment troughs. The peak seaward velocities occur at the top of the bed boundary layer, and the undertow is incorporated in addition to the rip currents. In the oblique incident case, the longshore currents are dominant, which result in a meandering flow along the depth contour, and the undertow is weaker compared to that in the normal incident condition.

rip currents, irregular bathymetry, 3-D numerical modelling

Introduction

The wave-induced currents generated in coastal regions are directly responsible for sediment transport and morphology evolutions. Therefore an accurate prediction of these currents is essential for coastal engineering applications. Wave-induced current phenomena were widely studied, including the wave setup[1,2], undertow[3]and longshore currents[4,5]. However, only the single-sloped plane beaches were often considered in the experiments. In real coastal areas, the underwater bathymetry is usually irregular (e.g., with rhythmic shorelines, sinuous cusps or gapped sand bars). As the waves propagate towards the shoreline, under the multiple processes of shoaling, breaking, refraction and diffraction, a more complex circulative flow could form, which is called the rip currents, which can be very intense with the highest velocity exceeding 2 m/s in sand bar gaps or cusp embayment troughs, to induce the strong seaward sediment transport, and subsequently to affect the coastal morpho-logy.

The structure of the rip currents is extremely complicated. Under different incident wave and bathymetry conditions, the flow patterns see variable features. Due to its complex nature, there were relatively few detailed explorations of the rip currents. Haller et al.[6]measured the rip currents on an artificial barred beach in laboratory using a wave tank and ADV, Peng and Zou[7]measured the rip currents on a similar sandbar beach using video-tracked drifters. But in those experiments only the horizontal velocity distribution was measured, leaving the vertical profiles aside. In order to investigate the detailed characteristics of the rip currents, Borthwick and Foote[8]installed a tri-cuspate beach in the UK Coastal Research Facility (UKCRF), and measured the 3-D structures of the current field.

As for numerical simulations, Bai et al.[9]modeled the rip currents under a rhythmic bathymetry using the quasi-3D SHORECIRC model, Rogers et al.[10]modeled the UKCRF experiments using the Godunovtype 2DH numerical model, Fang et al.[11]applied the Boussinesq equations to model the rip currents over the barred beach. However, as indicated by Borthwickand Foote[8], the rip currents have a significant 3-D nature, 2-D models would be inadequate, and a fully 3-D model should be used.

Xie[12]established a fully process-based 3-D wave-induced current model, which was satisfactorily validated using a series of experimental datasets. In this paper, the Xie[12]model was applied to simulate the flow structures of the rip currents including both the horizontal flow pattern and the vertical profile.

1. Model description

1.1 Hydrodynamic model

The governing equations of the hydrodynamic model are in the Reynolds form simplified from the original Navier-Stokes equations (see Eq.(1)-Eq.(4)). The contributions of the depth-varying residual momentum, the surface roller stresses and the turbulent mixings are included. Cartesian coordinates are used in the horizontal directions and the terrain-following sigma coordinate is used in the vertical direction.

where the vertical sigma coordinate σ=(z-η)/D ranges from σ=-1 at the bottom to σ=0 at the surface, t is the time,x and y are the horizontal coordinates,η is the free surface, U and V are the velocity components in x and ydirections, respectively,ω is the velocity component in σ coordinate, D is the water depth, g is the gravity acceleration, p is the pressure,M is the depth-varying residual momentum, R is the depth-varying roller momentum,KMcand AMcare the vertical and horizontal mixing coefficients combining waves and currents, respectively, ρ is the seawater density.

1.1.1 Wave-induced residual momentum

The formulation proposed by Lin and Zhang[13]is applied for the vertical distribution of the wave-induced residual momentum (see Eq.(5))

where E is the wave energy, n is the wave energy transfer rate, k is the wave number,δ is the Kronecker symbol, i and j represent the x,y directions, respectively.

1.1.2 Surface roller evolution

Based on the energy balance, Xie[12]derived an evolution model of the breaking-induced surface rollers, as expressed by Eq.(6). The model considers multiple factors including the roller energy transfer, the roller density, the bottom dissipation and the bed slope.

where Cg=Cn is the wave group celerity,αis the roller energy transfer factor, T is the wave period, ER=ρARC/2T is the roller energy, KR= 3(0.3+2.4s)/8, and s is the bed slope, ARis the roller area, n =n(cosθ,sinθ) is the wave vector, andρRis the roller density.

If the wave parameters are given, Eq.(6) can be solved by using an iteration algorithm from the brea-king point to the shoreline with an offshore boundary condition AR=0. The vertical profile of the roller momentumis expressed as an exponential function proposed by Haas and Warner[14], as expressed by Eq.(7). Note that because the depth integral of Rzshould be unity, it should be pre-normalized to Rznfollowing Eq.(8).

OnceAR(σ)is solved,ER(σ)can be calculated explicitly. The corresponding stresses in the governing equations caused by the roller can be determined, as in Eqs.(9)-(11).

1.1.3 Bottom shear stress

The wave-current combined bottom shear stress τcwis determined by Eq.(12). where τc=ρCDuc2is the bed shear stress by current only, CD=[κ-1ln(h+zb)/z0]-2is the drag coefficient, in which κ=0.4 is the Von Karman constant, h is the bed elevation, zbis the elevation of the first grid point above the bottom, z0is the roughness height, ucis the current velocityat the grid point nearest the bed, τw=0.5ρfwuw2is the shear stress due to waves only,uw=Hπ/Tsinh(kD) is the nearbottom wave orbitalvelocity, fwis the wave friction factor, B, P, Q are empirical coefficients.

The value of coefficient B is determined by Eq.(13) with analogous expressions for P and Q.

1.1.4 Turbulence mixing

The wave-current combined turbulent mixing

coefficients can be expres sed as in Eqs.(14)-(15).

where A and K represent the horizontal and the vertical turbulentmixing coefficients, respectively, and thesubscriptsM and W represent current and waves, respectively.

The horizontal mixing coefficient AMfor currents only is given by Eq.(16).

where ΔxandΔy are the horizontal grid steps, Csis an empirical factor. The current-induced vertica l mixingcoefficientKMis solved by using a Mellor-Yamada closure model.

Using the linear wave theory, Xie[12]derived the horizontal mixing coefficientAWfor waves only, as in Eq.(17)

The vertical mixing coefficient AWfor waves onlyis expressed as in Eq.(18)

where b is a calibration coefficient.

1.2 Wave model

The combined refraction/diffraction wavemodel (REF/DIF) is used as the wave driver for simulating monochromatic incident waves. The REF/DIF model is based on the parabolic mild-slope equation, and it can involve many processes, e.g. shoaling, refraction, energy dissipation, and irregular bottom bathymetry.

1.3 Wave-current interaction

The flow pattern of the rip currents is extremely complex, and the strong opposing currents affect the wavepropagation significantly. Therefore, the mutual interaction of the waves and currents should be considered in the simulation. In this paper, the REF/DIF procedure and the hydrodynamic procedure are coupled together through an iterative algorithm. After the current field reaches a stable state, the U and V fields feed back to the wave solver, and consequently the new wave parameters are calculated forthe preparation of the wave-related stresses, which are then in-corporated into the hydrodynamic equations for the solution in the next time step.

1.4 Solution technique

A finite difference method and a time-splitting technique are applied to solve the governing equations. The horizontal terms are treated explicitly, and the vertical terms are treated implicitly by using a doublesweep scheme. The arrangements of the variables follow the staggered C-grid system. The OGCM approach proposed by Oey[15]is used to model the inundation.

Fig.1 Bathymetry and the observation stations in the UKCRF experiment

Table 1 Incident wave parameters for the normal and oblique incident cases

2. Experimental cases of UKCRF

Borthwick and Foote[8]carried out laboratory studies of the 3-D structure of the rip currents over a tri-cuspate beach using the UKCRF. The wave basin hasthe plan dimensions of 27 m cross-shore by 36 m alongshore, and with the still water depth at the paddles of 0.5 m. The bed slope is 1:20. The experimental bathymetry and the locations of the observation profiles are shown in Fig.1. In the experiments, 2 wave observation sections and 7 velocity observation profiles are arranged both in the embayment and on the cusp horn.

In this paper, the normal incident case (Case B) and the oblique incident case (Case C) are considered, and the related p arameters are shown in Table 1.

Table 2 Input parameters in the numerical simulation

Fig.2 Arrangement of the σ layers

3. Numerical modelling of the rip currents

3.1 Model parameters

Table 2 shows the input parameters in the numerical simulations. In order to better describe the nearbed distributions of the current speed, the varying sigma discretization is used, where the spacing in the upper water column is selected asΔσ=0.1, and the near-bottom spacing is Δσ=0.01. The detailed arrangement for the sigma layers is illustrated in Fig.2.

Fig.3 Comparisons between the modeled and the measured wave heights for Case B

In order to estimate the model errors, two indices are applied, which are the root mean square (rms) errorand the correlation coefficient (COR). The former reflects the deviation between the measured and the simulated values, and the latter represents the linear correlation of two datasets. They are expressed aswhere Nis the total number of measurement points, I refers to a measurement point, meis the measurement value,om is the modeled value,N andrepresent the algebraic mean of the measured and the modeled values, respectively.

Fig.4 Comparis ons between the modeled and the measured current velocities for Case B

3.2 Normal incident case

Figure 3 shows the comparisonsbetween the modeled and the observed wave heights of two representative sections for the normal incident case (Case B). It indicates that the model can describe the wave propagation, including the shoaling and breaking processes. The rms errors are confined within the range of 0.01 m-0.02 m, and the correlations between the two datasets are satisfactory (88%-97%).

The comparisons with the observed velocities are shown in Fig.4. Theoretically, the V-velocities should be near-zero because of the symmetric nature of the bathymetry and the normal incident wave condition. However, the observed V-velocities are scattered. In fact, in the experiment, the rip current has unstable features and a trivial perturbation could lead to a deflection of the current direction. With above considerations in mind, in the comparisons, only the U-component is selected in the evaluation. As for the U-velocities, the rms errors are in the range of 0.03 m/s-0.06 m/s, and especially, the simulated velocities at sections P6 and P7 are larger than those observed. One reason is that the gradients of the observed wave heights are greater than those simulated in this area, which induces higher velocities (see Fig.3). The correlations between the two datasets seem not satisfactory, in which the lowest value is COR=12% (section P4). That is because the unstable nature of the rip currents makes the measurement data extremely scattered and their vertical variations are not smooth enough to infer the distribution trends clearly as compared to the simulation values. However, the comparisons do show that the simulated velocity profile of the rip currents captures the major distribution trend.

Additionally, it can be observed that both the measured and the simulated results indicate that in the embayment (P1-P4), the peak seadirected velocities do not occur at the bottom, but at some distance from the bed, say 0.8z/D-0.95z/D. That is because in the surfzone, the undertow also contributes to theflow structure. It could be imposed on the rip currents and make the maximum velocity occur at the top of the bed boundary layer.

Fig.5Planar distribution of the depth-averaged current velocities for Case B

The horizontal distributions of the rip currents are extremely complex due to the involvement of many processes e.g. the shoaling, breaking, refraction and diffraction. Under the impact of irregular bathymetry, the wave heights differ in bothxand ydirections while propagating onshore, consequently, residual momentum gradients are formed in each direction.

The flow field of depth-averaged rip currents for t he normal incident case is illustrated in Fig.5. It shows that pairs of counter-rotating circulation cells occur in each embayment. On the cusp horns, there is an onshore flow which fans out, divides and then feeds into the longshore currents that meet to form seaward rip currents at the embayment troughs. The rip currents are restricted in a relatively narrow zone, and then flow offshore with a large velocity (the maximum depth-averaged velocity of 0.3 m/s). The rip currents reach a short distance offshore in front of the breaking line, and die away at the rip heads. The comparison shows that the distribution of the simulated flow field agrees with that observed in the experiment.

Fig.6Comparison between the modeled and the measuredwave heights for Case C

3.3 Oblique incident case

Figure 6 gives comparisons between the modeled and the observed wave heights for the oblique incident case (Case C). Similar to the normal incidentcase, it isshown that the model can effectively describe the wave propagations with the rms errors in the range of 0.01m-0.02 m, and the correlations in the range of 86%-95%.

Fig.7 Comparis on between the modeled and the measured current velocities for Case C

The comparisons with the observed velocities are shown in Fig.7. Unlike the normal incidentcase, in the oblique incident case, both the U andV components are significant, hence, the rms errors and the correlations for each direction are given. It is estimated t hat t he r ms errors are in th e ran ge of 0.01 m /s-0.09m/sforthe U-velocities,andinthe rangeof 0.02 m/s-0.09m/s for theV-velocities. The maximum error occurs at P6, where the model overestimates the velocity magnitude. The correlations for most datasets aresatisfactory, with the worst correlation at P7 (–6% forV-velocity) because the curvature of the measured data is opposite to that in the model. Generally, the simulated velocity profile of the rip currents captures the major distribution trendsfor both the magnitude and the distribution characteristics.

The flow field of the depth-averaged rip currents for the oblique incident case is illustratedin Fig.8. It is indicated that the nearshore flow pattern is significantly different from that for the normal incident case. Because the incident waves are oblique, the gradientof the alongshore residual momentum contributes most to the nearshore currents. As a result, both the cell-like circulation structure and the rip currents are smoothened due to the strong longshore currents, and the flow is meandering along the bed contours. The maximum depth-averaged velocity is 0.48 m/s. The undertow is weaker compared to that in the normal incident condition. The comparison shows that the distribution of the simulated flow filed also agrees with that observed in the experiment.

Fig.8planar distribution of the depth-averaged current velocities for Case C.

To summarize, the flow structure of the rip currents under irregular bathymetry is extremely complex due to the multiple coastal processes and the wavecurrent interactions. All these factors make the laboratory measurements and the simulations very difficult. However, the comparisons with data of differentwave incident cases (normal and oblique) show that the established process-based 3-D numerical model can capture the major characteristics of the rip current field effectively for both the horizontal layout and the vertical profile. Generally, the model could provide someproper hydrodynamic information for further investigation of the morphodynamics in coastal areas.

4. Conclusions

(1) Using the process-based 3-D wave-induced current model, the rip current structures under irregular bathymetry were simulated. In the model, many processes including the depth-varying wave residual momentum, the surface rollers, the wave turbulent mixing and the wave-current interactions are considered.

(2) The comparisons with the laboratory measurement datasets indicate that the model can effectively describe the horizontal distribution and the vertical profile of the rip currents for both normal incident and oblique incident cases.

(3) The rip currents under irregular bathymetry show various characteristics under different wave incident conditions. For the normal incident case, pairs of counter-rotating circulation cells form, and offshore jet flows occur in the embayment troughs. The undertow contributes to the flow in the embayment and makes the peak seaward velocities occur at the top of the boundary layer. In the oblique incident case, the longshore currents are dominant, which results in a meandering flow along the depth contour, and the unde rtow is weaker compared to that in the normal inciden t condition.

Acknowledgements

This work was supported by the Central Public Institute Foundation of Tianjin Research Institute for Water Transport Engineering, Ministry of Transport (Grant No. TKS100102).

[1] HSU T., JOHN R. C. and WENG W. et al. Wave setup and setdown generated by obliquely incident waves[J]. Coastal Engineering, 2006, 53(10): 865-877.

[2] WEBER J. E., CHRISTENSEN K. H. and DENAMIEL C. Wave-induced setup of the mean surface over a sloping beach[J]. Continental Shelf Research, 2009, 29(11-12): 1448-1453.

[3] KURIYAMA Y., ITO Y. and YANAGISHIMA S. Cross-shore variation of long-term average longshore current velocity in the nearshore zone[J]. Continental Shelf Research, 2008, 28(3): 491-502.

[4] KURIYAMA Y., NAKATSUKASA T. A one-dimensional model for undertow and longshore current on a barred beach[J]. Coastal Engineering, 2000, 40(1): 39- 58.

[5] REN Chun-ping, ZOU Zhi-li and QIU Da-hong. Experimentalstudy of the instabilities of alongshore currents on plane beaches[J]. Coastal Engineering, 2012, 59(1): 72-89.

[6] HALLER M. C., DALRYMPLE R. A. and SVENDSEN I. A. Experimental study of nearshore dynamics on a barred beach with rip channels[J]. Journal of Geophysical Research, 2002, 107(14): 1-21.

[7]PENG Shi, ZOU Zhi-li. Experimental measurement of rip currents with video-tracked drifters[J]. Chinese Journal of Hydrodynamics, 2011, 26(6): 645-651(in Chinese).

[8]BORTHWICK A. G. L., FOOTE Y. L. M. Wave-induced currents at a tri-cuspate beach in the UKCRF[J]. Water and Maritime Engineering, 2002, 154(4): 251-263.

[9]BAI Zhi-gang, ZHANG Zhi-xian and CHEN Zhi-chun. A Quasi-3D nearshore circulation model applied in ripcurrent research[J]. Port and Waterway Engineering, 2007, (3): 12-17(in Chinese).

[10]ROGERS B. D., ALISTAIR G. L. and TAYLOR P. H. GODUNOV-type model of wave-induced nearshore currents at a multi-cusped beach in the UKCRF[C]. 28th International Conference of Coastal Engineering, Cardiff, Wales, UK, 2002, 760-771.

[11] FANG Ke-zhao, ZOU Zhi-li and LIU Zhong-bo. Numerical simulation of rip current generated on a barred beach[J]. Chinese Journal of Hydrodynamics, 2011, 26(4): 479-486(in Chinese).

[12] XIE Ming-xiao. Establishment, validation and discussions of a three dimensional wave-induced current model[J]. Ocean Modelling, 2011, 38(3-4): 230-243.

[13] LIN Peng-zhi, ZHANG Dan. The depth-dependent radiation stresses and their effect on coastal currents[C]. Proceedings of the 6th International Conference of Hydrodynamics: Hydrodynamics VI Theory and Applications. Perth, Australia, 2004, 247-253.

[14] HAAS K. A., WARNER J. C. Comparing a quasi-3D to a full 3D nearshore circulation model: SHORECIRC and ROMS[J]. Ocean Modeling, 2009, 26(1-2): 91- 103.

[15] OEY L. Y. An OGCM with movable land-sea boundaries[J]. Ocean Modelling, 2006, 13(2): 176-195.

10.1016/S1001-6058(11)60314-4

* Biography: XIE Ming-xiao (1982-), Male, Ph. D.