# Boundary conditions for 2D Boussinesq-type wave-current interaction equations.

Introduction

Boussinesq  developed the original formulation of the governing equations for a free surface flow, which included the effects of surface waves but in which the vertical dimension was eliminated. The formulation was in terms of the bottom velocity and was restricted to simulating waves moving over bathymetry with a flat bottom. The governing equations consist of one continuity equation and two momentum equations (in x and y directions). The governing equations were then called as Boussinesq equations.

Peregrine  developed two new formulations in two horizontal dimensions for the case of varying depth in terms of (i) the depth-averaged velocity vector and (ii) the velocity vector at still water level. The first formulation became known as the standard form of Boussinesq-type equations.

In 1993, Nwogu  developed a new approach in the derivation of a novel set of Boussinesq-type equations. The resulting equations are capable of simulating wave propagating over arbitrary bathymetry in terms of the horizontal velocity at an arbitrary level (Z = [Z.sub.1]) below still water level (-h [less than or equal to] [Z.sub.1] [less than or equal to] 0) in which h is the local still water depth. The Boussinesq-type wave equations by Nwogu were solved numerically in one dimension by Nwogu and in two dimensions by Wei and Kirby .

Derivation of two-dimensional (2D) boundary conditions for the governing equations of Nwogu was also discussed by Wei and Kirby.

In 1998, Chen et al.  extended the Boussinesq-type wave equations of Nwogu by incorporating a steady ambient current. As a result, the corresponding equations were capable of simulating wave-current interaction. The continuity equation can be expressed as

[[eta].sub.t] + [nabla]*[(h + [eta]) u] + [PI]= 0 (1)

and the momentum equation is

[u.sub.t] + g[nabla][eta] + (u * [nabla])u + [[LAMBDA].sup.t] + [[LAMBDA].sup.s] = 0 (2)

where

g = gravitational acceleration

[eta] = free surface elevation

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

u = horizontal velocity vector, u(x,y,z,t) = (u,v)

u = horizontal velocity in the x direction, u(x,y,z,t)

v = horizontal velocity in the y direction, v(x,y,z,t)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

The one dimensional (1D) form of these equations had been solved numerically by Chen et al. . However, the 2D one has not been solved yet. To do so, the suitable 2D boundary conditions need to be derived first. This is the objectives of the present study. Discussion regarding the numerical solution algorithm is not presented here, but can be found in Mera .

Boundary conditions for waves only case

The set of 2D boundary conditions for waves only case is discussed first. The other cases of currents only and waves plus currents will be considered in subsequent sub-section.

Incoming wave boundary conditions

The free surface elevation [eta] at the incoming wave boundary can be varied sinusoidally as

