Aiki et al. Progress in Earth and Planetary Science (2017) 4:11
DOI 10.1186/s40645-017-0121-1

Progress in Earth andPlanetary Science

RESEARCH ARTICLE

Open AccessTowards a seamlessly diagnosable
expression for the energy flux associated withboth equatorial and mid-latitude waves
Hidenori Aiki1,2* , Richard J. Greatbatch3,4 and Martin Claus3
Abstract
For mid-latitude Rossby waves (RWs) in the atmosphere, the expression for the energy flux for use in a model
diagnosis, and without relying on a Fourier analysis or a ray theory, has previously been derived using
quasi-geostrophic equations and is singular at the equator. By investigating the analytical solution of both equatorial
and mid-latitude waves, the authors derive an exact universal expression for the energy flux which is able to indicate
the direction of the group velocity at all latitudes for linear shallow water waves. This is achieved by introducing a
streamfunction as given by the inversion equation of Ertel’s potential vorticity, a novel aspect for considering theenergy flux. For ease of diagnosis from a model, an approximate version of the universal expression is explored andillustrated for a forced/dissipative equatorial basin mode simulated by a single-layer oceanic model that includes both
mid-latitude RWs and equatorial waves. Equatorial Kelvin Waves (KWs) propagate eastward along the equator, arepartially redirected poleward as coastal KWs at the eastern boundary of the basin, and then shed mid-latitude RWs
that propagate westward into the basin interior. The connection of the equatorial and coastal waveguides has beensuccessfully illustrated by the approximate expression of the group-velocity-based energy flux of the present study.
This will allow for tropical-extratropical interactions in oceanic and atmospheric model outputs to be diagnosed interms of an energy cycle in a future study.
Keywords: Group velocity, Model diagnosis, Equatorial Rossby waves, Equatorial mixed Rossby-gravity waves,Equatorial inertia-gravity waves, Equatorial Kelvin waves, Coastal Kelvin waves, Mid-latitude Rossby waves,Mid-latitude inertia-gravity waves, Tropical-extratropical interactions

Introduction
A feature of many phenomena in the equatorial oceans
is the role played by equatorial Kelvin waves (KWs),
examples being El Niño Southern Oscillation (ENSO;
Philander 1989) and the so-called Atlantic Niño (Merle
1980). KWs propagate along the equator and are partially redirected into coastal KWs at the eastern boundary,where they can influence off-equatorial latitudes (e.g.,
Lübbecke et al. 2010) as well as excite extratropical Rossby
waves (RWs) that subsequently propagate into the oceaninterior (McPhaden and Ripa 1990; Isachsen et al. 2007).
*Correspondence: aiki@nagoya-u.jp
Institute for Space-Earth Environmental Research, Nagoya University, Nagoya
City, 464-8601 Aichi, Japan
2
Application Laboratory, Japan Agency for Marine-Earth Science andTechnology, Yokohama, Japan
Full list of author information is available at the end of the article
1

A striking example of this behavior is the equatorial basinmode (Cane and Moore 1981). For the gravest basin mode,
the time scale is set by the time taken for an equatorial
KW to propagate across the basin and for the reflected
gravest long Rossby wave to return to the western boundary (that is 4L/c where L is the basin width and c is thephase propagation speed for KWs). In addition to waves
that are trapped on the equator, equatorial basin modes
also feature coastal KWs that propagate along the easternboundary and extratropical RWs that are excited by theseKWs and refocus on the equator, as described by Schopf
et al. (1981). There is growing evidence that equatorial
basin modes play an important role in equatorial ocean
dynamics. For example, basin modes have been associatedwith the equatorial deep jets (Johnson and Zhang 2003;
Brandt et al. 2011; Claus et al. 2016) and with the semiannual (Thierry et al. 2004) and annual cycles (Brandt© The Author(s). 2017 Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0
International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, andreproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to theCreative Commons license, and indicate if changes were made.

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

et al. 2016) in the equatorial Atlantic. However, the energy
cycle associated with equatorial basin modes has received
little attention and is an important factor when considering the forced/dissipative basin modes that one canrelate to observations. A particularly interesting example is the upward energy propagation associated with the
Atlantic equatorial deep jets (Johnson and Zhang 2003;
Brandt et al. 2011; Mathiessen et al. 2015). Yet, the detailedenergy cycle associated with the jets remains largely
unknown.
One way to approach the energy flux is to use ray theory.
However, ray theory is linked to the dispersion relationof a single type of wave and is not suitable for investigating the sequential connection of different types ofwaves that are associated with a basin mode. Likewise,a Fourier analysis is not suitable for the investigation ofwaves near the coastal boundaries of the ocean. In fact, it
is only for mid-latitude inertia-gravity waves (IGWs) thatthe flux of wave energy has been diagnosed from oceanic
model output (Cummins and Oey 1997; Niwa and Hibiya
2004; Furuichi et al. 2008). On the other hand, in theatmospheric literature, the model diagnosis of pseudomomentum (or wave activity) flux has been more popularthan the model diagnosis of the energy flux (Hoskins et al.
1983; Plumb 1986; Takaya and Nakamura 1997; Nakamura
and Solomon 2011).
Here, we seek a general expression that can be used
to diagnose the energy flux associated with linear shallow water waves at all latitudes from model output. This
manuscript is organized as follows. First we provide thetheoretical background. Then, we present an analyticalinvestigation that leads to a general expression for theenergy flux that can indicate the exact profile of thegroup velocity times wave energy for both equatorial andmid-latitude waves. The utility of the universal expression of energy flux as a model diagnostic is illustrated
for a forced/dissipative equatorial basin mode simulated
by a single-layer model. The model diagnosis is achievedby introducing an inversion for the linearized versionof Ertel’s potential vorticity. This is a novel aspect forconsidering the energy flux in the presence of a coastal
waveguide that connects the equatorial and mid-latitude
regions.

Theoretical backgroundWe use the shallow water equations for a single vertical normal mode (Gill 1982) appropriate to linear waves
in a rotating frame of reference and in the absence ofa mean flow. Let an arbitrary variable with an associated physical dimension be expressed by A∗ , and let
Cartesian-horizontal coordinates be labelled by the setof independent variables x∗ , y∗ , t ∗ , where each of x∗ , y∗increases eastward and northward, respectively, and u∗ , v∗
are the corresponding horizontal components of velocity

Page 2 of 18(a list of variables is given in Table 1)1 . The equations may
then be written as∂p∗∂u∗− f ∗ v∗ + ∗ = 0,
∗
∂t
∂x
∗
∂v∗
∂p
∗ ∗
+
f
u
+
= 0,
∂t ∗
∂y∗


∗
∂v∗
∂p∗
∗ 2 ∂u+c+ ∗ = 0,
∂t ∗
∂x∗
∂y

(1a)
(1b)
(1c)

where f ∗ = f0∗ + β ∗ y∗ is the Coriolis parameter, p∗ =p∗ (x∗ , y∗ , t ∗ ) corresponds to the pressure2 or geopotential,
and c∗ is a uniform constant representing the propagation speed of nonrotating gravity waves for a given mode.
Manipulation of (1a)–(1c) yields a prognostic equation forthe linearized version of Ertel’s potential vorticity (hereafter EPV and symbolized as q∗ ) to read∂u∗f∗ ∗∂  ∂v∗−−
p + v∗ β ∗ = 0,(2)∂t ∗ ∂x∗
∂y∗
c∗ 2



≡q∗

which is applicable to waves at all latitudes, such as midlatitude RWs, mid-latitude IGWs, and equatorial waves[i.e., equatorial RWs and IGWs, equatorial Rossby-gravity
waves (RGWs, i.e., Yanai waves), and equatorial KWs;
Matsuno 1966; Yanai and Maruyama 1966], understanding f0∗ = 0 for an equatorial β-plane and β ∗ = 0 for a
mid-latitude f -plane. Both mid-latitude IGWs (i.e., β ∗ =Table 1 List of symbols, where A∗ and A are arbitrary quantities
written dimensionally or non-dimensionally, respectivelyf ∗ = f0∗ + β ∗ y∗

Coriolis parameterc∗

Speed of long gravity wave

x, y, tCartesian coordinates wherein x and y increase
eastward and northwarda, bHorizontal vector with eastward and northward
components a and b

V = u, v

Horizontal velocity vector∇ ≡ ∂x , ∂y 

Horizontal gradient operatorp

Pressureq ≡ vx − uy − yp

Linearized Ertel’s potential vorticity: q∗ ≡ vx∗∗ − u∗y∗ −(f ∗ /c∗2 )p∗ϕ

Solution of ∇ 2 ϕ − y2 ϕ − 3ϕtt = q, see (16) & (17a)ϕ app

Solution of ∇ 2 ϕ app − y2 ϕ app = q, see (26a) & (18a)(u2 + v2 + p2 )/2

Wave energy: (u∗2 + v∗2 + p∗2 /c∗2 )/2θ = kx − ωtWave phasek

Zonal wavenumberω

Wave frequencyH(n)Hermite polynomial, see endnote 1nMeridional mode number of free equatorial waves

A

Phase average of A

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 3 of 180) and equatorial KWs (i.e., v∗ = 0) are characterized byq∗ = 0, as noted in Table 2.
On the other hand, a prognostic equation for wave
energy may be derived from (1a)–(1c) asp∗ 2∂ 1∗ 2 + v∗ 2 +u∂t ∗ 2c∗ 2+ ∇ ∗ · u∗ p∗ , v∗ p∗  = 0, (3)∂∂where ∇ ∗ ≡  ∂x∗ , ∂y∗  and the overbar symbol representsa phase-average operator (i.e., for a sinusoidal wave, A∗ =0 for A∗ = u∗ , v∗ , and p∗ ) or a low-pass time filter (for thisreason, we retain the local time derivative in (3) to allow
for slow time variations in the general case).
For mid-latitude IGWs in the ocean and atmosphere, thegroup velocity vector points in the same direction as theenergy flux vector in (3):V∗ p∗ = u∗ p∗ , v∗ p∗ ,(4)

a property that has been exploited by Cummins and Oey(1997), Niwa and Hibiya (2004), and Furuichi et al. (2008)
for a model diagnosis. However, for mid-latitude RWs, thevector in (4) does not point in the direction of the groupvelocity of the waves (Longuet-Higgins 1964; Masuda
1978; Cai and Huang 2013). In order to retrieve the correct
direction for the energy flux associated with mid-latitude
RWs, Orlanski and Sheldon (1993, hereafter OS93) havesuggested to modify (3), without affecting the horizontal
divergence of the energy flux, asp∗ 2∂ 1
u∗ 2 + v∗ 2 + ∗ 2∗∂t 2c∇ ∗ · u ∗ p∗ +

