Modeling of on Shore Propagation of Random Water Waves

Abstract

Two numerical models are investigated to model random water waves (RWWs) transformation due to mild depth variation. Modelling of steady on-shore propagation of small-amplitude RWWs is based on superposition principle of waves of different heights and directions. Each component is simulated through either the parabolic model (PM) or the elliptic model (EM). PM simulates weak refraction, diffraction, shoaling, and wave breaking. EM simulates strong refraction, diffraction, and shoaling. Both models neglect wave reflection. Comparison between PM and EM, in test cases that are experimentally measured, proved that both models give good results for unidirectional and narrow-directional RWW. However, EM is more accurate in modelling broad-directional RWWs.

1. Introduction

In studying many coastal engineering problems, it is essential to have accurate information on wave conditions in the area of interest. These wave conditions, including the wave height and the dominant wave direction, are usually obtained through a wave transformation model that transfers the wave characteristics from the location where the wave data are collected to the site of concern. Because of depth variation, coastal currents, artificial structures, and geological features, waves change their propagation direction and speed and redistribute their energy along wave crests as waves approach the shore. Inside the surf zone, where breaking is an important feature, waves have more severe transformation.

Researchers have developed the mild-slope equation to simulate shoaling, refraction, reflection, diffraction, and breaking of regular water waves, when the rate of change of depth bottom slope up to 1 : 3 [1] and the current is small within wave length. It is a two-dimensional partial differential equation of elliptic type that can be solved as a boundary-value problem using specified appropriate boundary conditions. Its computational requirements are too large relative to the ray tracing method, since the entire domain must be solved simultaneously and the grid size must be small enough to allow eight to ten nodes within a wavelength. Ebersole presented a modified method as an efficient way to solve the elliptic form of the mild-slope equation [2].

Situation becomes more complex when modeling of random water waves (RWWs) transformation. The extended mild-slope equation of Suh et al. [3] and Lee et el. [4] are compared analytically and numerically to determine their applicability to random wave transformation. The phenomenon of breaking of RWW was treated also in many works; see, for example, [5], where two numerical formulations of breaking phenomenon were implemented in a numerical model for RWW propagation. A low-frequency spectrum in a harbour excited by short and random incident waves was modeled in [6].

An alternative numerical scheme, the parabolic approximation, has been developed and applied to the mild-slope equation to reduce the computational effort. However, the parabolic approximation has two disadvantages compared to the elliptic formulation; it assumes weak refraction and no wave reflection [7].

In the following sections, we will introduce both the parabolic (PM) and elliptic models (EM) for single wave and extend them to model RWW propagation. The incident RWW is decomposed into a spectrum of multiple waves. Each wave, has its different height and direction, is propagated by PM or EM, and then superposition is imposed. In Section 4, both models are compared to experimental data of Vincent and Briggs [8] to investigate their limitations and accuracy.

2. The Parabolic Model

Application of the parabolic equation method to investigate wave transformation needs that the propagation directions of all concerned components of the wave field to be confined to some narrow band of directions centered about the dominant propagation direction. The allowed directional bandwidth is limited by the maximum allowed error in the principal direction.

2.1. Derivation of the Parabolic Governing Equation

Kirby started with the mild-slope equation [9]𝜕²Φ𝜕𝑡²−∇.𝑐𝑐_𝑔∇Φ+𝜔²−𝑘²𝑐𝑐_𝑔+𝑖𝜎𝜀2Φ=0,(1)𝑐_𝑔 = group velocity; 𝑐 = wave celerity; 𝑘 = wave number; 𝜎= wave frequency;𝜀 breaking coefficient; ∇=(𝜕/𝜕𝑥,𝜕/𝜕𝑦).

