#  Beyond Linear Decomposition: a Nonlinear Eigenspace Decomposition for a
Moist Atmosphere with Clouds

Antoine Remond-Tiedrez Department of Mathematics, University of
Wisconsin–Madison, Madison, WI, USA Leslie M. Smith Department of Mathematics
and Department of Atmospheric and Oceanic Sciences, University of
Wisconsin–Madison, Madison, WI, USA Samuel N. Stechmann 22footnotemark: 2

###### Abstract

A linear decomposition of states underpins many classical systems. This is the
case of the Helmholtz decomposition, used to split vector fields into
divergence-free and potential components, and of the dry Boussinesq system in
atmospheric dynamics, where identifying the slow and fast components of the
flow can be viewed as a decomposition. The dry Boussinesq system incorporates
two leading ingredients of mid-latitude atmospheric motion: rotation and
stratification. In both cases the leading order dynamics are linear so we can
rely on an eigendecomposition to decompose states.

Here we study the extension of dry Boussinesq to incorporate another important
ingredient in the atmosphere: moisture and clouds. The key challenge with this
system is that nonlinearities are present at leading order due to phase
boundaries at cloud edge. Therefore standard tools of linear algebra, relying
on eigenvalues and eigenvectors, are not applicable. The question we address
in this paper is this: in spite of the nonlinearities, can we find a
decomposition for this moist Boussinesq system?

We identify such a decomposition adapted to the nonlinear balances arising
from water phase boundaries. This decomposition combines perspectives from
partial differential equations (PDEs), the geometry, and the conserved energy.
Moreover it sheds light on two aspects of previous work. First, this
decomposition shows that the nonlinear elliptic PDE used for potential
vorticity and moisture inversion can be used outside the limiting system where
it was first derived. Second, we are able to rigorously justify, and interpret
geometrically, an existing numerical method for this elliptic PDE. This
decomposition may be important in applications because, like its linear
counterparts, it may be used to analyze observational data. Moreover, by
contrast with previous decompositions, it may be used even in the presence of
clouds.

MSC codes. Primary: 86A10, 76U60; Secondary: 37N10, 53Z05, 65N99.

Keywords: Potential vorticity inversion, clouds, wave-vortical decomposition,
inertia-gravity waves, water vapor, moisture, geophysical fluid dynamics,
phase interface.

###### Contents

  1. 1 Introduction
  2. 2 Background and motivation
    1. 2.1 Helmholtz decomposition
    2. 2.2 Dry Boussinesq and its decomposition
    3. 2.3 Moist Boussinesq
    4. 2.4 Literature on atmospheric dynamics with moisture and clouds, and mathematical theory
    5. 2.5 Simplified notation
  3. 3 Decomposition
  4. 4 Decomposition and slow-fast dynamics
  5. 5 Decomposition and energy
    1. 5.1 PDE-centric approach
    2. 5.2 Energy-centric approach
    3. 5.3 Metric-less approach
    4. 5.4 Summary
  6. 6 Iterative methods for the decomposition and their geometry
    1. 6.1 Convergence guarantees
    2. 6.2 Geometry
  7. 7 Summary and comparison with Helmholtz and the dry case
    1. 7.1 Helmholtz
    2. 7.2 Dry Boussinesq
    3. 7.3 Moist Boussinesq
  8. 8 Conclusion
    1. 8.1 Looking forward
  9. A Appendix
    1. A.1 Helmholtz decomposition
    2. A.2 Useful computations, identities, and facts
    3. A.3 Tools from convex analysis
    4. A.4 Exact nonlinear waves

Note to the reader: in Section 1 below we provide a “shortest path” to the
result, providing just enough background to state and motivate the result. The
main results and conclusions are summarized in Sections 7 and 8. For a more
detailed discussion of the motivation and context we refer the reader to
Section 2.

##  1 Introduction

Decomposing a vector field, or a state vector, into several constitutive
components is ubiquitous in fluid mechanics and atmospheric science. For
example, the Helmholtz decomposition will split a vector field into its
potential and divergence-free components, and it is used throughout fluid
mechanics [Cho68, Tem69, MB01]. Another example is the vortical–wave
decomposition of a state vector which is used in atmospheric and oceanic
sciences to understand the fast rotation and strong stratification regime of
the Boussinesq equations, e.g., [Blu72, Bar95, Maj03]. Crucially, as discussed
in Section 2, both of these decompositions can be viewed from multiple
perspectives and share salient features such as reconstruction of the original
field via elliptic PDE, a rich geometric interpretation of the decomposition
(involving orthogonality and projections), and a connection to slow and fast
components of appropriate dynamical systems.

In this paper we consider a moist Boussinesq system from atmospheric dynamics,
which brings additional realism and complexity due to clouds. We ask the
following question: does this moist, cloudy Boussinesq system also admit such
a decomposition?

While these decompositions have applications in slow–fast dynamical systems,
they can be investigated without reference to the dynamical evolution, by
considering the associated operator. To see this from a general perspective,
consider a dynamical system written abstractly as

$\partial_{t}U+\frac{1}{\epsilon}\mathcal{D}_{0}U+\mathcal{D}_{1}U=0,$

where $\epsilon$ is a small parameter, $\mathcal{D}_{0}$ is the leading-order
operator, and $\mathcal{D}_{1}$ is the next-order contribution to the operator
for the dynamical system. To leading order with respect to $\epsilon$, it is
the operator $\mathcal{D}_{0}$ that contains the essential information.
Understanding the dynamics can be tantamount to understanding the operator
$\mathcal{D}_{0}$. Hence, in the remainder of the paper, we will seldom refer
to the actual dynamical evolution, and instead the main object of interest is
the operator $\mathcal{D}_{0}$.

In past cases, such as the Helmholtz decomposition or Boussinesq equations,
the leading-order operator is linear ($\mathcal{D}_{0}=\mathcal{L}$), and the
decomposition involves a linear operator and linear eigenmodes. On the other
hand, for an atmosphere with clouds, additional nonlinearity arises from phase
changes of water, and $\mathcal{D}_{0}=\mathcal{N}$ is a nonlinear operator.
Hence a substantial challenge in the present paper is from nonlinearity, from
the distinction

$\mathcal{D}_{0}=\mathcal{L}\quad\mbox{versus}\quad\mathcal{D}_{0}=\mathcal{N}.$
(1.1)

For a nonlinear operator, it is not clear if eigenmodes can be identified.
Consequently, a decomposition

$U=U_{1}+U_{2}+\cdots+U_{N_{\text{eig}}}$ (1.2)

into $N_{\text{eig}}$ eigenmodes is likely to be impossible.

Despite the challenges of nonlinearity in (1.1)–(1.2), here we forge ahead
with a decomposition in the spirit of (1.2) as the main aim and main theme of
the paper.

In Section 2 we describe the background and main questions in more detail.
Sections 3, 4, 5, and 6 describe main results, and Sections 7 and 8 present a
summary of main results and conclusions, including links to the main theorems
for easy reference. We point the reader who may skip the background in Section
2 to note Section 2.5 where the simplified notation used in the main body of
the paper is recorded.

##  2 Background and motivation

In this section we first expound on known decompositions (namely the Helmholtz
decomposition and the vortical–wave decomposition of the dry Boussinesq
system) in order to set the stage for the decomposition of moist Boussinesq
states studied in this paper. We also motivate in more detail the context in
which the dry and moist Boussinesq equations arise.

We start with the Helmholtz decomposition for two reasons. First, we present
it here for folks familiar with fluid dynamics but not necessarily familiar
with geophysical and atmospheric fluid dynamics. Second, we present it here
because, as we will detail later, the Helmholtz decomposition plays a central
role in the other decompositions, in both the dry and the moist cases.

###  2.1 Helmholtz decomposition

The Helmholtz decomposition of a vector field $u$ splits it into two
components:

$u=\nabla p+\sigma$

where $\sigma$ is divergence-free. One of the main reasons for the ubiquity of
this decomposition, especially in the realm of fluid mechanics, is that this
single decomposition can be viewed through a variety of lenses. Here we recall
how it can be viewed geometrically, viewed in relation to elliptic PDE, and
viewed in relation to slow-fast dynamical systems. We emphasize these three
perspectives because these are precisely the perspectives that coincide with
each other in the dry case but diverge in the moist case.

We can view the Helmholtz decomposition geometrically: it decomposes the space
of three-dimensional vector fields into two orthogonal pieces, namely
potential vector fields and divergence-free vector fields, and thus each
component of the decomposition can be understood as a projection onto the
corresponding subspace. The Helmholtz decomposition is also related to
elliptic PDE since the scalar potential $p$ solves

$\Delta p=\nabla\cdot u$

while, if we write the solenoidal component as $\sigma=\nabla\times\psi+c$,
the vector potential $\psi$ solves

$\nabla\times(\nabla\times\psi)=\nabla\times u$

with $c$ given by $c=\fint u$. Here we use the notation $\fint f$ to denote
the average over the spatial domain of any scalar or vector function $f$.

Finally we can view the Helmholtz decomposition as a slow-fast decomposition
for an appropriate dynamical system. The relevant system is the low Mach
number limit of the compressible Euler equations (for an ideal isentropic gas)
[KM81, KM82, Maj84]. In this setting the slow dynamics correspond to
incompressible flow whereas the fast dynamics correspond to the acoustic wave
equation. Crucially: the solenoidal component of the Helmholtz decomposition
is precisely the slow piece and the potential component is precisely the fast
piece. We postpone a detailed discussion of this slow-fast dynamics to Section
A.1: this low Mach number limit is well-known, so we only provide details
later once it helps us compare and contrast these classical results
(pertaining to Helmholtz and dry Boussinesq) with the results we have obtained
in this paper.

###  2.2 Dry Boussinesq and its decomposition

Now we turn our attention to the dry Boussinesq system and its decomposition.
As we describe this decomposition we will highlight the properties that it
shares with the Helmholtz decomposition discussed above.

In oceanic and atmospheric application problems, one may need a decomposition
of observational or model data. The latter is often motivated by, or closely
related to, normal-mode analysis of a PDE model such as the dry Boussinesq
equations, e.g. [LM92, BCF14, LCBA20, ZNV+23]. Furthermore, a normal-mode
decomposition is often the starting point for numerical computations [Lei80,
Bar95, SW02, HMRP16, ELS20, VZ21]. For the Boussinesq equations, the normal-
mode analysis is also a wave-vortical decomposition, as further explained
below. Although many numerical computations focused on constant rotation and
buoyancy frequencies in periodic domains [Bar95, SW02, HMRP16], other analyses
also address more realistic PDE systems, such as variable-stratification
Boussinesq dynamics with rigid surface and bottom boundaries [ELS20] and the
shallow water equations on a sphere [Kas78, VZ21].

Let’s stop and discuss this terminology for a second: why is the decomposition
of dry Boussinesq states called a vortical-wave decomposition? It is called
this because the decomposition has two components. The first component is
characterised by the potential vorticity. This is why it is referred to as the
“vortical mode”. (In the dry case this corresponds to an honest-to-goodness
eigenmode. This is no longer true in the nonlinear moist case, which is why we
use the term “component” instead of “mode” in that case.) The second component
corresponds to the fast piece of the decomposition. It can then be shown that
the dynamics obeyed by this fast piece are that of linear oscillations. This
is why this component is known as the “wave mode”. Together the vortical mode
and the wave mode thus constitute the vortical-wave decomposition.

The dry Boussinesq system and its fast rotation and strong stratification
limit is one of the success stories of atmospheric science. In the associated
parameter regime the dry Boussinesq system may be written in terms of
conservation laws:

$\displaystyle\partial_{t}u+u\cdot\nabla u+\frac{1}{\varepsilon}(e_{3}\times
u+\nabla p-\theta e_{3})=0\text{ and }$ (2.1a)
$\displaystyle\partial_{t}\theta+u\cdot\nabla\theta+\frac{1}{\varepsilon}u_{3}=0$
(2.1b)

subject the incompressibility assumption $\nabla\cdot u=0$. Note that for
simplicity physical parameters of order one have be set to unity. Here the
unknowns are the velocity vector field $u$ and the potential temperature
$\theta$ (a scalar field). The first equation corresponds to conservation of
linear momentum while the second equation corresponds to the conservation of
potential temperature. Of course, as written above (2.1b) makes it look like
$\theta$ is _not_ conserved. This is because there is a large background
potential temperature gradient (this is a manifestation of the strong
stratification assumption) and so the quantity $\theta$ used above is actually
the _anomaly_ of the potential temperature about that background gradient.

Here the small parameter indicates that we are particularly interested in the
behaviour of this system as the stratification is strong and the rotation is
fast. In the limit we formally obtain from (2.1b) that the vertical component
of the velocity vanishes and we obtain from the (2.1a) the geostrophic balance
$u_{h}=\nabla^{\perp}_{h}p$ and the hydrostatic balance
$\theta=\partial_{3}p$. Here $v_{h}=(v_{1},\,v_{2},\,0)$ denotes the
horizontal part of a vector $v$ and $v_{h}^{\perp}$ denotes its
$\frac{\pi}{2}$–rotation $v_{h}^{\perp}=(-v_{2},\,v_{1},\,0)$. Note that the
velocity is horizontal since its vertical component vanishes, yet it still
depends on all three spatial coordinates (in this sense the flow is _not_ two-
dimensional). The discussion is formal here, but rigorous theorems are
available (see [EM96] and [ME98] for balanced and general initial data,
respectively).

In particular, the slow dynamics are characterised by the potential vorticity
$PV=\nabla_{h}^{\perp}\cdot u_{h}+\partial_{3}\theta$. This follows from two
observations which hold in the limit: (1) we may compute from (2.1a)–(2.1b)
that $PV$ is purely advected by the flow and (2) the geostrophic and
hydrostatic balances tell us that $p$ solves $\Delta p=PV$. This is the crux
of why the limiting Boussinesq dynamics on the slow timescale, known as the
(dry) quasi-geostrophic equations, are a success story: everything comes down
to a single conserved quantity, the potential vorticity, from which all other
information about the system may be recovered (by way of inverting a
Laplacian). The potential vorticity not only helps us understand the limit, it
also allows us to decompose states of the original dry Boussinesq system into
slow and fast components.

We now discuss this decomposition, also known as the vortical-wave
decomposition. This decomposition splits the state-vector $(u,\,\theta)$ into
two pieces:

$\begin{pmatrix}u\\\ \theta\end{pmatrix}=\begin{pmatrix}\nabla^{\perp}_{h}p\\\
\partial_{3}p\end{pmatrix}+\begin{pmatrix}\sigma_{h}^{\perp}+we_{3}\\\
\sigma_{3}\end{pmatrix}$

where $\sigma$ is divergence-free. The first component is known as the
_vortical mode_ and the second component is known as the _wave mode_. As in
the case of the Helmholtz decomposition above, this decomposition has several
equivalent interpretations. We can view this decomposition geometrically as
projecting into two $L^{2}$–orthogonal subspaces: the first piece is the set
of balanced states satisfying both geostrophic and hydrostatic balance and the
second piece is the set of states with vanishing potential vorticity. This is
similar to the Helmholtz decomposition where the first piece, the potential
component, satisfies the “balance” $u=\nabla p$ while the second piece belongs
to the set of states with vanishing divergence. Note that $L^{2}$ is used
here, not by default, but because it is precisely the energy conserved by the
dry Boussinesq system (2.1a)–(2.1b). We can view this decomposition in
relation to elliptic PDE since $p$ solves

$\Delta p=PV.$

Meanwhile, since $\sigma$ is divergence-free, we can write it in terms of a
vector potential as $\sigma=\nabla\times\psi+c$ for some constant $c$ where
$\psi$ solves

$\nabla\times(\nabla\times\psi)=\begin{pmatrix}j\\\
\partial_{3}w\end{pmatrix}$

for $j=\partial_{3}u_{h}-\nabla_{h}^{\perp}\theta$ denoting the thermal wind
imbalance111 A thermal wind refers to a velocity field which satisfies both
geostrophic and hydrostatic balance. This terminology is warranted by the fact
that such velocity fields are fully determined by the pressure through
geostrophic balance, which in turn is defined by the “thermal” variable,
namely the potential temperature, through hydrostatic balance. Here the
thermal wind _imbalance_ vanishes when both the geostrophic an hydrostatic
balances hold. This imbalance can therefore be thought of as measuring how far
a given state $(u,\,\theta)$ is from satisfying both of these balances. ,
$c=-\fint u_{h}^{\perp}+\fint\theta e_{3}$, and finally $w=u_{3}$. Finally we
can view this decomposition as a slow-fast decomposition for (2.1a)–(2.1b).
Here we describe only cursorily how this may be done, postponing a more
detailed discussion to Section 4. Indeed, the vortical mode corresponds
precisely to solutions of the limiting quasi-geostrophic system while the wave
mode corresponds to the component of the solution with time derivative
proportional to $1/\varepsilon$ (and hence fast). We discuss in more detail
why this decomposition works as a slow-fast decomposition in Section 4, where
we recall how the fast-wave averaging framework allows us to make this idea
precise in order to then use this framework to discuss the slow-fast
decomposition in the moist case.

###  2.3 Moist Boussinesq

We have seen in Sections 2.1 and 2.2 above that both the Helmholtz
decomposition and the vortical-wave decomposition possess several
interpretations, be they geometric, related to elliptic PDE, or about the
slow-fast nature of appropriate dynamical systems. In this section we
introduce the moist Boussinesq system, keeping in mind that our goal is to
find a decomposition for that system which possesses properties akin to those
discussed above.

First, we motivate the moist Boussinesq system itself. Its _dry_ counterpart
is potent for analytical, computational, and predictive purposes because it
focuses on two physical features which dominate mid-latitude atmospheric
dynamics: stratification and rotation (i.e. the Coriolis effect). The _moist_
Boussinesq system introduces another important physical feature: moisture.
Even without incorporating all of the complexities of the water cycle,
moisture plays an important role in atmospheric dynamics because of the role
it plays in the energy budget: as moist air cools, condensing water releases
latent heat energy (and vice-versa as air warms and liquid water evaporates).
When incorporating moisture into the dry system from (2.1a)–(2.1b), we obtain

