# New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
chap_misc.tex in NEMO/trunk/doc/latex/NEMO/subfiles – NEMO

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

Last change on this file since 10601 was 10601, checked in by acc, 4 years ago

Remove some redundant variables from in_out_manager.F90. Update namctl and text in manual. Make whitespace correction to OFF, SAO and SAS versions of nemogcm.F90

File size: 20.4 KB
Line
1\documentclass[../main/NEMO_manual]{subfiles}
2
3\begin{document}
4% ================================================================
5% Chapter --- Miscellaneous Topics
6% ================================================================
7\chapter{Miscellaneous Topics}
8\label{chap:MISC}
9
10\minitoc
11
12\newpage
13
14% ================================================================
15% Representation of Unresolved Straits
16% ================================================================
17\section{Representation of unresolved straits}
18\label{sec:MISC_strait}
19
20In climate modeling, it often occurs that a crucial connections between water masses is broken as
21the grid mesh is too coarse to resolve narrow straits.
22For example, coarse grid spacing typically closes off the Mediterranean from the Atlantic at
23the Strait of Gibraltar.
24In this case, it is important for climate models to include the effects of salty water entering the Atlantic from
25the Mediterranean.
26Likewise, it is important for the Mediterranean to replenish its supply of water from the Atlantic to
27balance the net evaporation occurring over the Mediterranean region.
28This problem occurs even in eddy permitting simulations.
29For example, in ORCA 1/4\deg several straits of the Indonesian archipelago (Ombai, Lombok...)
30are much narrow than even a single ocean grid-point.
31
32We describe briefly here the three methods that can be used in \NEMO to handle such improperly resolved straits.
33The first two consist of opening the strait by hand while ensuring that the mass exchanges through
34the strait are not too large by either artificially reducing the surface of the strait grid-cells or,
35locally increasing the lateral friction.
36In the third one, the strait is closed but exchanges of mass, heat and salt across the land are allowed.
37Note that such modifications are so specific to a given configuration that no attempt has been made to
38set them in a generic way.
39However, examples of how they can be set up is given in the ORCA 2\deg and 0.5\deg configurations.
40For example, for details of implementation in ORCA2, search: \texttt{IF( cp\_cfg == "orca" .AND. jp\_cfg == 2 )}
41
42% -------------------------------------------------------------------------------------------------------------
44% -------------------------------------------------------------------------------------------------------------
46\label{subsec:MISC_strait_hand}
47
48$\bullet$ reduced scale factor in the cross-strait direction to a value in better agreement with
49the true mean width of the strait (\autoref{fig:MISC_strait_hand}).
50This technique is sometime called "partially open face" or "partially closed cells".
51The key issue here is only to reduce the faces of $T$-cell
52(\ie change the value of the horizontal scale factors at $u$- or $v$-point) but not the volume of the $T$-cell.
53Indeed, reducing the volume of strait $T$-cell can easily produce a numerical instability at
54that grid point that would require a reduction of the model time step.
55The changes associated with strait management are done in \mdl{domhgr},
56just after the definition or reading of the horizontal scale factors.
57
58$\bullet$ increase of the viscous boundary layer thickness by local increase of the fmask value at the coast
59(\autoref{fig:MISC_strait_hand}).
60This is done in \mdl{dommsk} together with the setting of the coastal value of fmask (see  \autoref{sec:LBC_coast}).
61
62%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
63\begin{figure}[!tbp]
64  \begin{center}
65    \includegraphics[width=0.80\textwidth]{Fig_Gibraltar}
66    \includegraphics[width=0.80\textwidth]{Fig_Gibraltar2}
67    \caption{
68      \protect\label{fig:MISC_strait_hand}
69      Example of the Gibraltar strait defined in a $1^{\circ} \times 1^{\circ}$ mesh.
70      \textit{Top}: using partially open cells.
71      The meridional scale factor at $v$-point is reduced on both sides of the strait to account for
72      the real width of the strait (about 20 km).
73      Note that the scale factors of the strait $T$-point remains unchanged.
74      \textit{Bottom}: using viscous boundary layers.
75      The four fmask parameters along the strait coastlines are set to a value larger than 4,
76      \ie "strong" no-slip case (see \autoref{fig:LBC_shlat}) creating a large viscous boundary layer that
77      allows a reduced transport through the strait.
78    }
79  \end{center}
80\end{figure}
81%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
82
83%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
84\begin{figure}[!tbp]
85  \begin{center}
87    \caption{
89      Example of mask fields for the closea module. \textit{Left}: a
91      field. In this example, if ln\_closea is set to .true., the mean
92      freshwater flux over each of the American Great Lakes will be
93      set to zero, and the total residual for all the lakes, if
94      negative, will be put into the St Laurence Seaway in the area
95      shown.
96    }
97  \end{center}
98\end{figure}
99%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
100
101% ================================================================
102% Closed seas
103% ================================================================
104\section{Closed seas (\protect\mdl{closea})}
105\label{sec:MISC_closea}
106
107Some configurations include inland seas and lakes as ocean
108points. This is particularly the case for configurations that are
109coupled to an atmosphere model where one might want to include inland
110seas and lakes as ocean model points in order to provide a better
111bottom boundary condition for the atmosphere. However there is no
112route for freshwater to run off from the lakes to the ocean and this
113can lead to large drifts in the sea surface height over the lakes. The
114closea module provides options to either fill in closed seas and lakes
115at run time, or to set the net surface freshwater flux for each lake
116to zero and put the residual flux into the ocean.
117
118Prior to NEMO 4 the locations of inland seas and lakes was set via
119hardcoded indices for various ORCA configurations. From NEMO 4 onwards
120the inland seas and lakes are defined using mask fields in the
121domain configuration file. The options are as follows.
122
123\begin{enumerate}
124\item{{\bf No closea\_mask'' field is included in domain configuration
125  file.} In this case the closea module does nothing.}
126
127\item{{\bf A field called closea\_mask is included in the domain
128configuration file and ln\_closea=.false. in namelist namcfg.} In this
129case the inland seas defined by the closea\_mask field are filled in
130(turned to land points) at run time. That is every point in
131closea\_mask that is nonzero is set to be a land point.}
132
133\item{{\bf A field called closea\_mask is included in the domain
134configuration file and ln\_closea=.true. in namelist namcfg.} Each
135inland sea or group of inland seas is set to a positive integer value
137for an example). The net surface flux over each inland sea or group of
138inland seas is set to zero each timestep and the residual flux is
139distributed over the global ocean (ie. all ocean points where
141
143included in the domain configuration file and ln\_closea=.true. in
144namelist namcfg.} This option works as for option 3, except that if
145the net surface flux over an inland sea is negative (net
146precipitation) it is put into the ocean at specified runoff points. A
147net positive surface flux (net evaporation) is still spread over the
148global ocean. The mapping from inland seas to runoff points is defined
149by the closea\_mask\_rnf field. Each mapping is defined by a positive
150integer value for the inland sea(s) and the corresponding runoff
151points. An example is given in Figure
152\ref{fig:closea_mask_example}. If no mapping is provided for a
153particular inland sea then the residual is spread over the global
154ocean.}
155
157included in the domain configuration file and ln\_closea=.true. in
158namelist namcfg.} This option works the same as option 4 except that
159the nonzero net surface flux is sent to the ocean at the specified
160runoff points regardless of whether it is positive or negative. The
161mapping from inland seas to runoff points in this case is defined by
163\end{enumerate}
164
165There is a python routine to create the closea\_mask fields and append
166them to the domain configuration file in the utils/tools/DOMAINcfg directory.
167
168% ================================================================
169% Sub-Domain Functionality
170% ================================================================
171\section{Sub-domain functionality}
172\label{sec:MISC_zoom}
173
174\subsection{Simple subsetting of input files via NetCDF attributes}
175
176The extended grids for use with the under-shelf ice cavities will result in redundant rows around Antarctica if
177the ice cavities are not active.
178A simple mechanism for subsetting input files associated with the extended domains has been implemented to
179avoid the need to maintain different sets of input fields for use with or without active ice cavities.
180The existing 'zoom' options are overly complex for this task and marked for deletion anyway.
181This alternative subsetting operates for the j-direction only and works by optionally looking for and
182using a global file attribute (named: \np{open\_ocean\_jstart}) to determine the starting j-row for input.
183The use of this option is best explained with an example:
184consider an ORCA1 configuration using the extended grid bathymetry and coordinate files:
185\vspace{-10pt}
186\ifile{eORCA1\_bathymetry\_v2}
187\ifile{eORCA1\_coordinates}
188\noindent These files define a horizontal domain of 362x332.
189Assuming the first row with open ocean wet points in the non-isf bathymetry for this set is row 42
190(\fortran indexing) then the formally correct setting for \np{open\_ocean\_jstart} is 41.
191Using this value as the first row to be read will result in a 362x292 domain which is the same size as
192the original ORCA1 domain.
193Thus the extended coordinates and bathymetry files can be used with all the original input files for ORCA1 if
194the ice cavities are not active (\np{ln\_isfcav = .false.}).
195Full instructions for achieving this are:
196
197\noindent Add the new attribute to any input files requiring a j-row offset, i.e:
198\vspace{-10pt}
199\begin{cmds}
200ncatted  -a open_ocean_jstart,global,a,d,41 eORCA1_coordinates.nc
201ncatted  -a open_ocean_jstart,global,a,d,41 eORCA1_bathymetry_v2.nc
202\end{cmds}
203
204\noindent Add the logical switch to \ngn{namcfg} in the configuration namelist and set true:
205%--------------------------------------------namcfg--------------------------------------------------------
206
207\nlst{namcfg}
208%--------------------------------------------------------------------------------------------------------------
209
210\noindent Note the j-size of the global domain is the (extended j-size minus \np{open\_ocean\_jstart} + 1 ) and
211this must match the size of all datasets other than bathymetry and coordinates currently.
212However the option can be extended to any global, 2D and 3D, netcdf, input field by adding the:
213\vspace{-10pt}
214\begin{forlines}
215lrowattr=ln_use_jattr
216\end{forlines}
217optional argument to the appropriate \np{iom\_get} call and the \np{open\_ocean\_jstart} attribute to
218the corresponding input files.
219It remains the users responsibility to set \np{jpjdta} and \np{jpjglo} values in
220the \np{namelist\_cfg} file according to their needs.
221
222%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
223\begin{figure}[!ht]
224  \begin{center}
225    \includegraphics[width=0.90\textwidth]{Fig_LBC_zoom}
226    \caption{
227      \protect\label{fig:LBC_zoom}
228      Position of a model domain compared to the data input domain when the zoom functionality is used.
229    }
230  \end{center}
231\end{figure}
232%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
233
234
235% ================================================================
236% Accuracy and Reproducibility
237% ================================================================
238\section{Accuracy and reproducibility (\protect\mdl{lib\_fortran})}
239\label{sec:MISC_fortran}
240
241\subsection{Issues with intrinsinc SIGN function (\protect\key{nosignedzero})}
242\label{subsec:MISC_sign}
243
244The SIGN(A, B) is the \fortran intrinsic function delivers the magnitude of A with the sign of B.
245For example, SIGN(-3.0,2.0) has the value 3.0.
246The problematic case is when the second argument is zero, because, on platforms that support IEEE arithmetic,
247zero is actually a signed number.
248There is a positive zero and a negative zero.
249
250In \fninety, the processor was required always to deliver a positive result for SIGN(A, B) if B was zero.
251Nevertheless, in \fninety, the processor is allowed to do the correct thing and deliver ABS(A) when
252B is a positive zero and -ABS(A) when B is a negative zero.
253This change in the specification becomes apparent only when B is of type real, and is zero,
254and the processor is capable of distinguishing between positive and negative zero,
255and B is negative real zero.
256Then SIGN delivers a negative result where, under \fninety rules, it used to return a positive result.
257This change may be especially sensitive for the ice model,
258so we overwrite the intrinsinc function with our own function simply performing :   \\
259\verb?   IF( B >= 0.e0 ) THEN   ;   SIGN(A,B) = ABS(A)  ?    \\
260\verb?   ELSE                   ;   SIGN(A,B) =-ABS(A)     ?  \\
261\verb?   ENDIF    ? \\
262This feature can be found in \mdl{lib\_fortran} module and is effective when \key{nosignedzero} is defined.
263We use a CPP key as the overwritting of a intrinsic function can present performance issues with
264some computers/compilers.
265
266
267\subsection{MPP reproducibility}
268\label{subsec:MISC_glosum}
269
270The numerical reproducibility of simulations on distributed memory parallel computers is a critical issue.
271In particular, within NEMO global summation of distributed arrays is most susceptible to rounding errors,
272and their propagation and accumulation cause uncertainty in final simulation reproducibility on
273different numbers of processors.
274To avoid so, based on \citet{He_Ding_JSC01} review of different technics,
275we use a so called self-compensated summation method.
276The idea is to estimate the roundoff error, store it in a buffer, and then add it back in the next addition.
277
278Suppose we need to calculate $b = a_1 + a_2 + a_3$.
279The following algorithm will allow to split the sum in two
280($sum_1 = a_{1} + a_{2}$ and $b = sum_2 = sum_1 + a_3$) with exactly the same rounding errors as
281the sum performed all at once.
282\begin{align*}
283   sum_1 \ \  &= a_1 + a_2 \\
284   error_1     &= a_2 + ( a_1 - sum_1 ) \\
285   sum_2 \ \  &= sum_1 + a_3 + error_1 \\
286   error_2     &= a_3 + error_1 + ( sum_1 - sum_2 ) \\
287   b \qquad \ &= sum_2 \\
288\end{align*}
289An example of this feature can be found in \mdl{lib\_fortran} module.
290It is systematicallt used in glob\_sum function (summation over the entire basin excluding duplicated rows and
291columns due to cyclic or north fold boundary condition as well as overlap MPP areas).
292The self-compensated summation method should be used in all summation in i- and/or j-direction.
293See \mdl{closea} module for an example.
294Note also that this implementation may be sensitive to the optimization level.
295
296\subsection{MPP scalability}
297\label{subsec:MISC_mppsca}
298
299The default method of communicating values across the north-fold in distributed memory applications (\key{mpp\_mpi})
300uses a \textsc{MPI\_ALLGATHER} function to exchange values from each processing region in
301the northern row with every other processing region in the northern row.
302This enables a global width array containing the top 4 rows to be collated on every northern row processor and then
303folded with a simple algorithm.
304Although conceptually simple, this "All to All" communication will hamper performance scalability for
305large numbers of northern row processors.
306From version 3.4 onwards an alternative method is available which only performs direct "Peer to Peer" communications
307between each processor and its immediate "neighbours" across the fold line.
308This is achieved by using the default \textsc{MPI\_ALLGATHER} method during initialisation to
309help identify the "active" neighbours.
310Stored lists of these neighbours are then used in all subsequent north-fold exchanges to
311restrict exchanges to those between associated regions.
312The collated global width array for each region is thus only partially filled but is guaranteed to
313be set at all the locations actually required by each individual for the fold operation.
314This alternative method should give identical results to the default \textsc{ALLGATHER} method and
315is recommended for large values of \np{jpni}.
316The new method is activated by setting \np{ln\_nnogather} to be true ({\bf nammpp}).
317The reproducibility of results using the two methods should be confirmed for each new,
318non-reference configuration.
319
320% ================================================================
321% Model optimisation, Control Print and Benchmark
322% ================================================================
323\section{Model optimisation, control print and benchmark}
324\label{sec:MISC_opt}
325%--------------------------------------------namctl-------------------------------------------------------
326
327\nlst{namctl}
328%--------------------------------------------------------------------------------------------------------------
329
330Options are defined through the  \ngn{namctl} namelist variables.
331
332\subsection{Vector optimisation}
333
334\key{vectopt\_loop} enables the internal loops to collapse.
335This is very a very efficient way to increase the length of vector calculations and thus
336to speed up the model on vector computers.
337
338% Add here also one word on NPROMA technique that has been found useless, since compiler have made significant progress during the last decade.
339
340% Add also one word on NEC specific optimisation (Novercheck option for example)
341
342\subsection{Control print}
343
344The \np{ln\_ctl} switch was originally used as a debugging option in two modes:
345
346\begin{enumerate}
347\item{\np{ln\_ctl}: compute and print the trends averaged over the interior domain in all TRA, DYN, LDF and
348ZDF modules.
349This option is very helpful when diagnosing the origin of an undesired change in model results. }
350
351\item{also \np{ln\_ctl} but using the nictl and njctl namelist parameters to check the source of differences between
352mono and multi processor runs.}
353\end{enumerate}
354
355However, in recent versions it has also been used to force all processors to assume the
356reporting role. Thus when \np{ln\_ctl} is true all processors produce their own versions
357of files such as: ocean.output, layout.dat, etc.  All such files, beyond the the normal
358reporting processor (narea == 1), are named with a \_XXXX extension to their name, where
359XXXX is a 4-digit area number (with leading zeros, if required). Other reporting files
360such as run.stat (and its netCDF counterpart: run.stat.nc) and tracer.stat contain global
361information and are only ever produced by the reporting master (narea == 1). For version
3624.0 a start has been made to return \np{ln\_ctl} to its original function by introducing
363a new control structure which allows finer control over which files are produced. This
364feature is still evolving but it does already allow the user to: select individually the
365production of run.stat and tracer.stat files and to toggle the production of other files
366on processors other than the reporting master. These other reporters can be a simple
367subset of processors as defined by a minimum, maximum and incremental processor number.
368
369Note, that production of the run.stat and tracer.stat files require global communications.
370For run.stat, these are global min and max operations to find metrics such as the gloabl
371maximum velocity. For tracer.stat these are global sums of tracer fields. To improve model
372performance these operations are disabled by default and, where necessary, any use of the
373global values have been replaced with local calculations. For example, checks on the CFL
374criterion are now done on the local domain and only reported if a breach is detected.
375
376Experienced users may wish to still monitor this information as a check on model progress.
377If so, the best compromise will be to activate the files with:
378
379\begin{verbatim}
380     sn_cfctl%l_config = .TRUE.
381       sn_cfctl%l_runstat = .TRUE.
382       sn_cfctl%l_trcstat = .TRUE.
383\end{verbatim}
384
385and to use the new time increment setting to ensure the values are collected and reported
386at a suitably long interval. For example:
387
388\begin{verbatim}
389       sn_cfctl%ptimincr  = 25
390\end{verbatim}
391
392will carry out the global communications and write the information every 25 timesteps. This
393increment also applies to the time.step file which is otherwise updated every timestep.
394
395% ================================================================
396\biblio
397
398\pindex
399
400\end{document}
Note: See TracBrowser for help on using the repository browser.