The last term in (1) is dissipation function𝜀 to model frictional dissipation [10] or wave breaking [9]. Introducing the harmonic time representation of the wave potential as Φ=𝜑𝑒^{−𝑖𝜔𝑡},(2) where𝜔, angular frequency,𝑑𝜔𝑑𝑘=𝑐^𝑔𝑐,(3)₂𝑘²=𝜔²=𝑔𝑘tanh𝑘ℎ.(4) Then, (1) can be simplified as𝑐𝑐_𝑔𝜑_𝑥_𝑥+𝑐𝑐_𝑔𝜑_𝑦_𝑦+𝑘²𝑐𝑐_𝑔−𝑖𝜔2𝜀𝜑=0.(5) Rewrite (5) following the splitting method used by Kirby [7]𝜙_𝑥𝑥+𝑐𝑐^𝑔_𝑥𝑐𝑐_𝑔𝜙_𝑥+𝑘²1+𝑘𝑀²𝑐𝑐_𝑔𝜙=0,(6) where 𝑀𝜙=𝑐𝑐_𝑔𝜙_𝑦_𝑦−𝑖𝜎𝜀2𝜙.(7) The goal of the parabolic approximation is to split the reduced elliptic equation (6) for𝜙 into parabolic equations for a forward scattered wave 𝜙⁺ and a backward scattered wave𝜙⁻, where𝜙=𝜙⁺+𝜙⁻.(8) According to the parabolic assumption of no reflection the reflected wave 𝜙⁻ is neglected, and (6) is shown by Kirby and Dalrymple [9] to be parabolic equation of 𝜙⁺ because it has a first derivative in the longitudinal direction 𝑥𝜎𝑐_𝑔𝜙_𝑥+𝜎2𝑐^𝑔𝑥𝜙−𝑖𝑘𝜎𝑐_𝑔+34𝑀^𝜙−14𝑀^𝜙𝑘^𝑥𝑘²+𝑐^𝑔𝑥+12𝑘𝑐_𝑔4𝑘𝑀^𝜙_𝑥+𝜎𝜀2𝜙=0,(9) where 𝑀^𝜙=(𝑐𝑐_𝑔𝜙_𝑦)_𝑦.

The amplitude form of the parabolic equation is derived from (9) by making the substitution 𝜑=−𝑖𝑔2𝐴𝜎𝑒^{−𝑖−𝑘},(10) where 𝑘=(1/𝐵)∫^𝐵0𝑘(𝑥,𝑦)𝑑𝑦, where 𝐵 is the domain width 𝜕𝐴𝜕𝑥−𝑖𝑘−𝑘−𝑐𝑐^𝑔𝑘_𝑥+𝑐𝑐_𝑔2𝑐𝑐_𝑔𝑘−𝑐𝜀_𝑔𝐴_𝑦4𝑐𝑐_𝑔𝑘𝑖𝑘𝑘−3−𝑘𝑘^𝑥2−𝑐𝑐^𝑔𝑘_𝑥2𝑐𝑐_𝑔𝑘²+𝑐𝑐^𝑔_𝑦𝑥4𝑐𝑐_𝑔𝑘²𝐴^𝑦+𝑐𝑐^𝑔_𝑦4𝑐𝑐_𝑔𝑘²𝐴_𝑥𝑦+𝑐𝑐^𝑔_𝑦4𝑐𝑐_𝑔𝑘𝑖𝑘𝑘−3−𝑘𝑘^𝑥2−𝑐𝑐^𝑔𝑘_𝑥+𝑐𝑐_𝑔2𝑐𝑐_𝑔𝑘²_𝑥4𝑐𝑐_𝑔𝑘²𝐴^𝑦𝑦+14𝑘²𝐴_𝑦𝑦𝑥=0.(11) The finite difference scheme for (11) follows directly using the Crank-Nicolson method for performing an implicit update for each row in 𝑥-direction. Denote 𝑥 positions by “i” superscripts and 𝑦 positions “𝑗” subscripts. The computations proceed row by row by updating values of 𝐴 from the known “𝑖” row to the unknown “𝑖+1” row. The difference scheme will be 𝐶𝑃1^𝑖𝑗𝐴^{𝑖+1𝑗−1}+𝐶𝑃2^𝑖𝑗𝐴^𝑗𝑖+1+𝐶𝑃3^𝑗𝑖1𝐴^𝑖+1𝑗+1=𝐶1^𝑖𝑗𝐴^𝑖𝑗−1+𝐶2^𝑖𝑗𝐴^𝑖𝑗+𝐶3^𝑗𝑖1𝐴^𝑖𝑗+1.(12)