[eta] = 1/2 [H.sub.i] sin (k * x - [omega]t) (3) k * x = (k cos x + (k sin [[theta].sub.i]) y

In which [H.sub.i], is the incoming wave height, k, the wave number vector, x, the horizontal spatial vector, [omega], the angular frequency = 2[pi]/T, T, the wave period, and t is time.

For a locally constant depth, the continuity equation (Eq. 1) reduces to

[[eta].sub.t] + (h + [eta])([nabla]*u) + u*[nabla][eta] + [([alpha] + 1/3) [h.sup.3] + [alpha][h.sup.2][eta] - 1/2 h[[eta].sup.2] - 1/6 [[eta].sup.3]] [nabla]*[[nabla]([nabla]*u)] = 0 (4)

where

[alpha] = 1/2 [([Z.sub.2]).sup.2] + [Z.sub.2] 0.5 [less than or equal to] [alpha] [less than or equal to] 0

[Z.sub.2] = [Z.sub.1]/h -1 [less than or equal to] [Z.sub.2] [less than or equal to] 0

A periodic, small amplitude wave can be defined as

[eta] = [[eta].sub.a] exp[i(k * x) - [omega]t], u = [u.sub.a] exp[i(k * x) - [omega]t] (5)

where:

[[eta].sub.a] = amplitude of the water surface elevation;

[u.sub.a] = amplitude of the horizontal velocity vector;

i = [square root of -1]

[theta]i = incoming wave angle.

Substituting Equation 5 into Equation 4 gives the horizontal velocity vector at the incoming wave boundary, i.e.

u = [omega][eta]/k{h - [k.sup.2] [([alpha] / 1/3)[h.sup.3] + [alpha][h.sup.2][eta] - 1/2 h[[eta].sup.2] - 1/6 [[eta].sup.3]]} (6)

Outgoing wave boundary conditions

The boundary condition for [eta] is considered first. At the outgoing wave boundary, the Sommerfeld radiation condition  is used to allow the passage and egress of the wave energy, that is

[[eta].sub.t] + C * [nabla][eta] = 0 (7)

where the wave celerity C = |C|cos [theta]i + |C|sin [theta]i. The boundary condition for the velocity components is considered next. The depth-integrated continuity equation can be expressed in terms of the depth-averaged horizontal velocity vector [bar.u] as

[[eta].sub.t] + [nabla] * [(h + [eta])[bar.u]] = 0 (8)

To eliminate [[eta].sub.t], substitute Equation 8 into the Sommerfeld radiation condition (Eq. 7) to give

[nabla] * [(h + [eta])u ] = C * [nabla][eta] (9)

For a locally constant depth, Equation 9 may be integrated over the fluid domain to obtain the horizontal velocity vector

[bar.u] = C [eta]/h + [eta] (10)

Having solved for [bar.u] in Equation 10 is used to determined u.

Reflecting wave boundary conditions

The kinematic boundary condition at an impermeable wall can be stated as

u * n = 0 x [member of] [partial derivative][OMEGA] (11)

where n is an outward normal vector, [OMEGA] is the fluid domain, [partial derivative][OMEGA] is the boundary and x is a position in the boundary. Consider, for example, the case of an impermeable wall being parallel to the x-axis. Equation 11 is a boundary condition and can be written as

v = 0 x [member of] [partial derivative][OMEGA] (12)

The slope and the curvature of v normal to the impermeable wall is assumed to be zero and expressed respectively as

[partial derivative]V/[partial derivative]Y and [[partial derivative].sup.2]v/[partial derivative][y.sup.2] = 0 x [member of] [partial derivative][OMEGA] (13)

or can be written symply as

[v.sub.y] = 0 and [v.sub.yy] = 0 x [member of] [partial derivative][OMEGA] (14)

in which the subscript y denotes partial differentiation with respect to the y direction

The continuity equation (Eq. 1) can be expressed in terms of the volume flux vector Q as

[[eta].sub.t] + [nabla] * Q = 0 (15)

where

Q = (h + [eta])(u + [GAMMA]) - 1/6 ([h.sup.3] + [[eta].sup.3])[nabla]([nabla]*u) + 1/2 ([h.sup.2] - [[eta].sup.2])[nabla] [[nabla]*(hu)]

The kinematic boundary condition in terms of the volume flux vector at an impermeable wall as

Q * n = 0 x [member of] [partial derivative][OMEGA] (16)

For the case of the impermeable wall being parallel to the x-axis, the volume flux in the y-direction at the boundary becomes zero or

[Q.sup.y] = 0 x [member of] [partial derivative][OMEGA] (17)

i.e.

(h + [eta])v + [1/2 [([Z.sub.1]).sup.2] (h + [eta]) - 1/6([h.sup.3] + [[eta].sup.3])]([u.sub.yy] + [v.sub.yy]) + [[z.sub.1](h + [eta]) + 1/2([h.sup.2] + [[eta].sup.2])][[(hu).sub.xy] + [(hv).sub.yy]] = 0 x [member of] [partial derivative][OMEGA] (18)

where the subscripts x and y denote partial differentiation with respect to the x and y directions respectively. Substituting Equation 18 into Equation 1 gives a reflecting wave boundary condition for calculating the free surface elevation at the boundary wall as set out below

[[eta].sub.t] + [(h + [eta])u]x + [[PI].sup.x] = 0 x [member of] [partial derivative][OMEGA] (19)

where

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

For a locally constant depth, the horizontal velocity in the x-direction may be obtained by substituting Equations 12 and 13 into Equation 18 giving

[u.sub.xy] = 0 x [member of] [partial derivative][OMEGA] (20)

For a boundary being parallel to the x-axis, the boundary conditions are Equations 19, 20 and 12 for [eta], u and v respectively. The derivation of the present boundary conditions for waves only case follows Mera , who derived a set of boundary conditions for the Boussinesq-type wave equations of Nwogu.

Boundary conditions for currents only case

Inflow Boundary Conditions

At the upstream end, the depth-averaged velocity is specified but the boundary condition also needs to involve [eta]. Inflow currents are bound to flow in the xdirection. One way of linking [bar.u] and [eta] at the upstream end is to combine the Sommerfeld radiation condition (Eq. 7) and the continuity equation (Eq. 8) to give

2[[eta].sub.t] + [(h + [eta])[bar.u]]x + C[[eta].sub.x] = 0 x [member of] [partial derivative][OMEGA] (21)

and

v = 0 x [member of] [partial derivative][OMEGA] (22)

Outflow Boundary Conditions

At the downstream end, the free surface elevation can be predicted using Equation 7 in the x direction or

[[eta].sub.t] + C[[eta].sub.x] = 0 (23)

Then, the horizontal velocities are predicted using Equation 16 in the x direction or

[bar.u] = C [eta]/ h + [eta] (24)

and in the y direction

[bar.v] = 0 (25)

No-flow Boundary Conditions

At the boundaries, which are paralel to the flow, the boundary conditions are equivalent to the reflecting wave boundary conditions.

Boundary conditions for wave-current interaction cases

The governing equations considered here were derived based on a steady ambient current. In the model tests, the following procedure is adopted.

(1) Model is run with current only.

(2) The results from the model settle down to a steady state.

(3) After the steady state is reached, a sinusoidally varying surface elevation is imposed at the inflow or outflow boundary. This results in a wave train propagating into the computational domain.

Please note that breaking waves are not included in the existing governing equations and in the present boundary conditions considered in this paper. But, Purwanto  discussed a Quasi-Equilibrium Turbulent Energy (QETE) model with boundary conditions for breaking and non-breaking waves. The QETE model was intended for simulation of (density) current only case. The density current is a type of current that occurs when fluid flow enters a fluid body of different density.

Numerical set-up

The numerical set-up consists of a basin with a submerged shoal. The basin is 18 m long and 18.2 m wide (Figure 1). Side walls are at y = 0 and 18.2 m. The centre of the shoal was located at (x,y) = (13,9.22) m with the perimeter given by

[(x - 13).sup.2] + [(y - 9.22).sup.2] = [(2.57).sup.2] (26)

[FIGURE 1 OMITTED]

Current only case

In the first test, a flat water surface (initial value of [eta] = 0) and a steady velocity of 0.10 m/s is imposed at the Southern boundary (x = 0 m). The computation is carried out with a mesh of [DELTA]x = [DELTA]y = 0.1 m and [DELTA]t = 0.02 s. The imposed current flows from x = 0 m to x = 18 m, and reaches a steady state condition after about t = 65 s. Figures 2 shows a perspective view of the free surface elevation (upside down) at t = 65 s predicted by a numerical model based on the existing governing equations of Chen et al.  and the present boundary conditions. This figure indicates that the present boundary conditions go well with the existing governing equations.

[FIGURE 2 OMITTED]

In the second test, a flat water surface (initial value of [eta] = 0) and a constant inflowing velocity of 0.10 m/s is imposed in the opposite direction to that in first test, i.e. at the Northern boundary (x = 18 m) instead of at the Southern boundary. This leads to a steady current flowing from x = 18 m to x = 0 m of the basin (not presented here).

Waves and opposing current

Once the currents in the basin reach a steady state (after about t = 65 s) in the first test, the free surface elevation at the southern boundary (x = 18 m) is sinusoidally varied with time to generate an incident wave. The incoming wave specifications and the grid resolution remain the same as is used in the previous tests i.e. [DELTA]x = [DELTA]y = 0.1 m and [DELTA]t = 0.02 s. A wave train with a period of 1.0 s and a wave height of 0.0118 m comes prependicularly to the fluid domain. At the incoming wave boundary (x = 18 m), the ambient current is allowed to pass through, leaving the flow domain. A perspective view of the predicted surface elevation at t = 20 s and at t = 40 s are shown in Figures 3 and 4 respectively. From these figures can be seen that the present boundary conditions are still suitable to the existing governing equations.

[FIGURE 3 OMITTED]

[FIGURE 4 OMITTED]

Waves and current in same direction

On top of the steady current field (second test), a sinusoidal wave train with a period of 1.0 s and a wave height of 0.0118 m is imposed at x = 18 m. The incoming wave period and wave height are same as is used in the last test. Perspective views of the predicted free surface elevation at t = 20 s and at t = 40 s are shown in Figures 5 and 6 respectively in which the present boundary conditions are still going well with the existing governing equations.

[FIGURE 5 OMITTED]

[FIGURE 6 OMITTED]

Conclusions

A set of two-dimensional boundary conditions for the existing Boussinesq-type wave-current interaction equations is developed. The boundary conditions consist of for waves only case and for currents only case. Two kinds of the present boundary conditions are then combined in simulation of wave-current interaction.

A numerical set-up consists of a basin of 18 m long and 18.2 m wide with a submerged shoal. The numerical model runs stably at least for 40 s in simulation time as shown by the numerical results. This indicates that the present boundary conditions are suitable to the existing governing equations.

Acknowledgements

The author would like to thanks to Dr. Bruce Cathers from the University of New South Wales for his comments regarding the present material.

Received 9 April 2010; revised 6 November; accepted 28 November 2010.

References

[1.] Boussinesq, J., Theorie des ondes et des ramous qui se propagent le long d'un canal rectangulaire horizontal, en communi quant au liquide contenu dans ce canal des vitesses sensiblement pareilles de la surface au fond. J. Math. Pures Appl. 2nd Series 17, 1872, pp. 55-108.

[2.] Peregrine, D.H., Long Waves on a Beach. J. Fluid Mechanics 27(4), 1967, pp. 815-827.

[3.] Nwogu, O., Alternative Form of Boussinesq Equations for Nearshore Wave Propagation. J. Waterway Port Coastal and Ocean Eng, ASCE 119(6), 1993, pp. 618-638.

[4.] Wei, G. and Kirby, J.T., Time-dependent Numerical Code for Extended Boussinesq Equations, J. Waterway Port Coastal and Ocean Eng, ASCE 121(5), 1995, pp. 251-261.

[5.] Chen, Q., Madsen, P.A., Schaffer, H.A., and Basco, D.R., Wave-current Interaction Based on an Enhanced Boussinesq Approach. Coastal Engineering 33, 1998, pp. 11-39.

[6.] Mera, M., Boussinesq-type Numerical Models, Ph.D Thesis at the University of New South Wales, Sydney, Australia, 2002.

[7.] Mera, M., 2D Horizontal Weakly Non-linear Boussinesq-type Numerical Wave Model, Jurnal Teknologi Kelautan, Vol. 11, No. 1, 2007, pp. 3748.

[8.] Purwanto, B. S., Wind Affected Density Current Profile in a Small Semi-Enclosed Water Body, Civil Engineering Dimension, Vol. 11, No. 1, 2009, pp. 58-68.

Note: Discussion is expected before June, 1st 2011, and will be published in the "Civil Engineering Dimension" volume 13, number 2, September 2011.

Mera, M. (1)

(1) Civil Engineering, University of Andalas, Padang, INDONESIA Email: masmera@ft.unand.ac.id
Author: Printer friendly Cite/link Email Feedback Mera, M. Civil Engineering Dimension Report 9INDO Mar 1, 2011 2659 Numerical simulation of threat-independent progressive collapse. Kriging-based Timoshenko Beam Element for Static and Free vibration analyses. Approximation Approximation theory Boundary value problems Wave equation Wave equations