∂
∂y∗

p∗ 2
2f ∗

+
, v∗ p∗ −

∂
∂x∗

p∗ 2
2f ∗

= 0,
(5)

where each of u∗ and v∗ should be the sum of thegeostrophic and ageostrophic components and f ∗ = f0∗ +Table 2 Characteristics of different waves at various latitudes

β ∗ y∗ is understood. The energy flux vector in (5) consistsof two terms,V∗ p∗ − ∇ ∗ ×[ p∗ 2 /(2f ∗ )] z,(6)where V∗ p∗ is as in the gravity wave literature (i.e., V∗ isthe sum of the geostrophic and ageostrophic componentsof velocity). The second term in (6) is the additional rotational component required to reproduce the direction of
the group velocity of mid-latitude RWs (z is the upwardvertical unit vector). In Longuet-Higgins (1964), the second term of (6) has been expressed as −∇ ∗ ×[ f ∗ ψ ∗2 /2] z
where ψ ∗ is a streamfunction based on the assumptionof horizontally nondivergent velocity. This assumptionis hardly used in modern oceanography owing to thesmallness of the deformation radius. In quasi-geostrophic
theory, ψ ∗ = p∗ /f ∗ from which the connection with
(6) is clear.
The question naturally arises as to whether or not it
is possible to find a general expression for the additionalrotational flux, R∗ , that holds for waves at all latitudes
and is such that the corresponding energy flux V∗ p∗ + R∗
always points in the direction of the group velocity andthus constitutes a general expression for the energy flux
associated with waves at all latitudes. This is the main subject of the present study. In this study, we focus on wave
types for which the group velocity has been well formulated in the literature/textbook, as listed in Table 2. Of
particular interest is the energy flux associated with equatorial RWs given that the expression in (6) is singular at
the equator. The assumption of horizontally nondivergent
velocity in Longuet-Higgins (1964) is also inappropriate
for equatorial regions. In the next section, by investigating the analytical solution of equatorial waves, we derive
an exact universal expression for the rotational flux which,after being added to V∗ p∗ , is able to indicate the directionof the group velocity for linear waves at all latitudes.

Analytical investigation
We begin by revisiting analytical expressions for the profile of the energy flux associated with equatorial waves.
This investigation allows us to derive an expression for
the energy flux that points in the direction of the groupvelocity for waves at all latitudes.

V∗ p∗ parallel toq∗ = 0group velocity(ϕ app∗ = 0)

Equatorial Rossby wave

No

No

Equatorial mixed
Rossby-gravity wave

Depends on frequency No

Equatorial inertia-gravity
wave

Roughly yesNoEquatorial Kelvin wave

Yes

Yes

Coastal Kelvin wave

Yes

Yes

Mid-latitude Rossby wave

No

No (ϕ app∗  p∗ /f ∗ )ut − yv + px = 0,

(7a)

Yes

vt + yu + py = 0,(7b)pt + ux + vy = 0,(7c)Mid-latitude inertia-gravity Yes
wave

Energy flux associated with equatorial waves

We assume linear waves in the absence of a mean flowon an equatorial β-plane. As in√Matsuno (1966) and Gill
(1982), we use a time scale 1/ c∗ β ∗ and a length scale√c∗ /β ∗ to nondimensionalize the equation system (1a)–(1c) to giveAiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 4 of 18

where symbols without an asterisk indicate nondimensionalized quantities and subscripts indicate partial differentiations. Manipulation of (7a)–(7c) yields prognostic
equations for EPV and wave energy in a nondimensionalized form to read,∂t (vx − uy − yp) + v = 0,




(8)

∂t (u2 + v2 + p2 )/2 + ∇ · up, vp = 0,

(9)

≡q

∂
∂ ∂
where ∂t ≡ ∂t
, ∇ ≡  ∂x, ∂y , and for A = u, v, or p, A = 0
for sinusoidally varying waves.
In what follows, we assume v ≡ 0 which is appropriate
for equatorial RWs, RGWs, and IGWs (i.e., waves otherthan equatorial KWs). Then, we consider zonally propagating free waves with a relationship v ∝ cos θ, u ∝ sin θ,and p ∝ sin θ where θ ≡ kx − ωt is wave phase with k andω being wavenumber and wave frequency, respectively.
Substitution of these relationships to (7a)–(7c), followedby some manipulation, yields a characteristic equation forthe meridional structure of v to read,vyy + (ω2 − k 2 − k/ω − y2 )v = 0.

(10)

Matsuno (1966) has derived a solution for (7a)–(7c) and(10) to yield,
v = A cos θ exp(−y2 /2)H (n) ,

(11a)

2

2

(11b)

2

2

(11c)

u = (ωyvθ − kvyθ )/(ω − k ),p = (kyvθ − ωvyθ )/(ω − k ),group velocity times the meridional integral of the wave
energy (Philander 1989):
 +∞−∞up dy = (∂ω/∂k)

 +∞
−∞

(u2 + v2 + p2 )/2 dy.
(14a)

It should be noted that the identity (14a) does not hold if
it is evaluated without the meridional integral:
up = (∂ω/∂k)(u2 + v2 + p2 )/2.

(14b)

For low-frequency equatorial waves (with ω < 1—
see Fig. 1—, i.e., all equatorial RWs and westward
propagating RGWs), the meridional profiles of up and(∂ω/∂k)(u2 + v2 + p2 )/2 are shown by the dashed green
and solid black lines, respectively, in Fig. 2. It is clear that,when compared at a given latitude, up is not equal to thegroup velocity times wave energy. In particular, the meridional profile of up is sign-indefinite for low-frequency
equatorial waves (Fig. 2). On the other hand, as shown by
the dashed green and solid black lines in Fig. 3 for highfrequency equatorial waves (with ω > 1— see Fig. 1—, i.e.,
all equatorial IGWs and eastward propagating RGWs), themeridional profile of up provides a much better approximation for the group velocity times wave energy. The solid
blue line, dashed orange line, and purple dots in Figs. 2
and 3 are explained later in the manuscript.

where A is wave amplitude and the symbol H (n) is the
Hermite polynomial with n being the meridional mode
number3 . The subscript θ represents partial differentiation in terms of the wave phase [i.e., vθ ≡ ∂v/∂θ =−A sin θ exp(−y2 /2)H (n) ].
Substitution of (11a) to (10) yields,
ω3 − (k 2 + 2n + 1)ω − k = 0,(12)which is a unified dispersion relation for equatorial RWs,
RGWs, and IGWs. Partial differentiation of (12) withrespect to wavenumber k yields a unified expression for
the group velocity of equatorial waves,
2kω + 12ω2 + ω/k
∂ω
=
=
,
2
2
∂k
3ω − (k + 2n + 1)
2ω3 /k + 1(13)

where 2ω3 /k in the denominator has often been ignoredin previous studies when focusing on low-frequency equatorial waves (e.g., equatorial RWs; Gill 1982).
We now investigate the energy flux associated with
(7a)–(7c). It is known that, for zonally propagating equatorial waves, the meridional integral of up is equal to theFig. 1 Dispersion relation of free equatorial waves for a given
meridional number n, as given by Eq. (12). Axes have been√
non-nondimensionalized by β ∗ c∗ for wave frequency ω and by√ ∗ ∗β /c for zonal wavenumber k. The red dots locate the parameters
used to produce the 10 panels in Figs. 2 and 3

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 5 of 18

Fig. 2 Meridional profiles of the zonal component of phase-averaged energy flux associated with equatorial waves; the solid black line is(∂ω/∂k)(u2 + v2 + p2 )/2, the dashed green line is up, the dashed orange line is up + (pϕ app /2 + utt ϕ app )y , the solid blue line is up + (pϕ app /2)y , andthe purple dots are (ω/k)(E − v2 ) = (ω/k)(u2 − v2 + p2 )/2. All panels are for low-frequency equatorial waves with ω < 1: a westward propagating
RGWs, b short and c long RWs in the 1st meridional mode, and d short and e long RWs in the 2nd meridional mode. The associated values of
meridional-mode number n, zonal wavenumber k, wave frequency ω, and group velocity ∂ω/∂k are noted in each panel. For each of a–e, the wave∞
amplitude A in (11a) has been set to normalize the meridional integral of wave energy: −∞ (u2 + v2 + p2 )/2dy = 1

a

b

c

d

e

Fig. 3 Same as Fig. 2 except for high-frequency equatorial waves with ω > 1: a westward and b eastward propagating IGWs in the 2nd meridionalmode, c westward and d eastward propagating IGWs in the 1st meridional mode, and e eastward propagating RGWs

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 6 of 18

Identification of the additional rotational flux associatedwith equatorial waves

It is useful to derive the analytical expression for the difference between the left and right hand sides of (14b). A
first step for identifying the difference is to decompose
the zonal component of up into two parts, one that determines the meridional integral and one that does not affectit, as follows:
up


 

2

= y2 vv(ωk) − yvy v ω2 + k 2 + vy vy (ωk) / ω2 − k 2= vyy v(ωk) + vv ω3 k − ωk 3 − k 2


 

2
− (yvv/2)y − vv/2 ω2 + k 2 + vy vy (ωk) / ω2 − k 2
 
  
= vv 2ω3 k − 2ωk 3 − k 2 + ω2 + vy v y (2ωk)
  2 − (yvv)y ω2 + k 2 / 2 ω2 − k 2 = vv (2ωk + 1) / 2 ω2 − k 2  
2 + vy v(2ωk) − yvv ω2 + k 2 y / 2 ω2 − k 2,(14c)where the first equality has been derived using (11b)–(11c) and sin θ sin θ = cos θ cos θ and the second equality
has been derived using (10). Note that it is the secondof the two terms whose meridional integral is zero (noting that v and yv go to zero at large distances from theequator).
We now decompose the wave energy4 into two parts,
one that determines the meridional integral and one doesnot. We then have
u2 + v2 + p2 /2= vv/2 + y2 vv + vy vy ω2 + k 2
2 

  

− yvy v (4kω) / 2 ω2 − k 2  
= vv/2 + y2 vv − vyy v + vy v y ω2 + k 2
  
2 
− (yvv)y (2kω) + (vv) (2kω) / 2 ω2 − k 2 = vv k 4 − 2k 2 ω2 + ω4 + ω4 − k 4 + kω − k 3 /ω
  

2 
  
+ vy v y ω2 + k 2 − (yvv)y (2kω) / 2 ω2 − k 2


  
= vv 2ω2 + k/ω / 2 ω2 − k 2

