Abstract

Coupled-mode equations describing the propagation and scattering of internal waves over large-amplitude arbitrary topography in a two-dimensional stratified fluid are derived. They consist of a simple set of ordinary differential equations describing the evolution of modal amplitudes, based on an orthogonality condition that allows one to distinguish leftward- and rightward-propagating modes. The coupling terms expressing exchange of energy between modes are given in an analytical form using perturbation theory. This allows the derivation of a weak-topography approximate solution, generalizing previous linear solutions for a barotropic forcing that were described in 2002 by Llewellyn Smith and Young . In addition, the orthogonality condition derived is valid for a different set of eigenmodes defined on a sloping bottom, which shows a better convergence rate when compared with the standard set of modes. The work presented here provides a useful and simple framework for the investigation of internal wave propagation in an inhomogeneous ocean, along with theoretical insight.

1. Introduction

Ocean internal waves and internal tides have been identified as important constituents of the oceanic circulation (e.g., Wunsch and Ferrari 2004). They play an important role in the dissipation of energy and the mixing of water masses, with a major impact, as yet not fully understood, on the meridional overturning circulation (MOC; Kunze 2017). For instance, internal tides generated by barotropic tidal currents over topographic features are thought to account for 1 TW of energy dissipation in the deep ocean (Egbert and Ray 2001) of the 2.1 TW required to support the MOC (Munk and Wunsch 1998; Wunsch and Ferrari 2004).

One of the difficulties in understanding tidal-driven mixing and dissipation lies in their propagation, thus transporting and redistributing energy over large distances, with large variations in the distance traveled and hence large uncertainties in the localization of the ultimate breaking events, if any. Low-mode (modes 1 and 2) internal waves can propagate over very long distances, while higher modes generally dissipate and break locally. In addition, internal tides interact with their environment (topography, varying stratification and mesoscale currents), resulting in scattering and reflection, thus adding variability to their propagation (Garrett and Kunze 2007).

The conversion of the barotropic tide has been investigated by many authors (Cox and Sandstrom 1962; Baines 1973; Bell 1975a,b), and theoretical estimates, based on a small-topographic-amplitude approximation among others, have been formulated by Llewellyn Smith and Young (2002) and St. Laurent and Garrett (2002). These have been used to elaborate global maps of the estimated internal tide generation (Nycander 2005; Melet et al. 2013; Falahat et al. 2014), with reasonable accuracy at least in the deep ocean (Garrett and Kunze 2007). The subsequent propagation and scattering of the internal tide are less well documented, especially with regard to the theoretical aspects. The propagation and diffraction of internal wave beams were investigated by Rainville and Pinkel (2006), and the scattering was studied—again in the limit of weak topography—by Müller and Xu (1992), but questions remain, particularly in regions of large-amplitude topography such as coastal shelves (Kelly et al. 2012; Klymak et al. 2016).

Theoretical models addressing propagation for finite-amplitude topography can be classified into two groups, often formulated in a two-dimensional domain. The first uses the method of characteristics and conformal mapping to build solutions for isolated topography (Müller and Liu 2000; Bühler and Holmes-Cerfon 2011; Maas 2011). It is applicable to horizontally homogeneous stratification and a limited class of vertical profiles, although these can be general enough to apply to the ocean. Another approach is based on a vertical mode decomposition of the solution, which may be convenient because these modes propagate independently in a horizontally homogeneous medium over a flat bottom, with known horizontal wavenumbers and group velocities. In addition, such a decomposition is particularly relevant in the context of surface (barotropic)-internal-tide decomposition, and the generation of internal tides, where modal decomposition is a natural choice and has been used either for theoretical predictions (Llewellyn Smith and Young 2002) or for diagnosing observations or numerical simulations (e.g., Zhao et al. 2019). It was used by Griffiths and Grimshaw (2007) to investigate the generation of internal tide at continental shelves, a formulation reused later on by Griffiths (2011) to propose a global model of the internal tides. Maugé and Gerkema (2008) used a similar formulation in a more complex (and three-dimensional) domain, including background flow and weakly nonlinear effects, yet under the approximation of small topographic height. Kelly et al. (2012) derived an energy equation to diagnose the internal tide scattering on a continental slope, and Shimizu (2011) proposed a slightly different formulation based on the linearized shallow water equations with multiple layers. An intermediate method, introduced by Pétrélis et al. (2006) in the context of the barotropic tide conversion and applied more recently to the scattering of internal tide (Mathur et al. 2014), is based on representing the solution using a Green’s function and solving the resulting integral equation along the bottom boundary.

Using a vertical-mode decomposition in a horizontally inhomogeneous medium results in nonconstant modal amplitudes as well as nonzero horizontal derivatives of the eigenfunctions. The latter have to be evaluated numerically to be able to solve the equations that determine the former. Kelly et al. (2013) introduced a semianalytical two-dimensional model that confines the coupling at the interface of horizontal steps with homogeneous medium properties (equivalent depth), by matching the solutions between adjacent steps. This approach is a generalization of earlier work by Sjöberg and Stigebrandt (1992) where each step’s contribution was added independently (thus discarding backward reflection and interferences). Within this formulation, horizontal derivatives of the modes are replaced by integrals representing the product of mode functions between adjacent steps.

In this paper, we present a method based on modal decomposition that relies on perturbation theory to formulate a set of “coupled-mode equations” describing the evolution of modal amplitudes. In this approach, quantities F(x, y, z) are represented in the form

 
F(x,y,z)=mAm(x,y)Fm(z;x,y),
(1)

where the modes Fm are defined by one-dimensional eigenproblems with “frozen” properties; that is, horizontal variables enter only parametrically so that the eigenproblems are local with respect to x and y. The quantities Am(x) are range-dependent modal amplitudes the horizontal evolution of which is obtained. This method is known as the “reference waveguide method” (Brekhovskikh and Godin 1992) or “coupled-mode theory” (Desaubies and Dysthe 1995, and references therein). The result is a set of coupled-mode equations which contains coupling coefficients accounting for exchanges between modes. Obtaining the coupled-mode equations requires extensive integration by parts to account for the possibly different boundary conditions satisfied by F and Fm, as will be shown.

The general approach works in three dimensions with variable bathymetry, rotation, compressibility, stratification, and background velocity fields. In this paper we limit ourselves to the two-dimensional case for an incompressible ocean with rotation and vertical stratification N(z) so that the background state is at rest, so as to be able to focus on the details of the method and compare with previous results. We derive a set of eigenmodes and corresponding orthogonality condition that make explicit the distinction between counterpropagating modes, allowing us to project the initial linearized equation onto the basis of vertical modes in a concise way. In addition, the orthogonality condition enables us to use a different set of vertical modes that match the bottom boundary condition with finite-amplitude topography, thereby exhibiting more rapid convergence with n compared to the standard modes.

The use of perturbation theory allows one to recast the horizontal derivatives of modes as a function of the values of the modes themselves and the gradient of the topography. This manipulation is exact provided the topography is a differentiable function of x. Using a WKBJ approximation for the vertical modes, we can further derive analytical expression for the coupling between modes for both amplitude and energy flux. Injecting these coupling terms in the coupled-mode equations, using a weak-topography approximation, leads to analytical solution for the scattering of an incoming wave with arbitrary modal constituents.

The paper is organized as follows. Section 2 presents the linearized equations and the vertical modes used for projection. Some of the properties of these modes that are useful for the subsequent calculations, in particular results from perturbation theory, are described. Then, in section 3, the coupled-mode equations are derived and discussed in the context of energy exchange between modes. Section 4 examines two examples of internal wave scattering over topography, including solutions based on the coupled-mode equations. In section 5, we introduce a different set of vertical modes that have the boundary slope embedded in the bottom boundary condition, and discuss changes to the results that follow from this new set. Conclusions are given in section 6, and the  appendix contains some calculations whose results are needed in the main text.