2.2. Parabolic Modeling of Random Waves Propagation

It is assumed that the water surface is composed of multiple components of waves. Each wave has angular frequency 𝜔_𝑚 and direction𝜃_𝑛. The refraction, shoaling, and diffraction of discrete wave components are assumed to be governed by the parabolic model of Kirby and Dalrymple [9]. Rewrite (11), the governing equation of complex wave amplitude, by replacing 𝐴 by𝐴_𝑚𝑛. Indices 𝑚 and 𝑛 will be used to represent frequency and direction, respectively [11]𝜕𝐴_𝑚𝑛𝜕𝑥−𝑖𝑘^𝑚−𝑘_𝑚−𝑐𝑐^𝑔𝑘_𝑚_𝑥2𝑐𝑐_𝑔𝑘_𝑚−𝜀𝑐_𝑔𝐴^𝑚𝑛+𝑐𝑐^𝑔_𝑚_𝑦4𝑐𝑐_𝑔𝑘_𝑚𝑖𝑘_𝑚𝑘_𝑚−3−𝑘^𝑚𝑥𝑘^2𝑚−𝑐𝑐^𝑔𝑘_𝑚_𝑥+𝑐𝑐_𝑔2𝑐𝑐_𝑔𝑘²_𝑚_𝑚_𝑦𝑥4𝑐𝑐_𝑔𝑘²_𝑚𝐴^𝑚𝑛_𝑦+𝑐𝑐^𝑔_𝑚_𝑦4𝑐𝑐_𝑔𝑘²_𝑚𝐴_𝑚𝑛_𝑥𝑦+𝑐𝑐^𝑔_𝑚_𝑦4𝑐𝑐_𝑔𝑘𝑐𝑐_𝑔_𝑚𝑖𝑘_𝑚𝑘_𝑚−3−𝑘^𝑚𝑥𝑘^2𝑚−(𝑐𝑐^𝑔𝑘)_𝑚_𝑥+𝑐𝑐_𝑔2𝑐𝑐_𝑔𝑘²_𝑚_𝑚_𝑥4𝑐𝑐_𝑔𝑘²_𝑚𝐴^𝑚𝑛_𝑦𝑦+14𝑘^2𝑚𝐴_𝑚𝑛_𝑦𝑦𝑥=0,(13) where 𝑐_𝑔𝑚=group velocity 𝑐_𝑚 = wave celerity and 𝑘_𝑚= wave number of the 𝑚th component.

The discretization process of the directional spectrum results in wave components of amplitude 𝐴_𝑚𝑛with an associated frequency 𝑓_𝑚and an angle of incidence𝜃_𝑛. The transformed spectrum can be evaluated at any grid point by the superposition of the different wave components. Assuming a Rayleigh distribution of the wave heights and using the computed information about spectral components at locatio{𝑥,𝑦}, the significant wave height can be computed as𝐻_𝑠(𝑥,𝑦)=⎛⎜⎜⎝8^𝑁𝑓_𝑁𝑚=1^𝜃_𝑛=1||𝐴_𝑚𝑛(𝑥,𝑦)||²⎞⎟⎟⎠^1/2,(14) where 𝑁_𝑓 and 𝑁_𝜃= number of discretizations in frequency and direction, respectively.

3. The Elliptic Model

The mild-slope equation (1) has been used to study various kinds of combined refraction, diffraction, shoaling, and reflection phenomena. It was reported that the solutions from the mild-slope equation agreed excellently with the experimental data for waves scattered by a submerged shoal [12]. It was also found that the mild-slope equation could produce accurate solution even the bottom slope is large as 45° [10].

To decrease computational effort, a new numerical method to solve the boundary-value problem of the mild-slope equation is needed. Ebersole [2] presented this method as an efficient way to solve the elliptic form of the mild-slope equation. Balas and Inan [13] used the same technique of Ebersole [2] to solve a field wave transformation problem.

3.1. Derivation of the Elliptic Governing Equation

