Abstract

Holmes et al. (2019) have proposed a new theoretical framework for studying ocean heat uptake in potential temperature coordinates. One important step in their derivations requires understanding the temporal changes of the volume of water V with temperature greater than some value, which they write as the sum of two terms. The first one is due to the surface freshwater fluxes and is well defined, but the second one—attributed to the volume fluxes through the lower boundary of the domain—is given no explicit expression. What the authors mean exactly is unclear, however, because in the incompressible Boussinesq approximation, the use of a divergenceless velocity field implies that the sum of the volume fluxes through any kind of control volume must integrate to zero at all times. In this comment, we provide two alternative explicit mathematical expressions linking the volume change of Holmes et al. (2019) to the diabatic sources and sinks of heat that clarify their result. By contrasting Holmes et al.’s (2019) approach with that for a fully compressible ocean, it is concluded that the volume considered by Holmes et al. (2019) is best interpreted as a proxy for the Boussinesq mass M0 = ρ0V, where ρ0 is the reference Boussinesq density. If V were truly meant to represent volume rather than a proxy for the Boussinesq mass, the Boussinesq expression for dV/dt would have to be regarded as inaccurate because of its neglect of the volume changes resulting from mean density changes.

1. Introduction

Holmes et al. (2019) have recently published an interesting study of ocean heat uptake in potential temperature coordinates. A key part of the exercise involves deriving a theoretical expression for the time variations of the volume of water associated with different temperature classes. In Holmes et al. (2019), the issue is addressed in the context of the incompressible Boussinesq primitive equations. As is well known, such a system of equations conserves volume rather than mass. Alternatively, such a system can also be regarded as conserving the Boussinesq mass ρ0dV, but not the “true” mass ρdV, where ρ0 and ρ are the reference Boussinesq and fluid densities, respectively. As a result, one of its peculiar properties is that the volume flux integrated over the boundary V of any control volume V vanishes identically at all times, raising the question of how the control volume V varies with time. Mathematically,

 
VvndS=VvdV=0,
(1)

where n is the outward unit normal vector and v the 3D velocity.

In the particular case where V represents the volume of all water warmer than a given temperature θ, Holmes et al. (2019) suggest that the time variation of V should be given by the formula

 
Vt=G+Js,
(2)

[their Eq. (3)]. As shown in Fig. 1, Holmes et al. (2019) refer to G as “the volume flux across the Θ isotherm, or the water-mass transformation” and to JS as the surface volume flux. They provide the following mathematical expression for JS:

 
Js(θ,t)=θ(x,y,0,t)>θPER(x,y,t)dA,
(3)

[their Eq. (4)], where PER is described as “the net volume flux per unit area (m s−1) into the ocean associated with precipitation, evaporation, river runoff, and ice melt.” However, they did not provide any explicit mathematical expression for G.

Fig. 1.

(left) Schematic showing the volume of water V(θ,t) with temperature larger than θ, the surface volume flux Js, and the boundary volume flux G. The surface boundary of volume V is decomposed into two parts: Sθ and Aθ. Here, Sθ is where θ=constant and is shown in blue, while Aθ is at the ocean surface and is shown in red. (right) Adiabatic rearrangement of all parcels in the physical space so that θ surfaces in the reference space are flat. The reference depth zr is the depth associated with every θ value so that V(θ,t)=V^(zr). The term θ can be written as a function of zr and t: θ(x,y,z,t)=θr(zr,t). Thus, in the physical space, iso-θ surfaces are iso-zr surfaces.

Fig. 1.

(left) Schematic showing the volume of water V(θ,t) with temperature larger than θ, the surface volume flux Js, and the boundary volume flux G. The surface boundary of volume V is decomposed into two parts: Sθ and Aθ. Here, Sθ is where θ=constant and is shown in blue, while Aθ is at the ocean surface and is shown in red. (right) Adiabatic rearrangement of all parcels in the physical space so that θ surfaces in the reference space are flat. The reference depth zr is the depth associated with every θ value so that V(θ,t)=V^(zr). The term θ can be written as a function of zr and t: θ(x,y,z,t)=θr(zr,t). Thus, in the physical space, iso-θ surfaces are iso-zr surfaces.

Based on the fact that they refer to G as being due to the volume flux across the isotherm, Holmes et al. (2019) may have originally assumed (as happened to us and a few other colleagues) that the time variation of V should be governed by the following equation:

 
Vt=VvndS(WRONG!),
(4)

only to realize that the constraint (1) would yield the physically implausible result V/t=0, in clear contradiction with our physical intuition that V should, for instance, increase in an ocean experiencing net warming and decrease in an ocean experiencing net cooling.