2. Problem formulation and normal mode decomposition

a. Governing equations and boundary conditions

We consider the linearized equations for a stratified fluid in rotation under the Boussinesq approximation. We restrict ourselves to perturbations that are independent of one of the horizontal coordinates, thus considering a two-dimensional domain bounded by a rigid lid at z = 0 and bottom at z = −H(x). The equations are formulated in the frequency domain:

 
iωv+fez×v+p=bez,
(2)
 
iωb+N2w=0,
(3)
 
xu+zw=0,and
(4)

where N(z) is the buoyancy frequency and standard notation is used for the variables (e.g., Vallis 2006). The following boundary conditions hold:

 
w=0atz=0andudHdx+w=0atz=H(x).
(5)

The horizontal boundary conditions depend on the details of the configuration investigated, and will be specified later on.

The problem can be recast in terms of the velocities u and w, the variables present in the vertical boundary conditions. The equations become

 
zu+αxw=0and
(6)
 
xu+zw=0,
(7)

with α(z) = (N2ω2)/(ω2f2) being the square of the inverse characteristic slope, which has no horizontal variations in the problem considered.

We introduce the “state vector” Ψ = (u, w)T, which allows us to rewrite the problem in the following matrix form:

 
Ψx=LΨz,
(8)

with

 
L=(01α10).
(9)

b. Vertical normal modes

To apply a modal decomposition to the linear problem (8), we first need to define a complete set of modes. We consider the range-independent linear equations [i.e., with frozen properties with respect to x so that H(x) = H is treated as a constant] and apply a Fourier transform along the x direction. We obtain the following eigenvalue problem, with ϕn and φn being the vertical functions for u and w, respectively:

 
ϕn=iknαφnand
(10)
 
φn=iknϕn,
(11)

with eigenvalue kn. Here and in the rest of the paper, prime denotes a vertical derivative. This set of equations is usually formulated for one variable or the other; that is,

 
(ϕnα)+kn2ϕn=0andφn+kn2αφn=0.
(12)

As above we introduce the state vector for the normal modes Φn = (ϕn, φn)T that will be used below.

The boundary conditions corresponding to the rigid lid and flat bottom are φ(0) = φ(−H) = 0, implying ϕ′(0) = ϕ′(−H) = 0. With this choice of boundary conditions, both equations in (12) lead to an eigenproblem with real eigenvalue kn2. Then the eigenfunctions ϕn and φn are associated with the following orthogonality conditions:

 
H0ϕnϕmdz=0.5jn(δn,m+δn,m)and
(13)
 
H0αφnφmdz=0.5jn(δn,m+δn,m),
(14)

where jn is a real arbitrary normalization constant and δn,m is the usual Kronecker symbol. For simplicity, and without loss of generality, we will set jn = 1 in the rest of the paper.

The right-hand sides in the orthogonality conditions (13) and (14) reflect the fact that counterpropagating modes with identical vertical structure are degenerate. However, one can take a combination of the two conditions to distinguish between the modes. We thus define modes corresponding to the positive and negative roots by

 
k±n=±kn2.

With this notation, positive and negative roots are associated with modes propagating toward increasing and decreasing x, respectively, and a change of sign of φ. Unlike the usual convention, “mode” in this paper will thus refer to the state vector Φn; that is, counter-propagating modes have different structures, although they are composed of similar functions (ϕ±n, φ±n). A “generalized” orthogonality condition that makes the distinction between right-going and left-going modes can be formed by taking H0[(ϕnφm)(ϕmφn)]dz. Using (10) and (11), this gives

 
(ϕnφmϕmφn)H0=iH0[(knkm)ϕnϕmα(knkm)φnφm]dz.
(15)

By virtue of the boundary conditions, the LHS vanishes. Anticipating further calculation, we introduce the matrix

 
P=(100α)
(16)

and the notation H0dz for conciseness. The orthogonality condition is thus

 
ΦnPΦm=H0(ϕnϕmαφnφm)dz=δn,m,
(17)

where Φn=(ϕn*,φn*)=(ϕn,φn) is the conjugate transpose of Φn = (ϕn, φn)T, since choosing jn > 0 implies that the vertical modes ϕn are real while φn are imaginary.

c. Perturbation analysis and normal mode properties

When using modal analysis in an inhomogeneous domain, one uses locally defined normal modes to project the equations. This procedure results in terms that contain horizontal derivatives of the eigenmodes (see the next section and the  appendix), for which no analytical expression exists in general, thus requiring numerical evaluation. Here, we use perturbation theory (Courant and Hilbert 1966) to express these derivatives as a function of the eigenmodes themselves and the horizontal derivative of inhomogeneities of the medium. Interestingly, this method has been routinely used in the framework of ocean and atmospheric acoustics (e.g., Brekhovskikh and Godin 1992), but has not been used for internal wave propagation previously, to our knowledge.

The perturbation theory is applied as follows: consider the eigenproblem (10) and (11), in two media with slightly different properties: Hn = H, Hm = H + δH and αn = α, αm = α + δα. At this stage, the physical origin of the perturbations (e.g., horizontal inhomogeneities) does not matter. We have introduced here a perturbation in α for generality and further computation, although we will not use it to include variations in N2 because these are not addressed in this paper and would require a different set of linear equations. Perturbations of α may instead be caused by variations in ω, and this is used below to obtain an expression for the modal group velocity. As a result of these perturbations in α, we have kmkm + δkm, ϕmϕm + δϕm, and φmφm + δφm. Taking the first variation of (15) gives

 
(ϕnδφmδϕmφn)H0=iδkmδn,m+i(knkm)×H0(ϕnδϕmαφnδφm)dz+ikmH0φnφmδαdz.
(18)

The LHS vanishes at the upper boundary but not at the bottom. Expanding the lower boundary condition φm(−Hm) = 0 leads to δφm(H)=φm(H)δH=ikmϕm(H)δH. Hence

 
(ϕnδφmδϕmφn)|H=ϕnϕm|H×(ikmδH).
(19)

Since N2 does not vary with x in this paper, we first take δα = 0 and obtain the following equation:

 
ϕnϕm|HkmδH=δkmδn,m(knkm)×H0(ϕnδϕmαφnδφm)dz.
(20)

We can introduce the physical origin of the perturbation [here, δH = (dH/dx)δx] and take the limit of infinitely small perturbation to construct the horizontal derivative of the eigenmodes. This calculation is exact provided that H(x) is continuous and differentiable. Setting n = m in (20) gives the following expression for the horizontal derivative of the eigenvalues:

 
knx=ϕn2|HkndHdx,
(21)

Furthermore, combining (17) and (20) for nm leads to an expression for the horizontal derivative of the eigenvectors:

 
ΦnPΦmx=ϕnϕmknkm|HkmdHdx,nm.
(22)

Note that, while we restrict ourselves to variations of the topography along a single dimension here, the use of perturbation theory in 3D configurations (although using a different orthogonality condition, since right-/left-propagating mode separation would not hold) and for variations of other physical quantities, such as properties of the medium that would enter the eigenproblem as parameters, is straightforward.

This method also allows one to compute the modal group velocity: for this purpose, we repeat the above procedure with variable ω (and, thus, variable α) and δH = 0, to obtain

 
Cgnωkn=(ω2f2)2ωkn(α+1)φn2dz,
(23)

where the result

 
δα=(1+α)2ωω2f2δω
(24)

