Edge waves generated by atmospheric pressure disturbances moving along a shoreline on a sloping beach

Coastal Engineering 85 (2014) 43–59

Contents lists available at ScienceDirect

Coastal Engineering journal homepage: www.elsevier.com/locate/coastaleng

Edge waves generated by atmospheric pressure disturbances moving along a shoreline on a sloping beach Seung-Nam Seo a,b, Philip L.-F. Liu c,d a

Coastal Development & Ocean Energy Research Division, Korea Institute of Ocean Science and Technology, Ansan 426-744, Republic of Korea Ocean Science & Technology School, Korea Maritime University, Ansan 426-744, Republic of Korea School of Civil and Environmental Engineering, Cornell University, Ithaca 14853, USA d Institute of Hydrological and Ocean Sciences, National Central University, Jhongli, Taiwan b c

a r t i c l e

i n f o

Article history: Received 2 August 2013 Received in revised form 28 November 2013 Accepted 3 December 2013 Available online 3 January 2014 Keywords: Edge waves generated by moving atmospheric disturbances Method of contour Integration Numerical integration Resonance conditions

a b s t r a c t Edge waves generated by moving atmospheric disturbances parallel to the shoreline are investigated. Following a standard transformation method, an analytical expression of the surface elevation is derived, which consists of an infinite number of modes. Each mode is expressed as the sum of three singular integrals. Using the contour integration method, these singular integrals are converted to regular integrals, which are evaluated by numerical integration methods. The numerical results of two atmospheric pressure distributions studied by Greenspan (1956) are presented, and the resonance conditions are discussed. © 2013 Elsevier B.V. All rights reserved.

1. Introduction It has been reported that moving nearshore geophysical disturbances such as hurricanes and landslides can generate edge waves, which have special features as dispersive waves despite being derived from shallow water equations (Eckart, 1951; Whitham, 1979). These edge waves propagate along the shoreline in a packet form of timevarying finite extent and are essentially confined within a distance from the shore. The crest lines of an edge wave generated by a landslide moving in the cross-shore direction show a complicated pattern (Sammarco and Renzi, 2008; Seo and Liu, 2013). In contrast, edge wave patterns generated by moving atmospheric pressure parallel to the shoreline are much simpler with their crest lines perpendicular to the shoreline (An et al., 2012). For a long wave forced by a moving atmospheric pressure disturbance traveling in a constant depth, unbounded amplification occurs as the result of Proudman resonance (Proudman, 1953). For edge waves generated by a moving atmospheric pressure, resurgent waves typically have amplitudes of the order of 1 meter, wavelengths of a few hundred kilometers, and periods of several hours. Yankovsky (2009) demonstrated that hurricane landfall can generate large-scale edge waves when a storm approaches a coastline at nearly a right angle. Liu et al. (2002) presented an analytical solution of edge waves generated by an atmospheric pressure event in the vicinity of a circular island. Using the model, they identified the peaks in the observed 0378-3839/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.coastaleng.2013.12.002

spectra of sea level around the Belearic Islands in western Mediterranean as the fundamental modes of the edge waves. Monserrat et al. (2006) illustrated that atmospheric disturbances often generate meteotsunamis in the tsunami frequency band, and the wave can be amplified due to various resonance mechanisms. Vennell (2010) showed that when a storm moving slowly at a speed less than the shallow water wave celerity crosses a coast near a critical angle, a subcritical resonance occurs, and generates a large reflected wave traveling along the coast. The Earth's rotational effect on the resonance of the scattered waves generated by storms was investigated by Thiebaut and Vennell (2011). In the present study, we aim to give a plausible description of the amplification factor and resonance conditions for an edge wave generated by moving atmospheric pressure disturbances. Munk et al. (1956) summarized the data collected from four hurricanes that generated edge waves. The observed wave periods were in the range of 5.5 to 8 h, and the observed durations, defined by the total time during which wave motions have been experienced, were given as 16 to 30 h. Hence, the number of waves within a typical edge wave packet is at most 5. Since edge waves are dispersive and generate new crests during propagation, evidently the observed waves appeared rather shortly after wave generation. Accordingly the asymptotic methods adopted by Greenspan (1956), which are valid for a large time, may not provide accurate description to the solution. Greenspan (1956) obtained the surface elevation, consisting of an infinite number of modes, from linearized shallow water equations with

44

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

moving atmospheric pressure forcing. Each mode was expressed as a sum of three singular integrals. One of these integrals had two singular points and produced a quasi-steady wave behind the pressure center. The sum of the other two integrals yielded a complementary wave canceling the quasi-steady wave to make an edge wave in a finite packet form. In the rear part of the packet, the wave amplitudes gradually diminish to zero. To approximately evaluate the integrals for a large time, Greenspan (1956) used the method of steepest descent. When the pole in the integrand is very close to the saddle point, he pointed out that the steepest descent method failed. To address this, Greenspan evaluated only part of the integrand at the saddle point to approximately obtain the complementary wave. Even though Greenspan's solution is only an asymptotically approximate solution, it revealed the main features of edge wave propagation, and was in qualitative agreement with the observed data. In Greenspan's solution, the quasi-steady wave was also represented by an infinite series of sine functions; each of which had a constant amplitude and wavelength along the shoreline, and decayed exponentially in the offshore direction. The complementary wave is not only the most difficult one to evaluate, but also provides a decisive factor to the shape and extent of edge waves in a packet. Surprisingly, Greenspan (1956) approximated waves in a packet as a cut-off quasi-steady wave restricted in the interval from the moving pressure center back to half of its route. Thus, in his solution edge wave has constant amplitude in the packet and an abrupt water level change arises at the end of the packet unless the surface elevation at that point becomes zero. To improve upon this solution for edge waves, a better solution for the complementary wave is needed. For edge waves generated by a landslide moving in the on–off shore direction, the solution for surface elevation is expressed by regular integrals. Using numerical integration methods, Seo and Liu (2013) presented accurate results of surface elevation from the beginning of wave generation to a time of interest. In the present study, we shall employ the same methods to evaluate regular integrals derived from a moving pressure disturbance, which enables us to produce more accurate results on the edge wave packet. Closely following Greenspan (1956), we derive an analytical expression for surface elevation under moving atmospheric pressure disturbances. To convert the obtained singular integrals to regular integrals, the method of contour integration is used. We note that Greenspan (1956) applied the contour integration only to the first integral. While Greenspan used the method of steepest descent to asymptotically evaluate the other two integrals, we employ the contour integration method exclusively to all integrals. We also present a direct way to obtain the complementary wave by introducing an alternative integral effective in the region of nearly zero amplitude behind moving pressure center. Accordingly, our method and its contours are significantly different from Greenspan's approach. Furthermore, since we obtain our solutions by numerical integration, better wave profiles in a packet are obtained with reasonable accuracy and physically acceptable shape from the beginning of wave propagation. Numerical results are presented for two different atmospheric pressure disturbances and discussed to illustrate the main features of this study on edge wave propagation and its resonance. 2. Linear edge wave by a moving atmospheric pressure Cartesian coordinates are chosen with e x pointing out to sea, e y coinciding with a straight shoreline and et denoting time. Linear edge waves on a constant sloping beach can be derived from linear shallow water equations with surface elevation ζe and averaged horizontal veloce over a depth e ity u h (Greenspan, 1956).

∂ζe e  e e h ¼ 0; þ∇ u ∂et

ð2:1Þ

   e ∂u 1 ee  e ζe e Pa e x; e y; et −g ∇ x; e y; et ; ¼− ∇ ρ ∂et

ð2:2Þ

e denotes an atmospheric pressure distribution being applied on where P a sea surface, ρ the density of water and g the gravitational acceleration. To get dimensionless solution of Eqs. (2.1) and (2.2) over a constant sloping beach e hðe xÞ ¼ αe x, physical variables are scaled by the pressure disturbance size a and typical edge wave height ζ0. e x x¼ ; a

e y y¼ ; a

ζ¼

ζe ; ζ0

pffiffiffiffiffiffi gα t ¼ pffiffiffi et; a

Pa ≡

e e P P a ¼ a ; P 0 ρgζ 0

pffiffiffi ζ g e ¼ p0ffiffiffiffiffiffi u: u aα ð2:3Þ

Eliminating the velocity vector in Eqs. (2.1) and (2.2), and substituting the dimensionless variables above into the resulting equation, the non-dimensional equation of surface elevation becomes x

∂2 ζ ∂ζ ∂2 ζ ∂2 ζ þ − þ x ∂x2 ∂x ∂y2 ∂t 2 ! 2 ∂P a ∂ P ∂2 P ¼− þ x 2a þ x 2a ≡ −qðx; y; t Þ: ∂x ∂x ∂y

ð2:4Þ

