source: NEMO/branches/2020/ticket_2444/doc/latex/NEMO/subfiles/apdx_DOMAINcfg.tex @ 12769

Last change on this file since 12769 was 12769, checked in by mathiot, 8 months ago

#2444: update of the documentation (prior to changes suggested by Dave)

File size: 36.9 KB
Line 
1\documentclass[../main/NEMO_manual]{subfiles}
2
3\begin{document}
4
5\chapter{A brief guide to the DOMAINcfg tool}
6\label{apdx:DOMCFG}
7
8%    {\em 4.0} & {\em Andrew Coward} & {\em Created at v4.0 from materials removed from chap\_DOM that are still relevant to the \forcode{DOMAINcfg} tool and which illustrate and explain the choices to be made by the user when setting up new domains }  \\
9
10\thispagestyle{plain}
11
12\chaptertoc
13
14\paragraph{Changes record} ~\\
15
16{\footnotesize
17  \begin{tabularx}{\textwidth}{l||X|X}
18    Release & Author(s) & Modifications \\
19    \hline
20    {\em   4.0} & {\em ...} & {\em ...} \\
21    {\em   3.6} & {\em ...} & {\em ...} \\
22    {\em   3.4} & {\em ...} & {\em ...} \\
23    {\em <=3.4} & {\em ...} & {\em ...}
24  \end{tabularx}
25}
26
27\clearpage
28
29This appendix briefly describes some of the options available in the
30\forcode{DOMAINcfg} tool mentioned in \autoref{chap:DOM}.
31
32This tool will evolve into an independent utility with its own documentation but its
33current manifestation is mostly a wrapper for \NEMO\ \forcode{DOM} modules more aligned to
34those in the previous versions of \NEMO. These versions allowed the user to define some
35horizontal and vertical grids through additional namelist parameters. Explanations of
36these parameters are retained here for reference pending better documentation for
37\forcode{DOMAINcfg}. Please note that the namelist blocks named in this appendix refer to
38those read by \forcode{DOMAINcfg} via its own \forcode{namelist_ref} and
39\forcode{namelist_cfg} files. Although, due to their origins, these namelists share names
40with those used by \NEMO, they are not interchangeable and should be considered independent
41of those described elsewhere in this manual.
42
43%% =================================================================================================
44\section{Choice of horizontal grid}
45\label{sec:DOMCFG_hor}
46
47\begin{listing}
48  \nlst{namdom_domcfg}
49  \caption{\forcode{&namdom_domcfg}}
50  \label{lst:namdom_domcfg}
51\end{listing}
52
53The user has three options available in defining a horizontal grid, which involve the
54namelist variable \np{jphgr_mesh}{jphgr\_mesh} of the \nam{dom}{dom} (\texttt{DOMAINcfg} variant only)
55namelist.
56
57\begin{description}
58 \item [{\np{jphgr_mesh}{jphgr\_mesh}=0}]  The most general curvilinear orthogonal grids.
59  The coordinates and their first derivatives with respect to $i$ and $j$ are provided
60  in a input file (\ifile{coordinates}), read in \rou{hgr\_read} subroutine of the domhgr module.
61  This is now the only option available within \NEMO\ itself from v4.0 onwards.
62\item [{\np{jphgr_mesh}{jphgr\_mesh}=1 to 5}] A few simple analytical grids are provided (see below).
63  For other analytical grids, the \mdl{domhgr} module (\texttt{DOMAINcfg} variant) must be
64  modified by the user. In most cases, modifying the \mdl{usrdef\_hgr} module of \NEMO\ is
65  a better alternative since this is designed to allow simple analytical domains to be
66  configured and used without the need for external data files.
67\end{description}
68
69There are two simple cases of geographical grids on the sphere. With
70\np{jphgr_mesh}{jphgr\_mesh}=1, the grid (expressed in degrees) is regular in space,
71with grid sizes specified by parameters \np{ppe1_deg}{ppe1\_deg} and \np{ppe2_deg}{ppe2\_deg},
72respectively. Such a geographical grid can be very anisotropic at high latitudes
73because of the convergence of meridians (the zonal scale factors $e_1$
74become much smaller than the meridional scale factors $e_2$). The Mercator
75grid (\np{jphgr_mesh}{jphgr\_mesh}=4) avoids this anisotropy by refining the meridional scale
76factors in the same way as the zonal ones. In this case, meridional scale factors
77and latitudes are calculated analytically using the formulae appropriate for
78a Mercator projection, based on \np{ppe1_deg}{ppe1\_deg} which is a reference grid spacing
79at the equator (this applies even when the geographical equator is situated outside
80the model domain).
81
82In these two cases (\np{jphgr_mesh}{jphgr\_mesh}=1 or 4), the grid position is defined by the
83longitude and latitude of the south-westernmost point (\np{ppglamt0}
84and \np{ppgphi0}{ppgphi0}). Note that for the Mercator grid the user need only provide
85an approximate starting latitude: the real latitude will be recalculated analytically,
86in order to ensure that the equator corresponds to line passing through $t$-
87and $u$-points.
88
89Rectangular grids ignoring the spherical geometry are defined with
90\np{jphgr_mesh}{jphgr\_mesh} = 2, 3, 5. The domain is either an $f$-plane (\np{jphgr_mesh}{jphgr\_mesh} = 2,
91Coriolis factor is constant) or a beta-plane (\np{jphgr_mesh}{jphgr\_mesh} = 3, the Coriolis factor
92is linear in the $j$-direction). The grid size is uniform in meter in each direction,
93and given by the parameters \np{ppe1_m}{ppe1\_m} and \np{ppe2_m}{ppe2\_m} respectively.
94The zonal grid coordinate (\textit{glam} arrays) is in kilometers, starting at zero
95with the first $t$-point. The meridional coordinate (gphi. arrays) is in kilometers,
96and the second $t$-point corresponds to coordinate $gphit=0$. The input
97variable \np{ppglam0}{ppglam0} is ignored. \np{ppgphi0}{ppgphi0} is used to set the reference
98latitude for computation of the Coriolis parameter. In the case of the beta plane,
99\np{ppgphi0}{ppgphi0} corresponds to the center of the domain. Finally, the special case
100\np{jphgr_mesh}{jphgr\_mesh}=5 corresponds to a beta plane in a rotated domain for the
101GYRE configuration, representing a classical mid-latitude double gyre system.
102The rotation allows us to maximize the jet length relative to the gyre areas
103(and the number of grid points).
104
105%% =================================================================================================
106\section{Vertical grid}
107\label{sec:DOMCFG_vert}
108
109%% =================================================================================================
110\subsection{Vertical reference coordinate}
111\label{sec:DOMCFG_zref}
112
113\begin{figure}[!tb]
114  \centering
115  \includegraphics[width=0.66\textwidth]{DOMCFG_zgr}
116  \caption[DOMAINcfg: default vertical mesh for ORCA2]{
117    Default vertical mesh for ORCA2: 30 ocean levels (L30).
118    Vertical level functions for (a) T-point depth and (b) the associated scale factor for
119    the $z$-coordinate case.}
120  \label{fig:DOMCFG_zgr}
121\end{figure}
122
123The reference coordinate transformation $z_0(k)$ defines the arrays $gdept_0$ and
124$gdepw_0$ for $t$- and $w$-points, respectively. See \autoref{sec:DOMCFG_sco} for the
125S-coordinate options.  As indicated on \autoref{fig:DOM_index_vert} \jp{jpk} is the number of
126$w$-levels.  $gdepw_0(1)$ is the ocean surface.  There are at most \jp{jpk}-1 $t$-points
127inside the ocean, the additional $t$-point at $jk = jpk$ is below the sea floor and is not
128used.  The vertical location of $w$- and $t$-levels is defined from the analytic
129expression of the depth $z_0(k)$ whose analytical derivative with respect to $k$ provides
130the vertical scale factors.  The user must provide the analytical expression of both $z_0$
131and its first derivative with respect to $k$.  This is done in routine \mdl{domzgr}
132through statement functions, using parameters provided in the \nam{dom}{dom} namelist
133(\texttt{DOMAINcfg} variant).
134
135It is possible to define a simple regular vertical grid by giving zero stretching
136(\np[=0]{ppacr}{ppacr}).  In that case, the parameters \jp{jpk} (number of $w$-levels)
137and \np{pphmax}{pphmax} (total ocean depth in meters) fully define the grid.
138
139For climate-related studies it is often desirable to concentrate the vertical resolution
140near the ocean surface.  The following function is proposed as a standard for a
141$z$-coordinate (with either full or partial steps):
142\begin{gather}
143  \label{eq:DOMCFG_zgr_ana_1}
144    z_0  (k) = h_{sur} - h_0 \; k - \; h_1 \; \log  \big[ \cosh ((k - h_{th}) / h_{cr}) \big] \\
145    e_3^0(k) = \lt|    - h_0      -    h_1 \; \tanh \big[        (k - h_{th}) / h_{cr}  \big] \rt|
146\end{gather}
147
148where $k = 1$ to \jp{jpk} for $w$-levels and $k = 1$ to $k = 1$ for $t-$levels.  Such an
149expression allows us to define a nearly uniform vertical location of levels at the ocean
150top and bottom with a smooth hyperbolic tangent transition in between (\autoref{fig:DOMCFG_zgr}).
151
152A double hyperbolic tangent version (\np[=.true.]{ldbletanh}{ldbletanh}) is also available
153which permits finer control and is used, typically, to obtain a well resolved upper ocean
154without compromising on resolution at depth using a moderate number of levels.
155
156\begin{gather}
157  \label{eq:DOMCFG_zgr_ana_1b}
158    \begin{split}
159    z_0  (k) = h_{sur} - h_0 \; k &- \; h_1 \; \log  \big[ \cosh ((k - h_{th}) / h_{cr}) \big] \\
160                             \;   &- \; h2_1 \; \log  \big[ \cosh ((k - h2_{th}) / h2_{cr}) \big]
161    \end{split}
162\end{gather}
163\begin{gather}
164    \begin{split}
165    e_3^0(k) = \big|    - h_0    &-   h_1 \; \tanh \big[       (k - h_{th})  / h_{cr}   \big]  \\
166                                 &-  h2_1 \; \tanh \big[       (k - h2_{th}) / h2_{cr}  \big] \big|
167    \end{split}
168\end{gather}
169
170If the ice shelf cavities are opened (\np[=.true.]{ln_isfcav}{ln\_isfcav}), the definition
171of $z_0$ is the same.  However, definition of $e_3^0$ at $t$- and $w$-points is
172respectively changed to:
173\begin{equation}
174  \label{eq:DOMCFG_zgr_ana_2}
175  \begin{split}
176    e_3^T(k) &= z_W (k + 1) - z_W (k    ) \\
177    e_3^W(k) &= z_T (k    ) - z_T (k - 1)
178  \end{split}
179\end{equation}
180
181This formulation decreases the self-generated circulation into the ice shelf cavity
182(which can, in extreme case, leads to numerical instability). This is now the recommended formulation for all configurations using v4.0 onwards. The analytical derivation of thicknesses is maintained for backwards compatibility.
183
184The most used vertical grid for ORCA2 has $10~m$ ($500~m$) resolution in the surface
185(bottom) layers and a depth which varies from 0 at the sea surface to a minimum of
186$-5000~m$.  This leads to the following conditions:
187
188\begin{equation}
189  \label{eq:DOMCFG_zgr_coef}
190  \begin{array}{ll}
191    e_3 (1   + 1/2) =  10. & z(1  ) =     0. \\
192    e_3 (jpk - 1/2) = 500. & z(jpk) = -5000.
193  \end{array}
194\end{equation}
195
196With the choice of the stretching $h_{cr} = 3$ and the number of levels \jp{jpk}~$= 31$,
197the four coefficients $h_{sur}$, $h_0$, $h_1$, and $h_{th}$ in
198\autoref{eq:DOMCFG_zgr_ana_2} have been determined such that \autoref{eq:DOMCFG_zgr_coef}
199is satisfied, through an optimisation procedure using a bisection method.
200For the first standard ORCA2 vertical grid this led to the following values:
201$h_{sur} = 4762.96$, $h_0 = 255.58, h_1 = 245.5813$, and $h_{th} = 21.43336$.
202The resulting depths and scale factors as a function of the model levels are shown in
203\autoref{fig:DOMCFG_zgr} and given in \autoref{tab:DOMCFG_orca_zgr}.
204Those values correspond to the parameters \np{ppsur}{ppsur}, \np{ppa0}{ppa0}, \np{ppa1}{ppa1}, \np{ppkth}{ppkth} in \nam{cfg}{cfg} namelist.
205
206Rather than entering parameters $h_{sur}$, $h_0$, and $h_1$ directly, it is possible to
207recalculate them.  In that case the user sets \np{ppsur}{ppsur}~$=$~\np{ppa0}{ppa0}~$=$~\np{ppa1}{ppa1}~$=
208999999$., in \nam{cfg}{cfg} namelist, and specifies instead the four following parameters:
209\begin{itemize}
210\item \np{ppacr}{ppacr}~$= h_{cr}$: stretching factor (nondimensional).
211  The larger \np{ppacr}{ppacr}, the smaller the stretching.
212  Values from $3$ to $10$ are usual.
213\item \np{ppkth}{ppkth}~$= h_{th}$: is approximately the model level at which maximum stretching occurs
214  (nondimensional, usually of order 1/2 or 2/3 of \jp{jpk})
215\item \np{ppdzmin}{ppdzmin}: minimum thickness for the top layer (in meters).
216\item \np{pphmax}{pphmax}: total depth of the ocean (meters).
217\end{itemize}
218
219As an example, for the $45$ layers used in the DRAKKAR configuration those parameters are:
220\jp{jpk}~$= 46$, \np{ppacr}{ppacr}~$= 9$, \np{ppkth}{ppkth}~$= 23.563$, \np{ppdzmin}{ppdzmin}~$= 6~m$,
221\np{pphmax}{pphmax}~$= 5750~m$.
222
223\begin{table}
224  \centering
225  \begin{tabular}{c||r|r|r|r}
226    \hline
227    \textbf{LEVEL} & \textbf{gdept\_1d} & \textbf{gdepw\_1d} & \textbf{e3t\_1d } & \textbf{e3w\_1d} \\
228    \hline
229    1              & \textbf{     5.00} &               0.00 & \textbf{   10.00} &            10.00 \\
230    \hline
231    2              & \textbf{    15.00} &              10.00 & \textbf{   10.00} &            10.00 \\
232    \hline
233    3              & \textbf{    25.00} &              20.00 & \textbf{   10.00} &            10.00 \\
234    \hline
235    4              & \textbf{    35.01} &              30.00 & \textbf{   10.01} &            10.00 \\
236    \hline
237    5              & \textbf{    45.01} &              40.01 & \textbf{   10.01} &            10.01 \\
238    \hline
239    6              & \textbf{    55.03} &              50.02 & \textbf{   10.02} &            10.02 \\
240    \hline
241    7              & \textbf{    65.06} &              60.04 & \textbf{   10.04} &            10.03 \\
242    \hline
243    8              & \textbf{    75.13} &              70.09 & \textbf{   10.09} &            10.06 \\
244    \hline
245    9              & \textbf{    85.25} &              80.18 & \textbf{   10.17} &            10.12 \\
246    \hline
247    10             & \textbf{    95.49} &              90.35 & \textbf{   10.33} &            10.24 \\
248    \hline
249    11             & \textbf{   105.97} &             100.69 & \textbf{   10.65} &            10.47 \\
250    \hline
251    12             & \textbf{   116.90} &             111.36 & \textbf{   11.27} &            10.91 \\
252    \hline
253    13             & \textbf{   128.70} &             122.65 & \textbf{   12.47} &            11.77 \\
254    \hline
255    14             & \textbf{   142.20} &             135.16 & \textbf{   14.78} &            13.43 \\
256    \hline
257    15             & \textbf{   158.96} &             150.03 & \textbf{   19.23} &            16.65 \\
258    \hline
259    16             & \textbf{   181.96} &             169.42 & \textbf{   27.66} &            22.78 \\
260    \hline
261    17             & \textbf{   216.65} &             197.37 & \textbf{   43.26} &            34.30 \\
262    \hline
263    18             & \textbf{   272.48} &             241.13 & \textbf{   70.88} &            55.21 \\
264    \hline
265    19             & \textbf{   364.30} &             312.74 & \textbf{  116.11} &            90.99 \\
266    \hline
267    20             & \textbf{   511.53} &             429.72 & \textbf{  181.55} &           146.43 \\
268    \hline
269    21             & \textbf{   732.20} &             611.89 & \textbf{  261.03} &           220.35 \\
270    \hline
271    22             & \textbf{  1033.22} &             872.87 & \textbf{  339.39} &           301.42 \\
272    \hline
273    23             & \textbf{  1405.70} &            1211.59 & \textbf{  402.26} &           373.31 \\
274    \hline
275    24             & \textbf{  1830.89} &            1612.98 & \textbf{  444.87} &           426.00 \\
276    \hline
277    25             & \textbf{  2289.77} &            2057.13 & \textbf{  470.55} &           459.47 \\
278    \hline
279    26             & \textbf{  2768.24} &            2527.22 & \textbf{  484.95} &           478.83 \\
280    \hline
281    27             & \textbf{  3257.48} &            3011.90 & \textbf{  492.70} &           489.44 \\
282    \hline
283    28             & \textbf{  3752.44} &            3504.46 & \textbf{  496.78} &           495.07 \\
284    \hline
285    29             & \textbf{  4250.40} &            4001.16 & \textbf{  498.90} &           498.02 \\
286    \hline
287    30             & \textbf{  4749.91} &            4500.02 & \textbf{  500.00} &           499.54 \\
288    \hline
289    31             & \textbf{  5250.23} &            5000.00 & \textbf{  500.56} &           500.33 \\
290    \hline
291  \end{tabular}
292  \caption[Default vertical mesh in $z$-coordinate for 30 layers ORCA2 configuration]{
293    Default vertical mesh in $z$-coordinate for 30 layers ORCA2 configuration as
294    computed from \autoref{eq:DOMCFG_zgr_ana_2} using
295    the coefficients given in \autoref{eq:DOMCFG_zgr_coef}}
296  \label{tab:DOMCFG_orca_zgr}
297\end{table}
298%%%YY
299%% % -------------------------------------------------------------------------------------------------------------
300%% %        Meter Bathymetry
301%% % -------------------------------------------------------------------------------------------------------------
302%% =================================================================================================
303\subsection{Model bathymetry}
304\label{subsec:DOMCFG_bathy}
305
306Three options are possible for defining the bathymetry, according to the namelist variable
307\np{nn_bathy}{nn\_bathy} (found in \nam{dom}{dom} namelist (\texttt{DOMAINCFG} variant) ):
308\begin{description}
309\item [{\np[=0]{nn_bathy}{nn\_bathy}}]: a flat-bottom domain is defined.
310  The total depth $z_w (jpk)$ is given by the coordinate transformation.
311  The domain can either be a closed basin or a periodic channel depending on the parameter \np{jperio}{jperio}.
312\item [{\np[=-1]{nn_bathy}{nn\_bathy}}]: a domain with a bump of topography one third of the domain width at the central latitude.
313  This is meant for the "EEL-R5" configuration, a periodic or open boundary channel with a seamount.
314\item [{\np[=1]{nn_bathy}{nn\_bathy}}]: read a bathymetry and ice shelf draft (if needed).
315  The \ifile{bathy\_meter} file (Netcdf format) provides the ocean depth (positive, in meters) at
316  each grid point of the model grid.
317  The bathymetry is usually built by interpolating a standard bathymetry product (\eg\ ETOPO2) onto
318  the horizontal ocean mesh.
319  Defining the bathymetry also defines the coastline: where the bathymetry is zero,
320  no wet levels are defined (all levels are masked).
321\end{description}
322
323%% =================================================================================================
324\subsection{Choice of vertical grid}
325\label{sec:DOMCFG_vgrd}
326
327After reading the bathymetry, the algorithm for vertical grid definition differs between the different options:
328\begin{description}
329\item [\forcode{ln_zco = .true.}] set a reference coordinate transformation $z_0(k)$, and set $z(i,j,k,t) = z_0(k)$ where $z_0(k)$ is the closest match to the depth at $(i,j)$.
330\item [\forcode{ln_zps = .true.}] set a reference coordinate transformation $z_0(k)$, and calculate the thickness of the deepest level at
331  each $(i,j)$ point using the bathymetry, to obtain the final three-dimensional depth and scale factor arrays.
332\item [\forcode{ln_sco = .true.}] smooth the bathymetry to fulfill the hydrostatic consistency criteria and
333  set the three-dimensional transformation.
334\item [\forcode{s-z and s-zps}] smooth the bathymetry to fulfill the hydrostatic consistency criteria and
335  set the three-dimensional transformation $z(i,j,k)$,
336  and possibly introduce masking of extra land points to better fit the original bathymetry file.
337\end{description}
338
339%% =================================================================================================
340\subsubsection[$Z$-coordinate with uniform thickness levels (\forcode{ln_zco})]{$Z$-coordinate with uniform thickness levels (\protect\np{ln_zco}{ln\_zco})}
341\label{subsec:DOMCFG_zco}
342
343With this option the model topography can be fully described by the reference vertical
344coordinate and a 2D integer field giving the number of wet levels at each location
345(\forcode{bathy_level}). The resulting match to the real topography is likely to be poor
346though (especially with thick, deep levels) and slopes poorly represented. This option is
347rarely used in modern simulations but it can be useful for testing purposes.
348
349%% =================================================================================================
350\subsubsection[$Z$-coordinate with partial step (\forcode{ln_zps})]{$Z$-coordinate with partial step (\protect\np{ln_zps}{ln\_zps})}
351\label{subsec:DOMCFG_zps}
352
353In $z$-coordinate partial step, the depths of the model levels are defined by the
354reference analytical function $z_0(k)$ as described in \autoref{sec:DOMCFG_zref},
355\textit{except} in the bottom layer.  The thickness of the bottom layer is allowed to vary
356as a function of geographical location $(\lambda,\varphi)$ to allow a better
357representation of the bathymetry, especially in the case of small slopes (where the
358bathymetry varies by less than one level thickness from one grid point to the next).  The
359reference layer thicknesses $e_{3t}^0$ have been defined in the absence of bathymetry.
360With partial steps, layers from 1 to \jp{jpk}-2 can have a thickness smaller than
361$e_{3t}(jk)$.
362
363The model deepest layer (\jp{jpk}-1) is allowed to have either a smaller or larger
364thickness than $e_{3t}(jpk)$: the maximum thickness allowed is $2*e_{3t}(jpk - 1)$.
365
366This has to be kept in mind when specifying values in \nam{dom}{dom} namelist
367(\texttt{DOMAINCFG} variant), such as the maximum depth \np{pphmax}{pphmax} in partial steps.
368
369For example, with \np{pphmax}{pphmax}~$= 5750~m$ for the DRAKKAR 45 layer grid, the maximum ocean
370depth allowed is actually $6000~m$ (the default thickness $e_{3t}(jpk - 1)$ being
371$250~m$).  Two variables in the namdom namelist are used to define the partial step
372vertical grid.  The mimimum water thickness (in meters) allowed for a cell partially
373filled with bathymetry at level jk is the minimum of \np{rn_e3zps_min}{rn\_e3zps\_min} (thickness in
374meters, usually $20~m$) or $e_{3t}(jk)*$\np{rn_e3zps_rat}{rn\_e3zps\_rat} (a fraction, usually 10\%, of
375the default thickness $e_{3t}(jk)$).
376
377%% =================================================================================================
378\subsubsection[$S$-coordinate (\forcode{ln_sco})]{$S$-coordinate (\protect\np{ln_sco}{ln\_sco})}
379\label{sec:DOMCFG_sco}
380\begin{listing}
381  \nlst{namzgr_sco_domcfg}
382  \caption{\forcode{&namzgr_sco_domcfg}}
383  \label{lst:namzgr_sco_domcfg}
384\end{listing}
385Options are defined in \nam{zgr_sco}{zgr\_sco} (\texttt{DOMAINcfg} only).
386In $s$-coordinate (\np[=.true.]{ln_sco}{ln\_sco}), the depth and thickness of the model levels are defined from
387the product of a depth field and either a stretching function or its derivative, respectively:
388
389\begin{align*}
390  % \label{eq:DOMCFG_sco_ana}
391  z(k)   &= h(i,j) \; z_0 (k) \\
392  e_3(k) &= h(i,j) \; z_0'(k)
393\end{align*}
394
395where $h$ is the depth of the last $w$-level ($z_0(k)$) defined at the $t$-point location in the horizontal and
396$z_0(k)$ is a function which varies from $0$ at the sea surface to $1$ at the ocean bottom.
397The depth field $h$ is not necessary the ocean depth,
398since a mixed step-like and bottom-following representation of the topography can be used
399(\autoref{fig:DOM_z_zps_s_sps}) or an envelop bathymetry can be defined (\autoref{fig:DOM_z_zps_s_sps}).
400The namelist parameter \np{rn_rmax}{rn\_rmax} determines the slope at which
401the terrain-following coordinate intersects the sea bed and becomes a pseudo z-coordinate.
402The coordinate can also be hybridised by specifying \np{rn_sbot_min}{rn\_sbot\_min} and \np{rn_sbot_max}{rn\_sbot\_max} as
403the minimum and maximum depths at which the terrain-following vertical coordinate is calculated.
404
405Options for stretching the coordinate are provided as examples,
406but care must be taken to ensure that the vertical stretch used is appropriate for the application.
407
408The original default \NEMO\ s-coordinate stretching is available if neither of the other options are specified as true
409(\np[=.false.]{ln_s_SH94}{ln\_s\_SH94} and \np[=.false.]{ln_s_SF12}{ln\_s\_SF12}).
410This uses a depth independent $\tanh$ function for the stretching \citep{madec.delecluse.ea_JPO96}:
411
412\[
413  z = s_{min} + C (s) (H - s_{min})
414  % \label{eq:DOMCFG_SH94_1}
415\]
416
417where $s_{min}$ is the depth at which the $s$-coordinate stretching starts and
418allows a $z$-coordinate to placed on top of the stretched coordinate,
419and $z$ is the depth (negative down from the asea surface).
420\begin{gather*}
421  s = - \frac{k}{n - 1} \quad \text{and} \quad 0 \leq k \leq n - 1
422  % \label{eq:DOMCFG_s}
423 \\
424 \label{eq:DOMCFG_sco_function}
425  C(s) = \frac{[\tanh(\theta \, (s + b)) - \tanh(\theta \, b)]}{2 \; \sinh(\theta)}
426\end{gather*}
427
428A stretching function,
429modified from the commonly used \citet{song.haidvogel_JCP94} stretching (\np[=.true.]{ln_s_SH94}{ln\_s\_SH94}),
430is also available and is more commonly used for shelf seas modelling:
431
432\[
433  C(s) =   (1 - b) \frac{\sinh(\theta s)}{\sinh(\theta)}
434         + b       \frac{\tanh \lt[ \theta \lt(s + \frac{1}{2} \rt) \rt] -   \tanh \lt( \frac{\theta}{2} \rt)}
435                        {                                                  2 \tanh \lt( \frac{\theta}{2} \rt)}
436 \label{eq:DOMCFG_SH94_2}
437\]
438
439\begin{figure}[!ht]
440  \centering
441  \includegraphics[width=0.66\textwidth]{DOMCFG_sco_function}
442  \caption[DOMAINcfg: examples of the stretching function applied to a seamount]{
443    Examples of the stretching function applied to a seamount;
444    from left to right: surface, surface and bottom, and bottom intensified resolutions}
445  \label{fig:DOMCFG_sco_function}
446\end{figure}
447
448where $H_c$ is the critical depth (\np{rn_hc}{rn\_hc}) at which the coordinate transitions from pure $\sigma$ to
449the stretched coordinate, and $\theta$ (\np{rn_theta}{rn\_theta}) and $b$ (\np{rn_bb}{rn\_bb}) are the surface and
450bottom control parameters such that $0 \leqslant \theta \leqslant 20$, and $0 \leqslant b \leqslant 1$.
451$b$ has been designed to allow surface and/or bottom increase of the vertical resolution
452(\autoref{fig:DOMCFG_sco_function}).
453
454Another example has been provided at version 3.5 (\np{ln_s_SF12}{ln\_s\_SF12}) that allows a fixed surface resolution in
455an analytical terrain-following stretching \citet{siddorn.furner_OM13}.
456In this case the a stretching function $\gamma$ is defined such that:
457
458\begin{equation}
459  z = - \gamma h \quad \text{with} \quad 0 \leq \gamma \leq 1
460  % \label{eq:DOMCFG_z}
461\end{equation}
462
463The function is defined with respect to $\sigma$, the unstretched terrain-following coordinate:
464
465\begin{gather*}
466  % \label{eq:DOMCFG_gamma_deriv}
467  \gamma =   A \lt( \sigma   - \frac{1}{2} (\sigma^2     + f (\sigma)) \rt)
468           + B \lt( \sigma^3 - f           (\sigma) \rt) + f (\sigma)       \\
469  \intertext{Where:}
470 \label{eq:DOMCFG_gamma}
471  f(\sigma) = (\alpha + 2) \sigma^{\alpha + 1} - (\alpha + 1) \sigma^{\alpha + 2}
472  \quad \text{and} \quad \sigma = \frac{k}{n - 1}
473\end{gather*}
474
475This gives an analytical stretching of $\sigma$ that is solvable in $A$ and $B$ as a function of
476the user prescribed stretching parameter $\alpha$ (\np{rn_alpha}{rn\_alpha}) that stretches towards
477the surface ($\alpha > 1.0$) or the bottom ($\alpha < 1.0$) and
478user prescribed surface (\np{rn_zs}{rn\_zs}) and bottom depths.
479The bottom cell depth in this example is given as a function of water depth:
480
481\[
482  % \label{eq:DOMCFG_zb}
483  Z_b = h a + b
484\]
485
486where the namelist parameters \np{rn_zb_a}{rn\_zb\_a} and \np{rn_zb_b}{rn\_zb\_b} are $a$ and $b$ respectively.
487
488\begin{figure}[!ht]
489  \centering
490  \includegraphics[width=0.66\textwidth]{DOMCFG_compare_coordinates_surface}
491  \caption[DOMAINcfg: comparison of $s$- and $z$-coordinate]{
492    A comparison of the \citet{song.haidvogel_JCP94} $S$-coordinate (solid lines),
493    a 50 level $Z$-coordinate (contoured surfaces) and
494    the \citet{siddorn.furner_OM13} $S$-coordinate (dashed lines) in the surface $100~m$ for
495    a idealised bathymetry that goes from $50~m$ to $5500~m$ depth.
496    For clarity every third coordinate surface is shown.}
497  \label{fig:DOMCFG_fig_compare_coordinates_surface}
498\end{figure}
499 % >>>>>>>>>>>>>>>>>>>>>>>>>>>>
500
501This gives a smooth analytical stretching in computational space that is constrained to
502given specified surface and bottom grid cell thicknesses in real space.
503This is not to be confused with the hybrid schemes that
504superimpose geopotential coordinates on terrain following coordinates thus
505creating a non-analytical vertical coordinate that
506therefore may suffer from large gradients in the vertical resolutions.
507This stretching is less straightforward to implement than the \citet{song.haidvogel_JCP94} stretching,
508but has the advantage of resolving diurnal processes in deep water and has generally flatter slopes.
509
510As with the \citet{song.haidvogel_JCP94} stretching the stretch is only applied at depths greater than
511the critical depth $h_c$.
512In this example two options are available in depths shallower than $h_c$,
513with pure sigma being applied if the \np{ln_sigcrit}{ln\_sigcrit} is true and pure z-coordinates if it is false
514(the z-coordinate being equal to the depths of the stretched coordinate at $h_c$).
515
516Minimising the horizontal slope of the vertical coordinate is important in terrain-following systems as
517large slopes lead to hydrostatic consistency.
518A hydrostatic consistency parameter diagnostic following \citet{haney_JPO91} has been implemented,
519and is output as part of the model mesh file at the start of the run.
520
521%% =================================================================================================
522\subsubsection[\zstar- or \sstar-coordinate (\forcode{ln_linssh})]{\zstar- or \sstar-coordinate (\protect\np{ln_linssh}{ln\_linssh})}
523\label{subsec:DOMCFG_zgr_star}
524
525This option is described in the Report by Levier \textit{et al.} (2007), available on the \NEMO\ web site.
526
527\section{Ice shelf cavity definition}
528\label{subsec:zgrisf}
529
530  If the under ice shelf seas are opened (\np{ln_isfcav}{ln\_isfcav}), the depth of the ice shelf/ocean interface has to be include in
531  the \ifile{isfdraft\_meter} file (Netcdf format). This file need to include the \ifile{isf\_draft} variable.
532  A positive value will me an ice shelf/ocean or ice shelf bedrock interface below the reference 0m ssh.
533  The exact shape of the ice shelf cavity (grounding line position and minimum thickness of the water column under an ice shelf, ...) can be specify in \nam{lst:namzgr_isf}.
534
535
536\begin{listing}
537  \nlst{namzgr_isf_domcfg}
538  \caption{\forcode{&namzgr_isf}}
539  \label{lst:namzgr_isf}
540\end{listing}
541
542   The options available to define the shape of the under ice shelf cavities are listed in \nam{namzgr_isf}{namzgr\_isf} (\texttt{DOMAINcfg} only, \autoref{lst:namzgr_isf}).
543
544   \subsection{Model ice shelf draft definition}
545   \label{subsec:zgrisf_isfd}
546
547   First of all, the tool make sure, the ice shelf draft ($h_{isf}$) is sensible and compatible with the bathymetry.
548   There are 3 compulsory steps to achieve this:
549
550   \begin{description}
551   \item{\np{rn_isfdep_min}{rn\_isfdep\_min}:} this is the minimum ice shelf draft. This is to make sure there is no ridiculous thin ice shelf. If \np{rn_isfdep_min}{rn\_isfdep\_min} is smaller than the surface level, \np{rn_isfdep_min}{rn\_isfdep\_min} is set to $e3t\_1d(1)$.
552   Where $h_{isf} < MAX(e3t\_1d(1),\np{rn_isfdep_min}{rn\_isfdep\_min}$), $h_{isf}$ is set to \np{rn_isfdep_min}{rn\_isfdep\_min}.
553
554   \item{\np{rn_glhw_min}{rn\_glhw\_min}:} This parameter is used to define the grounding line position.
555   Where the difference between the bathymetry and the ice shelf draft is smaller than \np{rn_glhw_min}{rn\_glhw\_min}, the cell are grounded (ie masked).
556   This step is needed to take into account possible small mismatch between ice shelf draft value and bathymetry value (sources are coming from different grid, different data processes, rounding error, ...).
557
558   \item{\np{rn_isfhw_min}{rn\_isfhw\_min}:} This parameter is minimum water column thickness in the cavity.
559   Where the water column thickness is lower than \np{rn_isfhw_min}{rn\_isfhw\_min}, the ice shelf draft is adjusted to match this criterion.
560   If for any reason, this adjustement break the minimum ice shelf draft allowed (\np{rn_isfdep_min}{rn\_isfdep\_min}), the cell is masked.
561   \end{description}
562
563   Once all these adjustements are made, if the water column thickness contains one cell wide channels, these channels can be closed using \np{ln_isfchannel}{ln\_isfchannel}
564 
565   \subsection{Model top level definition}
566   After the definition of the ice shelf draft, the tool defines the top level.
567   The compulsory criterion is that the water column needs at least 2 wet cells in the water column at U- and V-points.
568   To do so, if there one cell wide water column, the tools adjust the ice shelf draft to fillful the requierement.\\
569
570   The process is the following:
571   \begin{description}
572   \item{step 1:} The top level is defined in the same way as the bottom level is defined.
573   \item{step 2:} The isolated grid point in the bathymetry are filled (as it is done in a domain without ice shelf)
574   \item{step 3:} The tools make sure, the top level is above or equal to the bottom level
575   \item{step 4:} If the water column at a U- or V- point is one wet cell wide, the ice shelf draft is adjusted. So the actual top cell become fully open and the new
576   top cell thickness is set to the minimum cell thickness allowed (following the same logic as for the bottom partial cell). This step is iterated 4 times to ensure the condition is fullfill along the 4 sides of the cell.
577   \end{description}
578
579   In case of steep slope and shallow water column, it likely that 2 cells are disconnected (bathymetry above its neigbourg ice shelf draft).
580   The option \np{ln_isfconnect}{ln\_isfconnect} allow the tool to force the connection between these 2 cells.
581   Some limiters in meter or levels on the digging allowed by the tool are available (respectively, \np{rn_zisfmax}{rn\_zisfmax} or \np{rn_kisfmax}{rn\_kisfmax}).
582   This will prevent the formation of subglacial lakes at the expense of long vertical pipe to connect cells at very different levels.
583
584   \subsection{Subglacial lakes}
585   Despite careful setting of your ice shelf draft and bathymetry input file as well as setting described in \autoref{subsec:zgrisf_isfd}, some situation are unavoidable.
586   For exemple if you setup your ice shelf draft and bathymetry to do ocean/ice sheet coupling,
587   you may decide to fill the whole antarctic with a bathymetry and an ice shelf draft value (ice/bedrock interface depth when grounded).
588   If you do so, the subglacial lakes will show up (Vostock for example). An other possibility is with coarse vertical resolution, some ice shelves could be cut in 2 parts:
589   one connected to the main ocean and an other one closed which can be considered as a subglacial sea be the model.\\
590
591   The namelist option \np{ln_isfsubgl}{ln\_isfsubgl} allow you to remove theses subglacial lakes.
592   This may be useful for esthetical reason or for stability reasons:
593
594   \begin{description}
595   \item $\bullet$ In a subglacial lakes, in case of very weak circulation (often the case), the only heat flux is the conductive heat flux through the ice sheet.
596         This will lead to constant freezing until water reaches -20C.
597         This is one of the defitiency of the 3 equation melt formulation (for details on this formulation, see: \autoref{sec:isf}).
598   \item $\bullet$ In case of coupling with an ice sheet model,
599         the ssh in the subglacial lakes and the main ocean could be very different (ssh initial adjustement for example),
600         and so if for any reason both a connected at some point, the model is likely to fall over.\\
601   \end{description}
602
603\section{Closed sea definition}
604\label{sec:clocfg}
605
606\begin{listing}
607  \nlst{namclo_domcfg}
608  \caption{\forcode{&namclo}}
609  \label{lst:namclo}
610\end{listing}
611
612The options available to define the closed seas and how closed sea net fresh water input will be redistributed by NEMO are listed in \nam{clo} (\texttt{DOMAINcfg} only, \autoref{lst:namclo}).
613The individual definition of each closed sea is managed by \np{sn_lake}{sn\_lake}. In this fields the user needs to defined:\\
614   \begin{description}
615   \item $\bullet$    the name of the closed sea (print output purposes).
616   \item $\bullet$    the seed location to define the area of the closed sea (if seed on land because not present in this configuration, this closed sea will be ignored).\\
617   \item $\bullet$    the seed location for the target area.
618   \item $\bullet$    the type of target area ('local','coast' or 'global'). See point 6 for definition of these cases.
619   \item $\bullet$    the type of redistribution scheme for the net fresh water flux over the closed sea (as a runoff in a target area, as emp in a target area, as emp globally). For the runoff case, if the net fwf is negative, it will be redistribut globally.
620   \item $\bullet$    the radius of the target area (not used for the 'global' case). So the target defined by a 'local' target area of a radius of 100km, for example, correspond to all the wet points within this radius. The coastal case will return only the coastal point within the specifid radius.
621   \item $\bullet$    the target id. This target id is used to group multiple lakes into the same river ouflow (Great Lakes for example).
622   \end{description}
623
624The closed sea module defines a number of masks in the \ifile{domain\_cfg} output:
625   \begin{description}
626   \item[\textit{mask\_opensea}:] a mask of the main ocean without all the closed seas closed. This mask is defined by a flood filling algorithm with an initial seed (localisation defined by \np{rn_lon_opnsea}{rn\_lon\_opnsea} and \np{rn_lat_opnsea}{rn\_lat\_opnsea}).
627   \item[\textit{mask\_csglo}, \textit{mask\_csrnf}, \textit{mask\_csemp}:] a mask of all the closed seas defined in the namelist by \np{sn_lake}{sn\_lake} for each redistribution scheme. The total number of defined closed seas has to be defined in \np{nn_closea}{nn\_closea}.
628   \item[\textit{mask\_csgrpglo}, \textit{mask\_csgrprnf}, \textit{mask\_csgrpemp}:] a mask of all the closed seas and targets grouped by target id for each type of redistribution scheme.
629   \item[\textit{mask\_csundef}:] a mask of all the closed sea not defined in \np{sn_lake}{sn\_lake}. This will allows NEMO to mask them if needed or to inform the user of potential minor issues in its bathymetry.
630   \end{description}
631   
632\subinc{\input{../../global/epilogue}}
633
634\end{document}
Note: See TracBrowser for help on using the repository browser.