that comes from the definition of α has been used. Under the hydrostatic approximation, the integral in the denominator of (23) becomes αφn2dz=1/2 and one recovers the hydrostatic modal group velocity (Zhao 2017).

3. Coupled-mode equations

a. Derivation of the coupled-mode equations

We now derive the coupled-mode equations that describe the evolution of the modal amplitude. To do so, we assume that normal modes actually depend on the horizontal coordinates in a parametric way (e.g., Desaubies and Dysthe 1995; Griffiths and Grimshaw 2007) by expanding the variables as

 
Ψ(x,z)=mAm(x)Φm(z;x),
(25)

which is the appropriate version of (1). At each location, Φm are defined by the eigenproblem (10) and (11) with frozen (i.e., the local value of) H. The initial system of equations in (8) is projected onto the basis of vertical modes, with the help of the equalities for the horizontal variations of these modes [(21) and (22)] obtained as a result of perturbation theory developed above. This method is known as the “reference waveguide method” (Brekhovskikh and Godin 1992).

Projection of the equations onto mode n is carried out by multiplying (8) by ΦnP and integrating over [−H(x), 0]. Since Ψ and the Φms do not obey the same condition at z = −H, one must integrate by parts and insert the modal expansion at an intermediate stage of the procedure (e.g., Llewellyn Smith and Young 2002; Kelly 2016). The procedure is detailed in the  appendix. One finally obtains a set of first-order ordinary differential equations (ODEs), the so-called coupled-mode equations, in the form

 
dAndxiknAn=mgnmAm.
(26)

This system of equations describes the evolution of the modal amplitudes over range. The right hand side terms account for energy exchange between modes expressed by the coupling coefficients gnm. These are given by

 
gnm={ϕn22|HdHdx,n=mknknkmϕnϕm|HdHdx,nm.
(27)

These expressions require the use of the results in the  appendix along with (22). With the normalization used in the paper, ϕnϕm scales like H−1 and gnm scales like H−1dH/dx.

The barotropic mode, with n = 0, is unaffected by the internal modes, since k0 = 0 because of the rigid lid condition. Thus, it can act only as a forcing generating internal modes, while g00 ensures mass conservation, that is, that ⟨u0⟩ is constant with respect to x (noting that ϕ0=1/H and φ0 = 0).

The set of equations in (26) combined with the normal mode definition (10) and (11) and the coupling coefficients (27) are a formally exact reformulation of the initial problem (8) using modal decomposition, provided that the depth profile H(x) is continuous and differentiable. It has the advantage of reducing the horizontal dependence of the solution to a simple system of ODEs where the distinction between leftward and rightward propagating mode appears explicitly. A similar approach was used in Griffiths and Grimshaw (2007) and in Kelly et al. (2016) in a more general setting (time-domain formulation, and varying stratification in a three-dimensional domain and background currents), in the form of coupled shallow-water equations. A closer formulation in a two-dimensional domain is the coupling equations for linear tides (CELT) model (Kelly et al. 2013), in which the horizontal domain is divided into sections with horizontally constant properties (N, H, and f) and the solution is constructed by applying matching conditions at the interfaces.

The coupled-mode equations formulation derived in this paper is different in several ways: first, it is continuous with respect to x. Therefore, one cannot derive the coupling terms for discontinuous topography (since perturbation theory does not apply exactly). However, the ODEs in (26) can be solved separately on each side of the discontinuity and the full solution can be obtained by applying dynamical conditions (on the pressure and the velocity) at the interface, following, for example, Robinson (1969). When doing so, the coupled-mode equations yield a linear matrix problem, akin to stepwise formulations as in the CELT model.

In addition, the horizontal derivatives of the eigenfunctions and eigenvalues are recast in terms of the variations of H, thus providing an analytical expression for the mode coupling. Note that this was partially done in the derivation by Griffiths and Grimshaw (2007), in a formulation that was valid for variations of H only (through the chain rule), whereas perturbation theory is applicable to variations of other properties [e.g., N(z, x), f] as well. We mention, however, that this could imply a different set of initial linear equations and additional complexity: for instance, nonhomogeneous stratification is associated with geostrophic currents and additional interaction (advective) terms.

b. Analytical expression for mode coupling and energy exchange

We now give expressions for the horizontal energy flux and the energy exchange between modes based on the coupled-mode equations derived previously. The horizontal flux of energy I is given by

 
I=p*u+pu*=2Re(p*u).
(28)

To obtain an expression for the modal energy flux, we now derive an expression for the modal projection of the pressure. From the vertical component of (2) and (3), we have

 
pz=ω2f2iωαw.
(29)

The procedure is now similar to other derivations given in this paper. To summarize, we multiply both sides of the above equality by φn, expand w=mAmφm, use the second orthogonality condition (14), and integrate the left-hand side by parts, expanding the pressure as p=mpmϕm, with boundary terms vanishing. Then, using (10) and the first orthogonality condition (13), we obtain

 
pn+pn=ω2f2ωkn(AnAn).
(30)

Note that leftward and rightward modes are degenerate in the orthogonality conditions (13) and (14) and are thus combined in (30). However, inserting the modal decomposition for pressure and velocity from (28) into (30), we obtain (after reordering the terms in the sum and using kn = −kn)

 
I=ω2f2ωn>0(|An|2kn+|An|2kn),
(31)

in which one recognizes the contribution of fluxes associated with left- (n < 0) and right-going (n > 0) waves. From this expression we can define the modal energy flux for mode n as

 
In=ω2f2ω|An|2kn.
(32)

Note that the barotropic energy flux is infinite according to this expression, which is a well-known limitation of the rigid-lid approximation. As mentioned previously, however, the barotropic mode only satisfies mass conservation and is an adiabatic reservoir of energy flux for the internal tides (i.e., forcing) in this limit.

From (26) and (21), the modal energy flux amplitude in mode n satisfies the following equation:

 
dIndx=2ω2f2ωmngnmknRe(An*Am)=2ω2f2ωdHdxmnϕnϕm|HknkmRe(An*Am).
(33)

From this expression, one sees that the exchange of energy between modes depends on the variations of the depth of fluid and on the values of the normal modes at the bottom (because gnm does). Moreover, the (knkm)−1 dependence shows that adjacent and copropagating modes are more prone to exchange energy with each other than counterpropagating modes.

We can obtain more useful information about the strength of the coupling (and, thereby, mode scattering) by using analytical expressions for the eigenmodes. We thus use the WKBJ approximation in the eigenproblem (10) and (11), which gives (e.g., Gill 1982)

 
ϕn=βHcos[Hzmn(z)dz]andφn=iβαHsin[Hzmn(z)dz],
(34)

where mn=knα and β(z)=Hα/α. The eigenmode is associated with the eigenvalue

 
kn=nπα.
(35)

We recall, for clarity, that kn, mn, and β are also functions of x (whereas α is not).

Injecting this form into the expression for the coupling coefficient (27), we obtain

 
gnm=β(H)HdHdx×{1/2,n=mnnm,nm.
(36)

This expression can be substituted into (33) to provide an analytical expression for the modal energy exchange:

 
dIndx=2ω2f2ωdHdxmnα(H)π(nm)Re(An*Am).
(37)

The expressions (36) and (37) show that the strength of coupling (in range) scales linearly with the relative topographic variations and the inverse difference between mode numbers.

4. Solutions to the coupled-mode equations

The set of coupled equations in (26), even when analytical expressions such as (36) are available from exact (or approximate WKBJ) solutions for the vertical modes, has no general analytical solution. In this section, we first derive approximate solutions valid in the limit of weak topographic coupling. Then numerical solutions of the coupled-mode equations are presented. These are used to illustrate the method presented in this paper and also to provide a reference for evaluating the accuracy of the approximate solution. We should emphasize that none of the solutions (approximate or numerical) here relies on using WKBJ approximation for the vertical modes. The horizontal (i.e., coupled-mode equations) and vertical (normal mode) solutions are distinct problems. The WKBJ approximation is only necessary if a fully analytical solution is desired (this will be discussed explicitly in what follows).

We first recast the coupled-mode equations by using a “phase compensated” amplitude:

 
Bn(x)=An(x)exp[iθn(x)],
(38)

where we have introduced the phase function θn(x)=xkn(x)dx. Inserting this into (26) gives

 
dBndx=mgnmBmexp(iθnm),
(39)

where θnm=x[kn(x)km(x)]dx. Hence Bn remains constant if gnm = 0; that is, changes in the variables B are associated with scattering. In the following, we will solve this new set of ODEs instead of (26).

Specifying the horizontal boundary conditions leads to a boundary-value problem. One can specify the value of all modes, that is, incoming and outgoing at a physical boundary, or implement a radiation condition by specifying Bn = 0 at the entrance boundary for each mode: the left and right boundary for n > 0 and n < 0, respectively.

a. Approximate solutions for weak coupling

As shown in section 3b, the strength of the coupling scales like the relative topographic variations and the inverse mode number difference. This can be used to derive approximate solutions when the coupling is expected to remain weak.

1) Adiabatic approximation

The so-called adiabatic approximation describes the evolution of modes when coupling terms, and thereby mode scattering, are neglected (Milder 1969; Brekhovskikh and Godin 1992; Desaubies and Dysthe 1995). In this case, energy flux conservation is ensured by retaining the diagonal terms gnn. The solution to (39) is then

 
Bn(0)=Bn(x0)exp[x0xgnn(x)dx],
(40)

where the exponential term accounts for the normalization of the vertical modes to ensure energy flux conservation. This solution is equivalent to the WKBJ solution of the depth-separated linear equations for baroclinic motions (in two dimensions) of Rainville and Pinkel (2006). Using in addition the WKBJ approximation for the vertical modes introduced above, (34), and assuming that β(H) is nearly constant (this is exactly true for stratification profiles with α = [a exp(Hz) + b]2, with a and b constant), the adiabatic solution can be further simplified to

 
Bn(0)(x)Bn(x0)[H(x0)H(x)]β(H)/2.
(41)

Under the adiabatic approximation, the modal amplitudes adjust to variations of the depth of fluid so that the energy flux is conserved mode-wise. This solution is not of great interest in the context of oceanic internal wave scattering, since scattering is explicitly neglected.

2) Weak topography approximation