A survey of the literature reveals that the above ambiguity can actually be traced back to Walin [1982, Eq. (2.2)] where, as above, the time variation of the volume is linked to the volume fluxes through its boundaries. In this comment, we seek to clarify the physics of the time variations of V by deriving an explicit mathematical expression for the term G. In particular, we aim to show that while JS is indeed related to the volume flux vndS through the ocean surface, this is not the case of the term G. Indeed, G is found to be related to the water mass transformation due to the diabatic sources and sinks of heat and salt, as expected from physical intuition, and also stated by Holmes et al. (2019). The derivation of such an expression is nontrivial, which might explain why it does not appear to have been published in the water-mass conversion literature before.

2. Theory

a. Linking G to water mass transformations

To link G to water mass transformations, and to evaluate the validity and accuracy of the Boussinesq form of the results obtained, we seek expressions valid for a fully compressible ocean first. Thus, the conservation equations for mass and heat that we take as our starting point are

 
ρt+(ρv)=0,
(5)
 
t(ρcpθ)+(ρcpθv)=(ρFθ),
(6)

where cp is a constant heat capacity and Fθ is the heat flux vector. In keeping with standard modeling practice, Eq. (6) assumes our definition of heat (cpθ) to be exactly conservative. Alternative and more accurate treatments would either entail the use of the McDougall (2003) Conservative Temperature or retaining the nonconservative production/destruction of θ as proposed by Tailleux (2015). Note that Eq. (6) implies for the Lagrangian derivative of θ:

 
DθDt=θ˙=1ρcp(ρFθ).
(7)

To obtain results pertaining to the incompressible Boussinesq approximation, one may simply replace ρ by the reference Boussinesq density ρ0 in any expression obtained from Eqs. (5) or (6).

To address the problem for a fully compressible ocean, we find it necessary to consider the mass M(θ,t) of the water masses of potential temperature greater than θ in addition to their volume V(θ,t). At the top, these water masses are bounded by the ocean free surface of equation z=η(x,y,t) and at the bottom, by the isothermal surface θ=constant of equation z=h(x,y,θ,t), denoted by Sθ in the following. As a result, the unit normal vectors at the top and bottom are respectively given by

 
n=kzη1+|zη|2,atz=η(x,y,t),and
(8)
 
n=θ|θ|=k+zh1+|zh|2,atz=h(x,y,θ,t),
(9)

where z denotes the horizontal gradient. At the ocean surface, the boundary conditions for mass and heat are

 
ρs(ηt+uszηws)=ρf(PE+R),and
(10)
 
ρsFθndS=Qnetdxdy,
(11)

where ρs=ρ(T,S,pa) is the surface density of seawater, ρf=ρ(T,0,pa) is the density of freswhater, Qnet is the net heat flux entering the ocean, while vs=(us,υs,ws) denotes the surface value of the velocity field. At the bottom, differentiating z=h(x,y,θ,t) yields

 
ht+ubzh+wb=hθDθDt=(θz)1DθDt,
(12)

where vb=(ub,υb,wb) is the value of the velocity field along the isothermal surface θ=constant.

By definition, the volume and mass of the water masses of potential temperature greater than θ can be written as

 
V(θ,t)=Aθh(x,y,θ,t)η(x,y,t)dzdxdy,and
(13)
 
M(θ,t)=Aθh(x,y,θ,t)η(x,y,t)ρdzdxdy,
(14)

where Aθ denotes the part of the ocean surface area capping V(θ,t) at its top. Taking the time derivative of Eq. (14), making use of Eqs. (5), (10), and (12), yields after some manipulation:

 
Mt=Aθ(ρsηt+ρbht)dxdyVρvndS=Aθρs(ηt+uszηws)dxdy+Aθρb(ht+ubzh+wb)dxdy=Aθρf(PE+R)dxdy+Aθρb(θz)1θ˙dxdy.
(15)

Equation (15) is a key result stating that only two physical processes can change M(θ,t) with time, namely, surface freshwater fluxes or diabatic modifications of θ along the isothermal lower surface Sθ. To derive an expression for the temporal evolution of V(θ,t), it is useful to define the volume mean density ρ¯(θ,t)=M(θ,t)/V(θ,t), which upon time differentiation can be shown to imply

 
Vt=1ρ¯MtVρ¯ρ¯t.
(16)