Following Liu [1], write the wave velocity potential as Φ=Z(𝑧)𝜙=−𝑖𝑔cosh(𝑧+ℎ)cosh𝑘ℎ𝐴𝜔𝑒^𝑖𝑠,(15) where both 𝐴 and 𝑠 are real functions. Substituting (15) into (1) and multiplying the resulting equation by𝐴/𝜔, for monochromatic steady propagation, we can get||∇𝑠||²=𝑘²+∇²𝐴𝐴+∇𝑐𝑐^𝑔.∇𝐴∇.𝑐_𝑔𝑐𝑐_𝑔𝐴,(16)∇𝑠𝑘𝐴𝜔²=0.(17) We can write ∇𝑠=|∇𝑠|𝑠^∧ and 𝑘𝜔=𝑐,then (17) will take the form∇.𝑐𝑐_𝑔||∇𝑠||𝑠^∧𝐴²=0.(18) The wave energy propagates in the 𝑠^∧ direction, and, due to effects of diffraction, the wave energy flux is no longer conservative along wave rays. The curves tangential to the effective wave number vector∇𝑠 may be viewed as the “effective wave rays” for the mild-slope equation; see Liu [1].

We now define 𝜃, the angle of wave propagation as∇𝑠=||∇𝑠||cos𝜃𝑖+sin𝜃𝑗.(19) The eikonal equation and the transport equation can be recast as ||∇𝑠||²=𝑘²+1𝐻𝜕²𝐻𝜕𝑥²+𝜕²𝐻𝜕𝑦²+1𝑐𝑐_𝑔×𝜕𝐻𝜕𝑥𝜕𝑐𝑐^𝑔𝜕𝑥+𝜕𝜕𝐻𝑦𝜕𝑐𝑐^𝑔𝜕𝑦,(20)𝜕𝜕𝑥𝐻²𝑐𝑐_𝑔||∇𝑠||cos𝜃+𝜕𝜕𝑦𝐻²𝑐𝑐_𝑔||∇𝑠||sin𝜃=0,(21) where 𝐻 = wave height = double wave amplitude, 𝐴.

The wave number vector is irrotational, hence𝜕𝜕𝑥||∇𝑠||sin𝜃+𝜕𝜕𝑦||∇𝑠||cos𝜃=0.(22) Equations (20), and (21), and (4) constitute the coupled governing equations for𝐻,|∇𝑠|,𝜃.

3.2. Numerical Solution

A three-step approximate iterative scheme to solve these equations is developed [2](i)First, using the known depth, the linear dispersion relation is solved for the wave number, 𝑘, and then the wave celerity, 𝑐, and the group velocity, 𝑐_𝑔. Snell’s law is formulated as sin𝜃_𝑜𝑐_𝑜=sin𝜃𝑐,(23) where 𝑐_𝑜 = deep water wave celerity = 𝑔/2𝜋; 𝜃_𝑜 = initial angle of incidence. Through Snell’s law, the direction 𝜃 is obtained at all grid points. The refraction and shoaling coefficients are calculated as 𝐾_𝑟=cos𝜃⁰cos𝜃^1/2,𝐾_𝑠=1(1+2𝑘ℎ/sinh2𝑘ℎ)^1/2.(24) An initial approximation for the wave height at grid points is estimated by 𝐻=𝐻_𝑜𝐾_𝑟𝐾_𝑠.(25)(ii)Second, the refraction problem is solved. We can write (22) as 𝜕𝐴𝑠𝜕𝑥+𝜕𝜕𝐴𝑦𝑐=0,(26) where𝐴𝑠=|∇𝑠|sin𝜃;𝐴𝑐=|∇𝑠|cos𝜃. Express the 𝑥-derivative as forward and the 𝑦-derivative as central 𝐴𝑠_𝑖+1,𝑗=𝐴𝑠_𝑖,𝑗+𝑟𝐴𝑐_𝑖,𝑗+1−𝐴𝑐_{𝑖,𝑗−1},(27) where𝑟=Δ𝑥/2Δ𝑦. Let|∇𝑠|=𝑘, from (27) 𝜃 is then obtained as 𝜃_𝑖+1,𝑗=sin⁻¹||∇𝐴𝑠𝑠||_𝑖+1,𝑗.(28) To the refracted wave height, we can write (21) as 𝜕𝐻𝑐𝜕𝑥+𝜕𝜕𝐻𝑦𝑠=0,(29) where𝐻𝑐=𝐻²𝑐𝑐_𝑔|∇𝑠|cos𝜃;𝐻𝑠=𝐻²𝑐𝑐_𝑔|∇𝑠|sin𝜃. Express the 𝑥-derivative as forward and the 𝑦-derivative as central 𝐻𝑐_𝑖+1,𝑗=𝐻𝑐_𝑖,𝑗+𝑟𝐻𝑠_𝑖,𝑗+1−𝐻𝑠_{𝑖,𝑗−1}.(30) The relation|∇𝑠|=𝑘is still used; the wave height is then obtained as 𝐻_𝑖+1,𝑗=𝐻𝑐𝑐𝑐_𝑔||∇𝑠||cos𝜃^{1/2𝑖+1,𝑗}.(31)(iii)To solve for diffraction, (20) is solved to obtain the modified wave number due to diffraction. In (20), express the 𝑥-derivative as backward and the 𝑦-derivatives as central.