$\displaystyle\partial_{t}u+u\cdot\nabla u+\frac{1}{\varepsilon}(e_{3}\times
u+\nabla p-\theta e_{3})=0,$ (2.2a)
$\displaystyle\partial_{t}\theta+u\cdot\nabla\theta+\frac{\Gamma_{\theta}}{\varepsilon}u_{3}=C,$
(2.2b) $\displaystyle\partial_{t}q_{v}+u\cdot\nabla
q_{v}-\frac{\Gamma_{q}}{\varepsilon}u_{3}=-C\text{ and }$ (2.2c)
$\displaystyle\partial_{t}q_{l}+u\cdot\nabla q_{l}=C.$ (2.2d)

In comparison to the dry Boussinesq system, the moist Boussinesq system
includes additional evolution equations for the water vapor mixing ratio,
$q_{v}$, in (2.2c) and for the liquid water mixing ratio, $q_{l}$, in (2.2d).
It also includes a source/sink term $C$ for condensation and evaporation, and
associated heating and cooling. Condensation occurs for $C>0$ and represents a
phase transition from vapor to liquid, and evaporation occurs for $C<0$ and
represents a phase transition from liquid to vapor. The constants
$\Gamma_{\theta}$ and $\Gamma_{q}$ correspond to parameters encoding the
strength of the background vertical gradients in potential temperature and
water vapor, respectively. We also note that the precise definition of the
water content $q_{v}$ or $q_{l}$ as a mixing ratio is measured as kilograms of
total water (vapor plus liquid) per kilogram of dry air [KM06, HDMSS13]. While
additional cloud microphysics processes will not be considered here, the moist
Boussinesq system above is valuable for more complex scenarios as well, since
it provides the starting point for extensions that include rainfall, ice,
other precipitation, and other complexities [GS96, KM06, SS17, WSSM19,
WSS+20].

Second, we rewrite the moist Boussinesq system in a way that facilitates the
main goal here of a slow–fast state decomposition. In particular, note that
the source terms $C$ in (2.2b)–(2.2d) appear to potentially cause
complications. However, the system can be rewritten so that the source terms
do not explicitly appear, through a convenient change of variables. This
reformulation has been used in several past studies [Bre87, PS10, HDMSS13,
SS17, MSS19, KSS23] and has facilitated theoretical advances such as an energy
decomposition [SS17, MSS19] and conservation of potential vorticity [KSS23],
even in the presence of clouds and phase changes.

For the convenient way of rewriting (2.2a)–(2.2d), we transform to a different
set of thermodynamic variables that is conserved. In particular, define the
equivalent potential temperature, $\theta_{e}$, and total water mixing ratio,
$q_{t}$, as

$\displaystyle\theta_{e}=\theta+q_{v},$ (2.3a) $\displaystyle
q_{t}=q_{v}+q_{l}.$ (2.3b)

Note that the equivalent potential temperature $\theta_{e}$ may be understood
physically as the potential temperature that a parcel of air would have if it
were heated by converting all of its water vapour into liquid water. The
evolution equations of $\theta_{e}$ and $q_{t}$ can be found by taking
appropriate linear combinations of (2.2b)–(2.2d), and they take the form

$\displaystyle\partial_{t}\theta_{e}+u\cdot\nabla\theta_{e}+\frac{1}{\epsilon}u_{3}=0,$
(2.4a) $\displaystyle\partial_{t}q_{t}+u\cdot\nabla
q_{t}-\frac{1}{\epsilon}u_{3}=0,$ (2.4b)

where we have now set the background gradient parameters from (2.2b) and
(2.2c) to be $\Gamma_{\theta}=2$ and $\Gamma_{q}=1$ without loss of generality
within the stably stratified regime, and for simplicity. The important feature
of (2.4a)–(2.4b) is that the source term of condensation and evaporation, $C$,
has been eliminated. As a result, (2.4a)–(2.4b) show that
$\theta_{e}+(z/\epsilon)$ and $q_{t}-(z/\epsilon)$ (which combine the
anomalies $\theta_{e}$ and $q_{t}$ with their background profiles $z/\epsilon$
and $-z/\epsilon$, respectively) are conserved along fluid parcel
trajectories. Physically, the equivalent potential temperature is conserved
because losses of water vapor $q_{v}$ are compensated by gains in heat
$\theta$, as indicated by the definition in (2.3a); and the total water mixing
ratio, $q_{t}=q_{v}+q_{l}$, is conserved because losses of water vapor $q_{v}$
are compensated by gains in liquid water $q_{l}$.

To complete the rewriting, the buoyancy term $\theta/\epsilon$ from (2.2a)
must also be rewritten in terms of $\theta_{e}$ and $q_{t}$. To do so, the
following transformation is used:

$\displaystyle\theta=\theta_{e}-\min(q_{t},q_{vs}),$ (2.5a) $\displaystyle
q_{v}=\min(q_{t},q_{vs}),$ (2.5b) $\displaystyle q_{l}=\max(0,q_{t}-q_{vs}),$
(2.5c)

which is the reverse of the transformation in (2.3a)–(2.3b). The saturation
water vapor value $q_{vs}$ appears here for distinguishing between the
unsaturated phase ($q_{v}<q_{vs}$), when no cloud is present, and the
saturated phase ($q_{v}=q_{vs}$), when a cloud is present. Note that $q_{vs}$
is not appearing here for the first time. It is already implicit in the
condensation/evaporation source term $C$ in (2.2b)–(2.2d). The value of
$q_{vs}$ will be constant here for simplicity. For explicit expressions for
$\theta$ in each phase, one can rewrite (2.5a) as

