source: NEMO/trunk/doc/latex/NEMO/subfiles/chap_misc.tex @ 14257

Last change on this file since 14257 was 14257, checked in by nicolasmartin, 4 months ago

Overall review of LaTeX sources (not tested completely as of now):

  • Reworking global files: main document.tex, add glossary.tex, cosmetic changes…
  • Ignore missing namelists (namsbc_isf, namsbc_iscpl and namptr)
  • Removal of references for unused indices (\hfile, \ifile and \jp)
  • Update of .svnignore and svn:ignore properties accordingly
  • Split of manual abstract into a common NEMO abs for all and a specific one for each engine
  • Shrinking variables names used in the frontmatter
File size: 27.7 KB
5\chapter{Miscellaneous Topics}
10\paragraph{Changes record} ~\\
13  \begin{tabularx}{\textwidth}{l||X|X}
14    Release & Author(s) & Modifications \\
15    \hline
16    {\em   4.0} & {\em ...} & {\em ...} \\
17    {\em   3.6} & {\em ...} & {\em ...} \\
18    {\em   3.4} & {\em ...} & {\em ...} \\
19    {\em <=3.4} & {\em ...} & {\em ...}
20  \end{tabularx}
25%% =================================================================================================
26\section{Representation of unresolved straits}
29In climate modeling, it often occurs that a crucial connections between water masses is broken as
30the grid mesh is too coarse to resolve narrow straits.
31For example, coarse grid spacing typically closes off the Mediterranean from the Atlantic at
32the Strait of Gibraltar.
33In this case, it is important for climate models to include the effects of salty water entering the Atlantic from
34the Mediterranean.
35Likewise, it is important for the Mediterranean to replenish its supply of water from the Atlantic to
36balance the net evaporation occurring over the Mediterranean region.
37This problem occurs even in eddy permitting simulations.
38For example, in ORCA 1/4\deg\ several straits of the Indonesian archipelago (Ombai, Lombok...)
39are much narrow than even a single ocean grid-point.
41We describe briefly here the two methods that can be used in \NEMO\ to handle such
42improperly resolved straits. The methods consist of opening the strait while ensuring
43that the mass exchanges through the strait are not too large by either artificially
44reducing the cross-sectional area of the strait grid-cells or, locally increasing the
45lateral friction.
47%% =================================================================================================
48\subsection{Hand made geometry changes}
51The first method involves reducing the scale factor in the cross-strait direction to a
52value in better agreement with the true mean width of the strait
53(\autoref{fig:MISC_strait_hand}).  This technique is sometime called "partially open face"
54or "partially closed cells".  The key issue here is only to reduce the faces of $T$-cell
55(\ie\ change the value of the horizontal scale factors at $u$- or $v$-point) but not the
56volume of the $T$-cell.  Indeed, reducing the volume of strait $T$-cell can easily produce
57a numerical instability at that grid point which would require a reduction of the model
58time step.  Thus to instigate a local change in the width of a Strait requires two steps:
62\item Add \texttt{e1e2u} and \texttt{e1e2v} arrays to the \np{cn_domcfg}{cn\_domcfg} file. These 2D
63arrays should contain the products of the unaltered values of: $\texttt{e1u}*\texttt{e2u}$
64and $\texttt{e1u}*\texttt{e2v}$ respectively. That is the original surface areas of $u$-
65and $v$- cells respectively.  These areas are usually defined by the corresponding product
66within the \NEMO\ code but the presence of \texttt{e1e2u} and \texttt{e1e2v} in the
67\np{cn_domcfg}{cn\_domcfg} file will suppress this calculation and use the supplied fields instead.
68If the model domain is provided by user-supplied code in \mdl{usrdef\_hgr}, then this
69routine should also return \texttt{e1e2u} and \texttt{e1e2v} and set the integer return
70argument \texttt{ie1e2u\_v} to a non-zero value. Values other than 0 for this argument
71will suppress the calculation of the areas.
73\item Change values of \texttt{e2u} or \texttt{e1v} (either in the \np{cn_domcfg}{cn\_domcfg} file or
74via code in  \mdl{usrdef\_hgr}), whereever a Strait reduction is required. The choice of
75whether to alter \texttt{e2u} or \texttt{e1v} depends. respectively,  on whether the
76Strait in question is North-South orientated (\eg\ Gibraltar) or East-West orientated (\eg
81The second method is to increase the viscous boundary layer thickness by a local increase
82of the fmask value at the coast. This method can also be effective in wider passages.  The
83concept is illustarted in the second part of  \autoref{fig:MISC_strait_hand} and changes
84to specific locations can be coded in \mdl{usrdef\_fmask}. The \forcode{usr_def_fmask}
85routine is always called after \texttt{fmask} has been defined according to the choice of
86lateral boundary condition as discussed in \autoref{sec:LBC_coast}. The default version of
87\mdl{usrdef\_fmask} contains settings specific to ORCA2 and ORCA1 configurations. These are
88meant as examples only; it is up to the user to verify settings and provide alternatives
89for their own configurations. The default \forcode{usr_def_fmask} makes no changes to
90\texttt{fmask} for any other configuration.
93  \centering
94  \includegraphics[width=0.66\textwidth]{MISC_Gibraltar}
95  \includegraphics[width=0.66\textwidth]{MISC_Gibraltar2}
96  \caption[Two methods to defined the Gibraltar strait]{
97    Example of the Gibraltar strait defined in a 1\deg\ $\times$ 1\deg\ mesh.
98    \textit{Top}: using partially open cells.
99    The meridional scale factor at $v$-point is reduced on both sides of the strait to
100    account for the real width of the strait (about 20 km).
101    Note that the scale factors of the strait $T$-point remains unchanged.
102    \textit{Bottom}: using viscous boundary layers.
103    The four fmask parameters along the strait coastlines are set to a value larger than 4,
104    \ie\ "strong" no-slip case (see \autoref{fig:LBC_shlat}) creating a large viscous boundary layer
105    that allows a reduced transport through the strait.}
106  \label{fig:MISC_strait_hand}
110  \centering
111  \includegraphics[width=0.66\textwidth]{MISC_closea_mask_example}
112  \caption[Mask fields for the \protect\mdl{closea} module]{
113    Example of mask fields for the \protect\mdl{closea} module.
114    \textit{Left}: a closea\_mask field;
115    \textit{Right}: a closea\_mask\_rnf field.
116    In this example, if \protect\np{ln_closea}{ln\_closea} is set to \forcode{.true.},
117    the mean freshwater flux over each of the American Great Lakes will be set to zero,
118    and the total residual for all the lakes, if negative, will be put into
119    the St Laurence Seaway in the area shown.}
120  \label{fig:MISC_closea_mask_example}
123%% =================================================================================================
124\section[Closed seas (\textit{closea.F90})]{Closed seas (\protect\mdl{closea})}
127Some configurations include inland seas and lakes as ocean
128points. This is particularly the case for configurations that are
129coupled to an atmosphere model where one might want to include inland
130seas and lakes as ocean model points in order to provide a better
131bottom boundary condition for the atmosphere. However there is no
132route for freshwater to run off from the lakes to the ocean and this
133can lead to large drifts in the sea surface height over the lakes. The
134closea module provides options to either fill in closed seas and lakes
135at run time, or to set the net surface freshwater flux for each lake
136to zero and put the residual flux into the ocean.
138Prior to \NEMO\ 4 the locations of inland seas and lakes was set via
139hardcoded indices for various ORCA configurations. From \NEMO\ 4 onwards
140the inland seas and lakes are defined using mask fields in the
141domain configuration file. The options are as follows.
144\item {{\bfseries No ``closea\_mask'' field is included in domain configuration
145  file.} In this case the closea module does nothing.}
147\item {{\bfseries A field called closea\_mask is included in the domain
148configuration file and ln\_closea=.false. in namelist namcfg.} In this
149case the inland seas defined by the closea\_mask field are filled in
150(turned to land points) at run time. That is every point in
151closea\_mask that is nonzero is set to be a land point.}
153\item {{\bfseries A field called closea\_mask is included in the domain
154configuration file and ln\_closea=.true. in namelist namcfg.} Each
155inland sea or group of inland seas is set to a positive integer value
156in the closea\_mask field (see \autoref{fig:MISC_closea_mask_example}
157for an example). The net surface flux over each inland sea or group of
158inland seas is set to zero each timestep and the residual flux is
159distributed over the global ocean (ie. all ocean points where
160closea\_mask is zero).}
162\item {{\bfseries Fields called closea\_mask and closea\_mask\_rnf are
163included in the domain configuration file and ln\_closea=.true. in
164namelist namcfg.} This option works as for option 3, except that if
165the net surface flux over an inland sea is negative (net
166precipitation) it is put into the ocean at specified runoff points. A
167net positive surface flux (net evaporation) is still spread over the
168global ocean. The mapping from inland seas to runoff points is defined
169by the closea\_mask\_rnf field. Each mapping is defined by a positive
170integer value for the inland sea(s) and the corresponding runoff
171points. An example is given in
172\autoref{fig:MISC_closea_mask_example}. If no mapping is provided for a
173particular inland sea then the residual is spread over the global
176\item {{\bfseries Fields called closea\_mask and closea\_mask\_emp are
177included in the domain configuration file and ln\_closea=.true. in
178namelist namcfg.} This option works the same as option 4 except that
179the nonzero net surface flux is sent to the ocean at the specified
180runoff points regardless of whether it is positive or negative. The
181mapping from inland seas to runoff points in this case is defined by
182the closea\_mask\_emp field.}
185There is a python routine to create the closea\_mask fields and append
186them to the domain configuration file in the utils/tools/DOMAINcfg directory.
188%% =================================================================================================
189\section{Sub-domain functionality}
192%% =================================================================================================
193\subsection{Simple subsetting of input files via NetCDF attributes}
195The extended grids for use with the under-shelf ice cavities will result in redundant rows
196around Antarctica if the ice cavities are not active.  A simple mechanism for subsetting
197input files associated with the extended domains has been implemented to avoid the need to
198maintain different sets of input fields for use with or without active ice cavities.  This
199subsetting operates for the j-direction only and works by optionally looking for and using
200a global file attribute (named: \np{open_ocean_jstart}{open\_ocean\_jstart}) to determine the starting j-row
201for input.  The use of this option is best explained with an example:
204\noindent Consider an ORCA1
205configuration using the extended grid domain configuration file: \textit{eORCA1\}
206This file define a horizontal domain of 362x332.  The first row with
207open ocean wet points in the non-isf bathymetry for this set is row 42 (\fortran\ indexing)
208then the formally correct setting for \np{open_ocean_jstart}{open\_ocean\_jstart} is 41.  Using this value as
209the first row to be read will result in a 362x292 domain which is the same size as the
210original ORCA1 domain.  Thus the extended domain configuration file can be used with all
211the original input files for ORCA1 if the ice cavities are not active (\np{ln\_isfcav =
212.false.}).  Full instructions for achieving this are:
215\item Add the new attribute to any input files requiring a j-row offset, i.e:
217ncatted  -a open_ocean_jstart,global,a,d,41
220\item Add the logical switch \np{ln_use_jattr}{ln\_use\_jattr} to \nam{cfg}{cfg} in the configuration
221namelist (if it is not already there) and set \forcode{.true.}
224\noindent Note that with this option, the j-size of the global domain is (extended
225j-size minus \np{open_ocean_jstart}{open\_ocean\_jstart} + 1 ) and this must match the \texttt{jpjglo} value
226for the configuration. This means an alternative version of \textit{eORCA1\} must
227be created for when \np{ln_use_jattr}{ln\_use\_jattr} is active. The \texttt{ncap2} tool provides a
228convenient way of achieving this:
231ncap2 -s 'jpjglo=292'
234The domain configuration file is unique in this respect since it also contains the value of \texttt{jpjglo}
235that is read and used by the model.
236Any other global, 2D and 3D, netcdf, input field can be prepared for use in a reduced domain by adding the
237\texttt{open\_ocean\_jstart} attribute to the file's global attributes.
238In particular this is true for any field that is read by \NEMO\ using the following optional argument to
239the appropriate call to \np{iom_get}{iom\_get}.
245Currently, only the domain configuration variables make use of this optional argument so
246this facility is of little practical use except for tests where no other external input
247files are needed or you wish to use an extended domain configuration with inputs from
248earlier, non-extended configurations. Alternatively, it should be possible to exclude
249empty rows for extended domain, forced ocean runs using interpolation on the fly, by
250adding the optional argument to \texttt{iom\_get} calls for the weights and initial
251conditions. Experimenting with this remains an exercise for the user.
253%% =================================================================================================
254\section[Accuracy and reproducibility (\textit{lib\_fortran.F90})]{Accuracy and reproducibility (\protect\mdl{lib\_fortran})}
257%% =================================================================================================
258\subsection[Issues with intrinsinc SIGN function (\texttt{\textbf{key\_nosignedzero}})]{Issues with intrinsinc SIGN function (\protect\key{nosignedzero})}
261The SIGN(A, B) is the \fortran\ intrinsic function delivers the magnitude of A with the sign of B.
262For example, SIGN(-3.0,2.0) has the value 3.0.
263The problematic case is when the second argument is zero, because, on platforms that support IEEE arithmetic,
264zero is actually a signed number.
265There is a positive zero and a negative zero.
267In \fninety, the processor was required always to deliver a positive result for SIGN(A, B) if B was zero.
268Nevertheless, in \fninety, the processor is allowed to do the correct thing and deliver ABS(A) when
269B is a positive zero and -ABS(A) when B is a negative zero.
270This change in the specification becomes apparent only when B is of type real, and is zero,
271and the processor is capable of distinguishing between positive and negative zero,
272and B is negative real zero.
273Then SIGN delivers a negative result where, under \fninety\ rules, it used to return a positive result.
274This change may be especially sensitive for the ice model,
275so we overwrite the intrinsinc function with our own function simply performing :   \\
276\verb?   IF( B >= 0.e0 ) THEN   ;   SIGN(A,B) = ABS(A)  ?    \\
277\verb?   ELSE                   ;   SIGN(A,B) =-ABS(A)     ?  \\
278\verb?   ENDIF    ? \\
279This feature can be found in \mdl{lib\_fortran} module and is effective when \key{nosignedzero} is defined.
280We use a CPP key as the overwritting of a intrinsic function can present performance issues with
281some computers/compilers.
283%% =================================================================================================
284\subsection{MPP reproducibility}
287The numerical reproducibility of simulations on distributed memory parallel computers is a critical issue.
288In particular, within \NEMO\ global summation of distributed arrays is most susceptible to rounding errors,
289and their propagation and accumulation cause uncertainty in final simulation reproducibility on
290different numbers of processors.
291To avoid so, based on \citet{he.ding_JS01} review of different technics,
292we use a so called self-compensated summation method.
293The idea is to estimate the roundoff error, store it in a buffer, and then add it back in the next addition.
295Suppose we need to calculate $b = a_1 + a_2 + a_3$.
296The following algorithm will allow to split the sum in two
297($sum_1 = a_{1} + a_{2}$ and $b = sum_2 = sum_1 + a_3$) with exactly the same rounding errors as
298the sum performed all at once.
300   sum_1 \ \  &= a_1 + a_2 \\
301   error_1     &= a_2 + ( a_1 - sum_1 ) \\
302   sum_2 \ \  &= sum_1 + a_3 + error_1 \\
303   error_2     &= a_3 + error_1 + ( sum_1 - sum_2 ) \\
304   b \qquad \ &= sum_2 \\
306An example of this feature can be found in \mdl{lib\_fortran} module.
307It is systematicallt used in glob\_sum function (summation over the entire basin excluding duplicated rows and
308columns due to cyclic or north fold boundary condition as well as overlap MPP areas).
309The self-compensated summation method should be used in all summation in i- and/or j-direction.
310See \mdl{closea} module for an example.
311Note also that this implementation may be sensitive to the optimization level.
313%% =================================================================================================
314\subsection{MPP scalability}
317The default method of communicating values across the north-fold in distributed memory applications (\key{mpp\_mpi})
318uses a \textsc{MPI\_ALLGATHER} function to exchange values from each processing region in
319the northern row with every other processing region in the northern row.
320This enables a global width array containing the top 4 rows to be collated on every northern row processor and then
321folded with a simple algorithm.
322Although conceptually simple, this "All to All" communication will hamper performance scalability for
323large numbers of northern row processors.
324From version 3.4 onwards an alternative method is available which only performs direct "Peer to Peer" communications
325between each processor and its immediate "neighbours" across the fold line.
326This is achieved by using the default \textsc{MPI\_ALLGATHER} method during initialisation to
327help identify the "active" neighbours.
328Stored lists of these neighbours are then used in all subsequent north-fold exchanges to
329restrict exchanges to those between associated regions.
330The collated global width array for each region is thus only partially filled but is guaranteed to
331be set at all the locations actually required by each individual for the fold operation.
332This alternative method should give identical results to the default \textsc{ALLGATHER} method and
333is recommended for large values of \np{jpni}{jpni}.
334The new method is activated by setting \np{ln_nnogather}{ln\_nnogather} to be true (\nam{mpp}{mpp}).
335The reproducibility of results using the two methods should be confirmed for each new,
336non-reference configuration.
338%% =================================================================================================
339\section{Model optimisation, control print and benchmark}
343  \nlst{namctl}
344  \caption{\forcode{&namctl}}
345  \label{lst:namctl}
348Options are defined through the  \nam{ctl}{ctl} namelist variables.
350%% =================================================================================================
351\subsection{Vector optimisation}
353\key{vectopt\_loop} enables the internal loops to collapse.
354This is very a very efficient way to increase the length of vector calculations and thus
355to speed up the model on vector computers.
357% Add here also one word on NPROMA technique that has been found useless, since compiler have made significant progress during the last decade.
359% Add also one word on NEC specific optimisation (Novercheck option for example)
361%% =================================================================================================
362\subsection{Status and debugging information output}
365NEMO can produce a range of text information output either: in the main output
366file (ocean.output) written by the normal reporting processor (narea == 1) or various
367specialist output files (e.g. layout.dat, run.stat, tracer.stat etc.). Some, for example
368run.stat and tracer.stat, contain globally collected values for which a single file is
369sufficient. Others, however, contain information that could, potentially, be different
370for each processing region. For computational efficiency, the default volume of text
371information produced is reduced to just a few files from the narea=1 processor.
373When more information is required for monitoring or debugging purposes, the various
374forms of output can be selected via the \np{sn_cfctl}{sn\_cfctl} structure. As well as simple
375on-off switches this structure also allows selection of a range of processors for
376individual reporting (where appropriate) and a time-increment option to restrict
377globally collected values to specified time-step increments.
379Most options within the structure are influenced by the top-level switches shown here
380with their default settings:
383   sn_cfctl%l_allon  = .FALSE.    ! IF T activate all options. If F deactivate all unless l_config is T
384     sn_cfctl%l_config = .TRUE.     ! IF .true. then control which reports are written with the following
387The first switch is a convenience option which can be used to switch on and off all
388sub-options. However, if it is false then switching off all sub-options is only done
389if \forcode{sn_cfctl%l\_config} is also false. Specifically, the logic is:
392  IF ( sn_cfctl%l_allon ) THEN
393    set all suboptions .TRUE.
394    and set procmin, procmax and procincr so that all regions are selected ([0,10000000,1], respectively)
395  ELSEIF ( sn_cfctl%l_config ) THEN
396    honour individual settings of the suboptions from the namelist
397  ELSE
398    set all suboptions .FALSE.
399  ENDIF
402Details of the suboptions follow but first an explanation of the stand-alone option:
403\forcode{sn_cfctl%l_glochk}.  This option modifies the action of the early warning checks
404carried out in \texttt{stpctl.F90}. These checks detect probable numerical instabilites
405by searching for excessive sea surface heights or velocities and salinity values
406outside a sensible physical range. If breaches are detected then the default behaviour
407is to locate and report the local indices of the grid-point in breach. These indices
408are included in the error message that precedes the model shutdown. When true,
409\forcode{sn_cfctl%l_glochk} modifies this action by performing a global location of
410the various minimum and maximum values and the global indices are reported. This has
411some value in locating the most severe error in cases where the first detected error
412may not be the worst culprit.
414\subsubsection{Control print suboptions}
416The options that can be individually selected fall into three categories:
418\begin{enumerate} \item{Time step progress information} This category includes
419\texttt{run.stat} and \texttt{tracer.stat} files which record certain physical and
420passive tracer metrics (respectively). Typical contents of \texttt{run.stat} include
421global maximums of ssh, velocity; and global minimums and maximums of temperature
422and salinity.  A netCDF version of \texttt{run.stat} (\texttt{}) is also
423produced with the same time-series data and this can easily be expanded to include
424extra monitoring information.  \texttt{tracer.stat} contains the volume-weighted
425average tracer value for each passive tracer. Collecting these metrics involves
426global communications and will impact on model efficiency so both these options are
427disabled by default by setting the respective options, \forcode{sn\_cfctl%runstat} and
428\forcode{sn\_cfctl%trcstat} to false. A compromise can be made by activating either or
429both of these options and setting the \forcode{sn\_cfctl%timincr} entry to an integer
430value greater than one. This increment determines the time-step frequency at which
431the global metrics are collected and reported.  This increment also applies to the
432time.step file which is otherwise updated every timestep.
433\item{One-time configuration information/progress logs}
435Some run-time configuration information and limited progress information is always
436produced by the first ocean process. This includes the \texttt{ocean.output} file
437which reports on all the namelist options read by the model and remains open to catch
438any warning or error messages generated during execution. A \texttt{layout.dat}
439file is also produced which details the MPI-decomposition used by the model. The
440suboptions: \forcode{sn\_cfctl%oceout} and \forcode{sn\_cfctl%layout} can be used
441to activate the creation of these files by all ocean processes.  For example,
442when \forcode{sn\_cfctl%oceout} is true all processors produce their own version of
443\texttt{ocean.output}.  All files, beyond the the normal reporting processor (narea == 1), are
444named with a \_XXXX extension to their name, where XXXX is a 4-digit area number (with
445leading zeros, if required). This is useful as a debugging aid since all processes can
446report their local conditions. Note though that these files are buffered on most UNIX
447systems so bug-hunting efforts using this facility should also utilise the \fortran:
449\forline|CALL FLUSH(numout)|
451statement after any additional write statements to ensure that file contents reflect
452the last model state. Associated with the \forcode{sn\_cfctl%oceout} option is the
453additional \forcode{sn\_cfctl%oasout} suboption. This does not activate its own output
454file but rather activates the writing of addition information regarding the OASIS
455configuration when coupling via oasis and the sbccpl routine. This information is
456written to any active \texttt{ocean.output} files.
457\item{Control sums of trends for debugging}
459NEMO includes an option for debugging reproducibility differences between
460a MPP and mono-processor runs.  This is somewhat dated and clearly only
461useful for this purpose when dealing with configurations that can be run
462on a single processor. The full details can be found in this report: \href{
464control print option in NEMO} The switches to activate production of the control sums
465of trends for either the physics or passive tracers are the \forcode{sn\_cfctl%prtctl}
466and \forcode{sn\_cfctl%prttrc} suboptions, respectively. Although, perhaps, of limited use for its
467original intention, the ability to produce these control sums of trends in specific
468areas provides another tool for diagnosing model behaviour.  If only the output from a
469select few regions is required then additional options are available to activate options
470for only a simple subset of processing regions. These are: \forcode{sn\_cfctl%procmin},
471\forcode{sn\_cfctl%procmax} and \forcode{sn\_cfctl%procincr} which can be used to specify
472the minimum and maximum active areas and the increment. The default values are set
473such that all regions will be active. Note this subsetting can also be used to limit
474which additional \texttt{ocean.output} and \texttt{layout.dat} files are produced if
475those suboptions are active.
480   sn_cfctl%l_glochk  = .false. ! Range sanity checks are local (F) or global (T). Set T for debugging only
481   sn_cfctl%l_allon   = .false. ! IF T activate all options. If F deactivate all unless l_config is T
482   sn_cfctl%l_config  = .true.  ! IF .true. then control which reports are written with the following
483   sn_cfctl%l_runstat = .false. ! switches and which areas produce reports with the proc integer settings.
484   sn_cfctl%l_trcstat = .false. ! The default settings for the proc integers should ensure
485   sn_cfctl%l_oceout  = .false. ! that  all areas report.
486   sn_cfctl%l_layout  = .false. !
487   sn_cfctl%l_prtctl  = .false. !
488   sn_cfctl%l_prttrc  = .false. !
489   sn_cfctl%l_oasout  = .false. !
490   sn_cfctl%procmin   = 0       ! Minimum area number for reporting [default:0]
491   sn_cfctl%procmax   = 1000000 ! Maximum area number for reporting [default:1000000]
492   sn_cfctl%procincr  = 1       ! Increment for optional subsetting of areas [default:1]
493   sn_cfctl%ptimincr  = 1       ! Timestep increment for writing time step progress info
Note: See TracBrowser for help on using the repository browser.