Let us introduce the following asymptotic expansion: Bn=jεjBn(j), where the small parameter ε is a nondimensional measure of the topographic amplitude δH/H. Since gnm scales with H−1dH/dx, the coupled-mode equations in (39) then become

 
dBn(j)dx=mgnmBn(j1)exp[iθnm(x)],
(42)

where the zeroth-order amplitude Bn(0) is constant and set by the boundary conditions. Note that the previous adiabatic solution can be seen as a zeroth-order approximation, in which the diagonal terms in gnn are explicitly retained.

Integrating (42) to first order in the expansion gives

 
Bn(1)=mBm(0)x0xgnmexp(iθnm)dx.
(43)

To go further, we can also expand the eigenvalues and eigenfunctions in powers of ε, and retain only their leading-order, x-independent part [e.g., if HH0 + εδH(x), then ϕ, φ, and k are function of H0 only at leading order]. The first-order solution then becomes

 
Bn(1)=Bn(0)ϕn022|H[H(x)H(x0)]mnBm(0)k0nk0nk0mϕ0nϕ0m|Hx0xdHdxei(k0nk0m)xdx,
(44)

where the 0 subscript refers to quantities evaluated at H = H0. This expression provides an approximate semianalytical solution for the mode amplitude evolution over arbitrary topography. The first term corresponds to the adiabatic correction to the modal amplitude, while the second accounts for energy transfer between modes. The solution obtained is a generalization (in a one-dimensional horizontal domain) of the linear, weak topography solution of Llewellyn Smith and Young (2002) for the conversion of the barotropic tide. Here, any incident barotropic or baroclinic mode with amplitude Bn(0) can act as a “forcing” tide, while the amplitude of the scattered waves Bn(1) are first order. It is worth mentioning that this first-order solution is accurate as long as the amplitude of Bn(1) remains small, so that the small-parameter expansion remains valid. Since the small parameter is εδH/H, this depends on the horizontal extension of the topographic irregularities (as the amplitude of the first-order solution can grow and accumulate over large distances) in addition to their magnitude. Using the WKBJ approximation for the vertical modes, as in section 1, it is possible to derive analytical solutions for a variety of profiles H(x) and boundary conditions.

b. Numerical solution: Method and example

For large-amplitude topography, the approximate solution derived above is likely to be inaccurate, and one is left with the full boundary value problem (39) to be solved numerically. To illustrate the method described in this paper, we consider an idealized configuration in which a first baroclinic mode, coming from the left, encounters a ridge. The latter is modeled by a cosine bump with

 
H(x)=H0ΔH2[1+cos(πxx0LH)]
(45)

for |xx0| < Lh and H(x) = H0 otherwise, where H0 is the depth away from the ridge and 2LH the ridge width. We use two different values of the ridge amplitude while keeping the far-field depth H0 = 3 and slope criticality max(dH/dx)(α)1/2 = 0.8 constant. The width of the topography changes accordingly. We take α = 1 (the results can be rescaled to realistic values by changing the aspect ratio of the domain and the amplitude of the velocities accordingly). A numerical solution is obtained and compared with the approximate solution derived above. It will be further compared with an exact 2D numerical solution in the next section.

