The equations of motion for most or all major ocean modeling software are implemented using a "Generalized Vertical Coordinate" (GVC). Generalized vertical coordinates contain a "fixed z" coordinate as a limiting case, but generalize to vertical coordinates that

• Are fixed in time but conform to topographic variations (σ coordinates);
• Vary in time in a diagnostic manner (for example, z-star coordinates, which vary with the sea surface displacement);
• Vary in time according to some prognostic equation (fully Lagrangian isopycnal coordinates, semi-Lagrangian hybrid coordinates / z-tilde coordinates that vary only with high frequency motion, and adaptive methods that prescribe GVC dynamics to obtain favorable properties like concentration in regions of strong stratification).

A fully general GVC typically also requires a "Lagrangian remapping" step to avoid extreme grid distortions in regions of persistent vertical velocities.

The implementation of GVC is likely a major refactoring of HydrostaticFreeSurfaceModel because it will change the equations of motion and could even potentially change the nature of its prognostic variables. For time-dependent GVC, tracer equations must be reformulated in terms of the "thickness-weighted" tracer concentration, which in our case means the tracer concentration normalized by the local grid spacing. This could mean either using the thickness-weighted tracer conservation as a state variable, or it could mean rewriting the time-stepping algorithm so that unweighted tracer conservation can be updated according to the conservation of thickness-weighted tracer.

A preliminary roadmap towards GVC in HydrostaticFreeSurfaceModel is

1. Introduce AbstractVerticalCoordinate and refactor the HydrostaticFreeSurfaceModel to integrate thickness-weighted equations. When using ZCoordinate, the resulting model produces identical results to the current implementation of HydrostaticFreeSurfaceModel. I think we will also want a rectilinear grid to use for testing to avoid modifying RegularRectilinearGrid, such as a refactored or reimplemented three-dimensionally StretchedRectilinearGrid from introducing a stretched grid in any or all directions #1532 which accepts the use of a generalized vertical coordinates.

2. Introduce ZStarCoordinate (or maybe FreeSurfaceWeightedVerticalCoordinate). "z-star" vertical coordinate is a relatively simple GVC, with a small, diagnostic time-dependence. A successful implementation will require integrating thickness-weighted equations, correct vertical derivatives, and ensuring correct horizontal pressure gradients.

3. Implement remapping as an alternative to vertical advection and test using ZCoordinate and ZStarCoordinate.

Of course, steps 3 might evolve depending on our experience during step 2.

An important choice that's part of step 1. is whether to use thickness-weighted variables as state variables rather than unweighted variables (eg hu rather than u, and hc rather than c. I suspect it makes sense to use thickness-weighted state variables because this simplifies the underlying implementation (for example, we only need hⁿ rather than both hⁿ⁺¹ and hⁿ to evolve a tracer field, and we can evaluate the prognostic equations in the same order that we do now), and because it is easy to provide abstractions for users to compute, inspect, and output unweighted variables u = hu / h. The pros and cons of this important choice should be discussed.

The discussion on #1549 is related (I decided to start a new issue because GVC are more general than terrain-following coordinates).

Some references

• MOM6 description (Adcroft et al 2019)
• Primer on “remap” methods (Griffies et al, 2019)
• “z-tilde” GVC used by the NEMO ocean model (Leclair and Madec 2011)
• Evaluation of various GVC in MPAS-Ocean (Petersen et al 2015)
• GVC in a nonhydrostatic compressible atmospheric model ENDgame (Kavcic and Thuburn, 2018)
• Correct calculation of pressure gradients in a finite volume model (Adcroft et al 2008)
  An “adaptive”, semi-Lagrangian GVC (Burchard and Beckers 2004) See also a 3D description (Hofmeister et al 2010) and realistic test case (Hofmeister et al 2011)
  “z star” GVC in MITgcm (Adcroft and Campin 2004)