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_DIA.tex in branches/2012/dev_UKMO_2012/DOC/TexFiles/Chapters – NEMO

source: branches/2012/dev_UKMO_2012/DOC/TexFiles/Chapters/Chap_DIA.tex @ 3602

Last change on this file since 3602 was 3602, checked in by rfurner, 11 years ago

Changes from dev_r3452_UKMO2_DIADCT, revision 3452 to 3530 merged in

File size: 54.1 KB
Line 
1% ================================================================
2% Chapter Ñ I/O & Diagnostics
3% ================================================================
4\chapter{Ouput and Diagnostics (IOM, DIA, TRD, FLO)}
5\label{DIA}
6\minitoc
7
8\newpage
9$\ $\newline    % force a new ligne
10
11% ================================================================
12%       Old Model Output
13% ================================================================
14\section{Old Model Output (default or \key{dimgout})}
15\label{DIA_io_old}
16
17The model outputs are of three types: the restart file, the output listing,
18and the output file(s). The restart file is used internally by the code when
19the user wants to start the model with initial conditions defined by a
20previous simulation. It contains all the information that is necessary in
21order for there to be no changes in the model results (even at the computer
22precision) between a run performed with several restarts and the same run
23performed in one step. It should be noted that this requires that the restart file
24contain two consecutive time steps for all the prognostic variables, and
25that it is saved in the same binary format as the one used by the computer
26that is to read it (in particular, 32 bits binary IEEE format must not be used for
27this file). The output listing and file(s) are predefined but should be checked
28and eventually adapted to the user's needs. The output listing is stored in
29the $ocean.output$ file. The information is printed from within the code on the
30logical unit $numout$. To locate these prints, use the UNIX command
31"\textit{grep -i numout}" in the source code directory.
32
33In the standard configuration, the user will find the model results in
34NetCDF files containing mean values (or instantaneous values if
35\key{diainstant} is defined) for every time-step where output is demanded.
36These outputs are defined in the \mdl{diawri} module.
37When defining \key{dimgout}, the output are written in DIMG format,
38an IEEE output format.
39
40Since version 3.2, an I/O server has been added which provides more
41flexibility in the choice of the fields to be output as well as how the
42writing work is distributed over the processors in massively parallel
43computing. It is presented in next section.
44
45
46%\gmcomment{                    % start of gmcomment
47
48% ================================================================
49% Diagnostics
50% ================================================================
51\section{Standard model Output (IOM)}
52\label{DIA_iom}
53
54
55Since version 3.2, iom\_put is the NEMO output interface. It was designed to be simple to use,
56flexible and efficient. Two main functionalities are covered by iom\_put:
57(1) the control of the output files through an external xml file defined by the user ;
58(2) the distribution (or not) of all task related to output files on dedicated processors.
59The first functionality allows the user to specify, without touching anything into the code,
60the way he want to output data: \\
61- choice of output frequencies that can be different for each file (including real months and years) \\
62- choice of file contents: decide which data will be written in which file (the same data can be
63outputted in different files)  \\
64- possibility to extract a subdomain (for example all TAO-PIRATA-RAMA moorings are already defined)  \\
65- choice of the temporal operation to perform: mean, instantaneous, min, max  \\
66- extremely large choice of data available   \\
67- redefine variables name and long\_name  \\
68In addition, iom\_put allows the user to output any variable (scalar, 2D or 3D) in the code
69in a very easy way. All details of iom\_put functionalities are listed in the following subsections.
70An example of the iodef.xml file that control the outputs can be found here:
71NEMOGCM/CONFIG/ORCA2\_LIM/EXP00/iodef.xml
72
73The second functionality targets outputs performances when running on a very large number of processes.
74The idea is to dedicate N specific processes to write the outputs, where N is defined by the user.
75In the current version, this functionality is technically working however, its performance are usually poor
76(for known reasons). Users can therefore test this functionality but they must be aware that expected
77performance improvement will not be achieved before the release 3.4.
78An example of xmlio\_server.def NEMOGCM/CONFIG/ORCA2\_LIM/EXP00/xmlio\_server.def
79 
80
81\subsection{Basic knowledge}
82
83
84\subsubsection{ XML basic rules}
85
86XML tags begin with the less-than character ("$<$") and end with the greater-than character (''$>$'').
87You use tags to mark the start and end of elements, which are the logical units of information
88in an XML document. In addition to marking the beginning of an element, XML start tags also
89provide a place to specify attributes. An attribute specifies a single property for an element,
90using a name/value pair, for example: $<$a b="x" c="y" b="z"$>$ ... $<$/a$>$.
91See \href{http://www.xmlnews.org/docs/xml-basics.html}{here} for more details.
92
93\subsubsection{Structure of the xml file used in NEMO}
94
95The xml file is split into 3 parts:
96\begin{description}
97\item[field definition]: define all variables that can be output \\
98all lines in between the following two tags\\
99\verb?   <field\_definition ...> ?  \\
100\verb?   </field\_definition ...> ?
101\item[file definition]: define the netcdf files to be created and the variables they will contain \\
102all lines in between the following two tags \\
103\verb?   <field\_definition> ?  \\
104\verb?   </field\_definition> ?
105\item[axis and grid definitions]: define the horizontal and vertical grids \\
106all lines in between the following two set of two tags\\
107\verb?   <axis\_definition ...> ?  \\
108\verb?   </axis\_definition ...> ?
109and \\
110\verb?   <grid\_definition ...> ?  \\
111\verb?   </grid\_definition ...> ?
112\end{description}
113
114\subsubsection{Inheritance and group }
115
116 Xml extensively uses the concept of inheritance. \\
117\\
118example 1: \\
119\vspace{-30pt}
120\begin{alltt}  {{\scriptsize   
121\begin{verbatim}
122   <field_definition operation="ave(X)" >
123      <field id="sst"                    />   <!-- averaged      sst -->
124      <field id="sss" operation="inst(X)"/>   <!-- instantaneous sss -->
125   </field_definition>
126\end{verbatim}
127}}\end{alltt} 
128
129The field ''sst'' which is part (or a child) of the field\_definition will inherit the value ''ave(X)''
130of the attribute ''operation'' from its parent ''field definition''. Note that a child can overwrite
131the attribute definition inherited from its parents. In the example above, the field ''sss'' will
132therefore output instantaneous values instead of average values.
133
134example 2: Use (or overwrite) attributes value of a field when listing the variables included in a file
135\vspace{-20pt}
136\begin{alltt}  {{\scriptsize
137\begin{verbatim}
138   <field_definition>
139      <field id="sst" description="sea surface temperature" />   
140      <field id="sss" description="sea surface salinity"    /> 
141   </field_definition>     
142
143   <file_definition>
144      <file id="file_1" />   
145            <field ref="sst"                              />  <!-- default def -->
146            <field ref="sss" description="my description" />  <!-- overwrite   -->
147      </file>   
148   </file_definition>
149\end{verbatim}
150}}\end{alltt} 
151
152With the help of the inheritance, the concept of group allow to define a set of attributes
153for several fields or files.
154
155example 3, group of fields: define a group ''T\_grid\_variables'' identified with the name
156''grid\_T''. By default variables of this group have no vertical axis but, following inheritance
157rules, ''axis\_ref'' can be redefined for the field ''toce'' that is a 3D variable.
158\vspace{-30pt}
159\begin{alltt}  {{\scriptsize
160\begin{verbatim}
161   <field_definition>
162      <group id="grid_T" axis_ref="none" grid_ref="T_grid_variables">
163            <field id="sst"/> 
164            <field id="sss"/> 
165            <field id="toce" axis_ref="deptht"/>  <!-- overwrite axis def -->
166      </group>
167   </field_definition>
168\end{verbatim}
169}}\end{alltt} 
170
171example 4, group of files: define a group of file with the attribute output\_freq equal to 432000 (5 days)
172\vspace{-30pt}
173\begin{alltt}  {{\scriptsize
174\begin{verbatim}
175   <file_definition>
176      <group id="5d" output_freq="432000">    <!-- 5d files -->
177         <file id="5d_grid_T" name="auto">   <!-- T grid file -->
178         ...
179         </file>
180         <file id="5d_grid_U" name="auto">   <!-- U grid file -->
181         ...
182         </file>
183      </group>
184   </file_definition>
185\end{verbatim}
186}}\end{alltt} 
187
188\subsubsection{Control of the xml attributes from NEMO}
189
190The values of some attributes are automatically defined by NEMO (and any definition
191given in the xml file is overwritten). By convention, these attributes are defined to ''auto''
192(for string) or ''0000'' (for integer) in the xml file (but this is not necessary).
193
194Here is the list of these attributes: \\
195\\
196\begin{tabular}{|l|c|c|c|}
197   \hline
198 \multicolumn{2}{|c|}{tag ids affected by automatic           }  & name      & attribute value \\
199  \multicolumn{2}{|c|}{definition of some of their attributes }  & attribute  &        \\
200   \hline
201   \hline
202    \multicolumn{2}{|c|}{field\_definition} & freq\_op & \np{rn\_rdt} \\
203   \hline
204    \multicolumn{2}{|c|}{SBC}                  & freq\_op & \np{rn\_rdt} $\times$ \np{nn\_fsbc} \\
205   \hline
206    1h, 2h, 3h, 4h, 6h, 12h & \_grid\_T, \_grid\_U,  &  name & filename defined by   \\
207    1d, 3d, 5d                     & \_grid\_V, \_grid\_W, &            & a call to rou{dia\_nam}  \\
208    1m, 2m, 3m, 4m, 6m    & \_icemod, \_ptrc\_T,  &            & following NEMO \\
209    1y, 2y, 5y, 10y               &  \_diad\_T, \_scalar   &            & nomenclature \\
210   \hline
211    \multicolumn{2}{|c|}{EqT, EqU, EqW} & jbegin, ni,      & according to the grid    \\
212    \multicolumn{2}{|c|}{                         } & name\_suffix &                                      \\
213   \hline
214    \multicolumn{2}{|c|}{TAO, RAMA and PIRATA moorings} & ibegin, jbegin,      & according to the grid    \\
215    \multicolumn{2}{|c|}{                                                       } & name\_suffix &                                      \\
216   \hline
217\end{tabular}
218
219
220\subsection{ Detailed functionalities }
221
222\subsubsection{Tag list}
223
224\begin{description}
225
226\item[context]: define the model using the xml file. Id is the only attribute accepted.
227Its value must be ''nemo'' or ''n\_nemo'' for the nth AGRIF zoom. Child of simulation tag.
228
229\item[field]: define the field to be output. Accepted attributes are axis\_ref, description, enable,
230freq\_op, grid\_ref, id (if child of field\_definition), level, operation, name, ref (if child of file),
231unit, zoom\_ref. Child of field\_definition, file or group of fields tag.
232
233\item[field\_definition]: definition of the part of the xml file corresponding to the field definition.
234Accept the same attributes as field tag. Child of context tag.
235
236\item[group]: define a group of file or field. Accept the same attributes as file or field.
237
238\item[file]: define the output file's characteristics. Accepted attributes are description, enable,
239output\_freq, output\_level, id, name, name\_suffix. Child of file\_definition or group of files tag.
240
241\item[file\_definition]: definition of the part of the xml file corresponding to the file definition.
242Accept the same attributes as file tag. Child of context tag.
243
244\item[axis]: definition of the vertical axis. Accepted attributes are description, id, positive, size, unit.
245Child of axis\_definition tag.
246
247\item[axis\_definition]: definition of the part of the xml file corresponding to the vertical axis definition.
248Accept the same attributes as axis tag. Child of context tag
249
250\item[grid]: definition of the horizontal grid. Accepted attributes are description and id.
251Child of axis\_definition tag.
252
253\item[grid\_definition]: definition of the part of the xml file corresponding to the horizontal grid definition.
254Accept the same attributes as grid tag. Child of context tag
255
256\item[zoom]: definition of a subdomain of an horizontal grid. Accepted attributes are description, id,
257i/jbegin, ni/j. Child of grid tag.
258
259\end{description}
260
261
262\subsubsection{Attributes list}
263
264Applied to a tag or a group of tags.
265
266% table to be added ?
267Another table, perhaps?
268
269%%%%
270
271Attribute
272Applied to?
273Definition
274Comment
275axis\_ref
276field
277String defining the vertical axis of the variable. It refers to the id of the vertical axis defined in the axis tag.
278Use ''non'' if the variable has no vertical axis
279
280%%%%%%
281
282\begin{description}
283
284\item[axis\_ref]: field attribute. String defining the vertical axis of the variable.
285It refers to the id of the vertical axis defined in the axis tag.
286Use ''none'' if the variable has no vertical axis
287
288\item[description]: this attribute can be applied to all tags but it is used only with the field tag.
289In this case, the value of description will be used to define, in the output netcdf file,
290the attributes long\_name and standard\_name of the variable.
291
292\item[enabled]: field and file attribute. Logical to switch on/off the output of a field or a file.
293
294\item[freq\_op]: field attribute (automatically defined, see part 1.4 (''control of the xml attributes'')).
295An integer defining the frequency in seconds at which NEMO is calling iom\_put for this variable.
296It corresponds to the model time step (rn\_rdt in the namelist) except for the variables computed
297at the frequency of the surface boundary condition (rn\_rdt ? nn\_fsbc in the namelist).   
298
299\item[grid\_ref]: field attribute. String defining the horizontal grid of the variable.
300It refers to the id of the grid tag.
301
302\item[ibegin]: zoom attribute. Integer defining the zoom starting point along x direction.
303Automatically defined for TAO/RAMA/PIRATA moorings (see part 1.4). 
304
305\item[id]: exists for all tag. This is a string defining the name to a specific tag that will be used
306later to refer to this tag. Tags of the same category must have different ids.
307
308\item[jbegin]: zoom attribute. Integer defining the zoom starting point along y direction.
309Automatically defined for TAO/RAMA/PIRATA moorings and equatorial section (see part 1.4).
310
311\item[level]: field attribute. Integer from 0 to 10 defining the output priority of a field.
312See output\_level attribute definition
313
314\item[operation]: field attribute. String defining the type of temporal operation to perform on a variable.
315Possible choices are ''ave(X)'' for temporal mean, ''inst(X)'' for instantaneous, ''t\_min(X)'' for temporal min
316and ''t\_max(X)'' for temporal max.
317
318\item[output\_freq]: file attribute. Integer defining the operation frequency in seconds.
319For example 86400 for daily mean.
320
321\item[output\_level]: file attribute. Integer from 0 to 10 defining the output priority of variables in a file:
322all variables listed in the file with a level smaller or equal to output\_level will be output.
323Other variables won't be output even if they are listed in the file. 
324
325\item[positive]: axis attribute (always .FALSE.). Logical defining the vertical axis convention used
326in \NEMO (positive downward). Define the attribute positive of the variable in the netcdf output file.
327
328\item[prec]: field attribute. Integer defining the output precision.
329Not implemented, we always output real4 arrays.
330
331\item[name]: field or file attribute. String defining the name of a variable or a file.
332If the name of a file is undefined, its id is used as a name. 2 files must have different names.
333Files with specific ids will have their name automatically defined (see part 1.4).
334Note that is name will be automatically completed by the cpu number (if needed) and ''.nc''
335
336\item[name\_suffix]: file attribute. String defining a suffix to be inserted after the name
337and before the cpu number and the ''.nc'' termination. Files with specific ids have an
338automatic definition of their suffix (see part 1.4).
339
340\item[ni]: zoom attribute. Integer defining the zoom extent along x direction.
341Automatically defined for equatorial sections (see part 1.4). 
342
343\item[nj]: zoom attribute. Integer defining the zoom extent along x direction.
344
345\item[ref]: field attribute. String referring to the id of the field we want to add in a file.
346
347\item[size]: axis attribute. use unknown...
348
349\item[unit]: field attribute. String defining the unit of a variable and the associated
350attribute in the netcdf output file.
351
352\item[zoom\_ref]: field attribute. String defining the subdomain of data on which
353the file should be written (to ouput data only in a limited area).
354It refers to the id of a zoom defined in the zoom tag.
355\end{description}
356
357
358\subsection{IO\_SERVER}
359
360\subsubsection{Attached or detached mode?}
361
362Iom\_put is based on the io\_server developed by Yann Meurdesoif from IPSL
363(see \href{http://forge.ipsl.jussieu.fr/ioserver/browser}{here} for the source code or
364see its copy in NEMOGCM/EXTERNAL directory).
365This server can be used in ''attached mode'' (as a library) or in ''detached mode''
366(as an external executable on n cpus). In attached mode, each cpu of NEMO will output
367its own subdomain. In detached mode, the io\_server will gather data from NEMO
368and output them split over n files with n the number of cpu dedicated to the io\_server.
369
370\subsubsection{Control the io\_server: the namelist file xmlio\_server.def}
371
372%
373%Again, a small table might be more readable?
374%Name
375%Type
376%Description
377%Comment
378%Using_server
379%Logical
380%Switch to use the server in attached or detached mode
381%(.TRUE. corresponding to detached mode).
382
383The control of the use of the io\_server is done through the namelist file of the io\_server
384called xmlio\_server.def.
385
386\textbf{using\_server}: logical, switch to use the server in attached or detached mode
387(.TRUE. corresponding to detached mode).
388
389\textbf{using\_oasis}: logical, set to .TRUE. if NEMO is used in coupled mode.
390
391\textbf{client\_id} = ''oceanx'' : character, used only in coupled mode.
392Specify the id used in OASIS to refer to NEMO. The same id must be used to refer to NEMO
393in the \$NBMODEL part of OASIS namcouple in the call of prim\_init\_comp\_proto in cpl\_oasis3f90
394
395\textbf{server\_id} = ''ionemo'' : character, used only in coupled mode.
396Specify the id used in OASIS to refer to the IO\_SERVER when used in detached mode.
397Use the same id to refer to the io\_server in the \$NBMODEL part of OASIS namcouple.
398
399\textbf{global\_mpi\_buffer\_size}: integer; define the size in Mb of the MPI buffer used by the io\_server.
400
401\subsubsection{Number of cpu used by the io\_server in detached mode}
402
403The number of cpu used by the io\_server is specified only when launching the model.
404Here is an example of 2 cpus for the io\_server and 6 cpu for opa using mpirun:
405
406\texttt{ -p 2 -e ./ioserver}
407
408\texttt{ -p 6 -e ./opa }
409
410
411\subsection{Practical issues}
412
413\subsubsection{Add your own outputs}
414
415It is very easy to add you own outputs with iom\_put. 4 points must be followed.
416\begin{description}
417\item[1-] in NEMO code, add a \\
418\texttt{      CALL iom\_put( 'identifier', array ) } \\
419where you want to output a 2D or 3D array.
420
421\item[2-] don't forget to add \\
422\texttt{   USE iom            ! I/O manager library }  \\
423in the list of used modules in the upper part of your module.
424
425\item[3-] in the file\_definition part of the xml file, add the definition of your variable using the same identifier you used in the f90 code.
426\vspace{-20pt}
427\begin{alltt}  {{\scriptsize
428\begin{verbatim}
429   <field_definition>
430      ...
431      <field id="identifier" description="blabla" />   
432      ...
433   </field_definition>
434\end{verbatim}
435}}\end{alltt} 
436attributes axis\_ref and grid\_ref must be consistent with the size of the array to pass to iom\_put.
437if your array is computed within the surface module each nn\_fsbc time\_step,
438add the field definition within the group defined with the id ''SBC'': $<$group id=''SBC''...$>$
439
440\item[4-] add your field in one of the output files   \\
441\vspace{-20pt}
442\begin{alltt}  {{\scriptsize
443\begin{verbatim}
444   <file id="file_1" .../>   
445      ...
446      <field ref="identifier" />   
447      ...
448   </file>   
449\end{verbatim}
450}}\end{alltt} 
451
452\end{description}
453
454\subsubsection{Several time axes in the output file}
455
456If your output file contains variables with different operations (see operation definition),
457IOIPSL will create one specific time axis for each operation. Note that inst(X) will have
458a time axis corresponding to the end each output period whereas all other operators
459will have a time axis centred in the middle of the output periods.
460
461\subsubsection{Error/bug messages from IOIPSL}
462
463If you get the following error in the standard output file:
464\vspace{-20pt}
465\begin{alltt}  {{\scriptsize
466\begin{verbatim}
467FATAL ERROR FROM ROUTINE flio_dom_set
468 --> too many domains simultaneously defined
469 --> please unset useless domains
470 --> by calling flio_dom_unset
471\end{verbatim}
472}}\end{alltt} 
473
474You must increase the value of dom\_max\_nb in fliocom.f90 (multiply it by 10 for example).
475
476If you mix, in the same file, variables with different freq\_op (see definition above),
477like for example variables from the surface module with other variables,
478IOIPSL will print in the standard output file warning messages saying there may be a bug.
479\vspace{-20pt}
480\begin{alltt}  {{\scriptsize
481\begin{verbatim}
482WARNING FROM ROUTINE histvar_seq   
483 --> There were 10 errors in the learned sequence of variables 
484 --> for file   4
485 --> This looks like a bug, please report it.
486\end{verbatim}
487}}\end{alltt} 
488
489Don't worry, there is no bug, everything is properly working!
490
491 %    }      %end  \gmcomment
492
493
494% ================================================================
495%       NetCDF4 support
496% ================================================================
497\section{NetCDF4 Support (\key{netcdf4})}
498\label{DIA_iom}
499
500Since version 3.3, support for NetCDF4 chunking and (loss-less) compression has
501been included.  These options build on the standard NetCDF output and allow
502the user control over the size of the chunks via namelist settings. Chunking
503and compression can lead to significant reductions in file sizes for a small
504runtime overhead. For a fuller discussion on chunking and other performance
505issues the reader is referred to the NetCDF4 documentation found
506\href{http://www.unidata.ucar.edu/software/netcdf/docs/netcdf.html#Chunking}{here}.
507
508The new features are only available when the code has been linked with a
509NetCDF4 library (version 4.1 onwards, recommended) which has been built
510with HDF5 support (version 1.8.4 onwards, recommended). Datasets created
511with chunking and compression are not backwards compatible with NetCDF3
512"classic" format but most analysis codes can be relinked simply with the
513new libraries and will then read both NetCDF3 and NetCDF4 files. NEMO
514executables linked with NetCDF4 libraries can be made to produce NetCDF3
515files by setting the \np{ln\_nc4zip} logical to false in the \textit{namnc4} 
516namelist:
517
518%------------------------------------------namnc4----------------------------------------------------
519\namdisplay{namnc4}
520%-------------------------------------------------------------------------------------------------------------
521
522If \key{netcdf4} has not been defined, these namelist parameters are not read.
523In this case, \np{ln\_nc4zip} is set false and dummy routines for a few
524NetCDF4-specific functions are defined. These functions will not be used but
525need to be included so that compilation is possible with NetCDF3 libraries.
526
527When using NetCDF4 libraries, \key{netcdf4} should be defined even if the
528intention is to create only NetCDF3-compatible files. This is necessary to
529avoid duplication between the dummy routines and the actual routines present
530in the library. Most compilers will fail at compile time when faced with
531such duplication. Thus when linking with NetCDF4 libraries the user must
532define \key{netcdf4} and control the type of NetCDF file produced via the
533namelist parameter.
534
535Chunking and compression is applied only to 4D fields and there is no
536advantage in chunking across more than one time dimension since previously
537written chunks would have to be read back and decompressed before being
538added to. Therefore, user control over chunk sizes is provided only for the
539three space dimensions. The user sets an approximate number of chunks along
540each spatial axis. The actual size of the chunks will depend on global domain
541size for mono-processors or, more likely, the local processor domain size for
542distributed processing. The derived values are subject to practical minimum
543values (to avoid wastefully small chunk sizes) and cannot be greater than the
544domain size in any dimension. The algorithm used is:
545
546\begin{alltt}  {{\scriptsize 
547\begin{verbatim}
548     ichunksz(1) = MIN( idomain_size,MAX( (idomain_size-1)/nn_nchunks_i + 1 ,16 ) )
549     ichunksz(2) = MIN( jdomain_size,MAX( (jdomain_size-1)/nn_nchunks_j + 1 ,16 ) )
550     ichunksz(3) = MIN( kdomain_size,MAX( (kdomain_size-1)/nn_nchunks_k + 1 , 1 ) )
551     ichunksz(4) = 1
552\end{verbatim}
553}}\end{alltt} 
554
555\noindent As an example, setting:
556\vspace{-20pt}
557\begin{alltt}  {{\scriptsize
558\begin{verbatim}
559     nn_nchunks_i=4, nn_nchunks_j=4 and nn_nchunks_k=31
560\end{verbatim}
561}}\end{alltt} \vspace{-10pt}
562
563\noindent for a standard ORCA2\_LIM configuration gives chunksizes of {\small\tt 46x38x1}
564respectively in the mono-processor case (i.e. global domain of {\small\tt 182x149x31}).
565An illustration of the potential space savings that NetCDF4 chunking and compression
566provides is given in table \ref{Tab_NC4} which compares the results of two short
567runs of the ORCA2\_LIM reference configuration with a 4x2 mpi partitioning. Note
568the variation in the compression ratio achieved which reflects chiefly the dry to wet
569volume ratio of each processing region.
570
571%------------------------------------------TABLE----------------------------------------------------
572\begin{table}  \begin{tabular}{lrrr}
573Filename & NetCDF3 & NetCDF4 & Reduction\\
574         &filesize & filesize & \% \\
575         &(KB)     & (KB)     & \\
576ORCA2\_restart\_0000.nc & 16420 & 8860 & 47\%\\
577ORCA2\_restart\_0001.nc & 16064 & 11456 & 29\%\\
578ORCA2\_restart\_0002.nc & 16064 & 9744 & 40\%\\
579ORCA2\_restart\_0003.nc & 16420 & 9404 & 43\%\\
580ORCA2\_restart\_0004.nc & 16200 & 5844 & 64\%\\
581ORCA2\_restart\_0005.nc & 15848 & 8172 & 49\%\\
582ORCA2\_restart\_0006.nc & 15848 & 8012 & 50\%\\
583ORCA2\_restart\_0007.nc & 16200 & 5148 & 69\%\\
584ORCA2\_2d\_grid\_T\_0000.nc & 2200 & 1504 & 32\%\\
585ORCA2\_2d\_grid\_T\_0001.nc & 2200 & 1748 & 21\%\\
586ORCA2\_2d\_grid\_T\_0002.nc & 2200 & 1592 & 28\%\\
587ORCA2\_2d\_grid\_T\_0003.nc & 2200 & 1540 & 30\%\\
588ORCA2\_2d\_grid\_T\_0004.nc & 2200 & 1204 & 46\%\\
589ORCA2\_2d\_grid\_T\_0005.nc & 2200 & 1444 & 35\%\\
590ORCA2\_2d\_grid\_T\_0006.nc & 2200 & 1428 & 36\%\\
591ORCA2\_2d\_grid\_T\_0007.nc & 2200 & 1148 & 48\%\\
592             ...            &  ... &  ... & ..  \\
593ORCA2\_2d\_grid\_W\_0000.nc & 4416 & 2240 & 50\%\\
594ORCA2\_2d\_grid\_W\_0001.nc & 4416 & 2924 & 34\%\\
595ORCA2\_2d\_grid\_W\_0002.nc & 4416 & 2512 & 44\%\\
596ORCA2\_2d\_grid\_W\_0003.nc & 4416 & 2368 & 47\%\\
597ORCA2\_2d\_grid\_W\_0004.nc & 4416 & 1432 & 68\%\\
598ORCA2\_2d\_grid\_W\_0005.nc & 4416 & 1972 & 56\%\\
599ORCA2\_2d\_grid\_W\_0006.nc & 4416 & 2028 & 55\%\\
600ORCA2\_2d\_grid\_W\_0007.nc & 4416 & 1368 & 70\%\\
601\end{tabular}
602\caption{   \label{Tab_NC4} 
603Filesize comparison between NetCDF3 and NetCDF4 with chunking and compression}
604\end{table}
605%----------------------------------------------------------------------------------------------------
606
607When \key{iomput} is activated with \key{netcdf4} chunking and
608compression parameters for fields produced via \np{iom\_put} calls are
609set via an equivalent and identically named namelist to \textit{namnc4} 
610in \np{xmlio\_server.def}. Typically this namelist serves the mean files
611whilst the \np{ namnc4} in the main namelist file continues to serve the
612restart files. This duplication is unfortunate but appropriate since, if
613using io\_servers, the domain sizes of the individual files produced by the
614io\_server processes may be different to those produced by the invidual
615processing regions and different chunking choices may be desired.
616 
617
618% -------------------------------------------------------------------------------------------------------------
619%       Tracer/Dynamics Trends
620% -------------------------------------------------------------------------------------------------------------
621\section[Tracer/Dynamics Trends (TRD)]
622                  {Tracer/Dynamics Trends  (\key{trdtra}, \key{trddyn},    \\ 
623                                                             \key{trddvor}, \key{trdmld})}
624\label{DIA_trd}
625
626%------------------------------------------namtrd----------------------------------------------------
627\namdisplay{namtrd} 
628%-------------------------------------------------------------------------------------------------------------
629
630When \key{trddyn} and/or \key{trddyn} CPP variables are defined, each
631trend of the dynamics and/or temperature and salinity time evolution equations
632is stored in three-dimensional arrays just after their computation ($i.e.$ at the end
633of each $dyn\cdots.F90$ and/or $tra\cdots.F90$ routines). These trends are then
634used in \mdl{trdmod} (see TRD directory) every \textit{nn\_trd } time-steps.
635
636What is done depends on the CPP keys defined:
637\begin{description}
638\item[\key{trddyn}, \key{trdtra}] : a check of the basin averaged properties of the momentum
639and/or tracer equations is performed ;
640\item[\key{trdvor}] : a vertical summation of the moment tendencies is performed,
641then the curl is computed to obtain the barotropic vorticity tendencies which are output ;
642\item[\key{trdmld}] : output of the tracer tendencies averaged vertically 
643either over the mixed layer (\np{nn\_ctls}=0),
644or       over a fixed number of model levels (\np{nn\_ctls}$>$1 provides the number of level),
645or       over a spatially varying but temporally fixed number of levels (typically the base
646of the winter mixed layer) read in \ifile{ctlsurf\_idx} (\np{nn\_ctls}=1) ;
647\end{description}
648
649The units in the output file can be changed using the \np{nn\_ucf} namelist parameter.
650For example, in case of salinity tendency the units are given by PSU/s/\np{nn\_ucf}.
651Setting \np{nn\_ucf}=86400 ($i.e.$ the number of second in a day) provides the tendencies in PSU/d.
652
653When \key{trdmld} is defined, two time averaging procedure are proposed.
654Setting \np{ln\_trdmld\_instant} to \textit{true}, a simple time averaging is performed,
655so that the resulting tendency is the contribution to the change of a quantity between
656the two instantaneous values taken at the extremities of the time averaging period.
657Setting \np{ln\_trdmld\_instant} to \textit{false}, a double time averaging is performed,
658so that the resulting tendency is the contribution to the change of a quantity between
659two \textit{time mean} values. The later option requires the use of an extra file, \ifile{restart\_mld} 
660(\np{ln\_trdmld\_restart}=true), to restart a run.
661
662
663Note that the mixed layer tendency diagnostic can also be used on biogeochemical models
664via the \key{trdtrc} and \key{trdmld\_trc} CPP keys.
665
666% -------------------------------------------------------------------------------------------------------------
667%       On-line Floats trajectories
668% -------------------------------------------------------------------------------------------------------------
669\section{On-line Floats trajectories (FLO) (\key{floats})}
670\label{FLO}
671%--------------------------------------------namflo-------------------------------------------------------
672\namdisplay{namflo} 
673%--------------------------------------------------------------------------------------------------------------
674
675The on-line computation of floats advected either by the three dimensional velocity
676field or constraint to remain at a given depth ($w = 0$ in the computation) have been
677introduced in the system during the CLIPPER project. The algorithm used is based
678either on the work of \cite{Blanke_Raynaud_JPO97} (default option), or on a $4^th$
679Runge-Hutta algorithm (\np{ln\_flork4}=true). Note that the \cite{Blanke_Raynaud_JPO97} 
680algorithm have the advantage of providing trajectories which are consistent with the
681numeric of the code, so that the trajectories never intercept the bathymetry.
682
683\subsubsection{ Input data: initial coordinates }
684
685Initial coordinates can be given with Ariane Tools convention ( IJK coordinates ,(\np{ln\_ariane}=true) )
686or with longitude and latitude.
687
688
689In case of Ariane convention, input filename is \np{"init\_float\_ariane"}. Its format is:
690
691\texttt{ I J K nisobfl itrash itrash }
692
693\noindent with:
694
695 - I,J,K  : indexes of initial position
696
697 - nisobfl: 0 for an isobar float, 1 for a float following the w velocity 
698
699 - itrash : set to zero; it is a dummy variable to respect Ariane Tools convention
700
701 - itrash :set to zero; it is a dummy variable to respect Ariane Tools convention
702
703\noindent Example:\\
704\noindent \texttt{ 100.00000  90.00000  -1.50000 1.00000   0.00000}\\
705\texttt{ 102.00000  90.00000  -1.50000 1.00000   0.00000}\\
706\texttt{ 104.00000  90.00000  -1.50000 1.00000   0.00000}\\
707\texttt{ 106.00000  90.00000  -1.50000 1.00000   0.00000}\\
708\texttt{ 108.00000  90.00000  -1.50000 1.00000   0.00000}\\
709
710
711In the other case ( longitude and latitude ), input filename is \np{"init\_float"}. Its format is:
712
713\texttt{ Long Lat depth nisobfl ngrpfl itrash}
714
715\noindent with:
716
717 - Long, Lat, depth  : Longitude, latitude, depth
718
719 - nisobfl: 0 for an isobar float, 1 for a float following the w velocity
720
721 - ngrpfl : number to identify searcher group
722
723 - itrash :set to 1; it is a dummy variable.
724
725\noindent Example:
726
727\noindent\texttt{  20.0 0.0 0.0 0 1 1 }\\
728\texttt{ -21.0 0.0 0.0 0 1 1 }\\
729\texttt{ -22.0 0.0 0.0 0 1 1 }\\
730\texttt{ -23.0 0.0 0.0 0 1 1 }\\
731\texttt{ -24.0 0.0 0.0 0 1 1 }\\
732
733\np{jpnfl} is the total number of floats during the run.
734When initial positions are read in a restart file ( \np{ln\_rstflo= .TRUE.} ),  \np{jpnflnewflo}
735can be added in the initialization file.
736
737\subsubsection{ Output data }
738
739\np{nn\_writefl} is the frequency of writing in float output file and \np{nn\_stockfl} 
740is the frequency of creation of the float restart file.
741
742Output data can be written in ascii files (\np{ln\_flo\_ascii = .TRUE.} ). In that case,
743output filename is \np{is trajec\_float}.
744
745Another possiblity of writing format is Netcdf (\np{ln\_flo\_ascii = .FALSE.} ). There are 2 possibilities:
746
747 - if (\key{iomput}) is used, outputs are selected in  \np{iodef.xml}. Here it is an example of specification
748   to put in files description section:
749
750\vspace{-30pt}
751\begin{alltt}  {{\scriptsize
752\begin{verbatim}
753
754     <group id="1d\_grid\_T" name="auto" description="ocean T grid variables" >   }
755       <file id="floats"  description="floats variables"> }\\
756           <field ref="traj\_lon"   name="floats\_longitude"   freq\_op="86400" />}
757           <field ref="traj\_lat"   name="floats\_latitude"    freq\_op="86400" />}
758           <field ref="traj\_dep"   name="floats\_depth"       freq\_op="86400" />}
759           <field ref="traj\_temp"  name="floats\_temperature" freq\_op="86400" />}
760           <field ref="traj\_salt"  name="floats\_salinity"    freq\_op="86400" />}
761           <field ref="traj\_dens"  name="floats\_density"     freq\_op="86400" />}
762           <field ref="traj\_group" name="floats\_group"       freq\_op="86400" />}
763       </file>}
764  </group>}
765
766\end{verbatim}
767}}\end{alltt}
768
769
770 -  if (\key{iomput}) is not used, a file called \np{trajec\_float.nc} will be created by IOIPSL library.
771
772
773
774See also \href{http://stockage.univ-brest.fr/~grima/Ariane/}{here} the web site describing
775the off-line use of this marvellous diagnostic tool.
776
777
778% -------------------------------------------------------------------------------------------------------------
779%       Harmonic analysis of tidal constituents
780% -------------------------------------------------------------------------------------------------------------
781\section{Harmonic analysis of tidal constituents (\key{diaharm}) }
782\label{DIA_diag_harm}
783
784A module is available to compute the amplitude and phase for tidal waves.
785This diagnostic is actived with \key{diaharm}.
786
787%------------------------------------------namdia_harm----------------------------------------------------
788\namdisplay{namdia_harm}
789%----------------------------------------------------------------------------------------------------------
790
791Concerning the on-line Harmonic analysis, some parameters are available in namelist:
792
793- \texttt{nit000\_han} is the first time step used for harmonic analysis
794
795- \texttt{nitend\_han} is the last time step used for harmonic analysis
796
797- \texttt{nstep\_han} is the time step frequency for harmonic analysis
798
799- \texttt{nb\_ana} is the number of harmonics to analyse
800
801- \texttt{tname} is an array with names of tidal constituents to analyse
802
803\texttt{nit000\_han} and \texttt{nitend\_han} must be between \texttt{nit000} and \texttt{nitend} of the simulation.
804The restart capability is not implemented.
805
806The Harmonic analysis solve this equation:
807\begin{equation}
808h_{i} - A_{0} + \sum^{nb\_ana}_{j=1}[A_{j}cos(\nu_{j}t_{j}-\phi_{j})] = e_{i}
809\end{equation}
810
811With $A_{j}$,$\nu_{j}$,$\phi_{j}$, the amplitude, frequency and phase for each wave and $e_{i}$ the error.
812$h_{i}$ is the sea level for the time $t_{i}$ and $A_{0}$ is the mean sea level. \\
813We can rewrite this equation:
814\begin{equation}
815h_{i} - A_{0} + \sum^{nb\_ana}_{j=1}[C_{j}cos(\nu_{j}t_{j})+S_{j}sin(\nu_{j}t_{j})] = e_{i}
816\end{equation}
817with $A_{j}=\sqrt{C^{2}_{j}+S^{2}_{j}}$ et $\phi_{j}=arctan(S_{j}/C_{j})$.
818
819We obtain in output $C_{j}$ and $S_{j}$ for each tidal wave.
820
821% -------------------------------------------------------------------------------------------------------------
822%       Sections transports
823% -------------------------------------------------------------------------------------------------------------
824\section{Transports across sections (\key{diadct}) }
825\label{DIA_diag_dct}
826
827A module is available to compute the transport of volume, heat and salt through sections. This diagnostic
828is actived with \key{diadct}.
829
830Each section is defined by the coordinates of its 2 extremities. The pathways between them are contructed
831using tools which can be found in  \texttt{NEMOGCM/TOOLS/SECTIONS\_DIADCT} and are written in a binary file
832 \texttt{section\_ijglobal.diadct\_ORCA2\_LIM} which is later read in by NEMO to compute on-line transports.
833
834The on-line transports module creates three output ascii files:
835
836- \texttt{volume\_transport} for volume transports (  unit: $10^{6} m^{3} s^{-1}$ )
837
838- \texttt{heat\_transport}   for heat transports   (  unit: $10^{15} W $ )
839
840- \texttt{salt\_transport}   for salt transports   (  unit: $10^{9}Kg s^{-1}$ )\\
841
842
843Namelist parameters control how frequently the flows are summed and the time scales over which
844 they are averaged, as well as the level of output for debugging:
845
846%------------------------------------------namdct----------------------------------------------------
847\namdisplay{namdct}
848%-------------------------------------------------------------------------------------------------------------
849
850\texttt{nn\_dct}: frequency of instantaneous transports computing
851
852\texttt{nn\_dctwri}: frequency of writing ( mean of instantaneous transports )
853
854\texttt{nn\_debug}: debugging of the section
855
856\subsubsection{ To create a binary file containing the pathway of each section }
857
858In \texttt{NEMOGCM/TOOLS/SECTIONS\_DIADCT/run}, the file \texttt{ {list\_sections.ascii\_global}}
859contains a list of all the sections that are to be computed (this list of sections is based on MERSEA project metrics).
860
861Another file is available for the GYRE configuration (\texttt{ {list\_sections.ascii\_GYRE}}).
862
863Each section is defined by:
864
865\noindent \texttt{ long1 lat1 long2 lat2 nclass (ok/no)strpond (no)ice section\_name }\\
866with:
867
868- \texttt{long1 lat1} , coordinates of the first extremity of the section;
869
870- \texttt{long2 lat2} , coordinates of the second extremity of the section;
871
872- \texttt{nclass} the number of bounds of your classes (e.g. 3 bounds for 2 classes);
873
874- \texttt{okstrpond} to compute heat and salt transport, \texttt{nostrpond} if no;
875
876- \texttt{ice}  to compute surface and volume ice transports, \texttt{noice} if no. \\
877
878
879\noindent The results of the computing of transports, and the directions of positive
880 and negative flow do not depend on the order of the 2 extremities in this file.\\ 
881
882
883\noindent If nclass =/ 0,the next lines contain the class type and the nclass bounds:
884
885\texttt{long1 lat1 long2 lat2 nclass (ok/no)strpond (no)ice section\_name}
886
887\texttt{classtype}
888
889\texttt{zbound1}
890
891\texttt{zbound2}
892
893\texttt{.}
894
895\texttt{.}
896
897\texttt{nclass-1}
898
899\texttt{nclass}
900
901\noindent where \texttt{classtype} can be:
902
903- \texttt{zsal} for salinity classes
904
905- \texttt{ztem} for temperature classes
906
907- \texttt{zlay} for depth classes
908
909- \texttt{zsigi} for insitu density classes
910
911- \texttt{zsigp} for potential density classes \\
912
913 
914The script \texttt{job.ksh} computes the pathway for each section and creates a binary file
915\texttt{section\_ijglobal.diadct\_ORCA2\_LIM} which is read by NEMO. \\
916
917It is possible to use this tools for new configuations: \texttt{job.ksh} has to be updated
918with the coordinates file name and path. \\
919
920
921Examples of two sections, the ACC\_Drake\_Passage with no classes, and the
922 ATL\_Cuba\_Florida with 4 temperature clases (5 class bounds), are shown:
923
924\noindent \texttt{ -68.    -54.5   -60.    -64.7  00 okstrpond noice ACC\_Drake\_Passage}
925
926\noindent \texttt{ -80.5    22.5   -80.5    25.5  05 nostrpond noice ATL\_Cuba\_Florida}
927
928\noindent \texttt{ztem}
929
930\noindent \texttt{-2.0}
931
932\noindent \texttt{ 4.5}
933
934\noindent \texttt{ 7.0}
935
936\noindent \texttt{12.0}
937
938\noindent \texttt{40.0}
939
940
941\subsubsection{ To read the output files }
942
943The output format is :
944 
945{\small\texttt{date, time-step number, section number, section name, section slope coefficient, class number,
946class name, class bound 1 , classe bound2, transport\_direction1 ,  transport\_direction2, transport\_total}}\\
947
948
949For sections with classes, the first \texttt{nclass-1} lines correspond to the transport for each class
950and the last line corresponds to the total transport summed over all classes. For sections with no classes, class number
951\texttt{ 1 } corresponds to \texttt{ total class } and this class is called  \texttt{N}, meaning \texttt{none}.\\
952
953
954\noindent \texttt{ transport\_direction1 } is the positive part of the transport ( \texttt{ > = 0 } ).
955
956\noindent \texttt{ transport\_direction2 } is the negative part of the transport ( \texttt{ < = 0 } ).\\
957
958
959\noindent The \texttt{section slope coefficient} gives information about the significance of transports signs and direction:\\
960
961
962
963\begin{tabular}{|c|c|c|c|p{4cm}|}
964\hline
965section slope coefficient & section type & direction 1 & direction 2 & total transport \\ \hline
9660.    &  horizontal & northward & southward & postive: northward  \\ \hline
9671000. &  vertical   & eastward  & westward  & postive: eastward  \\ \hline               
968\texttt{=/0, =/ 1000.}   &  diagonal   & eastward  & westward  & postive: eastward  \\ \hline               
969\end{tabular}
970
971
972
973% -------------------------------------------------------------------------------------------------------------
974%       Other Diagnostics
975% -------------------------------------------------------------------------------------------------------------
976\section{Other Diagnostics (\key{diahth}, \key{diaar5})}
977\label{DIA_diag_others}
978
979
980Aside from the standard model variables, other diagnostics can be computed
981on-line. The available ready-to-add diagnostics routines can be found in directory DIA.
982Among the available diagnostics the following ones are obtained when defining
983the \key{diahth} CPP key:
984
985- the mixed layer depth (based on a density criterion, \citet{de_Boyer_Montegut_al_JGR04}) (\mdl{diahth})
986
987- the turbocline depth (based on a turbulent mixing coefficient criterion) (\mdl{diahth})
988
989- the depth of the 20\deg C isotherm (\mdl{diahth})
990
991- the depth of the thermocline (maximum of the vertical temperature gradient) (\mdl{diahth})
992
993The poleward heat and salt transports, their advective and diffusive component, and
994the meriodional stream function can be computed on-line in \mdl{diaptr} by setting
995\np{ln\_diaptr} to true (see the \textit{namptr} namelist below). 
996When \np{ln\_subbas}~=~true, transports and stream function are computed
997for the Atlantic, Indian, Pacific and Indo-Pacific Oceans (defined north of 30\deg S)
998as well as for the World Ocean. The sub-basin decomposition requires an input file
999(\ifile{subbasins}) which contains three 2D mask arrays, the Indo-Pacific mask
1000been deduced from the sum of the Indian and Pacific mask (Fig~\ref{Fig_mask_subasins}).
1001
1002%------------------------------------------namptr----------------------------------------------------
1003\namdisplay{namptr} 
1004%-------------------------------------------------------------------------------------------------------------
1005%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
1006\begin{figure}[!t]     \begin{center}
1007\includegraphics[width=1.0\textwidth]{./TexFiles/Figures/Fig_mask_subasins.pdf}
1008\caption{   \label{Fig_mask_subasins}
1009Decomposition of the World Ocean (here ORCA2) into sub-basin used in to compute
1010the heat and salt transports as well as the meridional stream-function: Atlantic basin (red),
1011Pacific basin (green), Indian basin (bleue), Indo-Pacific basin (bleue+green).
1012Note that semi-enclosed seas (Red, Med and Baltic seas) as well as Hudson Bay
1013are removed from the sub-basins. Note also that the Arctic Ocean has been split
1014into Atlantic and Pacific basins along the North fold line.  }
1015\end{center}   \end{figure} 
1016%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
1017
1018In addition, a series of diagnostics has been added in the \mdl{diaar5}.
1019They corresponds to outputs that are required for AR5 simulations
1020(see Section \ref{MISC_steric} below for one of them).
1021Activating those outputs requires to define the \key{diaar5} CPP key.
1022\\
1023\\
1024
1025
1026
1027% ================================================================
1028% Steric effect in sea surface height
1029% ================================================================
1030\section{Diagnosing the Steric effect in sea surface height}
1031\label{DIA_steric}
1032
1033
1034Changes in steric sea level are caused when changes in the density of the water
1035column imply an expansion or contraction of the column. It is essentially produced
1036through surface heating/cooling and to a lesser extent through non-linear effects of
1037the equation of state (cabbeling, thermobaricity...).
1038Non-Boussinesq models contain all ocean effects within the ocean acting
1039on the sea level. In particular, they include the steric effect. In contrast,
1040Boussinesq models, such as \NEMO, conserve volume, rather than mass,
1041and so do not properly represent expansion or contraction. The steric effect is
1042therefore not explicitely represented.
1043This approximation does not represent a serious error with respect to the flow field
1044calculated by the model \citep{Greatbatch_JGR94}, but extra attention is required
1045when investigating sea level, as steric changes are an important
1046contribution to local changes in sea level on seasonal and climatic time scales.
1047This is especially true for investigation into sea level rise due to global warming.
1048
1049Fortunately, the steric contribution to the sea level consists of a spatially uniform
1050component that can be diagnosed by considering the mass budget of the world
1051ocean \citep{Greatbatch_JGR94}.
1052In order to better understand how global mean sea level evolves and thus how
1053the steric sea level can be diagnosed, we compare, in the following, the
1054non-Boussinesq and Boussinesq cases.
1055
1056Let denote
1057$\mathcal{M}$ the total mass of liquid seawater ($\mathcal{M}=\int_D \rho dv$),
1058$\mathcal{V}$ the total volume of seawater ($\mathcal{V}=\int_D dv$),
1059$\mathcal{A}$ the total surface of the ocean ($\mathcal{A}=\int_S ds$),
1060$\bar{\rho}$ the global mean seawater (\textit{in situ}) density ($\bar{\rho}= 1/\mathcal{V} \int_D \rho \,dv$), and
1061$\bar{\eta}$ the global mean sea level ($\bar{\eta}=1/\mathcal{A}\int_S \eta \,ds$).
1062
1063A non-Boussinesq fluid conserves mass. It satisfies the following relations:
1064\begin{equation} \label{Eq_MV_nBq} 
1065\begin{split} 
1066\mathcal{M} &\mathcal{V}  \;\bar{\rho}       \\
1067\mathcal{V} &\mathcal{A}  \;\bar{\eta} 
1068\end{split}
1069\end{equation}
1070Temporal changes in total mass is obtained from the density conservation equation :
1071\begin{equation}  \label{Eq_Co_nBq}
1072\frac{1}{e_3} \partial_t ( e_3\,\rho) + \nabla( \rho \, \textbf{U} ) = \left. \frac{\textit{emp}}{e_3}\right|_\textit{surface}
1073\end{equation}
1074where $\rho$ is the \textit{in situ} density, and \textit{emp} the surface mass
1075exchanges with the other media of the Earth system (atmosphere, sea-ice, land).
1076Its global averaged leads to the total mass change
1077\begin{equation}  \label{Eq_Mass_nBq}
1078\partial_t \mathcal{M} = \mathcal{A} \;\overline{\textit{emp}}
1079\end{equation}
1080where $\overline{\textit{emp}}=\int_S \textit{emp}\,ds$ is the net mass flux
1081through the ocean surface.
1082Bringing \eqref{Eq_Mass_nBq} and the time derivative of \eqref{Eq_MV_nBq} 
1083together leads to the evolution equation of the mean sea level
1084\begin{equation} \label{Eq_ssh_nBq}
1085  \partial_t \bar{\eta} =  \frac{\overline{\textit{emp}}}{ \bar{\rho}} 
1086               - \frac{\mathcal{V}}{\mathcal{A}}  \;\frac{\partial_t \bar{\rho} }{\bar{\rho}}
1087\end{equation}
1088The first term in equation \eqref{Eq_ssh_nBq} alters sea level by adding or
1089subtracting mass from the ocean.
1090The second term arises from temporal changes in the global mean
1091density; $i.e.$ from steric effects.
1092
1093In a Boussinesq fluid, $\rho$ is replaced by $\rho_o$ in all the equation except when $\rho$ 
1094appears multiplied by the gravity ($i.e.$ in the hydrostatic balance of the primitive Equations).
1095In particular, the mass conservation equation, \eqref{Eq_Co_nBq}, degenerates into
1096the incompressibility equation:
1097\begin{equation}  \label{Eq_Co_Bq}
1098\frac{1}{e_3} \partial_t ( e_3 ) + \nabla( \textbf{U} ) =  \left. \frac{\textit{emp}}{\rho_o \,e_3}\right|_ \textit{surface}
1099\end{equation}
1100and the global average of this equation now gives the temporal change of the total volume,
1101\begin{equation}  \label{Eq_V_Bq}
1102  \partial_t \mathcal{V} =   \mathcal{A} \;\frac{\overline{\textit{emp}}}{\rho_o} 
1103\end{equation}
1104Only the volume is conserved, not mass, or, more precisely, the mass which is conserved is the
1105Boussinesq mass, $\mathcal{M}_o = \rho_o \mathcal{V}$. The total volume (or equivalently 
1106the global mean sea level) is altered only by net volume fluxes across the ocean surface, 
1107not by changes in mean mass of the ocean: the steric effect is missing in a Boussinesq fluid.
1108 
1109Nevertheless, following \citep{Greatbatch_JGR94}, the steric effect on the volume can be
1110diagnosed by considering the mass budget of the ocean.
1111The apparent changes in $\mathcal{M}$, mass of the ocean, which are not induced by surface
1112mass flux must be compensated by a spatially uniform change in the mean sea level due to
1113expansion/contraction of the ocean \citep{Greatbatch_JGR94}. In others words, the Boussinesq
1114mass, $\mathcal{M}_o$, can be related to $\mathcal{M}$, the  total mass of the ocean seen
1115by the Boussinesq model, via the steric contribution to the sea level, $\eta_s$, a spatially
1116uniform variable, as follows:
1117\begin{equation}  \label{Eq_M_Bq}
1118   \mathcal{M}_o  =  \mathcal{M} + \rho_o \,\eta_s \,\mathcal{A} 
1119\end{equation}
1120Any change in $\mathcal{M}$ which cannot be explained by the net mass flux through
1121the ocean surface is converted into a mean change in sea level. Introducing the total density
1122anomaly, $\mathcal{D}= \int_D d_a \,dv$, where $d_a= (\rho -\rho_o ) / \rho_o$ 
1123is the density anomaly used in \NEMO (cf. \S\ref{TRA_eos}) in \eqref{Eq_M_Bq}
1124leads to a very simple form for the steric height:
1125\begin{equation}  \label{Eq_steric_Bq}
1126   \eta_s = - \frac{1}{\mathcal{A}} \mathcal{D} 
1127\end{equation}
1128
1129The above formulation of the steric height of a Boussinesq ocean requires four remarks.
1130First, one can be tempted to define $\rho_o$ as the initial value of $\mathcal{M}/\mathcal{V}$,
1131$i.e.$ set $\mathcal{D}_{t=0}=0$, so that the initial steric height is zero. We do not
1132recommend that. Indeed, in this case $\rho_o$ depends on the initial state of the ocean.
1133Since $\rho_o$ has a direct effect on the dynamics of the ocean (it appears in the pressure
1134gradient term of the momentum equation) it is definitively not a good idea when
1135inter-comparing experiments.
1136We better recommend to fixe once for all $\rho_o$ to $1035\;Kg\,m^{-3}$. This value is a
1137sensible choice for the reference density used in a Boussinesq ocean climate model since,
1138with the exception of only a small percentage of the ocean, density in the World Ocean
1139varies by no more than 2$\%$ from this value (\cite{Gill1982}, page 47).
1140
1141Second, we have assumed here that the total ocean surface, $\mathcal{A}$, does not
1142change when the sea level is changing as it is the case in all global ocean GCMs
1143(wetting and drying of grid point is not allowed).
1144 
1145Third, the discretisation of \eqref{Eq_steric_Bq} depends on the type of free surface
1146which is considered. In the non linear free surface case, $i.e.$ \key{vvl} defined, it is
1147given by
1148\begin{equation}  \label{Eq_discrete_steric_Bq}
1149   \eta_s =  - \frac{ \sum_{i,\,j,\,k} d_a\; e_{1t} e_{2t} e_{3t} }
1150                  { \sum_{i,\,j,\,k}         e_{1t} e_{2t} e_{3t} } 
1151\end{equation}
1152whereas in the linear free surface, the volume above the \textit{z=0} surface must be explicitly taken
1153into account to better approximate the total ocean mass and thus the steric sea level:
1154\begin{equation}  \label{Eq_discrete_steric_Bq}
1155   \eta_s =  - \frac{ \sum_{i,\,j,\,k} d_a\; e_{1t}e_{2t}e_{3t} + \sum_{i,\,j} d_a\; e_{1t}e_{2t} \eta }
1156                     {\sum_{i,\,j,\,k} e_{1t}e_{2t}e_{3t} + \sum_{i,\,j}           e_{1t}e_{2t} \eta } 
1157\end{equation}
1158
1159The fourth and last remark concerns the effective sea level and the presence of sea-ice.
1160In the real ocean, sea ice (and snow above it)  depresses the liquid seawater through
1161its mass loading. This depression is a result of the mass of sea ice/snow system acting
1162on the liquid ocean. There is, however, no dynamical effect associated with these depressions
1163in the liquid ocean sea level, so that there are no associated ocean currents. Hence, the
1164dynamically relevant sea level is the effective sea level, $i.e.$ the sea level as if sea ice
1165(and snow) were converted to liquid seawater \citep{Campin_al_OM08}. However,
1166in the current version of \NEMO the sea-ice is levitating above the ocean without
1167mass exchanges between ice and ocean. Therefore the model effective sea level
1168is always given by $\eta + \eta_s$, whether or not there is sea ice present.
1169
1170In AR5 outputs, the thermosteric sea level is demanded. It is steric sea level due to
1171changes in ocean density arising just from changes in temperature. It is given by:
1172\begin{equation}  \label{Eq_thermosteric_Bq}
1173   \eta_s = - \frac{1}{\mathcal{A}} \int_D d_a(T,S_o,p_o) \,dv
1174\end{equation}
1175where $S_o$ and $p_o$ are the initial salinity and pressure, respectively.
1176
1177Both steric and thermosteric sea level are computed in \mdl{diaar5} which needs
1178the \key{diaar5} defined to be called.
1179
1180% ================================================================
1181
1182
1183
1184
1185
1186
1187
1188
1189
1190
1191
1192
Note: See TracBrowser for help on using the repository browser.