$\theta=\left\\{\begin{array}[]{ll}\theta_{e}-q_{t}&\mbox{if}\quad
q_{t}<q_{vs},\\\\[4.0pt] \theta_{e}-q_{vs}&\mbox{if}\quad q_{t}\geq
q_{vs},\end{array}\right.$ (2.6)

which can be used for another expression of the buoyancy in terms of the
variables $\theta_{e}$ and $q_{t}$.

The end result of the rewriting is the moist Boussinesq system in the
following form:

$\displaystyle\partial_{t}u+u\cdot\nabla u+\frac{1}{\varepsilon}(e_{3}\times
u+\nabla p-(\theta_{e}-{\min}_{0}\,q_{t})e_{3})=0,$ (2.7a)
$\displaystyle\partial_{t}\theta_{e}+u\cdot\nabla\theta_{e}=-\frac{1}{\varepsilon}u_{3}\text{
and }$ (2.7b) $\displaystyle\partial_{t}q_{t}+u\cdot\nabla
q_{t}=\frac{1}{\varepsilon}u_{3},$ (2.7c)

subject to $\nabla\cdot u=0$, and following from (2.2a), (2.4a), (2.4b), and
(2.5a). Some additional specifications are as follows. The notation
$\min_{0}q_{t}=\min(0,\,q_{t})$ has been introduced in (2.7a), and the
saturation water vapor parameter $q_{vs}$ has been set to zero without loss of
generality, and for simplicity. Throughout this paper we work on the three-
dimensional torus (i.e. with triply-periodic boundary conditions). Also recall
that we set the background gradient parameters from (2.2b) and (2.2c) to be
$\Gamma_{\theta}=2$ and $\Gamma_{q}=1$ without loss of generality in the
stably stratified regime, and for simplicity. One can now see the value in
rewriting the moist Boussinesq system: while the $(\theta,q_{v},q_{l})$
formulation in (2.2a)–(2.2d) is helpful for seeing the additional cloud
physics that is added to the dry system, the $(\theta_{e},q_{t})$ formulation
in (2.7a)–(2.7c) is helpful for writing the moist system in a way that is
mathematically similar to the dry system. The main differences between the
moist and dry systems are that the moist system has two thermodynamic
variables ($\theta_{e}$ and $q_{t}$) instead of one ($\theta$), and the
buoyancy term in the moist case involves the nonlinear expression
$\theta_{e}-{\min}_{0}\,q_{t}$.

What makes the moist Boussinesq system particularly interesting is that
moisture introduces _phase boundaries_. In concrete terms: here we have
assumed without loss of generality that water saturation occurs at $q_{t}=0$
and so different physics are in play depending on the sign of $q_{t}$, i.e.
depending on whether water is above or below saturation (this is discussed in
more detail in Section 2). We therefore refer to the set
$\left\\{q_{t}=0\right\\}$ as the _phase boundary_ since it delineates the
separation between water in vapour form (outside of a cloud) and water in
vapour-plus-liquid form (inside of a cloud).

Finally, to end this subsection, we review some properties of the moist
Boussinesq system that are potentially relevant to the aim here—i.e., to the
search for a state decomposition.

As $\varepsilon\to 0$, solutions of the moist Boussinesq system (2.7a)–(2.7c)
are forced to satisfy two _balances_ , namely geostrophic balance
$u=u_{h}=\nabla_{h}^{\perp}p$ and hydrostatic balance
$\theta_{e}-\min_{0}q_{t}=\partial_{3}p$, which come from the fast rotation
and strong stratification, respectively. This warrants the following
definition.

######  Definition 2.1 (State space and balanced set).

We denote by $\mathbb{L}_{\sigma}^{2}\vcentcolon={(L^{2})}^{3}_{\sigma}\times
L^{2}\times L^{2}$ the state space of the moist Boussinesq system
(2.7a)–(2.7c). Here $u\in{(L^{2})}^{3}_{\sigma}$ means that $u\in
L^{2}(\mathbb{T}^{3};\mathbb{R}^{3})$ is divergence-free. We then define the
_balanced set_

$\mathcal{B}\vcentcolon=\left\\{(u,\,\theta_{e},\,q_{t})\in\mathbb{L}^{2}_{\sigma}:u_{h}=\nabla_{h}^{\perp}p,\,u_{3}=0,\,\text{and
}\theta_{e}-{\min}_{0}\,q_{t}=\partial_{3}p\text{ for some
}p\in\mathring{H}^{1}\right\\}$ (2.8)

where $\mathring{H}^{1}$ is the set of scalar fields in $H^{1}$ with vanishing
averages.

We note that a state belongs to the balanced set precisely when it satisfies
both geostrophic and hydrostatic balance.

Proceeding once again by analogy with the dry case, one may now identify _two_
slow quantities: the _equivalent_ potential vorticity
${PV}_{e}=\nabla_{h}^{\perp}\cdot u_{h}+\partial_{3}\theta_{e}$, defined now
in terms of the _equivalent_ potential temperature, and the moist variable
$M=\theta_{e}+q_{t}$ [SS17]. Similar to the dry case, now $p$ and $M$ fully
determine a balanced state by way of the ${PV}_{e}$-and-$M$ inversion

$\Delta p+\frac{1}{2}\partial_{3}{\min}_{0}\,(M-\partial_{3}p)={PV}_{e},$
(2.9)

which is obtained as before by inserting the balances into the definition of
${PV}_{e}$. This is a nonlinear elliptic PDE whose well-posedness was
established recently [RTSS24]. Using this PDE for ${PV}_{e}$-and-$M$
inversion, one can then phrase the limiting moist (and potentially
_precipitating_) quasi-geostrophic system as a nonlinear transport equation
for ${PV}_{e}$ and $M$ [SS17].

It is important to note that the nonlinear elliptic PDE in (2.9) only has
meaning for the _balanced states_ from (2.8), at this point in the story. This
is in contrast to the cases of the Helmholtz decomposition and dry Boussinesq
system from Sections 2.1 and 2.2, respectively, for which an elliptic PDE is
part of a state decomposition for _general_ states. One aim of the present
paper is to investigate whether the nonlinear ${PV}_{e}$-and-$M$ inversion in
(2.9) may have additional meaning as part of a state decomposition for
_general_ states.

The question is now this: can we find a decomposition for the moist Boussinesq
system which has the same favourable properties as the Helmholtz and
vortical–wave decompositions?

###  2.4 Literature on atmospheric dynamics with moisture and clouds, and
mathematical theory

In this section we review some literature on topics related to the main focus
of the paper, i.e., state decomposition for a moist atmosphere with clouds. We
mention both atmospheric sciences literature and mathematics literature.

Decompositions have been useful for numerous applications in atmospheric and
science, and they are commonly based on decomposing the dry portions
$(u,\,\theta)$ of the atmospheric state. Applications include analysis of
observational or model data [HMR85, ŽKT+15], data assimilation [Lor81, Dal93,
DB99, ZGK04, Ban21], and nonlinear wave interactions [GG99, MB03, BM04,
RDTM08, FKN09]. Other references are also described at the beginning of
Section 2.2. One aim here is to extend the ideas of dry decompositions to the
realm of a moist atmosphere with clouds.

For a moist atmosphere with clouds, decompositions can potentially provide
valuable insight into the behavior of rainfall and precipitation [KMS12]. Past
work has commonly used linearisation, and either dry or moist eigenmodes. By
adapting or modifying the dry linear eigenmodes, some studies have analyzed
observational and model data for understanding moisture-coupled phenomena such
as convectively coupled equatorial waves (CCEWs) [YHS03, GK12, OS16, MC18,
KGK+22], tropical intraseasonal oscillations such as the Madden–Julian
Oscillation (MJO) [ŽF15], and the Walker circulation [SO14]. Moist eigenmodes,
obtained by linearisation of moist dynamics, have also been used for analyzing
observational and model data for the MJO [OS15, SM15], and for asymptotic
analysis of weakly nonlinear wave interactions [CMS15, CMS16]. Balanced and
unbalanced moisture decompositions have been investigated with a type of
linearisation that applies with clouds present [WSSM19, WSS+20]. Applications
to data assimilation have also investigated linear eigenmodes that decompose
moisture into its balanced and unbalanced components [OWSS24].

Here we aim for a decomposition that does not rely on linearisation, so that
it may potentially include a more faithful representation of the nonlinearity
and the physics associated with clouds and latent heating.

The PDE community has recently provided rigorous treatment of moist
atmospheric dynamics with phase changes. A model with water vapor and liquid
water was analyzed in [CZT12, CZFTT13, BCZT14], where the velocity field is
prescribed. Then followed results for two-phase flow with dynamical velocity
evolving via the primitive equations [CZHK+15, LM20, TW15, TW16]. An
additional category of water was included in [CHT+18, TL22, HKLT17, HKLT20],
such that liquid water is divided into cloud water (that does not fall) and
precipitating water. The work on warm-rain microphysics flow may be separated
by the microphysics model determining how water is partitioned between vapor,
cloud water and rain water, as well as by whether the velocity is prescribed
or governed by the primitive equations. The case of prescribed velocity is
considered in [CHT+18] for Grabowski microphysics [Gra98], and in [HKLT17] for
Klein-Majda microphysics [KM06]. For primitive-equation velocity, theorems for
the Grabowski and Klein-Majda models are found, respectively, in [TL22] and
[HKLT20]. Incorporation of ice into water microphysics is discussed in
[CJTT21].

We note that the aforementioned rigorous studies of atmospheric flows with
phase changes treat the velocity as either prescribed, or as governed by the
primitive equations. Following [CT07], the barotropic-baroclinic decomposition
of the velocity is the basis for analysis of the primitive equations and its
extensions. Recently, the current authors obtained existence, uniqueness and
regularity results [RTSS24] for the nonlinear elliptic PDE (2.9) underlying
Boussinesq velocity dynamics in the limit of asymptotically fast rotation rate
and buoyancy frequencies [SS17, ZSS22]. Here, we show that the same nonlinear
elliptic PDE is the foundation for decomposition of the Boussinesq equations
more generally, in the regime of finite parameter values.

Pertaining to the dry Boussinesq evolution equations in the limit of
asymptotically fast rotation rate and buoyancy frequency, a series of rigorous
results established the decoupling between slow, balanced motions and fast,
unbalanced waves [EM96, BANZ97, EM98, ME98, Maj03]. For moist atmospheric
dynamics with phase changes, first steps to use the fast-wave averaging
framework [EM96, EM98, ME98, Maj03] appear in [ZSS21b, ZSS21a]. The latter
references present formal asymptotic analysis [ZSS21b] and supporting
numerical analyses [ZSS21a]. The inclusion of boundary layers and the diabatic
layer have been investigated in formal asymptotic analysis as well [KSPK22,
BHK23]. We note that rigorous fast-wave-averaging analysis for multi-phase
atmospheric flows remains open. In Section 4, we explain the path forward for
rigorous fast-wave-averaging analysis using the decomposition defined in
Section 3. To date, the evidence suggests that nonlinear oscillations lead to
coupling between fast and slow flow components, possibly even in the limiting
dynamics when distinct phases co-exist (see Figure 13 in the appendix and
[ZSS21a]).

###  2.5 Simplified notation

To simplify the notation used throughout the rest of the paper, we will from
now on use $\theta$ to denote the equivalent potential temperature
$\theta_{e}$ (which is not an issue since the potential temperature, denoted
$\theta$ in the text above, will not make any appearances in the sequel).
Similarly we will write $PV$ for the equivalent potential temperature
${PV}_{e}$ and $q$ to mean $q_{t}$, the total water content.

##  3 Decomposition

As discussed in Sections 1 and 2, we seek a decomposition adapted to the
leading-order operator of the moist Boussinesq system. We can see from
(2.7a)–(2.7c) that this operator is

$\mathcal{N}(u,\,\theta,\,q)=\begin{pmatrix}\mathbb{P}_{L}\left(u_{h}^{\perp}-(\theta-{\min}_{0}\,q\right)e_{3}\\\
u_{3}\\\ -u_{3}\end{pmatrix}$ (3.1)

where we note again that $\theta_{e}$ and $q_{t}$ are written as $\theta$ and
$q$ for notational simplicity, as specified in Section 2.5. Here
$\mathbb{P}_{L}$ denotes the Leray projector onto the space of divergence-free
vector fields. On the three-dimensional torus it takes the simple form
$\mathbb{P}_{L}u=-\nabla\times(\Delta^{-1}\nabla\times u)$.

Given that the operator $\mathcal{N}$ is nonlinear, it is likely impossible to
define and characterize a full set of the nonlinear eigenfunctions and
eigenvalues. Hence it is likely impossible to define a full eigenfunction
decomposition such as in (1.2). An alternative route to a decomposition must
be sought.

Here we propose to use the notions of null set and image as concepts that may
allow us to extend from the linear to the nonlinear realm. We propose a
decomposition that will split any state in the state space
$\mathbb{L}^{2}_{\sigma}$ into two pieces: one piece residing in the null set
of $\mathcal{N}$ and another piece residing in the image of $\mathcal{N}$. In
other words we will perform the decomposition

$\mathbb{L}^{2}_{\sigma}=\left\\{\mathcal{N}=0\right\\}+\operatorname{im}\mathcal{N}.$
(3.2)

If the operator $\mathcal{N}$ were _linear_ , then this decomposition would be
a (coarse) eigendecomposition, into the eigenspace of the zero eigenvalue and
the eigenspace of all nonzero eigenvalues. Such eigenspace ideas are helpful
as motivation for pursuing the null set and image as fundamental concepts for
defining a decomposition here. However, the operator $\mathcal{N}$ is
nonlinear, so the full meaning of (3.2) is not immediately clear. We aim to
demonstrate in this section that (3.2) is able to characterize the desired
decomposition in the nonlinear realm.

Note that in this section we will only present and discuss this decomposition
itself. We will not discus the various interpretations of this decomposition.
This will be done in later sections. In Section 4 we discuss how this
decomposition is related to slow-fast dynamics, in Section 5 we discuss how
this decomposition relates to the conserved energy and projections, and
finally in Section 6 we discuss iterative methods used to compute this
decomposition in practice and their geometric interpretation.

Since we wish to perform the decomposition foreshadowed in (3.2) we must begin
by understanding both the null set and the image of $\mathcal{N}$ as well as
possible. We first turn our attention to the null set. As shown in the result
below, this null set is precisely the set of states consisting of states in
both hydrostatic and geostrophic balance. This is why the null set of
$\mathcal{N}$ will be referred to, from now on, as the _balanced set_.

######  Proposition 3.1 (Alternate characterisations of the balanced set).

For any state $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ the following are
equivalent.

  1. 1.

$\mathcal{X}$ is balanced, i.e., $\mathcal{N}(u,\,\theta,\,q)=0$.

  2. 2.

$\mathcal{X}=(u,\,\theta,\,q)$ satisfies

$\partial_{3}u_{h}-\nabla_{h}^{\perp}(\theta-{\min}_{0}\,q)=0,\,u_{3}=0,\,\fint
u_{h}=0,\text{ and }\fint\theta-{\min}_{0}\,q=0.$ (3.3)

  3. 3.

$\mathcal{X}$ may be parameterised as $\mathcal{X}=\Phi(p,\,M)$ for some
unique $p\in\mathring{H}^{1}$ and $M\in L^{2}$, where

$\Phi(p,\,M)=\begin{pmatrix}\nabla^{\perp}_{h}p\\\
\partial_{3}p+\frac{1}{2}{\min}_{0}\,(M-\partial_{3}p)\\\
M-\partial_{3}p-\frac{1}{2}{\min}_{0}\,(M-\partial_{3}p)\end{pmatrix}.$

Note that this proposition characterizes the balanced set in three ways:

  * •

it is the null set of $\mathcal{N}$,

  * •

it is characterised as the null set of the three _measurements_ $j$, $a$, and
$w$ given by

$\displaystyle
j(u,\,\theta,\,q)=\partial_{3}u_{h}^{\perp}-\nabla_{h}^{\perp}(\theta-{\min}_{0}\,q),\,w(u,\,\theta,\,q)=u_{3}$
$\displaystyle\text{ and }a(u,\,\theta,\,q)=-\fint
u_{h}^{\perp}+\fint(\theta-{\min}_{0}\,q)e_{3},$

and

  * •

it may be parameterised via the _coordinates_ $p$ and $M$.

The measurement $j$ is commonly referred to as the _thermal wind imbalance_
since it quantifies how far the state is from obeying both geostrophic and
hydrostatic balanced. Further down we will see that the image of $\mathcal{N}$
also admits three comparable descriptions. The measurement $a$ depends on some
spatial averages of the state $(u,\,\theta,\,q)$. Why it depends on these
particular averages (e.g. it depends on the horizontal velocity average but
not the vertical velocity average) and why these averages are combined in this
way to constitute $a$ is explained in Remark 3.2.

###### Proof of Proposition 3.1.

It will be particularly convenient to consider the new unknown

$v\vcentcolon=-u_{h}^{\perp}+(\theta-{\min}_{0}\,q)e_{3}.$

Indeed, we see that using this unknown we may rephrase items 1, 2, and 3 in a
much simpler way. The first item may now be written as

$v=\nabla p\text{ and }u_{3}=0\text{ for some }p\in\mathring{H}^{1}$

while the second item may be written, using Lemma A.5, as

$\nabla\times v=0,\,\fint v=0,\text{ and }u_{3}=0.$

It then follows directly from the Helmholtz decomposition (see Corollary A.2)
that items 1 and 2 are equivalent. Finally we note that, by virtue of Lemma
A.6, item 3 is equivalent to

$\left\\{\begin{aligned} &u=\nabla_{h}^{\perp}p,\\\ &\theta+q=M,\text{ and
}\\\ &\theta-{\min}_{0}\,q=\partial_{3}p,\end{aligned}\right.$

which implies that $v=\nabla p$ and that $u_{3}=0$, from which item 1 follows.
∎

######  Remark 3.2 (The “good unknown”).

The so-called “good-unknown” $-u_{h}^{\perp}+(\theta-{\min}_{0}\,q)e_{3}$ will
appear throughout this paper. This is because this three-dimensional vector
field, which combines the horizontal velocity and the buoyancy, turns out to
be the one that needs to be decomposed according to a three-dimensional
Helmholtz decomposition. This contrasts with the approach commonly used in the
literature where instead of a three-dimensional Helmholtz decomposition of
this good unknown (or of its dry counterpart, the good unknown
$-u_{h}^{\perp}+\theta e_{3}$), the decomposition carried out is a two-
dimensional Helmholtz decomposition of only the horizontal velocity. It is
precisely because the typical decomposition is two-dimensional (even though
the problem is three-dimensional) that vertically sheared horizontal flows
often require special treatment. When an honest-to-goodness three-dimensional
decomposition is used, these vertically sheared horizontal flows naturally
accounted for (see Remark 3.10 where this is discussed in more detail).

Note that the role played by the good unknown in the good decomposition is
discussed in more detail in Remark 3.8. For now, in order to illustrate why it
is such a convenient unknown to work with, we record a couple of computations
from the dry case. In the dry case, the good unknown is

$-u_{h}^{\perp}+\theta e_{3}.$

We then observe from Lemma A.5 that its divergence is precisely the potential
vorticity since

$\nabla\cdot(-u_{h}^{\perp}+\theta e_{3})=-\nabla_{h}^{\perp}\cdot
u_{h}+\partial_{3}\theta=-PV$

while its curl is precisely the thermal wind imbalance since

$\nabla\times(-u_{h}^{\perp}+\theta
e_{3})=\partial_{3}u_{h}-\nabla_{h}^{\perp}\theta.$

The physical relevance of thermal wind imbalance is that it vanishes precisely
when the (dry) state $(u,\,\theta)$ is in both hydrostatic and geostrophic
balance (provided that both $u_{h}$ and $\theta$ have vanishing averages).
Once again: Remark 3.8 discusses in more detail the precise role that the
_moist_ good unknown $-u_{h}^{\perp}+(\theta-{\min}_{0}\,q)e_{3}$ plays in the
decomposition, but already the computations above, from the dry case, indicate
that essential quantities may readily be computed from the (dry) good unknown.

Now that we understand the null set quite well we turn our attention to the
image of $\mathcal{N}$. This set will actually be called the _wave set_. This
is because the dynamics induced on that set by the moist Boussinesq system are
oscillatory (in the dry case, the dynamics induced on the dry analog of the
wave set by the dry Boussinesq system are precisely the dynamics of linear
oscillations, i.e. waves).

######  Proposition 3.3 (Characterisation of the image of $\mathcal{N}$).

For any state $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ the following are
equivalent.

  1. 1.

$\mathcal{X}$ belongs to the image of $\mathcal{N}$.

  2. 2.

$\mathcal{X}=(u,\,\theta,\,q)$ satisfies

$PV=\nabla_{h}^{\perp}\cdot u_{h}+\partial_{3}\theta=0\text{ and
}M=\theta+q=0.$ (3.4)

  3. 3.

$\mathcal{X}$ may be parameterised as $\mathcal{X}=\Psi(\sigma,\,w)$ for some
unique $\sigma\in H^{1}_{\sigma}$ and $w\in L^{2}$ satisfying
$\nabla_{h}^{\perp}\cdot\sigma_{h}=\partial_{3}w$, where

$\Psi(\sigma,\,w)=\begin{pmatrix}\sigma_{h}^{\perp}+we_{3}\\\ \sigma_{3}\\\
-\sigma_{3}\end{pmatrix}.$

In either case we say that $\mathcal{X}\in\mathcal{W}$ for the _wave set_
$\mathcal{W}$. Why this set warrants the name “wave set” will be discussed in
Section 4 below when we discuss slow-fast decompositions in more detail.

Note that, as was done for the balanced set above, Proposition 3.3
characterises the wave set in three ways:

  * •

it is the image of $\mathcal{N}$,

  * •

it is characterised as the null set of the two measurements $PV$ and $M$, and

  * •

it may be parameterised via the coordinates $\sigma$ and $w$.

###### Proof.

First we suppose that item 1 holds, i.e. that

$\left\\{\begin{aligned}
&u=\mathbb{P}_{L}\left(v_{h}^{\perp}+(\phi-{\min}_{0}\,r)e_{3}\right),\\\
&\theta=v_{3},\text{ and }\\\ &q=-v_{3}\end{aligned}\right.$

for some $(v,\,\phi,\,r)\in\mathbb{L}^{2}_{\sigma}$. This means that
$u=v_{h}^{\perp}+(\phi-{\min}_{0}\,q)e_{3}+\nabla\pi$ for some
$\pi\in\mathring{H}^{1}$ and so

$\nabla_{h}^{\perp}\cdot
u_{h}+\partial_{3}\theta=\nabla_{h}^{\perp}\cdot(v_{h}^{\perp}+\nabla_{h}\pi)+\partial_{3}v_{3}=\nabla\cdot
v=0$

while clearly $\theta+q=0$. Thus item 2 holds.

Now suppose that item 2 holds. We may introduce the new unknown
$v=-u_{h}^{\perp}+\theta e_{3}$ such that now item 2 reads

$\nabla\cdot v=0\text{ and }\theta+q=0$

(by virtue of Lemma A.5). In other words $v=\sigma$ for some divergence-free
$\sigma$ while $\theta=-q$, and so

$\left\\{\begin{aligned} &u=\sigma_{h}^{\perp}+u_{3}e_{3},\\\
&\theta=\sigma_{3},\text{ and }\\\ &q=-\sigma_{3}.\end{aligned}\right.$

Therefore, for $w\vcentcolon=u_{3}$, we see that item 3 holds as desired since
$\nabla_{h}^{\perp}\cdot\sigma_{h}-\partial_{3}w=-\nabla\cdot u=0$.

Finally we suppose that item 3 holds. Then item 1 follows immediately from
choosing $v=\sigma$, $\phi=w$, and $r=0$ since $\sigma_{h}^{\perp}+we_{3}$ is
divergence-free and so the Leray projection of this vector field is equal to
itself. ∎

At this stage we have a good understanding of both the balanced set (which is
the null set of $\mathcal{N}$) and of the wave set (which is the image of
$\mathcal{N}$). We now seek to leverage this understanding to produce the
decomposition formally foreshadowed in (3.2).

In order to do this in a manner that will be useful later, especially when it
comes to dynamics, we need to discuss dynamics briefly. Recall that, at
leading order, the dynamics of interest are those governed by
$\partial_{t}+\mathcal{N}=0$. So how do the null set and image of
$\mathcal{N}$ relate to this dynamical system?

Well, note that Proposition 3.1 tells us that producing _balanced_ states is
easy: it suffices to specify $p$ and $M$ since then a balanced state may be
reconstructed from these two quantities. If the leading operator $\mathcal{N}$
were linear, then this would be enough since balanced solutions of
$\partial_{t}+\mathcal{N}$ would then be guaranteed to be slow, meaning that
they would have vanishing time derivatives. However, due to the presence of
moisture and phase boundaries, $\mathcal{N}$ is _not_ linear. This means that
we must _additionally_ impose that balanced solutions be slow. In more prosaic
terms: since any choice of the coordinates $p$ and $M$ produces a balanced
state, we must make sure that these coordinates are carefully chosen to also
be _slow_.

This is where our understanding of the image of $\mathcal{N}$, i.e. of the
wave set, comes in handy. Indeed: we have characterised this wave set as the
set of states whose $PV$ and $M$ vanish. This is another way of saying that
$PV$ and $M$ are precisely a complete set of _time-dependent_ slow
measurements! (Any constant measurement would be trivially slow.)

We have found our solution: to guarantee that our balanced component
determined by $p$ and $M$ is also slow, it must be reconstructed solely from
the slow measurements, $PV$ and $M$. If we think about this geometrically,
this is really a matter of _transversality_ : we need the wave set
$\mathcal{W}$ to be transverse to the balanced set $\mathcal{B}$ so that
balanced states in $\mathcal{B}$ may be fully characterised by the
measurements annihilating the wave set $\mathcal{W}$. This is taken care of in
the result below.

We can think about why the balanced set and the wave set need to be transverse
another way. Recall that we are after a decomposition of the form (3.2). This
decomposition is only unique if its two components are transverse, otherwise
the same state could be decomposed in several different ways.

######  Lemma 3.4 ($PV$ and $M$ characterise balanced states).

Given $PV\in H^{-1}$ and $M\in L^{2}$, where $H^{-1}$ denotes the dual of
$\mathring{H}^{1}$, there exists a unique $\mathcal{X}\in\mathcal{B}$ for
which $PV(\mathcal{X})=PV$ and $\mathcal{M}(\mathcal{X})=M$.

###### Proof.

By virtue of Proposition 3.1 it suffices to show that there exist a unique
$p\in\mathring{H}^{1}$ and a unique $\widetilde{M}\in L^{2}$ for which
$(PV\circ\Phi)(p,\,\widetilde{M})=PV$ and
$(\mathcal{M}\circ\Phi)(p,\,\widetilde{M})=M$. The latter identity tells us
that $\widetilde{M}=M$, and since

$(PV\circ\Phi)(p,\,M)=\Delta
p+\frac{1}{2}\partial_{3}{\min}_{0}\,(M-\partial_{3}p)$

the result follows from the uniqueness of a solution to nonlinear $PV$-and-$M$
inversion as proved in [RTSS24]. ∎

Now that we have verified that the wave set and the balanced set are
transverse, we are ready to obtain the decomposition that this entire section
is building towards.

######  Theorem 3.5 (Pre-decomposition).

We have the decomposition $\mathbb{L}^{2}_{\sigma}=\mathcal{B}+\mathcal{W}$ in
the sense that, for any $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ there exist
unique $\mathcal{X}_{B}\in\mathcal{B}$ and $\mathcal{X}_{W}\in\mathcal{W}$,
called the _balanced component_ of $\mathcal{X}$ and the _wave component_ of
$\mathcal{X}$ respectively, such that

  1. 1.

$\mathcal{X}=\mathcal{X}_{B}+\mathcal{X}_{W}$ and

  2. 2.

$PV(\mathcal{X})=PV(\mathcal{X}_{B})$ and
$\mathcal{M}(\mathcal{X})=\mathcal{M}(\mathcal{X}_{B})$ for
$\mathcal{M}(u,\,\theta,\,q)=\theta+q$ and
$PV(u,\,\theta,\,q)=\nabla_{h}^{\perp}\cdot u_{h}+\partial_{3}\theta$.

Moreover $\mathcal{X}_{B}$ and $\mathcal{X}_{W}$ may be computed explicitly as
follows.

  * •

$\mathcal{X}_{B}=\Phi(p,\,M)$ for $\Phi$ as in Proposition 3.1, where
$M\vcentcolon=\theta+q$ and where $p\in\mathring{H}^{1}$ is the unique
solution of

$\Delta p+\frac{1}{2}\partial_{3}{\min}_{0}\,(M-\partial_{3}p)=PV.$ (3.5)

  * •

$\mathcal{X}_{W}=\mathcal{X}-\mathcal{X}_{B}$, and, in light of Proposition
3.3, we also have $\mathcal{X}_{W}=\Psi(\sigma,\,w)$ where $w=u_{3}$ and

$\displaystyle\sigma$
$\displaystyle=-(u_{h}^{\perp}-u_{B,\,h}^{\perp})+(\theta-\theta_{B})e_{3}$
$\displaystyle=-(u_{h}^{\perp}-\nabla_{h}p)+(\theta-\partial_{3}p-\frac{1}{2}{\min}_{0}\,(M-\partial_{3}p))e_{3}.$

In particular $\sigma$ satisfies

$\nabla\times\left(A_{B}^{-1}\sigma\right)=\partial_{3}u_{h}-\nabla_{h}^{\perp}(\theta-
H_{B}q)+(\partial_{3}w)e_{3},$

where $H_{B}\vcentcolon=\mathds{1}(q_{B}<0)=\mathds{1}(M<\partial_{3}p)$ and
$A_{B}\vcentcolon=I-\frac{1}{2}H_{B}e_{3}\otimes e_{3}$, such that
$A_{B}^{-1}=I+H_{B}e_{3}\otimes e_{3}$, subject to $\nabla\cdot\sigma=0$ and
$\fint A_{B}^{-1}\sigma=-\fint u_{h}^{\perp}+\fint(\theta-H_{B}q)e_{3}$.

###### Proof.

The uniqueness of this decomposition follows from Lemma 3.4 and the fact that,
as per Proposition 3.3, elements of the wave set $\mathcal{W}$ have vanishing
$PV$ and $M$. The existence of the decomposition follows from the existence of
solutions to nonlinear $PV$-and-$M$ inversion since finding $\mathcal{X}_{B}$
is, by virtue of Proposition 3.1, equivalent to finding $p$ and $M$ and since
$p$ may be chosen to be the solution of nonlinear $PV$-and-$M$ inversion.

The latter characterisation of $\sigma$ follows from Lemmas A.5 and A.7:

$\nabla\times\left(A_{B}^{-1}\sigma\right)=\nabla\times\left[-u_{h}^{\perp}+(\theta-
H_{B}q)e_{3}\right]+\nabla\times\nabla
p=\partial_{3}u_{h}-\nabla_{h}^{\perp}(\theta-H_{B}q)+(\partial_{3}w)e_{3}$

and

$\fint A_{B}^{-1}\sigma=-\fint u_{h}^{\perp}+\fint(\theta-H_{B}q)e_{3}.\qed$

It is interesting to see $\nabla\times\left(A_{B}^{-1}\sigma\right)$ arise in
the pre-decomposition above, where matrix $A_{B}^{-1}$ depends on the balanced
moisture, $q_{B}$. It would be difficult to guess the form of
$\nabla\times\left(A_{B}^{-1}\sigma\right)$, with its coefficient matrix that
depends on balanced moisture, $q_{B}$. Here, the appearance of $A_{B}^{-1}$
follows from the decomposition into nullspace and image, and from seeking a
curl-like operator for determining the unbalanced coordinate $\sigma$. Also
note that $\nabla\times\left(A_{B}^{-1}\sigma\right)$ involves a linear
differential operator, once $q_{B}$ is known, but it is a nonlinear
transformation from $(u,\,\theta,\,q)$ to $\sigma$.

As an immediate corollary of this decomposition we make the following
observation: since both the balanced set and the wave set could be
parameterised and since any state can be decomposed into a balanced and a wave
piece, this means that we have obtained a global parametrisation of the entire
state space.

######  Corollary 3.6.

Any state $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ can be uniquely
parameterised by

$\left\\{\begin{aligned} u&=\nabla_{h}^{\perp}p+\sigma_{h}^{\perp}+we_{3},\\\
\theta&=\partial_{3}p+\frac{1}{2}{\min}_{0}\,(M-\partial_{3}p)+\sigma_{3},\text{
and }\\\
q&=M-\partial_{3}p-\frac{1}{2}{\min}_{0}\,(M-\partial_{3}p)-\sigma_{3}\end{aligned}\right.$

for $p\in\mathring{H}^{1}$, $M,\,w\in L^{2}$, and
$\sigma\in{\left(L^{2}\right)}^{3}_{\sigma}$ satisfying
$\nabla_{h}^{\perp}\cdot\sigma_{h}=\partial_{3}w$.

###### Proof.

Theorem 3.5 tells us that any $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ may be
uniquely decomposed as

$\mathcal{X}=\mathcal{X}_{B}+\mathcal{X}_{W}\in\mathcal{B}+\mathcal{W}.$

Then Propositions 3.1 and 3.3 tell us that $\mathcal{X}_{B}=\Phi(p,\,M)$ and
$\mathcal{X}_{W}=\Psi(\sigma,\,w)$ for some uniquely determined $p$, $M$,
$\sigma$, and $w$. So indeed

$\mathcal{X}=\Phi(p,\,M)+\Psi(\sigma,\,w)$

as desired. ∎

We do not even have to stop there: both the balanced set and the wave set had
two _dual_ descriptions. They were described using a parametrisation, which
produced Corollary 3.6 above. But they were also described using the
measurements that annihilate them. This suggests that using the full set of
measurements, namely both those characterising the balanced set and those
characterising the wave set, would yield a complete description of the state
space. We verify this in the result below.

######  Proposition 3.7 (Global nonlinear change of coordinates).

Consider the _measurement space_

$\mathfrak{M}\vcentcolon=H^{-1}\times
L^{2}\times{\left(H^{-1}\right)}^{2}\times L^{2}\times\mathbb{R}^{3}.$

The map $\mathcal{M}:\mathbb{L}^{2}_{\sigma}\to\mathfrak{M}$ given by

$\mathcal{M}(u,\,\theta,\,q)=\begin{pmatrix}\nabla_{h}^{\perp}\cdot
u_{h}+\partial_{3}\theta\\\ \theta+q\\\
\partial_{3}u_{h}-\nabla_{h}^{\perp}(\theta-{\min}_{0}\,q)\\\ u_{3}\\\ -\fint
u_{h}^{\perp}+\fint(\theta-{\min}_{0}\,q)e_{3}\end{pmatrix}=\begin{pmatrix}PV\\\
M\\\ j\\\ w\\\ a\end{pmatrix}$

is invertible.

###### Proof.

To verify the invertibility of this change of coordinates we will use the
coordinates $\pi$, $M$, $\zeta$, and $w$ defined implicitly by

$\left\\{\begin{aligned} u&=\nabla_{h}^{\perp}\pi+\zeta_{h}^{\perp}+we_{3},\\\
\theta&=\partial_{3}\pi+\zeta_{3}+\frac{1}{2}{\min}_{0}\,(M-\partial_{3}\pi-\zeta_{3}),\text{
and }\\\
q&=M-\partial_{3}\pi-\zeta_{3}-\frac{1}{2}{\min}_{0}\,(M-\partial_{3}\pi-\zeta_{3})\end{aligned}\right.$
(3.6)

where $\nabla_{h}^{\perp}\cdot\zeta_{h}=\partial_{3}w$ and
$\nabla\cdot\zeta=0$. Note that the coordinates $(\pi,\,M,\,\zeta,\,w)$ differ
from the coordinates $(p,\,M,\,\sigma,\,w)$ used in Corollary 3.6. The key
advantage of the coordinates used here is that

$\operatorname{sign}q=\operatorname{sign}(M-\partial_{3}\pi-\zeta_{3}),$ (3.7)

which is why they are particularly useful here when it comes to inverting
$\mathcal{M}$.

However these coordinates have the less advantageous property that $\pi$
depends on $PV$, $M$, _and_ $\zeta$, as we will see below. This is different
from $p$ (from Corollary 3.6), which only depends on $PV$ and $M$. That is the
reason why these $\pi$ and $\zeta$ coordinates are _not_ used for the
decomposition (instead of $p$ and $\sigma$, respectively): the resulting
balanced component computed from $\pi$ and $M$ would _not_ be slow.

We close this parenthesis regarding the choice of coordinates and go back to
the matter of invertibility. In light of (3.7), Lemma A.6 and Lemma A.7 tell
us that (3.6) is equivalent to

$\left\\{\begin{aligned}
&-u_{h}^{\perp}+(\theta-{\min}_{0}\,q)e_{3}=\nabla\pi+\zeta,\\\
&u_{3}=w,\text{ and }\\\ &\theta+q=M.\end{aligned}\right.$ (3.8)

In particular it follows from the uniqueness of the Helmholtz decomposition
that these coordinates are uniquely determined by $u$, $\theta$, and $q$.

So now it only remains to show that, for any tuple $(PV,\,M,\,j,\,w,\,a)$,
corresponding coordinates $(\pi,\,M,\,\zeta,\,w)$ exist. Plugging (3.6) into
the definition of the measurement map $\mathcal{M}$ we see, using (3.8), that

$\left\\{\begin{aligned}
PV&=\nabla\pi+\frac{1}{2}\partial_{3}{\min}_{0}\,(M-\partial_{3}\pi-\zeta_{3}),\\\
M&=M,\\\ j&=\partial_{3}\zeta_{h}^{\perp}-\nabla_{h}^{\perp}\zeta_{3},\\\
w&=w,\text{ and }\\\ a&=\fint\zeta.\end{aligned}\right.$

Since $\nabla_{h}^{\perp}\cdot\zeta_{h}=\partial_{3}w$ we see using Lemma A.4
that

$\nabla\times\zeta=j+(\partial_{3}w)e_{3}.$

We are thus guaranteed a (unique) solution $\zeta$ of

$\left\\{\begin{aligned} \nabla\times\zeta&=j+(\partial_{3}w)e_{3},\\\
\nabla\cdot\zeta&=0,\text{ and }\\\ \fint\zeta&=a.\end{aligned}\right.$ (3.9)

The solvability of nonlinear $PV$-and-$M$ inversion [RTSS24] then guarantees a
unique solution $\pi$ in $\mathring{H}^{1}$ of

$\Delta\pi+\frac{1}{2}\partial_{3}{\min}_{0}\,\left[\left(M-\zeta_{3}\right)-\partial_{3}\pi\right]=PV.$
(3.10)

In other words, if we define

$\widetilde{\mathfrak{C}}\vcentcolon=\left\\{(\pi,\,M,\,\zeta,\,w)\in\mathring{H}^{1}\times
L^{2}\times{\left(L^{2}\right)}_{\sigma}^{2}\times
L^{2}:\nabla_{h}^{\perp}\cdot\zeta_{h}=\partial_{3}w\right\\}$

to be the space where the coordinates $(\pi,\,M,\,\zeta,\,w)$ reside and
denote by $(u,\,\theta,\,q)=\mathcal{C}(\pi,\,M,\,\zeta,\,w)$ the map defined
by (3.6), we have shown that for any $(PV,\,M,\,j,\,w,\,a)\in\mathfrak{M}$
there exists a unique $(\pi,\,M,\,\zeta,\,w)\in\widetilde{\mathfrak{C}}$ such
that
$(\mathcal{M}\circ\mathcal{C})(\pi,\,M,\,\zeta,\,w)=(PV,\,M,\,j,\,w,\,a)$.
This means that we can take
$(u,\,\theta,\,q)=\mathcal{C}(\pi,\,M,\,\zeta,\,w)$ to satisfy
$\mathcal{M}(u,\,\theta,\,q)=(PV,\,M,\,j,\,w,\,a)$. Moreover, since we know
that this state $(u,\,\theta,q)$ is uniquely determined by its coordinates
$(\pi,\,M,\,\zeta,\,w)$ this proves that there is a _unique_ state
$(u,\,\theta,\,q)$ satisfying
$\mathcal{M}(u,\,\theta,\,q)=(PV,\,M,\,j,\,w,\,a)$. This proves that
$\mathcal{M}$ is invertible, as desired. Moreover the inverse may be computed
by solving, in order, (3.9) then (3.10) and then plugging the solutions into
(3.6). ∎

######  Remark 3.8 (The global change of coordinates in the dry case).

The proof of Proposition 3.7 above is a bit of a mouthful. In order to shed
some light on the underlying ideas, without having to worry about the details
of the proof itself, we discuss here what an analog global change of
coordinates would look like in the dry case.

In that case we could break down the change of coordinates as a composition of
invertible pieces, thus illustrating more clearly what the moving pieces are.
Indeed, in the dry case we have the following chain of inversions:

$\displaystyle\begin{pmatrix}u\\\
\theta\end{pmatrix}\leftrightarrow\begin{pmatrix}u_{h}\\\ u_{3}\\\
\theta\end{pmatrix}\leftrightarrow\begin{pmatrix}-u_{h}^{\perp}+\theta
e_{3}\\\ u_{3}\end{pmatrix}=\vcentcolon\begin{pmatrix}g\\\
w\end{pmatrix}\overset{(1)}{\leftrightarrow}\begin{pmatrix}\nabla\cdot g\\\
\nabla\times g\\\ \fint g\\\ w\end{pmatrix}$
$\displaystyle\overset{(2)}{=}\begin{pmatrix}\nabla_{h}\cdot
u_{h}+\partial_{3}\theta\\\
\partial_{3}u_{h}-\nabla_{h}^{\perp}\theta+(\partial_{3}w)e_{3}\\\ -\fint
u_{h}^{\perp}+\fint\theta e_{3}\\\ w\end{pmatrix}$
$\displaystyle=\vcentcolon\begin{pmatrix}PV\\\ j+(\partial_{3}w)e_{3}\\\ a\\\
w\end{pmatrix}\leftrightarrow\begin{pmatrix}PV\\\ j\\\ a\\\ w\end{pmatrix}$

where (1) follows from the Helmholtz decomposition (recorded in Theorem A.1)
and (2) follows from Lemma A.5 and the fact that the velocity field $u$ is
divergence-free.

It is therefore tempting to assume that a similar chain of invertible maps
produces the global nonlinear change of coordinates for the _moist_ case
recorded in Proposition 3.7. This is not the case. Of course, we would need to
perform the following inversions initially in the moist case, by virtue of
Lemma A.6:

$\begin{pmatrix}u\\\ \theta\\\
q\end{pmatrix}\leftrightarrow\begin{pmatrix}u_{h}\\\ u_{3}\\\ \theta+q\\\
\theta-{\min}_{0}\,q\end{pmatrix}=\vcentcolon\begin{pmatrix}u_{h}\\\ w\\\ M\\\
b\end{pmatrix}.$

This is not the source of the issue. The source of the issue is that, in the
moist case, the potential vorticity $PV$ and the thermal wind imbalance $j$
come from _different_ good unknowns. Indeed:

$PV=\nabla\cdot\left(-u_{h}^{\perp}+\theta e_{3}\right)\text{ while
}j={\left[\nabla\times\left(-u_{h}^{\perp}+\left(\theta-{\min}_{0}\,q\right)e_{3}\right)\right]}_{h}.$

The good unknowns

$-u_{h}^{\perp}+\theta e_{3}\text{ and
}-u_{h}^{\perp}+\left(\theta-{\min}_{0}\,q\right)e_{3}$

are not the same. Note that this is _not_ a technical issue. This comes from
the dynamics in the moist case being fundamentally different. In particular,
the operator in the moist case loses the skew-symmetry present in the dry
case, and this is precisely the source of the discrepancy between the good
unknowns used for $PV$ and $j$. Indeed, the $(\theta,\,q)$–dependence in the
equation governing $\partial_{t}u$ is of the form

$\partial_{t}u=\dots+(\theta-{\min}_{0}\,q)e_{3}=\dots+(\theta-qH)e_{3}$

for $H\vcentcolon=\mathds{1}(q<0)$ while the $u$–dependence in the equation
governing $\partial_{t}\theta$ and $\partial_{t}q$ is of the form

$\displaystyle\partial_{t}\theta=\dots-u_{3}\text{ and }$
$\displaystyle\partial_{t}q=\dots+u_{3}.$

In matrix form this reads

$\partial_{t}\begin{pmatrix}u\\\ \theta\\\
q\end{pmatrix}=\begin{pmatrix}\dots&e_{3}&-He_{3}\\\ -e_{3}^{T}&\dots&\dots\\\
e_{3}^{T}&\dots&\dots\end{pmatrix}\begin{pmatrix}u\\\ \theta\\\
q\end{pmatrix}.$

This means that the good unknown for $PV$ comes from the first column, while
the good unknown for $j$ comes from the first row. Said another way: $PV$
comes from the _image_ of the operator, whereas $j$ comes from the _null set_
of the operator. Since that matrix is _not_ skew-symmetric, i.e. since the
operator is _not_ skew-adjoint, it follows that $PV$ and $j$ come from
_different_ good unknowns.

######  Remark 3.9.

To conclude this section, we note that we now have _three_ descriptions of the
state space available to us:

  * •

the “primitive” description in terms of the unknowns $u$, $\theta$, and $q$,

  * •

the description in terms of the coordinates $p$, $M$, $\sigma$, and $w$ (see
Corollary 3.6 ), and

  * •

the description in terms of the measurements $PV$, $M$, $j$, $a$, and $w$ (see
Proposition 3.7).

Crucially: the latter two descriptions are compatible with the decomposition!
We can summarize this with the following diagram, where we give a name to the
maps converting between each of these descriptions.

${(u,\,\theta,\,q)\in\mathbb{L}^{2}_{\sigma}}$${(p,\,M,\,\sigma,\,w)\in\mathfrak{C}}$${(PV,\,M,\,j,\,a,\,w)\in\mathfrak{M}}$$\scriptstyle{\mathcal{M}}$$\scriptstyle{\mathcal{S}}$$\scriptstyle{\mathcal{C}}$

Recall that $\mathbb{L}^{2}_{\sigma}$ is the _state space_ introduced in
Definition 2.1 and that $\mathfrak{M}$ is the _measurement space_ defined in
Proposition 3.7. Here $\mathfrak{C}$ is the _coordinate space_ defined as
follows:

$\mathfrak{C}=\left\\{(p,\,M,\,\sigma,\,w)\in\mathring{H}^{1}\times
L^{2}\times{\left(L^{2}\right)}_{\sigma}^{2}\times
L^{2}:\nabla_{h}^{\perp}\cdot\sigma_{h}=\partial_{3}w\right\\}$

The state-to-measurement map $\mathcal{M}$ is then as defined in Proposition
3.7, where it is shown to be invertible, the coordinate-to-state map
$\mathcal{S}$ corresponds to the parametrisation recorded in Corollary 3.6,
which is invertible since Theorem 3.5 tells us how to compute its inverse, and
the measurement-to-coordinate map $\mathcal{C}$ is well-defined by virtue of
the invertibility of the other two maps.

######  Remark 3.10 (Vertically Sheared Horizontal Flows).

This remark discusses vertically sheared horizontal flows, which are
particular solutions to the Boussinesq system characterized by vanishing
thermodynamic variables (namely $\theta$ and $q$ in the moist case) and by a
purely $z$–dependent and horizontal velocity $u=u_{h}(z)$.

The purpose of this remark is two-fold. First it shows what the global change
of coordinates of Corollary 3.6 and Proposition 3.7 looks like in a specific
case, namely when it is applied to these vertically sheared horizontal flow.
Second it shows one benefit of the decomposition presented here: these flows
are accounted for _naturally_.

First, we recall that combining Corollary 3.6 and Proposition 3.7 gives us
_three_ descriptions of any given state.

  * •

The description in physical space in terms of $(u,\,\theta,\,q)$.

  * •

The description in coordinate space in terms of $(p,\,M,\,\sigma,\,w)$.

  * •

The description in measurement space in terms of $(PV,M,\,j,\,w,\,a)$.

Let us now see what vertically sheared horizontal flows look like in each of
these three descriptions.

  * •

In physical space a vertically sheared horizontal flow corresponds to any
state $(u,\,\theta,\,q)$ where $\theta$ and $q$ vanish and where $u=u_{h}(z)$,
i.e. $u_{3}=0$ and $u$ only depends on $z$.

  * •

In coordinate space a vertically sheared horizontal flow corresponds to any
state described by $(p,\,M,\,\sigma,\,w)$ where $p$, $M$, and $w$ vanish and
where $\sigma=\sigma_{h}(z)$.

  * •

In measurement space a vertically sheared horizontal flow corresponds to any
any state described by $(PV,\,M,\,j,\,w,\,a)$ where $PV$, $M$, and $w$ vanish,
where $j=j(z)$, and where $a=(a_{1},\,a_{2},\,0)$ is any purely horizontal
vector in $\mathbb{R}^{3}$.

In other words, since a vertically sheared horizontal flows is fully
characterized by a map $\gamma:\mathbb{T}^{1}\to\mathbb{R}^{2}$, here are
three descriptions of the same such flow.

Physical space Coordinate space Measurement space $\displaystyle
u=\gamma_{h}(z)$ $\displaystyle p=0$ $\displaystyle PV=0$
$\displaystyle\theta=0$ $\displaystyle M=0$ $\displaystyle M=0$ $\displaystyle
q=0$ $\displaystyle\sigma=-\gamma^{\perp}_{h}(z)$ $\displaystyle
j=\gamma^{\prime}(z)$ $\displaystyle w=0$ $\displaystyle w=0$ $\displaystyle
a=-\fint\gamma_{h}^{\perp}(z)dz$

Second, we note that as described above our decomposition accounts for these
vertically sheared horizontal flows very naturally. This comes from the fact
that our decomposition leverages a _three-dimensional_ Helmholtz decomposition
of the _mixed vector_ fields $-u^{\perp}_{h}+\theta e_{3}$ (for $PV$) and
$-u_{h}^{\perp}+(\theta-{\min}_{0}\,q)e_{3}$ (for $j$). This is different from
what is often done where a _two-dimensional_ Helmholtz decomposition is
performed _solely on_ the horizontal velocity $u_{h}$. This two-dimensional
Helmholtz decomposition then leaves a two-dimensional average term of the form
$u_{h}(z)$ to be determined, which is precisely the term accounting for these
vertically sheared horizontal flows. In other words: our decomposition
accounts for these intrinsically, instead of having to treat them as a
separate add-on after the fact.

######  Remark 3.11.

It is important to note that the decomposition discussed in this section is
_always_ valid, even when the state in $\mathbb{L}^{2}_{\sigma}$ thus
decomposed is _not_ a solution of the moist Boussinesq dynamical system. In
other words: the decomposition of Theorem 3.5 is valid for any state in the
state space $\mathbb{L}^{2}_{\sigma}$, similarly to how the Helmholtz
decomposition is valid for any three-dimensional vector field.

However, as discussed in the next section, the decomposition we introduced in
this context does also have additional meaning in the context of the
_dynamics_ of the moist Boussinesq system (2.7a)–(2.7c).

##  4 Decomposition and slow-fast dynamics

In Section 3 above we have introduced a preliminary decomposition of the
leading-order dynamics. In particular we have discussed how to then make sure
that the balanced component is slow. The key step in doing so was ensuring
that the pressure, which acts as a streamfunction for the balanced component
of the velocity, be slow. This differs from previous work [ZSS21b, ZSS21a,
WSSM19, WSS+20] where the streamfunction constructed was approximately slow
but not precisely slow.

A question remains: “did we extract _all_ of the slow components of a state in
this manner?”. In the dry case, that question may be answered in the
affirmative since the fast dynamics are that of linear oscillations which do
not feedback into the slow component. In the moist case, however, the fast
oscillations are nonlinear and so they _may_ feed back into the slow
component. To make this idea precise we appeal to the framework of fast-wave
averaging. Fast-wave averaging is a technique used to identify the limiting
behaviour of certain dynamical systems which admit two timescales. The
prototypical case in which two timescales would be present is the following
dynamical system in $\mathbb{R}^{n_{1}+n_{2}}$:

$\partial_{t}\begin{pmatrix}x_{1}\\\
x_{2}\end{pmatrix}=\begin{pmatrix}f_{1}(x_{1})\\\
\frac{1}{\varepsilon}f_{2}(x_{2})\end{pmatrix}$

for $(x_{1},\,x_{2})\in\mathbb{R}^{n_{1}+n_{2}}$, where $\varepsilon>0$ is a
small parameter. In this case we would say that the dynamics of $x_{1}$ are
slow since they are of order one, whereas the dynamics of $x_{2}$ are fast
since they are of order $1/\varepsilon$. In general, dynamical systems will
not admit an _exact_ two timescale decomposition as in the prototypical
example above. In general, that decomposition is only approximately valid for
finite $\varepsilon$ and we are then interested in the limiting behaviour of
the slow and fast components as $\varepsilon\to 0$.

Where does fast-wave averaging fit in all this? Well: possessing two
timescales is a feature a dynamical system may or may not manifest, whereas
fast-wave averaging is a technique used to identify the limiting behaviour of
the slow and fast timescales (for such dynamical systems). Fast-wave averaging
happens in four steps.

  1. 1.

Assume a two-timescale asymptotic expansion of the unknowns (which requires a
sub-linear growth condition).

  2. 2.

Plug this expansion into the equations of motion to identify the resulting
equations at leading-order $1/\varepsilon$ and at next-order.

  3. 3.

Solve the leading-order equations over the fast timescale in terms of a yet-
to-be-determined “background profile” at the slow timescale.

  4. 4.

Determine this background profile at the slow timescale by using a solvability
condition for the next-order equations.

We will now run through these four steps for the moist Boussinesq system
(2.7a)–(2.7c). At each step we will highlight the differences with the dry
case. This is because the dry case is well-understood however in the moist
case we only provide here a sketch of the fast-wave averaging process, leaving
a rigorous treatment for future work.

First we assume a two-timescale expansion of the unknowns. This is done
exactly as in the dry case and we assume that, for
$\mathcal{X}=(u,\,\theta,\,q)$,

$\mathcal{X}=\mathcal{X}^{0}(t,\,\tau)|_{\tau=t/\varepsilon}+\varepsilon\mathcal{X}^{1}(t,\,\tau)|_{\tau=t/\varepsilon}+O(\varepsilon^{2})$
(4.1)

for some $\mathcal{X}^{0}$ and $\mathcal{X}^{1}$ to be determined where $t$ is
the slow timescale and $\tau=t/\varepsilon$ is the fast timescale, subject to
the _sub-linear growth condition_

$\mathcal{X}^{1}(t,\,\tau)=o(\tau)\text{ in $L^{2}$ uniformly in $t$ as
}\tau\to\infty,$ (4.2)

or equivalently for every $\delta>0$ there exists $\tau_{*}>0$ such that if
$\tau>\tau_{*}$ then

$\frac{{\left\lvert\left\lvert\mathcal{X}^{1}(t,\,\tau)\right\rvert\right\rvert}_{L^{2}}}{\tau}<\delta\text{
for all }t.$

Note that this sub-linearity condition is imposed for the following reason. If
the next-order term $\mathcal{X}^{1}$ is allowed to grow linearly in $\tau$
then we would have that

$\varepsilon\mathcal{X}^{1}(t,\,\tau)|_{\tau=t/\varepsilon}\sim\varepsilon\left(\mathcal{X}^{1}(t,\,0)+\frac{t}{\varepsilon}\right)=\varepsilon\mathcal{X}^{1}(t,\,0)+t,$

i.e. the term $t$ of order one would appear in an expression that ought to be
of order $\varepsilon$.

The next step is to insert the two-timescale assumption (4.1) into the moist
Boussinesq system (2.7a)–(2.7c) and read off the leading-order and next-order
equations that come out. In order to carry this out we record here two
preliminary computations. First we note that

$\partial_{t}\left[f(t,\,\tau)\middle|_{\tau=t/\varepsilon}\right]=\left[(\partial_{t}f)(t,\,\tau)+\frac{1}{\varepsilon}(\partial_{\tau}f)(t,\,\tau)\middle]\right|_{\tau=t/\varepsilon}.$

Second we note that, introducing the notation
$H\vcentcolon=\mathds{1}(q^{0}+\varepsilon q^{1}<0)$ and
$H^{0}\vcentcolon=\mathds{1}(q^{0}<0)$,

${\min}_{0}\,(q_{0}+\varepsilon q^{1})=(q^{0}+\varepsilon
q^{1})H=\underbrace{q^{0}H^{0}}_{{\min}_{0}\,q^{0}}+q^{0}(H-H^{0})+\varepsilon
q^{1}H$

where

$\displaystyle q^{0}(H-H^{0})$
$\displaystyle=q^{0}\left[\mathds{1}(q^{0}+\varepsilon
q^{1}<0)-\mathds{1}(q^{0}<0)\right]$
$\displaystyle=q^{0}\left[\mathds{1}(q^{0}+\varepsilon
q^{1}<0,\,q^{0}\geqslant 0)-\mathds{1}(q^{0}<0,\,q^{0}+\varepsilon
q^{1}\geqslant 0\right]$ $\displaystyle=q^{0}\left[\mathds{1}(0\leqslant
q^{0}<-\varepsilon q^{1})-\mathds{1}(-\varepsilon q^{1}\leqslant
q^{0}<0)\right]$
$\displaystyle=\varepsilon\cdot\frac{q^{0}}{\varepsilon}\left[\mathds{1}\left(0\leqslant\frac{q^{0}}{\varepsilon}<-q^{1}\right)-\mathds{1}\left(-q^{1}\leqslant\frac{q^{0}}{\varepsilon}<0\right)\right]$
$\displaystyle=\vcentcolon\varepsilon R(q^{0};\,q^{1}).$

In particular we note that $R(q^{0};\,q^{1})\in O(1)$ provided that $q^{1}$ is
bounded.

Using these two observations we plug (4.1) into (2.7a)–(2.7c) and read off the
leading-order equations, at order $1/\varepsilon$, to be

$\displaystyle\partial_{\tau}u^{0}+e_{3}\times u^{0}+\nabla
p^{0}-\theta^{0}e_{3}+{\min}_{0}\,q^{0}e_{3}=0,$ (4.3a)
$\displaystyle\partial_{\tau}\theta^{0}+u_{3}^{0}=0,\text{ and }$ (4.3b)
$\displaystyle\partial_{\tau}q^{0}-u_{3}^{0}=0$ (4.3c)

while the next-order equations, at order one (i.e. order $\varepsilon^{0}$) to
be

$\displaystyle\partial_{\tau}u^{1}+e_{3}\times u^{1}+\nabla
p^{1}-\theta^{1}e_{3}+q^{1}He_{3}=-\partial_{t}u^{0}-u^{0}\cdot\nabla
u^{0}+R(q^{0};\,q^{1}),$ (4.4a)
$\displaystyle\partial_{\tau}\theta^{1}+u_{3}^{1}=-\partial_{t}\theta^{0}-u^{0}\cdot\nabla\theta^{0},\text{
and }$ (4.4b)
$\displaystyle\partial_{\tau}q^{1}-u_{3}^{1}=-\partial_{t}q^{0}-u^{0}\cdot\nabla
q^{0}.$ (4.4c)

We can also write these equations more compactly as

$\displaystyle\partial_{\tau}\mathcal{X}^{0}+\mathcal{N}(\mathcal{X}^{0})\text{
and }$ (4.5a)
$\displaystyle\partial_{\tau}\mathcal{X}^{1}+\mathcal{L}_{H}(\mathcal{X}^{1})=-\partial_{t}\mathcal{X}^{0}-u^{0}\cdot\nabla\mathcal{X}^{0}+\mathcal{R}(\mathcal{X}^{0};\,\mathcal{X}^{1}),$
(4.5b)

for $\mathcal{N}$ as defined in (3.1), $\mathcal{L}_{H}$ its formal
linearisation given by

$\mathcal{L}_{H}\begin{pmatrix}u\\\ \theta\\\
q\end{pmatrix}=\begin{pmatrix}\mathbb{P}_{L}\left[e_{3}\times u-(\theta-
qH)e_{3}\right]\\\ u_{3}\\\ -u_{3}\end{pmatrix},$ (4.6)

and the remainder term $\mathcal{R}$ given by
$\mathcal{R}(u^{0},\,\theta^{0},\,q^{0};\,u^{1},\,\theta^{1},\,q^{1})=(R(q^{0};\,q^{1}),\,0,\,0)$.

Crucially: we see that the dynamics at leading-order are governed by
$\partial_{\tau}+\mathcal{N}=0$ for $\mathcal{N}$ as introduced in (3.1). The
sets $\mathcal{B}$ and $\mathcal{W}$ characterised in Propositions 3.1 and 3.3
thus play an essential role since they come from that same operator
$\mathcal{N}$. Indeed, as discussed in Proposition A.9, a solution
$\mathcal{X}$ of $\partial_{\tau}\mathcal{X}+\mathcal{N}(\mathcal{X})=0$ may
be written as

$\mathcal{X}(\tau)=\mathcal{X}_{B}(0)+e^{-\tau\mathcal{L}_{W,\,H_{B}}}\mathcal{X}_{W}(0).$
(4.7)

Note that here we abused notation a little bit in order to be reminiscent of
what happens in the dry case. In the dry case the leading order dynamics are
$\partial_{\tau}+\mathcal{L}$ for $\mathcal{L}$ given by

$\mathcal{L}(u,\,\theta)=\begin{pmatrix}\mathbb{P}_{L}\left(u_{h}^{\perp}-\theta
e_{3}\right)\\\ u_{3}\end{pmatrix},$ (4.8)

as per (2.1a)–(2.1b). Since the operator $\mathcal{L}$ is linear and constant,
it follows that the solution map of $\partial_{\tau}+\mathcal{L}$ is precisely
given by the semigroup $e^{-t\mathcal{L}}$. The operator
$\mathcal{L}_{W,\,H_{B}}$ in the moist case is linear but with _variable_
coefficients (namely through its dependence on $H_{B}$), yet we use the same
notation to represent the solution map, keeping it mind that it is slightly
abusive notation since the solution map is _not_ a semigroup.

We are thus ready to carry out the third step of the fast-wave averaging
process. We deduce from (4.7) that the solution of the leading-order equations
(4.5a) is given by

$\underline{\mathcal{X}}^{0}(t,\,\tau)=\underline{\mathcal{X}}^{0}_{B}(t)+e^{-\tau\mathcal{L}_{W,\,H_{B}}}\underline{\mathcal{X}}_{W}^{0}(t)$

where the background slow profile
$\underline{\mathcal{X}}^{0}(t)=\underline{\mathcal{X}}_{B}^{0}(t)+\underline{\mathcal{X}}_{W}^{0}(t)$
is to be determined (this is done in the next step of the fast-wave averaging
process). This is where the first essential difference with the dry case
manifests itself. In the dry case, the leading-order solution takes the form

$\underline{\mathcal{X}}^{0}(t,\,\tau)=\underline{\mathcal{X}}^{0}_{B}(t)+e^{-\tau\mathcal{L}_{W}}\underline{\mathcal{X}}^{0}_{W}(t),$
(4.9)

where $\mathcal{L}_{W}$ denotes the restriction of $\mathcal{L}$ from (4.8)
onto its image. Crucially: since the operator $\mathcal{L}_{W}$ is skew-
adjoint in $L^{2}$ the oscillations in (4.9) are linear oscillations with
constant frequencies. Here, in the moist case, _the frequencies depend on the
phase_.

To see why this matters, suppose for now that we have completed the fourth and
last step of the fast-wave averaging process. That is, suppose we have
determined the dynamics of the slow part $\underline{\mathcal{X}}^{0}(t)$.
Since we are interested in the limiting behaviour of solutions, we then wish
to study the limiting behaviour of $\mathcal{X}^{0}(t,\,\tau)$ as
$\varepsilon\to 0$. In the _dry_ case this means studying the limiting
behaviour of

$e^{-\tau\mathcal{L}_{W}}\underline{\mathcal{X}}_{W}^{0}(t).$

Crucially: since the dry oscillations have constant frequencies, these
oscillations have vanishing time averages, and so the weak limit of this term
is simply zero!

In the moist case we must by contrast study the limiting behaviour of

$e^{-\tau\mathcal{L}_{W,\,H_{B}}}\underline{\mathcal{X}}_{W}^{0}(t).$

But now, as shown in an example first recorded in [ZSS21a] (see Corollary
A.22, and the simple example it arises from in Proposition A.19), phase-
dependent oscillations do _not_ have vanishing time averages. This wave term
therefore converges weakly to a _non-zero_ time average!

This tells us exactly how to augment our “pre-decomposition” above, in Theorem
3.5, to now make sure we include _all_ of the slow components: we take the
slow balanced piece from the previous section and add to it the time-average
of the wave piece.

This produces the following slow-fast decomposition, where for any space $S$,
$L^{1}S$ is used to denote $L^{1}([0,\,\infty);\,S)$.

######  Theorem 4.1 (Slow-fast decomposition).

We have the decomposition

$L^{1}\mathbb{L}^{2}_{\sigma}=L^{1}\mathcal{B}+\overline{\mathcal{W}}+\widetilde{\mathcal{W}}$

where

$\overline{\mathcal{W}}\vcentcolon=\left\\{\mathcal{X}\in
L^{1}\mathcal{W}:\mathcal{X}\text{ is independent of time
$\tau$}\right\\}\text{ and
}\widetilde{\mathcal{W}}\vcentcolon=\left\\{\mathcal{X}\in
L^{1}\mathcal{W}:\lim_{\tau\to\infty}\fint_{0}^{\tau}\mathcal{X}(\tau^{\prime})d\tau^{\prime}=0\right\\},$

in the sense that, for every $\mathcal{X}\in L^{1}\mathbb{L}^{2}_{\sigma}$
there exist unique $\mathcal{X}_{B}\in L^{1}\mathcal{B}$,
$\overline{\mathcal{X}}_{W}\in\overline{\mathcal{W}}$, and
$\widetilde{\mathcal{X}}_{W}\in\widetilde{\mathcal{W}}$ such that

  1. 1.

$\mathcal{X}=\mathcal{X}_{B}+\overline{\mathcal{X}}_{W}+\widetilde{\mathcal{X}}_{W}$
and

  2. 2.

$PV(\mathcal{X})=PV(\mathcal{X}_{B})$ and
$\mathcal{M}(\mathcal{X})=\mathcal{M}(\mathcal{X}_{B})$ for $PV$ and
$\mathcal{M}$ as in Theorem 3.5.

Moreover $\mathcal{X}_{B}(\tau)$ and $\mathcal{X}_{W}(\tau)$ may be computed
explicitly at every instant $\tau$ as in Theorem 3.5 such that then
$\overline{\mathcal{X}}_{W}\vcentcolon=\lim_{\tau\to\infty}\fint_{0}^{\tau}\mathcal{X}_{W}(\tau^{\prime})d\tau^{\prime}$
and
$\widetilde{\mathcal{X}}\vcentcolon=\mathcal{X}-\overline{\mathcal{X}}_{W}$.

###### Proof.

This follows immediately from applying Theorem 3.5 to $\mathcal{X}(\tau)$ at
every instant $\tau\geqslant 0$ and then constructing
$\overline{\mathcal{X}}_{W}$ and $\widetilde{\mathcal{X}}_{W}$ as indicated. ∎

Note that the technical aspects of this theorem were already taken care of in
the proof of the “pre-decomposition” of Theorem 3.5. So there is nothing new
to prove here. Nonetheless, as explained in this section in the lead-up to
Theorem 4.1, this statement now has an _interpretation_ that Theorem 3.5 did
not. Once viewed through the lens of fast-wave averaging, we can indeed
interpret Theorem 4.1 as a _slow-fast_ decomposition:

$\mathcal{X}=\underbrace{\mathcal{X}_{B}+\overline{\mathcal{X}}_{W}}_{\mathcal{X}_{\text{slow}}}+\underbrace{\widetilde{\mathcal{X}}_{W}}_{\mathcal{X}_{\text{fast}}}.$

In other words: the fast-wave averaging framework now, at least formally,
ensures that this decomposition accounts for the entirety of the slow
components of the solution!

Finally, for completeness we record the fourth and final step of fast-wave
averaging. The idea is that we can solve (4.5b) for $\mathcal{X}^{1}$, but we
are not particularly interested in what that solution looks like. What we care
about is making sure that this next-order solution satisfies the sub-linear
growth condition (4.2). Ultimately, this comes down to a condition on the
right-hand side forcing terms of (4.5b). Since these forcing terms depend on
the leading-order solution $\mathcal{X}^{0}$, enforcing the sub-linear growth
condition on $\mathcal{X}^{1}$ turns into a condition on $\mathcal{X}^{0}$.

In the dry case, this condition takes the form

$\fint_{0}^{\tau}\partial_{t}\underline{\mathcal{X}}^{0}+e^{\tau^{\prime}\mathcal{L}}\left(e^{-\tau^{\prime}\mathcal{L}}\bar{u}^{0}\cdot\nabla
e^{-\tau^{\prime}\mathcal{L}}\underline{\mathcal{X}}_{0}\right)d\tau^{\prime}=0$

since we may use the fact that the semigroup $e^{-\tau\mathcal{L}}$ is norm-
preserving. In the moist case however the solution map
$e^{-\tau\mathcal{L}_{W,\,H_{B}}}$ is _not_ norm-preserving (see Proposition
A.12 ) and so we cannot perform the same trick. Identifying what the
solvability condition looks like in the moist case is therefore left open for
future work.

######  Remark 4.2 (Slow-fast interpretation in different limiting regimes).

In this section we have then described how, once in the limiting regime where
both the Froude and Rossby numbers approach zero, the decomposition introduced
in this paper has an additional interpretation as a slow-fast decomposition.
Nonetheless, this same decomposition may still be of use in other limiting
regimes, provided the interpretation of each of its components as slow or fast
is updated appropriately.

Consider for example the case where the Froude number is still taken to
approach zero but the Rossby number is kept constant. In that case the
vertically sheared horizontal flows discussed in Remark 3.10 are no longer
fast, but now become slow!

##  5 Decomposition and energy

In the dry case, the decomposition of states into their slow and fast
components has the additional properties that it is an orthogonal
decomposition, where the inner product used is the one induced by the
conserved energy.

In the moist case, which is the focus of this paper, we therefore ask the
following question: how does our decomposition of Theorem 3.5 relate to the
energy conserved by the moist Boussinesq system? The energy in question is
given by

$E=\frac{1}{2}\int{\lvert u\rvert}^{2}+\theta^{2}+{\min}_{0}^{2}\,q+M^{2},$
(5.1)

where an arbitrary coefficient of the $M^{2}$ term has been set to 1 for
simplicity. For a demonstration that this energy is conserved, see Proposition
A.10.

In the dry case, since the decomposition comes down to linear projections, and
since the energy gives rise to an inner product, the relation between the
decomposition and the energy can be phrased as a simple question: are these
projections orthogonal with respect to this inner product? In the moist case,
both the decomposition and the energy are _nonlinear_ and so it is not clear
how one would go about even _asking_ if the two are related. We offer several
options in Sections 5.1-5.3, summarizing and contrasting each option, along
with their benefits and inconveniences, in Section 5.4.

Note that, in each of the three cases discussed below, the strategy is the
same. Namely we seek to identify a way to quantify how far apart two states
$(u_{1},\,\theta_{1},\,q_{1})$ and $(u_{2},\,\theta_{2},\,q_{2})$ are in such
a manner that

  1. 1.

this agrees with the decomposition of Theorem 3.5, and gives it further
meaning as a projection, in the sense that the closest balanced state to any
given state must be its balanced component and

  2. 2.

this agrees with the energy (when the two are comparable).

Finally, note that the authors believe that the quantification provided in
Section 5.3 is really the right one to use for describing the decomposition
from a projection perspective. However, it may appear a bit odd at first. To
better motivate it we therefore discuss, in Sections 5.1 and 5.2, alternative
approaches which are less conclusive but help us understand why the _good_
approach, discussed in Section 5.3, must take the form it does. Also, the
alternative approaches in Sections 5.1 and 5.2 have value in their own right,
for other applications beyond the focus here on a nonlinear eigenspace
decomposition.

Before we dive into the details of each of these three approaches, we must
elucidate the following question: what does it mean for two states to be
comparable using the energy?

######  Remark 5.1 (Comparing states using the energy).

Recall that the energy of a state $\mathcal{X}=(u,\,\theta,\,q)$ is

$E(\mathcal{X})=\int{\lvert
u\rvert}^{2}+\theta^{2}+{\min}_{0}^{2}\,q+{(\theta+q)}^{2}.$

The question is: given two states $\mathcal{X}_{1}$ and $\mathcal{X}_{2}$,
_when_ can we use the energy to quantify the discrepancy between
$\mathcal{X}_{1}$ and $\mathcal{X}_{2}$?

A naive attempt would be to use the energy of their _difference_ , i.e.

$E(\mathcal{X}_{1}-\mathcal{X}_{2})=\int{\lvert
u_{1}-u_{2}\rvert}^{2}+{(\theta_{1}-\theta_{2})}^{2}+{\min}_{0}^{2}\,(q_{1}-q_{2})+{((\theta_{1}+q_{1})-(\theta_{2}+q_{2}))}^{2}.$
(5.2)

This is not physically meaningful. Why? Because of the nonlinearity. When
computing $E(\mathcal{X}_{1}-\mathcal{X}_{2})$ the nonlinear term
${\min}_{0}^{2}\,(q_{1}-q_{2})$ only turns on when $q_{1}-q_{2}>0$. However
this nonlinearity is not supposed to care about the sign of $q_{1}-q_{2}$.
Instead the nonlinearity is there to account for the fact that, depending on
the signs of $q_{1}$ and $q_{2}$, the energetic weights of $q_{1}$ and $q_{2}$
change.

What is the takeaway? The only time we can directly compare two states
$\mathcal{X}_{1}$ and $\mathcal{X}_{2}$ using the energy is when the signs of
$q_{1}$ and $q_{2}$ agree everywhere in the spatial domain. The sign of $q$ is
a mathematical proxy used to tell us which _phase_ the state is in at any
point in the spatial domain. We may thus rephrase the observation above as
follows.

The only time we can directly compare two states $\mathcal{X}_{1}$ and
$\mathcal{X}_{2}$

is when their phases coincide everywhere in the spatial domain.

It will therefore be convenient to encode the spatial location of a state’s
phases. This is done by using the indicator function

$H\vcentcolon=\mathds{1}(q<0).$

Note that since the water content $q=q(t,\,x)$ is a function of $t$ and $x$,
then the indicator $H$ also depends on $t$ and $x$. As above, we will often
use notation which suppresses the $(t,x)$-dependence in an effort to improve
readability.

So let us consider two states $\mathcal{X}_{1}$ and $\mathcal{X}_{2}$ whose
phases agree, i.e. for which $H_{1}=H_{2}=\vcentcolon H$. On one hand

$\displaystyle E(\mathcal{X}_{1})$ $\displaystyle=\int{\lvert
u_{1}\rvert}^{2}+\theta_{1}^{2}+{\min}_{0}^{2}\,q_{1}+{(\theta_{1}+q_{1})}^{2}$
$\displaystyle=\int{\lvert
u_{1}\rvert}^{2}+\theta_{1}^{2}+q_{1}^{2}H+{(\theta_{1}+q_{1})}^{2}$

while on the other hand

$\displaystyle E(\mathcal{X}_{2})$ $\displaystyle=\int{\lvert
u_{2}\rvert}^{2}+\theta_{2}^{2}+{\min}_{0}^{2}\,q_{2}+{(\theta_{2}+q_{2})}^{2}$
$\displaystyle=\int{\lvert
u_{2}\rvert}^{2}+\theta_{2}^{2}+q_{2}^{2}H+{(\theta_{2}+q_{2})}^{2}.$

A natural way to compare $\mathcal{X}_{1}$ and $\mathcal{X}_{2}$ using the
energy is then to treat the indicator function $H$ as _fixed_ in the
expressions above for both $E(\mathcal{X}_{1})$ and $E(\mathcal{X}_{2})$ and
to consider

$\int{\lvert
u_{1}-u_{2}\rvert}^{2}+{(\theta_{1}-\theta_{2})}^{2}+{(q_{1}-q_{2})}^{2}H+{((\theta_{1}+q_{1})-(\theta_{2}+q_{2}))}^{2}.$
(5.3)

This is _not_ the same as $E(\mathcal{X}_{1}-\mathcal{X}_{2})$, as can be seen
by comparing (5.2) with (5.3). Nonetheless the quantity appearing in (5.3) is,
so far, the only way we have of using the energy to compare two states, as
long as their phases agree everywhere.

So we give this quantity a name.

######  Definition 5.2.

For _any_ indicator function $H$ and _any_ state
$\mathcal{X}=(u,\,\theta,\,q)\in\mathbb{L}^{2}_{\sigma}$ let

$E_{H}(\mathcal{X})\vcentcolon=\int{\lvert
u\rvert}^{2}+\theta^{2}+q^{2}H+{(\theta+q)}^{2}.$

Note that this notation may then be used even when (especially when!) the sign
of $q$ does not agree with $H$. That is precisely what happens in (5.3): $H$
encodes the signs of $q_{1}$ and $q_{2}$, _not_ the sign of $q_{1}-q_{2}$
(that was precisely the naive nonphysical approach we sought to avoid!)

Here we have motivated this “fixed $H$” energy by talking about restricting
our attention to a subset of states whose phases agree. There is another
reason to look at this energy: it is the quadratic approximation of the
conserved energy about a state whose phase are encoded in $H$. So this energy
also appears in situations where linearisation are used, since linearising the
dynamical equations is essentially the same as working with the quadratic
approximation of the conserved energy. This interpretation of the “fixed $H$”
energy will be more prevalent in later sections, such as Section 6.2 where
geometrical notions relying on tangent spaces (a.k.a. linearisations) are
discussed. Linearisations and the “fixed $H$” energy are also relevant for
applications such as adjoint models, sensitivity analysis, and variational
data assimilation, which rely on tangent linear approximations to the dynamics
[LDT86, PD97, Err97, ER99, Mah99, BBP+01, AZD08, DARR14].

###  5.1 PDE-centric approach

We first discuss a PDE-centric approach to relating the decomposition and the
energy. Why start here? There are two reasons to begin with a PDE-centric
approach.

The first reason to take such an approach is that elliptic PDEs underpin the
decomposition in both the case of the Helmholtz decomposition and in the dry
case. Since a nonlinear elliptic PDE is one of the few structures that have
already been identified to be present in the moist case, this PDE therefore
presents a natural starting point.

The second reason to take a PDE-centric approach is that a similar approach
was fruitful in identifying the slow–fast decomposition. Indeed, the slow–fast
decomposition was obtained by recognizing the important role played by the
slow measurements $PV$ and $M$, and then inverting the nonlinear $PV$-and-$M$
inversion PDE (3.5) in order to obtain $p$. The variables $p$ and $M$ were
then used to characterise the balanced part of a state.

Equipped with Propositions 3.1 and 3.3 we can actually go a step further: the
slow measurements $PV$ and $M$ characterise the balanced component of a state
while the fast measurements $j$, $w$, and $a$ characterise its wave component.
Indeed, we have the following.

  * •

On one hand Proposition 3.1 tells us that the measurements

$\displaystyle
j(u,\,\theta,\,q)=\partial_{3}u_{h}-\nabla_{h}^{\perp}(\theta-{\min}_{0}\,q),\,w(u,\,\theta,\,q)=u_{3},$
$\displaystyle\text{ and }a(u,\,\theta,\,q)=-\fint
u_{h}^{\perp}+\fint(\theta-{\min}_{0}\,\theta)e_{3}$

are transversal to the balanced set $\mathcal{B}$, and so we may expect them
to characterize elements of the wave set $\mathcal{W}$.

  * •

On the other hand Proposition 3.3 tells us that the measurements

$PV(u,\,\theta,\,q)=\nabla_{h}^{\perp}\cdot u_{h}+\partial_{3}\theta\text{ and
}\\\ M(u,\,\theta,\,q)=\theta+q$

are transversal to the wave set $\mathcal{W}$, and so we may similarly expect
them to characterize elements of the balanced set $\mathcal{B}$.

Both of these assertions turn out to be true, and thus give rise to the global
change of coordinates recorded in Proposition 3.7 which allows us to translate
between a state $(u,\,\theta,\,q)$ and the tuple of measurements
$(PV,\,M,\,j,\,w,\,a)$.

We will now use this change of variable from Proposition 3.7 to define a
distance on the state space $\mathbb{L}^{2}_{\sigma}$. We will do this by
using the change of variables to pull-back a metric on the measurement space
$\mathcal{M}$ and thus produce a metric on the state space
$\mathbb{L}^{2}_{\sigma}$ (see Proposition A.13). The execution of this
strategy requires a few steps. In order to have a clearer idea of where we are
headed, we record below what would happen if we tried doing the same thing in
the _dry_ case. In the dry case we can write the energy in terms of the
measurements (namely $PV$, $j$, $w$, and $a$). In the language that will be
used in this section, we can therefore view the energy, and the $L^{2}$ metric
it corresponds to, as the pullback of a metric on the (dry) measurement space!
Since we will ultimately, in this section, define the distance of interest in
terms of the _coordinates_ (think $p$, $M$, $\sigma$, and $w$), we also record
in Proposition 5.3 below how the energy may be written in a third way, in
terms of the coordinates (which are $p$, $\sigma$, and $w$ in the dry case).

######  Proposition 5.3 (Parseval identity in the dry case).

For any $u\in{\left(L^{2}\right)}_{\sigma}^{3}$ and $\theta\in L^{2}$ we have
that

$\underbrace{\int{\lvert
u\rvert}^{2}+\theta^{2}}_{E}=\int{\lvert\nabla\Delta^{-1}PV\rvert}^{2}+{\left\lvert{(\nabla\times)}^{-1}\left(j+(\partial_{3}w)e_{3},\,a\right)\right\rvert}^{2}=\underbrace{\int{\lvert\nabla
p\rvert}^{2}+{\lvert\sigma\rvert}^{2}+w^{2}}_{\sim d^{2}}$

where $p$ solves

$\Delta p=\underbrace{\nabla_{h}^{\perp}\cdot u_{h}+\partial_{3}\theta}_{PV}$

while $\sigma$ solves

$\nabla\times\sigma=\underbrace{\partial_{3}u_{h}-\nabla_{h}^{\perp}\theta}_{j}+(\partial_{3}w)e_{3}\text{
subject to }\fint\sigma=\underbrace{\fint-u_{h}^{\perp}+\fint\theta
e_{3}}_{a},$

and where $w\vcentcolon=u_{3}$. Here we write
$\sigma={(\nabla\times)}^{-1}(\omega,\,a)$ to mean that $\sigma$ solves

$\nabla\times\sigma=\omega\text{ subject to }\nabla\cdot\sigma=0\text{ and
}\fint\sigma=a.$

###### Proof.

We recall that the decomposition in the dry case is

$\begin{pmatrix}u\\\ \theta\end{pmatrix}=\begin{pmatrix}\nabla_{h}^{\perp}p\\\
\partial_{3}p\end{pmatrix}+\begin{pmatrix}\sigma_{h}^{\perp}+we_{3}\\\
\sigma_{3}\end{pmatrix}$

where $\sigma$ is divergence-free and satisfies
$\nabla_{h}^{\perp}\cdot\sigma_{h}=\partial_{3}w$. We may then compute that

$PV=\nabla_{h}^{\perp}\cdot u_{h}+\partial_{3}\theta=\Delta p$

while Lemma A.4 tells us that

$\nabla\times\sigma=\partial_{3}\sigma_{h}^{\perp}-\nabla_{h}^{\perp}\sigma_{3}+(\nabla_{h}^{\perp}\cdot\sigma_{h})e_{3}=\partial_{3}u_{h}-\nabla_{h}^{\perp}\theta+(\partial_{3}w)e_{3},$

with moreover

$\fint\sigma=\fint\left(-u_{h}^{\perp}+\theta e_{3}\right),$

as claimed. Finally we may compute that, since $\sigma$ is divergence-free and
hence $\int\sigma\cdot\nabla p=0$,

$\displaystyle E=\frac{1}{2}\int{\lvert u\rvert}^{2}+\theta^{2}$
$\displaystyle=\frac{1}{2}\int{\lvert u_{h}\rvert}^{2}+u_{3}^{2}+\theta^{2}$
$\displaystyle=\frac{1}{2}\int{\lvert\nabla_{h}p+\sigma_{h}\rvert}^{2}+w^{2}+{(\partial_{3}p+\sigma_{3})}^{2}$
$\displaystyle=\frac{1}{2}\int{\lvert\nabla
p\rvert}^{2}+{\lvert\sigma\rvert}^{2}+w^{2},$

as desired. ∎

The benefit of the approach carried out in this section, where we pull back a
metric from measurement space onto state space, is that it will immediately
produce a metric which agrees with the decomposition. The drawback of this
approach is that it does not, inherently, relate to the conserved energy in
any way. To see why a metric produced in this way will automatically agree
with the decomposition relies on the following observation: if two states have
the same $PV$ and $M$, then they must have the same balanced component. This
is proved in Lemma 5.4 below. Therefore projecting a state onto the balanced
set is the same as keeping its $PV$ and $M$ measurements fixed and sending its
remaining measurements, namely $j$, $w$, and $a$, to zero. This is why a
metric pulled back from the measurement space $\mathcal{M}$ will necessarily
agree with the decomposition. As mentioned above, this relies on the following
result.

######  Lemma 5.4 ($PV$ and $M$ characterize balanced components).

Given any $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ the unique balanced state in
$\mathcal{B}$ whose $PV$ and $M$ agree with the $PV$ and $M$ of $\mathcal{X}$
is precisely the balanced component $\mathcal{X}_{B}$.

###### Proof.

Let $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ and let
$\mathcal{Y}\in\mathcal{B}$ be a balanced state whose $PV$ and $M$ agrees with
those of $\mathcal{X}$, i.e.

$PV(\mathcal{X})=PV(\mathcal{Y})\text{ and }M(\mathcal{X})=M(\mathcal{Y}).$

We know from Proposition 3.1 that since $\mathcal{Y}$ is balanced it is
uniquely characterised by $p\in\mathring{H}^{1}$ and $M\in L^{2}$ via
$\mathcal{Y}=\Phi(p,\,M)$. Moreover we may compute that

$M(\Phi(p,\,M))=M\text{ and }PV(\Phi(p,\,M))=\Delta
p+\frac{1}{2}\partial_{3}{\min}_{0}\,(M-\partial_{3}p).$

By the well-posedness of nonlinear $PV$-and-$M$ inversion, we know that there
is a unique $p$ solving

$\Delta
p+\frac{1}{2}\partial_{3}{\min}_{0}\,(M-\partial_{3}p)=PV(\mathcal{Y})=PV(\mathcal{X}),$

and Proposition 3.1 tells us that the resulting state is precisely
$\mathcal{X}_{B}$. This verifies that indeed $\mathcal{Y}=\mathcal{X}_{B}$. ∎

Using Lemma 5.4 above we could then immediately deduce that a distance on
$\mathbb{L}^{2}_{\sigma}$ defined via

$\displaystyle
d\left(\mathcal{X}_{1},\,\mathcal{X}_{2}\right)^{2}\vcentcolon={\left\lvert\left\lvert\mathcal{M}(\mathcal{X}_{1})-\mathcal{M}(\mathcal{X}_{2})\right\rvert\right\rvert}_{\mathfrak{M}}^{2}$
$\displaystyle={\left\lvert\left\lvert{PV}_{1}-{PV}_{2}\right\rvert\right\rvert}_{H^{-1}}^{2}+{\left\lvert\left\lvert
M_{1}-M_{2}\right\rvert\right\rvert}_{L^{2}}^{2}+{\left\lvert\left\lvert
j_{1}-j_{2}\right\rvert\right\rvert}_{H^{-1}}^{2}+{\left\lvert\left\lvert
w_{1}-w_{2}\right\rvert\right\rvert}_{L^{2}}^{2}+{\lvert
a_{1}-a_{2}\rvert}^{2}$ (5.4)

would automatically agree with the decomposition. However this distance has no
chance of agreeing with the energy! This is because its dimensions, or units,
are all wrong. Indeed: this naive distance above treats $PV$ and $M$ on the
same footing even though they do not have the same dimensions. In order to
produce a metric on the state space $\mathbb{L}^{2}_{\sigma}$ pulled back from
the measurement space $\mathcal{M}$ which has a chance of agreeing with the
energy we must therefore _dimensionalize_ the naive distance from (5.4).

In order to perform this dimensionalisation we will use the fact that the
“fixed $H$” energy introduced in Definition 5.2 may sometimes be written not
in terms of the state $(u,\,\theta,\,q)$ but in terms of its coordinates
$(p,\,M,\,\sigma,\,w)$. To see this, we first define the following.

######  Definition 5.5.

Fix an indicator function $H$. We define, for any
$(p,\,M,\,\sigma,\,w)\in\mathcal{C}$,

$\widetilde{E}_{H}(p,\,M,\,\sigma,\,w)\vcentcolon=\frac{1}{2}\int\underbrace{A_{H}\nabla
p\cdot\nabla p}_{=\vcentcolon
Q_{H}^{+}(p)}+\left(1+\frac{H}{2}\right)M^{2}+\underbrace{A_{H}^{-1}\sigma\cdot\sigma}_{=\vcentcolon
Q_{H}^{-}(\sigma)}+w^{2},$

where $A_{H}=I-\frac{H}{2}e_{3}\otimes e_{3}$ and so
$A_{H}^{-1}=I+He_{3}\otimes e_{3}$.

######  Remark 5.6 (The elliptic matrix $A_{H}$ and its inverse).

The matrix $A_{H}$ introduced in Definition 5.5 above will be encountered many
times in the sequel. Where does it come from? The simplest explanation is that
it comes from the linearisation of nonlinear $PV$-and-$M$ inversion. Indeed,
upon linearising $PV$-and-$M$ inversion

$\Delta p+\frac{1}{2}\partial_{3}{\min}_{0}\,(M-\partial_{3}p)=PV$

we obtain, upon differentiating,

$\Delta\partial_{\alpha}\pi+\frac{1}{2}\mathds{1}(M<\partial_{3}p)(\partial_{\alpha}M-\partial_{\alpha}\partial_{3}\pi)=\partial_{\alpha}PV.$

In other words $\partial_{\alpha}p$ satisfies, for
$H\vcentcolon=\mathds{1}(M<\partial_{3}p)$,

$\nabla\cdot(A_{H}\nabla\partial_{\alpha}p)=\partial_{\alpha}PV-\frac{1}{2}H\partial_{\alpha}M.$

In particular we may readily verify that the identity

$\left(1+H\right)\left(1-\frac{H}{2}\right)=1$

holds for any $H$ satisfying $H^{2}=H$ (which is the case if $H$ is an
indicator function), and so the inverse of $A_{H}$ is indeed given by

$A_{H}^{-1}=I+He_{3}\otimes e_{3}.$

In general the form of the fixed Heaviside energy provided in Definition 5.5
does not agree with the form of the fixed Heaviside energy provided in
Definition 5.2. Sometimes, they do.

######  Proposition 5.7 (Rewriting the energy).

Let $E_{H}$ be as introduced in Definition 5.2 and $\widetilde{E}_{H}$ be as
in Definition 5.5. Let $(u,\,\theta,q)\in\mathbb{L}^{2}_{\sigma}$ be
corresponding coordinates $(p,\,M,\,\sigma,\,w)$ as per Corollary 3.6. In
general we have that

$E_{H}(u,\,\theta,\,q)=\frac{1}{2}\int{\lvert
u_{h}\rvert}^{2}+u_{3}^{2}+\left(1-\frac{H}{2}\right){(\theta-
qH)}^{2}+\left(1+\frac{H}{2}\right){(\theta+q)}^{2}.$

Moreover, if $q$ and its balanced component $q_{B}$ share the same interface,
i.e. if

$H\vcentcolon=\mathds{1}(q<0)=\mathds{1}(q_{B}<0),$

where note that $\mathds{1}(q_{B}<0)=\mathds{1}(M<\partial_{3}p)$, then

$E_{H}(u,\,\theta,\,q)=\widetilde{E}_{H}(p,\,M,\,\sigma,\,w).$

###### Proof.

The first identity follows from an immediate computation since $H^{2}=H$ and
hence

$\theta^{2}+q^{2}H=\left(1-\frac{H}{2}\right){(\theta-
qH)}^{2}+\frac{H}{2}{(\theta+q)}^{2}.$

To deduce the second identity we first note that Corollary 3.6 and Lemma A.7
tell us that, if
$H=H_{B}\vcentcolon=\mathds{1}(q_{B}<0)=\mathds{1}(M<\partial_{3}p)$, then

$\left\\{\begin{aligned} &u_{h}=\nabla_{h}^{\perp}p+\sigma_{h}^{\perp},\\\
&u_{3}=w,\\\ &\theta-qH=\partial_{3}p+(1+H)\sigma_{3},\text{ and }\\\
&\theta+q=M.\end{aligned}\right.$

Therefore

$\displaystyle\int{\lvert
u_{h}\rvert}^{2}+u_{3}^{2}+\left(1-\frac{H}{2}\right){(\theta-
qH)}^{2}+\left(1+\frac{H}{2}\right){(\theta+q)}^{2}$
$\displaystyle=\int{\lvert\nabla_{h}p+\sigma_{h}\rvert}^{2}+w^{2}+\left(1-\frac{H}{2}\right){\left(\partial_{3}p+(1+H)\sigma_{3}\right)}^{2}+\left(1+\frac{H}{2}\right)M^{2}.$

Finally, for $A_{H}$ and $A_{H}^{-1}$ as in Definition 5.5 we compute that,
since $A_{H}$ and $A_{H}^{-1}$ are symmetric and since $\sigma$ is divergence-
free,

$\displaystyle\int{\lvert\nabla_{h}p+\sigma_{h}\rvert}^{2}+\left(1-\frac{H}{2}\right){\left(\partial_{3}p+(1+H)\sigma\right)}^{2}$
$\displaystyle=\int A_{H}(\nabla p+A_{H}^{-1}\sigma)\cdot(\nabla
p+A_{H}^{-1}\sigma)$ $\displaystyle=\int A_{H}\nabla p\cdot\nabla
p+2\underbrace{\int\nabla p\cdot\sigma}_{=0}+\int
A_{H}^{-1}\sigma\cdot\sigma.$

The second identity follows. ∎

This is a good point to take a stock. Recall that we seek to define a metric
pulled back from the measurement space $\mathcal{M}$ onto the state space
$\mathbb{L}^{2}_{\sigma}$. Such a metric will automatically agree with the
decomposition but, in order to make it agree with the energy, we must
dimensionalize it appropriately. The identity recorded in Proposition 5.7
above tells us how to do this. The last remaining piece required to define our
pullback metric is thus to identify how the measurements ($PV$, $M$, $j$, $w$,
and $a$) will determine the coordinates ($p$, $M$, $\sigma$, and $w$). Note
that such a map is already implicit in Corollary 3.6, however we will do
something slightly different: since the “fixed $H$” energy already uses a
fixed Heaviside, we will leverage this to define a map taking measurements to
coordinates in a _linear_ way. This linearity will make the ensuing pullback
distance easier to work with.

######  Lemma 5.8 (Definition and invertibility of the fixed Heaviside
parametrisation $\mathcal{C}_{H}$).

Fix an indicator function $H$ and define
$A_{H}\vcentcolon=I-\frac{H}{2}e_{3}\otimes e_{3}$, such that
$A_{H}^{-1}=I+He_{3}\otimes e_{3}$ (since $(1-H/2)(1+H)\equiv 1$). We define
$\mathcal{C}_{H}:\mathfrak{M}\to\mathfrak{C}$, where $\mathfrak{C}$ is the
_coordinate space_ defined in (3.9). Given
$(PV,\,M,\,j,\,w,\,a)\in\mathfrak{M}$, let $M\vcentcolon=M$, $w\vcentcolon=w$,
and let $p$ and $\sigma$ solve

$\nabla\cdot(A_{H}\nabla p)=PV-\frac{1}{2}\partial_{3}(HM)$ (5.5)

and

$\nabla\times(A_{H}^{-1}\sigma)=j+(\partial_{3}w)e_{3}\text{ subject to
}\nabla\cdot\sigma=0\text{ and }\fint A_{H}^{-1}\sigma=a.$ (5.6)

The map $\mathcal{C}_{H}$ is well-defined and invertible.

Note that if
$H=H_{B}\vcentcolon=\mathds{1}(q_{B}<0)=\mathds{1}(M<\partial_{3}p)$ then
$\mathcal{C}_{H}$ is precisely the map implicit in Corollary 3.6 which assigns
the coordinates $p$, $M$, $\sigma$, and $w$ to a state $(u,\,\theta,\,q)$
through its measurements $PV$, $M$, $j$, $w$, and $a$.

###### Proof of Lemma 5.8.

It suffices to show that both (5.5) and (5.6) have a unique solution once
$PV$, $M$, $j$, $w$, and $a$ are fixed. The well-posedness of (5.5) follows
from the fact that $A_{H}$ is uniformly elliptic since
$A_{H}(x)\geqslant\frac{1}{2}I$ for all $x$. The well-posedness of (5.6)
follows from the fact that $A_{H}^{-1}$ is uniformly elliptic, since
$A_{H}^{-1}(x)\geqslant I$ for all $x$, and so Proposition A.8 applies to
(5.6). ∎

######  Remark 5.9.

In the definition of $\mathcal{C}_{H}$ we are linearising the coordinates, but
_not_ the measurements. This is because it is the full nonlinear measurements
which are well-adapted to the balanced $\mathcal{B}$ and to the wave set
$\mathcal{W}$, so changing those would mean that the resulting distance (in
Definition 5.10 below) would fail to be compatible with the energy. It should
not yet be clear why using linearised versions of the coordinates is
particularly helpful. This will be discussed in Section 5.3 below.

We are now ready to define a dimensionalized version of the naive pullback
metric of (5.4).

######  Definition 5.10 (Dimensionalized Parseval distance).

Fix an indicator function $H$. Define
$d_{H}:\mathbb{L}^{2}_{\sigma}\times\mathbb{L}^{2}_{\sigma}\to[0,\,\infty)$ as
follows:

$\frac{1}{2}{d_{H}(\mathcal{X}_{1},\,\mathcal{X}_{2})}^{2}=\widetilde{E}_{H}\left(\mathcal{C}_{H}\left[\mathcal{M}(\mathcal{X}_{1})-\mathcal{M}(\mathcal{X}_{2})\right]\right).$

In other words, given
$\mathcal{X}_{1},\,\mathcal{X}_{2}\in\mathbb{L}^{2}_{\sigma}$, let

$M\vcentcolon=M(\mathcal{X}_{1})-M(\mathcal{X}_{2}),\,w\vcentcolon=w(\mathcal{X}_{1})-w(\mathcal{X}_{2}),\text{
and }a\vcentcolon=a(\mathcal{X}_{1})-a(\mathcal{X}_{2}),$

let $p$ solve, for
$PV\vcentcolon={PV}(\mathcal{X}_{1})-{PV}(\mathcal{X}_{2})$,

$\nabla\cdot(A_{H}\nabla p)=PV-\frac{1}{2}\partial_{3}(HM)$

and let $\sigma$ solve, for
$j\vcentcolon=j(\mathcal{X}_{1})-j(\mathcal{X}_{2})$,

$\nabla\times(A_{H}^{-1}\sigma)=j+(\partial_{3}w)e_{3}\text{ subject to }\fint
A_{H}^{-1}\sigma=a\text{ and }\nabla\cdot\sigma=0.$

Then

$\frac{1}{2}{d_{H}(\mathcal{X}_{1},\,\mathcal{X}_{2})}^{2}=\int
Q_{H}^{+}(p)+\left(1+\frac{H}{2}\right)M^{2}+Q_{H}^{-}(\sigma)+w^{2}$

for $Q_{H}^{\pm}$ as introduced in Definition 5.5.

######  Proposition 5.11.

For any indicator function $H$ the function $d_{H}$ defined in Definition 5.10
is a metric which agrees with the extraction of a balanced component. Moreover
this distance is compatible with the energy when
$H=H_{1}=H_{B,\,1}=H_{2}=H_{B,\,2}$, i.e. in that case

$\displaystyle\frac{1}{2}{d_{H}\left(\mathcal{X}_{1},\,\mathcal{X}_{2}\right)}^{2}$
$\displaystyle=E_{H}(\mathcal{X}_{1}-\mathcal{X}_{2})$
$\displaystyle=\int{\lvert
u_{1}-u_{2}\rvert}^{2}+{(\theta_{1}-\theta_{2})}^{2}+{(q_{1}-q_{2})}^{2}H+{(M_{1}-M_{2})}^{2}$

where $H=\mathds{1}(q_{1}<0)=\mathds{1}(q_{2}<0)$.

###### Proof.

First we show that $d_{H}$ is a metric. Since $\widetilde{E}_{H}$ is positive-
definite and viewed here as taking arguments in $\mathfrak{C}$ it is a
quadratic form which induces an inner product, and hence a metric which we
denote $d_{E}$, on $\mathfrak{C}$. Since $\mathcal{C}_{H}$ is linear and since
both $\mathcal{C}_{H}$ and $\mathcal{M}$ are invertible (see Proposition 3.7
and Lemma 5.8 ) and hence injective it follows from Proposition A.13 that
$d_{H}$ is the pullback metric of $d_{E}$ under
$\mathcal{C}_{H}\circ\mathcal{M}$. So $d_{H}$ is indeed a metric on
$\mathbb{L}^{2}_{\sigma}$.

Now we show that $d_{H}$ conditionally agrees with the energy. Suppose that
$H=H_{1}=H_{B,\,1}=H_{2}=H_{B,\,2}$. The argument is easier to follow if we
rewrite $d_{H}$ as follows, using the fact that $\mathcal{C}_{H}$ is linear:

$\frac{1}{2}{d_{H}(\mathcal{X}_{1},\,\mathcal{X}_{2})}^{2}=\widetilde{E}_{H}\left[\mathcal{C}_{H}\left(\mathcal{M}\left(\mathcal{X}_{1}\right)\right)-\mathcal{C}_{H}\left(\mathcal{M}\left(\mathcal{X}_{2}\right)\right)\right].$

Let
$(p_{1},\,M_{1},\,\sigma_{1},\,w_{1})\vcentcolon=\mathcal{C}_{H}\left(\mathcal{M}\left(\mathcal{X}_{1}\right)\right)$.
Since $H=H_{B,\,1}=\mathds{1}(M_{1}<\partial_{3}p_{1})$ we see that

$\displaystyle\nabla\cdot\left(A_{H}\nabla
p_{1}\right)={PV}_{1}-\frac{1}{2}\partial_{3}(HM_{1})$
$\displaystyle\iff\Delta
p_{1}+\frac{1}{2}\partial_{3}{\min}_{0}\,(M_{1}-\partial_{3}p_{1})={PV}_{1}.$
(5.7)

Similarly, since $H_{1}=H_{B,\,1}$ we have that

$\displaystyle j_{1}$
$\displaystyle=\partial_{3}u_{h,\,1}-\nabla_{h}^{\perp}(\theta_{1}-{\min}_{0}\,q_{1})$
$\displaystyle=\partial_{3}u_{h,\,1}-\nabla_{h}^{\perp}(\theta_{1}-H_{B,\,1}q_{1})$

and so $\sigma_{1}$ solves

$\nabla\times(A_{H}^{-1}\sigma_{1})=\partial_{3}u_{h,\,1}-\nabla_{h}^{\perp}(\theta_{1}-H_{B,\,1}q_{1})+(\partial_{3}w_{1})e_{3}.$
(5.8)

We deduce from (5.7) and (5.8) that $p_{1}$ and $\sigma_{1}$ are precisely
characterized as in Theorem 3.5 by

$\mathcal{X}_{1}=\Phi(p_{1},\,M_{1})+\Psi(\sigma_{1},\,w_{1}).$

In other words: the approximate inversion $\mathcal{C}_{H}\circ\mathcal{M}$ is
_exact_ when $H=H_{1}=H_{B,\,1}$. Moreover Proposition 5.7 tells us that
$\widetilde{E}_{H}$ agrees with $E_{H}$ and with the nonlinear energy when
$H=H_{1}=H_{B,\,1}$ and so we conclude that indeed if
$H=H_{1}=H_{B,\,1}=H_{2}=H_{B,\,2}$ then

$\displaystyle\frac{1}{2}{d_{H}(\mathcal{X}_{1},\,\mathcal{X}_{2})}^{2}$
$\displaystyle=E_{H}(\mathcal{X}_{1}-\mathcal{X}_{2})$
$\displaystyle=\int{\lvert
u_{1}-u_{2}\rvert}^{2}+{(\theta_{1}-\theta_{2})}^{2}+{(q_{1}-q_{2})}^{2}H+{(M_{1}-M_{2})}^{2},$

where $H=\mathds{1}(q_{1}<0)=\mathds{1}(q_{2}<0)$.

Finally we show that $d_{H}$ agrees with the extraction of a balanced
component. We fix $\mathcal{X}\in\mathbb{L}^{2}_{\sigma}$ and let
$\mathcal{Y}\in\mathcal{B}$. Since $\mathcal{Y}\in\mathcal{B}$ we know from
Proposition 3.1 that

$j(\mathcal{Y})=0,\,w(\mathcal{Y})=0,\text{ and }a(\mathcal{Y})=0.$

Let us define

$M\vcentcolon=M(\mathcal{X})-M(\mathcal{Y})\text{ and
}PV\vcentcolon=PV(\mathcal{X})-PV(\mathcal{Y}),$ (5.9)

and let $p$ be the solution of

$\nabla\cdot(A_{H}\nabla p)=PV-\frac{1}{2}\partial_{3}(HM)$ (5.10)

while $\sigma$ is the solution of

$\nabla\times(A_{H}^{-1}\sigma)=j(\mathcal{X})+\left[\partial_{3}w(\mathcal{X})\right]e_{3}\text{
subject to }\fint A_{H}^{-1}\sigma=a(\mathcal{X})\text{ and
}\nabla\cdot\sigma=0.$

Then we may write the distance between $\mathcal{X}$ and $\mathcal{Y}$ via

$\frac{1}{2}{d_{H}(\mathcal{X},\,\mathcal{Y})}^{2}=\int
Q_{H}^{+}(p)+\left(1+\frac{H}{2}\right)M^{2}+Q_{H}^{-}(\sigma)+w^{2}$

for $Q_{H}^{\pm}$ as in Definition 5.5. Crucially: $\sigma$ and $w$ only
depend on $\mathcal{X}$! So we may simplify the minimization problem as
follows:

$\operatorname*{arg\,min}_{\mathcal{Y}\in\mathcal{B}}d_{H}(\mathcal{X},\,\mathcal{Y})=\operatorname*{arg\,min}_{\mathcal{Y}\in\mathcal{B}}\int
Q_{H}^{+}(p)+\left(1+\frac{H}{2}\right)M^{2}.$

In particular, since $Q_{H}^{+}(v)=A_{H}v\cdot v\geqslant\frac{1}{2}{\lvert
v\rvert}^{2}$, i.e. since $Q_{H}^{+}$ is positive-definite, the unique
minimizer corresponds to

$p=0\text{ and }M=0.$

Plugging this into (5.9), and then into (5.10), this means that
$M(\mathcal{X})=M(\mathcal{Y})$, and hence that
$PV(\mathcal{X})=PV(\mathcal{Y})$. In other words the unique minimizer is the
balanced state $\mathcal{Y}$ whose $PV$ and $M$ agree with $\mathcal{X}$. By
Lemma 5.4, this means precisely that $\mathcal{Y}=\mathcal{X}_{B}$, as
desired. ∎

This concludes our first attempt at quantifying how far apart two states are
in a way that agrees with

  1. 1.

the extraction of a balanced component and

  2. 2.

the conserved energy.

We see that the first of these two requirements was met by the distance
introduced in Definition 5.10 above. However the second requirement was only
met in the restricted setting where the fixed Heaviside agreed with the
Heavisides arising from both states being compared.

In some sense, that only the first requirement is satisfied is not surprising.
Indeed, the two distances introduced in this section are built on top of the
measurements from the global change of coordinates of Proposition 3.7 and
these measurements have been chosen precisely to be compatible with the
decomposition. However, these measurements have, a priori, little to do with
the conserved energy and so it is actually quite positive that the distance
built in this way is in any way comparable to the energy.

This motivates an alternate approach: what if, instead of using
decomposition–amenable measurements as a starting point, we used the energy
itself as a starting point to define a distance? This is the approach we take
in the following section.

###  5.2 Energy-centric approach

Recall that we are after a way to measure how different two states are.
Moreover we wish to do so in a manner which is compatible with both the
balanced-unbalanced decomposition underpinning Theorem 3.5 _and_ the conserved
energy. In the previous section we approached this problem with the first
requirement in mind, and eventually discovered that this approach failed to
(generically) satisfy the second requirement. In this section we proceed the
other way around: we enforce the second requirement by putting the conserved
energy front and center in the construction of a new distance. The question
will then be to identify whether or not this metric agrees with the
decomposition.

Also, recall the difficulty that arises in trying to define a distance or
metric or inner product that agrees with the energy in (5.1): the nonlinear
switch from the minimum function or Heaviside function. This nonlinear switch
is activated based on the total water $q$, and it is not clear how to define
the nonlinear switch when two different total water variables, $q_{1}$ and
$q_{2}$, are involved in quantifying the distance between two different
states. This difficulty and others are discussed in more detail in the text
preceding Section 5.1.

A key idea is helpful for overcoming the difficulties of the nonlinear switch:
in defining a metric, begin by working _pointwise in space_. Then, depending
on which phase that point belongs to, the conserved energy (5.1) determines a
natural inner product to use.

More precisely, this section proceeds as follows. We begin by defining a
metric _pointwise_ which measures the discrepancy between two states
$\mathcal{X}_{1}(x)$ and $\mathcal{X}_{2}(x)$ when evaluated at a specific
point $x$ of the domain $\mathbb{T}^{3}$. This pointwise metric is defined to
agree with the conserved energy density and so it is defined piecewise, taking
different forms depending on the sign of $q(x)$. We then record the properties
of this pointwise metric: it is well-defined, a metric, and agrees with the
energy unconditionally. We also record how to compute this distance explicitly
in a special case. We then make the key observation of this section: this
pointwise metric does _not_ agree with the extraction of a balanced component,
i.e. it does not agree with the decomposition. Finally we conclude this
section by defining the induced metric on the state space
$\mathbb{L}^{2}_{\sigma}$ obtained by integrating the pointwise metric over
the domain and explain why we do not expect this metric to agree with the
decomposition either.

As mentioned above, a key idea in the definition of this metric is to work
_pointwise in space_. In particular, depending on which phase that point
belongs to, the conserved energy (5.1) determines a natural inner product to
use. We introduce appropriate notation with this in mind.

######  Definition 5.12.

For any $S_{1}=(u_{1},\,\theta_{1},\,q_{1})\in\mathbb{R}^{5}$ and
$S_{2}=(u_{2},\,\theta_{2},\,q_{2})\in\mathbb{R}^{5}$, where
$u_{1},\,u_{2}\in\mathbb{R}^{3}$, we define the inner products

${\langle S_{1},\,S_{2}\rangle}_{u}\vcentcolon=u_{1}\cdot
u_{2}+\theta_{1}\theta_{2}+q_{1}q_{2}+(\theta_{1}+q_{1})(\theta_{2}+q_{2})$

and

${\langle S_{1},\,S_{2}\rangle}_{s}\vcentcolon=u_{1}\cdot
u_{2}+\theta_{1}\theta_{2}+(\theta_{1}+q_{1})(\theta_{2}+q_{2}).$

We denote by ${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{u}$
and ${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{s}$ the
corresponding norms.

The inner products introduced in Definition 5.12 are chosen precisely so that
the ensuing norms
${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{u}$ and
${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{s}$ agree with the
conserved energy density. Indeed, when $q<0$ we have that

${\left\lvert\left\lvert(u,\,\theta,\,q)\right\rvert\right\rvert}_{u}^{2}={\lvert
u\rvert}^{2}+\theta^{2}+q^{2}+{(\theta+q)}^{2}$

and when $q\geqslant 0$ we have that

${\left\lvert\left\lvert(u,\,\theta,\,q)\right\rvert\right\rvert}_{s}^{2}={\lvert
u\rvert}^{2}+\theta^{2}+{(\theta+q)}^{2},$

such that indeed they agree with the conserved energy density as recorded in
(5.1).

$q=0$$\theta$$q$$b=0$$b=-1$$b=1$$b=2$$M=3$$M=2$$M=1$$M=0$$M=-1$$M=-2$$\widehat{M}_{s}$$\widehat{b}_{s}$$\widehat{M}_{u}$$\widehat{b}_{u}$
Figure 1: The vectors $\left\\{\widehat{M}_{s},\,\widehat{b}_{s}\right\\}$ and
$\left\\{\widehat{M}_{u},\,\widehat{b}_{u}\right\\}$ form orthonormal bases
with respect to the inner products ${\langle\,\cdot\,,\,\cdot\,\rangle}_{s}$
and ${\langle\,\cdot\,,\,\cdot\,\rangle}_{u}$, respectively, as they were
introduced in Definition 5.12, namely $\widehat{M}_{s}=(0,\,1)$,
$\widehat{b}_{s}=(1,\,-1)$, $\widehat{M}_{u}=\frac{1}{\sqrt{6}}(1,\,1)$, and
$\widehat{b}_{u}=\frac{1}{\sqrt{2}}(1,\,-1)$. Even though, especially in the
saturated region $\left\\{q>0\right\\}$, labelling these vectors with the
moist variable $M$ and the buoyancy $b$ may at first appear perplexing, this
is justified by the fact that their corresponding _dual_ vectors are precisely
(multiples of) $M$ and $b$. Indeed:
${\langle\widehat{M}_{s},\,(\theta,\,q)\rangle}_{s}=\theta+q=M$,
${\langle\widehat{b}_{s},\,(\theta,\,q)\rangle}_{s}=\theta=b_{s}$,
${\langle\widehat{M}_{u},\,(\theta,\,q)\rangle}_{u}=\sqrt{\frac{3}{2}}(\theta+q)=\sqrt{\frac{3}{2}}M$,
and
${\langle\widehat{b}_{u},\,(\theta,\,q)\rangle}_{u}=\frac{1}{\sqrt{2}}\theta-q=\frac{1}{\sqrt{2}}b_{u}$.

Next we wish to define the _pointwise_ metric. As alluded to at the start of
this section, this metric will be defined _piecewise_ (depending on the sign
of $q$) and so we introduce notation which describes each of these pieces.

######  Definition 5.13.

We consider the following subsets of $\mathbb{R}^{5}$.

  * •

$\mathcal{U}\vcentcolon=\left\\{(u,\,\theta,\,q)\in\mathbb{R}^{5}:q<0\right\\}$
is the set of unsaturated states,

  * •

$\mathcal{S}\vcentcolon=\left\\{(u,\,\theta,\,q)\in\mathbb{R}^{5}:q\geqslant
0\right\\}$ is the set of saturated states, and

  * •

$\mathcal{I}\vcentcolon=\left\\{(u,\,\theta,\,q)\in\mathbb{R}^{5}:q=0\right\\}$
is the state-space interface between unsaturated states and saturated states.

We have the necessary notation in hand to define the following _Snell-type_
metric. We refer to this metric as being of Snell-type since shortest paths
between states in different phases exhibit the piecewise affine behaviour
observed by light passing across a boundary (see also Figure 2).

######  Definition 5.14 (Snell-type metric, pointwise).

We define $d_{\text{Snell}}:\mathbb{R}^{5}\to\mathbb{R}^{5}\to[0,\,\infty)$ as
follows.

  * •

If $S_{1},\,S_{2}\in\mathcal{U}$ then
$d_{\text{Snell}}(S_{1},\,S_{2})\vcentcolon={\left\lvert\left\lvert
S_{1}-S_{2}\right\rvert\right\rvert}_{u}$.

  * •

If $S_{1},\,S_{2}\in\mathcal{S}$ then
$d_{\text{Snell}}(S_{1},\,S_{2})\vcentcolon={\left\lvert\left\lvert
S_{1}-S_{2}\right\rvert\right\rvert}_{s}$.

  * •

If $S_{1}\in\mathcal{S}$ and $S_{2}\in\mathcal{U}$ then

$d_{\text{Snell}}\vcentcolon=\min_{I\in\mathcal{I}}\;{\left\lvert\left\lvert
S_{1}-I\right\rvert\right\rvert}_{s}+{\left\lvert\left\lvert
S_{2}-I\right\rvert\right\rvert}_{u}.$

  * •

If $S_{1}\in\mathcal{U}$ and $S_{2}\in\mathcal{S}$ then
$d_{\text{Snell}}(S_{1},\,S_{2})\vcentcolon=d_{\text{Snell}}(S_{2},\,S_{1})$.

$q=0$$\theta=-3,\,q=1$$\theta=1,\,q=-1$$\theta_{*}=-4+\frac{3}{2}\sqrt{2}\approx-1.4$
Figure 2: As introduced in Definition 5.14 we see that here that the “Snell-
type” metric warrants its name. Shown here is the geodesic between two points
in the $(\theta,\,q)$ plane which live on different sides of the state-space
interface $\left\\{q=0\right\\}$. As is familiar from Snell’s refraction law
of optics, the geodesic intersects the interface at an angle. Note that we may
use Lemma 5.17 in order to compute the intersection point and the value of
$\theta_{*}$ exactly.

First we verify that this distance function is well-defined.

######  Lemma 5.15 (The Snell-type distance is well-defined).

For any $S_{1}\in\mathcal{S}$ and any $S_{2}\in\mathcal{U}$ there is a unique
minimizer of

$\min_{I\in\mathcal{I}}{\left\lvert\left\lvert
S_{1}-I\right\rvert\right\rvert}_{s}+{\left\lvert\left\lvert
S_{2}-I\right\rvert\right\rvert}_{u}.$

###### Proof.

To verify that $d_{\text{Snell}}$ is well-defined it suffices to verify that,
for any $S_{1}\in\mathcal{S}$ and any $S_{2}\in\mathcal{U}$, there exists a
unique minimizer of

$\min_{I\in\mathcal{I}}{\left\lvert\left\lvert
S_{1}-I\right\rvert\right\rvert}_{s}+{\left\lvert\left\lvert
S_{2}-I\right\rvert\right\rvert}_{u}.$ (5.11)

To do so, fix $S\in\mathbb{R}^{5}\setminus\mathcal{I}$ and consider
$f:\mathbb{R}^{4}\to\mathbb{R}^{5}$ defined by, for every
$I=(u_{h},\,w,\,\theta,\,0)\in\mathcal{I}$,

$f(u_{h},\,w,\,\theta)={\left\lvert\left\lvert S-I\right\rvert\right\rvert}$

where ${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}$ denotes
${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{s}$ or
${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{u}$. To verify
that (5.11) has a unique minimizer it suffices to show that $f$ is strictly
convex. Since $f$ is smooth (as $S$ is not in $\mathcal{I}$, and $f$ and its
derivatives would blow up precisely when $S=I$, but $I\in\mathcal{I}$) it
suffices to verify that

$\nabla^{2}f>0.$

So we compute:

$\partial_{i}f=\partial_{i}{\left\lvert\left\lvert
S-I\right\rvert\right\rvert}=\frac{\langle
I-S,\,\partial_{i}I\rangle}{{\left\lvert\left\lvert
S-I\right\rvert\right\rvert}}$

and hence

$\partial_{i}\partial_{j}f=\partial_{j}\left(\frac{\langle
I-S,\,\partial_{i}I\rangle}{{\left\lvert\left\lvert
S-I\right\rvert\right\rvert}}\right)=\frac{\langle\partial_{i}I,\,\partial_{j}I\rangle{{\left\lvert\left\lvert
S-I\right\rvert\right\rvert}}^{2}-\langle\partial_{i}I,\,S-I\rangle\langle\partial_{j}I,\,S-I\rangle}{{\left\lvert\left\lvert
S-I\right\rvert\right\rvert}^{3}}.$

To make the remaining computation easier to follow, let us write

$x^{i}\vcentcolon=\partial_{i}I\text{ and }z\vcentcolon=S-I.$

Then

${\left\lvert\left\lvert z\right\rvert\right\rvert}^{3}\nabla^{2}f=\langle
x^{i},\,x^{j}\rangle{\left\lvert\left\lvert
z\right\rvert\right\rvert}^{2}-\langle x^{i},\,z\rangle\langle
x^{j},\,z\rangle=\vcentcolon A$

and so it suffices to show that $A>0$ (since ${\left\lvert\left\lvert
z\right\rvert\right\rvert}\neq 0$ as long as $S\notin\mathcal{I}$). To
simplify further, let us introduce

$D\vcentcolon=\text{diag}\,{\left({\left\lvert\left\lvert
x^{i}\right\rvert\right\rvert}\right)}_{i=1}^{4}.$

Then

$A>0\iff D^{T}AD>0$

where, since ${\left(x^{i}\right)}_{i=1}^{4}$ forms an orthogonal set with
respect to both
${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{s}$ and
${\left\lvert\left\lvert\,\cdot\,\right\rvert\right\rvert}_{u}$, we have that
$D^{T}AD=I{\left\lvert\left\lvert z\right\rvert\right\rvert}^{2}-B$ for

$B_{ij}\vcentcolon=\langle\hat{x}^{i},\,z\rangle\langle\hat{x}^{j},\,z\rangle.$

To conclude the proof it thus suffices to prove that
$B<{\left\lvert\left\lvert z\right\rvert\right\rvert}^{2}I$. First we show
that $B\leqslant{\left\lvert\left\lvert z\right\rvert\right\rvert}^{2}I$.