Equation (16) states that in a fully compressible ocean, the volume V can change either as the result of an addition/subtraction of mass or due to a change in the mean density ρ¯. In the context of sea level change arising from global warming, the first effect is generally associated with land ice melting and the second effect to thermal expansion, both being known to contribute O(1) mm yr−1 to the globally averaged sea level. This means that the two terms in Eq. (16) are often of comparable importance, and hence that it is in general not possible to justify neglecting the second term.

Now, replacing ρ by ρ0 in the above expressions yields the following Boussinesq limits:

 
MtMt|Boussinesq=ρ0[JS+Aθ(θz)1θ˙dxdyG],and
(17)
 
Vt1ρ0Mt|Boussinesq=JS+Aθ(θz)1θ˙dxdyG,
(18)

where JS is defined as in Holmes et al. (2019). Based on the above considerations, it is clear that only Eq. (17) remains accurate in the Boussinesq limit, since it is, in general, inaccurate to neglect volume changes due to mean density changes. In fact, it is well accepted that Boussinesq ocean models that aim to predict sea level change must somehow account for the latter effect, often using a procedure based on that proposed by Greatbatch (1994), for instance. At this stage, it is important to point out that although Holmes et al. (2019) interpret their equation for V/t as pertaining to the volume V(θ,t), it is equally possible to interpret it as an equation pertaining to the Boussinesq mass ρ0V instead. Doing so seems more logical, since it seems obvious that Holmes et al.’s (2019) framework physically relies on combining the mass and heat budgets, rather than the volume and heat budgets. It is such an interpretation that is assumed in the following, which allows us to stop worrying about the possible impacts of neglecting the volume changes due to mean density changes on ocean heat uptake. As a result, the expression for G in either Eq. (17) or (18) represents the desired expression explicitly linking G to water mass transformations, as postulated but not demonstrated by Holmes et al. (2019).

b. Link to volume-integrated diabatic processes

Although Eq. (18) explicitly links G to water mass transformations, it is arguably impractical for diagnosing G from ocean model outputs owing to its dependence on the boundary values of Dθ/Dt along the isothermal surface θ=constant. This is why in the following we seek a more practical alternative expression linking G to volume-integrated diabatic effects instead. For simplicity, we restrict our discussion to the case of an incompressible Boussinesq ocean, as in Holmes et al. (2019). To that end, we use a pdf approach that physically amounts to sorting the ocean according to potential temperature. The underlying idea is to define a reference potential temperature profile θr(z,t) and reference depth zr=zr(θ,t) through the following relations:

 
V(θ,t)=Aθh(x,y,θ,t)η(x,y,t)dzdxdy=zr(θ,t)0A(z)dz,and
(19)
 
Aθh(x,y,θ,t)η(x,y,t)θ(x,y,z,t)dzdxdy=zr(θ,t)0A(z)θr(z,t)dz,
(20)

where A(z) denotes the area of the ocean at the depth z. By construction, θr(z,t) and zr(θ,t) satisfy θr(zr,t)=θ at all times. Such a property defines a one-to-one relation between θ and zr, allowing one to regard any function F(θ,t) of θ and time t alternatively as a function of zr and t through the identity F(θ,t)=F(θr(zr,t),t)=F^(zr,t). For instance, V(θ,t)=V^(zr) or h(x,y,θ,t)=h^(x,y,zr,t), with the hat being used to denote the zr-based representation. An advantage of the zr representation is that V^(zr) is independent of time at fixed zr, which implies

 
V^t|zr=Vθθrt(zr,t)+Vt=0,
(21)

or alternatively,

 
Vt|θ=Vθθrt(zr,t).
(22)

Equation (22) shows that an alternative approach to deriving an expression for V/t is via deriving an expression for θr/t. This is achieved here by differentiating Eq. (20) with time at fixed zr; after making use of the Green theorem and of the boundary conditions for heat and mass, one eventually arrives at

 
zr0A(z)θrt(z,t)dz=tAθh^(x,y,zr,t)η(x,y,t)θdzdxdy
 
=Aθθs(EPR)dxdyPr(zr,t)θSθvndS+AθQnetρ0cpdxdyH(zr,t)SθFθcpndS,
(23)

where θs=θ(x,y,z=η,t) is the surface value of θ. As stated in the introduction, the use of divergenceless velocity field v=0 imposes at each time

 
SθvndSJS(θ,t)=0.
(24)

In Eq. (20), the last term represents the diathermal heat flux due to the parameterized mixing processes through the isothermal surface Sθ. Similarly as in Hochet et al. (2019), we find it useful to represent it as downgradient diffusion in terms of an effective diffusivity Keff such that

 
SθFθndS=KeffcpA(zr)θrz.
(25)