We seek edge wave solutions of Eq. (2.4) subject to conditions that wave motions diminish at a great distance from the wave generation area and a pressure disturbance is suddenly set in motion at t = 0. Since the solutions are well documented in literatures (Greenspan, 1956; Sammarco and Renzi, 2008; Seo and Liu, 2013; Whitham, 1979), we present, without derivation, the surface elevation by a moving pressure with a constant speed V along the positive y direction. ∞ Z 1 X ζ ðx; y; t Þ ¼ − Re π n¼0

9 8 −jkjx > = < e " Ln ð2jkjxÞQ n ðjkjÞ × #> −ik ðy−VtÞ −i ðk y−ωn t Þ −i ðk yþωn t Þ dk: e e e > > k − þ j j −∞ : ; k2 V 2 −ω2n 2ωn ðkV−ωn Þ 2ωn ðkV þ ωn Þ ∞

ð2:5Þ In Eq. (2.5), only the real part, denoted by Re, has physical meaning. The angular frequency ωn and the forcing function Qn(|k|) are, respectively, defined by 2

ωn ¼ jkjð2n þ 1Þ;

n ¼ 0; 1; 2; ⋯;

ð2:6Þ

and Z Q n ðjkjÞ ¼

∞ 0

e

−jkjσ

^ðσ ; kÞ dσ ; Ln ð2jkjσ Þ q

^ðx; kÞeik Vt ¼ q

Z

∞ −∞

iky

qðx; y; t Þe dy:

ð2:7Þ ^ denotes spatial Here |k| is used for including negative values of k, and q component of the Fourier transform of q shown in Eq. (2.4). When one of the edge wave mode speeds is close to the moving pressure speed, ωn/kn → V, the integrand becomes large but bounded, as will be shown later. Under the condition, resurgent wave occurs as pointed out by Monserrat et al. (2006) and Vennell (2010). In the present study, the same Fourier transform pair as in Greenspan (1956) is used to correctly produce the resurgent waves behind moving pressure center. And following the standard notation (Abramowitz and Stegun, 1972; Mei, 1989), Laguerre polynomial Ln is defined as Ln ðxÞ ¼

n X ð−1Þm n!xm : ðn−mÞ!m!m! m¼0

ð2:8Þ

We note that the definition in Eq. (2.8) is slightly different from that in Greenspan (1956). It can be readily shown that the surface elevation is finite at the shoreline (x = 0) and vanishes away from the shore.

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

For some cases, it is more convenient to express the solution in terms of the velocity potential ϕ, defined by u = ∇ϕ. As presented in Greenspan (1956), the differential equation for ϕ has the same form as Eq. (2.4) except the forcing term. x

∂2 ϕ ∂ϕ ∂2 ϕ ∂2 ϕ ∂P a þ − 2 ¼ ≡ −qðx; y; t Þ: þ x ∂t ∂x2 ∂x ∂y2 ∂t

ð2:9Þ

After solving ϕ from Eq. (2.9) subject to the imposed conditions, surface elevation can be computed by ∂ϕ ζ ¼ −P a − : ∂t

45

pressure has a hump-shaped cross-section in the alongshore direction and a singular point at the center so that we limit the location of the center on land, x0 b 0 as in the previous studies. As shown by Greenspan (1956), the special pressure in Eq. (3.1) produces only one mode so that it gives a significant simplification in analysis. Applying Fourier transform on Eq. (3.1) and substituting the resulting equation into Eq. (2.5) yields the surface elevation as ζ ðx; y; t Þ ¼

1 2

Z 2

∞

e

−jkj ðx−x0 Þ

−∞ −ik ðy−Vt Þ

k e ×4 0 jkj−k0

ð2:10Þ

 pffiffiffiffi   pffiffiffiffi  3 −i k y− jkj t −i k yþ jkj t jkje jkje þ  pffiffiffiffiffiffiffiffiffiffiffiffiffi þ  pffiffiffiffiffiffiffiffiffiffiffiffiffi5dk: 2 jkj−k jkj=k0 2 jkj þ k jkj=k0

ð3:2Þ

