1 | \documentclass[../main/NEMO_manual]{subfiles} |
---|
2 | |
---|
3 | \begin{document} |
---|
4 | |
---|
5 | \chapter{Space Domain (DOM)} |
---|
6 | \label{chap:DOM} |
---|
7 | |
---|
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 | \chaptertoc |
---|
18 | |
---|
19 | \paragraph{Changes record} ~\\ |
---|
20 | |
---|
21 | {\footnotesize |
---|
22 | \begin{tabularx}{0.8\textwidth}{l||X|X} |
---|
23 | Release & |
---|
24 | Author(s) & |
---|
25 | Modifications \\ |
---|
26 | \hline |
---|
27 | {\em 4.0 } & |
---|
28 | {\em Simon M\"{u}ller \& Andrew Coward \newline \newline |
---|
29 | Simona Flavoni and Tim Graham } & |
---|
30 | {\em Compatibility changes: many options moved to external domain configuration tools |
---|
31 | (see \autoref{apdx:DOMCFG}). \newline |
---|
32 | Updates } \\ |
---|
33 | {\em 3.6 } & |
---|
34 | {\em Rachid Benshila, Christian \'{E}th\'{e}, Pierre Mathiot and Gurvan Madec } & |
---|
35 | {\em Updates } \\ |
---|
36 | {\em $\leq$ 3.4 } & |
---|
37 | {\em Gurvan Madec and S\'{e}bastien Masson } & |
---|
38 | {\em First version } |
---|
39 | \end{tabularx} |
---|
40 | } |
---|
41 | |
---|
42 | \clearpage |
---|
43 | |
---|
44 | Having defined the continuous equations in \autoref{chap:MB} and |
---|
45 | chosen a time discretisation \autoref{chap:TD}, |
---|
46 | we need to choose a grid for spatial discretisation and related numerical algorithms. |
---|
47 | In the present chapter, we provide a general description of the staggered grid used in \NEMO, |
---|
48 | and other relevant information about the DOM (DOMain) source code modules. |
---|
49 | |
---|
50 | %% ================================================================================================= |
---|
51 | \section{Fundamentals of the discretisation} |
---|
52 | \label{sec:DOM_basics} |
---|
53 | |
---|
54 | %% ================================================================================================= |
---|
55 | \subsection{Arrangement of variables} |
---|
56 | \label{subsec:DOM_cell} |
---|
57 | |
---|
58 | \begin{figure} |
---|
59 | \centering |
---|
60 | \includegraphics[width=0.33\textwidth]{DOM_cell} |
---|
61 | \caption[Arrangement of variables in the unit cell of space domain]{ |
---|
62 | Arrangement of variables in the unit cell of space domain. |
---|
63 | $t$ indicates scalar points where |
---|
64 | temperature, salinity, density, pressure and horizontal divergence are defined. |
---|
65 | $(u,v,w)$ indicates vector points, and $f$ indicates vorticity points where |
---|
66 | both relative and planetary vorticities are defined.} |
---|
67 | \label{fig:DOM_cell} |
---|
68 | \end{figure} |
---|
69 | |
---|
70 | The numerical techniques used to solve the Primitive Equations in this model are based on |
---|
71 | the traditional, centred second-order finite difference approximation. |
---|
72 | Special attention has been given to the homogeneity of the solution in the three spatial directions. |
---|
73 | The arrangement of variables is the same in all directions. |
---|
74 | It consists of cells centred on scalar points ($t$, $S$, $p$, $\rho$) with |
---|
75 | vector points $(u, v, w)$ defined in the centre of each face of the cells (\autoref{fig:DOM_cell}). |
---|
76 | This is the generalisation to three dimensions of the well-known ``C'' grid in |
---|
77 | Arakawa's classification \citep{mesinger.arakawa_bk76}. |
---|
78 | The relative and planetary vorticity, $\zeta$ and $f$, are defined in the centre of each |
---|
79 | vertical edge and the barotropic stream function $\psi$ is defined at horizontal points overlying |
---|
80 | the $\zeta$ and $f$-points. |
---|
81 | |
---|
82 | The ocean mesh (\ie\ the position of all the scalar and vector points) is defined by |
---|
83 | the transformation that gives $(\lambda,\varphi,z)$ as a function of $(i,j,k)$. |
---|
84 | The grid-points are located at integer or integer and a half value of $(i,j,k)$ as indicated on |
---|
85 | \autoref{tab:DOM_cell}. |
---|
86 | In all the following, |
---|
87 | subscripts $u$, $v$, $w$, $f$, $uw$, $vw$ or $fw$ indicate the position of the grid-point where |
---|
88 | the scale factors are defined. |
---|
89 | Each scale factor is defined as the local analytical value provided by \autoref{eq:MB_scale_factors}. |
---|
90 | As a result, the mesh on which partial derivatives $\pd[]{\lambda}$, $\pd[]{\varphi}$ and |
---|
91 | $\pd[]{z}$ are evaluated is a uniform mesh with a grid size of unity. |
---|
92 | Discrete partial derivatives are formulated by |
---|
93 | the traditional, centred second order finite difference approximation while |
---|
94 | the scale factors are chosen equal to their local analytical value. |
---|
95 | An important point here is that the partial derivative of the scale factors must be evaluated by |
---|
96 | centred finite difference approximation, not from their analytical expression. |
---|
97 | This preserves the symmetry of the discrete set of equations and |
---|
98 | therefore satisfies many of the continuous properties (see \autoref{apdx:INVARIANTS}). |
---|
99 | A similar, related remark can be made about the domain size: |
---|
100 | when needed, an area, volume, or the total ocean depth must be evaluated as |
---|
101 | the product or sum of the relevant scale factors (see \autoref{eq:DOM_bar} in the next section). |
---|
102 | |
---|
103 | \begin{table} |
---|
104 | \centering |
---|
105 | \begin{tabular}{|l|l|l|l|} |
---|
106 | \hline |
---|
107 | t & $i $ & $j $ & $k $ \\ |
---|
108 | \hline |
---|
109 | u & $i + 1/2$ & $j $ & $k $ \\ |
---|
110 | \hline |
---|
111 | v & $i $ & $j + 1/2$ & $k $ \\ |
---|
112 | \hline |
---|
113 | w & $i $ & $j $ & $k + 1/2$ \\ |
---|
114 | \hline |
---|
115 | f & $i + 1/2$ & $j + 1/2$ & $k $ \\ |
---|
116 | \hline |
---|
117 | uw & $i + 1/2$ & $j $ & $k + 1/2$ \\ |
---|
118 | \hline |
---|
119 | vw & $i $ & $j + 1/2$ & $k + 1/2$ \\ |
---|
120 | \hline |
---|
121 | fw & $i + 1/2$ & $j + 1/2$ & $k + 1/2$ \\ |
---|
122 | \hline |
---|
123 | \end{tabular} |
---|
124 | \caption[Location of grid-points]{ |
---|
125 | Location of grid-points as a function of integer or |
---|
126 | integer and a half value of the column, line or level. |
---|
127 | This indexing is only used for the writing of the semi-discrete equations. |
---|
128 | In the code, the indexing uses integer values only and |
---|
129 | is positive downwards in the vertical with $k=1$ at the surface. |
---|
130 | (see \autoref{subsec:DOM_Num_Index})} |
---|
131 | \label{tab:DOM_cell} |
---|
132 | \end{table} |
---|
133 | |
---|
134 | Note that the definition of the scale factors |
---|
135 | (\ie\ as the analytical first derivative of the transformation that |
---|
136 | results in $(\lambda,\varphi,z)$ as a function of $(i,j,k)$) |
---|
137 | is specific to the \NEMO\ model \citep{marti.madec.ea_JGR92}. |
---|
138 | As an example, a scale factor in the $i$ direction is defined locally at a $t$-point, |
---|
139 | whereas many other models on a C grid choose to define such a scale factor as |
---|
140 | the distance between the $u$-points on each side of the $t$-point. |
---|
141 | Relying on an analytical transformation has two advantages: |
---|
142 | firstly, there is no ambiguity in the scale factors appearing in the discrete equations, |
---|
143 | since they are first introduced in the continuous equations; |
---|
144 | secondly, analytical transformations encourage good practice by |
---|
145 | the definition of smoothly varying grids |
---|
146 | (rather than allowing the user to set arbitrary jumps in thickness between adjacent layers) |
---|
147 | \citep{treguier.dukowicz.ea_JGR96}. |
---|
148 | An example of the effect of such a choice is shown in \autoref{fig:DOM_zgr_e3}. |
---|
149 | \begin{figure} |
---|
150 | \centering |
---|
151 | \includegraphics[width=0.5\textwidth]{DOM_zgr_e3} |
---|
152 | \caption[Comparison of grid-point position, vertical grid-size and scale factors]{ |
---|
153 | Comparison of (a) traditional definitions of grid-point position and grid-size in the vertical, |
---|
154 | and (b) analytically derived grid-point position and scale factors. |
---|
155 | For both grids here, the same $w$-point depth has been chosen but |
---|
156 | in (a) the $t$-points are set half way between $w$-points while |
---|
157 | in (b) they are defined from an analytical function: |
---|
158 | $z(k) = 5 \, (k - 1/2)^3 - 45 \, (k - 1/2)^2 + 140 \, (k - 1/2) - 150$. |
---|
159 | Note the resulting difference between the value of the grid-size $\Delta_k$ and |
---|
160 | those of the scale factor $e_k$.} |
---|
161 | \label{fig:DOM_zgr_e3} |
---|
162 | \end{figure} |
---|
163 | |
---|
164 | %% ================================================================================================= |
---|
165 | \subsection{Discrete operators} |
---|
166 | \label{subsec:DOM_operators} |
---|
167 | |
---|
168 | Given the values of a variable $q$ at adjacent points, |
---|
169 | the differencing and averaging operators at the midpoint between them are: |
---|
170 | \begin{alignat*}{2} |
---|
171 | % \label{eq:DOM_di_mi} |
---|
172 | \delta_i [q] &= & &q (i + 1/2) - q (i - 1/2) \\ |
---|
173 | \overline q^{\, i} &= &\big\{ &q (i + 1/2) + q (i - 1/2) \big\} / 2 |
---|
174 | \end{alignat*} |
---|
175 | |
---|
176 | Similar operators are defined with respect to $i + 1/2$, $j$, $j + 1/2$, $k$, and $k + 1/2$. |
---|
177 | Following \autoref{eq:MB_grad} and \autoref{eq:MB_lap}, |
---|
178 | the gradient of a variable $q$ defined at a $t$-point has |
---|
179 | its three components defined at $u$-, $v$- and $w$-points while |
---|
180 | its Laplacian is defined at the $t$-point. |
---|
181 | These operators have the following discrete forms in the curvilinear $s$-coordinates system: |
---|
182 | \begin{gather*} |
---|
183 | % \label{eq:DOM_grad} |
---|
184 | \nabla q \equiv \frac{1}{e_{1u}} \delta_{i + 1/2} [q] \; \, \vect i |
---|
185 | + \frac{1}{e_{2v}} \delta_{j + 1/2} [q] \; \, \vect j |
---|
186 | + \frac{1}{e_{3w}} \delta_{k + 1/2} [q] \; \, \vect k \\ |
---|
187 | % \label{eq:DOM_lap} |
---|
188 | \Delta q \equiv \frac{1}{e_{1t} \, e_{2t} \, e_{3t}} |
---|
189 | \; \lt[ \delta_i \lt( \frac{e_{2u} \, e_{3u}}{e_{1u}} \; \delta_{i + 1/2} [q] \rt) |
---|
190 | + \delta_j \lt( \frac{e_{1v} \, e_{3v}}{e_{2v}} \; \delta_{j + 1/2} [q] \rt) \; \rt] |
---|
191 | + \frac{1}{e_{3t}} |
---|
192 | \delta_k \lt[ \frac{1 }{e_{3w}} \; \delta_{k + 1/2} [q] \rt] |
---|
193 | \end{gather*} |
---|
194 | |
---|
195 | Following \autoref{eq:MB_curl} and \autoref{eq:MB_div}, |
---|
196 | a vector $\vect A = (a_1,a_2,a_3)$ defined at vector points $(u,v,w)$ has |
---|
197 | its three curl components defined at $vw$-, $uw$, and $f$-points, and |
---|
198 | its divergence defined at $t$-points: |
---|
199 | \begin{multline*} |
---|
200 | % \label{eq:DOM_curl} |
---|
201 | \nabla \times \vect A \equiv \frac{1}{e_{2v} \, e_{3vw}} |
---|
202 | \Big[ \delta_{j + 1/2} (e_{3w} \, a_3) |
---|
203 | - \delta_{k + 1/2} (e_{2v} \, a_2) \Big] \vect i \\ |
---|
204 | + \frac{1}{e_{2u} \, e_{3uw}} |
---|
205 | \Big[ \delta_{k + 1/2} (e_{1u} \, a_1) |
---|
206 | - \delta_{i + 1/2} (e_{3w} \, a_3) \Big] \vect j \\ |
---|
207 | + \frac{1}{e_{1f} \, e_{2f}} |
---|
208 | \Big[ \delta_{i + 1/2} (e_{2v} \, a_2) |
---|
209 | - \delta_{j + 1/2} (e_{1u} \, a_1) \Big] \vect k |
---|
210 | \end{multline*} |
---|
211 | \[ |
---|
212 | % \label{eq:DOM_div} |
---|
213 | \nabla \cdot \vect A \equiv \frac{1}{e_{1t} \, e_{2t} \, e_{3t}} |
---|
214 | \Big[ \delta_i (e_{2u} \, e_{3u} \, a_1) + \delta_j (e_{1v} \, e_{3v} \, a_2) \Big] |
---|
215 | + \frac{1}{e_{3t}} \delta_k (a_3) |
---|
216 | \] |
---|
217 | |
---|
218 | The vertical average over the whole water column is denoted by an overbar and |
---|
219 | is for a masked field $q$ (\ie\ a quantity that is equal to zero inside solid areas): |
---|
220 | \begin{equation} |
---|
221 | \label{eq:DOM_bar} |
---|
222 | \bar q = \frac{1}{H} \int_{k^b}^{k^o} q \; e_{3q} \, dk \equiv \frac{1}{H_q} \sum \limits_k q \; e_{3q} |
---|
223 | \end{equation} |
---|
224 | where $H_q$ is the ocean depth, which is the masked sum of the vertical scale factors at $q$ points, |
---|
225 | $k^b$ and $k^o$ are the bottom and surface $k$-indices, |
---|
226 | and the symbol $\sum \limits_k$ refers to a summation over all grid points of the same type in |
---|
227 | the direction indicated by the subscript (here $k$). |
---|
228 | |
---|
229 | In continuous form, the following properties are satisfied: |
---|
230 | \begin{gather} |
---|
231 | \label{eq:DOM_curl_grad} |
---|
232 | \nabla \times \nabla q = \vect 0 \\ |
---|
233 | \label{eq:DOM_div_curl} |
---|
234 | \nabla \cdot (\nabla \times \vect A) = 0 |
---|
235 | \end{gather} |
---|
236 | |
---|
237 | It is straightforward to demonstrate that these properties are verified locally in discrete form as |
---|
238 | soon as the scalar $q$ is taken at $t$-points and the vector $\vect A$ has its components defined at |
---|
239 | vector points $(u,v,w)$. |
---|
240 | |
---|
241 | Let $a$ and $b$ be two fields defined on the mesh, with a value of zero inside continental areas. |
---|
242 | It can be shown that the differencing operators ($\delta_i$, $\delta_j$ and |
---|
243 | $\delta_k$) are skew-symmetric linear operators, |
---|
244 | and further that the averaging operators ($\overline{\cdots}^{\, i}$, $\overline{\cdots}^{\, j}$ and |
---|
245 | $\overline{\cdots}^{\, k}$) are symmetric linear operators, \ie |
---|
246 | \begin{alignat}{5} |
---|
247 | \label{eq:DOM_di_adj} |
---|
248 | &\sum \limits_i a_i \; \delta_i [b] &\equiv &- &&\sum \limits_i \delta _{ i + 1/2} [a] &b_{i + 1/2} \\ |
---|
249 | \label{eq:DOM_mi_adj} |
---|
250 | &\sum \limits_i a_i \; \overline b^{\, i} &\equiv & &&\sum \limits_i \overline a ^{\, i + 1/2} &b_{i + 1/2} |
---|
251 | \end{alignat} |
---|
252 | |
---|
253 | In other words, |
---|
254 | the adjoint of the differencing and averaging operators are $\delta_i^* = \delta_{i + 1/2}$ and |
---|
255 | $(\overline{\cdots}^{\, i})^* = \overline{\cdots}^{\, i + 1/2}$, respectively. |
---|
256 | These two properties will be used extensively in the \autoref{apdx:INVARIANTS} to |
---|
257 | demonstrate integral conservative properties of the discrete formulation chosen. |
---|
258 | |
---|
259 | %% ================================================================================================= |
---|
260 | \subsection{Numerical indexing} |
---|
261 | \label{subsec:DOM_Num_Index} |
---|
262 | |
---|
263 | \begin{figure} |
---|
264 | \centering |
---|
265 | \includegraphics[width=0.33\textwidth]{DOM_index_hor} |
---|
266 | \caption[Horizontal integer indexing]{ |
---|
267 | Horizontal integer indexing used in the \fortran\ code. |
---|
268 | The dashed area indicates the cell in which |
---|
269 | variables contained in arrays have the same $i$- and $j$-indices} |
---|
270 | \label{fig:DOM_index_hor} |
---|
271 | \end{figure} |
---|
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 |
---|
275 | the use of integer values for $t$-points only while |
---|
276 | all the other points involve integer and a half values. |
---|
277 | Therefore, a specific integer indexing has been defined for points other than $t$-points |
---|
278 | (\ie\ velocity and vorticity grid-points). |
---|
279 | Furthermore, the direction of the vertical indexing has been reversed and |
---|
280 | the surface level set at $k = 1$. |
---|
281 | |
---|
282 | %% ================================================================================================= |
---|
283 | \subsubsection{Horizontal indexing} |
---|
284 | \label{subsec:DOM_Num_Index_hor} |
---|
285 | |
---|
286 | The indexing in the horizontal plane has been chosen as shown in \autoref{fig:DOM_index_hor}. |
---|
287 | For an increasing $i$ index ($j$ index), |
---|
288 | the $t$-point and the eastward $u$-point (northward $v$-point) have the same index |
---|
289 | (see the dashed area in \autoref{fig:DOM_index_hor}). |
---|
290 | A $t$-point and its nearest north-east $f$-point have the same $i$-and $j$-indices. |
---|
291 | |
---|
292 | %% ================================================================================================= |
---|
293 | \subsubsection{Vertical indexing} |
---|
294 | \label{subsec:DOM_Num_Index_vertical} |
---|
295 | |
---|
296 | In the vertical, the chosen indexing requires special attention since |
---|
297 | the direction of the $k$-axis in the \fortran\ code is the reverse of |
---|
298 | that used in the semi-discrete equations and given in \autoref{subsec:DOM_cell}. |
---|
299 | The sea surface corresponds to the $w$-level $k = 1$, |
---|
300 | which is the same index as the $t$-level just below (\autoref{fig:DOM_index_vert}). |
---|
301 | The last $w$-level ($k = jpk$) either corresponds to or is below the ocean floor while |
---|
302 | the last $t$-level is always outside the ocean domain (\autoref{fig:DOM_index_vert}). |
---|
303 | Note that a $w$-point and the directly underlaying $t$-point have a common $k$ index |
---|
304 | (\ie\ $t$-points and their nearest $w$-point neighbour in negative index direction), |
---|
305 | in contrast to the indexing on the horizontal plane where |
---|
306 | the $t$-point has the same index as the nearest velocity points in |
---|
307 | the positive direction of the respective horizontal axis index |
---|
308 | (compare the dashed area in \autoref{fig:DOM_index_hor} and \autoref{fig:DOM_index_vert}). |
---|
309 | Since the scale factors are chosen to be strictly positive, |
---|
310 | a \textit{minus sign} is included in the \fortran\ implementations of |
---|
311 | \textit{all the vertical derivatives} of the discrete equations given in this manual in order to |
---|
312 | accommodate the opposing vertical index directions in implementation and documentation. |
---|
313 | |
---|
314 | \begin{figure} |
---|
315 | \centering |
---|
316 | \includegraphics[width=0.33\textwidth]{DOM_index_vert} |
---|
317 | \caption[Vertical integer indexing]{ |
---|
318 | Vertical integer indexing used in the \fortran\ code. |
---|
319 | Note that the $k$-axis is oriented downward. |
---|
320 | The dashed area indicates the cell in which |
---|
321 | variables contained in arrays have a common $k$-index.} |
---|
322 | \label{fig:DOM_index_vert} |
---|
323 | \end{figure} |
---|
324 | |
---|
325 | %% ================================================================================================= |
---|
326 | \section{Spatial domain configuration} |
---|
327 | \label{subsec:DOM_config} |
---|
328 | |
---|
329 | Two typical methods are available to specify the spatial domain configuration; |
---|
330 | they can be selected using parameter \np{ln_read_cfg}{ln\_read\_cfg} parameter in |
---|
331 | namelist \nam{cfg}{cfg}. |
---|
332 | |
---|
333 | If \np{ln_read_cfg}{ln\_read\_cfg} is set to \forcode{.true.}, |
---|
334 | the domain-specific parameters and fields are read from a NetCDF input file, |
---|
335 | whose name (without its .nc suffix) can be specified as |
---|
336 | the value of the \np{cn_domcfg}{cn\_domcfg} parameter in namelist \nam{cfg}{cfg}. |
---|
337 | |
---|
338 | If \np{ln_read_cfg}{ln\_read\_cfg} is set to \forcode{.false.}, |
---|
339 | the domain-specific parameters and fields can be provided (\eg\ analytically computed) by |
---|
340 | subroutines \mdl{usrdef\_hgr} and \mdl{usrdef\_zgr}. |
---|
341 | These subroutines can be supplied in the \path{MY_SRC} directory of the configuration, |
---|
342 | and default versions that configure the spatial domain for the GYRE reference configuration are |
---|
343 | present in the \path{./src/OCE/USR} directory. |
---|
344 | |
---|
345 | In version 4.0 there are no longer any options for reading complex bathymetries and |
---|
346 | performing a vertical discretisation at run-time. |
---|
347 | Whilst it is occasionally convenient to have a common bathymetry file and, for example, |
---|
348 | to run similar models with and without partial bottom boxes and/or sigma-coordinates, |
---|
349 | supporting such choices leads to overly complex code. |
---|
350 | Worse still is the difficulty of ensuring the model configurations intended to be identical are |
---|
351 | indeed so when the model domain itself can be altered by runtime selections. |
---|
352 | The code previously used to perform vertical discretisation has been incorporated into |
---|
353 | an external tool (\path{./tools/DOMAINcfg}) which is briefly described in \autoref{apdx:DOMCFG}. |
---|
354 | |
---|
355 | The next subsections summarise the parameter and fields related to |
---|
356 | the configuration of the whole model domain. |
---|
357 | These represent the minimum information that must be provided either via |
---|
358 | the \np{cn_domcfg}{cn\_domcfg} file or |
---|
359 | set by code inserted into user-supplied versions of the \texttt{usrdef\_*} subroutines. |
---|
360 | The requirements are presented in three sections: |
---|
361 | the domain size (\autoref{subsec:DOM_size}), the horizontal mesh (\autoref{subsec:DOM_hgr}), |
---|
362 | and the vertical grid (\autoref{subsec:DOM_zgr}). |
---|
363 | |
---|
364 | %% ================================================================================================= |
---|
365 | \subsection{Domain size} |
---|
366 | \label{subsec:DOM_size} |
---|
367 | |
---|
368 | The total size of the computational domain is set by the parameters \texttt{jpiglo}, \texttt{jpjglo} and |
---|
369 | \texttt{jpkglo} for the $i$, $j$ and $k$ directions, respectively. |
---|
370 | Note, that the variables \texttt{jpi} and \texttt{jpj} refer to |
---|
371 | the size of each processor subdomain when the code is run in parallel using domain decomposition |
---|
372 | (\key{mpp\_mpi} defined, see \autoref{sec:LBC_mpp}). |
---|
373 | |
---|
374 | The name of the configuration is set through parameter \np{cn_cfg}{cn\_cfg}, |
---|
375 | and the nominal resolution through parameter \np{nn_cfg}{nn\_cfg} |
---|
376 | (unless in the input file both of variables \texttt{ORCA} and \texttt{ORCA\_index} are present, |
---|
377 | in which case \np{cn_cfg}{cn\_cfg} and \np{nn_cfg}{nn\_cfg} are set from these values accordingly). |
---|
378 | |
---|
379 | The global lateral boundary condition type is selected from 8 options using parameters \texttt{l\_Iperio}, \texttt{l\_Jperio}, \texttt{l\_NFold} and \texttt{c\_NFtype}. |
---|
380 | See \autoref{sec:LBC_jperio} for details on the available options and |
---|
381 | the corresponding values for \texttt{l\_Iperio}, \texttt{l\_Jperio}, \texttt{l\_NFold} and \texttt{c\_NFtype}. |
---|
382 | |
---|
383 | %% ================================================================================================= |
---|
384 | \subsection[Horizontal grid mesh (\textit{domhgr.F90}]{Horizontal grid mesh (\protect\mdl{domhgr})} |
---|
385 | \label{subsec:DOM_hgr} |
---|
386 | |
---|
387 | %% ================================================================================================= |
---|
388 | \subsubsection{Required fields} |
---|
389 | \label{sec:DOM_hgr_fields} |
---|
390 | |
---|
391 | The explicit specification of a range of mesh-related fields are required for |
---|
392 | the definition of a configuration. |
---|
393 | These include: |
---|
394 | |
---|
395 | \begin{clines} |
---|
396 | integer Ni0glo, NjOglo, jpkglo /* global domain sizes (without MPI halos) */ |
---|
397 | logical l\_Iperio, l\_Jperio /* lateral global domain b.c.: i- j-periodicity */ |
---|
398 | logical l\_NFold /* lateral global domain b.c.: North Pole folding */ |
---|
399 | char(1) c\_NFtype /* type of North pole Folding: T or F point */ |
---|
400 | real glamt, glamu, glamv, glamf /* geographic longitude (t,u,v and f points respectively) */ |
---|
401 | real gphit, gphiu, gphiv, gphif /* geographic latitude */ |
---|
402 | real e1t, e1u, e1v, e1f /* horizontal scale factors */ |
---|
403 | real 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 | %% ================================================================================================= |
---|
414 | \subsubsection{Optional fields} |
---|
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. |
---|
440 | |
---|
441 | \smallskip |
---|
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 | \label{subsec:DOM_zgr} |
---|
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: |
---|
465 | \begin{enumerate} |
---|
466 | \item the bathymetry given in meters; |
---|
467 | \item the number of levels of the model (\texttt{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. |
---|
472 | \end{enumerate} |
---|
473 | |
---|
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} |
---|
490 | \end{figure} |
---|
491 | |
---|
492 | The choice of a vertical coordinate is made when setting up the configuration; |
---|
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. |
---|
522 | |
---|
523 | \medskip |
---|
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): |
---|
526 | |
---|
527 | \begin{itemize} |
---|
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}). |
---|
531 | \end{itemize} |
---|
532 | |
---|
533 | Additionally, hybrid combinations of the three main coordinates are available: |
---|
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 | %% ================================================================================================= |
---|
556 | \subsubsection{Required fields} |
---|
557 | \label{sec:DOM_zgr_fields} |
---|
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} |
---|
564 | int ln_zco, ln_zps, ln_sco /* flags for z-coord, z-coord with partial steps and s-coord */ |
---|
565 | int ln_isfcav /* flag for ice shelf cavities */ |
---|
566 | double e3t_1d, e3w_1d /* reference vertical scale factors at T and W points */ |
---|
567 | double e3t_0, e3u_0, e3v_0, e3f_0, e3w_0 /* vertical scale factors 3D coordinate at T,U,V,F and W points */ |
---|
568 | double e3uw_0, e3vw_0 /* vertical scale factors 3D coordinate at UW and VW points */ |
---|
569 | int bottom_level, top_level /* last wet T-points, 1st wet T-points (for ice shelf cavities) */ |
---|
570 | /* For reference: */ |
---|
571 | float bathy_metry /* bathymetry used in setting top and bottom levels */ |
---|
572 | \end{clines} |
---|
573 | |
---|
574 | This set of vertical metrics is sufficient to describe the initial depth and thickness of |
---|
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 \texttt{bottom\_level} |
---|
578 | (and, possibly, \texttt{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 \texttt{bottom\_level} and \texttt{top\_level} 2D arrays define |
---|
588 | the \texttt{bottom\_level} and top wet levels in each grid column. |
---|
589 | Without ice cavities, \texttt{top\_level} is essentially a land mask (0 on land; 1 everywhere else). |
---|
590 | With ice cavities, \texttt{top\_level} determines the first wet point below the overlying ice shelf. |
---|
591 | |
---|
592 | %% ================================================================================================= |
---|
593 | \subsubsection{Level bathymetry and mask} |
---|
594 | \label{subsec:DOM_msk} |
---|
595 | |
---|
596 | From \texttt{top\_level} and \texttt{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*} |
---|
610 | |
---|
611 | Note that, without ice shelves cavities, |
---|
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) |
---|
616 | exactly in the same way as for the bottom boundary. |
---|
617 | |
---|
618 | %% The specification of closed lateral boundaries requires that at least |
---|
619 | %% the first and last rows and columns of the \textit{mbathy} array are set to zero. |
---|
620 | %% In the particular case of an east-west cyclical boundary condition, \textit{mbathy} has its last column equal to |
---|
621 | %% the second one and its first column equal to the last but one (and so too the mask arrays) |
---|
622 | %% (see \autoref{fig:LBC_jperio}). |
---|
623 | |
---|
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 | %% ================================================================================================= |
---|
649 | \subsection{Output grid files} |
---|
650 | \label{subsec:DOM_meshmask} |
---|
651 | |
---|
652 | 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 |
---|
659 | checking or confirmation is required. |
---|
660 | |
---|
661 | Alternatively, all the arrays relating to a particular ocean model configuration |
---|
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. |
---|
681 | |
---|
682 | \begin{description} |
---|
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}). |
---|
695 | \end{description} |
---|
696 | |
---|
697 | \subinc{\input{../../global/epilogue}} |
---|
698 | |
---|
699 | \end{document} |
---|