Then, we can get||∇𝑠||^2𝑖,𝑗=𝑘^2𝑖,𝑗+1𝐻_𝑖,𝑗𝐻^{𝑖−2,𝑗}−2𝐻_{𝑖−1,𝑗}+𝐻_𝑖,𝑗(Δ𝑥)²+𝐻^𝑖,𝑗+1−2𝐻_𝑖,𝑗+𝐻_𝑖+1,𝑗+1(Δ𝑦)²𝑐𝑐_𝑔_𝑖,𝑗𝐻_𝑖,𝑗×𝐻^𝑖,𝑗−𝐻_{𝑖−1,𝑗}(Δ𝑥)𝑐𝑐^𝑔_𝑖,𝑗−𝑐𝑐_𝑔_{𝑖−1,𝑗}+𝐻(Δ𝑥)_𝑖,𝑗+1−𝐻_{𝑖,𝑗−1}2(Δ𝑦)𝑐𝑐^𝑔_𝑖,𝑗+1−𝑐𝑐_𝑔_{𝑖,𝑗−1}2(Δ𝑦).(32)

Now, the modified wave number is inserted in (27) and solved for angle of propagation. (29) is then resolved for the wave height, and (32) is solved for the modified wave number. This algorithm is repeated until the error in the wave height is less than 1%. This solution method ignores reflection, since it marches toward the shoreline. This procedure takes advantage of the elliptic form of the mild-slope equation of full refraction.

3.3. Elliptic Modeling of Random Waves Propagation

The method adopted, for modeling of random waves propagation, is based on spectral calculation method of Goda [14], which assumes linear behavior between different components of the directional spectrum. The incident directional spectrum is discretized into components. Each component has its own frequency and direction. The solution of each component is carried through the elliptic model of Ebersole [2]. The transformed spectrum can be evaluated at any grid point by the superposition of the different wave components. Assuming a Rayleigh distribution of the wave heights and using the computed information about spectral components at location{𝑥,𝑦}, the significant wave height can be computed as𝐻_𝑠(𝑥,𝑦)=⎛⎜⎜⎝2^𝑁𝑓_𝑁𝑚=1^𝜃_𝑛=1||𝐻_𝑚𝑛(𝑥,𝑦)||²⎞⎟⎟⎠^1/2,(33) where 𝑁_𝑓 and 𝑁_𝜃 = number of discretizations in frequency and direction, respectively.

4. Comparison between the Parabolic and Elliptic Models

To investigate the limits of (PM) and (EM) in modeling the propagation of RWW, both are tested versus the experimental data of Vincent and Briggs [8]. The test includes only nonbreaking series. Nonbreaking series will enable us to compare their abilities in modeling of refraction, diffraction, and shoaling of RWW over a submerged shoal.

4.1. Experiment Setup

The basin is approximately 35 m wide and 29 m long. The basin is flat of uniform depth of 45.72 cm, except a shoal. The shoal center is located at 𝑥 = 6.01 m and 𝑦 = 13.72 m. The elliptical shoal has a major radius of 3.96 m, minor radius of 3.05 m, and height of 30.48 cm at the center. The shoal perimeter is given by3𝑥.0′5²+3𝑦.9′6²=1,(34) where 𝑥^=(𝑥−6.01)m;𝑦′=(𝑦−13.72)m.

