Changeset 13463 for NEMO/branches/2019/dev_r11351_fldread_with_XIOS/doc/latex/NEMO/subfiles/chap_DOM.tex
- Timestamp:
- 2020-09-14T17:40:34+02:00 (4 years ago)
- Location:
- NEMO/branches/2019/dev_r11351_fldread_with_XIOS
- Files:
-
- 6 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11351_fldread_with_XIOS
- Property svn:externals
-
old new 3 3 ^/utils/build/mk@HEAD mk 4 4 ^/utils/tools@HEAD tools 5 ^/vendors/AGRIF/dev @HEADext/AGRIF5 ^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS ext/AGRIF 6 6 ^/vendors/FCM@HEAD ext/FCM 7 7 ^/vendors/IOIPSL@HEAD ext/IOIPSL 8 9 # SETTE 10 ^/utils/CI/sette@13382 sette
-
- Property svn:externals
-
NEMO/branches/2019/dev_r11351_fldread_with_XIOS/doc
-
Property
svn:externals
set to
^/utils/badges badges
^/utils/logos logos
-
Property
svn:externals
set to
-
NEMO/branches/2019/dev_r11351_fldread_with_XIOS/doc/latex
- Property svn:ignore deleted
-
NEMO/branches/2019/dev_r11351_fldread_with_XIOS/doc/latex/NEMO
-
Property
svn:externals
set to
^/utils/figures/NEMO figures
-
Property
svn:externals
set to
-
NEMO/branches/2019/dev_r11351_fldread_with_XIOS/doc/latex/NEMO/subfiles
- Property svn:ignore
-
old new 1 *.aux 2 *.bbl 3 *.blg 4 *.dvi 5 *.fdb* 6 *.fls 7 *.idx 1 *.ind 8 2 *.ilg 9 *.ind10 *.log11 *.maf12 *.mtc*13 *.out14 *.pdf15 *.toc16 _minted-*
-
- Property svn:ignore
-
NEMO/branches/2019/dev_r11351_fldread_with_XIOS/doc/latex/NEMO/subfiles/chap_DOM.tex
r11330 r13463 2 2 3 3 \begin{document} 4 % ================================================================ 5 % Chapter 2 ——— Space and Time Domain (DOM) 6 % ================================================================ 4 7 5 \chapter{Space Domain (DOM)} 8 6 \label{chap:DOM} 9 7 10 \minitoc 11 12 % Missing things: 13 % - istate: description of the initial state ==> this has to be put elsewhere.. 14 % perhaps in MISC ? By the way the initialisation of T S and dynamics 15 % should be put outside of DOM routine (better with TRC staff and off-line 16 % tracers) 17 % -geo2ocean: how to switch from geographic to mesh coordinate 18 % - domclo: closed sea and lakes.... management of closea sea area : specific to global configuration, both forced and coupled 19 20 \vfill 21 \begin{figure}[b] 22 \subsubsection*{Changes record} 23 \begin{tabular}{m{0.08\linewidth}||m{0.32\linewidth}|m{0.6\linewidth}} 24 Release & Author(s) & Modifications \\ 25 \hline 26 {\em 4.0} & {\em Simon M{\"u}ller \& Andrew Coward} & {\em Compatibility changes for v4.0. Major simplication has moved many of the options to external domain configuration tools. For now this information has been retained in \autoref{apdx:DOMAINcfg} } \\ 27 {\em 3.x} & {\em Sebastien Masson, Gurvan Madec \& Rashid Benshila } & {\em } \\ 28 \end{tabular} 29 \end{figure} 30 31 \newpage 32 33 Having defined the continuous equations in \autoref{chap:PE} and chosen a time discretization \autoref{chap:STP}, 34 we need to choose a grid for spatial discretization and related numerical algorithms. 8 % Missing things 9 % - istate: description of the initial state ==> this has to be put elsewhere.. 10 % perhaps in MISC ? By the way the initialisation of T S and dynamics 11 % should be put outside of DOM routine (better with TRC staff and off-line 12 % tracers) 13 % - geo2ocean: how to switch from geographic to mesh coordinate 14 % - domclo: closed sea and lakes.... 15 % management of closea sea area: specific to global cfg, both forced and coupled 16 17 \thispagestyle{plain} 18 19 \chaptertoc 20 21 \paragraph{Changes record} ~\\ 22 23 {\footnotesize 24 \begin{tabularx}{0.8\textwidth}{l||X|X} 25 Release & 26 Author(s) & 27 Modifications \\ 28 \hline 29 {\em 4.0 } & 30 {\em Simon M\"{u}ller \& Andrew Coward \newline \newline 31 Simona Flavoni and Tim Graham } & 32 {\em Compatibility changes: many options moved to external domain configuration tools 33 (see \autoref{apdx:DOMCFG}). \newline 34 Updates } \\ 35 {\em 3.6 } & 36 {\em Rachid Benshila, Christian \'{E}th\'{e}, Pierre Mathiot and Gurvan Madec } & 37 {\em Updates } \\ 38 {\em $\leq$ 3.4 } & 39 {\em Gurvan Madec and S\'{e}bastien Masson } & 40 {\em First version } 41 \end{tabularx} 42 } 43 44 \clearpage 45 46 Having defined the continuous equations in \autoref{chap:MB} and 47 chosen a time discretisation \autoref{chap:TD}, 48 we need to choose a grid for spatial discretisation and related numerical algorithms. 35 49 In the present chapter, we provide a general description of the staggered grid used in \NEMO, 36 and other relevant information about the DOM (DOMain) source-code modules . 37 38 % ================================================================ 39 % Fundamentals of the Discretisation 40 % ================================================================ 50 and other relevant information about the DOM (DOMain) source code modules. 51 52 %% ================================================================================================= 41 53 \section{Fundamentals of the discretisation} 42 54 \label{sec:DOM_basics} 43 55 44 % ------------------------------------------------------------------------------------------------------------- 45 % Arrangement of Variables 46 % ------------------------------------------------------------------------------------------------------------- 56 %% ================================================================================================= 47 57 \subsection{Arrangement of variables} 48 58 \label{subsec:DOM_cell} 49 59 50 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 51 \begin{figure}[!tb] 52 \begin{center} 53 \includegraphics[width=\textwidth]{Fig_cell} 54 \caption{ 55 \protect\label{fig:cell} 56 Arrangement of variables. 57 $t$ indicates scalar points where temperature, salinity, density, pressure and 58 horizontal divergence are defined. 59 $(u,v,w)$ indicates vector points, and $f$ indicates vorticity points where both relative and 60 planetary vorticities are defined. 61 } 62 \end{center} 60 \begin{figure} 61 \centering 62 \includegraphics[width=0.33\textwidth]{DOM_cell} 63 \caption[Arrangement of variables in the unit cell of space domain]{ 64 Arrangement of variables in the unit cell of space domain. 65 $t$ indicates scalar points where 66 temperature, salinity, density, pressure and horizontal divergence are defined. 67 $(u,v,w)$ indicates vector points, and $f$ indicates vorticity points where 68 both relative and planetary vorticities are defined.} 69 \label{fig:DOM_cell} 63 70 \end{figure} 64 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 65 66 The numerical techniques used to solve the Primitive Equations in this model are based on the traditional, 67 centred second-order finite difference approximation. 71 72 The numerical techniques used to solve the Primitive Equations in this model are based on 73 the traditional, centred second-order finite difference approximation. 68 74 Special attention has been given to the homogeneity of the solution in the three spatial directions. 69 75 The arrangement of variables is the same in all directions. 70 It consists of cells centred on scalar points ($t$, $S$, $p$, $\rho$) with vector points $(u, v, w)$ defined in 71 the centre of each face of the cells (\autoref{fig:cell}). 72 This is the generalisation to three dimensions of the well-known ``C'' grid in Arakawa's classification 73 \citep{mesinger.arakawa_bk76}. 74 The relative and planetary vorticity, $\zeta$ and $f$, are defined in the centre of each vertical edge and 75 the barotropic stream function $\psi$ is defined at horizontal points overlying the $\zeta$ and $f$-points. 76 77 The ocean mesh (\ie the position of all the scalar and vector points) is defined by the transformation that 78 gives $(\lambda,\varphi,z)$ as a function of $(i,j,k)$. 79 The grid-points are located at integer or integer and a half value of $(i,j,k)$ as indicated on \autoref{tab:cell}. 80 In all the following, subscripts $u$, $v$, $w$, $f$, $uw$, $vw$ or $fw$ indicate the position of 81 the grid-point where the scale factors are defined. 82 Each scale factor is defined as the local analytical value provided by \autoref{eq:scale_factors}. 76 It consists of cells centred on scalar points ($t$, $S$, $p$, $\rho$) with 77 vector points $(u, v, w)$ defined in the centre of each face of the cells (\autoref{fig:DOM_cell}). 78 This is the generalisation to three dimensions of the well-known ``C'' grid in 79 Arakawa's classification \citep{mesinger.arakawa_bk76}. 80 The relative and planetary vorticity, $\zeta$ and $f$, are defined in the centre of each 81 vertical edge and the barotropic stream function $\psi$ is defined at horizontal points overlying 82 the $\zeta$ and $f$-points. 83 84 The ocean mesh (\ie\ the position of all the scalar and vector points) is defined by 85 the transformation that gives $(\lambda,\varphi,z)$ as a function of $(i,j,k)$. 86 The grid-points are located at integer or integer and a half value of $(i,j,k)$ as indicated on 87 \autoref{tab:DOM_cell}. 88 In all the following, 89 subscripts $u$, $v$, $w$, $f$, $uw$, $vw$ or $fw$ indicate the position of the grid-point where 90 the scale factors are defined. 91 Each scale factor is defined as the local analytical value provided by \autoref{eq:MB_scale_factors}. 83 92 As a result, the mesh on which partial derivatives $\pd[]{\lambda}$, $\pd[]{\varphi}$ and 84 93 $\pd[]{z}$ are evaluated is a uniform mesh with a grid size of unity. 85 Discrete partial derivatives are formulated by the traditional, centred second order finite difference approximation 86 while the scale factors are chosen equal to their local analytical value. 94 Discrete partial derivatives are formulated by 95 the traditional, centred second order finite difference approximation while 96 the scale factors are chosen equal to their local analytical value. 87 97 An important point here is that the partial derivative of the scale factors must be evaluated by 88 98 centred finite difference approximation, not from their analytical expression. 89 This preserves the symmetry of the discrete set of equations and therefore satisfies many of90 the continuous properties (see \autoref{apdx:C}).99 This preserves the symmetry of the discrete set of equations and 100 therefore satisfies many of the continuous properties (see \autoref{apdx:INVARIANTS}). 91 101 A similar, related remark can be made about the domain size: 92 when needed, an area, volume, or the total ocean depth must be evaluated as the product or sum of the relevant scale factors 93 (see \autoref{eq:DOM_bar} in the next section). 94 95 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 96 \begin{table}[!tb] 97 \begin{center} 98 \begin{tabular}{|p{46pt}|p{56pt}|p{56pt}|p{56pt}|} 99 \hline 100 t & $i $ & $j $ & $k $ \\ 101 \hline 102 u & $i + 1/2$ & $j $ & $k $ \\ 103 \hline 104 v & $i $ & $j + 1/2$ & $k $ \\ 105 \hline 106 w & $i $ & $j $ & $k + 1/2$ \\ 107 \hline 108 f & $i + 1/2$ & $j + 1/2$ & $k $ \\ 109 \hline 110 uw & $i + 1/2$ & $j $ & $k + 1/2$ \\ 111 \hline 112 vw & $i $ & $j + 1/2$ & $k + 1/2$ \\ 113 \hline 114 fw & $i + 1/2$ & $j + 1/2$ & $k + 1/2$ \\ 115 \hline 116 \end{tabular} 117 \caption{ 118 \protect\label{tab:cell} 119 Location of grid-points as a function of integer or integer and a half value of the column, line or level. 120 This indexing is only used for the writing of the semi -discrete equations. 121 In the code, the indexing uses integer values only and is positive downwards in the vertical with $k=1$ at the surface. 122 (see \autoref{subsec:DOM_Num_Index}) 123 } 124 \end{center} 102 when needed, an area, volume, or the total ocean depth must be evaluated as 103 the product or sum of the relevant scale factors (see \autoref{eq:DOM_bar} in the next section). 104 105 \begin{table} 106 \centering 107 \begin{tabular}{|l|l|l|l|} 108 \hline 109 t & $i $ & $j $ & $k $ \\ 110 \hline 111 u & $i + 1/2$ & $j $ & $k $ \\ 112 \hline 113 v & $i $ & $j + 1/2$ & $k $ \\ 114 \hline 115 w & $i $ & $j $ & $k + 1/2$ \\ 116 \hline 117 f & $i + 1/2$ & $j + 1/2$ & $k $ \\ 118 \hline 119 uw & $i + 1/2$ & $j $ & $k + 1/2$ \\ 120 \hline 121 vw & $i $ & $j + 1/2$ & $k + 1/2$ \\ 122 \hline 123 fw & $i + 1/2$ & $j + 1/2$ & $k + 1/2$ \\ 124 \hline 125 \end{tabular} 126 \caption[Location of grid-points]{ 127 Location of grid-points as a function of integer or 128 integer and a half value of the column, line or level. 129 This indexing is only used for the writing of the semi-discrete equations. 130 In the code, the indexing uses integer values only and 131 is positive downwards in the vertical with $k=1$ at the surface. 132 (see \autoref{subsec:DOM_Num_Index})} 133 \label{tab:DOM_cell} 125 134 \end{table} 126 %>>>>>>>>>>>>>>>>>>>>>>>>>>>>127 135 128 136 Note that the definition of the scale factors 129 (\ie as the analytical first derivative of the transformation that137 (\ie\ as the analytical first derivative of the transformation that 130 138 results in $(\lambda,\varphi,z)$ as a function of $(i,j,k)$) 131 is specific to the \NEMO model \citep{marti.madec.ea_JGR92}.139 is specific to the \NEMO\ model \citep{marti.madec.ea_JGR92}. 132 140 As an example, a scale factor in the $i$ direction is defined locally at a $t$-point, 133 141 whereas many other models on a C grid choose to define such a scale factor as … … 136 144 firstly, there is no ambiguity in the scale factors appearing in the discrete equations, 137 145 since they are first introduced in the continuous equations; 138 secondly, analytical transformations encourage good practice by the definition of smoothly varying grids 139 (rather than allowing the user to set arbitrary jumps in thickness between adjacent layers) \citep{treguier.dukowicz.ea_JGR96}. 140 An example of the effect of such a choice is shown in \autoref{fig:zgr_e3}. 141 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 142 \begin{figure}[!t] 143 \begin{center} 144 \includegraphics[width=\textwidth]{Fig_zgr_e3} 145 \caption{ 146 \protect\label{fig:zgr_e3} 147 Comparison of (a) traditional definitions of grid-point position and grid-size in the vertical, 148 and (b) analytically derived grid-point position and scale factors. 149 For both grids here, the same $w$-point depth has been chosen but 150 in (a) the $t$-points are set half way between $w$-points while 151 in (b) they are defined from an analytical function: 152 $z(k) = 5 \, (k - 1/2)^3 - 45 \, (k - 1/2)^2 + 140 \, (k - 1/2) - 150$. 153 Note the resulting difference between the value of the grid-size $\Delta_k$ and 154 those of the scale factor $e_k$. 155 } 156 \end{center} 146 secondly, analytical transformations encourage good practice by 147 the definition of smoothly varying grids 148 (rather than allowing the user to set arbitrary jumps in thickness between adjacent layers) 149 \citep{treguier.dukowicz.ea_JGR96}. 150 An example of the effect of such a choice is shown in \autoref{fig:DOM_zgr_e3}. 151 \begin{figure} 152 \centering 153 \includegraphics[width=0.5\textwidth]{DOM_zgr_e3} 154 \caption[Comparison of grid-point position, vertical grid-size and scale factors]{ 155 Comparison of (a) traditional definitions of grid-point position and grid-size in the vertical, 156 and (b) analytically derived grid-point position and scale factors. 157 For both grids here, the same $w$-point depth has been chosen but 158 in (a) the $t$-points are set half way between $w$-points while 159 in (b) they are defined from an analytical function: 160 $z(k) = 5 \, (k - 1/2)^3 - 45 \, (k - 1/2)^2 + 140 \, (k - 1/2) - 150$. 161 Note the resulting difference between the value of the grid-size $\Delta_k$ and 162 those of the scale factor $e_k$.} 163 \label{fig:DOM_zgr_e3} 157 164 \end{figure} 158 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 159 160 % ------------------------------------------------------------------------------------------------------------- 161 % Vector Invariant Formulation 162 % ------------------------------------------------------------------------------------------------------------- 165 166 %% ================================================================================================= 163 167 \subsection{Discrete operators} 164 168 \label{subsec:DOM_operators} 165 169 166 Given the values of a variable $q$ at adjacent points, the differencing and averaging operators at167 the midpoint between them are:170 Given the values of a variable $q$ at adjacent points, 171 the differencing and averaging operators at the midpoint between them are: 168 172 \begin{alignat*}{2} 169 % \label{eq: di_mi}173 % \label{eq:DOM_di_mi} 170 174 \delta_i [q] &= & &q (i + 1/2) - q (i - 1/2) \\ 171 175 \overline q^{\, i} &= &\big\{ &q (i + 1/2) + q (i - 1/2) \big\} / 2 … … 173 177 174 178 Similar operators are defined with respect to $i + 1/2$, $j$, $j + 1/2$, $k$, and $k + 1/2$. 175 Following \autoref{eq:PE_grad} and \autoref{eq:PE_lap}, the gradient of a variable $q$ defined at 176 a $t$-point has its three components defined at $u$-, $v$- and $w$-points while 179 Following \autoref{eq:MB_grad} and \autoref{eq:MB_lap}, 180 the gradient of a variable $q$ defined at a $t$-point has 181 its three components defined at $u$-, $v$- and $w$-points while 177 182 its Laplacian is defined at the $t$-point. 178 183 These operators have the following discrete forms in the curvilinear $s$-coordinates system: 179 \ [184 \begin{gather*} 180 185 % \label{eq:DOM_grad} 181 186 \nabla q \equiv \frac{1}{e_{1u}} \delta_{i + 1/2} [q] \; \, \vect i 182 187 + \frac{1}{e_{2v}} \delta_{j + 1/2} [q] \; \, \vect j 183 + \frac{1}{e_{3w}} \delta_{k + 1/2} [q] \; \, \vect k 184 \] 185 \begin{multline*} 188 + \frac{1}{e_{3w}} \delta_{k + 1/2} [q] \; \, \vect k \\ 186 189 % \label{eq:DOM_lap} 187 190 \Delta q \equiv \frac{1}{e_{1t} \, e_{2t} \, e_{3t}} 188 191 \; \lt[ \delta_i \lt( \frac{e_{2u} \, e_{3u}}{e_{1u}} \; \delta_{i + 1/2} [q] \rt) 189 + \delta_j \lt( \frac{e_{1v} \, e_{3v}}{e_{2v}} \; \delta_{j + 1/2} [q] \rt) \; \rt] \\192 + \delta_j \lt( \frac{e_{1v} \, e_{3v}}{e_{2v}} \; \delta_{j + 1/2} [q] \rt) \; \rt] 190 193 + \frac{1}{e_{3t}} 191 194 \delta_k \lt[ \frac{1 }{e_{3w}} \; \delta_{k + 1/2} [q] \rt] 192 \end{multline*} 193 194 Following \autoref{eq:PE_curl} and \autoref{eq:PE_div}, a vector $\vect A = (a_1,a_2,a_3)$ defined at 195 vector points $(u,v,w)$ has its three curl components defined at $vw$-, $uw$, and $f$-points, and 195 \end{gather*} 196 197 Following \autoref{eq:MB_curl} and \autoref{eq:MB_div}, 198 a vector $\vect A = (a_1,a_2,a_3)$ defined at vector points $(u,v,w)$ has 199 its three curl components defined at $vw$-, $uw$, and $f$-points, and 196 200 its divergence defined at $t$-points: 197 \begin{multline }201 \begin{multline*} 198 202 % \label{eq:DOM_curl} 199 203 \nabla \times \vect A \equiv \frac{1}{e_{2v} \, e_{3vw}} … … 206 210 \Big[ \delta_{i + 1/2} (e_{2v} \, a_2) 207 211 - \delta_{j + 1/2} (e_{1u} \, a_1) \Big] \vect k 208 \end{multline }209 \ begin{equation}212 \end{multline*} 213 \[ 210 214 % \label{eq:DOM_div} 211 215 \nabla \cdot \vect A \equiv \frac{1}{e_{1t} \, e_{2t} \, e_{3t}} 212 216 \Big[ \delta_i (e_{2u} \, e_{3u} \, a_1) + \delta_j (e_{1v} \, e_{3v} \, a_2) \Big] 213 217 + \frac{1}{e_{3t}} \delta_k (a_3) 214 \ end{equation}215 216 The vertical average over the whole water column is denoted by an overbar and is for217 a masked field $q$ (\iea quantity that is equal to zero inside solid areas):218 \] 219 220 The vertical average over the whole water column is denoted by an overbar and 221 is for a masked field $q$ (\ie\ a quantity that is equal to zero inside solid areas): 218 222 \begin{equation} 219 223 \label{eq:DOM_bar} … … 221 225 \end{equation} 222 226 where $H_q$ is the ocean depth, which is the masked sum of the vertical scale factors at $q$ points, 223 $k^b$ and $k^o$ are the bottom and surface $k$-indices, and the symbol $\sum \limits_k$ refers to a summation over 224 all grid points of the same type in the direction indicated by the subscript (here $k$). 227 $k^b$ and $k^o$ are the bottom and surface $k$-indices, 228 and the symbol $\sum \limits_k$ refers to a summation over all grid points of the same type in 229 the direction indicated by the subscript (here $k$). 225 230 226 231 In continuous form, the following properties are satisfied: … … 232 237 \end{gather} 233 238 234 It is straightforward to demonstrate that these properties are verified locally in discrete form as soon as235 the scalar $q$ is taken at $t$-points and the vector $\vect A$ has its components defined at239 It is straightforward to demonstrate that these properties are verified locally in discrete form as 240 soon as the scalar $q$ is taken at $t$-points and the vector $\vect A$ has its components defined at 236 241 vector points $(u,v,w)$. 237 242 238 243 Let $a$ and $b$ be two fields defined on the mesh, with a value of zero inside continental areas. 239 It can be shown that the differencing operators ($\delta_i$, $\delta_j$ and $\delta_k$) 240 are skew-symmetric linear operators, and further that the averaging operators $\overline{\cdots}^{\, i}$, 241 $\overline{\cdots}^{\, j}$ and $\overline{\cdots}^{\, k}$) are symmetric linear operators, \ie 242 \begin{alignat}{4} 244 It can be shown that the differencing operators ($\delta_i$, $\delta_j$ and 245 $\delta_k$) are skew-symmetric linear operators, 246 and further that the averaging operators ($\overline{\cdots}^{\, i}$, $\overline{\cdots}^{\, j}$ and 247 $\overline{\cdots}^{\, k}$) are symmetric linear operators, \ie 248 \begin{alignat}{5} 243 249 \label{eq:DOM_di_adj} 244 250 &\sum \limits_i a_i \; \delta_i [b] &\equiv &- &&\sum \limits_i \delta _{ i + 1/2} [a] &b_{i + 1/2} \\ … … 247 253 \end{alignat} 248 254 249 In other words, the adjoint of the differencing and averaging operators are $\delta_i^* = \delta_{i + 1/2}$ and 255 In other words, 256 the adjoint of the differencing and averaging operators are $\delta_i^* = \delta_{i + 1/2}$ and 250 257 $(\overline{\cdots}^{\, i})^* = \overline{\cdots}^{\, i + 1/2}$, respectively. 251 These two properties will be used extensively in the \autoref{apdx: C} to258 These two properties will be used extensively in the \autoref{apdx:INVARIANTS} to 252 259 demonstrate integral conservative properties of the discrete formulation chosen. 253 260 254 % ------------------------------------------------------------------------------------------------------------- 255 % Numerical Indexing 256 % ------------------------------------------------------------------------------------------------------------- 261 %% ================================================================================================= 257 262 \subsection{Numerical indexing} 258 263 \label{subsec:DOM_Num_Index} 259 264 260 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 261 \begin{figure}[!tb] 262 \begin{center} 263 \includegraphics[width=\textwidth]{Fig_index_hor} 264 \caption{ 265 \protect\label{fig:index_hor} 266 Horizontal integer indexing used in the \fortran code. 267 The dashed area indicates the cell in which variables contained in arrays have the same $i$- and $j$-indices 268 } 269 \end{center} 265 \begin{figure} 266 \centering 267 \includegraphics[width=0.33\textwidth]{DOM_index_hor} 268 \caption[Horizontal integer indexing]{ 269 Horizontal integer indexing used in the \fortran\ code. 270 The dashed area indicates the cell in which 271 variables contained in arrays have the same $i$- and $j$-indices} 272 \label{fig:DOM_index_hor} 270 273 \end{figure} 271 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 272 273 The array representation used in the \fortran code requires an integer indexing. 274 However, the analytical definition of the mesh (see \autoref{subsec:DOM_cell}) is associated with the use of 275 integer values for $t$-points only whileall the other points involve integer and a half values.274 275 The array representation used in the \fortran\ code requires an integer indexing. 276 However, the analytical definition of the mesh (see \autoref{subsec:DOM_cell}) is associated with 277 the use of integer values for $t$-points only while 278 all the other points involve integer and a half values. 276 279 Therefore, a specific integer indexing has been defined for points other than $t$-points 277 (\ie velocity and vorticity grid-points). 278 Furthermore, the direction of the vertical indexing has been reversed and the surface level set at $k = 1$. 279 280 % ----------------------------------- 281 % Horizontal Indexing 282 % ----------------------------------- 280 (\ie\ velocity and vorticity grid-points). 281 Furthermore, the direction of the vertical indexing has been reversed and 282 the surface level set at $k = 1$. 283 284 %% ================================================================================================= 283 285 \subsubsection{Horizontal indexing} 284 286 \label{subsec:DOM_Num_Index_hor} 285 287 286 The indexing in the horizontal plane has been chosen as shown in \autoref{fig: index_hor}.288 The indexing in the horizontal plane has been chosen as shown in \autoref{fig:DOM_index_hor}. 287 289 For an increasing $i$ index ($j$ index), 288 290 the $t$-point and the eastward $u$-point (northward $v$-point) have the same index 289 (see the dashed area in \autoref{fig:index_hor}). 290 A $t$-point and its nearest northeast $f$-point have the same $i$-and $j$-indices. 291 292 % ----------------------------------- 293 % Vertical indexing 294 % ----------------------------------- 291 (see the dashed area in \autoref{fig:DOM_index_hor}). 292 A $t$-point and its nearest north-east $f$-point have the same $i$-and $j$-indices. 293 294 %% ================================================================================================= 295 295 \subsubsection{Vertical indexing} 296 296 \label{subsec:DOM_Num_Index_vertical} 297 297 298 In the vertical, the chosen indexing requires special attention since the direction of the $k$-axis in299 the \fortran code is the reverse of that used in the semi -discrete equations and300 given in \autoref{subsec:DOM_cell}.301 The sea surface corresponds to the $w$-level $k = 1$, which is the same index as the $t$-level just below302 (\autoref{fig:index_vert}).298 In the vertical, the chosen indexing requires special attention since 299 the direction of the $k$-axis in the \fortran\ code is the reverse of 300 that used in the semi-discrete equations and given in \autoref{subsec:DOM_cell}. 301 The sea surface corresponds to the $w$-level $k = 1$, 302 which is the same index as the $t$-level just below (\autoref{fig:DOM_index_vert}). 303 303 The last $w$-level ($k = jpk$) either corresponds to or is below the ocean floor while 304 the last $t$-level is always outside the ocean domain (\autoref{fig:index_vert}). 305 Note that a $w$-point and the directly underlaying $t$-point have a common $k$ index (\ie $t$-points and their 306 nearest $w$-point neighbour in negative index direction), in contrast to the indexing on the horizontal plane where 307 the $t$-point has the same index as the nearest velocity points in the positive direction of the respective horizontal axis index 308 (compare the dashed area in \autoref{fig:index_hor} and \autoref{fig:index_vert}). 304 the last $t$-level is always outside the ocean domain (\autoref{fig:DOM_index_vert}). 305 Note that a $w$-point and the directly underlaying $t$-point have a common $k$ index 306 (\ie\ $t$-points and their nearest $w$-point neighbour in negative index direction), 307 in contrast to the indexing on the horizontal plane where 308 the $t$-point has the same index as the nearest velocity points in 309 the positive direction of the respective horizontal axis index 310 (compare the dashed area in \autoref{fig:DOM_index_hor} and \autoref{fig:DOM_index_vert}). 309 311 Since the scale factors are chosen to be strictly positive, 310 a \textit{minus sign} is included in the \fortran implementations of \textit{all the vertical derivatives} of 311 the discrete equations given in this manual in order to accommodate the opposing vertical index directions in implementation and documentation. 312 313 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 314 \begin{figure}[!pt] 315 \begin{center} 316 \includegraphics[width=\textwidth]{Fig_index_vert} 317 \caption{ 318 \protect\label{fig:index_vert} 319 Vertical integer indexing used in the \fortran code. 320 Note that the $k$-axis is oriented downward. 321 The dashed area indicates the cell in which variables contained in arrays have a common $k$-index. 322 } 323 \end{center} 312 a \textit{minus sign} is included in the \fortran\ implementations of 313 \textit{all the vertical derivatives} of the discrete equations given in this manual in order to 314 accommodate the opposing vertical index directions in implementation and documentation. 315 316 \begin{figure} 317 \centering 318 \includegraphics[width=0.33\textwidth]{DOM_index_vert} 319 \caption[Vertical integer indexing]{ 320 Vertical integer indexing used in the \fortran\ code. 321 Note that the $k$-axis is oriented downward. 322 The dashed area indicates the cell in which 323 variables contained in arrays have a common $k$-index.} 324 \label{fig:DOM_index_vert} 324 325 \end{figure} 325 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 326 327 % ------------------------------------------------------------------------------------------------------------- 328 % Domain configuration 329 % ------------------------------------------------------------------------------------------------------------- 326 327 %% ================================================================================================= 330 328 \section{Spatial domain configuration} 331 329 \label{subsec:DOM_config} 332 330 333 \nlst{namcfg} 334 335 Two typical methods are available to specify the spatial domain 336 configuration; they can be selected using parameter \np{ln\_read\_cfg} 337 parameter in namelist \ngn{namcfg}. 338 339 If \np{ln\_read\_cfg} is set to \forcode{.true.}, the domain-specific parameters 340 and fields are read from a netCDF input file, whose name (without its .nc 341 suffix) can be specified as the value of the \np{cn\_domcfg} parameter in 342 namelist \ngn{namcfg}. 343 344 If \np{ln\_read\_cfg} is set to \forcode{.false.}, the domain-specific 345 parameters and fields can be provided (\eg analytically computed) by subroutines 346 \mdl{usrdef\_hgr} and \mdl{usrdef\_zgr}. These subroutines can be supplied in 347 the \path{MY_SRC} directory of the configuration, and default versions that 348 configure the spatial domain for the GYRE reference configuration are present in 349 the \path{src/OCE/USR} directory. 350 351 In version 4.0 there are no longer any options for reading complex bathmetries and 352 performing a vertical discretization at run-time. Whilst it is occasionally convenient 353 to have a common bathymetry file and, for example, to run similar models with and 354 without partial bottom boxes and/or sigma-coordinates, supporting such choices leads to 355 overly complex code. Worse still is the difficulty of ensuring the model configurations 356 intended to be identical are indeed so when the model domain itself can be altered by runtime 357 selections. The code previously used to perform vertical discretization has be incorporated 358 into an external tool (\path{tools/DOMAINcfg}) which is briefly described in \autoref{apdx:DOMAINcfg}. 359 360 The next subsections summarise the parameter and fields related to the 361 configuration of the whole model domain. These represent the minimum information 362 that must be provided either via the \np{cn\_domcfg} file or set by code 363 inserted into user-supplied versions of the \mdl{usrdef\_*} subroutines. The 364 requirements are presented in three sections: the domain size 365 (\autoref{subsec:DOM_size}), the horizontal mesh 366 (\autoref{subsec:DOM_hgr}), and the vertical grid 367 (\autoref{subsec:DOM_zgr}). 368 369 % ----------------------------------- 370 % Domain Size 371 % ----------------------------------- 331 Two typical methods are available to specify the spatial domain configuration; 332 they can be selected using parameter \np{ln_read_cfg}{ln\_read\_cfg} parameter in 333 namelist \nam{cfg}{cfg}. 334 335 If \np{ln_read_cfg}{ln\_read\_cfg} is set to \forcode{.true.}, 336 the domain-specific parameters and fields are read from a NetCDF input file, 337 whose name (without its .nc suffix) can be specified as 338 the value of the \np{cn_domcfg}{cn\_domcfg} parameter in namelist \nam{cfg}{cfg}. 339 340 If \np{ln_read_cfg}{ln\_read\_cfg} is set to \forcode{.false.}, 341 the domain-specific parameters and fields can be provided (\eg\ analytically computed) by 342 subroutines \mdl{usrdef\_hgr} and \mdl{usrdef\_zgr}. 343 These subroutines can be supplied in the \path{MY_SRC} directory of the configuration, 344 and default versions that configure the spatial domain for the GYRE reference configuration are 345 present in the \path{./src/OCE/USR} directory. 346 347 In version 4.0 there are no longer any options for reading complex bathymetries and 348 performing a vertical discretisation at run-time. 349 Whilst it is occasionally convenient to have a common bathymetry file and, for example, 350 to run similar models with and without partial bottom boxes and/or sigma-coordinates, 351 supporting such choices leads to overly complex code. 352 Worse still is the difficulty of ensuring the model configurations intended to be identical are 353 indeed so when the model domain itself can be altered by runtime selections. 354 The code previously used to perform vertical discretisation has been incorporated into 355 an external tool (\path{./tools/DOMAINcfg}) which is briefly described in \autoref{apdx:DOMCFG}. 356 357 The next subsections summarise the parameter and fields related to 358 the configuration of the whole model domain. 359 These represent the minimum information that must be provided either via 360 the \np{cn_domcfg}{cn\_domcfg} file or 361 set by code inserted into user-supplied versions of the \texttt{usrdef\_*} subroutines. 362 The requirements are presented in three sections: 363 the domain size (\autoref{subsec:DOM_size}), the horizontal mesh (\autoref{subsec:DOM_hgr}), 364 and the vertical grid (\autoref{subsec:DOM_zgr}). 365 366 %% ================================================================================================= 372 367 \subsection{Domain size} 373 368 \label{subsec:DOM_size} 374 369 375 The total size of the computational domain is set by the parameters 376 \np{jpiglo}, \np{jpjglo} and \np{jpkglo} for the $i$, $j$ and $k$ 377 directions, respectively. Note, that the variables \forcode{jpi} and \forcode{jpj} 378 refer to the size of each processor subdomain when the code is run in 379 parallel using domain decomposition (\key{mpp\_mpi} defined, see 380 \autoref{sec:LBC_mpp}). 381 382 The name of the configuration is set through parameter \np{cn\_cfg}, 383 and the nominal resolution through parameter \np{nn\_cfg} (unless in 384 the input file both of variables \forcode{ORCA} and \forcode{ORCA_index} 385 are present, in which case \np{cn\_cfg} and \np{nn\_cfg} are set from these 386 values accordingly). 387 388 The global lateral boundary condition type is selected from 8 options 389 using parameter \np{jperio}. See \autoref{sec:LBC_jperio} for 390 details on the available options and the corresponding values for 391 \np{jperio}. 392 393 % ================================================================ 394 % Domain: Horizontal Grid (mesh) 395 % ================================================================ 396 \subsection{Horizontal grid mesh (\protect\mdl{domhgr})} 370 The total size of the computational domain is set by the parameters \jp{jpiglo}, \jp{jpjglo} and 371 \jp{jpkglo} for the $i$, $j$ and $k$ directions, respectively. 372 Note, that the variables \texttt{jpi} and \texttt{jpj} refer to 373 the size of each processor subdomain when the code is run in parallel using domain decomposition 374 (\key{mpp\_mpi} defined, see \autoref{sec:LBC_mpp}). 375 376 The name of the configuration is set through parameter \np{cn_cfg}{cn\_cfg}, 377 and the nominal resolution through parameter \np{nn_cfg}{nn\_cfg} 378 (unless in the input file both of variables \texttt{ORCA} and \texttt{ORCA\_index} are present, 379 in which case \np{cn_cfg}{cn\_cfg} and \np{nn_cfg}{nn\_cfg} are set from these values accordingly). 380 381 The global lateral boundary condition type is selected from 8 options using parameter \jp{jperio}. 382 See \autoref{sec:LBC_jperio} for details on the available options and 383 the corresponding values for \jp{jperio}. 384 385 %% ================================================================================================= 386 \subsection[Horizontal grid mesh (\textit{domhgr.F90}]{Horizontal grid mesh (\protect\mdl{domhgr})} 397 387 \label{subsec:DOM_hgr} 398 388 399 % ================================================================ 400 % Domain: List of hgr-related fields needed 401 % ================================================================ 389 %% ================================================================================================= 402 390 \subsubsection{Required fields} 403 391 \label{sec:DOM_hgr_fields} 404 The explicit specification of a range of mesh-related fields are required for the definition of a configuration. These include: 405 406 \begin{Verbatim}[fontsize=\tiny] 407 int jpiglo, jpjglo, jpkglo /* global domain sizes */ 408 int jperio /* lateral global domain b.c. */ 409 double glamt, glamu, glamv, glamf /* geographic longitude (t,u,v and f points respectively) */ 410 double gphit, gphiu, gphiv, gphif /* geographic latitude */ 411 double e1t, e1u, e1v, e1f /* horizontal scale factors */ 412 double e2t, e2u, e2v, e2f /* horizontal scale factors */ 413 \end{Verbatim} 414 415 The values of the geographic longitude and latitude arrays at indices $i,j$ correspond to the analytical expressions of the longitude $\lambda$ and latitude $\varphi$ as a function of $(i,j)$, evaluated at the values as specified in Table \autoref{tab:cell} for the respective grid-point position. The calculation of the values of the horizontal scale factor arrays in general additionally involves partial derivatives of $\lambda$ and $\varphi$ with respect to $i$ and $j$, evaluated for the same arguments as $\lambda$ and $\varphi$. 416 392 393 The explicit specification of a range of mesh-related fields are required for 394 the definition of a configuration. 395 These include: 396 397 \begin{clines} 398 int jpiglo, jpjglo, jpkglo /* global domain sizes */ 399 int jperio /* lateral global domain b.c. */ 400 double glamt, glamu, glamv, glamf /* geographic longitude (t,u,v and f points respectively) */ 401 double gphit, gphiu, gphiv, gphif /* geographic latitude */ 402 double e1t, e1u, e1v, e1f /* horizontal scale factors */ 403 double e2t, e2u, e2v, e2f /* horizontal scale factors */ 404 \end{clines} 405 406 The values of the geographic longitude and latitude arrays at indices $i,j$ correspond to 407 the analytical expressions of the longitude $\lambda$ and latitude $\varphi$ as a function of $(i,j)$, 408 evaluated at the values as specified in \autoref{tab:DOM_cell} for the respective grid-point position. 409 The calculation of the values of the horizontal scale factor arrays in general additionally involves 410 partial derivatives of $\lambda$ and $\varphi$ with respect to $i$ and $j$, 411 evaluated for the same arguments as $\lambda$ and $\varphi$. 412 413 %% ================================================================================================= 417 414 \subsubsection{Optional fields} 418 \begin{Verbatim}[fontsize=\tiny] 419 /* Optional: */ 420 int ORCA, ORCA_index /* configuration name, configuration resolution */ 421 double e1e2u, e1e2v /* U and V surfaces (if grid size reduction in some straits) */ 422 double ff_f, ff_t /* Coriolis parameter (if not on the sphere) */ 423 \end{Verbatim} 424 425 NEMO can support the local reduction of key strait widths by altering individual values of 426 e2u or e1v at the appropriate locations. This is particularly useful for locations such as 427 Gibraltar or Indonesian Throughflow pinch-points (see \autoref{sec:MISC_strait} for 428 illustrated examples). The key is to reduce the faces of $T$-cell (\ie change the value of 429 the horizontal scale factors at $u$- or $v$-point) but not the volume of the cells. Doing 430 otherwise can lead to numerical instability issues. In normal operation the surface areas 431 are computed from $\texttt{e1u} * \texttt{e2u}$ and $\texttt{e1v} * \texttt{e2v}$ but in 432 cases where a gridsize reduction is required, the unaltered surface areas at $u$ and $v$ 433 grid points (\texttt{e1e2u} and \texttt{e1e2v}, respectively) must be read or pre-computed 434 in \mdl{usrdef\_hgr}. If these arrays are present in the \np{cn\_domcfg} file they are 435 read and the internal computation is suppressed. Versions of \mdl{usrdef\_hgr} which set 436 their own values of \texttt{e1e2u} and \texttt{e1e2v} should set the surface-area 437 computation flag: \texttt{ie1e2u\_v} to a non-zero value to suppress their re-computation. 415 416 \begin{clines} 417 /* Optional: */ 418 int ORCA, ORCA_index /* configuration name, configuration resolution */ 419 double e1e2u, e1e2v /* U and V surfaces (if grid size reduction in some straits) */ 420 double ff_f, ff_t /* Coriolis parameter (if not on the sphere) */ 421 \end{clines} 422 423 \NEMO\ can support the local reduction of key strait widths by 424 altering individual values of e2u or e1v at the appropriate locations. 425 This is particularly useful for locations such as Gibraltar or Indonesian Throughflow pinch-points 426 (see \autoref{sec:MISC_strait} for illustrated examples). 427 The key is to reduce the faces of $T$-cell 428 (\ie\ change the value of the horizontal scale factors at $u$- or $v$-point) but 429 not the volume of the cells. 430 Doing otherwise can lead to numerical instability issues. 431 In normal operation the surface areas are computed from $e1u * e2u$ and $e1v * e2v$ but 432 in cases where a gridsize reduction is required, 433 the unaltered surface areas at $u$ and $v$ grid points 434 (\texttt{e1e2u} and \texttt{e1e2v}, respectively) must be read or pre-computed in \mdl{usrdef\_hgr}. 435 If these arrays are present in the \np{cn_domcfg}{cn\_domcfg} file they are read and 436 the internal computation is suppressed. 437 Versions of \mdl{usrdef\_hgr} which set their own values of \texttt{e1e2u} and \texttt{e1e2v} should 438 set the surface-area computation flag: 439 \texttt{ie1e2u\_v} to a non-zero value to suppress their re-computation. 438 440 439 441 \smallskip 440 Similar logic applies to the other optional fields: \texttt{ff\_f} and \texttt{ff\_t} 441 which can be used to provide the Coriolis parameter at F- and T-points respectively if the 442 mesh is not on a sphere. If present these fields will be read and used and the normal 443 calculation ($2*\Omega*\sin(\varphi)$) suppressed. Versions of \mdl{usrdef\_hgr} which set 444 their own values of \texttt{ff\_f} and \texttt{ff\_t} should set the Coriolis computation 445 flag: \texttt{iff} to a non-zero value to suppress their re-computation. 446 447 Note that longitudes, latitudes, and scale factors at $w$ points are exactly 448 equal to those of $t$ points, thus no specific arrays are defined at $w$ points. 449 450 451 % ================================================================ 452 % Domain: Vertical Grid (domzgr) 453 % ================================================================ 454 \subsection[Vertical grid (\textit{domzgr.F90})] 455 {Vertical grid (\protect\mdl{domzgr})} 442 Similar logic applies to the other optional fields: 443 \texttt{ff\_f} and \texttt{ff\_t} which can be used to 444 provide the Coriolis parameter at F- and T-points respectively if the mesh is not on a sphere. 445 If present these fields will be read and used and 446 the normal calculation ($2 * \Omega * \sin(\varphi)$) suppressed. 447 Versions of \mdl{usrdef\_hgr} which set their own values of \texttt{ff\_f} and \texttt{ff\_t} should 448 set the Coriolis computation flag: 449 \texttt{iff} to a non-zero value to suppress their re-computation. 450 451 Note that longitudes, latitudes, and scale factors at $w$ points are exactly equal to 452 those of $t$ points, thus no specific arrays are defined at $w$ points. 453 454 %% ================================================================================================= 455 \subsection[Vertical grid (\textit{domzgr.F90})]{Vertical grid (\protect\mdl{domzgr})} 456 456 \label{subsec:DOM_zgr} 457 %-----------------------------------------namdom------------------------------------------- 458 \nlst{namdom} 459 %------------------------------------------------------------------------------------------------------------- 460 461 In the vertical, the model mesh is determined by four things: 457 458 \begin{listing} 459 \nlst{namdom} 460 \caption{\forcode{&namdom}} 461 \label{lst:namdom} 462 \end{listing} 463 464 In the vertical, the model mesh is determined by four things: 462 465 \begin{enumerate} 463 \item the bathymetry given in meters; 464 \item the number of levels of the model (\jp{jpk}); 465 \item the analytical transformation $z(i,j,k)$ and the vertical scale factors (derivatives of the transformation); and 466 \item the masking system, \ie the number of wet model levels at each 467 $(i,j)$ location of the horizontal grid. 466 \item the bathymetry given in meters; 467 \item the number of levels of the model (\jp{jpk}); 468 \item the analytical transformation $z(i,j,k)$ and the vertical scale factors 469 (derivatives of the transformation); and 470 \item the masking system, 471 \ie\ the number of wet model levels at each $(i,j)$ location of the horizontal grid. 468 472 \end{enumerate} 469 473 470 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 471 \begin{figure}[!tb] 472 \ begin{center}473 \includegraphics[width=\textwidth]{Fig_z_zps_s_sps}474 \caption{475 \protect\label{fig:z_zps_s_sps}476 The ocean bottom as seen by the model:477 (a) $z$-coordinate with full step,478 (b) $z$-coordinate with partial step,479 (c) $s$-coordinate: terrain following representation,480 (d) hybrid $s-z$ coordinate,481 (e) hybrid $s-z$ coordinate with partial step, and482 ( f) same as (e) but in the non-linear free surface (\protect\np{ln\_linssh}\forcode{ = .false.}).483 Note that the non-linear free surface can be used with any of the 5 coordinates (a) to (e).484 485 \ end{center}474 \begin{figure} 475 \centering 476 \includegraphics[width=0.5\textwidth]{DOM_z_zps_s_sps} 477 \caption[Ocean bottom regarding coordinate systems ($z$, $s$ and hybrid $s-z$)]{ 478 The ocean bottom as seen by the model: 479 \begin{enumerate*}[label=(\textit{\alph*})] 480 \item $z$-coordinate with full step, 481 \item $z$-coordinate with partial step, 482 \item $s$-coordinate: terrain following representation, 483 \item hybrid $s-z$ coordinate, 484 \item hybrid $s-z$ coordinate with partial step, and 485 \item same as (e) but in the non-linear free surface 486 (\protect\np[=.false.]{ln_linssh}{ln\_linssh}). 487 \end{enumerate*} 488 Note that the non-linear free surface can be used with any of the 5 coordinates (a) to (e).} 489 \label{fig:DOM_z_zps_s_sps} 486 490 \end{figure} 487 %>>>>>>>>>>>>>>>>>>>>>>>>>>>>488 491 489 492 The choice of a vertical coordinate is made when setting up the configuration; 490 it is not intended to be an option which can be changed in the middle of an 491 experiment. The one exception to this statement being the choice of linear or 492 non-linear free surface. In v4.0 the linear free surface option is implemented 493 as a special case of the non-linear free surface. This is computationally 494 wasteful since it uses the structures for time-varying 3D metrics for fields 495 that (in the linear free surface case) are fixed. However, the linear 496 free-surface is rarely used and implementing it this way means a single configuration 497 file can support both options. 498 499 By default a non-linear free surface is used (\np{ln\_linssh} set to \forcode{ = 500 .false.} in \ngn{namdom}): the coordinate follow the time-variation of the free 501 surface so that the transformation is time dependent: $z(i,j,k,t)$ 502 (\eg \autoref{fig:z_zps_s_sps}f). When a linear free surface is assumed 503 (\np{ln\_linssh} set to \forcode{ = .true.} in \ngn{namdom}), the vertical 504 coordinates are fixed in time, but the seawater can move up and down across the 505 $z_0$ surface (in other words, the top of the ocean in not a rigid lid). 506 507 Note that settings: \np{ln\_zco}, \np{ln\_zps}, \np{ln\_sco} and \np{ln\_isfcav} mentioned 508 in the following sections appear to be namelist options but they are no longer truly 509 namelist options for NEMO. Their value is written to and read from the domain configuration file 510 and they should be treated as fixed parameters for a particular configuration. They are 511 namelist options for the \forcode{DOMAINcfg} tool that can be used to build the 512 configuration file and serve both to provide a record of the choices made whilst building the 513 configuration and to trigger appropriate code blocks within NEMO. 514 These values should not be altered in the \np{cn\_domcfg} file. 493 it is not intended to be an option which can be changed in the middle of an experiment. 494 The one exception to this statement being the choice of linear or non-linear free surface. 495 In v4.0 the linear free surface option is implemented as 496 a special case of the non-linear free surface. 497 This is computationally wasteful since it uses the structures for time-varying 3D metrics 498 for fields that (in the linear free surface case) are fixed. 499 However, the linear free-surface is rarely used and 500 implementing it this way means a single configuration file can support both options. 501 502 By default a non-linear free surface is used 503 (\np{ln_linssh}{ln\_linssh} set to \forcode{=.false.} in \nam{dom}{dom}): 504 the coordinate follow the time-variation of the free surface so that 505 the transformation is time dependent: $z(i,j,k,t)$ (\eg\ \autoref{fig:DOM_z_zps_s_sps}f). 506 When a linear free surface is assumed 507 (\np{ln_linssh}{ln\_linssh} set to \forcode{=.true.} in \nam{dom}{dom}), 508 the vertical coordinates are fixed in time, but 509 the seawater can move up and down across the $z_0$ surface 510 (in other words, the top of the ocean in not a rigid lid). 511 512 Note that settings: 513 \np{ln_zco}{ln\_zco}, \np{ln_zps}{ln\_zps}, \np{ln_sco}{ln\_sco} and \np{ln_isfcav}{ln\_isfcav} 514 mentioned in the following sections appear to be namelist options but 515 they are no longer truly namelist options for \NEMO. 516 Their value is written to and read from the domain configuration file and 517 they should be treated as fixed parameters for a particular configuration. 518 They are namelist options for the \texttt{DOMAINcfg} tool that can be used to 519 build the configuration file and serve both to provide a record of the choices made whilst 520 building the configuration and to trigger appropriate code blocks within \NEMO. 521 These values should not be altered in the \np{cn_domcfg}{cn\_domcfg} file. 515 522 516 523 \medskip 517 The decision on these choices must be made when the \np{cn \_domcfg} file is constructed.518 Three main choices are offered (\autoref{fig: z_zps_s_sps}a-c):524 The decision on these choices must be made when the \np{cn_domcfg}{cn\_domcfg} file is constructed. 525 Three main choices are offered (\autoref{fig:DOM_z_zps_s_sps}a-c): 519 526 520 527 \begin{itemize} 521 \item $z$-coordinate with full step bathymetry (\np {ln\_zco}\forcode{ = .true.}),522 \item $z$-coordinate with partial step ($zps$) bathymetry (\np {ln\_zps}\forcode{ = .true.}),523 \item Generalized, $s$-coordinate (\np {ln\_sco}\forcode{ = .true.}).528 \item $z$-coordinate with full step bathymetry (\np[=.true.]{ln_zco}{ln\_zco}), 529 \item $z$-coordinate with partial step ($zps$) bathymetry (\np[=.true.]{ln_zps}{ln\_zps}), 530 \item Generalized, $s$-coordinate (\np[=.true.]{ln_sco}{ln\_sco}). 524 531 \end{itemize} 525 532 526 533 Additionally, hybrid combinations of the three main coordinates are available: 527 $s-z$ or $s-zps$ coordinate (\autoref{fig:z_zps_s_sps}d and \autoref{fig:z_zps_s_sps}e). 528 529 A further choice related to vertical coordinate concerns the presence (or not) of ocean 530 cavities beneath ice shelves within the model domain. A setting of \np{ln\_isfcav} as 531 \forcode{.true.} indicates that the domain contains ocean cavities, otherwise the top, 532 wet layer of the ocean will always be at the ocean surface. This option is currently only 533 available for $z$- or $zps$-coordinates. In the latter case, partial steps are also applied 534 at the ocean/ice shelf interface. 535 536 Within the model, the arrays describing the grid point depths and vertical scale factors 537 are three set of three dimensional arrays $(i,j,k)$ defined at \textit{before}, 538 \textit{now} and \textit{after} time step. The time at which they are defined is 539 indicated by a suffix: $\_b$, $\_n$, or $\_a$, respectively. They are updated at each 540 model time step. The initial fixed reference coordinate system is held in variable names 541 with a $\_0$ suffix. When the linear free surface option is used 542 (\np{ln\_linssh}\forcode{ = .true.}), \textit{before}, \textit{now} and \textit{after} 543 arrays are initially set to their reference counterpart and remain fixed. 544 534 $s-z$ or $s-zps$ coordinate (\autoref{fig:DOM_z_zps_s_sps}d and \autoref{fig:DOM_z_zps_s_sps}e). 535 536 A further choice related to vertical coordinate concerns 537 the presence (or not) of ocean cavities beneath ice shelves within the model domain. 538 A setting of \np{ln_isfcav}{ln\_isfcav} as \forcode{.true.} indicates that 539 the domain contains ocean cavities, 540 otherwise the top, wet layer of the ocean will always be at the ocean surface. 541 This option is currently only available for $z$- or $zps$-coordinates. 542 In the latter case, partial steps are also applied at the ocean/ice shelf interface. 543 544 Within the model, 545 the arrays describing the grid point depths and vertical scale factors are 546 three set of three dimensional arrays $(i,j,k)$ defined at 547 \textit{before}, \textit{now} and \textit{after} time step. 548 The time at which they are defined is indicated by a suffix: $\_b$, $\_n$, or $\_a$, respectively. 549 They are updated at each model time step. 550 The initial fixed reference coordinate system is held in variable names with a $\_0$ suffix. 551 When the linear free surface option is used (\np[=.true.]{ln_linssh}{ln\_linssh}), 552 \textit{before}, \textit{now} and \textit{after} arrays are initially set to 553 their reference counterpart and remain fixed. 554 555 %% ================================================================================================= 545 556 \subsubsection{Required fields} 546 557 \label{sec:DOM_zgr_fields} 547 The explicit specification of a range of fields related to the vertical grid are required for the definition of a configuration. These include: 548 549 \begin{Verbatim}[fontsize=\tiny] 558 559 The explicit specification of a range of fields related to the vertical grid are required for 560 the definition of a configuration. 561 These include: 562 563 \begin{clines} 550 564 int ln_zco, ln_zps, ln_sco /* flags for z-coord, z-coord with partial steps and s-coord */ 551 565 int ln_isfcav /* flag for ice shelf cavities */ … … 556 570 /* For reference: */ 557 571 float bathy_metry /* bathymetry used in setting top and bottom levels */ 558 \end{ Verbatim}572 \end{clines} 559 573 560 574 This set of vertical metrics is sufficient to describe the initial depth and thickness of 561 every gridcell in the model regardless of the choice of vertical coordinate. With constant 562 z-levels, e3 metrics will be uniform across each horizontal level. In the partial step 563 case each e3 at the \np{bottom\_level} (and, possibly, \np{top\_level} if ice cavities are 564 present) may vary from its horizontal neighbours. And, in s-coordinates, variations can 565 occur throughout the water column. With the non-linear free-surface, all the coordinates 566 behave more like the s-coordinate in that variations occurr throughout the water column 567 with displacements related to the sea surface height. These variations are typically much 568 smaller than those arising from bottom fitted coordinates. The values for vertical metrics 569 supplied in the domain configuration file can be considered as those arising from a flat 570 sea surface with zero elevation. 571 572 The \np{bottom\_level} and \np{top\_level} 2D arrays define the \np{bottom\_level} and top 573 wet levels in each grid column. Without ice cavities, \np{top\_level} is essentially a land 574 mask (0 on land; 1 everywhere else). With ice cavities, \np{top\_level} determines the 575 first wet point below the overlying ice shelf. 576 577 578 579 % ------------------------------------------------------------------------------------------------------------- 580 % level bathymetry and mask 581 % ------------------------------------------------------------------------------------------------------------- 575 every gridcell in the model regardless of the choice of vertical coordinate. 576 With constant z-levels, e3 metrics will be uniform across each horizontal level. 577 In the partial step case each e3 at the \jp{bottom\_level} 578 (and, possibly, \jp{top\_level} if ice cavities are present) 579 may vary from its horizontal neighbours. 580 And, in s-coordinates, variations can occur throughout the water column. 581 With the non-linear free-surface, all the coordinates behave more like the s-coordinate in that 582 variations occur throughout the water column with displacements related to the sea surface height. 583 These variations are typically much smaller than those arising from bottom fitted coordinates. 584 The values for vertical metrics supplied in the domain configuration file can be considered as 585 those arising from a flat sea surface with zero elevation. 586 587 The \jp{bottom\_level} and \jp{top\_level} 2D arrays define 588 the \jp{bottom\_level} and top wet levels in each grid column. 589 Without ice cavities, \jp{top\_level} is essentially a land mask (0 on land; 1 everywhere else). 590 With ice cavities, \jp{top\_level} determines the first wet point below the overlying ice shelf. 591 592 %% ================================================================================================= 582 593 \subsubsection{Level bathymetry and mask} 583 594 \label{subsec:DOM_msk} 584 595 585 586 From \np{top\_level} and \np{bottom\_level} fields, the mask fields are defined as follows: 587 \begin{alignat*}{2} 588 tmask(i,j,k) &= & & 589 \begin{cases} 590 0 &\text{if $ k < top\_level(i,j)$} \\ 591 1 &\text{if $bottom\_level(i,j) \leq k \leq top\_level(i,j)$} \\ 592 0 &\text{if $ k > bottom\_level(i,j)$} 593 \end{cases} 594 \\ 595 umask(i,j,k) &= & &tmask(i,j,k) * tmask(i + 1,j, k) \\ 596 vmask(i,j,k) &= & &tmask(i,j,k) * tmask(i ,j + 1,k) \\ 597 fmask(i,j,k) &= & &tmask(i,j,k) * tmask(i + 1,j, k) \\ 598 & &* &tmask(i,j,k) * tmask(i + 1,j, k) \\ 599 wmask(i,j,k) &= & &tmask(i,j,k) * tmask(i ,j,k - 1) \\ 600 \text{with~} wmask(i,j,1) &= & &tmask(i,j,1) 601 \end{alignat*} 596 From \jp{top\_level} and \jp{bottom\_level} fields, the mask fields are defined as follows: 597 \begin{align*} 598 tmask(i,j,k) &= 599 \begin{cases} 600 0 &\text{if $ k < top\_level(i,j)$} \\ 601 1 &\text{if $ bottom\_level(i,j) \leq k \leq top\_level(i,j)$} \\ 602 0 &\text{if $k > bottom\_level(i,j) $} 603 \end{cases} \\ 604 umask(i,j,k) &= tmask(i,j,k) * tmask(i + 1,j, k) \\ 605 vmask(i,j,k) &= tmask(i,j,k) * tmask(i ,j + 1,k) \\ 606 fmask(i,j,k) &= tmask(i,j,k) * tmask(i + 1,j, k) * tmask(i,j,k) * tmask(i + 1,j, k) \\ 607 wmask(i,j,k) &= tmask(i,j,k) * tmask(i ,j,k - 1) \\ 608 \text{with~} wmask(i,j,1) &= tmask(i,j,1) 609 \end{align*} 602 610 603 611 Note that, without ice shelves cavities, 604 masks at $t-$ and $w-$points are identical with the numerical indexing used (\autoref{subsec:DOM_Num_Index}). 605 Nevertheless, $wmask$ are required with ocean cavities to deal with the top boundary (ice shelf/ocean interface) 612 masks at $t-$ and $w-$points are identical with the numerical indexing used 613 (\autoref{subsec:DOM_Num_Index}). 614 Nevertheless, 615 $wmask$ are required with ocean cavities to deal with the top boundary (ice shelf/ocean interface) 606 616 exactly in the same way as for the bottom boundary. 607 617 … … 612 622 %% (see \autoref{fig:LBC_jperio}). 613 623 614 615 % -------------------------------------------------------------------------------------------------616 % Closed seas 617 %------------------------------------------------------------------------------------------------- 618 \subsection{Closed seas} \label{subsec:DOM_closea} 619 620 When a global ocean is coupled to an atmospheric model it is better to represent all large 621 water bodies (\eg great lakes, Caspian sea...) even if the model resolution does not allow 622 their communication with the rest of the ocean. This is unnecessary when the ocean is 623 forced by fixed atmospheric conditions, so these seas can be removed from the ocean 624 domain. The user has the option to set the bathymetry in closed seas to zero (see 625 \autoref{sec:MISC_closea}) and to optionally decide on the fate of any freshwater 626 imbalance over the area. The options are explained in \autoref{sec:MISC_closea} but it 627 should be noted here that a successful use of these options requires appropriate mask 628 fields to be present in the domain configuration file. Among the possibilities are: 629 630 \begin{Verbatim}[fontsize=\tiny] 631 int closea_mask /* non-zero values in closed sea areas for optional masking */ 632 int closea_mask_rnf /* non-zero values in closed sea areas with runoff locations (precip only) */ 633 int closea_mask_emp /* non-zero values in closed sea areas with runoff locations (total emp)*/634 \end{Verbatim} 635 636 % ------------------------------------------------------------------------------------------------------------- 637 % Grid files 638 % -------------------------------------------------------------------------------------------------------------624 % Closed seas 625 %% ================================================================================================= 626 \subsection{Closed seas} 627 \label{subsec:DOM_closea} 628 629 When a global ocean is coupled to an atmospheric model it is better to 630 represent all large water bodies (\eg\ Great Lakes, Caspian sea, \dots) even if 631 the model resolution does not allow their communication with the rest of the ocean. 632 This is unnecessary when the ocean is forced by fixed atmospheric conditions, 633 so these seas can be removed from the ocean domain. 634 The user has the option to 635 set the bathymetry in closed seas to zero (see \autoref{sec:MISC_closea}) and to 636 optionally decide on the fate of any freshwater imbalance over the area. 637 The options are explained in \autoref{sec:MISC_closea} but 638 it should be noted here that a successful use of these options requires 639 appropriate mask fields to be present in the domain configuration file. 640 Among the possibilities are: 641 642 \begin{clines} 643 int closea_mask /* non-zero values in closed sea areas for optional masking */ 644 int closea_mask_rnf /* non-zero values in closed sea areas with runoff locations (precip only) */ 645 int closea_mask_emp /* non-zero values in closed sea areas with runoff locations (total emp) */ 646 \end{clines} 647 648 %% ================================================================================================= 639 649 \subsection{Output grid files} 640 650 \label{subsec:DOM_meshmask} 641 651 642 \nlst{namcfg} 643 644 Most of the arrays relating to a particular ocean model configuration dicussed in this 645 chapter (grid-point position, scale factors) can be saved in a file if namelist parameter 646 \np{ln\_write\_cfg} (namelist \ngn{namcfg}) is set to \forcode{.true.}; the output 647 filename is set thorugh parameter \np{cn\_domcfg\_out}. This is only really useful 648 ifthe fields are computed in subroutines \mdl{usrdef\_hgr} or \mdl{usrdef\_zgr} and652 Most of the arrays relating to a particular ocean model configuration discussed in this chapter 653 (grid-point position, scale factors) can be saved in a file if 654 namelist parameter \np{ln_write_cfg}{ln\_write\_cfg} (namelist \nam{cfg}{cfg}) is set to 655 \forcode{.true.}; 656 the output filename is set through parameter \np{cn_domcfg_out}{cn\_domcfg\_out}. 657 This is only really useful if 658 the fields are computed in subroutines \mdl{usrdef\_hgr} or \mdl{usrdef\_zgr} and 649 659 checking or confirmation is required. 650 660 651 \nlst{namdom}652 653 661 Alternatively, all the arrays relating to a particular ocean model configuration 654 (grid-point position, scale factors, depths and masks) can be saved in a file called 655 \texttt{mesh\_mask} if namelist parameter \np{ln\_meshmask} (namelist \ngn{namdom}) is set 656 to \forcode{.true.}. This file contains additional fields that can be useful for 657 post-processing applications 658 659 % ================================================================ 660 % Domain: Initial State (dtatsd & istate) 661 % ================================================================ 662 \section[Initial state (\textit{istate.F90} and \textit{dtatsd.F90})] 663 {Initial state (\protect\mdl{istate} and \protect\mdl{dtatsd})} 664 \label{sec:DTA_tsd} 665 %-----------------------------------------namtsd------------------------------------------- 666 \nlst{namtsd} 667 %------------------------------------------------------------------------------------------ 668 669 Basic initial state options are defined in \ngn{namtsd}. By default, the ocean starts 670 from rest (the velocity field is set to zero) and the initialization of temperature and 671 salinity fields is controlled through the \np{ln\_tsd\_init} namelist parameter. 662 (grid-point position, scale factors, depths and masks) can be saved in 663 a file called \texttt{mesh\_mask} if 664 namelist parameter \np{ln_meshmask}{ln\_meshmask} (namelist \nam{dom}{dom}) is set to 665 \forcode{.true.}. 666 This file contains additional fields that can be useful for post-processing applications. 667 668 %% ================================================================================================= 669 \section[Initial state (\textit{istate.F90} and \textit{dtatsd.F90})]{Initial state (\protect\mdl{istate} and \protect\mdl{dtatsd})} 670 \label{sec:DOM_DTA_tsd} 671 672 \begin{listing} 673 \nlst{namtsd} 674 \caption{\forcode{&namtsd}} 675 \label{lst:namtsd} 676 \end{listing} 677 678 Basic initial state options are defined in \nam{tsd}{tsd}. 679 By default, the ocean starts from rest (the velocity field is set to zero) and 680 the initialization of temperature and salinity fields is controlled through the \np{ln_tsd_init}{ln\_tsd\_init} namelist parameter. 672 681 673 682 \begin{description} 674 \item[\np{ln\_tsd\_init}\forcode{= .true.}] 675 Use T and S input files that can be given on the model grid itself or on their native 676 input data grids. In the latter case, the data will be interpolated on-the-fly both in 677 the horizontal and the vertical to the model grid (see \autoref{subsec:SBC_iof}). The 678 information relating to the input files are specified in the \np{sn\_tem} and 679 \np{sn\_sal} structures. The computation is done in the \mdl{dtatsd} module. 680 \item[\np{ln\_tsd\_init}\forcode{= .false.}] 681 Initial values for T and S are set via a user supplied \rou{usr\_def\_istate} routine 682 contained in \mdl{userdef\_istate}. The default version sets horizontally uniform T and 683 profiles as used in the GYRE configuration (see \autoref{sec:CFG_gyre}). 683 \item [{\np[=.true.]{ln_tsd_init}{ln\_tsd\_init}}] Use T and S input files that can be given on 684 the model grid itself or on their native input data grids. 685 In the latter case, 686 the data will be interpolated on-the-fly both in the horizontal and the vertical to the model grid 687 (see \autoref{subsec:SBC_iof}). 688 The information relating to the input files are specified in 689 the \np{sn_tem}{sn\_tem} and \np{sn_sal}{sn\_sal} structures. 690 The computation is done in the \mdl{dtatsd} module. 691 \item [{\np[=.false.]{ln_tsd_init}{ln\_tsd\_init}}] Initial values for T and S are set via 692 a user supplied \rou{usr\_def\_istate} routine contained in \mdl{userdef\_istate}. 693 The default version sets horizontally uniform T and profiles as used in the GYRE configuration 694 (see \autoref{sec:CFGS_gyre}). 684 695 \end{description} 685 696 686 \biblio 687 688 \pindex 697 \subinc{\input{../../global/epilogue}} 689 698 690 699 \end{document}
Note: See TracChangeset
for help on using the changeset viewer.