Inconsistencies in evolution equation for density variables and misleading water mass transformations in density space

[hdrake]

The water mass transformation approach implemented here (following Groeskamp et al 2019 and described in more detail in Drake et al. in prep) requires, as a starting point, a locally-valid conservation equation for a tracer $\lambda$. For applications to ocean circulation, an insightful choice of $\lambda$ is seawater density because the gravitational force opposes motion across density surfaces and constrains the flow to largely follow surfaces of constant density.

The dependence of in-situ density on pressure means that in-situ density is not conserved for adiabatic displacement across isobars, motivating the definition of potential density (referenced to a given pressure). However, potential density surfaces are only approximately neutral within a few hundred dbar of their reference pressure and therefore can not be used for studies of the global ocean overturning circulation that spans the whole vertical water column. The current implementation (circa 07b2e6f and #40) uses potential density and thus a limitation of the results is that some unknown fraction of the apparent diapycnal water mass transformations that occur far from the reference pressure are in fact along neutral surfaces (and vice versa). It also seems that there may be a conceptual error in the current implementation in that the thermal expansion and haline contraction coefficients that are being used are referenced to the local pressure, whereas a fixed reference pressure is used for the density variable itself.

While the non-linear equation of state makes it impossible to define a density variable that is everywhere neutral, approximately neutral surfaces can be constructed that stitch together locally neutral tangent plans to minimize non-neutrality (Stanley et al. 2021; @geoffstanley's neutralocean python package). One can also write (and evaluate) an approximate conservation equation for this such a neutral density variable, which only requires multiplying in-situ density tendencies by an O(1) correction factor $b \equiv |\nabla \gamma^{n}|/ |\nabla \rho|$. I argue that this should be the favored approach for density-space calculations in xwmt (and dependent packages).

[hdrake]

After conversations with Jan Zika and Geoff Stanley, I am convinced that using locally referenced thermal expansion and haline contraction coefficients for water mass transformations in potential density coordinates is incorrect. I have corrected this in hdrake@0a25a53 but this means that my answers in #40 are now substantially changed relative to @jetesdal's methods in both the main branch 07b2e6f and the 2D+3D merged branch 797cd88.