2 
  
 
,
+ vy v ω2 + k 2 − yvv(2kω) y / 2 ω2 − k 2(14d)where the first equality has been derived using (11b)–(11c) and sin θ sin θ = cos θ cos θ , and the third equality
has been derived using (10). As before, it is the secondof the two terms whose meridional integral is zero. Using(14c)–(14d), we now obtain an analytical expression forthe difference between the right and left hand sides of(14b) to yield(∂ω/∂k)(u2 + v2 + p2 )/2 − up(vy v)y(2ω2 + ω/k)(ω2 + k 2 )−2ωk=2(ω2 − k 2 )22ω3 /k + 1(yvv)y(2ω2 + ω/k)2kω
2
2
+
k
)
−
(ω
−
2(ω2 − k 2 )22ω3 /k + 1(vy v)y (2ω4 + 2ω2 k 2 + ω3 /k + ωk − 4ω4 − 2ωk)2(ω2 − k 2 )2 (2ω3 /k + 1)
(yvv)y (4ω3 k + 2ω2 − 2ω5 /k − 2ω3 k − ω2 − k 2 )−2(ω2 − k 2 )2 (2ω3 /k + 1)
(vyθ vθ )y (ω/k − 2ω2 ) − (yvθ vθ )y (1 − 2ω3 /k)=2(ω2 − k 2 )(2ω3 /k + 1)=

[ (ωvyθ − kyvθ )vθ ]y +[ (−kvyθ + ωyvθ )vθ ]y (2ω2 )2k(ω2 − k 2 )(2ω3 /k + 1)−(pvθ )y − (2utt vθ )y
,
(14e)
=
2k(1 + 2ω3 /k)=

where the first and second equalities have been derived
using (13), the third equality has been derived usingcos θ cos θ = sin θ sin θ , and the last equality has been
derived using (11b)–(11c). The last line of (14e) has been
written as the meridional gradient of scalar quantities.
Thus, the meridional integral of (14e) vanishes for equatorial waves (with a meridionally decaying structure) andis consistent with (14a).
Using (14e), we can now rewrite the zonal componentof the group velocity times wave energy as(∂ω/∂k)(u2 + v2 + p2 )/2 = up + (pϕ/2 + utt ϕ)y ,3

ϕ ≡ −vθ /(k + 2ω ),

(15a)
(15b)

where the scalar quantity ϕ has been introduced. We haveconfirmed that, as long as ϕ is set by (15b), the meridional profile of the zonal energy flux, up + (pϕ/2 + utt ϕ)y ,in (15a) is precisely identical to (∂ω/∂k)(u2 + v2 + p2 )/2
for all types of equatorial waves in Figs. 2 and 3. Namely,
all solid black lines in Figs. 2 and 3 may be drawn using
either expression. As far as we know, (15a) and (15b) havenot been mentioned in previous studies and therefore
constitute a new result.
Inversion equations for Ertel’s potential vorticity

The definition of ϕ, as given by (15b), is based on a Fourier
expansion. However, we have found that (15b) may be
rewritten into an expression which contains none of θ, k,and ω to read
∇ 2 ϕ − y2 ϕ − 3ϕtt = −vθ /ω
= q,

(16)

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 7 of 18

where ∇ 2 ≡ ∂xx + ∂yy is understood, the first line has been
derived using (10), and the second line has been derived
using (8) [i.e., qt = −ωqθ = −v and thus −ωqθθ =ωq = −vθ ]. The new Eq. (16) of EPV is the cornerstone of the present study, because it suggests a possibility
for the scalar quantity ϕ to be estimated without using a
Fourier analysis. This feature is important for identifying
the direction of the energy flux of waves in the presence ofcoastal boundaries.
To summarize, in order to reproduce the profile of
the group velocity times wave energy without relying on
a Fourier analysis, we have obtained a new expression
for the energy flux that has turned out to be associatedwith the streamfunction Eq. (16). Equation (16) may be
rewritten into a dimensional form as∇ ∗2 ϕ ∗ − (f ∗ /c∗ )2 ϕ ∗ − (3/c∗ 2 )ϕt∗∗ t∗ = q∗ , (17a)where ∇ ∗ ≡ ∂x∗ , ∂y∗  and q∗ = v∗x∗ − u∗y∗ − (f ∗ /c∗2 )p∗ .
The exact profile of the group velocity times wave energy
may be reproduced by the right hand side of (15a) and ishere rewritten into a vector and dimensional form asV∗ p∗ − ∇ ∗ ×[ (p∗ ϕ ∗ )/2 + (u∗t∗ t∗ ϕ ∗ )/β ∗ ] z.

(17b)

The additional rotational flux in (17b) corrects the profile of the energy flux, without affecting the divergenceof the energy flux. The quantity ϕ ∗ in (17b) is the solution of the accurate streamfunction Eq. (17a) associatedwith EPV in a dimensional form. We note in passingthat for zonally propagating equatorial waves, as given by(11a)–(11c), v∗ p∗ vanishes owing to the phase relationshipbetween v∗ and p∗ [see (11a) and (11c)] and the meridional
component of the additional rotational flux, −(p∗ ϕ ∗ /2 +u∗t∗ t∗ ϕ ∗ /β ∗ )x∗ , also vanishes.
Equatorial KWs

So far, we have not investigated the energy flux of equatorial KWs. Since KWs are gravity waves, V∗ p∗ becomesequal to the group velocity times wave energy. Namely, theadditional rotational flux is absent. KWs are also characterized by q∗ = 0; hence, the EPV equation (17a) yields
ϕ ∗ = 0. The result is that, in the case of KWs, the expression for the energy flux, as given by (17b) reduces to V∗ p∗ ,which is consistent with the nature of gravity waves.
Boundary conditions and the connection to mid-latitude
regions

Consider a basin with closed zonal boundaries (i.e., theeastern and western coastlines of a basin of arbitrary
shape). It is clear that the flux V∗ p∗ in (17b) has nocomponent normal to the zonal boundaries. Hence, theadditional rotational flux in (17b) should also have nocomponent crossing the closed boundaries. This require-ment is satisfied in the present study by solving (17a) witha boundary condition ofϕ ∗ = 0.

(17c)

In a general situation in the ocean, waves propagating
eastward along the equatorial waveguide are partiallyredirected poleward as KWs along the eastern boundary
where they can shed RWs that then propagate westward
into the ocean interior (Cane and Moore 1981; Philander
1989; Chelton and Schlax 1996; Isachsen et al. 2007).
We now investigate whether or not the set of (17a)and (17b) is applicable to off-equatorial regions where
small-amplitude perturbations are characterized by either
mid-latitude RWs or IGWs. For perturbations associatedwith mid-latitude RWs, the solution ϕ ∗ of (17a) corresponds to the geostrophic streamfunction for which ϕ ∗ p∗ /f ∗ is a reasonable approximation in an interior region(i.e., far from coastal boundaries), noting that ∇ ∗ 2 ϕ ∗ corresponds to v∗x∗ − u∗y∗ . The result is that the energy fluxin (17b) automatically reduces to the expression of OS93
for mid-latitude RWs5 . On the other hand, if perturbations associated with mid-latitude IGWs are given, theinversion Eq. (17a) of EPV, which equals zero, yields, withϕ ∗ = 0 on the boundaries, ϕ ∗ = 0 everywhere. Thus, theenergy flux in (17b) automatically reduces to V∗ p∗ which
represents the group velocity of mid-latitude IGWs times
wave energy. We conclude that the set of (17a) and (17b)
can represent the exact profile of the group velocity times
wave energy associated with both mid-latitude IGWs andRWs, which may be reconfirmed using almost the sameprocedure as in the “Identification of the additional rotational flux associated with equatorial waves” section. See
Appendix 1 for details.

Methods/Experimental
The rest of this manuscript presents an example illustrating the diagnosis of the energy flux from a model. To be
useful for our discussion, the exact universal expression
for both equatorial and mid-latitude waves, as given by
the set of (17a) and (17b), is hereafter referred to as thelevel-0 energy flux. In practice, the level-0 expression of
the energy flux is not straightforward to compute frommodel output, since the second-order time derivative termin (17a) makes it difficult to solve for ϕ ∗ .
For the present study, we investigate the consequenceof artificially removing the second-order time derivative
term from (17a) to give∇ ∗2 ϕ app∗ − (f ∗ /c∗ )2 ϕ app∗ = q∗ ,(18a)which may be justified at least for low-frequency waves(e.g., both equatorial and mid-latitude RWs) based onscale analysis. The superscript of ϕ app∗ indicates that thesolution of (18a) may be regarded as an approximation forAiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 8 of 18

the solution ϕ ∗ of the accurate streamfunction Eq. (17a)associated with EPV. Then, we replace ϕ ∗ in (17b) withϕ app∗ to readV∗ p∗ − ∇ ∗ ×[ (p∗ ϕ app∗ )/2 + (u∗t∗ t∗ ϕ app∗ )/β ∗ ] z,
(18b)
which is diagnosable6 from model output and is referred
to as the level-1 expression of the energy flux in thepresent study. As shown by the dashed orange lines inFig. 2, the level-1 expression provides a nice approximation for the group-velocity-based energy flux of lowfrequency equatorial waves, but not for high-frequency
equatorial waves in Fig. 3. Next, with the form of the additional rotational flux −∇ ∗ ×[ p∗ 2 /(2f ∗ )] z in (6) in mind,
we investigate the consequence of simplifying (18b) asV∗ p∗ − ∇ ∗ × (p∗ ϕ app∗ /2)z,(18c)

which we refer to as the level-2 expression for the energy
flux. As shown by the solid blue lines in Figs. 2 and 3,
the level-2 expression provides an approximation for thegroup-velocity-based energy flux of both low- and highfrequency equatorial waves, although there can be someerror. Further discussion of the level-2 approximation is
given in Appendices 2 and 3 where it is noted that thelevel-2 approximation is comparable in accuracy to thepseudomomentum (or wave-activity) flux used in previous studies (Randel and Williamson 1990; Brunet and
Haynes 1996; Fukutomi and Yasunari 2002; Wakata andKitaya 2002; Kawatani et al. 2010).
We now contrast both the level-1 and level-2 energy
fluxes with the expressions in previous studies, given by(6) and (4), using a solution from a linear shallow watermodel. This illustrates the potential of the expression
given by (18b) and (18c) for use as a model diagnostic(see Table 3). Suitable for this purpose is an equatorial basin mode solution since it is associated with both
equatorial and coastal waveguides as well as the radiation of mid-latitude RWs into the basin interior. Furthermore, as noted in the Introduction section, the equatorial

basin mode, first studied by Cane and Moore (1981), hasrecently attracted attention because of its importance inthe dynamics of the equatorial Atlantic Ocean. Indeed, theannual cycle, the semi-annual cycle, and the interannual
variability associated with the Atlantic equatorial deep jets
(Brandt et al. 2011) all appear to be resonant excitations ofequatorial basin modes [see Brandt et al. (2016) and Clauset al. (2016) for more details].
Model set-up

To illustrate the importance of dissipation for explaining
the observed cross-equatorial width of the equatorial deep
jets, Greatbatch et al. (2012, hereafter G12) have simulated
a forced/dissipative basin mode solution using a singlelayer reduced-gravity linear model. The model is set up
in spherical coordinates, with a rectangular domain in latitude/longitude space of roughly the same width as the
Atlantic Ocean at the equator (that is 55◦ in longitude)
and reaching to 10◦ N/S on either side of the equator7 .
All lateral boundaries are closed. In both G12 and Clauset al. (2014, hereafter C14), the model has been forced by
an idealized oscillatory forcing with a period of 4.5 years
in the zonal momentum equation to mimic the forcing of
the jets, together with a lateral mixing of momentum thatprovides dissipation. [See Ascani et al. (2015) for a discussion on the forcing of the equatorial deep jets, the detailsof which are not important here]. It should be noted that
4.5 years is roughly the time taken for an equatorial KW
and the reflected long gravest equatorial RW, to travel
across the basin for the vertical mode that is closest to
resonance. As noted in G12 and C14, the (westward) propagation speed of equatorial long RWs is three times lessthan the (eastward) propagation speed of equatorial KWs[see the dispersion relation (12)].
Our model has been set up as in G12 and C14. The gravity wave speed is set equal to c∗ = 0.17 m/s [see the upper
panel in Fig.
√ 4 of C14]. The equatorial deformation radius
becomes c∗ /β ∗ = 87 km, with a consequence that disturbances further than a few degrees from the equator in
our model experiment may be regarded as mid-latitude

Table 3 List of energy flux vectors and EPV-based streamfunctions in dimensional form and their location in the text and figures
Approx.

Energy flux vector

Equation

Level-0(15a), (17b), (23a)Solid black–(18b), (26b)Dashed orange–

Level-2V∗ p∗ − ∇ ∗ ×[ (p∗ ϕ ∗ )/2 + (u∗t∗ t∗ ϕ ∗ )/β ∗ ] zV∗ p∗ − ∇ ∗ ×[ (p∗ ϕ app∗ )/2 + (u∗t∗ t∗ ϕ app∗ )/β ∗ ] zV∗ p∗ − ∇ ∗ × (p∗ ϕ app∗ )/2z

(18c), (26c)

Solid blue

(c)

QG

V∗ p∗ − ∇ ∗ ×[ p∗2 /(2f ∗ )] z(5), (6)–(b)

f -PlaneV∗ p∗(3), (4)Dashed green(a)Level-1

Definition of EPV-based streamfunctions∇ ∗2 ϕ ∗ − (f ∗ /c∗ )2 ϕ ∗ − (3/c∗2 )ϕt∗∗ t∗ = q∗ ,(15b), (16), (17a), (23c)∇ ∗2 ϕ app∗ − (f ∗ /c∗ )2 ϕ app∗ = q∗ ,(18a), (26a)Figs. 2 & 3

Figs. 6 & 7

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 9 of 18

condition of no motion and no pressure anomaly, themodel has been integrated for 20 cycles (i.e., 90 years)
using the oscillatory forcing which is sufficient for a steady
oscillatory state to be reached. Since the model code is
fully non-linear, we have set the amplitude of the forcing to a small value, 1.0 × 10−10 m/s2 to ensure thatlinear dynamics prevails. Indeed, the magnitude of velocity associated with the gravest basin mode may be scaled
as 10−10 m/s2 × 4.5 years/(2π) = 0.0023 m/s, which
results in a Froude number of (0.0023 m/s)/c∗ = 0.014(nondimensional). These parameters are summarized inTable 4. Below, we show results from an experiment whichcorresponds to the “full” case in G12. In particular, theoscillatory zonal forcing is spatially uniform and acts over
the whole model domain. All the model results shown
below are averages over the last model cycle.

Table 4 Parameters in the model experiment of the present
study

Equatorial inertial periodc∗ = 0.17 m/s√ ∗ ∗c /β = 87 km√
2π/ c∗ β ∗ = 37 days

Forcing periodT ∗ = 4.5 years

Forcing amplitude

10−10 m/s2

Long gravity wave speed
Equatorial deformation radius

Forcing area

Full domain

Domain size55◦ (zonal) × 20◦ (meridional)

Forcing Froude number(0.0023 m/s)/c∗ = 0.014

Horizontal resolution

0.1◦ (zonal) × 0.1◦ (meridional)

Lateral eddy viscosity

10 m2 /s

RWs, even though they are part of the equatorial basinmode resonance. As in G12, our model has been formulated in a spherical coordinate system with a grid spacingof 0.1◦ in both longitude and latitude. The coefficient8 ofeddy viscosity has been set to 10 m2 /s. From an initial

Results and discussion
At each time step of the model output, we have calculated the EPV-based streamfunction ϕ app∗ (contours in

a

b

c

d

e

f

g

h

Fig. 4 Left panels are the snapshots of the EPV-based streamfunction ϕ app∗ (solid and dotted line contours for positive and negative values,
respectively, with an interval of 60 m2 /s) and thickness anomaly (color shading, normalized). Right panels are the snapshots of the zonal componentof velocity u∗ (color shading, m/s). These snapshots have been adapted from the 20th cycle of the model experiment, wherein the period of the√oscillatory forcing is T ∗ = 4.5 years; the reference layer thickness is H∗ = 100 m, and the equatorial deformation radius is c∗ /β ∗ = 87 km. See also
movie (Additional file 1)

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

the left panels of Fig. 4) by solving the spherical coordinate version of (18a) with the boundary condition of
ϕ app∗ = 0. The color shading in Fig. 4 shows the snapshots of thickness anomaly (left panels) and the zonal
component of velocity u∗ (right panels). The movie of themodel experiment is found in Additional file 1. RWs are
identified by the correlation (anticorrelation) between the
EPV-based streamfunction and thickness anomaly in thenorthern (southern) hemisphere. This follows from thecorrespondence between the EPV-based streamfunction
and the geostrophic streamfunction for the case of midlatitude RWs, as noted earlier. As noted in G12 and C14,
the (westward) propagation speed of equatorial long RWs
is three times smaller than the (eastward) propagation
speed of equatorial KWs [see the dispersion relation (12)].
It takes a three-quarter cycle (i.e., 3T ∗ /4) for equatorial
long RWs to travel westward from the eastern boundary
to the western boundary of the model domain (see redlines in Fig. 5a). After reflection at the western boundary,it takes only a quarter of a cycle (i.e., T ∗ /4) for equatorial KWs to travel eastward to the eastern boundaryof the model domain (see blue lines in Fig. 5a), wheresome disturbances are deflected poleward along the eastern boundary to be the source of mid-latitude RWs whichthen propagate westward (Fig. 5b).
In Fig. 6, the divergence of the horizontal energy flux,
given by ∇ ∗ · V∗ p∗ , is shown for the whole model domainusing color shading. Red indicates regions of a net energy
input, and blue indicates regions of a net dissipation. It isclear that the main region of energy input is in the central
part of the basin along the equator, where the strongest
zonal velocities are found, and that the main regionsof energy loss are associated with the RWs that radiateaway from the eastern boundary. Arrows in Fig. 6a show
the energy flux used in the gravity-wave literature, V∗ p∗ ,which is mostly westward along the equator and eastward
in the immediate off-equatorial region. This can be clearlyseen in Fig. 7a which shows a blow-up of the eastern equatorial region. Figures 6b and 7b show the energy flux givenby (6), which has been adapted from OS93, where only
regions more than 1◦ latitude away from the equator areplotted to avoid the singularity in the Coriolis parameter f ∗ at the equator. From these figures (especially theblow-up of the eastern equatorial region in Fig. 7a, b), it isclear that the energy flux is strongly reversed when compared to V∗ p∗ in the immediate off-equatorial region andis now strongly eastward in association with RWs that areradiated from the eastern boundary.
From Figs. 6c and 7c, it is clear that when the set ofEqs. (18a), (18c) and (17c) is used to estimate the energy
flux, the westward flux associated with the off-equatorial
RWs is part of a recirculation of energy in the eastern partof the basin (Fig. 7c) with eastward energy flux along theequator and westward energy flux off the equator. The

Page 10 of 18abFig. 5 a Hovmoller diagram at the equator for thickness anomaly(color shading, normalized) with the blue and red lines indicating thephase speeds of equatorial Kelvin waves and long Rossby waves,respectively. b Hovmoller diagram at 2◦ N for thickness anomaly (color
shading, normalized) and the EPV-based streamfunction (contours,
plotted as in the left panels of Fig. 4)

eastward flux along the equator in Figs. 6c and 7c is in the
opposite direction to the westward V∗ p∗ flux in Figs. 6aand 7a along the equator in the same region. This indicates the role of the rotational flux contribution in (18c)which counters the westward V∗ p∗ flux along the equator.
This westward flux is associated with the equatorial RWs
but represents an overestimation of the energy flux associated with these waves (see Fig. 2). When the rotational
flux is added, what emerges is the eastward flux associated with the KW which, in turn, leads to a poleward flux
arising from KWs propagating along the eastern boundary and, in turn, leads to the westward flux associatedwith the off-equatorial RWs that are excited at the easternboundary. Here, in terms of the transfer of wave energy,
the equatorial waveguide has been connected to the eastern coastal waveguide and, in turn, to the basin interior at
off-equatorial latitudes, which is at the heart of the present
study.
Finally, we note that the forcing period of T ∗ = 4.5 years
is muchlonger than the equatorial inertial period of√
2π/ c∗ β ∗ = 37 days. It can be said that the simulated
equatorial basin mode consists of low-frequency equatorial waves, as in Fig. 2, and mid-latitude RWs. We recall theAiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 11 of 18

a

b

c

Fig. 6 Comparisons of three expressions for the horizontal flux ofwave energy (arrows): a the pressure flux [see Eq. (3)], b the pressure
flux plus the additional rotational flux of Orlanski and Sheldon (1993)[see Eq. (4)], and c the pressure flux plus the additional rotational fluxof the present study [see Eq. (15)]. In both b and c, the additional
fluxes have been calculated using the rotation operator in thespherical coordinate system to be consistent with the model
formulation. Green contours in b and c show the distributions ofp∗2 /(2f ∗ ) and p∗ ϕ app∗ /2, respectively (solid and dotted lines indicate
positive and negative values, respectively, with an interval of
10−2 m4 /s3 ) . Color shading in all panels shows the horizontal
divergence of the time-averaged energy flux. Note that the additionalrotational flux in b and c has no influence on the divergence. All
quantities have been calculated from the output of the same
experiment with a time-average between t∗ = 19T ∗ and 20T ∗Fig. 7 Same as Fig. 6 except that this is an enlarged view of theeastern equatorial region of the model domain

we have learned that the latter quantity (not shown) isthree orders of magnitude smaller than the former. Thus,
we conclude that, in the diagnosis of the simulated basinmode, the expression of the energy flux, as given by (18c),
has provided a nice approximation for the group velocity
times wave energy.

Conclusions
small difference between the solid blue and dashed orange
lines in Fig. 2, the former and the latter of which may be
written as u∗ p∗ + (p∗ ϕ app∗ /2)y∗ and u∗ p∗ + (p∗ ϕ app∗ /2 +u∗t∗ t∗ ϕ app∗ /β ∗ )y∗ , respectively, in a dimensional form (seelevel-2 and level-1, respectively, in Table 3). Since arrows
in Figs. 6c and 7c have been plotted using the expression
which corresponds to the solid blue line in Fig. 2, we havechecked for any improvement by using the expression
which corresponds to the dashed orange lines in the same
figure. The checking has been done by comparing the
distribution of p∗ ϕ app∗ /2 and u∗t∗ t∗ ϕ app∗ /β ∗ , from which

In previous studies of the ocean, the energy flux of waves
in model output has been diagnosed using V∗ p∗ , whereV∗ is the horizontal component of velocity perturbationand p∗ corresponds to the pressure perturbation. This isappropriate for understanding the energy flux associatedwith mid-latitude inertia-gravity waves (IGWs). For midlatitude Rossby waves (RWs), however, the direction ofV∗ p∗ differs from the group velocity and hence the energy
flux, by a rotational vector flux with zero divergence. The
rotational flux to be added to V∗ p∗ for estimating thegroup velocity of mid-latitude RWs has previously beenAiki et al. Progress in Earth and Planetary Science (2017) 4:11

derived using quasi-geostrophic equations and is singularat the equator.
By investigating the analytical solution of both equatorial waves (“Analytical investigation” section) and midlatitude waves (Appendix 1), we have derived an exact
universal9 expression for the rotational flux which, after
being added to V∗ p∗ , is able to indicate the profile of
the group velocity times wave energy for linear waves
at all latitudes. This is what we call the level-0 expression of the energy flux. The level-0 energy flux is written
using the solution ϕ ∗ of (17a), previously unmentioned in
the literature, which we refer to as the accurate streamfunction associated with Ertel’s potential vorticity (EPV).
Equation (17a) is the cornerstone of the present study,
because it suggests a possibility for the energy flux to be
estimated (i) without using a Fourier analysis nor ray theory and (ii) in the presence of coastal boundaries, which
will allow for tropical-extratropical interactions in model
output to be diagnosed in terms of an energy cycle ina future study. Presently, the level-0 energy flux is notpractical for use as a model diagnostic, since the secondorder time derivative term in (17a) makes it difficult tosolve for ϕ ∗ . Thus, we hope that a future study is able todevelop a numerical algorithm to solve (17a) for ϕ ∗ . We
also note the need to extend the theory to a continuously
stratified ocean and also to test out the theory in the presence of a sheared mean flow, both of which topics await
a future study. This is a new step from the recent understanding of energetics in the atmosphere and ocean that
had been focused on, for example, the global mapping ofenergy conversion rates associated with various physicalprocesses (e.g., baroclinic and barotropic instabilities) andexternal forcing (Iwasaki 2001; Aiki and Richards 2008;
Zhai et al. 2012).
The potential of our analysis as a model diagnostic is
illustrated in the present study for a forced/dissipativeequatorial basin mode simulated by a single-layer model.
The model result includes both mid-latitude RWs (maintained by coastal KWs propagating poleward along theeastern boundary) and equatorial RWs (maintained by thereflection of equatorial KWs at the eastern boundary).
We have used approximate expressions for the energy
flux (what we call the level-1 and level-2 energy fluxes)
that is based on the inversion equation (18a) of EPV
and which is shown to be good approximations to thelevel-0 expression in the case of the model run being
considered. Since (18a) is seamlessly solvable at all latitudes with ϕ app∗ = 0 at coastlines, the source of the
westward energy flux of mid-latitude RWs in the model
output has been successfully illustrated in the present
study. To our knowledge, this is the first attempt to diagnose the energy cycle of a tropical-extratropical interaction associated with the connection of the equatorial andcoastal waveguides.

Page 12 of 18

Endnotes1
While the energy flux of waves at all latitudes is considered in the present study, the pseudomomentum (orwave-activity) flux of waves at all latitudes is consideredin Aiki et al. (2015, hereafter ATG15). Both the formulations of the present study and ATG15 may be reproduced
even if a spherical coordinate system is used. The useof a Cartesian horizontal coordinate system in both thepresent study and ATG15 is for the purpose of simplicity,which will allow for the results of the two studies to belinked in a future study. A related discussion appears inAppendix 3.
2What we call pressure, energy, and momentum inthe present study are actually dynamic pressure, energy
density, and momentum density, respectively, following
ATG15.
3
dH (n) /dy = 2nH (n−1) , H (n+1) = 2yH (n) − 2nH (n−1) ,H (0) = 1, H (1) = 2y, H (2) = 4y2 − 2, H (3) = 8y3 − 12y,
H (4) = 16y4 − 48y2 + 12.
4
The factor ∂ω/∂k to calculate the energy flux is addedin (14e).
5
The second term in the square brackets of (17b) vanishes as u∗t∗ t∗ ϕ ∗  (−p∗y∗ t∗ t∗ /f ∗ )(p∗ /f ∗ ) = 0 where thephase relationship of plane waves is understood.
6
We use the term “diagnosable” to indicate that thequantity is readily estimated from quantities in model
output without relying on a Fourier analysis.
7
In a related paper, Claus et al. (2014) also used thissolution to investigate the influence of the barotropic
mean flow on the Atlantic equatorial deep jets. The
Atlantic equatorial deep jets are resonant with the gravest
basin mode for a high-order baroclinic mode (typicallythe 15th vertical normal mode) and consist of verticallystacked zonal jets that oscillate at a given depth with aperiod of around 4.5 years.
8
This is lower than the value recommended by G12 for
capturing the observed width of the deep jets but is chosen here since it is not so large as to prevent focusing of
RWs on the equator. In the inviscid solution of Cane andMoore (1981), there is a singularity on the equator at thecenter of the basin due to RW focusing as described by
Schopf et al. (1981).
9
In the present manuscript, we have used the term
“exact” to refer to the level-0 expression, in contrast
to approximate expressions (i.e., level-1 and -2). Likewise, we have used the term “universal” to indicate theAiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 13 of 18

ability to handle all wave types in Table 2, for which thegroup velocity has been well formulated in the literature/textbook.
10 Although it is not in the list of wave types in Table 2,
IGWs on a mid-latitude β-plane may be characterized as
1 corresponds toα
1, δ 2 ≤ 1, γ 2 < 1 where α(19b). Thus, the net content in the square brackets on thelast line of (24c) becomes O(1). Given α in front of c∗ v∗ v∗
on the last line of (24c), we may justify (23d) for IGWs
on a mid-latitude β-plane. It can be said that the righthand side of (24c) becomes significantly nonzero when
the assumption of plane waves in the meridional direction
becomes inconsistent (Anderson and Gill 1979).
11 While the pseudomomentum flux itself (E ∗ − v∗ v∗ ) isdiagnosable from model output, the pseudomomentumflux-based expression of the energy flux (E∗ − v∗ v∗ )ω∗ /k ∗is not easily diagnosable from model output because ofmultiplication by the phase speed (see Appendix 3 fordetails).

Is the streamfunction Eq. (17a) associated with EPV
applicable to mid-latitude waves?

Manipulation of the shallow water equation system (1a)–(1c) yields a characteristic equation associated with themeridional component of velocity to readv∗t∗ t∗ t∗ − c∗2 v∗x∗ x∗ + v∗y∗ y∗ ∗ + f ∗2 v∗t∗ − β ∗ c∗2 v∗x∗ = 0,t(19a)which is applicable to both mid-latitude and equatorial
regions. In what follows, we consider plane waves oneither an f -plane or a mid-latitude β-plane (i.e., f ∗ =f0∗ + β ∗ y∗ and |f0∗ |  |β ∗ y∗ |) and thus assumef ∗2  f0∗2 .

(19b)Then, (19a) may be simplified asv∗t∗ t∗ t∗ − c∗2 v∗x∗ x∗ + v∗y∗ y∗ ∗ + f0∗2 v∗t∗ − β ∗ c∗2 v∗x∗ = 0.
t

(19c)
The Coriolis parameter f0∗ in (19c) is constant that allows
us to assume a horizontally monochromatic wave in acomplex form
(20a)where i is the unit imaginary number, A ∗ is wave amplitude, and θ = k ∗ x∗ + l∗ y∗ − ω∗ t ∗ is wave phase (k ∗and l∗ are the zonal and meridional components of awavenumber vector, respectively, and ω∗ is wave phase).
For simplicity, all A ∗ , k ∗ , l∗ , and ω∗ are assumed to beconstant. Substitution of (20a) to both (1a) and (1c) yields
a solution for u∗ and p∗ to read

(20b)
(20c)

where f ∗ = f0∗ + β ∗ y∗ . On the other hand, substitution of(20a) to (19c) yieldsω∗3 − c∗2 k ∗2 + l∗2 ω∗ − f0∗2 ω∗ − β ∗ c∗2 k ∗ = 0, (21)which is a universal expression for the dispersion relation of the various types of waves in mid-latitude regions.
For example, substitution of β ∗ = 0 to (21) yields
a classical dispersion relation for mid-latitude IGWs
(i.e., waves on an f -plane), and substitution of ω∗2c∗2 k ∗2 to (21) yields a classical dispersion relation for
mid-latitude RWs.
An expression for the zonal component of group velocity may be derived using (21) to read2c∗2 k ∗ ω∗ + β ∗ c∗2
∂ω∗


=
∗
∂k
3ω∗2 − c∗2 k ∗2 + l∗2 − f0∗2=

Appendix 1v∗ = A ∗ eiθ , u∗ = f ∗ ω∗ v∗θ + c∗2 k ∗ l∗ v∗ / ω∗2 − c∗2 k ∗2 ,p∗ = f ∗ k ∗ v∗θ + ω∗ l∗ v∗ c∗2 / ω∗2 − c∗2 k ∗2 ,2c∗2 ω∗2 k ∗ + β ∗ c∗2 ω∗.
2ω∗3 + β ∗ c∗2 k ∗(22a)

We now identify the content of (A∗ B∗ )y∗ in the following
equation:


∂ω∗ 1
p∗2
u∗ p∗ + A∗ B∗ y∗ =u∗2 + v∗2 + ∗2
∗
∂k 2
c

, (22b)where each of A∗ and B∗ are quantities associated withthe set of u∗ , v∗ , p∗ , c∗ , and f ∗ . A first step for investigating (22b) is to decompose u∗ p∗ into two parts: one thatis associated with the numerator of (22a) and one that is
written as the meridional derivative of a scalar quantity, as
follows:
v∗ v∗ f ∗2 + c∗2 l∗2 c∗2 ω∗ k ∗u ∗ p∗ =2ω∗2 − c∗2 k ∗2v∗ v∗ ω∗3 − c∗2 k ∗2 ω∗ − β ∗ c∗2 k ∗ c∗2 k ∗2ω∗2 − c∗2 k ∗2
 ∗f v∗ v∗ y∗ c∗4 k ∗2v∗ v∗ c∗2 ω∗ k ∗− =  ∗22
ω − c∗2 k ∗2ω∗2 − c∗2 k ∗2v∗ v∗ 2c∗2 ω∗ k ∗ + β ∗ c∗2=2 ω∗2 − c∗2 k ∗2
 ∗



f v∗ v∗ y∗ c∗2 ω∗2 + c∗2 k ∗2
−
,
(22c)
2

2 ω∗2 − c∗2 k ∗2
where the first equality has been derived using both(20b)–(20c) and the set of v∗ v∗ = v∗θ v∗θ and v∗θ v∗ = 0 andthe approximate equality in the middle has been derived
using both the dispersion relation (21) and (19b). Then,
we decompose the wave energy in (22b) into two parts,
one that is associated with the denominator of (22a) andAiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 14 of 18

one that is written as the meridional derivative of a scalarquantity. We then have12u∗2 + v∗2 +u∗ p∗ +(p∗ ϕ ∗ /2 + u∗t∗ t∗ ϕ ∗ /β ∗ )y∗

p∗2
c∗2

 ∗2
2
ω − c∗2 k ∗2=2
2 ω∗2 − c∗2 k ∗2v∗ v∗ ω∗2 f ∗2 + c∗4 k ∗2 l∗2 + k ∗2 f ∗2 c∗2 + ω∗2 l∗2 c∗2+2
2 ω∗2 − c∗2 k ∗2

2 


v∗ v∗ ω∗2 − c∗2 k ∗2 + f ∗2 + c∗2 l∗2 ω∗2 + c∗2 k ∗2=2
2 ω∗2 − c∗2 k ∗22v∗ v∗ ω∗2 − c∗2 k ∗22
2 ω∗2 − c∗2 k ∗2
 ∗2


v∗ v∗ ω − c∗2 k ∗2 − β ∗ c∗2 k ∗ /ω∗ ω∗2 + c∗2 k ∗2+2
2 ω∗2 − c∗2 k ∗2



 ∗
f v∗ v∗ y∗ c∗2 k ∗ ω∗2 + c∗2 k ∗2v∗ v∗ ω∗2−
=  ∗22
ω − c∗2 k ∗22ω∗ ω∗2 − c∗2 k ∗2f ∗ v∗ v∗ y∗ 2c∗2 ω∗ k ∗v∗ v∗ 2ω∗3 + β ∗ c∗2 k ∗

 −
=

2 ,
2ω∗ ω∗2 − c∗2 k ∗2
2 ω∗2 − c∗2 k ∗2

where the last equality has been derived using (20a)–(20c).
Equation (22e) may be rewritten as

=

v∗ v∗

ϕ∗ ≡
(22d)

=

−(f ∗ v∗θ v∗θ )y∗ c∗2 [ 1 − 2ω∗3 /(β ∗ c∗2 k ∗ )]2(ω∗2 − c∗2 k ∗2 )[ 2ω∗3 /(β ∗ c∗2 k ∗ ) + 1]−v∗θ,k ∗ + 2ω∗3 /(β ∗ c∗2 )(23b)

has been introduced. The definition of ϕ ∗ , as given by(23b), is based on a Fourier expansion and may be rewritten into an expression which contains none of θ, k ∗ , l∗ , andω∗ to read= q∗ ,(23c)

where the first equality has been derived using (20b) andthe second equality has been derived using (2) [i.e., qt∗∗ =∗ = ω∗ q∗ = −β ∗ v∗ ].
−ω∗ qθ∗ = −β ∗ v∗ and thus −ω∗ qθθθ
As far as we know, the set of (23a) and (23c) has not been
mentioned in previous studies for mid-latitude waves andhas turned out to be almost the same as the set of (17b)and (17a) that has been derived for equatorial waves.
We now consider the meridional flux of wave energy.
We would like to show thatv∗ p∗ − (p∗ ϕ ∗ /2 + u∗t∗ t∗ ϕ ∗ /β ∗ )x∗



0

=

−(f ∗ v∗ v∗ )y∗ c∗2  (4c∗2 ω∗3 k ∗2 + 2β ∗ c∗2 ω∗2 k ∗ )2(ω∗2 − c∗2 k ∗2 )22ω∗3 + β ∗ c∗2 k ∗(−2ω∗5 − 2c∗2 ω∗3 k ∗2 − β ∗ c∗2 ω∗2 k ∗ − β ∗ c∗4 k ∗3 ) +2ω∗3 + β ∗ c∗2 k ∗−(f ∗ v∗θ v∗θ )y∗ c∗2 (β ∗ c∗2 k − 2ω∗3 )2(ω∗2 − c∗2 k ∗2 )(2ω∗3 + β ∗ c∗2 k ∗ )(23a)∇ ∗2 ϕ ∗ − (f0∗ /c∗ )2 ϕ ∗ − (3/c∗2 )ϕt∗∗ t∗ = −β ∗ v∗θ /ω∗

p∗2
∂ω∗ 1
∗2 + v∗2 +− u ∗ p∗u∂k ∗ 2c∗2−(f ∗ v∗ v∗ )y∗ c∗2  (2c∗2 ω∗2 k ∗ + β ∗ c∗2 ω)2k ∗ ω∗2(ω∗2 − c∗2 k ∗2 )22ω∗3 + β ∗ c∗2 k ∗− (ω∗2 + c∗2 k ∗2 )

=

,

where

where the first equality has been derived using both(20b)–(20c) and the set of v∗ v∗ = v∗θ v∗θ and v∗θ v∗ = 0 andthe approximated equality in the middle has been derived
using both the dispersion relation (21) and (19b). The setof (22c) and (22d) allows us to identify the content of(A∗ B∗ )y∗ in (22b) to read

=

∂ω∗ 1
p∗2
∗2 + v∗2 +u∂k ∗ 2
c∗2

∂ω∗ 1
∂l∗ 2

u∗2 + v∗2 +

p∗2
c∗2

.

(23d)

It turns out that the second term on the left hand side,
associated with the additional rotational flux, vanishes
when evaluated using the analytical solution of waves [i.e.,(v∗ v∗ )x∗ = kv∗θ v∗ = 0], which is as in Longuet-Higgins(1964). This is attributed to the assumption of all A ∗ ,k ∗ , l∗ , and ω∗ being constant in particular in the zonal
direction. An expression for the meridional component ofgroup velocity may be derived from (21) to read2c∗2 l∗ ω∗
∂ω∗
=
∗
∂l
3ω∗2 − c∗2 (k ∗2 + l∗2 ) − f0∗2=−[ (f ∗ k ∗ v∗θ + ω∗ l∗ v∗ )c∗2 v∗θ ]y∗= ∗ ∗22k (ω − c∗2 k ∗2 )[ 2ω∗3 /(β ∗ c∗2 k ∗ ) + 1]2c∗2 l∗ ω∗2.
2ω∗3 + β ∗ c∗2 k ∗(24a)Then, we calculate the left hand side of (23d) using (20a)–(20b) as[ (f ∗ ω∗ v∗ + c∗2 k ∗ l∗ v∗ )v∗θ ]y∗ 2ω∗2 /β ∗+ ∗ ∗2 θ ∗2 ∗2
2k (ω − c k )[ 2ω∗3 /(β ∗ c∗2 k ∗ ) + 1]v∗ p∗ =−(p∗ v∗θ )y∗ − (2u∗t∗ t∗ v∗θ )y∗ /β ∗,=2k ∗ [ 1 + 2ω∗3 /(β ∗ c∗2 k ∗ )](22e)v∗ v∗ c∗2 ω∗ l∗,ω∗2 − c∗2 k ∗2

(24b)where v∗θ v∗ = 0 has been used. We now calculate the difference of the meridional component of the group velocity

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

Page 15 of 18

times wave energy and v∗ p∗ using the set of (22d), (24a),and (24b) to yield
∂ω∗ 1
∂l∗ 2

u∗2 + v∗2 +

p∗2
c∗2

− v∗ p∗(f v∗ v∗ )y∗2c∗4 k ∗ l∗ ω∗3(2ω∗3 + β ∗ c∗2 k ∗ ) (ω∗2 − c∗2 k ∗2 )2c∗ v∗ v∗2β ∗ c∗3 k ∗ l∗=−
(2 + β ∗ c∗2 k ∗ /ω∗3 ) ω∗4 (1 − c∗2 k ∗2 /ω∗2 )2c∗2 k ∗ l∗2δ 2
αc∗ v∗ v∗ ,=−(2 + αδ 2 γ )(1 − γ 2 )2 ω∗2
=−

(24c)

Approximate expressions for the energy flux

The exact profile of the group velocity times wave energy
is given by the set of (15a) and (16), which is what wecall the level-0 energy flux. Owing to the last term on
the left hand side of (16) that contains the second-orderpartial differentiation with respect to time, the procedureof inverting EPV, without using a Fourier analysis, is still
complicated.
Hence, we investigate the consequence of artificiallyremoving the second-order time derivative term from (16)as
∇ 2 ϕ app − y2 ϕ app = q,

where the last line has been written using the set of
nondimensional parameters. These are defined asα ≡ β ∗ c∗ /f0∗2 , δ ≡ f0∗ /ω∗ , γ ≡ c∗ k ∗ /ω∗ .

Appendix 2(24d)where the superscript of ϕ app indicates that the solution of(26a) may be regarded as an approximation for the solution ϕ of the accurate streamfunction Eq. (16) associatedwith EPV. We have calculated the meridional profiles ofup + (pϕ app /2 + utt ϕ app )y ,It can be said that the last line of (24c) represents the contribution of higher order terms in an asymptotic expansion based on α, δ, and γ . This contribution should not
be confused with the universal expression of the additional rotational flux which has already been clarified at(23a) and (23d). It should be also noted that the net content within the square brackets on the last line of (24c) is
nondimensional, for which we shall make scale analysis in
the next paragraph.
The quantity αc∗ v∗ v∗ on the last line of (24c) may
be interpreted as a reference for the magnitude of theenergy flux of mid-latitude RWs. Mid-latitude RWs may
be characterized as|αγ | =β ∗ c∗2 /f0∗2≥ 1, δ 2  1, γ 2  1.
|ω∗ /k ∗ |(25a)Thus, the net content within the square brackets on thelast line of (24c) approximates to zero, which justifies(23d) for mid-latitude RWs. On the other hand, for midlatitude IGWs, c∗ v∗ v∗ on the last line of (24c) represents areference for the magnitude of the energy flux. IGWs on
an f -plane may be characterized as
α = 0, δ 2 ≤ 1, γ 2 < 1.

(25b)Thus, the last line of (24c) vanishes, which justifies (23d)for IGWs on an f -plane10 .
To summarize, the streamfunction Eq. (17a) associatedwith EPV and the universal expression of the additionalrotational flux in (17b) applies to both mid-latitude andequatorial waves, in particular for wave types consideredin the present study, as listed in Table 2.

(26a)

(26b)

as shown by the dashed orange lines in Fig. 2 for lowfrequency equatorial waves (e.g., equatorial RWs) and in
Fig. 3 for high-frequency equatorial waves (e.g., equatorial
IGWs). Since this is an analytical investigation, we haveused ϕ app = −vθ /(k − ω3 ) which has been derived fromthe EPV inversion Eq. (26a) with the use of the characteristic Eq. (10). All panels in Fig. 2 show a nice agreement between the dashed orange line given by (26b) andthe solid black line, (∂ω/∂k)(u2 + v2 + p2 ). By contrast,
all panels in Fig. 3 show a finite disagreement between
the dashed orange line given by (26b), up + (pϕ app /2 +utt ϕ app )y , and the solid-black line, (∂ω/∂k)(u2 + v2 + p2 ).
It would be nice if there is a unified approximation
for the energy flux that is able to represent the profileof the group velocity times the energy of both low- andhigh-frequency equatorial waves. We have found that thisrequirement is roughly satisfied if (26b) is simplified as
up + (pϕ app /2)y ,(26c)where ϕ app = −vθ /(k − ω3 ) is the solution of (26a).
The profile of (26c) is shown by the solid blue lines
in Figs. 2 and 3 for low- and high-frequency equatorial waves, respectively. This expression provides what
we think is a potentially useful approximation for thegroup velocity times wave energy (the solid black lines)
for all types of equatorial waves, as we show in the“Methods/Experimental” section.
In the present study, (26b) and its vector and dimensional form (18b) are referred to as the level-1 energy flux.
Likewise, (26c) and its vector and dimensional form (18c)are referred to as the level-2 energy flux.
Why do we appreciate the level-2 energy flux regardless of the error? An expression for pseudomomentum (orAiki et al. Progress in Earth and Planetary Science (2017) 4:11

wave-activity) flux has long been used for the model diagnosis of the direction of the group velocity of waves in
the atmosphere (and also the ocean), including in lowlatitude regions (Ripa 1982; Hoskins et al. 1983; Plumb1986; Haynes 1988; Randel and Williamson 1990; Brunet
and Haynes 1996; Fukutomi and Yasunari 2002; Wakata
and Kitaya 2002; Kawatani et al. 2010). Using the analytical solution of equatorial waves, we have calculated the
profile of the traditional pseudomomentum flux11 times
the phase velocity of waves (see Appendix 3), as shownby the purple dots in Figs. 2 and 3. Interestingly, for lowfrequency waves, the profile of the pseudomomentumflux-based expression (the purple dots) is almost thesame as that of the level-2 energy flux (the blue solidline). On the other hand, for high-frequency waves,the profile of the pseudomomentum-flux-based expression (the purple dots) is similar to that of the level-1
energy flux (the orange dashed line) and quite different
from the exact, level-0 energy flux to which the level-2
energy flux is a better approximation. Thus, the level-2
energy flux is, in general, an improvement on the traditional model diagnosis of group velocity based on thepseudomomentum flux.
Concerning extension to mid-latitude waves, both thelevel-1 and level-2 energy fluxes satisfy all conditions
noted in the last paragraph of the “Boundary conditions
and the connection to mid-latitude regions” section. Notethat the inversion Eq. (18a) of EPV is seamlessly solvable
at all latitudes with the boundary condition of ϕ app∗ =0. To summarize, the set of (18a) and (18c) [together
with the boundary condition (17c)]—what we call thelevel-2 expression—originates from a trade-off between
mathematical exactness and practical accessibility. The
mathematical exactness for retrieving the group velocity
of equatorial waves times wave energy has been achievedby the set of (17a) and (17b)—what we call the level0 expression. However, its accessibility is harmed by thesecond-order time derivative term in the streamfunction
equation (16) associated with EPV. On the other hand,concerning the practical accessibility, the set of (18a) and(18c)—the level-2 expression—has the advantages that (i)it is seamlessly solvable at all latitudes and (ii) it provides
a unified expression for all types of waves with whichto estimate the direction of the group velocity. We havenoted, for equatorial waves, that the profile of the level-2
energy flux is somewhat better than that of the traditional
pseudomomentum flux. It should be also noted that theenergy flux given by (18c) satisfies the boundary condition of no flux through coastlines [using (17c)], an issue
not considered in previous studies for the pseudomomentum flux. With these requirements in mind, we hope that
future studies can lead to either an improved approximation or a numerical algorithm for the level-0 energy
flux.

Page 16 of 18

Appendix 3
Similarity between the level-2 energy flux of this study andthe pseudomomentum flux in previous studies

Ripa (1982) has derived a conservation equation forpseudomomentum (or wave activity) associated withageostrophic waves. His equation may be reproduced
using (1a)–(1c) as∂∂t ∗p∗ u ∗q∗ 2−+∇ ∗ · E∗ − v∗ v∗ , v∗ u∗  = 0,



2β ∗
c∗ 2
IB fluxIB pseudomomentum

(27a)
E∗ ≡

1
p∗ 2
u∗ 2 + v∗ 2 + ∗ 2
2
c

,

(27b)

where the prognostic quantity may be referred to as theimpulse-bolus (IB) pseudomomentum (Aiki et al. 2015,
hereafter ATG15) and E∗ is the wave energy. Note that
the IB pseudomomentum given here is the shallow waterversion of that given by Eq. (27a) in ATG15. It has been
known that the expression of the flux in (27a) can indicate the direction of the group velocity of different types ofwaves, in particular, mid-latitude RWs and IGWs (Hoskins
et al. 1983; Plumb 1986; Haynes 1988). Another nice feature of the IB pseudomomentum Eq. (27a) is that it does
not contain a singularity at the equator. In order to investigate the origin of these features, ATG15 have shown intheir Eq. (18a) an identity between the IB pseudomomentum and the classical energy-based (CE) pseudomomentum to read (again, written here for the shallow waterequations)E∗(ω∗ /k ∗ )  
CE pseudomomentump∗ u ∗
q∗ 2
−
2
2β ∗
c∗




=

IB pseudomomentum ∗ ∗
u q
∂

−

∂y∗

2β ∗

+

∂
∂x∗

 ∗ ∗
v q
,
2β ∗
(28a)

which may be derived from (1a)–(1c) of the present study.
Application of a low-pass temporal filter to (27b), andthen, understanding the phase relationship between v∗ =−qt∗∗ /β ∗ and q∗ yields
q∗ 2
∂
E∗
p∗ u ∗
=
−
− ∗
2
∗
∗
∗
∗
(ω /k )2β∂yc ∗ ∗
u q
.
2β ∗

(28b)

Substitution of (28b) to a low-pass time-filtered version of(28a) yields
∂ ∗
E +
∂t ∗
 ∗ ∗
∂
ω∗ ∗
u q
∗
∗
∗
∗
∗
= 0,
∇ · E −v v , v u + ∗
k∗
∂t
2β ∗

(29)

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

which is a prognostic equation for the wave energy
wherein the zonal component of the flux is proportional
to that in the IB pseudomomentum equation (27a).
It is easy to expect that the expression of the flux in (29)
can indicate the direction of the group velocity of midlatitude RWs and IGWs (Hoskins et al. 1983; Plumb 1986;
Haynes 1988). For equatorial waves, here, we investigate the meridional profile of (E∗ − v∗ v∗ )ω∗ /k ∗ as shownby the purple dots in Figs. 2 and 3 for low- and highfrequency waves, respectively. For low-frequency waves(Fig. 2), the meridional profile of (E∗ − v∗ v∗ )ω∗ /k ∗ (thepurple dots) is almost the same as that of the level-2
energy flux (the blue solid line), showing that the level2 energy flux and the IB flux are closely related. For
high-frequency waves (Fig. 3), the meridional profile of(E∗ − v∗ v∗ )ω∗ /k ∗ (the purple dots) is nearly the same asthat of the level-1 energy flux (the orange dashed line),
indicating that the level-2 energy flux is somewhat better
than the IB flux.
In fact, without relying on the level-0 expression, we
have arrived at the level-2 expression of the energy flux byextending the investigation of ATG15 concerning the algebraic structure of the IB flux (to be explained in a future
study). ATG15 have addressed the importance of a waveinduced scalar quantity and symbolized it as : it vanishes
for mid-latitude IGWs (i.e., waves with no perturbationof EPV) and becomes nonzero for mid-latitude RWs (i.e.,wave with a perturbation of EPV). Here, we suggest that
= (p∗ η∗ )y∗ /2 is closely linked to (p∗ ϕ app∗ )y∗ /2 in thepresent study (η∗ is meridional displacement). This is why
the level-2 expression for the energy flux in the present
study can indicate the direction of the group velocity ofdifferent types of waves, an issue we shall discuss in a
future study.
Note that the IB flux in (27a) has already been used
for the model diagnosis of waves in low-latitude regions(Randel and Williamson 1990; Brunet and Haynes 1996;
Fukutomi and Yasunari 2002; Wakata and Kitaya 2002;
Kawatani et al. 2010). We suggest that, despite the certain
inaccuracy associated with equatorial waves as comparedwith the level-0 expression, the level-2 expression of theenergy flux in the present study will be at least as usefulas the IB flux which has long been used in the atmospheric (and oceanic) literature. For oceanic applications,the level-2 energy flux brings two new advantages over the
IB flux: (i) the level-2 energy flux satisfies a no-normalflux boundary condition at coastlines, and (ii) the wave
energy is a sign-definite quantity while the IB pseudomomentum is not.
Overall, we address the balance of (i) model accessibility, (ii) unified treatment for different types of waves,(iii) mathematical accuracy, and (iv) boundary conditions at coastlines. With these requirements in mind,
we hope future studies can lead to either an improved

Page 17 of 18

approximation or a numerical algorithm for the level-0
energy flux, wherein the profile of the IB flux will providea reference for accuracy because the IB flux has long beenused in previous studies.

Additional file
Additional file 1: Movie of the model experiment. See the caption ofFig. 4 for details. (MP4 2365 kb)
Abbreviations
EPV: Ertel’s potential vorticity; IGW: Inertia gravity wave; KW: Kelvin wave; RGW:
Mixed Rossby-gravity wave; RW: Rossby wave
Acknowledgements
This manuscript has been improved by comments from two anonymous
reviewers. HA thanks Paal Erik Isachsen for the helpful discussions and RJG isgrateful to the GEOMAR for ongoing support.
Funding
This study was supported by JSPS KAKENHI Grant Numbers 26400474 and
15H02129 and also by the Deutsche Forschungsgemeinschaft as part of theSonderforschungsbereich 754 “Climate - Biogeochemistry Interactions in the
Tropical Ocean,” by the German Federal Ministry of Education and Research aspart of the cooperative project SACUS (03G0837A), and by the European
Union 7th Framework Programme (FP7 2007-2013) under grant agreement603521 PREFACE project.
Authors’ contributions
HA proposed the topic and performed the analytical investigation. RJG helped
write the manuscript. MC helped with the numerical investigation. All authors
read and approved the final manuscript.
Competing interests
The authors declare that they have no competing interest.

Publisher’s NoteSpringer Nature remains neutral with regard to jurisdictional claims inpublished maps and institutional affiliations.
Author details
1 Institute for Space-Earth Environmental Research, Nagoya University, Nagoya
City, 464-8601 Aichi, Japan. 2 Application Laboratory, Japan Agency for
Marine-Earth Science and Technology, Yokohama, Japan. 3 GEOMAR
Helmholtz-Zentrum für Ozeanforschung Kiel, Kiel, Germany. 4 Faculty ofMathematics and Natural Sciences, University of Kiel, Kiel, Germany.
Received: 28 June 2016 Accepted: 6 March 2017

References
Aiki H, Richards KJ (2008) Energetics of the global ocean: the role oflayer-thickness form drag. J Phys Oceanogr 38:1845–1869
Aiki H, Takaya K, Greatbatch RJ (2015) A divergence-form wave-inducedpressure inherent in the extension of the Eliassen-Palm theory to athree-dimensional framework for waves at all latitudes. J Atmos Sci
72:2822–2849
Anderson DLT, Gill AE (1979) Beta dispersion of inertial waves. J Geophys Res
84:1836–1842
Ascani F, Firing E, McCreary JP, Brandt P, Greatbatch RJ (2015) The deep
equatorial ocean circulation in wind-forced numerical solutions. J Phys
Oceanogr 45:1709–1734
Brandt P, Funk A, Hormann V, Dengler M, Greatbatch RJ (2011) Interannual
atmospheric variability forced by the deep equatorial Atlantic Ocean.
Nature 473:497–500
Brandt P, Claus M, Greatbatch RJ, Kopte R, Toole JM, Johns WE (2016) Annual
and semi-annual cycle of equatorial Atlantic circulation associated withbasin mode resonance. J Phys Oceanogr 46:3011–3029

Aiki et al. Progress in Earth and Planetary Science (2017) 4:11

Brunet G, Haynes PH (1996) Low-latitude reflection of Rossby wave trains. J
Atmos Sci 53:482–496
Cai M, Huang B (2013) A new look at the physics of Rossby waves: amechanical-Coriolis oscillation. J Atmos Sci 70:303–316
Cane MA, Moore DW (1981) A note on low-frequency equatorial basin modes.
J Phys Oceanogr 11:1794–1806
Chelton DB, Schlax MG (1996) Global observations of oceanic Rossby waves.
Science 272:234–238
Claus M, Greatbatch RJ, Brandt P (2014) Influence of the barotropic mean flowon the width and the structure of the Atlantic equatorial deep jets. J Phys
Oceanogr 44:2485–2497
Claus M, Greatbatch RJ, Brandt P, Toole J (2016) Forcing of the Atlantic
equatorial deep jets derived from observations. J Phys Oceanogr
46:3549–3562
Cummins PF, Oey LY (1997) Simulation of barotropic and baroclinic tides off
northern British Columbia. J Phys Oceanogr 27:762–781
Fukutomi Y, Yasunari T (2002) Tropical-extratropical interaction associatedwith the 10–25-day oscillation over the western Pacific during thenorthern summer. J Meteo Soc Japan 80:311–331
Furuichi N, Hibiya T, Niwa Y (2008) Model-predicted distribution ofwind-induced internal wave energy in the world’s oceans. J Geophys Res
113:C09034
Gill AE (1982) Atmosphere–ocean dynamics. Academic Press, London
Greatbatch RJ, Brandt P, Claus M, Didwischus S-H, Fu Y (2012) On the width ofthe equatorial deep jets. J Phys Oceanogr 42:1729–1740
Haynes PH (1988) Forced, dissipative generalizations of finite-amplitude
wave-activity conservation relations for zontal and nonzonal basic flows. J
Atmos Sci 45:2352–2362
Hoskins BJ, James IN, White GH (1983) The shape, propagation and mean-flow
interaction of large-scale weather systems. J Atmos Sci 40:1595–1612
Isachsen PE, LaCasce JJ, Pedlosky J (2007) Rossby wave instability and apparent
phase speeds in large ocean basins. J Phys Oceanogr 37:1177–1191
Iwasaki T (2001) Atmospheric energy cycle viewed from wave-mean-flow
interaction and Lagrangian mean circulation. J Atmos Sci 58:3036–3052
Johnson GC, Zhang D (2003) Structure of the Atlantic Ocean equatorial deep
jets. J Phys Oceanogr 33:600–609
Kawatani Y, Sato K, Dunkerton TJ, Watanabe S, Miyahara S, Takahashi M (2010)
The roles of equatorial trapped waves and internal inertia-gravity waves indriving the quasi-biennial oscillation. Part II: three-dimensional distributionof wave forcing. J Atmos Sci 67:981–997
Lübbecke JF, Böning CW, Keenlyside N, Xie S-P (2010) On the connection
between Benguela and equatorial Atlantic Ninos and the role of the South
Atlantic Anticyclone. J Geophys Res 115:C09015
Longuet-Higgins MS (1964) On group velocity and energy flux in planetary
wave motion. Deep-Sea Res 11:35–42
Masuda A (1978) Group velocity and energy transport by Rossby waves. J
Oceanogr Soc Jpn 34:1–7
Matsuno T (1966) Quasi-geostrophic motions in the equatorial area. J Meteo
Soc Japan 44:25–43
Matthiessen J-D, Greatbatch RJ, Brandt P, Claus M, Didwischus S-H (2015)
Influence of the equatorial deep jets on the north equatorial
countercurrent. Ocean Dyn 65:1095–1102
McPhaden MJ, Ripa P (1990) Wave-mean flow interactions in the equatorial
ocean. Annu Rev Fluid Mech 20:167–205
Merle J (1980) Annual and interannual variability of temperature in the easternequatorial Atlantic—the hypothesis of an Atlantic El Nino. Oceanol Acta
3:209–220
Nakamura N, Solomon A (2011) Finite-amplitude wave activity and mean flow
adjustments in the atmospheric general circulation. Part II: analysis in theisentropic coordinates. J Atmos Sci 68:2783–2799
Niwa Y, Hibiya T (2004) Three-dimensional numerical simulation of M2 internal
tides in the East China Sea. J Geophys Res 109:C04027
Orlanski I, Sheldon J (1993) A case of downstream baroclinic development
over western north America. Mon Wea Rev 121:2929–2950
Philander SGH (1989) El Nino, La Nina, and the Southern Oscillation. Academic
Press, London
Plumb RA (1986) Three-dimensional propagation of transientquasi-geostrophic eddies and its relationship with the eddy forcing of thetime mean flow. J Atmos Sci 43:1657–1678

Page 18 of 18

Randel WJ, Williamson DL (1990) A comparison of the climate simulated by
the NCAR community climate model (CCM1:R15) with ECMWF analysis. J
Climate 3:608–633
Ripa P (1982) Nonlinear wave-wave interactions in a one-layer reduced-gravity
model on the equatorial β plane. J Phys Oceanogr 12:97–111
Schopf PS, Anderson DLT, Smith R (1981) Beta-dispersion of low-frequency
Rossby waves. Dyn Atmos Oceans 5:187–214
Takaya K, Nakamura H (1997) A formulation of a wave activity flux for stationary
Rossby waves on a zonally varying basic flow. Geophys Res Lett
24:2985–2988
Thierry V, Treguier AM, Mercier H (2004) Numerical study of the annual andsemi-annual fluctuations in the deep equatorial Atlantic Ocean. Ocean
Model 6:1–30
Wakata Y, Kitaya S (2002) Annual variability of sea surface height and upper
layer thickness in the Pacific Ocean. J Oceanogr 58:439–450
Yanai M, Maruyama T (1966) Stratospheric wave disturbances propagatingover the equatorial pacific. J Meteo Soc Japan 44:291–294
Zhai X, Johnson HL, Marshall DP, Wunsch C (2012) On the wind power input tothe ocean general circulation. J Phys Oceanogr 42:1357–1365

Submit your manuscript to ajournal and beneﬁt from:
7 Convenient online submission
7 Rigorous peer review
7 Immediate publication on acceptance
7 Open access: articles freely available online
7 High visibility within the ﬁeld
7 Retaining the copyright to your article

Submit your next manuscript at 7 springeropen.com