The solution of Eq. (2.9) can be written as ∞ 1X π n¼0 (

Z

ϕ¼−

−jkjx

× e

∞ −∞

"

Ln ð2jkjxÞ Q n ðjkjÞ

#)

jkje−ik ðy−VtÞ jkje−i ðk y−ωn t Þ jkje−i ðk yþωn tÞ þ − 2ωn ðkV−ωn Þ 2ωn ðkV þ ωn Þ k2 V 2 −ω2n

dk:

ð2:11Þ Substituting Eq. (2.11) into Eq. (2.10) gives surface elevation for the potential formulation. ∞ i X ζ ¼ −P a þ Re π n¼0

Z

8 −jkjx 9 > < e " Ln ð2jkjxÞQ n ðjkjÞ × = #> −ik ðy−VtÞ −i ðk y−ωn t Þ −i ðk yþωn t Þ dk: kV e e e > > j k j − − −∞ : 2ðkV−ωn Þ 2ðkV þ ωn Þ ; k2 V 2 −ω2n ∞

ð2:12Þ Here Qn(|k|) is defined in Eq. (2.7) with the forcing function q given in Eq. (2.9). The elevation in Eq. (2.12) consists of two parts: static component produced by the direct effect of pressure distribution Pa on the surface and the dynamic component generating edge waves. The integrands in Eq. (2.12) decaying relatively slowly, compared to Eq. (2.5), may give difficulties in evaluating the integral if Qn does not have sufficiently fast decay. Although the expressions for the surface elevation given in Eqs. (2.5) and (2.12) appear to be very different, they must be identical because they are the solutions for the same problem. We will confirm this point for a special case in Appendix B, in which fundamental wave number k0 = 1/V2 is introduced. These solutions in general form can be reduced to the specific ones presented by Greenspan (1956). Surface elevations in both Eqs. (2.5) and (2.12) are expressed by a sum of three singular integrals. To evaluate these singular integrals, we will use the method of contour integration in the complex plane. One crucial point in performing contour integration essentially resides in the representation of |k| on the complex plane. Following Greenspan (1956), function k⁎ in the complex k plane is used and defined by k ¼

pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi k2 þ δ2 :

ð2:13Þ

In order to render this multivalued function k⁎ to be single valued, two branch cuts from ±iδ to ±i∞, respectively, are needed to limit the domain. Hence, k⁎ ≥ 0 on the entire real axis. Whenever it is necessary, we will take the limit of δ → 0 for evaluating the integrals correctly. 3. Hump-shape moving atmospheric pressure by Munk et al. (1956) The atmospheric pressure function used by Munk et al. (1956) was Pa ¼

ðx−x0 Þ ; ðx−x0 Þ2 þ ðy−VtÞ2

ð3:1Þ

which was also studied by Greenspan (1956) with x0 = −1. The pressure center is located at (x0,Vt), moving with a speed V. This special

Eq. (3.2) was also derived by Greenspan (1956) in dimensional form. Hereafter k denotes a complex variable in performing the contour integrations. And the contours are taken such that contributions from a large semi-circle and a small circle around the branch points tend to zero as their limiting processes are performed. Therefore, we only need to consider the contributions from indented semicircles around poles and branch lines, which will be shown shortly. 3.1. The first integral (quasi-steady wave) To evaluate the first integral in Eq. (3.2), we consider a complex integral with the closed contour shown in Fig. 1. −jkj ðx−x0 Þ

∮C e

e−ik ðy−VtÞ dk ¼ 0: k −k0

ð3:3Þ

According to Jordan Lemma, the upper half-circle with a branch cut must be taken for y − Vt b 0 and the lower half-circle for y − Vt N 0. For y = Vt, a special treatment must be employed. As shown in Appendix A, for y − Vt b 0 evaluating the contour integral in Eq. (3.3) gives Z

∞ −∞

−jkjðx−x0 Þ −ikðy−VtÞ

e

e −k ðx−x0 Þ dk ¼ 4π e 0 sin½k0 ðy−VtÞ j kj−k0 Z ∞ i r ðx−x0 Þ r ðy−Vt Þ e e dr: þ2Re r−ik0 0

ð3:4Þ

The first term on the right hand side of Eq. (3.4) is the contribution from both indented semicircles and the second from both edges of the branch cut. For y N Vt, the right contour in Fig. 1 must be used and it can be readily shown that there is no contribution from the indented semicircles. On the segment along both edges of the branch cut, the contribution becomes Z

Z e−jkjðx−x0 Þ e−ikðy−VtÞ dk ¼ 2Re j kj−k0 −∞ ∞

∞

−r ðy−VtÞ

e

0

eir ðx−x0 Þ dr: r−ik0

ð3:5Þ

Greenspan (1956) expressed the integrals in Eqs. (3.4) and (3.5) as confluent hyper-geometric functions, but in this study we evaluate them numerically. For y = Vt, the integrand is symmetric about k = 0 so that the integral can be reduced to Z

∞ −∞

−jkjðx−x0 Þ

e

jkj−k0

Z dk ¼ 2

∞ 0

e

−kðx−x0 Þ

k−k0

dk;

ðx−x0 N0Þ:

ð3:6Þ

46

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

Fig. 1. Contours for the first integral subject to different conditions.

This integral can be evaluated by using formula 3.352.4 of Gradshteyn and Ryzhik (1980).

3.2. The second and third integrals (complementary wave)

Z

In evaluating the second integral in Eq. (3.2), a complex integral is considered with the closed contour shown in the left panel of Fig. 2. As mentioned in Section 1, the present method of contour integration is significantly different from an asymptotic method employed by Greenspan (1956).

∞ 0

e−kðx−x0 Þ dk ¼ − exp½−k0 ðx−x0 ÞEi½k0 ðx−x0 Þ; k−k0

ð3:7Þ

where Ei denotes the exponential integral. Since Ei tends to −∞ as its argument goes to zero, it is necessary that x ≠ x0. This case was omitted in Greenspan (1956), which bridges surface elevations in front of the pressure center to those behind it. More importantly, the computed values at y = Vt from Eq. (3.7) can be used to check accuracy of numerical integration for Eqs. (3.4) and (3.5). In summary, the first integral can be evaluated by using Eqs. (3.4), (3.5) or (3.7). The contribution from the poles, Eq. (3.4), produces the quasi-steady wave, valid only for 0 ≤ y b Vt, and plays a dominant role. We point out again that the contribution from the poles arises only from an indented contour with an upper semicircle. Z k0 ∞ −jkjðx−x0 Þ e−i k ðy−VtÞ e dk 2 8 −∞ jkj−k0 Z ∞ ir ðx−x0 Þ > −k ðx−x0 Þ r ðy−Vt Þ e > > dr; 2π e 0 sin½k0 ðy−Vt Þ þ Re e > > < ð r−ik0 Þ 0 Z i r ð x−x Þ 0 ∞ ¼ k0 −r ðy−VtÞ e > > Re e dr; > > ð r−ik0 Þ > 0 : − exp½−k0 ðx−x0 ÞEi½k0 ðx−x0 Þ;

ζ1 ¼

−k ðx−x0 Þ

∮C e

 pffiffiffiffi  k e−i ky− k t  pffiffiffiffiffiffiffiffiffiffiffiffi dk ¼ 0: k −k k =k0

There is no contribution from the pole in both contours. In evaluating contour integrals along the branch lines in Fig. 2, a branch cut for square root is needed. Here it is taken on the non-negative real axis. As shown in Appendix A, the second integral in Eq. (3.2) becomes

ζ2 ¼

1 4

Z

∞

−jkjðx−x0 Þ

e

−∞

 pffiffiffiffi  Z −i ky− jkjt jkje 1  pffiffiffiffiffiffi  dk ¼ Rei 2 jkj−k jkj V

ybVt

y ¼ Vt

ð3:8Þ Integrals from the branch lines in Eq. (3.8) are regular and the exponential decay term helps to ensure fast convergence. Hence it can be easily computed by a standard method of numerical integration. In Appendix B, we show analytically that surface elevation for the potential formulation Eq. (2.12) gives the same result as in Eq. (3.8) for y = Vt.

∞ 0

−y r−t

e

 pffiffiffiffiffiffi r=2 e pffiffiffiffiffiffiffiffiffiffi iπ=4 dr: 1 þ r=k0 e

pffiffiffiffiffiffi

r=2 −i rðx−x0 Þ−t

ð3:10Þ

:

yN Vt

ð3:9Þ

Because the exponential term in the right hand side decays fast with increasing r, this regular integral is easily computed by a numerical integration method and gives very small contribution to the final result. For the third integral, the following contour integral is considered with the closed contour shown in the right panel of Fig. 2.

−k ðx−x0 Þ

∮C e

 pffiffiffiffi  k e−i kyþ k t  pffiffiffiffiffiffiffiffiffiffiffiffi dk ¼ 0: k þ k k =k0

Fig. 2. Contours for the second integral (left) and the third (right).

ð3:11Þ

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

47

In a similar way to the second integral Eq. (3.9), we get 1 ζ3 ¼ 4

Z

∞

−jkjðx−x0 Þ

e

−∞

 pffiffiffiffi  Z −i kyþ jkjt 1 jkje  pffiffiffiffiffiffi  dk ¼ Rei 2 jkj þ k jkj V

∞

 pffiffiffiffiffiffi r=2 e pffiffiffiffiffiffiffiffiffiffi dr: 1 þ r=k0 e−i3π=4

−y rþt

e

0

pffiffiffiffiffiffi

r=2 −i rðx−x0 Þþt

ð3:12Þ Since the exponential term in the integrand does not show monotonic decay with r, some difficulties in evaluation are readily anticipated which essentially depends on the relative magnitude between y and t in the exponential term. However, the third integral is eventually bounded by self-cancelation due to high oscillations of the integrand as r increases. Consequently, a special treatment is needed to address the difficulty. Evaluation of the second and third integrals is based on a branch cut for the square root taken on the non-negative real axis. As a result, Eq. (3.10) is not the complex conjugate to Eq. (3.12) under this branch cut. In contrast, when a branch cut for the square root is taken on non-positive real axis, we can show that the evaluation of both integrals becomes the complex conjugate to each other, which was also mentioned by Greenspan (1956). Here we reiterate again that Greenspan used the method of steepest descent to obtain the solution approximately. It can be also shown that the sum of two integrals, Eqs. (3.10) and (3.12), does not depend on the choice of branch cut for the square root. In comparison with Eq. (3.10), Eq. (3.12) plays a major role in canceling the quasi-steady pffiffiffiffiffiffiffiffiffiffiffiffi wave. To elucidate its role, the change of variables ξ ¼ r=2k0 is introduced into Eq. (3.12). Then, we have  pffiffiffiffiffiffi Z ∞ r=2 e pffiffiffiffiffiffiffiffiffiffi −i3π=4 Re i dr ¼ −4k0 1 þ r=k0 e 0 0

−k0 ξ ð2y ξ−Vt Þ ξe ξ cosðk0 ξ½2ðx−x0 Þξ þ Vt Þ− dξ: × 2 2 ð1−ξÞ sinðk0 ξ½2ðx−x0 Þξ þ Vt Þ ð1−ξÞ þ ξ Z

∞

−y rþt

pffiffiffiffiffiffi

r=2 −i r ðx−x0 Þþt

e

ð3:13Þ

Fig. 3. Contour for alternative integral.

its radius tends to zero. From the contour integral of Eq. (3.15), we can get Z

∞

−jkjðx−x0 Þ



e

−jkjðx−x0 Þ

−∞

pffiffiffiffi 

Z pffiffiffiffiffiffi dk ¼ Re jkj þ kV jkj

jkje

−i kyþ

jkjt



∞ −∞

e

−jkjðx−x0 Þ

e

pffiffiffiffi 

i kyþ

jkjt

pffiffiffiffiffiffiffiffiffiffi dk: 1 þ k= k0 jkj

ð3:14Þ To evaluate the proposed integral, we will consider a contour integral with the closed contour shown in Fig. 3.

−k ðx−x0 Þ

0 ¼ ∮C e

pffiffiffiffi  jkjt

 pffiffiffiffiffiffi r=2 e −k ðx−x0 Þ pffiffiffiffiffiffiffiffiffiffi −i3π=4 dr ¼ −4k0 πe 0 sin½k0 ðy−Vt Þ 1 þ r=k0 e 0  pffiffiffiffiffiffi pffiffiffiffiffiffi Z ∞ −y r−t r=2 −i r ðx−x0 Þ−t r=2 e e   þRe i dr: ð3:17Þ pffiffiffiffiffiffiffiffiffiffi 0 1 þ r=k0 eiπ=4 ∞

−y rþt

e

pffiffiffiffiffiffi

r=2 −i r ðx−x0 Þþt

Adding Eq. (3.10) and (3.17), finally we have −k ðx−x Þ

To circumvent difficulty in numerical evaluation for the 3rd integral of Eq. (3.12), we consider an equivalent integral. ∞

kyþ

where the first term on the right hand side of Eq. (3.16) comes from an integral along the indented small circle and the second from both branch lines. We used a branch cut for the square root on non-negative real axis. From Eqs. (3.14) and (3.16), Eq. (3.12) can be rewritten as Z

3.3. Alternative integral

Z Re



−k0 jx−x0 j −i k0 ðy−Vt Þ  e pffiffiffiffiffiffiffiffiffiffi dk ¼ −4πk0 i e −∞ 1 þ k= k0 jkj  pffiffiffiffiffiffi pffiffiffiffiffiffi Z ∞ −y r−t r=2 −i r ðx−x0 Þ−t r=2 e e   dr; ð3:16Þ þ2Re i pffiffiffiffiffiffiffiffiffiffi 0 1 þ r=k0 eiπ=4

e

Re i

The major contribution of Eq. (3.13) occurs at small ξ because of the highly oscillatory behavior of the integrand for large ξ. Moreover, the exponential term is greater than or equal to 1 for ξ ≤ Vt/2y. For y ≪ Vt, we have experienced numerical difficulties even though a specialized routine for highly oscillatory integral such as “Qwaf” in Quadpack (Piessens et al., 1983) was used. In the following section an alternative approach is introduced to overcome this problem.

ei

 pffiffiffiffi ei ykþt k  pffiffiffiffiffiffiffiffiffiffi dk: 1 þ k= k0 k

ð3:15Þ

It can be shown that the contribution from the large circle vanishes by Jordan Lemma as its radius goes to infinity and the contribution from the small circle around the branch point also becomes zero as

0 sin½k ðy−VtÞ ζ 2 þ ζ 3 ¼ −2k0 πe 0 pffiffiffiffiffiffi 0  pffiffiffiffiffiffi Z ∞ −yr−t r=2 −i rðx−x0 Þ−t r=2 e e pffiffiffiffiffiffiffiffiffiffi dr: þ Rei 1 þ r=k0 eiπ=4 0

ð3:18Þ

The second term on the right hand side of Eq. (3.18) is twice of the integral in Eq. (3.10) which is very small as mentioned before. Thus, Eq. (3.18) is essentially represented by the first term. Further, the first term cancels the quasi-steady wave shown in Eq. (3.8). As an edge wave packet passes by a location, the surface elevation returns back to the still water level. Thus, the computed surface elevation in the region from the origin to the trailing end of edge wave packet should be nearly the still water level. Since the sum of three integrals, Eqs. (3.8) and (3.18), is negligibly small over the range of y for a large t, the representation for the complementary wave is very effective for this region. So far, Jordan lemma is used to estimate an integral along a large semi-circle. When y = 0, Jordan lemma may not be applied so that this case should be treated separately. For completeness, we present a brief description of the complementary wave at y = 0 in Appendix C, where the wave is composed of both the same canceling term and very small correction of an integral.

48

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

4. Gaussian atmospheric pressure distribution Greenspan (1956) presented another example of Gaussian atmospheric pressure function, which is expressed by Pa ¼ e

−ðx−x0 Þ2 −ðy−VtÞ2

ð4:1Þ

Because of its Fourier transform, the potential function formulation has advantages to obtain surface elevation easily as pointed out by ^ obtained from Fourier transform of Greenspan (1956). Substituting q Eq. (4.1) into Eq. (2.7) gives pffiffiffi −k2 =4 Q n ðjkjÞ ¼ −ikV πe

Z

∞

e

−jkjσ

0

Ln ð2jkjσ Þe

2

−ðσ −x0 Þ

pffiffiffi −k2 =4 An ðjkjÞ dσ ≡ −ikV πe jkj

ð4:2Þ From Eqs. (4.1), (4.2) and (2.12) together with modal wave number defined by kn = k0(2n + 1), surface elevation can be expressed by −ðx−x0 Þ2 −ðy−Vt Þ2

ζ ðx; y; t Þ ¼ −e 1 ζ 1;n ¼ pffiffiffi π

Z

ζ 3;n

∞  X

 ζ 1;n þ ζ 2;n þ ζ 3;n ;

n¼0 ∞

2

−k =4−jkjx

e

−∞

1 ζ 2;n ¼ − pffiffiffi π

þ

Z

∞

Ln ð2jkjxÞAn ðjkjÞ

−k2 =4−jkjx

−ik ðy−VtÞ

jkje jkj−k n

−i k y−Vt

ke

dk; pffiffiffiffiffiffiffi jkjkn

 pffiffiffiffiffiffiffiffiffiffi dk; −∞ 2 k− jkjkn  pffiffiffiffiffiffiffi Z ∞ −i k yþVt jkjkn 2 1 ke −k =4−jkjx  ¼ − pffiffiffi e Ln ð2jkjxÞAn ðjkjÞ pffiffiffiffiffiffiffiffiffiffi dk: π −∞ 2 k þ jkjkn e

Ln ð2jkjxÞAn ðjkjÞ

ð4:3Þ

Z

∞ 0

−ν

e

−ðν=jkj−x0 Þ2

Ln ð2νÞe

dν:

ð4:4Þ

This function depends on the center position x0 only among parameters. When k is real, Eq. (4.4) can be computed by Gauss–Laguerre quadrature. It can be also expressed by a parabolic cylinder function (Gradshteyn and Ryzhik, 1980). For n = 0, using formula 7.4.2 of Abramowitz and Stegun (1972), Eq. (4.4) becomes A0 ðkÞ ¼

 pffiffiffi k k2 =4−kx0 k π e erfc −x0 : 2 2

4.1. Contour integration To evaluate the first integral in Eq. (4.3), we consider an integral in complex plane with the closed contour shown in Fig. 6. −k2 =4−k x

∮C e

ð4:5Þ

Ln ð2k xÞAn ðk Þ

k e−ikðy−VtÞ dk ¼ 0: k −kn

ð4:6Þ

The contours are taken to produce contributions only from the indented semicircles around poles and/or radial lines. This integral must be considered separately for three cases as discussed in Section 3.1. It is found, however, that the case y = Vt can be included into any of other two cases largely due to the strong decay in the exponential term −k2/4. Using similar procedures described in Section 3, the first integral is changed to 2 ζ 1;n ¼ pffiffiffi π

It is noted that Greenspan (1956) did not present the complete solutions for surface elevation as shown in Eq. (4.3). However, his dimensional solution for quasi-steady waves can be deduced from ζ1,n. It can be observed that the exponential part, exp(−k2/4), in Eq. (4.3) helps the integrands to decay quickly, which in turn accelerates fast convergence of the integrals. By changing of variables ν = |k|σ, the function An(|k|) can be rewritten as An ðjkjÞ ¼

The function An(k) with real argument is bounded as k → ∞. As shown in Fig. 4, the function depends on the mode number n. For a small argument and n N 0, the function shows oscillation. The complex valued function An with a complex argument shows a higher level of oscillation for a small argument. Fig. 5 shows its behavior along the radial line with an angle π/4 and the function is also bounded as r → ∞. Because of the enhanced oscillatory behavior for the case of complex argument, Gauss–Laguerre quadrature fails to yield accurate results. The numerical results in Fig. 5 are obtained by routines in Quadpack (Piessens et al., 1983).

( 2π e

−k2n =4−kn x

kn Ln ð2kn xÞAn ðkn Þ sin½kn ðy−VtÞ þ I; I;

y≤ Vt ; y NVt

ð4:7Þ

with Z I ¼ Re

∞

−r

e

ffi

xþjy−Vt j p 2

h

−iprffi x−jy−Vtjþ

e

2

ffi

i

r p 2 2

0

    iπ=4 iπ=4 An r e Ln 2xr e

iπ=4

re dr: r−kn e−iπ=4

ð4:8Þ

It can be shown that the quasi-steady wave in Eq. (4.7) is equal to that presented in dimensional form by Greenspan (1956). To evaluate the regular integral in Eq. (4.7), Gauss–Laguerre quadrature can be used. We have found that it gives accurate results only for locations away from the pressure center y = Vt. Close to the center, more accurate Quadpack routines are needed. The second integral can be written as 

−k2 =4−k x

∮C e

Ln ð2k xÞAn ðk Þ

pffiffiffiffiffiffiffi 

−i ky−

ke

k kn Vt

pffiffiffiffiffiffiffiffiffiffi k− k kn

dk ¼ 0;

ð4:9Þ

and its contour is shown on the left in Fig. 7. As noted before, there is no contribution from the pole in these contours. The branch cut of the

Fig. 4. Function An(k) with real argument for different moving center x0.

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

49

Fig. 5. Function An(reiπ/4) with complex argument for different x0: Real part of function (upper panel) and imaginary part (lower panel).

square root in the regular integrals is taken on the non-negative real axis. Then the integral becomes ζ 2;n

Z 1 ¼ pffiffiffi Re π

∞

−a1 r ib1

e

e

Vt

e

pffiffiffiffiffi

rkn ei5π=8

0

    iπ=4 iπ=4 An r e Ln 2x r e

ð4:10Þ

e−i3π=4 × pffiffiffiffiffiffiffiffiffiffi −iπ=8  dr 1 þ kn =r e

−k2 =4−k x

∮C e

 pffiffiffiffiffiffiffi ke−i kyþVt k kn Ln ð2k xÞAn ðk Þ  pffiffiffiffiffiffiffiffiffiffi dk ¼ 0: k þ k kn

ð4:12Þ

Applying method of contour integration, Eq. (4.12) is changed to Z 1 ζ 3;n ¼ pffiffiffi Re π

with

 xþy r r a1 ¼ pffiffiffi ; b1 ¼ − pffiffiffi x−y þ pffiffiffi : 2 2 2 2

For the third integral, a contour integral with the contour shown on the right in Fig. 7 is considered.

ð4:11Þ

With increasing distance from the origin, the exponential decay term in the integral pffiffiffiffiffiffiffidecreases monotonically. And the component of time exponent Vt rkn cosð5π=8Þb0 also helps fast convergence of the integral. Gauss–Laguerre quadrature can be used to compute Eq. (4.10), which gives very small contribution to the complementary wave.

∞

−a1 r ib1 Vt

e

e e

pffiffiffiffiffi

rkn e−i3π=8

0

    iπ=4 iπ=4 An re Ln 2xre

ð4:13Þ

e−i3π=4 × pffiffiffiffiffiffiffiffiffiffi −iπ=8  dr; 1− kn =r e where parameters a1 and b1 are given in Eq. (4.11). Because the timedependent exponential term in the integrand does not decay monotonically, numerical integration of this integral is non-trivial. Since the integral plays a major role as the complementary wave, a high level of numerical accuracy is required. Quadpack routines are again used to get numerical results of Eq. (4.13). However, they sometimes return

Fig. 6. Contours for y ≤ Vt (left) and y N Vt (right) of the first integral.

50

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

Fig. 7. Contours for the second integral (left) and the third (right).

signs of unreliable level of accuracy. Under this circumstance, a special treatment is devised to produce more accurate result. 4.2. Alternative integral Adopting a similar way used in Section 3.3, we devise a technique to get the complementary wave solution more directly. Noting that the third integral is equal to an alternative integral shown in the right hand side of the following equation, Z Re



pffiffiffiffiffiffiffi

Z ke−i kyþVt jkjkn Ln ð2jkjxÞAn ðjkjÞ  pffiffiffiffiffiffiffiffiffiffi dk ¼ Re −∞ k þ jkjkn  pffiffiffiffiffiffiffi k ei k yþVt jkjkn −k2 =4−jkjx ×e Ln ð2jkjxÞAn ðjkjÞ  pffiffiffiffiffiffiffiffiffiffi dk; k þ jkjkn ∞

−k2 =4−jkjx

e



−k2 =4−k x

∮C e

∞

−k2 =4−jkjx

−∞

¼ 4πe

−2Re

∞

ð4:14Þ

0

pffiffiffiffiffi

e e

rkn e

0

−i3π=8

−i3π=4     e iπ=4 iπ=4  An re Ln 2xre pffiffiffiffiffiffiffiffiffiffi −iπ=8  dr ¼ 1− kn =r e

2

−k =4−k x

n Ln ð2kn xÞAn ðkn Þkn sin½kn ðy−Vt Þ −4πe Z ∞ pffiffiffiffiffi i5π=8     −a1 r ib1 Vt rkn e iπ=4 iπ=4 e e e Ln 2xre þRe An re

−i3π=4

e pffiffiffiffiffiffiffiffiffiffi dr: 1 þ kn =r e−iπ=8

Substituting Eq. (4.17) into Eq. (4.13) and adding it to Eq. (4.10) gives the complementary wave effective in the region of nearly zero wave amplitude.

ζ 2;n þ ζ 3;n

8 2 > −2πe−k =4−kn x Ln ð2kn xÞAn ðkn Þkn sin½kn ðy−VtÞ 2 < Z pffiffiffiffiffi     ∞ ¼ pffiffiffi −a r ib Vt rkn ei5π=8 iπ=4 iπ=4 π> e 1 e 1e Ln 2xre An re : þRe 0

pffiffiffiffiffiffiffi k kn

9 > = : e pffiffiffiffiffiffiffiffiffiffi −iπ=8 dr > ; 1 þ kn =r e −i3π=4

ð4:18Þ

ð4:15Þ

 pffiffiffiffiffiffiffi kei kyþVt jkjkn Ln ð2jkjxÞAn ðjkjÞ  pffiffiffiffiffiffiffiffiffiffi dk; k þ jkjkn

Ln ð2kn xÞAn ðkn Þkn sin½kn ðy−Vt Þ pffiffiffiffiffi     −a r ib Vt rkn ei5π=8 iπ=4 iπ=4 e 1 e 1e Ln 2xr e An r e

−a1 r ib1 Vt

ð4:17Þ

Summing all integrals for surface elevation, Eqs. (4.7) and (4.18), the resurgent wave generated from poles in the first integral is canceled by the wave produced by the pole shown in Fig. 8. For y = 0, the second and third integrals may be treated separately by similar method described in Appendix C. However, it has been observed that the representation in Eq. (4.18) does not give rise to notable numerical problems. 5. Numerical results

−k2n =4−kn x

Z

e

−∞

ke Ln ð2k xÞAn ðk Þ  pffiffiffiffiffiffiffiffiffiffi dk ¼ 0: k þ k kn

e

∞

0

From the shape of indented semicircle, the pole produces a contribution. Evaluating Eq. (4.15) and taking a branch cut for the square root on the non-negative real axis, we have Z Re

Z Re

∞

we consider the following integral with the close contour shown in Fig. 8. i kyþVt

where parameters a1 and b1 are already given in Eq. (4.11). From Eqs. (4.13), (4.14) and (4.16), the major component of complementary wave can be written in another form, which is effective on locations where the wave amplitude tends nearly to zero.

−i3π=4

1þ

e pffiffiffiffiffiffiffiffiffiffi −iπ=8 dr kn =r e

ð4:16Þ

Fig. 8. Contour for the complementary wave of Gaussian pressure distribution.

In the present study, both Gauss–Laguerre quadrature and routines in Quadpack (Piessens et al., 1983) are used for evaluating regular integrals. When an integrand shows a high variation near the lower limit of integration, more quadrature points are needed in this rapidly varying interval. To this end, we employ a modified version of the Gauss– Laguerre quadrature described by Seo and Liu (2013). Another way to improve the accuracy is simply to increase the number of quadrature points. We use both techniques, and at least 150 terms in the quadrature are employed. Under these circumstances, quadruple precision is required for the calculation of the abscissa and weight. Since the Gauss–Laguerre quadrature did not produce a sufficient level of accuracy in many cases, the Quadpack routines were used for most evaluations. The leading-order solution, the quasi-steady wave in Eqs. (3.8) and (4.7), is expressed as a function of k0 and x0. The atmospheric pressure shape determines the forcing term, which largely affects the solution. ffi e is written as V pffiffiffiffiffiffiffiffi Because the dimensional speed V gαa , the bottom slope α and pressure size a also influence the solution. The dimensional

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

51

Fig. 9. Normalized moving speed for different bottom slopes, storm size and speed.

elevation ζe in Eq. (2.3) is given as (P0/ρg)ζ, with P0 being a pressure anomaly. Therefore, the surface elevation is affected by both the bottom slope and atmospheric properties such as shape, disturbance size and strength (a and P0), moving speed, and center position e x0 . The results in this section will be discussed with regard to the effect of these parameters.

For hurricanes observed on the east coast of the United States, the typical moving speed and storm size a were on the order of 60 km/h and 200 km, respectively. Using an average bottom slope α = 5 × 10−4, Greenspan (1956) obtained 3 b k0 b 4, which leads to 0.5 b V b 0.58. Under this condition, he showed that only the fundamental mode of the edge wave could be generated. As shown in Fig. 9, a smaller storm

Fig. 10. Snapshots of edge waves generated by a hump-shaped pressure distribution at different time frames with parameters x0 = −0.1 and V = 1.

52

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

has a higher normalized speed when other parameters remain the same, and the bottom slope affects V significantly. In order to compare the previous results, we illustrate edge waves generated by a moving high atmospheric pressure. Edge wave by a low pressure system has only sign difference compared to that by high pressure system. 5.1. Hump-shaped moving atmospheric pressure Fig. 10 shows the computed surface elevations at different time frames of an edge wave generated by a hump-shaped pressure distribution with parameters x0 = −0.1 and V = 1. Thus, the pressure center is located on land along (− 0.1, t). Early on, at t = 2, we can clearly see that the static effect of the pressure distribution on the sea surface is dominant. As time progresses, an edge wave packet is formed, which results from the interaction between the quasi-steady wave and the complementary wave. The extent of the edge wave packet increases with time, and trailing waves have decreasing amplitude away from the moving pressure center. It is also clearly shown that the largest amplitude appears at the shore, and a quick decrease in amplitude in the offshore direction confines the waves to near the shore. At the top of each panel, the computed extreme values are presented, and the subscript of ζ denotes the number of mode summations starting from 0. In order to examine evolution of the edge wave and its characteristics generated by a hump-shaped moving atmospheric pressure, extreme surface elevations and their positions along the shoreline are plotted in Fig. 11. In the left panel of Fig. 11, the positions of extreme surface elevation signifying the maximum (positive) and minimum (negative) elevations along the shoreline are shown for various moving pressure centers. In the computation, the same wave conditions as in Fig. 10 are used. Early on, the highest elevation in a wave packet appears ahead of the moving center. At slightly less than Vt = 4, the highest elevation shifts to the second wave crest in the packet. For a short time after the shift, the speed of the second wave crest is slightly slower than that of moving pressure. As time goes on, the highest wave speed is recovered to that of the moving pressure. On the other hand, the lowest elevation appears at the trough of the first wave in the packet for a while. After Vt = 23, the lowest elevation occurs at the second wave trough. In the right panel of Fig. 11, extreme surface elevations are plotted against the positions of the moving pressure center. The maximum elevation in a wave packet increases with the distance from the origin. After the maximum elevation attains its largest value at around Vt = 22.5, it begins to decrease. The minimum elevation in a given packet shows a more complicated pattern. It has the first lowest elevation at around Vt = 14.5, after which the minimum elevation increases until shifting to the second wave occurs.

Fig. 12. Surface elevation at the shoreline by a hump-shaped moving pressure with speed V = 1 and x0 = −0.1. The black lines denote the sum of all wave components and the red dotted lines are the quasi-steady wave created by poles in the first integral.

Surface elevation profiles along the shoreline (x = 0) are presented in Fig. 12 generated by a hump-shaped pressure distribution with parameters x0 = − 0.1 and V = 1. The black line denotes the edge wave elevation from all components, and the dotted one is that of the quasi-steady wave. The wave in Eq. (3.8) has an amplitude of 2πk0 exp[− k0(x − x0)], which shows a good approximation of the edge wave in an interval, shown in the dotted lines in Fig. 12, from the pressure center to a point behind it. According to Greenspan (1956), this point is given as 0.5 Vt and from this point to the origin (0 b y b 0.5 Vt) surface elevation returns back to the still water level. From Fig. 12 we can see considerable deviation from Greenspan's approximation. At a later time, it can be observed that most waves in a packet have approximately the same amplitude as the quasi-steady wave. The rear part of a packet which is heavily affected by the complementary wave has a significantly different shape from that of the front waves.

Fig. 11. Positions of extreme surface elevations and extreme values of surface elevation for edge waves generated by a hump-shaped pressure distribution with parameters x0 = −0.1 and V = 1.

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

Fig. 13. Amplitudes of quasi-steady wave at the shoreline for a hump-shaped moving pressure with different speeds and centers: x0 = −0.1, −0.2 and −1.

The complementary wave is mainly computed from Eq. (3.12) by using Quadpack routines. The numerical accuracy in the decreasing wave region becomes worse as time progresses. If a notification of unreliable level of accuracy is returned from the routines, an alternative representation of the complementary wave in Eq. (3.17) is used to replace Eq. (3.12). The elevation of the last wave in the wave packet abruptly drops to the still water level, which is caused by a switch from Eq. (3.12) to Eq. (3.17), and indicates difficulty in computing the complementary wave. In Fig. 13, the amplification factors of the quasi-steady wave at the shoreline are plotted for a hump-shaped moving pressure against varying speeds. In Eq. (3.8), the amplitude of the quasi-steady wave depends on both the moving speed and pressure center location on land. We have conducted several numerical experiments in examining the surface elevation profiles for slow-moving systems. Because of the factor k0 = 1/V2 in Eq. (3.8), slow-moving pressures, which roughly mean a larger system, create bigger waves. When the pressure center is at an inland point farther away from the shoreline, the amplification factor is reduced appreciably. Moreover, when x0 = −1, it can be anticipated that only static depression in the first integral occurs for V ≤ 0.3,

53

in which the amplification factor amounts to 0.001. For a given position x0, the maximum amplitude for a hump-shaped pressure ocpffiffiffiffiffiffiffiffiffiffi curs at V ¼ −x0, which can be obtained by finding zeroes of derivative of the quasi-steady wave amplitude at the shoreline with respect to k0. Thus, a distant moving disturbance from the shoreline produces the maximum amplitude at a higher speed. The case of x0 = − 1 was studied by previous researchers (Greenspan, 1956; Munk et al., 1956), but surface profiles were not presented. In Fig. 14, the surface elevation at the shoreline is plotted for different moving speeds V = 0.3 and 1, respectively, at t = 30. In order to clearly show the difference between the edge wave (black line), the quasi-steady wave (red dotted line) and the contribution from branch cut lines of the first integral (blue dashed line), the quasisteady wave is plotted over the same range even though the quasisteady wave is valid only behind the pressure center. It can be concluded from the result for V = 1 that in the rear of the wave packet, the quasi-steady wave should not be used as an approximation of the edge wave. And the complementary wave affects both the reducing amplitude and changing phase. We can also see that the extent of the wave packet presented by Greenspan (1956), Vt/2 b y b Vt, is only a crude approximation. For V = 0.3, evidently the static depression being plotted in the dashed line is produced by the regular integrals in Eq. (3.8), because the elevation of the quasi-steady wave is drawn in a straight line with nearly zero height. If a moving speed is slower than a critical speed, static depression is dominant with small wave packets developed in front and behind of the pressure center and it moves together with the atmospheric pressure. The critical speed may be defined by a speed to produce the amplitude of the fundamental mode in the quasi-steady wave less than 0.01 in the present scale. In Fig. 14, the static depression for V = 0.3 is virtually the same as the surface drop in this figure resolution by the given pressure distribution, according to the inverse barometric rule described by Munk et al. (1956). The surface elevation is calculated by the analytical expression in Eq. (3.8), but the Quadpack routine is used to evaluate the integrals in the equation. The smooth transition of the surface elevation around the pressure center demonstrates the accuracy of the present numerical methods. To examine evolution of surface profile generated by a slow-moving pressure, we plot the snapshots of surface profiles along the shoreline at different times in Fig. 15. As time progresses, small changes in the wave form can be observed. Behind the moving center, wave packet is not well developed. Changes in the front wave packet are rather perceptible. The static depression without change of form moves together with the atmospheric pressure.

Fig. 14. Surface elevation at the shoreline by a hump-shaped moving pressure with x0 = −1 and different speeds at t = 30: Edge wave (black lines), the quasi-steady wave (red dotted lines) and the other contribution of the first integral (blue dashed line).

54

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

Fig. 15. Snapshots of surface profiles along the shoreline by a hump-shaped atmospheric pressure with a slow-moving speed V = 0.34 and center at an inland line x0 = −1.

5.2. Gaussian pressure A Gaussian pressure distribution is regarded as a better model than a hump-shaped pressure distribution for a realistic atmospheric pressure system. In Eq. (4.3), each integral has a function An(|k|) with real argument, which can be computed by the Gauss–Laguerre quadrature. When the argument of the function is complex, more accurate Quadpack routines should be used. Fig. 16 shows the computed surface profiles for a Gaussian pressure with a very high speed V = 1 and centered at the shoreline, x0 = 0. Although this speed might be regarded as an upper limit of typical hurricane conditions, it is included to examine the effect of the second mode. In this example, the early phase of the propagation is also dominated by the effect of static depression. As time passes, the dynamic effect is reinforced to produce an edge wave packet. The second mode makes the crest strikingly peaked, and the finite extent of the wave packet is clearly seen, which is increasing with time. In Fig. 17 computed at t = 20 with a moving speed V = 1, we demonstrate the effect of a higher mode on the surface elevation at the shoreline. The black line denotes the edge wave from all three integrals, and the red dotted line with symbol △ shows the quasi-steady wave. The red dotted line is drawn where the quasi-steady wave appears as a good approximation to the edge wave. The second mode changes the surface profile significantly, and both the crest and trough of the edge wave become more peaked. At the first trough and crest behind the pressure center, the first and second modes are approximately in phase, as shown in the third panel, which enhances the elevation at the points. Away from the pressure center, the wave phase of the first mode is substantially modified by the complementary wave. Hence, the summation of two modes becomes irregular, as shown in the middle

panel. In this example, the last wave in the packet does not show an abrupt change in elevation, as found in the previous example. Each modal amplitude in the quasi-steady wave, evaluated at the shoreline, is plotted in Fig. 18 against the moving pressure speeds with different center locations. The amplitude has a sign carried over from An(|k|). For a moving system of V = 0.4, the normalized amplitude of the first mode amounts to 0.002. In the present velocity potential formulation, the static depression should be included in this leading-order amplitude, as in Eq. (2.12). Since the static depression is restricted to the neighborhood of the pressure center and the first crest of the quasi-steady wave is generally located away from the center, the contribution from the static depression at the point can be negligible. It can be expected from Fig. 18 that a slow-moving pressure less than or equal to V = 0.4 essentially produces a static depression (or protrusion) with small wave packet in front and behind of the pressure center, which moves together with its disturbance. An et al. (2012) also showed the similar result, as shown in their Fig. 4, using a long wave numerical model. For slow-moving systems with speed less than a critical speed, defined differently by the phase speed of the fundamental mode whose dimensionless wave length is set to 2, they create only a surface depression without any wave pattern. The amplitude of the first mode in the quasi-steady wave attains the highest value of 4.16 at a speed of V = 0.77. For a moving system of V ≤ 0.7, an edge wave may be created exclusively by the first mode, because all of the higher mode amplitudes are essentially negligible. In the neighborhood of V = 1, the second mode cannot be negligible in comparison with the first mode, as shown in Fig. 17. In Fig. 18, the magnitude of the quasi-steady wave amplitude is clearly affected by the moving center position. For V ≤ 1, only two modes are effectively generated. Interestingly, in the graph “Mode Sum,” the moving speed showing the maximum amplitude does not depend on its center position, and as the speed increases, all amplitudes of each mode move towards zero. For a small k0 (or equivalently large V), however, the exponential term of the quasi-steady wave amplitude in Eq. (4.7) becomes almost unity and An(|k|) shown in Fig. 4 has a positive small value, which produces a small positive modal amplitude. For a large V, the amplitude of the first mode is significantly smaller than that of a higher mode. To demonstrate large time behaviors of the edge wave propagation along the shoreline, we plot the surface elevations at different times for various speeds in Fig. 19. The surface profiles are computed for a time longer than the observed hurricane duration by Munk et al. (1956). The wave crests are created during propagation, and the moving pressure speed strongly affects both the wave amplitude and shape in a packet. In the figure, the black line denotes the surface elevation summed by the first two modes. For V b 0.8, the amplitude of the second mode is so small that the surface elevation is essentially formed by the first mode. For the slow-moving pressure in Fig. 19a, the effect of static depression is clearly seen near the moving center located at 0.5 t and the generated wave looks quite irregular compared to the others. Among these plots, the maximum height appears for the case of V = 0.77, and the elevation from x0 = 0.5 is higher than that from x0 = 0. In the computation of a Gaussian pressure distribution by An et al. (2012), e ¼ 50 m/s, and x0 = 0 were used, which is a = 200 km, α = 1/400, V equivalent to V = 0.71. The snapshots of the surface elevation in Fig. 2 of An et al. look very similar to those in Fig. 19b. In the diminishing amplitude zone, the integral for the complementary wave in Eq. (4.13) is switched to its alternative integral in Eq. (4.17), which occasionally brings small but sudden changes in elevation, especially after a long time. It should be noted that in the present method the numerical disturbances neither propagate in space nor in time. Even though they are produced, no notable problems on the surface profile can be found in Fig. 19. Based on the results presented in Figs. 18 and 19, resurgent edge waves are possibly created in the

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

55

Fig. 16. Snapshots of computed surface elevation summed up to the second mode for a Gaussian pressure with speed V = 1 and center at x0 = 0.

range of 0.7 ≤ V ≤ 1, the amplitude of which amounts to approximately 4 or larger. 6. Conclusions Edge waves generated by moving atmospheric pressures with a constant speed parallel to a straight coastline are analyzed. For a constant

sloping beach, an edge wave is composed of infinite modes, each mode of which is expressed by a sum of three singular integrals. To evaluate the integrals, the method of contour integration is applied to convert them to regular integrals. These regular integrals are accurately evaluated by numerical integration methods. The first integral has two singular points which generate a quasisteady wave behind the pressure center and an additional regular

Fig. 17. Surface elevation at the shoreline for a Gaussian pressure with a speed V = 1 and center at χ0 = 0. The black line denotes the edge wave for a given mode at t = 20, and the red lines with symbols is from the quasi-steady wave by two poles of the first integral. The right panel shows elevation summed up for the first two modes, while the left panel is for the first mode only.

56

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

When the wave is computed by the first integral instead of the quasi-steady wave, a better approximation of the edge wave can be achieved. To determine the time duration of an edge wave, the complementary wave must be computed accurately, and the quasi-steady wave approximation by Greenspan (1956) leads to an appreciable error. For a Gaussian atmospheric pressure, regarded as a more realistic model for natural atmospheric pressure distribution, the resurgent edge waves are possibly created in the range of 0.7 ≤ V ≤ 1.0, in which most typical hurricanes occur, and its amplitude amounts to approximately 4 or larger. Under this condition, at most two modes are generated. It can be observed from the expression of the surface elevation that the elevation is affected by both the bottom slope and atmospheric properties such as shape, disturbance size, pressure anomaly, moving speed, and center position. Acknowledgments

Fig. 18. Modal amplitudes of quasi-steady wave evaluated at the shoreline for Gaussian pressure distributions with different moving speeds.

This work was supported in part by the Korea government research project under grant PM57700, and by KIOST projects, PE98916 and PE98976. Appendix A. Contour integral for a hump-shape pressure When δ → 0, a point on the indented semi-circle centered at k0, shown in the left contour of Fig. 1, can be written by

integral. This regular integral makes a smooth transition of the surface elevation in the neighborhood of the pressure center, where the sea surface is also directly deformed by the pressure distribution according to the inverse barometer rule. The sum of the other two integrals, with one singular point in each integral, produces a complementary wave against the quasi-steady wave to generate an edge wave packet of a finite extent behind the pressure center. Once the edge wave packet is generated, the trailing waves occupy roughly the first half of the route of the moving pressure. The amplitude of the rear waves in the packet is gradually diminishing to zero, and the surface elevation from the end of the packet up to the origin returns back to the still water level due to the action of the complementary waves. The numerical evaluation experiences severe difficulties in the region of nearly zero amplitude. To circumvent this, an alternative integral is used to effectively generate the complementary wave in the region of diminishing amplitude, although this integral is not uniformly valid over the whole range. Even though numerical disturbances are produced in the computed elevation caused by using the alternative integral, the computed profiles are shown to have acceptable shape without any notable problems. In numerical experiments, it has clearly been shown that the edge wave propagates in a packet form when the moving pressure speed is faster than a critical speed, which is defined by the speed to produce amplitude of the fundamental mode in the quasisteady wave less than 0.01 in the present scale. For many cases, the quasi-steady wave can give a reasonable approximation of the waves only in the first half of a wave packet. Because of this distinctive feature, the amplitude of waves in the first half of a wave packet is almost constant after an edge wave packet has been fully formed. The edge waves in the present study are generated by a constant forcing on water from a moving atmospheric pressure. In contrast, the edge waves by a landslide movement on an infinite sloping beach are formed by the landslide forcing with a finite duration. When the landmass moves below a certain depth, the effect of it on the free surface motion essentially vanishes. After this, consequently, the amplitude of landslide edge wave decreases with time (Sammarco and Renzi, 2008; Seo and Liu, 2013).

iθ

k ¼ k0 þ εe ;

  2 k ¼ k0 þ ε cosθ þ O ε ; ε≪1 :

ðA:1Þ

The contribution from this small semicircle yields Z

−k ðx−x0 Þ −ikðy−VtÞ

e

lim

ε→0; δ→0

C iþ

e k −k0

−k0 ðx−x0 Þ −ik0 ðy−VtÞ

dk ¼ −2πie

e

:

ðA:2Þ

Similarly, the contribution from the small semicircle centered at − k0 gives the complex conjugate to Eq. (A.2). Summing up the contributions from the indented semicircles gives Z

e

−k ðx−x0 Þ −ikðy−VtÞ

e k −k0

Ci

dk ¼ −4πe

−k0 ðx−x0 Þ

sin½k0 ðy−VtÞ:

ðA:3Þ

On the segment along the right edge of the branch cut, a point can be expressed as k ¼ ir;

k ¼ ir;

ðA:4Þ

as δ→0;

and we have Z lim

R→∞; ε→0

C Sr

Z e−k ðx−x0 Þ e−ikðy−Vt Þ dk ¼ − k −k0

∞ 0

−irðx−x0 Þ

e

erðy−Vt Þ dr: r þ ik0

ðA:5Þ

And on the segment along the left edge of the branch cut, a point is given by k ¼ ir;

k ¼ −ir;

as δ→0:

ðA:6Þ

The integral along this contour becomes the complex conjugate to Eq. (A.5). Applying the residue theorem to Eq. (3.3) for y − Vt b 0, Eq. (3.4) can be obtained.

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

57

Fig. 19. Surface elevations along the shoreline at different times by Gaussian pressure with different moving speeds: The red dotted lines denote surface elevations by the first mode and the black lines by the first two modes.

For the contribution from the right edge of branch pffiffiffiffiffiffiffiffiffi cut of contour integral Eq. (3.9), a branch cut for square root, −i r , is taken on the non-negative real axis. On the right edge of branch cut, we have Z

−k ðx−x0 Þ

k e−i



yk−t

pffiffiffiffi k

pffiffiffiffiffiffiffiffiffiffiffiffi dk ¼ lim i R→∞; k −k k =k0 ε→0 pffiffiffiffiffiffi  pffiffiffiffiffiffi Z ∞ −y r−t r=2 i rðx−x0 Þ−t r=2 e e pffiffiffiffiffiffiffiffiffiffi dr: ¼i 1 þ r=k0 e−iπ=4 0

lim

R→∞; ε→0

e

C Sr

Z

R

e ε

i r ðx−x0 Þ

pffiffiffiffiffiffi ð−irÞe−y rþit −ir pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi dr ð−irÞ þ ir ð−irÞ=k0

ðA:7Þ

And it can be shown that contribution on the left edge of branch cut becomes complex conjugate to Eq. (A.7).

Appendix B. Surface elevation at the pressure center from velocity potential representation In order to show that the two different expressions, Eqs. (2.5) and (2.12), for a surface elevation at the pressure center are indeed identical, it is convenient to use a moving coordinate system (x′ = x, y′ = y − Vt, t′ = t) as in Greenspan (1956). Fourier transform of q in Eq. (2.9) becomes ^ðx′; kÞ ¼ V q

Z

∞ −∞

2 2 ∂P a iky′ −ðx′−xo Þ pffiffiffi −k =4 πe : e dy′ ¼ −ikVe ∂y′

ðB:1Þ

58

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59

Fig. 20. Contours for the second integral (left) and the third (right) at y = 0.

Contribution from the branch lines gives

Substituting Eq. (B.1) into Eq. (2.7) also leads to only one mode. Q 0 ðjkjÞ ¼ −

ikVπ jkjx0 e : 2jkj

ðB:2Þ

From Eqs. (B.2) and (2.12), surface elevation on the potential formulation can be written as Z ðx−x0 Þ 1 ∞ −jkj ðx−x0 Þ þ e ðx−x0 Þ2 þ ðy−VtÞ2 2 −∞ 2 3 −ik ðy−Vt Þ ke−i ðk y−ω0 tÞ ke−i ðk yþω0 tÞ 5 jkje  dk: ×4 þ pffiffiffiffiffiffiffiffiffiffi  − pffiffiffiffiffiffiffiffiffiffi jkj−k0 2 2 jkjk −k jkjk þ k

ðB:3Þ

ζ ¼−

0

pffiffiffiffi pffiffiffiffiffiffi  pffiffiffiffiffiffi Z ∞ −t r=2 −i rðx−x0 Þ−t r=2 k e−k ðx−x0 Þ eit k e e pffiffiffiffiffi dk ¼ −i pffiffiffiffiffiffiffiffiffiffi lim R→∞; k −kV k 1− r=k0 eiπ=4 CS 0 ε→0 pffiffiffiffiffiffi  pffiffiffiffiffiffi Z ∞ −t r=2 −i rðx−x0 Þ−t r=2 e e pffiffiffiffiffiffiffiffiffiffi × dr−i dr: 1 þ r=k0 eiπ=4 0 Z

From all contributions in Eq. (C.1), we get 1 4

Z

0

For y = Vt and x − x0 N 0, the sum of the static depression and the first integral, with the help of Eq. (3.7), becomes Z ∞ −kðx−x0 Þ Z − k ðx−x0 Þ 1 1 ∞ e jj 1 ke jkj þ þ ζ1 ¼ − dk ¼ − dk ðx−x0 Þ 2 −∞ jkj−k0 ðx−x0 Þ k−k0 0 Z ∞ −kðx−x0 Þ e dk ¼ −k0 exp½−k0 ðx−x0 ÞEi½k0 ðx−x0 Þ: ¼ k0 0 k−k0

ðB:4Þ Comparing Eqs. (B.4) and (3.8), it is evident that two different expressions for a surface elevation are the same.

For this case, we consider a contour integral that originated from the second integral Eq. (3.2) −k ðx−x0 Þ

0 ¼ ∮C e

pffiffiffiffi

it

pffiffiffiffiffi dk; k −kV k

−k ðx−x0 Þ

lim

ε→0

e Ci

it

k e

k −kV

pffiffiffiffiffi dk ¼ 4πi k0 e k

0

e−t

r=2 −i rðx−x0 Þ−t

pffiffiffiffiffiffi

e

pffiffiffiffiffiffiffiffiffiffi 1 þ r=k0 eiπ=4

r=2

9 > > = : > dr > ;

For the third integral in Eq. (3.2), we consider an integral. −k ðx−x0 Þ

0 ¼ ∮C e

k e−it

pffiffiffiffi k

pffiffiffiffiffi dk: k þ kV k

ðC:5Þ

The indented semicircle and branch lines contribute value. The regular integrals obtained from the contour integration becomes Z

pffiffiffiffi e−jkjðx−x0 Þ e−i jkjt 1  pffiffiffiffiffiffi  dk ¼ 4 −∞ jkj þ k jkj V 8 pffiffiffiffi −k0 ðx−x0 Þ −it k0 > > < 4πik0 e pffiffiffiffiffiffi e pffiffiffiffiffiffi Z Z × ∞ e−t r=2 ei r ðx−x0 Þ−t r=2 > > pffiffiffiffiffiffiffiffiffiffi −iπ=4 dr−i : −i 1 þ r=k0 e 0 ∞

9 > > pffiffiffiffiffiffi  pffiffiffiffiffiffi = : ∞ e−t r=2 ei rðx−x0 Þ−t r=2 > pffiffiffiffiffiffiffiffiffiffi −iπ=4 dr > ; 1− r=k e 0 0

In summary, adding Eqs. (C.4) and (C.6) gives the complementary wave for y = 0 which can be written, after simplifying the result, as −k ðx−x0 Þ

e 0 ζ 2 þ ζ 3 ¼ 2πk08
∞ 0

e−t

sinðk VtÞ pffiffiffiffiffiffi 0  pffiffiffiffiffiffi r=2 −i r ðx−x0 Þ−t r=2 e drg: 1−iðr=k0 Þ

ðC:7Þ

References

pffiffiffiffi k

∞



pffiffiffiffiffiffi

ðC:1Þ

and its contour is shown in Fig. 20 The contours are chosen such that contributions come from pole and branch cut lines. In evaluating regular integrals obtained from contour integration, a branch cut for the square root is taken on the non-positive real axis, which produces a decreasing value as time increases. The second integral on the indented semicircle becomes Z

−jkjðx−x0 Þ i

e 

ðC:6Þ

k

k e

pffiffiffiffi jkjt e 1 pffiffiffiffiffiffi  dk ¼ 4 −∞ jkj−k jkj V 8 pffiffiffiffi −k0 ðx−x0 Þ it k0 > > þ < −4πik0 epffiffiffiffiffiffi  e pffiffiffiffiffiffi Z ∞ −t r=2 −i rðx−x0 Þ−t r=2 Z × e e > > p ffiffiffiffiffiffiffiffiffi ffi i dr þ i : iπ=4 1− r=k0 e 0 ∞

ðC:4Þ

1 4

Appendix C. Complementary wave at y = 0 for a hump-shape pressure

ðC:3Þ

−k0 ðx−x0 Þ it

e

pffiffiffiffi k0

:

ðC:2Þ

Abramowitz, M., Stegun, I.A., 1972. Handbook of Mathematical Functions. Dover, New York. An, C., Liu, P.L.-F., Seo, S.N., 2012. Large-scale edge waves generated by a moving atmospheric pressure. Theor. Appl. Mech. Lett. 2, 042001.

S.-N. Seo, P.L.-F. Liu / Coastal Engineering 85 (2014) 43–59 Eckart, C., 1951. Surface waves in water of variable depth. Tech. Rep. 100. Scripps Inst. Ocean. Wave Report. Gradshteyn, I.S., Ryzhik, I.M., 1980. Table of Integrals, Series and Products. Academic Press, Inc., Orlando. Greenspan, H.P., 1956. The generation of edge waves by moving pressure distributions. J. Fluid Mech. 1, 574–590. Liu, P.L.-F., Monserrat, S., Marcos, M., 2002. Analytical simulation of edge waves observed around Balearic Islands. Geophys. Res. Lett. 29 (17), 1–4 (28). Mei, C.C., 1989. The Applied Dynamics of Ocean Surface Waves. World Scientific Publishing Co., Singapore. Monserrat, S., Vilibic, I., Rabinovich, A.B., 2006. Meteotsunamis: atmospherically induced destructive ocean waves in the tsunami frequency band. Nat. Hazards Earth Syst. Sci. 6, 1035–1051. Munk, W., Snodgrass, F., Carrier, G., 1956. Science 123, 127–132.

59

Piessens, R., De Doncker, E., Uberhuber, C., Kahaner, D., 1983. Quadpack —A Subroutine Package Automatic Integration. Springer-Verlag, Berlin. Proudman, J., 1953. Dynamical Oceanography. John Wiley & Sons, Inc., New York. Sammarco, P., Renzi, E., 2008. Landslide tsunami propagating along a plane beach. J. Fluid Mech. 598, 107–119. Seo, S.N., Liu, P.L.-F., 2013. Edge waves generated by the landslide on a sloping beach. Coast. Eng. 73, 133–150. Thiebaut, S., Vennell, R., 2011. Resonance and long waves generated by storms obliquely crossing shelf topography in a rotating ocean. J. Fluid Mech. 682, 261–288. Vennell, R., 2010. Resonance and trapping of topographic transient ocean waves generated by a moving atmospheric disturbance. J. Fluid Mech. 650, 427–442. Whitham, G.B., 1979. Lectures on Wave Propagation. Springer-Verlag, New York. Yankovsky, A.E., 2009. Large-scale edge waves generated by hurricane landfall. J. Geophys. Res. 114, C03014.