The water depth over the shoal is given byℎ=0.9144−.7621−𝑥^23.81+4𝑦.9^25^1/2m.(35) The inputs of the experiment for random waves are in the form of the directional wave spectrum as𝑆(𝑓,𝜃)=𝐺(𝑓)∗𝐷(𝑓,𝜃),(36) where 𝐺(𝑓) represent a one-dimensional frequency spectrum and 𝐷(𝑓,𝜃) is a directional spreading function that satisfies ^02𝜋𝐷(𝑓,𝜃)𝑑𝜃=1.(37)

The irregular wave frequency spectrum produced by the laboratory was TMA [15].

There are two frequency spectra, narrow and broad, which are paired with two different directional spreading, narrow and broad directional spectra narrow and broad.

4.2. Comparison of Results

The results are compared for the nonbreaking series of monochromatic and spectral distribution from case 1 up to case 7 including all varieties of random wave spectra. Test results are all at transect 4 which lies behind the shoal and shows both refraction and diffraction phenomena. For all test cases: wave period is 1.3 sec and representative wave height is 2.54 cm. They have one monochromatic incident wave M2 and two input unidirectional wave frequency spectra (narrow U3 and broad U4). However, they have for input directional wave spectra,

N: narrow directional spreading,

B: broad directional spreading,

3: narrow frequency spreading,

4: broad frequency spreading.

The input spectra for spectral test cases are all combinations of the shapes of Figures 1 and 2. In Figures 3, 4, 5, and 6, we have the output of comparison between the PM and EM as compared to experimental data at Section 4.

Figure 1 

Input frequency spectra.

Figure 2 

Input directional spectra.

(a)

(b)

Figure 3 

Comparison between PM and EM and experimental data for cases (U3 and N3).

(a)

(b)

Figure 4 

Comparison between PM and EM and experimental data for cases (U4 and N4).

(a)

(b)

Figure 5 

Comparison between PM and EM and experimental data for cases (B3 and B4).

Figure 6 

Comparison between PM and EM and experimental data for case (M2).

In general, the numerical models and experimental data show that the difference between the monochromatic wave (M2) and random wave cases (N3, B3, N4, and B4) is dramatic. The pattern associated with the monochromatic wave shows wave height amplification of about 2.5. The wide directional spectral waves (B3 and B4), in contrast, have no amplification greater than 1.2 (which is almost 50% less than the monochromatic case). On the other hand, the narrow directional spread cases (N3 and N4) reach amplification of 1.8 (which is almost 30% less than monochromatic case). The unidirectional cases (U3 and U4) give results close to the monochromatic case.

Comparison of the four directional spectral cases indicates that directional spreading is a more effective parameter than frequency spreading. The patterns of narrow spread cases (N3 and N4) are reasonably similar and those of broad spread cases (B3 and B4) are also reasonably similar.

The unidirectional cases are more like the monochromatic wave than the directional cases. Comparison between the numerical models and the experimental data show that the PM and EM give good results for the monochromatic wave, the unidirectional cases, and the narrow spread cases. In contrast, the EM gives good results for the broad spread cases, while PM does not. Also, the peak relative wave height is almost perfectly reached by the EM but the PM gives slight overestimation in the monochromatic, the narrow spread, and the unidirectional cases and a higher overestimation in the broad spread cases.

Finally, the monochromatic wave representation is good for the unidirectional spectra and adequate for the narrow directional spectra, but it is not accurate for the broad directional spectra.

5. Conclusions

(1) Directional spreading of random water waves is more significant than their energy spreading in the frequency space.

(2) The difference in results between PM and EM depends mainly on the directional spreading of the incident wave spectra. The results are close for monochromatic, unidirectional, and narrow directional spectra, while the results for broad directional spectra are different. EM can cover all incident directions, and hence, its results are more accurate than PM.

(3) The parabolic model is efficient in cases of weak refraction, because the less efficient elliptic model gives the same results.