By making use of Eqs. (24) and (25), Eq. (23) may thus be rewritten as

 
zr0A(z)θrt(z,t)dz=PrθJS+HKeffA(zr)θrz.
(26)

Our sought-for evolution equation for θr is then simply obtaining by differentiating Eq. (26) with respect to zr. After some rearrangement, this yields

 
θrt=1A(zr)(HPrθJS)zr+1A(zr)zr[KeffA(zr)θrzr].
(27)

Interestingly, note that Eq. (27) can be written as a classical vertical advection diffusion equation,

 
θrt+weffθrz=1A(zr)zr[KeffA(zr)θrzr],
(28)

by introducing the pseudo effective vertical velocity weff as follows:

 
weff(zr,t)=1A(zr)(HPrθJS)θ.
(29)

To make the connection with Holmes et al. (2019) results, we may use the fact that

 
A(zr)=Vzr=Vθθrzr,
(30)

combined with Eqs. (27) and (22) to show that

 
Vt|θ=(HPrθJS)θ+θ(KeffAθrzr),
(31)

thus implying for G,

 
G=Vt|θJS=Hθ+θ(KeffAθrzr)+[(Pr+θJS)θJS]=0.
(32)

Now, Eq. (32) is directly comparable with Eq. (14) of Holmes et al. (2019), namely,

 
GHolmes=1ρ0cp[Fθ+Mθ+Iθ].
(33)

In Holmes et al. (2019), the term F in Eq. (33) is linked to the total surface flux, which in our expression is linked to H, while both M and I are related to the parameterized and numerical mixing, which in our expression appears as an effective diffusive flux. The main advantage of Eq. (32) is that all of its terms are arguably more familiar and easier to diagnose from numerical model outputs than local values of θ˙ along isothermal surfaces Sθ for which we have less physical intuition.

3. Conclusions

In this comment, we derived two mathematically equivalent expressions [Eqs. (18) and (32)] for the term G entering Eq. (3) of Holmes et al. (2019) governing the time variations of V(θ,t), which hopefully can help clarify a confusion in the theory of water masses dating back to Walin (1982). However, from the consideration of the fully compressible case, we believe that the evolution equation for V is best interpreted as pertaining to the Boussinesq mass M0=ρ0V rather than the volume V itself, since the boundary conditions that enter the problem belong to the mass budget. If the expression for V/t discussed by Holmes et al. (2019) were truly for the volume itself, it would be arguably inaccurate owing to its neglect of the volume changes due to mean density changes, which is often of comparable importance to that due to mass changes.

Acknowledgments

This work was supported by the OUTCROP NERC Grant NE/R010536/1. The comments of a reviewer helped clarify and simplify the derivations of this paper.

REFERENCES

REFERENCES
Greatbatch
,
R. J.
,
1994
:
A note on the representation of steric sea level in models that conserve volume rather than mass
.
J. Geophys. Res.
,
99
,
12 767
12 771
, https://doi.org/10.1029/94JC00847.
Hochet
,
A.
,
R.
Tailleux
,
D.
Ferreira
, and
T.
Kuhlbrodt
,
2019
:
Isoneutral control of effective diapycnal mixing in numerical ocean models with neutral rotated diffusion tensors
.
Ocean Sci.
,
15
,
21
32
, https://doi.org/10.5194/os-15-21-2019.
Holmes
,
R. M.
,
J. D.
Zika
, and
M. H.
England
,
2019
:
Diathermal heat transport in a global ocean model
.
J. Phys. Oceanogr.
,
49
,
141
161
, https://doi.org/10.1175/JPO-D-18-0098.1.
McDougall
,
T. J.
,
2003
:
Potential enthalpy: A conservative oceanic variable for evaluating heat content and heat fluxes
.
J. Phys. Oceanogr.
,
33
,
945
963
, https://doi.org/10.1175/1520-0485(2003)033<0945:PEACOV>2.0.CO;2.
Tailleux
,
R.
,
2015
:
Observational and energetics constraints on the non-conservation of potential/conservative temperature and implications for ocean modelling
.
Ocean Modell.
,
88
,
26
37
, https://doi.org/10.1016/j.ocemod.2015.02.001.
Walin
,
G.
,
1982
:
On the relation between sea-surface heat flow and thermal circulation in the ocean
.
Tellus
,
34
,
187
195
, https://doi.org/10.3402/tellusa.v34i2.10801.

Footnotes

The original article that was the subject of this comment/reply can be found at http://journals.ametsoc.org/doi/abs/10.1175/JPO-D-18-0098.1.

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