Changeset 11608


Ignore:
Timestamp:
2019-09-27T12:24:49+02:00 (22 months ago)
Author:
nicolasmartin
Message:

Review of "Model Basics" chapter

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/doc/latex/NEMO/subfiles/chap_model_basics.tex

    r11598 r11608  
    1616    Release & Author(s) & Modifications \\ 
    1717    \hline 
    18     {\em   4.0} & {\em ...} & {\em ...} \\ 
    19     {\em   3.6} & {\em ...} & {\em ...} \\ 
    20     {\em   3.4} & {\em ...} & {\em ...} \\ 
    21     {\em <=3.4} & {\em ...} & {\em ...} 
     18    {\em   4.0} & {\em Mike Bell                       } & {\em Update       } \\ 
     19    {\em   3.6} & {\em Gurvan Madec                    } & {\em Minor changes} \\ 
     20    {\em <=3.4} & {\em Gurvan Madec and Steven Alderson} & {\em First version} \\ 
    2221  \end{tabularx} 
    2322} 
     
    3837plus the following additional assumptions made from scale considerations: 
    3938 
    40 \begin{enumerate} 
    41 \item \textit{spherical Earth approximation}: the geopotential surfaces are assumed to be oblate spheriods 
    42   that follow the Earth's bulge; these spheroids are approximated by spheres with 
    43   gravity locally vertical (parallel to the Earth's radius) and independent of latitude 
     39\begin{labeling}{Neglect of additional Coriolis terms} 
     40\item [\textit{Spherical Earth approximation}] The geopotential surfaces are assumed to 
     41  be oblate spheroids that follow the Earth's bulge; 
     42  these spheroids are approximated by spheres with gravity locally vertical 
     43  (parallel to the Earth's radius) and independent of latitude 
    4444  \citep[][section 2]{white.hoskins.ea_QJRMS05}. 
    45 \item \textit{thin-shell approximation}: the ocean depth is neglected compared to the earth's radius 
    46 \item \textit{turbulent closure hypothesis}: the turbulent fluxes 
     45\item [\textit{Thin-shell approximation}] The ocean depth is neglected compared to the earth's radius 
     46\item [\textit{Turbulent closure hypothesis}] The turbulent fluxes 
    4747  (which represent the effect of small scale processes on the large-scale) 
    4848  are expressed in terms of large-scale features 
    49 \item \textit{Boussinesq hypothesis}: density variations are neglected except in their contribution to 
    50   the buoyancy force 
     49\item [\textit{Boussinesq hypothesis}] Density variations are neglected except in 
     50  their contribution to the buoyancy force 
    5151  \begin{equation} 
    5252    \label{eq:MB_PE_eos} 
    5353    \rho = \rho \ (T,S,p) 
    5454  \end{equation} 
    55 \item \textit{Hydrostatic hypothesis}: the vertical momentum equation is reduced to a balance between 
    56   the vertical pressure gradient and the buoyancy force 
     55\item [\textit{Hydrostatic hypothesis}] The vertical momentum equation is reduced to 
     56  a balance between the vertical pressure gradient and the buoyancy force 
    5757  (this removes convective processes from the initial Navier-Stokes equations and so 
    5858  convective processes must be parameterized instead) 
     
    6161    \pd[p]{z} = - \rho \ g 
    6262  \end{equation} 
    63 \item \textit{Incompressibility hypothesis}: the three dimensional divergence of the velocity vector $\vect U$ 
    64   is assumed to be zero. 
     63\item [\textit{Incompressibility hypothesis}] The three dimensional divergence of 
     64  the velocity vector $\vect U$ is assumed to be zero. 
    6565  \begin{equation} 
    6666    \label{eq:MB_PE_continuity} 
    6767    \nabla \cdot \vect U = 0 
    6868  \end{equation} 
    69 \item \textit{Neglect of additional Coriolis terms}: the Coriolis terms that vary with the cosine of latitude are neglected. 
    70   These terms may be non-negligible where the Brunt-Vaisala frequency $N$ is small, either in the deep ocean or 
    71   in the sub-mesoscale motions of the mixed layer, or near the equator \citep[][section 1]{white.hoskins.ea_QJRMS05}. 
    72   They can be consistently included as part of the ocean dynamics \citep[][section 3(d)]{white.hoskins.ea_QJRMS05} and are 
    73   retained in the MIT ocean model. 
    74 \end{enumerate} 
     69\item [\textit{Neglect of additional Coriolis terms}] The Coriolis terms that vary with 
     70  the cosine of latitude are neglected. 
     71  These terms may be non-negligible where the Brunt-V\"{a}is\"{a}l\"{a} frequency $N$ is small, 
     72  either in the deep ocean or in the sub-mesoscale motions of the mixed layer, 
     73  or near the equator \citep[][section 1]{white.hoskins.ea_QJRMS05}. 
     74  They can be consistently included as part of the ocean dynamics 
     75  \citep[][section 3(d)]{white.hoskins.ea_QJRMS05} and are retained in the MIT ocean model. 
     76\end{labeling} 
    7577 
    7678Because the gravitational force is so dominant in the equations of large-scale motions, 
     
    7880$k$ is the local upward vector and $(i,j)$ are two vectors orthogonal to $k$, 
    7981\ie\ tangent to the geopotential surfaces. 
    80 Let us define the following variables: $\vect U$ the vector velocity, $\vect U = \vect U_h + w \, \vect k$ 
     82Let us define the following variables: 
     83$\vect U$ the vector velocity, $\vect U = \vect U_h + w \, \vect k$ 
    8184(the subscript $h$ denotes the local horizontal vector, \ie\ over the $(i,j)$ plane), 
    8285$T$ the potential temperature, $S$ the salinity, $\rho$ the \textit{in situ} density. 
     
    8689  \label{eq:MB_PE} 
    8790  \begin{gather} 
    88     \intertext{$-$ the momentum balance} 
     91    \shortintertext{$-$ the momentum balance} 
    8992    \label{eq:MB_PE_dyn} 
    90     \pd[\vect U_h]{t} = - \lt[ (\nabla \times \vect U) \times \vect U + \frac{1}{2} \nabla \lt( \vect U^2 \rt) \rt]_h 
    91                         - f \; k \times \vect U_h - \frac{1}{\rho_o} \nabla_h p 
    92                         + \vect D^{\vect U} + \vect F^{\vect U} \\ 
    93     \intertext{$-$ the heat and salt conservation equations} 
     93    \pd[\vect U_h]{t} = - \lt[ (\nabla \times \vect U) \times \vect U + \frac{1}{2} \nabla \lt( \vect U^2 \rt) \rt]_h - f \; k \times \vect U_h - \frac{1}{\rho_o} \nabla_h p + \vect D^{\vect U} + \vect F^{\vect U} 
     94    \shortintertext{$-$ the heat and salt conservation equations} 
    9495    \label{eq:MB_PE_tra_T} 
    9596    \pd[T]{t} = - \nabla \cdot (T \ \vect U) + D^T + F^T \\ 
     
    102103(\autoref{eq:MB_PE_eos}), $\rho_o$ is a reference density, $p$ the pressure, 
    103104$f = 2 \vect \Omega \cdot k$ is the Coriolis acceleration 
    104 (where $\vect \Omega$ is the Earth's angular velocity vector), and $g$ is the gravitational acceleration. 
     105(where $\vect \Omega$ is the Earth's angular velocity vector), 
     106and $g$ is the gravitational acceleration. 
    105107$\vect D^{\vect U}$, $D^T$ and $D^S$ are the parameterisations of small-scale physics for momentum, 
    106108temperature and salinity, and $\vect F^{\vect U}$, $F^T$ and $F^S$ surface forcing terms. 
    107 Their nature and formulation are discussed in \autoref{sec:MB_zdf_ldf} and \autoref{subsec:MB_boundary_condition}. 
     109Their nature and formulation are discussed in \autoref{sec:MB_zdf_ldf} and 
     110\autoref{subsec:MB_boundary_condition}. 
    108111 
    109112%% ================================================================================================= 
     
    116119where $H$ is the depth of the ocean bottom and $\eta$ is the height of the sea surface 
    117120(discretisation can introduce additional artificial ``side-wall'' boundaries). 
    118 Both $H$ and $\eta$ are referenced to a surface of constant geopotential (\ie\ a mean sea surface height) on which $z = 0$. 
    119 (\autoref{fig:MB_ocean_bc}). 
    120 Through these two boundaries, the ocean can exchange fluxes of heat, fresh water, salt, and momentum with 
    121 the solid earth, the continental margins, the sea ice and the atmosphere. 
    122 However, some of these fluxes are so weak that even on climatic time scales of thousands of years 
    123 they can be neglected. 
     121Both $H$ and $\eta$ are referenced to a surface of constant geopotential 
     122(\ie\ a mean sea surface height) on which $z = 0$ (\autoref{fig:MB_ocean_bc}). 
     123Through these two boundaries, the ocean can exchange fluxes of heat, fresh water, salt, 
     124and momentum with the solid earth, the continental margins, the sea ice and the atmosphere. 
     125However, some of these fluxes are so weak that 
     126even on climatic time scales of thousands of years they can be neglected. 
    124127In the following, we briefly review the fluxes exchanged at the interfaces between the ocean and 
    125128the other components of the earth system. 
    126129 
    127 \begin{figure}[!ht] 
     130\begin{figure} 
    128131  \centering 
    129132  \includegraphics[width=0.66\textwidth]{Fig_I_ocean_bc} 
     
    136139 
    137140\begin{description} 
    138 \item [Land - ocean interface:]  the major flux between continental margins and the ocean is a mass exchange of fresh water through river runoff. 
     141\item [Land - ocean] The major flux between continental margins and the ocean is a mass exchange of 
     142  fresh water through river runoff. 
    139143  Such an exchange modifies the sea surface salinity especially in the vicinity of major river mouths. 
    140   It can be neglected for short range integrations but has to be taken into account for long term integrations as 
     144  It can be neglected for short range integrations but 
     145  has to be taken into account for long term integrations as 
    141146  it influences the characteristics of water masses formed (especially at high latitudes). 
    142147  It is required in order to close the water cycle of the climate system. 
    143   It is usually specified as a fresh water flux at the air-sea interface in the vicinity of river mouths. 
    144 \item [Solid earth - ocean interface:]  heat and salt fluxes through the sea floor are small, except in special areas of little extent. 
    145   They are usually neglected in the model 
    146   \footnote{ 
     148  It is usually specified as a fresh water flux at the air-sea interface in 
     149  the vicinity of river mouths. 
     150\item [Solid earth - ocean] Heat and salt fluxes through the sea floor are small, 
     151  except in special areas of little extent. 
     152  They are usually neglected in the model \footnote{ 
    147153    In fact, it has been shown that the heat flux associated with the solid Earth cooling 
    148     (\ie\ the geothermal heating) is not negligible for the thermohaline circulation of the world ocean 
    149     (see \autoref{subsec:TRA_bbc}). 
     154    (\ie\ the geothermal heating) is not negligible for 
     155    the thermohaline circulation of the world ocean (see \autoref{subsec:TRA_bbc}). 
    150156  }. 
    151157  The boundary condition is thus set to no flux of heat and salt across solid boundaries. 
    152   For momentum, the situation is different. There is no flow across solid boundaries, 
    153   \ie\ the velocity normal to the ocean bottom and coastlines is zero (in other words, 
    154   the bottom velocity is parallel to solid boundaries). This kinematic boundary condition 
    155   can be expressed as: 
     158  For momentum, the situation is different. 
     159  There is no flow across solid boundaries, 
     160  \ie\ the velocity normal to the ocean bottom and coastlines is zero 
     161  (in other words, the bottom velocity is parallel to solid boundaries). 
     162  This kinematic boundary condition can be expressed as: 
    156163  \begin{equation} 
    157164    \label{eq:MB_w_bbc} 
     
    160167  In addition, the ocean exchanges momentum with the earth through frictional processes. 
    161168  Such momentum transfer occurs at small scales in a boundary layer. 
    162   It must be parameterized in terms of turbulent fluxes using bottom and/or lateral boundary conditions. 
     169  It must be parameterized in terms of turbulent fluxes using 
     170  bottom and/or lateral boundary conditions. 
    163171  Its specification depends on the nature of the physical parameterisation used for 
    164172  $\vect D^{\vect U}$ in \autoref{eq:MB_PE_dyn}. 
    165   It is discussed in \autoref{eq:MB_zdf}.% and Chap. III.6 to 9. 
    166 \item [Atmosphere - ocean interface:]  the kinematic surface condition plus the mass flux of fresh water PE (the precipitation minus evaporation budget) 
    167   leads to: 
     173  It is discussed in \autoref{eq:MB_zdf}. % and Chap. III.6 to 9. 
     174\item [Atmosphere - ocean] The kinematic surface condition plus the mass flux of fresh water PE 
     175  (the precipitation minus evaporation budget) leads to: 
    168176  \[ 
    169177    % \label{eq:MB_w_sbc} 
    170178    w = \pd[\eta]{t} + \lt. \vect U_h \rt|_{z = \eta} \cdot \nabla_h (\eta) + P - E 
    171179  \] 
    172   The dynamic boundary condition, neglecting the surface tension (which removes capillary waves from the system) 
    173   leads to the continuity of pressure across the interface $z = \eta$. 
     180  The dynamic boundary condition, neglecting the surface tension 
     181  (which removes capillary waves from the system) leads to 
     182  the continuity of pressure across the interface $z = \eta$. 
    174183  The atmosphere and ocean also exchange horizontal momentum (wind stress), and heat. 
    175 \item [Sea ice - ocean interface:]  the ocean and sea ice exchange heat, salt, fresh water and momentum. 
     184\item [Sea ice - ocean] The ocean and sea ice exchange heat, salt, fresh water and momentum. 
    176185  The sea surface temperature is constrained to be at the freezing point at the interface. 
    177186  Sea ice salinity is very low ($\sim4-6 \, psu$) compared to those of the ocean ($\sim34 \, psu$). 
    178   The cycle of freezing/melting is associated with fresh water and salt fluxes that cannot be neglected. 
     187  The cycle of freezing/melting is associated with fresh water and salt fluxes that 
     188  cannot be neglected. 
    179189\end{description} 
    180190 
     
    198208\] 
    199209Two strategies can be considered for the surface pressure term: 
    200 $(a)$ introduce of a  new variable $\eta$, the free-surface elevation, 
     210\begin{enumerate*}[label={(\alph*)}] 
     211\item introduce of a new variable $\eta$, the free-surface elevation, 
    201212for which a prognostic equation can be established and solved; 
    202 $(b)$ assume that the ocean surface is a rigid lid, 
     213\item assume that the ocean surface is a rigid lid, 
    203214on which the pressure (or its horizontal gradient) can be diagnosed. 
     215\end{enumerate*} 
    204216When the former strategy is used, one solution of the free-surface elevation consists of 
    205217the excitation of external gravity waves. 
     
    212224modifies certain other longwave dynamics (\eg\ barotropic Rossby or planetary waves). 
    213225The rigid-lid hypothesis is an obsolescent feature in modern OGCMs. 
    214 It has been available until the release 3.1 of \NEMO, and it has been removed in release 3.2 and followings. 
     226It has been available until the release 3.1 of \NEMO, 
     227and it has been removed in release 3.2 and followings. 
    215228Only the free surface formulation is now described in this document (see the next sub-section). 
    216229 
     
    221234In the free surface formulation, a variable $\eta$, the sea-surface height, 
    222235is introduced which describes the shape of the air-sea interface. 
    223 This variable is solution of a prognostic equation which is established by forming the vertical average of 
    224 the kinematic surface condition (\autoref{eq:MB_w_bbc}): 
     236This variable is solution of a prognostic equation which is established by 
     237forming the vertical average of the kinematic surface condition (\autoref{eq:MB_w_bbc}): 
    225238\begin{equation} 
    226239  \label{eq:MB_ssh} 
    227240  \pd[\eta]{t} = - D + P - E \quad \text{where} \quad D = \nabla \cdot \lt[ (H + \eta) \; \overline{U}_h \, \rt] 
    228241\end{equation} 
    229 and using (\autoref{eq:MB_PE_hydrostatic}) the surface pressure is given by: $p_s = \rho \, g \, \eta$. 
    230  
    231 Allowing the air-sea interface to move introduces the external gravity waves (EGWs) as 
     242and using (\autoref{eq:MB_PE_hydrostatic}) the surface pressure is given by: 
     243$p_s = \rho \, g \, \eta$. 
     244 
     245Allowing the air-sea interface to move introduces 
     246the \textbf{E}xternal \textbf{G}ravity \textbf{W}aves (EGWs) as 
    232247a class of solution of the primitive equations. 
    233248These waves are barotropic (\ie\ nearly independent of depth) and their phase speed is quite high. 
     
    236251Two choices can be made regarding the implementation of the free surface in the model, 
    237252depending on the physical processes of interest. 
    238  
    239 $\bullet$ If one is interested in EGWs, in particular the tides and their interaction with 
    240 the baroclinic structure of the ocean (internal waves) possibly in shallow seas, 
    241 then a non linear free surface is the most appropriate. 
    242 This means that no approximation is made in \autoref{eq:MB_ssh} and that 
    243 the variation of the ocean volume is fully taken into account. 
    244 Note that in order to study the fast time scales associated with EGWs it is necessary to 
    245 minimize time filtering effects 
    246 (use an explicit time scheme with very small time step, or a split-explicit scheme with reasonably small time step, 
    247 see \autoref{subsec:DYN_spg_exp} or \autoref{subsec:DYN_spg_ts}). 
    248  
    249 $\bullet$ If one is not interested in EGW but rather sees them as high frequency noise, 
    250 it is possible to apply an explicit filter to slow down the fastest waves while 
    251 not altering the slow barotropic Rossby waves. 
    252 If further, an approximative conservation of heat and salt contents is sufficient for the problem solved, 
    253 then it is sufficient to solve a linearized version of \autoref{eq:MB_ssh}, 
    254 which still allows to take into account freshwater fluxes applied at the ocean surface \citep{roullet.madec_JGR00}. 
    255 Nevertheless, with the linearization, an exact conservation of heat and salt contents is lost. 
    256  
    257 The filtering of EGWs in models with a free surface is usually a matter of discretisation of 
    258 the temporal derivatives, 
     253\begin{itemize} 
     254\item If one is interested in EGWs, in particular the tides and their interaction with 
     255  the baroclinic structure of the ocean (internal waves) possibly in shallow seas, 
     256  then a non linear free surface is the most appropriate. 
     257  This means that no approximation is made in \autoref{eq:MB_ssh} and that 
     258  the variation of the ocean volume is fully taken into account. 
     259  Note that in order to study the fast time scales associated with EGWs it is necessary to 
     260  minimize time filtering effects 
     261  (use an explicit time scheme with very small time step, 
     262  or a split-explicit scheme with reasonably small time step, 
     263  see \autoref{subsec:DYN_spg_exp} or \autoref{subsec:DYN_spg_ts}). 
     264\item If one is not interested in EGWs but rather sees them as high frequency noise, 
     265  it is possible to apply an explicit filter to slow down the fastest waves while 
     266  not altering the slow barotropic Rossby waves. 
     267  If further, an approximative conservation of heat and salt contents is sufficient for 
     268  the problem solved, then it is sufficient to solve a linearized version of \autoref{eq:MB_ssh}, 
     269  which still allows to take into account freshwater fluxes applied at the ocean surface 
     270  \citep{roullet.madec_JGR00}. 
     271  Nevertheless, with the linearization, an exact conservation of heat and salt contents is lost. 
     272\end{itemize} 
     273The filtering of EGWs in models with a free surface is usually 
     274a matter of discretisation of the temporal derivatives, 
    259275using a split-explicit method \citep{killworth.webb.ea_JPO91, zhang.endoh_JGR92} or 
    260 the implicit scheme \citep{dukowicz.smith_JGR94} or 
    261 the addition of a filtering force in the momentum equation \citep{roullet.madec_JGR00}. 
    262 With the present release, \NEMO\  offers the choice between 
    263 an explicit free surface (see \autoref{subsec:DYN_spg_exp}) or 
    264 a split-explicit scheme strongly inspired the one proposed by \citet{shchepetkin.mcwilliams_OM05} 
    265 (see \autoref{subsec:DYN_spg_ts}). 
    266  
    267 %% ================================================================================================= 
    268 \section{Curvilinear \textit{z-}coordinate system} 
     276the implicit scheme \citep{dukowicz.smith_JGR94} or the addition of a filtering force in 
     277the momentum equation \citep{roullet.madec_JGR00}. 
     278With the present release, \NEMO\ offers the choice between an explicit free surface 
     279(see \autoref{subsec:DYN_spg_exp}) or a split-explicit scheme strongly inspired the one proposed by 
     280\citet{shchepetkin.mcwilliams_OM05} (see \autoref{subsec:DYN_spg_ts}). 
     281 
     282%% ================================================================================================= 
     283\section{Curvilinear \textit{z}-coordinate system} 
    269284\label{sec:MB_zco} 
    270285 
     
    283298A solution consists of introducing an appropriate coordinate transformation that 
    284299shifts the singular point onto land \citep{madec.imbard_CD96, murray_JCP96}. 
    285 As a consequence, it is important to solve the primitive equations in various curvilinear coordinate systems. 
    286 An efficient way of introducing an appropriate coordinate transform can be found when using a tensorial formalism. 
     300As a consequence, 
     301it is important to solve the primitive equations in various curvilinear coordinate systems. 
     302An efficient way of introducing an appropriate coordinate transform can be found when 
     303using a tensorial formalism. 
    287304This formalism is suited to any multidimensional curvilinear coordinate system. 
    288 Ocean modellers mainly use three-dimensional orthogonal grids on the sphere (spherical earth approximation), 
    289 with preservation of the local vertical. Here we give the simplified equations for this particular case. 
    290 The general case is detailed by \citet{eiseman.stone_SR80} in their survey of the conservation laws of fluid dynamics. 
     305Ocean modellers mainly use three-dimensional orthogonal grids on the sphere 
     306(spherical earth approximation), with preservation of the local vertical. 
     307Here we give the simplified equations for this particular case. 
     308The general case is detailed by \citet{eiseman.stone_SR80} in 
     309their survey of the conservation laws of fluid dynamics. 
    291310 
    292311Let $(i,j,k)$ be a set of orthogonal curvilinear coordinates on 
     
    303322\begin{equation} 
    304323  \label{eq:MB_scale_factors} 
    305   \begin{aligned} 
    306     e_1 &= (a + z) \lt[ \lt( \pd[\lambda]{i} \cos \varphi \rt)^2 + \lt( \pd[\varphi]{i} \rt)^2 \rt]^{1/2} \\ 
    307     e_2 &= (a + z) \lt[ \lt( \pd[\lambda]{j} \cos \varphi \rt)^2 + \lt( \pd[\varphi]{j} \rt)^2 \rt]^{1/2} \\ 
    308     e_3 &= \lt( \pd[z]{k} \rt) 
    309   \end{aligned} 
     324    e_1 = (a + z) \lt[ \lt( \pd[\lambda]{i} \cos \varphi \rt)^2 + \lt( \pd[\varphi]{i} \rt)^2 \rt]^{1/2} \quad e_2 = (a + z) \lt[ \lt( \pd[\lambda]{j} \cos \varphi \rt)^2 + \lt( \pd[\varphi]{j} \rt)^2 \rt]^{1/2} \quad e_3 = \lt( \pd[z]{k} \rt) 
    310325\end{equation} 
    311326 
    312 \begin{figure}[!tb] 
     327\begin{figure} 
    313328  \centering 
    314   \includegraphics[width=0.66\textwidth]{Fig_I_earth_referential} 
     329  \includegraphics[width=0.33\textwidth]{Fig_I_earth_referential} 
    315330  \caption[Geographical and curvilinear coordinate systems]{ 
    316331    the geographical coordinate system $(\lambda,\varphi,z)$ and the curvilinear 
     
    328343\begin{subequations} 
    329344  % \label{eq:MB_discrete_operators} 
    330   \begin{gather} 
     345  \begin{align} 
    331346    \label{eq:MB_grad} 
    332     \nabla q =   \frac{1}{e_1} \pd[q]{i} \; \vect i 
    333                + \frac{1}{e_2} \pd[q]{j} \; \vect j 
    334                + \frac{1}{e_3} \pd[q]{k} \; \vect k \\ 
     347    \nabla q &=   \frac{1}{e_1} \pd[q]{i} \; \vect i + \frac{1}{e_2} \pd[q]{j} \; \vect j + \frac{1}{e_3} \pd[q]{k} \; \vect k \\ 
    335348    \label{eq:MB_div} 
    336     \nabla \cdot \vect A =   \frac{1}{e_1 \; e_2} \lt[ \pd[(e_2 \; a_1)]{\partial i} + \pd[(e_1 \; a_2)]{j} \rt] 
    337                            + \frac{1}{e_3} \lt[ \pd[a_3]{k} \rt] 
    338   \end{gather} 
    339   \begin{multline} 
     349    \nabla \cdot \vect A &=   \frac{1}{e_1 \; e_2} \lt[ \pd[(e_2 \; a_1)]{\partial i} + \pd[(e_1 \; a_2)]{j} \rt] + \frac{1}{e_3} \lt[ \pd[a_3]{k} \rt] \\ 
    340350    \label{eq:MB_curl} 
    341       \nabla \times \vect{A} =   \lt[ \frac{1}{e_2} \pd[a_3]{j} - \frac{1}{e_3} \pd[a_2]{k}   \rt] \vect i \\ 
    342                                + \lt[ \frac{1}{e_3} \pd[a_1]{k} - \frac{1}{e_1} \pd[a_3]{i}   \rt] \vect j \\ 
    343                                + \frac{1}{e_1 e_2} \lt[ \pd[(e_2 a_2)]{i} - \pd[(e_1 a_1)]{j} \rt] \vect k 
    344   \end{multline} 
    345   \begin{gather} 
     351      \nabla \times \vect{A} &=   \lt[ \frac{1}{e_2} \pd[a_3]{j} - \frac{1}{e_3} \pd[a_2]{k}   \rt] \vect i + \lt[ \frac{1}{e_3} \pd[a_1]{k} - \frac{1}{e_1} \pd[a_3]{i}   \rt] \vect j + \frac{1}{e_1 e_2} \lt[ \pd[(e_2 a_2)]{i} - \pd[(e_1 a_1)]{j} \rt] \vect k \\ 
    346352    \label{eq:MB_lap} 
    347     \Delta q = \nabla \cdot (\nabla q) \\ 
     353    \Delta q &= \nabla \cdot (\nabla q) \\ 
    348354    \label{eq:MB_lap_vector} 
    349     \Delta \vect A = \nabla (\nabla \cdot \vect A) - \nabla \times (\nabla \times \vect A) 
    350   \end{gather} 
     355    \Delta \vect A &= \nabla (\nabla \cdot \vect A) - \nabla \times (\nabla \times \vect A) 
     356  \end{align} 
    351357\end{subequations} 
    352 where $q$ is a scalar quantity and $\vect A = (a_1,a_2,a_3)$ a vector in the $(i,j,k)$ coordinates system. 
     358where $q$ is a scalar quantity and 
     359$\vect A = (a_1,a_2,a_3)$ a vector in the $(i,j,k)$ coordinates system. 
    353360 
    354361%% ================================================================================================= 
     
    357364 
    358365In order to express the Primitive Equations in tensorial formalism, 
    359 it is necessary to compute the horizontal component of the non-linear and viscous terms of the equation using 
    360 \autoref{eq:MB_grad}) to \autoref{eq:MB_lap_vector}. 
    361 Let us set $\vect U = (u,v,w) = \vect U_h + w \; \vect k $, the velocity in the $(i,j,k)$ coordinates system and 
    362 define the relative vorticity $\zeta$ and the divergence of the horizontal velocity field $\chi$, by: 
     366it is necessary to compute the horizontal component of the non-linear and viscous terms of 
     367the equation using \autoref{eq:MB_grad}) to \autoref{eq:MB_lap_vector}. 
     368Let us set $\vect U = (u,v,w) = \vect U_h + w \; \vect k $, 
     369the velocity in the $(i,j,k)$ coordinates system, 
     370and define the relative vorticity $\zeta$ and the divergence of the horizontal velocity field $\chi$, 
     371by: 
    363372\begin{gather} 
    364373  \label{eq:MB_curl_Uh} 
     
    371380$e_3$  is a function of the single variable $k$, 
    372381$NLT$ the nonlinear term of \autoref{eq:MB_PE_dyn} can be transformed as follows: 
    373 \begin{alignat*}{2} 
    374   &NLT &=   &\lt[ (\nabla \times {\vect U}) \times {\vect U} + \frac{1}{2} \nabla \lt( {\vect U}^2 \rt) \rt]_h \\ 
    375   &    &=   &\lt( 
    376     \begin{array}{*{20}c} 
    377                 \lt[ \frac{1}{e_3} \pd[u]{k} - \frac{1}{e_1} \pd[w]{i} \rt] w - \zeta \; v   \\ 
    378                 \zeta \; u - \lt[ \frac{1}{e_2} \pd[w]{j} - \frac{1}{e_3} \pd[v]{k} \rt] \ w 
    379     \end{array} 
    380                                                                                              \rt) 
    381           + \frac{1}{2} \lt( 
    382     \begin{array}{*{20}c} 
    383                              \frac{1}{e_1} \pd[(u^2 + v^2 + w^2)]{i} \\ 
    384                              \frac{1}{e_2} \pd[(u^2 + v^2 + w^2)]{j} 
    385     \end{array} 
    386                                                                      \rt) \\ 
    387   &    &=   &\lt( 
    388     \begin{array}{*{20}c} 
    389                   -\zeta \; v \\ 
    390                    \zeta \; u 
    391     \end{array} 
    392                               \rt) 
    393           + \frac{1}{2} \lt( 
    394     \begin{array}{*{20}c} 
    395                              \frac{1}{e_1} \pd[(u^2 + v^2)]{i} \\ 
    396                              \frac{1}{e_2} \pd[(u^2 + v^2)]{j} 
    397     \end{array} 
    398                                                                \rt) \\ 
    399   &    &  &+ \frac{1}{e_3} \lt( 
    400     \begin{array}{*{20}c} 
    401                                 w \; \pd[u]{k} \\ 
    402                                 w \; \pd[v]{k} 
    403     \end{array} 
    404                                                \rt) 
    405            - \lt( 
    406     \begin{array}{*{20}c} 
    407                   \frac{w}{e_1} \pd[w]{i} - \frac{1}{2 e_1} \pd[w^2]{i} \\ 
    408                   \frac{w}{e_2} \pd[w]{j} - \frac{1}{2 e_2} \pd[w^2]{j} 
    409     \end{array} 
    410                                                                         \rt) 
    411 \end{alignat*} 
    412 The last term of the right hand side is obviously zero, and thus the nonlinear term of 
    413 \autoref{eq:MB_PE_dyn} is written in the $(i,j,k)$ coordinate system: 
     382\begin{align*} 
     383  NLT &= \lt[ (\nabla \times {\vect U}) \times {\vect U} + \frac{1}{2} \nabla \lt( {\vect U}^2 \rt) \rt]_h \\ 
     384      &= \lt( 
     385        \begin{array}{*{20}c} 
     386          \lt[ \frac{1}{e_3} \pd[u]{k} - \frac{1}{e_1} \pd[w]{i} \rt] w - \zeta \; v   \\ 
     387          \zeta \; u - \lt[ \frac{1}{e_2} \pd[w]{j} - \frac{1}{e_3} \pd[v]{k} \rt] \ w 
     388        \end{array} 
     389  \rt) 
     390  + \frac{1}{2} \lt( 
     391  \begin{array}{*{20}c} 
     392    \frac{1}{e_1} \pd[(u^2 + v^2 + w^2)]{i} \\ 
     393    \frac{1}{e_2} \pd[(u^2 + v^2 + w^2)]{j} 
     394  \end{array} 
     395  \rt) \\ 
     396      &= \lt( 
     397        \begin{array}{*{20}c} 
     398          -\zeta \; v \\ 
     399          \zeta \; u 
     400        \end{array} 
     401  \rt) 
     402  + \frac{1}{2} \lt( 
     403  \begin{array}{*{20}c} 
     404    \frac{1}{e_1} \pd[(u^2 + v^2)]{i} \\ 
     405    \frac{1}{e_2} \pd[(u^2 + v^2)]{j} 
     406  \end{array} 
     407  \rt) 
     408  + \frac{1}{e_3} \lt( 
     409  \begin{array}{*{20}c} 
     410    w \; \pd[u]{k} \\ 
     411    w \; \pd[v]{k} 
     412  \end{array} 
     413  \rt) 
     414  - \lt( 
     415  \begin{array}{*{20}c} 
     416    \frac{w}{e_1} \pd[w]{i} - \frac{1}{2 e_1} \pd[w^2]{i} \\ 
     417    \frac{w}{e_2} \pd[w]{j} - \frac{1}{2 e_2} \pd[w^2]{j} 
     418  \end{array} 
     419  \rt) 
     420\end{align*} 
     421The last term of the right hand side is obviously zero, 
     422and thus the \textbf{N}on\textbf{L}inear \textbf{T}erm ($NLT$) of \autoref{eq:MB_PE_dyn} is written in 
     423the $(i,j,k)$ coordinate system: 
    414424\begin{equation} 
    415425  \label{eq:MB_vector_form} 
    416   NLT =   \zeta \; \vect k \times \vect U_h + \frac{1}{2} \nabla_h \lt( \vect U_h^2 \rt) 
    417         + \frac{1}{e_3} w \pd[\vect U_h]{k} 
     426  NLT = \zeta \; \vect k \times \vect U_h + \frac{1}{2} \nabla_h \lt( \vect U_h^2 \rt) + \frac{1}{e_3} w \pd[\vect U_h]{k} 
    418427\end{equation} 
    419428 
     
    421430For some purposes, it can be advantageous to write this term in the so-called flux form, 
    422431\ie\ to write it as the divergence of fluxes. 
    423 For example, the first component of \autoref{eq:MB_vector_form} (the $i$-component) is transformed as follows: 
    424 \begin{alignat*}{2} 
     432For example, 
     433the first component of \autoref{eq:MB_vector_form} (the $i$-component) is transformed as follows: 
     434\begin{alignat*}{3} 
    425435  &NLT_i &= &- \zeta \; v + \frac{1}{2 \; e_1} \pd[ (u^2 + v^2) ]{i} + \frac{1}{e_3} w \ \pd[u]{k} \\ 
    426   &      &=  &\frac{1}{e_1 \; e_2} \lt( -v \pd[(e_2 \, v)]{i} + v \pd[(e_1 \, u)]{j} \rt) 
    427             + \frac{1}{e_1 e_2} \lt( e_2 \; u \pd[u]{i} + e_2 \; v \pd[v]{i} \rt) \\ 
    428   &      & &+ \frac{1}{e_3} \lt( w \; \pd[u]{k} \rt) \\ 
    429   &      &=  &\frac{1}{e_1 \; e_2} \lt[ - \lt( v^2 \pd[e_2]{i} + e_2 \, v \pd[v]{i} \rt) 
    430                                      + \lt( \pd[ \lt( e_1 \, u \, v \rt)]{j} -         e_1 \, u \pd[v]{j} \rt) \rt. \\ 
    431   &      &                       &\lt. + \lt( \pd[ \lt( e_2 \, u \, u \rt)]{i} - u \pd[ \lt( e_2 u \rt)]{i} \rt) 
    432                                      + e_2 v \pd[v]{i}                                                         \rt] \\ 
     436  &      &= &\frac{1}{e_1 \; e_2} \lt( -v \pd[(e_2 \, v)]{i} + v \pd[(e_1 \, u)]{j} \rt) + \frac{1}{e_1 e_2} \lt( e_2 \; u \pd[u]{i} + e_2 \; v \pd[v]{i} \rt) + \frac{1}{e_3} \lt( w \; \pd[u]{k} \rt) \\ 
     437  &      &= &\frac{1}{e_1 \; e_2} \lt[ - \lt( v^2 \pd[e_2]{i} + e_2 \, v \pd[v]{i} \rt) + \lt( \pd[ \lt( e_1 \, u \, v \rt)]{j} - e_1 \, u \pd[v]{j} \rt) \rt. \lt. + \lt( \pd[ \lt( e_2 \, u \, u \rt)]{i} - u \pd[ \lt( e_2 u \rt)]{i} \rt) + e_2 v \pd[v]{i} \rt] \\ 
    433438  &      & &+ \frac{1}{e_3} \lt( \pd[(w \, u)]{k} - u \pd[w]{k} \rt) \\ 
    434   &      &=  &\frac{1}{e_1 \; e_2} \lt( \pd[(e_2 \, u \, u)]{i} + \pd[(e_1 \, u \, v)]{j} \rt) 
    435             + \frac{1}{e_3} \pd[(w \, u)]{k} \\ 
    436   &      & &+ \frac{1}{e_1 e_2} \lt[ - u \lt( \pd[(e_1 v)]{j} - v \, \pd[e_1]{j} \rt) 
    437                                   - u \pd[(e_2 u)]{i}                              \rt] 
    438             - \frac{1}{e_3} \pd[w]{k} u \\ 
     439  &      &= &\frac{1}{e_1 \; e_2} \lt( \pd[(e_2 \, u \, u)]{i} + \pd[(e_1 \, u \, v)]{j} \rt) + \frac{1}{e_3} \pd[(w \, u)]{k} + \frac{1}{e_1 e_2} \lt[ - u \lt( \pd[(e_1 v)]{j} - v \, \pd[e_1]{j} \rt) - u \pd[(e_2 u)]{i} \rt] - \frac{1}{e_3} \pd[w]{k} u \\ 
    439440  &      & &+ \frac{1}{e_1 e_2} \lt( - v^2 \pd[e_2]{i} \rt) \\ 
    440   &      &= &\nabla \cdot (\vect U \, u) - (\nabla \cdot \vect U) \ u 
    441             + \frac{1}{e_1 e_2} \lt( -v^2 \pd[e_2]{i} + u v \, \pd[e_1]{j} \rt) \\ 
    442   \intertext{as $\nabla \cdot {\vect U} \; = 0$ (incompressibility) it becomes:} 
     441  &      &= &\nabla \cdot (\vect U \, u) - (\nabla \cdot \vect U) \ u + \frac{1}{e_1 e_2} \lt( -v^2 \pd[e_2]{i} + u v \, \pd[e_1]{j} \rt) \\ 
     442  \shortintertext{as $\nabla \cdot {\vect U} \; = 0$ (incompressibility) it becomes:} 
    443443  &      &= &\, \nabla \cdot (\vect U \, u) + \frac{1}{e_1 e_2} \lt( v \; \pd[e_2]{i} - u \; \pd[e_1]{j} \rt) (-v) 
    444444\end{alignat*} 
     
    447447\begin{equation} 
    448448  \label{eq:MB_flux_form} 
    449   NLT =   \nabla \cdot \lt( 
    450     \begin{array}{*{20}c} 
    451                             \vect U \, u \\ 
    452                             \vect U \, v 
    453     \end{array} 
    454                                          \rt) 
    455         + \frac{1}{e_1 e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \vect k \times \vect U_h 
     449  NLT = \nabla \cdot \lt( 
     450  \begin{array}{*{20}c} 
     451    \vect U \, u \\ 
     452    \vect U \, v 
     453  \end{array} 
     454  \rt) 
     455  + \frac{1}{e_1 e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \vect k \times \vect U_h 
    456456\end{equation} 
    457457 
    458458The flux form has two terms, 
    459 the first one is expressed as the divergence of momentum fluxes (hence the flux form name given to this formulation) 
    460 and the second one is due to the curvilinear nature of the coordinate system used. 
    461 The latter is called the \textit{metric} term and can be viewed as a modification of the Coriolis parameter: 
     459the first one is expressed as the divergence of momentum fluxes 
     460(hence the flux form name given to this formulation) and 
     461the second one is due to the curvilinear nature of the coordinate system used. 
     462The latter is called the \textit{metric} term and can be viewed as 
     463a modification of the Coriolis parameter: 
    462464\[ 
    463465  % \label{eq:MB_cor+metric} 
     
    472474the following tensorial formalism: 
    473475 
    474 \begin{itemize} 
    475 \item \textbf{Vector invariant form of the momentum equations}: 
     476\begin{description} 
     477\item [Vector invariant form of the momentum equations] 
    476478  \begin{equation} 
    477479    \label{eq:MB_dyn_vect} 
    478     \begin{split} 
     480    \begin{gathered} 
    479481    % \label{eq:MB_dyn_vect_u} 
    480       \pd[u]{t} = &+ (\zeta + f) \, v - \frac{1}{2 e_1} \pd[]{i} (u^2 + v^2) 
    481                    - \frac{1}{e_3} w \pd[u]{k} - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) \\ 
    482                   &+ D_u^{\vect U} + F_u^{\vect U} \\ 
    483       \pd[v]{t} = &- (\zeta + f) \, u - \frac{1}{2 e_2} \pd[]{j} (u^2 + v^2) 
    484                    - \frac{1}{e_3} w \pd[v]{k} - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) \\ 
    485                   &+ D_v^{\vect U} + F_v^{\vect U} 
    486     \end{split} 
     482      \pd[u]{t} = + (\zeta + f) \, v - \frac{1}{2 e_1} \pd[]{i} (u^2 + v^2) - \frac{1}{e_3} w \pd[u]{k} - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) + D_u^{\vect U} + F_u^{\vect U} \\ 
     483      \pd[v]{t} = - (\zeta + f) \, u - \frac{1}{2 e_2} \pd[]{j} (u^2 + v^2) - \frac{1}{e_3} w \pd[v]{k} - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) + D_v^{\vect U} + F_v^{\vect U} 
     484    \end{gathered} 
    487485  \end{equation} 
    488 \item \textbf{flux form of the momentum equations}: 
     486\item [Flux form of the momentum equations] 
    489487  % \label{eq:MB_dyn_flux} 
    490488  \begin{multline*} 
    491489    % \label{eq:MB_dyn_flux_u} 
    492     \pd[u]{t} = + \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, v \\ 
    493                 - \frac{1}{e_1 \; e_2} \lt( \pd[(e_2 \, u \, u)]{i} + \pd[(e_1 \, v \, u)]{j} \rt) \\ 
    494                 - \frac{1}{e_3} \pd[(w \, u)]{k} - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) 
    495                 + D_u^{\vect U} + F_u^{\vect U} 
     490    \pd[u]{t} = + \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, v - \frac{1}{e_1 \; e_2} \lt( \pd[(e_2 \, u \, u)]{i} + \pd[(e_1 \, v \, u)]{j} \rt) \\ 
     491    - \frac{1}{e_3} \pd[(w \, u)]{k} - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) + D_u^{\vect U} + F_u^{\vect U} 
    496492  \end{multline*} 
    497493  \begin{multline*} 
    498494    % \label{eq:MB_dyn_flux_v} 
    499     \pd[v]{t} = - \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, u \\ 
    500                 - \frac{1}{e_1 \; e_2} \lt( \pd[(e_2 \, u \, v)]{i} + \pd[(e_1 \, v \, v)]{j} \rt) \\ 
    501                 - \frac{1}{e_3} \pd[(w \, v)]{k} - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) 
    502                 + D_v^{\vect U} + F_v^{\vect U} 
     495    \pd[v]{t} = - \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, u - \frac{1}{e_1 \; e_2} \lt( \pd[(e_2 \, u \, v)]{i} + \pd[(e_1 \, v \, v)]{j} \rt) \\ 
     496    - \frac{1}{e_3} \pd[(w \, v)]{k} - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) + D_v^{\vect U} + F_v^{\vect U} 
    503497  \end{multline*} 
    504   where $\zeta$, the relative vorticity, is given by \autoref{eq:MB_curl_Uh} and $p_s$, the surface pressure, 
    505   is given by: 
     498  where $\zeta$, the relative vorticity, is given by \autoref{eq:MB_curl_Uh} and 
     499  $p_s$, the surface pressure, is given by: 
    506500  \[ 
    507501  % \label{eq:MB_spg} 
     
    518512  \] 
    519513  where the divergence of the horizontal velocity, $\chi$ is given by \autoref{eq:MB_div_Uh}. 
    520  
    521 \item \textbf{tracer equations}: 
    522   \begin{equation} 
    523   \begin{split} 
    524     \pd[T]{t} = & - \frac{1}{e_1 e_2} \lt[ \pd[(e_2 T \, u)]{i} + \pd[(e_1 T \, v)]{j} \rt] 
    525                 - \frac{1}{e_3} \pd[(T \, w)]{k} + D^T + F^T \\ 
    526     \pd[S]{t} = & - \frac{1}{e_1 e_2} \lt[ \pd[(e_2 S \, u)]{i} + \pd[(e_1 S \, v)]{j} \rt] 
    527                 - \frac{1}{e_3} \pd[(S \, w)]{k} + D^S + F^S \\ 
    528     \rho = & \rho \big( T,S,z(k) \big) 
    529   \end{split} 
    530   \end{equation} 
    531 \end{itemize} 
    532  
    533 The expression of $\vect D^{U}$, $D^{S}$ and $D^{T}$ depends on the subgrid scale parameterisation used. 
     514\item [Tracer equations] 
     515  \begin{gather*} 
     516    \pd[T]{t} = - \frac{1}{e_1 e_2} \lt[ \pd[(e_2 T \, u)]{i} + \pd[(e_1 T \, v)]{j} \rt] - \frac{1}{e_3} \pd[(T \, w)]{k} + D^T + F^T \\ 
     517    \pd[S]{t} = - \frac{1}{e_1 e_2} \lt[ \pd[(e_2 S \, u)]{i} + \pd[(e_1 S \, v)]{j} \rt] - \frac{1}{e_3} \pd[(S \, w)]{k} + D^S + F^S \\ 
     518    \rho = \rho \big( T,S,z(k) \big) 
     519  \end{gather*} 
     520\end{description} 
     521 
     522The expression of $\vect D^{U}$, $D^{S}$ and $D^{T}$ depends on 
     523the subgrid scale parameterisation used. 
    534524It will be defined in \autoref{eq:MB_zdf}. 
    535525The nature and formulation of $\vect F^{\vect U}$, $F^T$ and $F^S$, the surface forcing terms, 
     
    545535varying from more than 6,000 meters in abyssal trenches to zero at the coast. 
    546536Last but not least, the ocean stratification exerts a strong barrier to vertical motions and mixing. 
    547 Therefore, in order to represent the ocean with respect to 
    548 the first point a space and time dependent vertical coordinate that follows the variation of the sea surface height 
     537Therefore, in order to represent the ocean with respect to the first point 
     538a space and time dependent vertical coordinate that follows the variation of the sea surface height 
    549539\eg\ an \zstar-coordinate; 
    550 for the second point, a space variation to fit the change of bottom topography 
     540for the second point, 
     541a space variation to fit the change of bottom topography 
    551542\eg\ a terrain-following or $\sigma$-coordinate; 
    552 and for the third point, one will be tempted to use a space and time dependent coordinate that 
    553 follows the isopycnal surfaces, \eg\ an isopycnic coordinate. 
    554  
    555 In order to satisfy two or more constraints one can even be tempted to mixed these coordinate systems, as in 
    556 HYCOM (mixture of $z$-coordinate at the surface, isopycnic coordinate in the ocean interior and $\sigma$ at 
    557 the ocean bottom) \citep{chassignet.smith.ea_JPO03} or 
    558 OPA (mixture of $z$-coordinate in vicinity the surface and steep topography areas and $\sigma$-coordinate elsewhere) 
    559 \citep{madec.delecluse.ea_JPO96} among others. 
     543and for the third point, 
     544one will be tempted to use a space and time dependent coordinate that follows the isopycnal surfaces, 
     545\eg\ an isopycnic coordinate. 
     546 
     547In order to satisfy two or more constraints one can even be tempted to mixed these coordinate systems, 
     548as in HYCOM (mixture of $z$-coordinate at the surface, isopycnic coordinate in the ocean interior and 
     549$\sigma$ at the ocean bottom) \citep{chassignet.smith.ea_JPO03} or 
     550OPA (mixture of $z$-coordinate in vicinity the surface and steep topography areas and 
     551$\sigma$-coordinate elsewhere) \citep{madec.delecluse.ea_JPO96} among others. 
    560552 
    561553In fact one is totally free to choose any space and time vertical coordinate by 
     
    565557  s = s(i,j,k,t) 
    566558\end{equation} 
    567 with the restriction that the above equation gives a single-valued monotonic relationship between $s$ and $k$, 
    568 when $i$, $j$ and $t$ are held fixed. 
    569 \autoref{eq:MB_s} is a transformation from the $(i,j,k,t)$ coordinate system with independent variables into 
     559with the restriction that the above equation gives a single-valued monotonic relationship between 
     560$s$ and $k$, when $i$, $j$ and $t$ are held fixed. 
     561\autoref{eq:MB_s} is a transformation from 
     562the $(i,j,k,t)$ coordinate system with independent variables into 
    570563the $(i,j,s,t)$ generalised coordinate system with $s$ depending on the other three variables through 
    571564\autoref{eq:MB_s}. 
    572565This so-called \textit{generalised vertical coordinate} \citep{kasahara_MWR74} is in fact 
    573 an Arbitrary Lagrangian--Eulerian (ALE) coordinate. 
    574 Indeed, one has a great deal of freedom in the choice of expression for $s$. The choice determines 
    575 which part of the vertical velocity (defined from a fixed referential) will cross the levels (Eulerian part) and 
    576 which part will be used to move them (Lagrangian part). 
    577 The coordinate is also sometime referenced as an adaptive coordinate \citep{hofmeister.burchard.ea_OM10}, 
    578 since the coordinate system is adapted in the course of the simulation. 
     566an \textbf{A}rbitrary \textbf{L}agrangian--\textbf{E}ulerian (ALE) coordinate. 
     567Indeed, one has a great deal of freedom in the choice of expression for $s$. 
     568The choice determines which part of the vertical velocity (defined from a fixed referential) 
     569will cross the levels (Eulerian part) and which part will be used to move them (Lagrangian part). 
     570The coordinate is also sometimes referenced as an adaptive coordinate 
     571\citep{hofmeister.burchard.ea_OM10}, since the coordinate system is adapted in 
     572the course of the simulation. 
    579573Its most often used implementation is via an ALE algorithm, 
    580574in which a pure lagrangian step is followed by regridding and remapping steps, 
     
    593587The horizontal velocity in ocean models measures motions in the horizontal plane, 
    594588perpendicular to the local gravitational field. 
    595 That is, horizontal velocity is mathematically the same regardless of the vertical coordinate, be it geopotential, 
    596 isopycnal, pressure, or terrain following. 
     589That is, horizontal velocity is mathematically the same regardless of the vertical coordinate, 
     590be it geopotential, isopycnal, pressure, or terrain following. 
    597591The key motivation for maintaining the same horizontal velocity component is that 
    598592the hydrostatic and geostrophic balances are dominant in the large-scale ocean. 
    599 Use of an alternative quasi -horizontal velocity, for example one oriented parallel to the generalized surface, 
     593Use of an alternative quasi -horizontal velocity, 
     594for example one oriented parallel to the generalized surface, 
    600595would lead to unacceptable numerical errors. 
    601596Correspondingly, the vertical direction is anti -parallel to the gravitational force in 
     
    604599the surface of a constant generalized vertical coordinate. 
    605600 
    606 It is the method used to measure transport across the generalized vertical coordinate surfaces which differs between 
    607 the vertical coordinate choices. 
    608 That is, computation of the dia-surface velocity component represents the fundamental distinction between 
     601It is the method used to measure transport across the generalized vertical coordinate surfaces which 
     602differs between the vertical coordinate choices. 
     603That is, 
     604computation of the dia-surface velocity component represents the fundamental distinction between 
    609605the various coordinates. 
    610 In some models, such as geopotential, pressure, and terrain following, this transport is typically diagnosed from 
    611 volume or mass conservation. 
    612 In other models, such as isopycnal layered models, this transport is prescribed based on assumptions about 
    613 the physical processes producing a flux across the layer interfaces. 
     606In some models, such as geopotential, pressure, and terrain following, 
     607this transport is typically diagnosed from volume or mass conservation. 
     608In other models, such as isopycnal layered models, 
     609this transport is prescribed based on assumptions about the physical processes producing 
     610a flux across the layer interfaces. 
    614611 
    615612In this section we first establish the PE in the generalised vertical $s$-coordinate, 
     
    620617\subsection{\textit{S}-coordinate formulation} 
    621618 
    622 Starting from the set of equations established in \autoref{sec:MB_zco} for the special case $k = z$ and 
    623 thus $e_3 = 1$, we introduce an arbitrary vertical coordinate $s = s(i,j,k,t)$, 
     619Starting from the set of equations established in \autoref{sec:MB_zco} for 
     620the special case $k = z$ and thus $e_3 = 1$, 
     621we introduce an arbitrary vertical coordinate $s = s(i,j,k,t)$, 
    624622which includes $z$-, \zstar- and $\sigma$-coordinates as special cases 
    625623($s = z$, $s = \zstar$, and $s = \sigma = z / H$ or $ = z / \lt( H + \eta \rt)$, resp.). 
    626624A formal derivation of the transformed equations is given in \autoref{apdx:SCOORD}. 
    627 Let us define the vertical scale factor by $e_3 = \partial_s z$  ($e_3$ is now a function of $(i,j,k,t)$ ), 
     625Let us define the vertical scale factor by $e_3 = \partial_s z$ 
     626($e_3$ is now a function of $(i,j,k,t)$ ), 
    628627and the slopes in the $(i,j)$ directions between $s$- and $z$-surfaces by: 
    629628\begin{equation} 
     
    632631  \sigma_2 = \frac{1}{e_2} \; \lt. \pd[z]{j} \rt|_s 
    633632\end{equation} 
    634 We also introduce $\omega$, a dia-surface velocity component, defined as the velocity 
    635 relative to the moving $s$-surfaces and normal to them: 
     633We also introduce $\omega$, a dia-surface velocity component, 
     634defined as the velocity relative to the moving $s$-surfaces and normal to them: 
    636635\[ 
    637636  % \label{eq:MB_sco_w} 
     
    642641(see \autoref{sec:SCOORD_momentum}): 
    643642 
    644 \begin{itemize} 
    645 \item \textbf{Vector invariant form of the momentum equation}: 
     643\begin{description} 
     644\item [Vector invariant form of the momentum equation] 
     645  \begin{gather*} 
     646    % \label{eq:MB_sco_u_vector} 
     647    \pd[u]{t} = + (\zeta + f) \, v - \frac{1}{2 \, e_1} \pd[]{i} (u^2 + v^2) - \frac{1}{e_3} \omega \pd[u]{k} - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) - g \frac{\rho}{\rho_o} \sigma_1 + D_u^{\vect U} + F_u^{\vect U} \\ 
     648    % \label{eq:MB_sco_v_vector} 
     649    \pd[v]{t} = - (\zeta + f) \, u - \frac{1}{2 \, e_2} \pd[]{j} (u^2 + v^2) - \frac{1}{e_3} \omega \pd[v]{k} - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) - g \frac{\rho}{\rho_o} \sigma_2 + D_v^{\vect U} + F_v^{\vect U} 
     650  \end{gather*} 
     651\item [Flux form of the momentum equation] 
    646652  \begin{multline*} 
    647   % \label{eq:MB_sco_u_vector} 
    648     \pd[u]{t} = + (\zeta + f) \, v - \frac{1}{2 \, e_1} \pd[]{i} (u^2 + v^2) - \frac{1}{e_3} \omega \pd[u]{k} \\ 
    649                 - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) - g \frac{\rho}{\rho_o} \sigma_1 
    650                 + D_u^{\vect U} + F_u^{\vect U} 
    651   \end{multline*} 
    652   \begin{multline*} 
    653   % \label{eq:MB_sco_v_vector} 
    654     \pd[v]{t} = - (\zeta + f) \, u - \frac{1}{2 \, e_2} \pd[]{j}(u^2 + v^2) - \frac{1}{e_3} \omega \pd[v]{k} \\ 
    655                 - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) - g \frac{\rho}{\rho_o} \sigma_2 
    656                 + D_v^{\vect U} + F_v^{\vect U} 
    657   \end{multline*} 
    658 \item \textbf{Flux form of the momentum equation}: 
    659   \begin{multline*} 
    660   % \label{eq:MB_sco_u_flux} 
    661     \frac{1}{e_3} \pd[(e_3 \, u)]{t} = + \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, v \\ 
    662                                        - \frac{1}{e_1 \; e_2 \; e_3} \lt( \pd[(e_2 \, e_3 \, u \, u)]{i} + \pd[(e_1 \, e_3 \, v \, u)]{j} \rt) \\ 
    663                                        - \frac{1}{e_3} \pd[(\omega \, u)]{k} 
    664                                        - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) 
    665                                        - g \frac{\rho}{\rho_o} \sigma_1 + D_u^{\vect U} + F_u^{\vect U} 
     653    % \label{eq:MB_sco_u_flux} 
     654    \frac{1}{e_3} \pd[(e_3 \, u)]{t} = + \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, v - \frac{1}{e_1 \; e_2 \; e_3} \lt( \pd[(e_2 \, e_3 \, u \, u)]{i} + \pd[(e_1 \, e_3 \, v \, u)]{j} \rt) \\ 
     655    - \frac{1}{e_3} \pd[(\omega \, u)]{k} - \frac{1}{e_1} \pd[]{i} \lt( \frac{p_s + p_h}{\rho_o} \rt) - g \frac{\rho}{\rho_o} \sigma_1 + D_u^{\vect U} + F_u^{\vect U} 
    666656  \end{multline*} 
    667657  \begin{multline*} 
    668658  % \label{eq:MB_sco_v_flux} 
    669     \frac{1}{e_3} \pd[(e_3 \, v)]{t} = - \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, u \\ 
    670                                        - \frac{1}{e_1 \; e_2 \; e_3} \lt( \pd[( e_2 \; e_3 \, u \, v)]{i} + \pd[(e_1 \; e_3 \, v \, v)]{j} \rt) \\ 
    671                                        - \frac{1}{e_3} \pd[(\omega \, v)]{k} 
    672                                        - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) 
    673                                        - g \frac{\rho}{\rho_o}\sigma_2 + D_v^{\vect U} + F_v^{\vect U} 
     659    \frac{1}{e_3} \pd[(e_3 \, v)]{t} = - \lt[ f + \frac{1}{e_1 \; e_2} \lt( v \pd[e_2]{i} - u \pd[e_1]{j} \rt) \rt] \, u - \frac{1}{e_1 \; e_2 \; e_3} \lt( \pd[( e_2 \; e_3 \, u \, v)]{i} + \pd[(e_1 \; e_3 \, v \, v)]{j} \rt) \\ 
     660    - \frac{1}{e_3} \pd[(\omega \, v)]{k} - \frac{1}{e_2} \pd[]{j} \lt( \frac{p_s + p_h}{\rho_o} \rt) - g \frac{\rho}{\rho_o}\sigma_2 + D_v^{\vect U} + F_v^{\vect U} 
    674661  \end{multline*} 
    675662  where the relative vorticity, $\zeta$, the surface pressure gradient, 
     
    682669    \chi = \frac{1}{e_1 e_2 e_3} \lt( \pd[(e_2 e_3 \, u)]{i} + \pd[(e_1 e_3 \, v)]{j} \rt) 
    683670  \] 
    684 \item \textit{tracer equations}: 
    685   \begin{multline*} 
    686   % \label{eq:MB_sco_t} 
    687     \frac{1}{e_3} \pd[(e_3 \, T)]{t} = - \frac{1}{e_1 e_2 e_3} \lt(   \pd[(e_2 e_3 \, u \, T)]{i} 
    688                                                                     + \pd[(e_1 e_3 \, v \, T)]{j} \rt) \\ 
    689                                        - \frac{1}{e_3} \pd[(T \, \omega)]{k} + D^T + F^S 
    690   \end{multline*} 
    691   \begin{multline} 
    692   % \label{eq:MB_sco_s} 
    693     \frac{1}{e_3} \pd[(e_3 \, S)]{t} = - \frac{1}{e_1 e_2 e_3} \lt(   \pd[(e_2 e_3 \, u \, S)]{i} 
    694                                                                     + \pd[(e_1 e_3 \, v \, S)]{j} \rt) \\ 
    695                                        - \frac{1}{e_3} \pd[(S \, \omega)]{k} + D^S + F^S 
    696   \end{multline} 
    697 \end{itemize} 
     671\item [Tracer equations] 
     672  \begin{gather*} 
     673    % \label{eq:MB_sco_t} 
     674    \frac{1}{e_3} \pd[(e_3 \, T)]{t} = - \frac{1}{e_1 e_2 e_3} \lt(   \pd[(e_2 e_3 \, u \, T)]{i} + \pd[(e_1 e_3 \, v \, T)]{j} \rt) - \frac{1}{e_3} \pd[(T \, \omega)]{k} + D^T + F^S \\ 
     675    % \label{eq:MB_sco_s} 
     676    \frac{1}{e_3} \pd[(e_3 \, S)]{t} = - \frac{1}{e_1 e_2 e_3} \lt(   \pd[(e_2 e_3 \, u \, S)]{i} + \pd[(e_1 e_3 \, v \, S)]{j} \rt) - \frac{1}{e_3} \pd[(S \, \omega)]{k} + D^S + F^S 
     677  \end{gather*} 
     678\end{description} 
    698679The equation of state has the same expression as in $z$-coordinate, 
    699680and similar expressions are used for mixing and forcing terms. 
     
    709690\label{subsec:MB_zco_star} 
    710691 
    711 \begin{figure}[!b] 
     692\begin{figure} 
    712693  \centering 
    713694  \includegraphics[width=0.66\textwidth]{Fig_z_zstar} 
    714695  \caption[Curvilinear z-coordinate systems (\{non-\}linear free-surface cases and re-scaled \zstar)]{ 
    715     (a) $z$-coordinate in linear free-surface case ; 
    716     (b) $z$-coordinate in non-linear free surface case ; 
    717     (c) re-scaled height coordinate 
    718     (become popular as the \zstar-coordinate \citep{adcroft.campin_OM04}).} 
     696    \begin{enumerate*}[label={(\alph*)}] 
     697    \item $z$-coordinate in linear free-surface case; 
     698    \item $z$-coordinate in non-linear free surface case; 
     699    \item re-scaled height coordinate 
     700      (become popular as the \zstar-coordinate \citep{adcroft.campin_OM04}). 
     701    \end{enumerate*} 
     702  } 
    719703  \label{fig:MB_z_zstar} 
    720704\end{figure} 
    721705 
    722 In this case, the free surface equation is nonlinear, and the variations of volume are fully taken into account. 
    723 These coordinates systems is presented in a report \citep{levier.treguier.ea_rpt07} available on the \NEMO\ web site. 
    724  
    725 The \zstar coordinate approach is an unapproximated, non-linear free surface implementation which allows one to 
    726 deal with large amplitude free-surface variations relative to the vertical resolution \citep{adcroft.campin_OM04}. 
    727 In the \zstar formulation, 
    728 the variation of the column thickness due to sea-surface undulations is not concentrated in the surface level, 
    729 as in the $z$-coordinate formulation, but is equally distributed over the full water column. 
     706In this case, the free surface equation is nonlinear, 
     707and the variations of volume are fully taken into account. 
     708These coordinates systems is presented in a report \citep{levier.treguier.ea_rpt07} available on 
     709the \NEMO\ web site. 
     710 
     711The \zstar\ coordinate approach is an unapproximated, non-linear free surface implementation which 
     712allows one to deal with large amplitude free-surface variations relative to the vertical resolution 
     713\citep{adcroft.campin_OM04}. 
     714In the \zstar\ formulation, the variation of the column thickness due to sea-surface undulations is 
     715not concentrated in the surface level, as in the $z$-coordinate formulation, 
     716but is equally distributed over the full water column. 
    730717Thus vertical levels naturally follow sea-surface variations, with a linear attenuation with depth, 
    731718as illustrated by \autoref{fig:MB_z_zstar}. 
    732 Note that with a flat bottom, such as in \autoref{fig:MB_z_zstar}, the bottom-following $z$ coordinate and \zstar are equivalent. 
     719Note that with a flat bottom, such as in \autoref{fig:MB_z_zstar}, 
     720the bottom-following $z$ coordinate and \zstar\ are equivalent. 
    733721The definition and modified oceanic equations for the rescaled vertical coordinate \zstar, 
    734722including the treatment of fresh-water flux at the surface, are detailed in Adcroft and Campin (2004). 
    735723The major points are summarized here. 
    736 The position (\zstar) and vertical discretization (\zstar) are expressed as: 
     724The position (\zstar) and vertical discretization ($\delta \zstar$) are expressed as: 
    737725\[ 
    738726  % \label{eq:MB_z-star} 
    739   H + \zstar = (H + z)  / r \quad \text{and}  \quad \delta \zstar 
    740               = \delta z / r \quad \text{with} \quad r 
    741               = \frac{H + \eta}{H} . 
     727  H + \zstar = (H + z)  / r \quad \text{and}  \quad \delta \zstar = \delta z / r \quad \text{with} \quad r = \frac{H + \eta}{H} 
    742728\] 
    743729Simple re-organisation of the above expressions gives 
    744730\[ 
    745731  % \label{eq:MB_zstar_2} 
    746   \zstar = H \lt( \frac{z - \eta}{H + \eta} \rt) . 
     732  \zstar = H \lt( \frac{z - \eta}{H + \eta} \rt) 
    747733\] 
    748734Since the vertical displacement of the free surface is incorporated in the vertical coordinate \zstar, 
    749 the upper and lower boundaries are at fixed  \zstar position, 
     735the upper and lower boundaries are at fixed \zstar\ position, 
    750736$\zstar = 0$ and $\zstar = -H$ respectively. 
    751737Also the divergence of the flow field is no longer zero as shown by the continuity equation: 
    752738\[ 
    753   \pd[r]{t} = \nabla_{\zstar} \cdot \lt( r \; \vect U_h \rt) + \pd[r \; w^*]{\zstar} = 0 . 
     739  \pd[r]{t} = \nabla_{\zstar} \cdot \lt( r \; \vect U_h \rt) + \pd[r \; w^*]{\zstar} = 0 
    754740\] 
    755 This \zstar coordinate is closely related to the "eta" coordinate used in many atmospheric models 
    756 (see Black (1994) for a review of eta coordinate atmospheric models). 
     741This \zstar\ coordinate is closely related to the $\eta$ coordinate used in many atmospheric models 
     742(see Black (1994) for a review of $\eta$ coordinate atmospheric models). 
    757743It was originally used in ocean models by Stacey et al. (1995) for studies of tides next to shelves, 
    758744and it has been recently promoted by Adcroft and Campin (2004) for global climate modelling. 
    759745 
    760 The surfaces of constant \zstar are quasi -horizontal. 
    761 Indeed, the \zstar coordinate reduces to $z$ when $\eta$ is zero. 
     746The surfaces of constant \zstar\ are quasi-horizontal. 
     747Indeed, the \zstar\ coordinate reduces to $z$ when $\eta$ is zero. 
    762748In general, when noting the large differences between 
    763749undulations of the bottom topography versus undulations in the surface height, 
    764 it is clear that surfaces constant \zstar are very similar to the depth surfaces. 
     750it is clear that surfaces constant \zstar\ are very similar to the depth surfaces. 
    765751These properties greatly reduce difficulties of computing the horizontal pressure gradient relative to 
    766752terrain following sigma models discussed in \autoref{subsec:MB_sco}. 
    767753Additionally, since $\zstar = z$ when $\eta = 0$, 
    768 no flow is spontaneously generated in an unforced ocean starting from rest, regardless the bottom topography. 
    769 This behaviour is in contrast to the case with "s"-models, where pressure gradient errors in the presence of 
    770 nontrivial topographic variations can generate nontrivial spontaneous flow from a resting state, 
     754no flow is spontaneously generated in an unforced ocean starting from rest, 
     755regardless the bottom topography. 
     756This behaviour is in contrast to the case with ``s''-models, 
     757where pressure gradient errors in the presence of nontrivial topographic variations can generate 
     758nontrivial spontaneous flow from a resting state, 
    771759depending on the sophistication of the pressure gradient solver. 
    772 The quasi -horizontal nature of the coordinate surfaces also facilitates the implementation of 
    773 neutral physics parameterizations in \zstar  models using the same techniques as in $z$-models 
     760The quasi-horizontal nature of the coordinate surfaces also facilitates the implementation of 
     761neutral physics parameterizations in \zstar\ models using the same techniques as in $z$-models 
    774762(see Chapters 13-16 of \cite{griffies_bk04}) for a discussion of neutral physics in $z$-models, 
    775763as well as \autoref{sec:LDF_slp} in this document for treatment in \NEMO). 
    776764 
    777 The range over which \zstar  varies is time independent $-H \leq \zstar \leq 0$. 
     765The range over which \zstar\ varies is time independent $-H \leq \zstar \leq 0$. 
    778766Hence, all cells remain nonvanishing, so long as the surface height maintains $\eta > -H$. 
    779 This is a minor constraint relative to that encountered on the surface height when using $s = z$ or $s = z - \eta$. 
    780  
    781 Because \zstar  has a time independent range, all grid cells have static increments ds, 
     767This is a minor constraint relative to that encountered on the surface height when 
     768using $s = z$ or $s = z - \eta$. 
     769 
     770Because \zstar\ has a time independent range, all grid cells have static increments ds, 
    782771and the sum of the vertical increments yields the time independent ocean depth. %k ds = H. 
    783 The \zstar coordinate is therefore invisible to undulations of the free surface, 
     772The \zstar\ coordinate is therefore invisible to undulations of the free surface, 
    784773since it moves along with the free surface. 
    785 This property means that no spurious vertical transport is induced across surfaces of constant \zstar  by 
    786 the motion of external gravity waves. 
     774This property means that no spurious vertical transport is induced across 
     775surfaces of constant \zstar\  by the motion of external gravity waves. 
    787776Such spurious transport can be a problem in z-models, especially those with tidal forcing. 
    788 Quite generally, the time independent range for the \zstar  coordinate is a very convenient property that 
    789 allows for a nearly arbitrary vertical resolution even in the presence of large amplitude fluctuations of 
    790 the surface height, again so long as $\eta > -H$. 
     777Quite generally, 
     778the time independent range for the \zstar\ coordinate is a very convenient property that 
     779allows for a nearly arbitrary vertical resolution even in 
     780the presence of large amplitude fluctuations of the surface height, again so long as $\eta > -H$. 
    791781%end MOM doc %%% 
    792782 
     
    795785\label{subsec:MB_sco} 
    796786 
    797 %% ================================================================================================= 
    798 \subsubsection{Introduction} 
    799  
    800787Several important aspects of the ocean circulation are influenced by bottom topography. 
    801 Of course, the most important is that bottom topography determines deep ocean sub-basins, barriers, sills and 
    802 channels that strongly constrain the path of water masses, but more subtle effects exist. 
    803 For example, the topographic $\beta$-effect is usually larger than the planetary one along continental slopes. 
     788Of course, the most important is that bottom topography determines deep ocean sub-basins, barriers, 
     789sills and channels that strongly constrain the path of water masses, but more subtle effects exist. 
     790For example, 
     791the topographic $\beta$-effect is usually larger than the planetary one along continental slopes. 
    804792Topographic Rossby waves can be excited and can interact with the mean current. 
    805793In the $z$-coordinate system presented in the previous section (\autoref{sec:MB_zco}), 
     
    809797large localized depth gradients associated with large localized vertical velocities. 
    810798The response to such a velocity field often leads to numerical dispersion effects. 
    811 One solution to strongly reduce this error is to use a partial step representation of bottom topography instead of 
    812 a full step one \cite{pacanowski.gnanadesikan_MWR98}. 
     799One solution to strongly reduce this error is to use a partial step representation of 
     800bottom topography instead of a full step one \cite{pacanowski.gnanadesikan_MWR98}. 
    813801Another solution is to introduce a terrain-following coordinate system (hereafter $s$-coordinate). 
    814802 
    815 The $s$-coordinate avoids the discretisation error in the depth field since the layers of 
    816 computation are gradually adjusted with depth to the ocean bottom. 
    817 Relatively small topographic features as well as  gentle, large-scale slopes of the sea floor in the deep ocean, 
    818 which would be ignored in typical $z$-model applications with the largest grid spacing at greatest depths, 
     803The $s$-coordinate avoids the discretisation error in the depth field since 
     804the layers of computation are gradually adjusted with depth to the ocean bottom. 
     805Relatively small topographic features as well as 
     806gentle, large-scale slopes of the sea floor in the deep ocean, 
     807which would be ignored in typical $z$-model applications with 
     808the largest grid spacing at greatest depths, 
    819809can easily be represented (with relatively low vertical resolution). 
    820 A terrain-following model (hereafter $s$-model) also facilitates the modelling of the boundary layer flows over 
    821 a large depth range, which in the framework of the $z$-model would require high vertical resolution over 
     810A terrain-following model (hereafter $s$-model) also facilitates 
     811the modelling of the boundary layer flows over a large depth range, 
     812which in the framework of the $z$-model would require high vertical resolution over 
    822813the whole depth range. 
    823 Moreover, with a $s$-coordinate it is possible, at least in principle, to have the bottom and the sea surface as 
     814Moreover, 
     815with a $s$-coordinate it is possible, at least in principle, to have the bottom and the sea surface as 
    824816the only boundaries of the domain (no more lateral boundary condition to specify). 
    825 Nevertheless, a $s$-coordinate also has its drawbacks. Perfectly adapted to a homogeneous ocean, 
     817Nevertheless, a $s$-coordinate also has its drawbacks. 
     818Perfectly adapted to a homogeneous ocean, 
    826819it has strong limitations as soon as stratification is introduced. 
    827820The main two problems come from the truncation error in the horizontal pressure gradient and 
    828821a possibly increased diapycnal diffusion. 
    829822The horizontal pressure force in $s$-coordinate consists of two terms (see \autoref{apdx:SCOORD}), 
    830  
    831823\begin{equation} 
    832824  \label{eq:MB_p_sco} 
     
    836828The second term in \autoref{eq:MB_p_sco} depends on the tilt of the coordinate surface and 
    837829leads to a truncation error that is not present in a $z$-model. 
    838 In the special case of a $\sigma$-coordinate (i.e. a depth-normalised coordinate system $\sigma = z/H$), 
    839 \citet{haney_JPO91} and \citet{beckmann.haidvogel_JPO93} have given estimates of the magnitude of this truncation error. 
    840 It depends on topographic slope, stratification, horizontal and vertical resolution, the equation of state, 
    841 and the finite difference scheme. 
     830In the special case of a $\sigma$-coordinate 
     831(i.e. a depth-normalised coordinate system $\sigma = z/H$), 
     832\citet{haney_JPO91} and \citet{beckmann.haidvogel_JPO93} have given estimates of 
     833the magnitude of this truncation error. 
     834It depends on topographic slope, stratification, horizontal and vertical resolution, 
     835the equation of state, and the finite difference scheme. 
    842836This error limits the possible topographic slopes that a model can handle at 
    843837a given horizontal and vertical resolution. 
    844838This is a severe restriction for large-scale applications using realistic bottom topography. 
    845 The large-scale slopes require high horizontal resolution, and the computational cost becomes prohibitive. 
     839The large-scale slopes require high horizontal resolution, 
     840and the computational cost becomes prohibitive. 
    846841This problem can be at least partially overcome by mixing $s$-coordinate and 
    847 step-like representation of bottom topography \citep{gerdes_JGR93*a,gerdes_JGR93*b,madec.delecluse.ea_JPO96}. 
     842step-like representation of bottom topography 
     843\citep{gerdes_JGR93*a,gerdes_JGR93*b,madec.delecluse.ea_JPO96}. 
    848844However, the definition of the model domain vertical coordinate becomes then a non-trivial thing for 
    849845a realistic bottom topography: 
    850 an envelope topography is defined in $s$-coordinate on which a full or 
    851 partial step bottom topography is then applied in order to adjust the model depth to the observed one 
    852 (see \autoref{subsec:DOM_zgr}. 
    853  
    854 For numerical reasons a minimum of diffusion is required along the coordinate surfaces of 
    855 any finite difference model. 
     846an envelope topography is defined in $s$-coordinate on which 
     847a full or partial step bottom topography is then applied in order to 
     848adjust the model depth to the observed one (see \autoref{subsec:DOM_zgr}). 
     849 
     850For numerical reasons a minimum of diffusion is required along 
     851the coordinate surfaces of any finite difference model. 
    856852It causes spurious diapycnal mixing when coordinate surfaces do not coincide with isoneutral surfaces. 
    857853This is the case for a $z$-model as well as for a $s$-model. 
    858 However, density varies more strongly on $s$-surfaces than on horizontal surfaces in regions of 
    859 large topographic slopes, implying larger diapycnal diffusion in a $s$-model than in a $z$-model. 
    860 Whereas such a diapycnal diffusion in a $z$-model tends to weaken horizontal density (pressure) gradients and thus 
    861 the horizontal circulation, it usually reinforces these gradients in a $s$-model, creating spurious circulation. 
    862 For example, imagine an isolated bump of topography in an ocean at rest with a horizontally uniform stratification. 
     854However, density varies more strongly on $s$-surfaces than on horizontal surfaces in 
     855regions of large topographic slopes, 
     856implying larger diapycnal diffusion in a $s$-model than in a $z$-model. 
     857Whereas such a diapycnal diffusion in a $z$-model tends to 
     858weaken horizontal density (pressure) gradients and thus the horizontal circulation, 
     859it usually reinforces these gradients in a $s$-model, creating spurious circulation. 
     860For example, imagine an isolated bump of topography in 
     861an ocean at rest with a horizontally uniform stratification. 
    863862Spurious diffusion along $s$-surfaces will induce a bump of isoneutral surfaces over the topography, 
    864863and thus will generate there a baroclinic eddy. 
    865864In contrast, the ocean will stay at rest in a $z$-model. 
    866 As for the truncation error, the problem can be reduced by introducing the terrain-following coordinate below 
    867 the strongly stratified portion of the water column (\ie\ the main thermocline) \citep{madec.delecluse.ea_JPO96}. 
    868 An alternate solution consists of rotating the lateral diffusive tensor to geopotential or to isoneutral surfaces 
    869 (see \autoref{subsec:MB_ldf}). 
     865As for the truncation error, the problem can be reduced by 
     866introducing the terrain-following coordinate below the strongly stratified portion of the water column 
     867(\ie\ the main thermocline) \citep{madec.delecluse.ea_JPO96}. 
     868An alternate solution consists of rotating the lateral diffusive tensor to 
     869geopotential or to isoneutral surfaces (see \autoref{subsec:MB_ldf}). 
    870870Unfortunately, the slope of isoneutral surfaces relative to the $s$-surfaces can very large, 
    871 strongly exceeding the stability limit of such a operator when it is discretized (see \autoref{chap:LDF}). 
    872  
    873 The $s$-coordinates introduced here \citep{lott.madec.ea_OM90,madec.delecluse.ea_JPO96} differ mainly in two aspects from 
    874 similar models: 
    875 it allows a representation of bottom topography with mixed full or partial step-like/terrain following topography; 
     871strongly exceeding the stability limit of such a operator when it is discretized 
     872(see \autoref{chap:LDF}). 
     873 
     874The $s$-coordinates introduced here \citep{lott.madec.ea_OM90,madec.delecluse.ea_JPO96} 
     875differ mainly in two aspects from similar models: 
     876it allows a representation of bottom topography with 
     877mixed full or partial step-like/terrain following topography; 
    876878It also offers a completely general transformation, $s=s(i,j,z)$ for the vertical coordinate. 
    877879 
    878880%% ================================================================================================= 
    879 \subsection{\texorpdfstring{Curvilinear \ztilde-coordinate}{}} 
     881\subsection{Curvilinear \ztilde-coordinate} 
    880882\label{subsec:MB_zco_tilde} 
    881883 
    882 The \ztilde -coordinate has been developed by \citet{leclair.madec_OM11}. 
     884The \ztilde-coordinate has been developed by \citet{leclair.madec_OM11}. 
    883885It is available in \NEMO\ since the version 3.4 and is more robust in version 4.0 than previously. 
    884886Nevertheless, it is currently not robust enough to be used in all possible configurations. 
     
    889891\label{sec:MB_zdf_ldf} 
    890892 
    891 The hydrostatic primitive equations describe the behaviour of a geophysical fluid at space and time scales larger than 
    892 a few kilometres in the horizontal, a few meters in the vertical and a few minutes. 
    893 They are usually solved at larger scales: the specified grid spacing and time step of the numerical model. 
     893The hydrostatic primitive equations describe the behaviour of a geophysical fluid at 
     894space and time scales larger than a few kilometres in the horizontal, 
     895a few meters in the vertical and a few minutes. 
     896They are usually solved at larger scales: the specified grid spacing and time step of 
     897the numerical model. 
    894898The effects of smaller scale motions (coming from the advective terms in the Navier-Stokes equations) 
    895899must be represented entirely in terms of large-scale patterns to close the equations. 
     
    902906small scale processes \textit{in fine} balance the surface input of kinetic energy and heat. 
    903907 
    904 The control exerted by gravity on the flow induces a strong anisotropy between the lateral and vertical motions. 
     908The control exerted by gravity on the flow induces a strong anisotropy between 
     909the lateral and vertical motions. 
    905910Therefore subgrid-scale physics \textbf{D}$^{\vect U}$, $D^{S}$ and $D^{T}$  in 
    906911\autoref{eq:MB_PE_dyn}, \autoref{eq:MB_PE_tra_T} and \autoref{eq:MB_PE_tra_S} are divided into 
    907 a lateral part \textbf{D}$^{l \vect U}$, $D^{l S}$ and $D^{l T}$ and 
     912a  lateral part \textbf{D}$^{l \vect U}$, $D^{l S}$ and $D^{l T}$ and 
    908913a vertical part \textbf{D}$^{v \vect U}$, $D^{v S}$ and $D^{v T}$. 
    909 The formulation of these terms and their underlying physics are briefly discussed in the next two subsections. 
     914The formulation of these terms and their underlying physics are briefly discussed in 
     915the next two subsections. 
    910916 
    911917%% ================================================================================================= 
     
    913919\label{subsec:MB_zdf} 
    914920 
    915 The model resolution is always larger than the scale at which the major sources of vertical turbulence occur 
    916 (shear instability, internal wave breaking...). 
     921The model resolution is always larger than the scale at which 
     922the major sources of vertical turbulence occur (shear instability, internal wave breaking...). 
    917923Turbulent motions are thus never explicitly solved, even partially, but always parameterized. 
    918 The vertical turbulent fluxes are assumed to depend linearly on the gradients of large-scale quantities 
    919 (for example, the turbulent heat flux is given by $\overline{T' w'} = -A^{v T} \partial_z \overline T$, 
     924The vertical turbulent fluxes are assumed to depend linearly on 
     925the gradients of large-scale quantities (for example, 
     926the turbulent heat flux is given by $\overline{T' w'} = -A^{v T} \partial_z \overline T$, 
    920927where $A^{v T}$ is an eddy coefficient). 
    921928This formulation is analogous to that of molecular diffusion and dissipation. 
     
    927934  \label{eq:MB_zdf} 
    928935  \begin{gathered} 
    929     \vect D^{v \vect U} = \pd[]{z} \lt( A^{vm} \pd[\vect U_h]{z} \rt) \ , \\ 
    930           D^{vT}       = \pd[]{z} \lt( A^{vT} \pd[T]{z}         \rt) \quad \text{and} \quad 
     936    \vect D^{v \vect U} = \pd[]{z} \lt( A^{vm} \pd[\vect U_h]{z} \rt) \text{,} \ 
     937          D^{vT}       = \pd[]{z} \lt( A^{vT} \pd[T]{z}         \rt) \ \text{and} \ 
    931938          D^{vS}       = \pd[]{z} \lt( A^{vT} \pd[S]{z}         \rt) 
    932939  \end{gathered} 
    933940\end{equation} 
    934 where $A^{vm}$ and $A^{vT}$ are the vertical eddy viscosity and diffusivity coefficients, respectively. 
     941where $A^{vm}$ and $A^{vT}$ are the vertical eddy viscosity and diffusivity coefficients, 
     942respectively. 
    935943At the sea surface and at the bottom, turbulent fluxes of momentum, heat and salt must be specified 
    936944(see \autoref{chap:SBC} and \autoref{chap:ZDF} and \autoref{sec:TRA_bbl}). 
    937945All the vertical physics is embedded in the specification of the eddy coefficients. 
    938946They can be assumed to be either constant, or function of the local fluid properties 
    939 (\eg\ Richardson number, Brunt-Vais\"{a}l\"{a} frequency, distance from the boundary ...), 
     947(\eg\ Richardson number, Brunt-V\"{a}is\"{a}l\"{a} frequency, distance from the boundary ...), 
    940948or computed from a turbulent closure model. 
    941949The choices available in \NEMO\ are discussed in \autoref{chap:ZDF}). 
     
    948956(which can be solved explicitly if the resolution is sufficient since 
    949957their underlying physics are included in the primitive equations), 
    950 and a sub mesoscale turbulence which is never explicitly solved even partially, but always parameterized. 
    951 The formulation of lateral eddy fluxes depends on whether the mesoscale is below or above the grid-spacing 
    952 (\ie\ the model is eddy-resolving or not). 
     958and a sub mesoscale turbulence which is never explicitly solved even partially, 
     959but always parameterized. 
     960The formulation of lateral eddy fluxes depends on whether 
     961the mesoscale is below or above the grid-spacing (\ie\ the model is eddy-resolving or not). 
    953962 
    954963In non-eddy-resolving configurations, the closure is similar to that used for the vertical physics. 
    955 The lateral turbulent fluxes are assumed to depend linearly on the lateral gradients of large-scale quantities. 
     964The lateral turbulent fluxes are assumed to depend linearly on 
     965the lateral gradients of large-scale quantities. 
    956966The resulting lateral diffusive and dissipative operators are of second order. 
    957 Observations show that lateral mixing induced by mesoscale turbulence tends to be along isopycnal surfaces 
     967Observations show that 
     968lateral mixing induced by mesoscale turbulence tends to be along isopycnal surfaces 
    958969(or more precisely neutral surfaces \cite{mcdougall_JPO87}) rather than across them. 
    959 As the slope of neutral surfaces is small in the ocean, a common approximation is to assume that 
    960 the `lateral' direction is the horizontal, \ie\ the lateral mixing is performed along geopotential surfaces. 
     970As the slope of neutral surfaces is small in the ocean, 
     971a common approximation is to assume that the ``lateral'' direction is the horizontal, 
     972\ie\ the lateral mixing is performed along geopotential surfaces. 
    961973This leads to a geopotential second order operator for lateral subgrid scale physics. 
    962 This assumption can be relaxed: the eddy-induced turbulent fluxes can be better approached by assuming that 
     974This assumption can be relaxed: 
     975the eddy-induced turbulent fluxes can be better approached by assuming that 
    963976they depend linearly on the gradients of large-scale quantities computed along neutral surfaces. 
    964977In such a case, the diffusive operator is an isoneutral second order operator and 
    965978it has components in the three space directions. 
    966 However, 
    967 both horizontal and isoneutral operators have no effect on mean (\ie\ large scale) potential energy whereas 
     979However, both horizontal and isoneutral operators have no effect on 
     980mean (\ie\ large scale) potential energy whereas 
    968981potential energy is a main source of turbulence (through baroclinic instabilities). 
    969982\citet{gent.mcwilliams_JPO90} proposed a parameterisation of mesoscale eddy-induced turbulence which 
    970983associates an eddy-induced velocity to the isoneutral diffusion. 
    971984Its mean effect is to reduce the mean potential energy of the ocean. 
    972 This leads to a formulation of lateral subgrid-scale physics made up of an isoneutral second order operator and 
    973 an eddy induced advective part. 
     985This leads to a formulation of lateral subgrid-scale physics made up of 
     986an isoneutral second order operator and an eddy induced advective part. 
    974987In all these lateral diffusive formulations, 
    975988the specification of the lateral eddy coefficients remains the problematic point as 
    976 there is no really satisfactory formulation of these coefficients as a function of large-scale features. 
     989there is no really satisfactory formulation of these coefficients as 
     990a function of large-scale features. 
    977991 
    978992In eddy-resolving configurations, a second order operator can be used, 
     
    983997not interfering with the resolved mesoscale activity. 
    984998Another approach is becoming more and more popular: 
    985 instead of specifying explicitly a sub-grid scale term in the momentum and tracer time evolution equations, 
     999instead of specifying explicitly a sub-grid scale term in 
     1000the momentum and tracer time evolution equations, 
    9861001one uses an advective scheme which is diffusive enough to maintain the model stability. 
    987 It must be emphasised that then, all the sub-grid scale physics is included in the formulation of 
    988 the advection scheme. 
     1002It must be emphasised that then, 
     1003all the sub-grid scale physics is included in the formulation of the advection scheme. 
    9891004 
    9901005All these parameterisations of subgrid scale physics have advantages and drawbacks. 
    991 They are not all available in \NEMO. For active tracers (temperature and salinity) the main ones are: 
     1006They are not all available in \NEMO. 
     1007For active tracers (temperature and salinity) the main ones are: 
    9921008Laplacian and bilaplacian operators acting along geopotential or iso-neutral surfaces, 
    9931009\citet{gent.mcwilliams_JPO90} parameterisation, and various slightly diffusive advection schemes. 
    994 For momentum, the main ones are: Laplacian and bilaplacian operators acting along geopotential surfaces, 
     1010For momentum, the main ones are: 
     1011Laplacian and bilaplacian operators acting along geopotential surfaces, 
    9951012and UBS advection schemes when flux form is chosen for the momentum advection. 
    9961013 
     
    10011018\begin{equation} 
    10021019  \label{eq:MB_iso_tensor} 
    1003   D^{lT} = \nabla \vect . \lt( A^{lT} \; \Re \; \nabla T \rt) \quad \text{with} \quad 
    1004   \Re = 
    1005     \begin{pmatrix} 
    1006       1    & 0    & -r_1          \\ 
    1007       0    & 1    & -r_2          \\ 
    1008       -r_1 & -r_2 & r_1^2 + r_2^2 \\ 
    1009     \end{pmatrix} 
     1020  D^{lT} = \nabla \vect . \lt( A^{lT} \; \Re \; \nabla T \rt) \quad \text{with} \quad \Re = 
     1021  \begin{pmatrix} 
     1022    1    & 0    & -r_1          \\ 
     1023    0    & 1    & -r_2          \\ 
     1024    -r_1 & -r_2 & r_1^2 + r_2^2 \\ 
     1025  \end{pmatrix} 
    10101026\end{equation} 
    10111027where $r_1$ and $r_2$ are the slopes between the surface along which the diffusive operator acts and 
     
    10141030the rotation between geopotential and $s$-surfaces, 
    10151031while it is only an approximation for the rotation between isoneutral and $z$- or $s$-surfaces. 
    1016 Indeed, in the latter case, two assumptions are made to simplify \autoref{eq:MB_iso_tensor} \citep{cox_OM87}. 
    1017 First, the horizontal contribution of the dianeutral mixing is neglected since the ratio between iso and 
    1018 dia-neutral diffusive coefficients is known to be several orders of magnitude smaller than unity. 
     1032Indeed, in the latter case, two assumptions are made to simplify \autoref{eq:MB_iso_tensor} 
     1033\citep{cox_OM87}. 
     1034First, the horizontal contribution of the dianeutral mixing is neglected since 
     1035the ratio between iso and dia-neutral diffusive coefficients is known to be 
     1036several orders of magnitude smaller than unity. 
    10191037Second, the two isoneutral directions of diffusion are assumed to be independent since 
    10201038the slopes are generally less than $10^{-2}$ in the ocean (see \autoref{apdx:DIFFOPERS}). 
     
    10271045they are equal to $\sigma_1$ and $\sigma_2$, respectively (see \autoref{eq:MB_sco_slope}). 
    10281046 
    1029 For \textit{isoneutral} diffusion $r_1$ and $r_2$ are the slopes between the isoneutral and computational surfaces. 
     1047For \textit{isoneutral} diffusion $r_1$ and $r_2$ are the slopes between 
     1048the isoneutral and computational surfaces. 
    10301049Therefore, they are different quantities, but have similar expressions in $z$- and $s$-coordinates. 
    10311050In $z$-coordinates: 
     
    10481067where $ \vect U^\ast = \lt( u^\ast,v^\ast,w^\ast \rt)$ is a non-divergent, 
    10491068eddy-induced transport velocity. This velocity field is defined by: 
    1050 \begin{gather} 
     1069\begin{gather*} 
    10511070  % \label{eq:MB_eiv} 
    1052   u^\ast =   \frac{1}{e_3}            \pd[]{k} \lt( A^{eiv} \;        \tilde{r}_1 \rt) \\ 
     1071  u^\ast =   \frac{1}{e_3}            \pd[]{k} \lt( A^{eiv} \;        \tilde{r}_1 \rt) \quad 
    10531072  v^\ast =   \frac{1}{e_3}            \pd[]{k} \lt( A^{eiv} \;        \tilde{r}_2 \rt) \\ 
    10541073  w^\ast = - \frac{1}{e_1 e_2} \lt[   \pd[]{i} \lt( A^{eiv} \; e_2 \, \tilde{r}_1 \rt) 
    10551074                                     + \pd[]{j} \lt( A^{eiv} \; e_1 \, \tilde{r}_2 \rt) \rt] 
    1056 \end{gather} 
     1075\end{gather*} 
    10571076where $A^{eiv}$ is the eddy induced velocity coefficient 
    10581077(or equivalently the isoneutral thickness diffusivity coefficient), 
    1059 and $\tilde r_1$ and $\tilde r_2$ are the slopes between isoneutral and \textit{geopotential} surfaces. 
    1060 Their values are thus independent of the vertical coordinate, but their expression depends on the coordinate: 
    1061 \begin{align} 
     1078and $\tilde r_1$ and $\tilde r_2$ are the slopes between 
     1079isoneutral and \textit{geopotential} surfaces. 
     1080Their values are thus independent of the vertical coordinate, 
     1081but their expression depends on the coordinate: 
     1082\begin{equation} 
    10621083  \label{eq:MB_slopes_eiv} 
    10631084  \tilde{r}_n = 
     
    10671088    \end{cases} 
    10681089  \quad \text{where~} n = 1, 2 
    1069 \end{align} 
     1090\end{equation} 
    10701091 
    10711092The normal component of the eddy induced velocity is zero at all the boundaries. 
    1072 This can be achieved in a model by tapering either the eddy coefficient or the slopes to zero in the vicinity of 
    1073 the boundaries. 
     1093This can be achieved in a model by tapering either the eddy coefficient or the slopes to zero in 
     1094the vicinity of the boundaries. 
    10741095The latter strategy is used in \NEMO\ (cf. \autoref{chap:LDF}). 
    10751096 
     
    11011122\end{align*} 
    11021123 
    1103 Such a formulation ensures a complete separation between the vorticity and horizontal divergence fields 
    1104 (see \autoref{apdx:INVARIANTS}). 
     1124Such a formulation ensures a complete separation between 
     1125the vorticity and horizontal divergence fields (see \autoref{apdx:INVARIANTS}). 
    11051126Unfortunately, it is only available in \textit{iso-level} direction. 
    11061127When a rotation is required 
    1107 (\ie\ geopotential diffusion in $s$-coordinates or isoneutral diffusion in both $z$- and $s$-coordinates), 
    1108 the $u$ and $v$-fields are considered as independent scalar fields, so that the diffusive operator is given by: 
     1128(\ie\ geopotential diffusion in $s$-coordinates or 
     1129isoneutral diffusion in both $z$- and $s$-coordinates), 
     1130the $u$ and $v$-fields are considered as independent scalar fields, 
     1131so that the diffusive operator is given by: 
    11091132\begin{gather*} 
    11101133  % \label{eq:MB_lapU_iso} 
     
    11221145\subsubsection{Lateral bilaplacian momentum diffusive operator} 
    11231146 
    1124 As for tracers, the bilaplacian order momentum diffusive operator is a re-entering Laplacian operator with 
     1147As for tracers, 
     1148the bilaplacian order momentum diffusive operator is a re-entering Laplacian operator with 
    11251149the harmonic eddy diffusion coefficient set to the square root of the biharmonic one. 
    11261150Nevertheless it is currently not available in the iso-neutral case. 
Note: See TracChangeset for help on using the changeset viewer.