Instead of solving the BVP (39) directly, we reduce the numerical complexity by turning the problem into a set of iterative initial value problems, using the fact that coupling terms between counterpropagating modes are weaker than for copropagating modes. The following set of equations is obtained:

 
dBn(i)dx={m0gnmBm(i)exp(iθnm)m<0gnmBm(i1)exp(iθnm),n>0m0gnmBm(i)exp(iθnm)m>0gnmBm(i1)exp(iθnm),n<0.
(46)

At each iteration, the above set of coupled ODEs is solved in both directions separately, with the contribution from counterpropagating modes taken from the previous iteration.

Boundary conditions are specified on the zeroth-order amplitude Bn(0) while Bn(i>0)=0 is imposed at the “entrance” boundary (i.e., left and right boundary for positive and negative n, respectively). Note that the barotropic amplitude is given by the adiabatic solution and thus B0(i>0)=0, because it is not altered by the baroclinic modes as shown in section 3. The numerical solution is computed using a finite number of modes N in each direction and a finite number of iterations. When reconstructing the physical solution, we only use the first 3N/4 modes because the highest modes resolved are inaccurate: indeed, they should exchange energy with higher modes, but these latter modes are not present due to truncation, so that the corresponding coupling terms are missing. This is similar to the aliasing effect due to nonlinearities in pseudospectral methods. Since coupling between modes is most important for closest modes, the lower modes are only weakly affected by this truncation effect.

The solution to (46) is computed using 50 vertical modes (in each direction) for a “moderate amplitude” ridge ΔH/H = 1/6 with B1(0)=1 on the left boundary, thus forcing the incident first baroclinic mode. Only a few iterations (two here) are necessary to obtain an accurate approximation to the solution, which is shown in Fig. 1. One sees the formation of beams emanating from the topography in both directions: 14% of the incoming flux (contained in the first baroclinic mode) is forward-scattered into higher copropagating modes, while 5% is reflected (i.e., backscattered). As can be seen in the right panel of Fig. 1, the weak-topography solution (44) provides a reasonable approximation to the numerical solution, predicting overall 15% of energy flux forward scattered and 4% reflected. This agreement is lost for higher values of the amplitude of the topography, as expected. The solution for ΔH/H = 1/2 is shown in Fig. 2. In this case, much more of the incoming energy flux is forward scattered into higher copropagating modes (80%) which thus carry most of the energy at the right boundary (nearly 4 times as much as the remaining first baroclinic mode). In contrast the amplitude of the reflected beam is much weaker and accounts for only 0.4% of the incoming energy flux. By comparison, the weak-topography approximate solution predicts 125% and 0.01% forward scattered and backscattered, respectively. It greatly overestimates the generation of copropagating low-order modes (mode 2) while failing to predict large-enough amplitudes for the higher counterpropagating modes.

Fig. 1.

Numerical solution for a moderate amplitude of topography ΔH/H = 1/6: (left) Reconstructed vertical velocity field w=nAnφn (imaginary part, i.e., w at t = π/ω). Absolute value of the amplitude |An| of modes propagating (top right) to the right and (top center) to the left, for the first 10 modes. Note that horizontal extent of the domain is zoomed over the bump and that the values of the y axis differ for positive vs negative modes. (bottom right) Norm of the modal amplitudes at the outgoing boundaries (i.e., left and right boundary for negative and positive modes, respectively; black bars) compared with the weak-topography approximate solution (gray bars).

Fig. 1.

Numerical solution for a moderate amplitude of topography ΔH/H = 1/6: (left) Reconstructed vertical velocity field w=nAnφn (imaginary part, i.e., w at t = π/ω). Absolute value of the amplitude |An| of modes propagating (top right) to the right and (top center) to the left, for the first 10 modes. Note that horizontal extent of the domain is zoomed over the bump and that the values of the y axis differ for positive vs negative modes. (bottom right) Norm of the modal amplitudes at the outgoing boundaries (i.e., left and right boundary for negative and positive modes, respectively; black bars) compared with the weak-topography approximate solution (gray bars).

Fig. 2.

As in Fig. 1, but for a large amplitude of topography ΔH/H = 1/2 (note that the different slope of the beam is due to the different horizontal length of the topography and the size of the domain shown).

Fig. 2.

As in Fig. 1, but for a large amplitude of topography ΔH/H = 1/2 (note that the different slope of the beam is due to the different horizontal length of the topography and the size of the domain shown).

5. Normal modes over a sloping bottom

a. Definition of the sloping-bottom modes

The set of normal modes ϕn and φn, defined in section 2b and used in the subsequent derivation and solution of the coupled-mode equations, forms a complete basis that is uniformly convergent over [−H, 0] (Jackson 1914, his theorem 1). Because the boundary conditions are different from those of the solution to (8), the convergence rate of the basis scales as 1/N, where N is the total number of modes. In this section, we introduce a different set of vertical modes with modified boundary conditions that incorporate the bottom slope. This decomposition can be used in place of the previous one, although the coupling coefficients are different. We discuss the implication and potential utility of this different decomposition, its range of validity and the corresponding modal coupling coefficients.

We consider the previous eigenproblem (10) and (11) with a modified bottom boundary condition:

 
γϕn+φn=0,z=H.
(47)

Here, one can recognize the “natural” boundary condition (5) by identifying γ with the local value of dH/dx. The corresponding eigenproblem is solved either numerically using a collocation method with Chebyshev polynomials or analytically for particular stratifications (e.g., constant α in the examples treated in this paper). The corresponding solution is shown in Fig. 3, compared with the standard, flat-bottom modes defined in section 2b. Note that in the limit of vanishing slope (γ → 0), both solutions are identical. Despite the different boundary condition, the bracket in the LHS of (15) still vanishes, so the solutions still obey the orthogonality condition (17) and define a complete basis suitable for modal expansion. However, the modified eigenproblem is not in standard Sturm–Liouville form. In particular, the eigenmodes as well as the eigenvalues now take complex values.

Fig. 3.

Comparison of flat-bottom (blue) and sloping-bottom vertical modes with γ = 0.2 (orange) and 0.8 (red) for mode numbers (left) 1 and (right) 3. The real and imaginary parts are represented by continuous and dashed lines, respectively. Note that here α = 1 so that the slope is critical for γ = 1.

Fig. 3.

Comparison of flat-bottom (blue) and sloping-bottom vertical modes with γ = 0.2 (orange) and 0.8 (red) for mode numbers (left) 1 and (right) 3. The real and imaginary parts are represented by continuous and dashed lines, respectively. Note that here α = 1 so that the slope is critical for γ = 1.

b. Properties of the modal expansion

This new basis may be useful for several reasons, which are discussed here, among with some noteworthy properties of the “sloping bottom” normal modes. First, the convergence of the corresponding series is expected to be more rapid with n than for the “flat bottom” set, because the boundary condition at the bottom is the same for the basis function and the solution. Figure 4 shows the rate of convergence computed on the basis of an exact solution in a wedge with α = 1 and two different slopes. The solution (shown in Fig. 4, upper-left panel) corresponds to a first baroclinic “wedge” mode propagating upslope in the apex (Wunsch 1969):

 
Ψwedge={exp[iqlog(α1xz)]/(α1xz)+exp[iqlog(α1x+z)]/(α1x+z)exp[iqlog(α1xz)]/(xαz)exp[iqlog(α1x+z)]/(x+αz)},
(48)

with q = 2/log[(α−1 + γ)/(α−1γ)] and γ being the bottom slope. The convergence of the series is measured as the normalized absolute value of the error for u and w given a number of modes N used in the projection (Fig. 4, upper-right panel). (The latter is performed in a standard way, numerically evaluating An=ΦnPΨwedge using Clenshaw–Curtis quadrature.) For the flat-bottom mode decomposition, w exhibits the highest error and lowest convergence rates: slightly slower that 1/N. On the contrary, the sloping-bottom mode decomposition exhibits higher convergence rate, close to 1/N2, especially for highest bottom slope, and a significant decrease of the error associated with w. This is further illustrated in Fig. 4 (lower-left panel) by the magnitude of the modal amplitude resulting from the projection on both type of eigenmodes, which shows that higher modes have much weaker amplitude for the sloping-bottom modes than for the flat-bottom modes.

Fig. 4.

Wedge solution and reconstruction using flat-bottom and sloping-bottom modes: Real part (in phase) of the (top left) horizontal (u) and (middle left) vertical (w) solution for an upslope first baroclinic wedge mode [(48), with α = 1 and γ = 0.5]. Realistic values can be restored by rescaling x and γ by α1/2 and α−1/2, respectively, and the aspect ratio in the figures is such that critical slope is a 1/1 (rise/distance) oblique line in the plots. (top right) Normalized absolute error |0Nanϕnudz|/|u|dz (and similar for w) as a function of the number of modes N used in the projection for flat-bottom (dashed) and sloping-bottom (continuous) modes, for two values of the bottom slope γ = 0.2 and 0.8. Also shown are the projection and reconstruction of the wedge solution using N = 10 modes. (bottom left) Norm of the modal amplitude for the flat-bottom (blue) and tilt-bottom (orange) mode projection, and (bottom center) u and (bottom right) w reconstructed fields (same convention) superimposed to the wedge mode (black). Note the error of the flat-bottom modal reconstruction near the bottom, especially on the vertical velocity (here shown at x = 4). These results are insensitive to the horizontal location, because the wedge solutions are self-similar. Note that in the bottom panels the wedge mode (black) is mostly obscured by the other solutions.

Fig. 4.

Wedge solution and reconstruction using flat-bottom and sloping-bottom modes: Real part (in phase) of the (top left) horizontal (u) and (middle left) vertical (w) solution for an upslope first baroclinic wedge mode [(48), with α = 1 and γ = 0.5]. Realistic values can be restored by rescaling x and γ by α1/2 and α−1/2, respectively, and the aspect ratio in the figures is such that critical slope is a 1/1 (rise/distance) oblique line in the plots. (top right) Normalized absolute error |0Nanϕnudz|/|u|dz (and similar for w) as a function of the number of modes N used in the projection for flat-bottom (dashed) and sloping-bottom (continuous) modes, for two values of the bottom slope γ = 0.2 and 0.8. Also shown are the projection and reconstruction of the wedge solution using N = 10 modes. (bottom left) Norm of the modal amplitude for the flat-bottom (blue) and tilt-bottom (orange) mode projection, and (bottom center) u and (bottom right) w reconstructed fields (same convention) superimposed to the wedge mode (black). Note the error of the flat-bottom modal reconstruction near the bottom, especially on the vertical velocity (here shown at x = 4). These results are insensitive to the horizontal location, because the wedge solutions are self-similar. Note that in the bottom panels the wedge mode (black) is mostly obscured by the other solutions.

While the orthogonality condition (17) is common to the two sets of basis function, it is worth noting that (13) and (14) are not verified by modes on a sloping bottom. Instead we have

 
H0ϕnϕmdz=H0αφnφmdz=iγϕnϕmkn+km|H,nm.
(49)

This means that the vertically integrated modal energy density (or horizontal flux), given by similar integrals, has no meaning within this decomposition. In other words, sloping-bottom modes are intrinsically coupled with each other, in the sense that one cannot define an energy per mode unless dH/dx = 0 (i.e., the bottom is flat).

To gain further insight in the structure of these modes, we use the WKBJ approximation. Then, incorporating the boundary condition at the top, the solution takes the form

 
φn=m(z)1/2sin[Hzm(z)dz],

with mkα, and ϕn follows from (11). Note that this solution is exact for some stratification profiles (Baines 1973). The quantities m and k take discrete value as a result of the boundary condition at the bottom, which becomes

 
tan(kH0αdz)=iγα(H).
(50)

This equality shows that the eigenvalue with the sloping bottom condition is equal to the flat-bottom eigenvalue plus an additional imaginary part. This imaginary part is associated with enhanced amplitude of the modal structure near the bottom (cf. Figs. 3 and 4, lower-right panels). This behavior is consistent with the physics of internal wave beam reflection over a slope, which is known to exhibit amplification of the energy density near the slope as the slope get higher, in particular when the angle becomes close to criticality.

A major caveat of the sloping-bottom boundary condition is that it can give rise to singularities. While it is beyond the scope of this manuscript to derive a complete spectral theory for the eigenproblem with sloping-bottom boundary condition, this singular behavior can be discussed using the above WKBJ solution. The relation (50) has an inverse hyperbolic tangent singularity [with ℑ(m) → ∞] for α(H)γ=1, which coincides with the well-known criticality condition of internal waves on a slope. This is associated with increasing confinement of the eigenmode structure, with very large amplitude in the vicinity of the bottom (due to the increasing value of the imaginary part of m) for α(H)γ1. Away from criticality, the WKBJ solution is well defined. This criticality issue is absent when using the standard modes, as the bottom slope does not enter their definition.

c. Coupling coefficients and coupled-mode equations for sloping-bottom modes

We now list the counterparts to the formulas established in the text in the presence of a bottom slope. The perturbation analysis now includes a change in the slope, γn = γ and γm = γ + δγ. We deduce from (47)

 
δφmφmδH+γδϕmγϕmδH+ϕmδγ=0
(51)

at z =−H. Some algebra and use of (10), (11), and (47) then lead to

 
ϕnδφmδϕmφn=ϕnϕm[ikm(1αγ2)δH+δγ],z=H.
(52)

Then (20)(22) become

 
ϕnϕm|H[km(1αγ2)δHiδγ]=δkmδn,m(knkm)H0(ϕnδϕmαφnδφm)dz,
(53)

and

 
knx=ϕn2|H[kn(1αγ2)dHdxidγdx],
(54)
 
ΦnPΦmx=ϕnϕmknkm|H[km(1αγ2)dHdxidγdx],nm.
(55)

The projection of the linear system of (8) onto the new vertical modes is carried out in the same fashion as over a flat bottom, considering a frozen value of γ = dH/dx in addition to H itself. In this case, however, the boundary conditions for expansion functions and independent variables are the same, so there is no need to use integration by parts. The counterpart to (A3) is

 
ΦnPΦnx=12{dHdx(ΦnPΦn)H}=dHdx(1αγ2)ϕn22|H.
(56)

After some calculation, one obtains the same coupled-mode equations that are in (26) as before, but with the following coupling coefficients:

 
gnm={ϕn22|H[1α(dHdx)2]dHdx,n=m1knkmϕnϕm|H{km[1α(dHdx)2]dHdxid2Hdx2},nm.
(57)

d. Comparison of the solutions using flat-bottom and sloping-bottom modes

A comparison of the solution to the coupled-mode equations and the direct two-dimensional reconstructed fields is shown in Fig. 5. It is further compared to a numerical solution of the 2D initial problem (8) using a collocation method with Chebyshev polynomials in both directions (with 80 collocations points in the vertical and 60 in the horizontal), after recasting the problem in a rectangular domain using a stretched vertical coordinate. The solution is computed with the same topography (with ΔH/H = 1/6) and in the same domain as in Fig. 1 but, for simplicity, we consider instead an initial value problem in which a first baroclinic mode is coming from the left and no waves exit the domain through the left boundary (i.e., An = 0, n ≠ 1). This implies that weak-amplitude (left-propagating) modes are entering the domain at the right boundary and vanish above the topography through interactions with each other and right-propagating modes. Both methods yield solutions very similar to the numerical solution. As previously shown in Fig. 4, the solution above the topography is better reproduced with sloping-bottom mode decomposition, especially near the bottom. Yet, the modal amplitudes at the right boundary (where dH/dx = 0) are very similar with both types of eigenmodes and compare very well with the amplitudes computed from the projection of the numerical solution onto the set of eigenmodes, up to n ≈ ±20. Thus, while using sloping-bottom modes has the advantage of better reconstructing the dynamical fields above nonflat topography, especially near the bottom, they do not provide significant improvement on the modal amplitudes.

Fig. 5.

Comparison of the resolution of the coupled-mode equation using flat-bottom and sloping-bottom modes with a full 2D numerical calculation: Real part of (top) u and (middle) w for the (left) 2D numerical solution and the modal solution using (center) flat-bottom and (right) sloping-bottom modes. Also shown are modal amplitude (absolute value) at the (bottom left) right boundary and vertical profiles of (bottom center) u and (bottom right) w (real value) at the location of maximum slope on the left side of the ridge.

Fig. 5.

Comparison of the resolution of the coupled-mode equation using flat-bottom and sloping-bottom modes with a full 2D numerical calculation: Real part of (top) u and (middle) w for the (left) 2D numerical solution and the modal solution using (center) flat-bottom and (right) sloping-bottom modes. Also shown are modal amplitude (absolute value) at the (bottom left) right boundary and vertical profiles of (bottom center) u and (bottom right) w (real value) at the location of maximum slope on the left side of the ridge.

6. Conclusions

We have derived a representation of linear internal waves propagating in a domain with nonuniform, arbitrary bathymetry, based on a normal mode expansion and coupled-mode equations. It constitutes a simpler, yet exact, reformulation of the initial set of equations, provided variations in range of the environment are smooth (i.e., differentiable). While we have restricted the analysis to variations of bottom depth, the method can be generalized to horizontal variations of any background fields, provided an orthogonality condition can be found and the background state is stationary.

This formulation relies on a complete basis of vertical normal modes. This basis is defined locally, and the propagation and interaction between modes are expressed using variational calculus. This representation has several advantages. The coupled-mode equations express the modal amplitude evolution in a straightforward manner, making explicit the distinction between modes propagating in opposite directions and displaying the coupling term between modes. Using perturbation theory, we derived an analytical expression for the latter, from which an equation for the energy exchange between modes, that is, wave scattering, follows. Our results show that the strength of the energy exchanges scales like the relative topography variations (H−1dH/dx) and the inverse mode number difference. Besides providing a useful framework for the analytical investigation of internal wave scattering, it also saves us from evaluating the horizontal gradients of modal structures numerically, and hence can reduce the numerical cost of the method.

We obtained an approximate solution for the scattering of incoming internal tide by isolated topography by solving the coupled-mode equations in the limit of weak topography. This solution is an extension to baroclinic incoming tide to the solution derived by Llewellyn Smith and Young (2002) in the context of barotropic tide conversion.

Last, we introduced a different vertical mode decomposition that includes the bottom slope in the boundary condition. The basis functions in this representation have the advantage of satisfying the same boundary conditions as the physical fields. While these vertical modes can exhibit a singularity that seems to coincide with slope criticality, the corresponding series converge faster than the usual vertical modes, and they provide a more accurate reconstruction of the dynamical fields especially near the bottom. A potentially interesting property of this decomposition is the intensification of the mode amplitude near the bottom, which matches with physical observations at near-critical slopes and can be related to processes like boundary streaming. In combination with an appropriate parameterization, it could also provide a way to describe bottom-intensified mixing with no need to capture the smallest vertical wavenumbers generated. However, the energy flux can no longer be decomposed into separate modes.

The coupled-mode formulation derived in this paper can be useful for analytical and numerical investigation of internal tide dynamics, as well as for diagnostics of observations and outputs from numerical models. Future work will address the generalization of the present formulation to nonlinearities, three-dimensional configurations and cases with a background geostrophic flow.

Acknowledgments

This work was supported by ONR Award N00014-13-1-0347. We are grateful to an anonymous reviewer and to Jonas Nycander for suggestions that helped to improve the paper.

APPENDIX

Series Expansion and Galerkin Integration

Let us consider two functions f(x, z) and g(x, z). We are concerned with series expansions of the form g=mAmgm in the vertical integrals ⟨fxg⟩ and ⟨fzg⟩. If g and gm do not satisfy the same boundary conditions, one must perform the integrations with care in a Galerkin fashion to avoid exchanging differentiation and summation. This is done using Leibniz’s rule or integrations by parts back and forth, with the series expansion performed at an intermediate stage. We find

 
fxg=x(fg)gxf=ddxfg(fgdzdx)z1z2gxf=ddx(mAmfgm)(fgdzdx)z1z2mAmgmxf=m(dAmdxfgm+Amfxgm)+[f(mAmgmg)dzdx]z1z2.
(A1)

For the terms explicitly involved in the derivation of the coupled-mode equations, this gives for the projection of the LHS of (8):

 
ΦnPΨx=dAndx+mAmΦnPΦmxdHdx[ΦnP(ΨmAmΦm)]|H.
(A2)

The second term on the RHS is given by (22) for nm. For n = m, we make use of Leibniz’s rule again to obtain

 
ΦnPΦnx=12[dHdx(ΦnPΦn)H]=dHdxϕn22|H.
(A3)

Last, for the projection of the RHS of (8), one has

 
ΦnPLΨz=mAmΦnPLzΦm(ΦnPLΨ)|H=mAmΦnPLΦmiknδn,m[ΦnPL(ΨmAmΦm)]|H.
(A4)

The first term in the RHS is computed using the eigenproblem equations in (10) and (11).

REFERENCES

REFERENCES
Baines
,
P. G.
,
1973
:
The generation of internal tides by flat-bump topography
.
Deep-Sea Res.
,
20
,
179
205
, https://doi.org/10.1016/0011-7471(73)90050-8.
Bell
,
T. H.
,
1975a
:
Lee waves in stratified flows with simple harmonic time dependence
.
J. Fluid Mech.
,
67
,
705
722
, https://doi.org/10.1017/S0022112075000560.
Bell
,
T. H.
,
1975b
:
Topographically generated internal waves in the open ocean
.
J. Geophys. Res.
,
80
,
320
327
, https://doi.org/10.1029/JC080i003p00320.
Brekhovskikh
,
L. M.
, and
O.
Godin
,
1992
:
Acoustics of Layered Media II: Point Sources and Bounded Beams. Springer Series on Wave Phenomena, Vol. 10, Springer, 531 pp
.
Bühler
,
O.
, and
M.
Holmes-Cerfon
,
2011
:
Decay of an internal tide due to random topography in the ocean
.
J. Fluid Mech.
,
678
,
271
293
, https://doi.org/10.1017/jfm.2011.115.
Courant
,
R.
, and
D.
Hilbert
,
1966
:
Methods of Mathematical Physics. Vol. 2, John Wiley and Sons, 856 pp
.
Cox
,
D.
, and
H.
Sandstrom
,
1962
:
Coupling of internal and surface waves in water of variable depth
.
J. Oceanogr. Soc. Japan
,
20
,
499
513
.
Desaubies
,
Y.
, and
K.
Dysthe
,
1995
:
Normal-mode propagation in slowly varying ocean waveguides
.
J. Acoust. Soc. Amer.
,
97
,
933
946
, https://doi.org/10.1121/1.412072.
Egbert
,
G. D.
, and
R. D.
Ray
,
2001
:
Estimates of M2 tidal energy dissipation from TOPEX/Poseidon altimeter data
.
J. Geophys. Res.
,
106
,
22 475
22 502
, https://doi.org/10.1029/2000JC000699.
Falahat
,
S.
,
J.
Nycander
,
F.
Roquet
, and
M.
Zarroug
,
2014
:
Global calculation of tidal energy conversion into vertical normal modes
.
J. Phys. Oceanogr.
,
44
,
3225
3244
, https://doi.org/10.1175/JPO-D-14-0002.1.
Garrett
,
C.
, and
E.
Kunze
,
2007
:
Internal tide generation in the deep ocean
.
Annu. Rev. Fluid Mech.
,
39
,
57
87
, https://doi.org/10.1146/annurev.fluid.39.050905.110227.
Gill
,
A. E.
,
1982
:
Atmosphere–Ocean Dynamics. International Geophysics Series, Vol. 30, Academic Press, 662 pp
.
Griffiths
,
S. D.
,
2011
:
Global modelling of internal tides. Second Norway-Scotland Internal Waves Symp., Edinburgh, United Kingdom, Royal Society, 4 pp.
Griffiths
,
S. D.
, and
R. H. J.
Grimshaw
,
2007
:
Internal tide generation at the continental shelf modeled using a modal decomposition: Two-dimensional results
.
J. Phys. Oceanogr.
,
37
,
428
451
, https://doi.org/10.1175/JPO3068.1.
Jackson
,
D.
,
1914
:
On the degree of convergence of Sturm–Liouville series
.
Trans. Amer. Math. Soc.
,
15
,
439
466
, https://doi.org/10.1090/S0002-9947-1914-1500989-4.
Kelly
,
S. M.
,
2016
:
The vertical mode decomposition of surface and internal tides in the presence of a free surface and arbitrary topography
.
J. Phys. Oceanogr.
,
46
,
3777
3788
, https://doi.org/10.1175/JPO-D-16-0131.1.
Kelly
,
S. M.
,
J. D.
Nash
,
K. I.
Martini
,
M. H.
Alford
, and
E.
Kunze
,
2012
:
The cascade of tidal energy from low to high modes on a continental slope
.
J. Phys. Oceanogr.
,
42
,
1217
1232
, https://doi.org/10.1175/JPO-D-11-0231.1.
Kelly
,
S. M.
,
N. L.
Jones
, and
J. D.
Nash
,
2013
:
A coupled model for Laplace’s tidal equations in a fluid with one horizontal dimension and variable depth
.
J. Phys. Oceanogr.
,
43
,
1780
1797
, https://doi.org/10.1175/JPO-D-12-0147.1.
Kelly
,
S. M.
,
P. F. J.
Lermusiaux
,
T. F.
Duda
, and
P. J.
Haley
Jr.
,
2016
:
A coupled-mode shallow-water model for tidal analysis: Internal tide reflection and refraction by the Gulf Stream
.
J. Phys. Oceanogr.
,
46
,
3661
3679
, https://doi.org/10.1175/JPO-D-16-0018.1.
Klymak
,
J. M.
,
H. L.
Simmons
,
D.
Braznikov
,
S.
Kelly
,
J. A.
MacKinnon
,
M. H.
Alford
,
R.
Pinkel
, and
J. D.
Nash
,
2016
:
Reflection of linear internal tides from realistic topography: The Tasman continental slope
.
J. Phys. Oceanogr.
,
46
,
3321
3337
, https://doi.org/10.1175/JPO-D-16-0061.1.
Kunze
,
E.
,
2017
:
The internal-wave-driven meridional overturning circulation
.
J. Phys. Oceanogr.
,
47
,
2673
2689
, https://doi.org/10.1175/JPO-D-16-0142.1.
Llewellyn Smith
,
S. G.
, and
W. R.
Young
,
2002
:
Conversion of the barotropic tide
.
J. Phys. Oceanogr.
,
32
,
1554
1566
, https://doi.org/10.1175/1520-0485(2002)032<1554:COTBT>2.0.CO;2.
Maas
,
L. R. M.
,
2011
:
Topographies lacking tidal conversion
.
J. Fluid Mech.
,
684
,
5
24
, https://doi.org/10.1017/jfm.2011.245.
Mathur
,
M.
,
G. S.
Carter
, and
T.
Peacock
,
2014
:
Topographic scattering of the low-mode internal tide in the deep ocean
.
J. Geophys. Res. Oceans
,
119
,
2165
2182
, https://doi.org/10.1002/2013JC009152.
Maugé
,
R.
, and
T.
Gerkema
,
2008
:
Generation of weakly nonlinear nonhydrostatic internal tides over large topography: A multi-modal approach
.
Nonlinear Processes Geophys.
,
15
,
233
244
, https://doi.org/10.5194/npg-15-233-2008.
Melet
,
A.
,
R.
Hallberg
,
S.
Legg
, and
K.
Polzin
,
2013
:
Sensitivity of the ocean state to the vertical distribution of internal-tide-driven mixing
.
J. Phys. Oceanogr.
,
43
,
602
615
, https://doi.org/10.1175/JPO-D-12-055.1.
Milder
,
D. M.
,
1969
:
Ray and wave invariants for SOFAR channel propagation
.
J. Acoust. Soc. Amer.
,
46
,
1259
1263
, https://doi.org/10.1121/1.1911850.
Müller
,
P.
, and
N.
Xu
,
1992
:
Scattering of oceanic internal gravity waves off random bottom topography
.
J. Phys. Oceanogr.
,
22
,
474
488
, https://doi.org/10.1175/1520-0485(1992)022<0474:SOOIGW>2.0.CO;2.
Müller
,
P.
, and
X.
Liu
,
2000
:
Scattering of internal waves at finite topography in two dimensions. Part I: Theory and case studies
.
J. Phys. Oceanogr.
,
30
,
532
549
, https://doi.org/10.1175/1520-0485(2000)030<0532:SOIWAF>2.0.CO;2.
Munk
,
W.
, and
C.
Wunsch
,
1998
:
Abyssal recipes II: energetics of tidal and wind mixing
.
Deep-Sea Res. I
,
45
,
1977
2010
, https://doi.org/10.1016/S0967-0637(98)00070-3.
Nycander
,
J.
,
2005
:
Generation of internal waves in the deep ocean by tides
.
J. Geophys. Res.
,
110
, https://doi.org/10.1029/2004JC002487.
Pétrélis
,
F.
,
S. G.
Llewellyn Smith
, and
W. R.
Young
,
2006
:
Tidal conversion at a submarine ridge
.
J. Phys. Oceanogr.
,
36
,
1053
1071
, https://doi.org/10.1175/JPO2879.1.
Rainville
,
L.
, and
R.
Pinkel
,
2006
:
Propagation of low-mode internal waves through the ocean
.
J. Phys. Oceanogr.
,
36
,
1220
1236
, https://doi.org/10.1175/JPO2889.1.
Robinson
,
R.
,
1969
:
The effects of a vertical barrier on internal waves
.
Deep-Sea Res. Oceanogr. Abstr.
,
16
,
421
429
, https://doi.org/10.1016/0011-7471(69)90030-8.
Shimizu
,
K.
,
2011
:
A theory of vertical modes in multilayer stratified fluids
.
J. Phys. Oceanogr.
,
41
,
1694
1707
, https://doi.org/10.1175/2011JPO4546.1.
Sjöberg
,
B.
, and
A.
Stigebrandt
,
1992
:
Computations of the geographical distribution of the energy flux to mixing processes via internal tides and the associated vertical circulation in the ocean
.
Deep-Sea Res.
,
39
,
269
291
, https://doi.org/10.1016/0198-0149(92)90109-7.
St. Laurent
,
L.
, and
C.
Garrett
,
2002
:
The role of internal tides in mixing the deep ocean
.
J. Phys. Oceanogr.
,
32
,
2882
2899
, https://doi.org/10.1175/1520-0485(2002)032<2882:TROITI>2.0.CO;2.
Vallis
,
G. K.
,
2006
:
Atmospheric and Oceanic Fluid Dynamics. Cambridge University Press, 745 pp
.
Wunsch
,
C.
,
1969
:
Progressive internal waves on slopes
.
J. Fluid Mech.
,
35
,
131
144
, https://doi.org/10.1017/S0022112069001005.
Wunsch
,
C.
, and
R.
Ferrari
,
2004
:
Vertical mixing, energy, and the general circulation of the oceans
.
Annu. Rev. Fluid Mech.
,
36
,
281
314
, https://doi.org/10.1146/annurev.fluid.36.050802.122121.
Zhao
,
Z.
,
2017
:
Propagation of the semidiurnal internal tide: Phase velocity versus group velocity
.
Geophys. Res. Lett.
,
44
,
11 942
11 950
, https://doi.org/10.1002/2017GL076008.
Zhao
,
Z.
,
J.
Wang
,
D.
Menemenlis
,
L.-L.
Fu
,
S.
Chen
, and
B.
Qiu
,
2019
:
Decomposition of the multimodal multidirectional M2 internal tide field
.
J. Atmos. Oceanic Technol.
, https://doi.org/10.1175/JTECH-D-19-0022.1.

Footnotes

© 2020 American Meteorological Society. For information regarding reuse of this content and general copyright information, consult the AMS Copyright Policy (www.ametsoc.org/PUBSReuseLicenses).