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.
Changeset 11954 for NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/doc/latex/NEMO/subfiles/chap_time_domain.tex – NEMO

Ignore:
Timestamp:
2019-11-22T17:15:18+01:00 (4 years ago)
Author:
acc
Message:

Branch 2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles. Merge in trunk changes up to 11943 in preparation for end of year merge

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/doc/latex/NEMO/subfiles/chap_time_domain.tex

    r11599 r11954  
    1313 
    1414{\footnotesize 
    15   \begin{tabularx}{\textwidth}{l||X|X} 
    16     Release & Author(s) & Modifications \\ 
     15  \begin{tabularx}{0.5\textwidth}{l||X|X} 
     16    Release          & Author(s)                                       & 
     17    Modifications                                                      \\ 
    1718    \hline 
    18     {\em   4.0} & {\em ...} & {\em ...} \\ 
    19     {\em   3.6} & {\em ...} & {\em ...} \\ 
    20     {\em   3.4} & {\em ...} & {\em ...} \\ 
    21     {\em <=3.4} & {\em ...} & {\em ...} 
     19    {\em        4.0} & {\em J\'{e}r\^{o}me Chanut \newline Tim Graham} & 
     20    {\em Review \newline Update                                      } \\ 
     21    {\em        3.6} & {\em Christian \'{E}th\'{e}                   } & 
     22    {\em Update                                                      } \\ 
     23    {\em $\leq$ 3.4} & {\em Gurvan Madec                             } & 
     24    {\em First version                                               } \\ 
    2225  \end{tabularx} 
    2326} 
     
    2629 
    2730% Missing things: 
    28 %  - daymod: definition of the time domain (nit000, nitend and the calendar) 
    29  
    30 \gmcomment{STEVEN :maybe a picture of the directory structure in the introduction which could be referred to here, 
    31   would help  ==> to be added} 
    32  
    33 Having defined the continuous equations in \autoref{chap:MB}, we need now to choose a time discretization, 
     31% - daymod: definition of the time domain (nit000, nitend and the calendar) 
     32 
     33\cmtgm{STEVEN :maybe a picture of the directory structure in the introduction which 
     34could be referred to here, would help  ==> to be added} 
     35 
     36Having defined the continuous equations in \autoref{chap:MB}, 
     37we need now to choose a time discretization, 
    3438a key feature of an ocean model as it exerts a strong influence on the structure of the computer code 
    3539(\ie\ on its flowchart). 
    36 In the present chapter, we provide a general description of the \NEMO\  time stepping strategy and 
     40In the present chapter, we provide a general description of the \NEMO\ time stepping strategy and 
    3741the consequences for the order in which the equations are solved. 
    3842 
     
    4751\end{equation} 
    4852where $x$ stands for $u$, $v$, $T$ or $S$; 
    49 RHS is the Right-Hand-Side of the corresponding time evolution equation; 
     53RHS is the \textbf{R}ight-\textbf{H}and-\textbf{S}ide of the corresponding time evolution equation; 
    5054$\rdt$ is the time step; 
    5155and the superscripts indicate the time at which a quantity is evaluated. 
    52 Each term of the RHS is evaluated at a specific time stepping depending on the physics with which it is associated. 
     56Each term of the RHS is evaluated at a specific time stepping depending on 
     57the physics with which it is associated. 
    5358 
    5459The choice of the time stepping used for this evaluation is discussed below as well as 
    5560the implications for starting or restarting a model simulation. 
    5661Note that the time stepping calculation is generally performed in a single operation. 
    57 With such a complex and nonlinear system of equations it would be dangerous to let a prognostic variable evolve in 
    58 time for each term separately. 
     62With such a complex and nonlinear system of equations it would be dangerous to 
     63let a prognostic variable evolve in time for each term separately. 
    5964 
    6065The three level scheme requires three arrays for each prognostic variable. 
     
    6267The third array, although referred to as $x_a$ (after) in the code, 
    6368is usually not the variable at the after time step; 
    64 but rather it is used to store the time derivative (RHS in \autoref{eq:TD}) prior to time-stepping the equation. 
    65 The time stepping itself is performed once at each time step where implicit vertical diffusion is computed, \ie\ in the \mdl{trazdf} and \mdl{dynzdf} modules. 
     69but rather it is used to store the time derivative (RHS in \autoref{eq:TD}) 
     70prior to time-stepping the equation. 
     71The time stepping itself is performed once at each time step where 
     72implicit vertical diffusion is computed, 
     73\ie\ in the \mdl{trazdf} and \mdl{dynzdf} modules. 
    6674 
    6775%% ================================================================================================= 
     
    6977\label{sec:TD_leap_frog} 
    7078 
    71 The time stepping used for processes other than diffusion is the well-known leapfrog scheme 
    72 \citep{mesinger.arakawa_bk76}. 
     79The time stepping used for processes other than diffusion is 
     80the well-known \textbf{L}eap\textbf{F}rog (LF) scheme \citep{mesinger.arakawa_bk76}. 
    7381This scheme is widely used for advection processes in low-viscosity fluids. 
    74 It is a time centred scheme, \ie\ the RHS in \autoref{eq:TD} is evaluated at time step $t$, the now time step. 
     82It is a time centred scheme, \ie\ the RHS in \autoref{eq:TD} is evaluated at 
     83time step $t$, the now time step. 
    7584It may be used for momentum and tracer advection, pressure gradient, and Coriolis terms, 
    7685but not for diffusion terms. 
    7786It is an efficient method that achieves second-order accuracy with 
    7887just one right hand side evaluation per time step. 
    79 Moreover, it does not artificially damp linear oscillatory motion nor does it produce instability by 
    80 amplifying the oscillations. 
     88Moreover, it does not artificially damp linear oscillatory motion 
     89nor does it produce instability by amplifying the oscillations. 
    8190These advantages are somewhat diminished by the large phase-speed error of the leapfrog scheme, 
    82 and the unsuitability of leapfrog differencing for the representation of diffusion and Rayleigh damping processes. 
     91and the unsuitability of leapfrog differencing for the representation of diffusion and 
     92Rayleigh damping processes. 
    8393However, the scheme allows the coexistence of a numerical and a physical mode due to 
    8494its leading third order dispersive error. 
    8595In other words a divergence of odd and even time steps may occur. 
    86 To prevent it, the leapfrog scheme is often used in association with a Robert-Asselin time filter 
    87 (hereafter the LF-RA scheme). 
    88 This filter, first designed by \citet{robert_JMSJ66} and more comprehensively studied by \citet{asselin_MWR72}, 
     96To prevent it, the leapfrog scheme is often used in association with 
     97a \textbf{R}obert-\textbf{A}sselin time filter (hereafter the LF-RA scheme). 
     98This filter, 
     99first designed by \citet{robert_JMSJ66} and more comprehensively studied by \citet{asselin_MWR72}, 
    89100is a kind of laplacian diffusion in time that mixes odd and even time steps: 
    90101\begin{equation} 
     
    99110However, the second order truncation error is proportional to $\gamma$, which is small compared to 1. 
    100111Therefore, the LF-RA is a quasi second order accurate scheme. 
    101 The LF-RA scheme is preferred to other time differencing schemes such as predictor corrector or trapezoidal schemes, 
    102 because the user has an explicit and simple control of the magnitude of the time diffusion of the scheme. 
    103 When used with the 2nd order space centred discretisation of the advection terms in 
     112The LF-RA scheme is preferred to other time differencing schemes such as 
     113predictor corrector or trapezoidal schemes, because the user has an explicit and simple control of 
     114the magnitude of the time diffusion of the scheme. 
     115When used with the 2$^nd$ order space centred discretisation of the advection terms in 
    104116the momentum and tracer equations, LF-RA avoids implicit numerical diffusion: 
    105 diffusion is set explicitly by the user through the Robert-Asselin 
    106 filter parameter and the viscosity and diffusion coefficients. 
     117diffusion is set explicitly by the user through the Robert-Asselin filter parameter and 
     118the viscosity and diffusion coefficients. 
    107119 
    108120%% ================================================================================================= 
     
    110122\label{sec:TD_forward_imp} 
    111123 
    112 The leapfrog differencing scheme is unsuitable for the representation of diffusion and damping processes. 
     124The leapfrog differencing scheme is unsuitable for 
     125the representation of diffusion and damping processes. 
    113126For a tendency $D_x$, representing a diffusion term or a restoring term to a tracer climatology 
    114127(when present, see \autoref{sec:TRA_dmp}), a forward time differencing scheme is used : 
     
    119132 
    120133This is diffusive in time and conditionally stable. 
    121 The conditions for stability of second and fourth order horizontal diffusion schemes are \citep{griffies_bk04}: 
     134The conditions for stability of second and fourth order horizontal diffusion schemes are 
     135\citep{griffies_bk04}: 
    122136\begin{equation} 
    123137  \label{eq:TD_euler_stability} 
     
    128142  \end{cases} 
    129143\end{equation} 
    130 where $e$ is the smallest grid size in the two horizontal directions and $A^h$ is the mixing coefficient. 
     144where $e$ is the smallest grid size in the two horizontal directions and 
     145$A^h$ is the mixing coefficient. 
    131146The linear constraint \autoref{eq:TD_euler_stability} is a necessary condition, but not sufficient. 
    132147If it is not satisfied, even mildly, then the model soon becomes wildly unstable. 
    133 The instability can be removed by either reducing the length of the time steps or reducing the mixing coefficient. 
     148The instability can be removed by either reducing the length of the time steps or 
     149reducing the mixing coefficient. 
    134150 
    135151For the vertical diffusion terms, a forward time differencing scheme can be used, 
    136 but usually the numerical stability condition imposes a strong constraint on the time step. To overcome the stability constraint, a 
    137 backward (or implicit) time differencing scheme is used. This scheme is unconditionally stable but diffusive and can be written as follows: 
     152but usually the numerical stability condition imposes a strong constraint on the time step. 
     153To overcome the stability constraint, a backward (or implicit) time differencing scheme is used. 
     154This scheme is unconditionally stable but diffusive and can be written as follows: 
    138155\begin{equation} 
    139156  \label{eq:TD_imp} 
     
    141158\end{equation} 
    142159 
    143 %%gm 
    144 %%gm   UPDATE the next paragraphs with time varying thickness ... 
    145 %%gm 
    146  
    147 This scheme is rather time consuming since it requires a matrix inversion. For example, the finite difference approximation of the temperature equation is: 
     160\cmtgm{UPDATE the next paragraphs with time varying thickness ...} 
     161 
     162This scheme is rather time consuming since it requires a matrix inversion. 
     163For example, the finite difference approximation of the temperature equation is: 
    148164\[ 
    149165  % \label{eq:TD_imp_zdf} 
     
    159175\end{equation} 
    160176where 
    161 \begin{align*} 
    162   c(k) &= A_w^{vT} (k) \, / \, e_{3w} (k)     \\ 
    163   d(k) &= e_{3t}   (k)       \, / \, (2 \rdt) + c_k + c_{k + 1}    \\ 
    164   b(k) &= e_{3t}   (k) \; \lt( T^{t - 1}(k) \, / \, (2 \rdt) + \text{RHS} \rt) 
    165 \end{align*} 
    166  
    167 \autoref{eq:TD_imp_mat} is a linear system of equations with an associated matrix which is tridiagonal. 
    168 Moreover, 
    169 $c(k)$ and $d(k)$ are positive and the diagonal term is greater than the sum of the two extra-diagonal terms, 
     177\[ 
     178  c(k) = A_w^{vT} (k) \, / \, e_{3w} (k) \text{,} \quad 
     179  d(k) = e_{3t}   (k)       \, / \, (2 \rdt) + c_k + c_{k + 1} \quad \text{and} \quad 
     180  b(k) = e_{3t}   (k) \; \lt( T^{t - 1}(k) \, / \, (2 \rdt) + \text{RHS} \rt) 
     181\] 
     182 
     183\autoref{eq:TD_imp_mat} is a linear system of equations with 
     184an associated matrix which is tridiagonal. 
     185Moreover, $c(k)$ and $d(k)$ are positive and 
     186the diagonal term is greater than the sum of the two extra-diagonal terms, 
    170187therefore a special adaptation of the Gauss elimination procedure is used to find the solution 
    171188(see for example \citet{richtmyer.morton_bk67}). 
     
    175192\label{sec:TD_spg_ts} 
    176193 
    177 The leapfrog environment supports a centred in time computation of the surface pressure, \ie\ evaluated 
    178 at \textit{now} time step. This refers to as the explicit free surface case in the code (\np[=.true.]{ln_dynspg_exp}{ln\_dynspg\_exp}). 
    179 This choice however imposes a strong constraint on the time step which should be small enough to resolve the propagation 
    180 of external gravity waves. As a matter of fact, one rather use in a realistic setup, a split-explicit free surface 
    181 (\np[=.true.]{ln_dynspg_ts}{ln\_dynspg\_ts}) in which barotropic and baroclinic dynamical equations are solved separately with ad-hoc 
    182 time steps. The use of the time-splitting (in combination with non-linear free surface) imposes some constraints on the design of 
    183 the overall flowchart, in particular to ensure exact tracer conservation (see \autoref{fig:TD_TimeStep_flowchart}). 
    184  
    185 Compared to the former use of the filtered free surface in \NEMO\ v3.6 (\citet{roullet.madec_JGR00}), the use of a split-explicit free surface is advantageous 
    186 on massively parallel computers. Indeed, no global computations are anymore required by the elliptic solver which saves a substantial amount of communication 
    187 time. Fast barotropic motions (such as tides) are also simulated with a better accuracy. 
    188  
    189 %\gmcomment{ 
    190 \begin{figure}[!t] 
     194The leapfrog environment supports a centred in time computation of the surface pressure, 
     195\ie\ evaluated at \textit{now} time step. 
     196This refers to as the explicit free surface case in the code 
     197(\np[=.true.]{ln_dynspg_exp}{ln\_dynspg\_exp}). 
     198This choice however imposes a strong constraint on the time step which 
     199should be small enough to resolve the propagation of external gravity waves. 
     200As a matter of fact, one rather use in a realistic setup, 
     201a split-explicit free surface (\np[=.true.]{ln_dynspg_ts}{ln\_dynspg\_ts}) in which 
     202barotropic and baroclinic dynamical equations are solved separately with ad-hoc time steps. 
     203The use of the time-splitting (in combination with non-linear free surface) imposes 
     204some constraints on the design of the overall flowchart, 
     205in particular to ensure exact tracer conservation (see \autoref{fig:TD_TimeStep_flowchart}). 
     206 
     207Compared to the former use of the filtered free surface in \NEMO\ v3.6 (\citet{roullet.madec_JGR00}), 
     208the use of a split-explicit free surface is advantageous on massively parallel computers. 
     209Indeed, no global computations are anymore required by the elliptic solver which 
     210saves a substantial amount of communication time. 
     211Fast barotropic motions (such as tides) are also simulated with a better accuracy. 
     212 
     213%\cmtgm{ 
     214\begin{figure} 
    191215  \centering 
    192   \includegraphics[width=0.66\textwidth]{Fig_TimeStepping_flowchart_v4} 
     216  \includegraphics[width=0.66\textwidth]{TD_TimeStepping_flowchart_v4} 
    193217  \caption[Leapfrog time stepping sequence with split-explicit free surface]{ 
    194218    Sketch of the leapfrog time stepping sequence in \NEMO\ with split-explicit free surface. 
    195     The latter combined with non-linear free surface requires the dynamical tendency being 
    196     updated prior tracers tendency to ensure conservation. 
     219    The latter combined with non-linear free surface requires 
     220    the dynamical tendency being updated prior tracers tendency to ensure conservation. 
    197221    Note the use of time integrated fluxes issued from the barotropic loop in 
    198222    subsequent calculations of tracer advection and in the continuity equation. 
     
    203227 
    204228%% ================================================================================================= 
    205 \section{Modified Leapfrog -- Asselin filter scheme} 
     229\section{Modified LeapFrog -- Robert Asselin filter scheme (LF-RA)} 
    206230\label{sec:TD_mLF} 
    207231 
    208 Significant changes have been introduced by \cite{leclair.madec_OM09} in the LF-RA scheme in order to 
    209 ensure tracer conservation and to allow the use of a much smaller value of the Asselin filter parameter. 
     232Significant changes have been introduced by \cite{leclair.madec_OM09} in 
     233the LF-RA scheme in order to ensure tracer conservation and to 
     234allow the use of a much smaller value of the Asselin filter parameter. 
    210235The modifications affect both the forcing and filtering treatments in the LF-RA scheme. 
    211236 
    212 In a classical LF-RA environment, the forcing term is centred in time, 
    213 \ie\ it is time-stepped over a $2 \rdt$ period: 
     237In a classical LF-RA environment, 
     238the forcing term is centred in time, \ie\ it is time-stepped over a $2 \rdt$ period: 
    214239$x^t = x^t + 2 \rdt Q^t$ where $Q$ is the forcing applied to $x$, 
    215 and the time filter is given by \autoref{eq:TD_asselin} so that $Q$ is redistributed over several time step. 
     240and the time filter is given by \autoref{eq:TD_asselin} so that 
     241$Q$ is redistributed over several time step. 
    216242In the modified LF-RA environment, these two formulations have been replaced by: 
    217243\begin{gather} 
     
    222248                    - \gamma \, \rdt \, \lt( Q^{t + \rdt / 2} - Q^{t - \rdt / 2} \rt) 
    223249\end{gather} 
    224 The change in the forcing formulation given by \autoref{eq:TD_forcing} (see \autoref{fig:TD_MLF_forcing}) 
    225 has a significant effect: 
    226 the forcing term no longer excites the divergence of odd and even time steps \citep{leclair.madec_OM09}. 
     250The change in the forcing formulation given by \autoref{eq:TD_forcing} 
     251(see \autoref{fig:TD_MLF_forcing}) has a significant effect: 
     252the forcing term no longer excites the divergence of odd and even time steps 
     253\citep{leclair.madec_OM09}. 
    227254% forcing seen by the model.... 
    228255This property improves the LF-RA scheme in two aspects. 
    229256First, the LF-RA can now ensure the local and global conservation of tracers. 
    230257Indeed, time filtering is no longer required on the forcing part. 
    231 The influence of the Asselin filter on the forcing is explicitly removed by adding a new term in the filter 
    232 (last term in \autoref{eq:TD_RA} compared to \autoref{eq:TD_asselin}). 
     258The influence of the Asselin filter on the forcing is explicitly removed by 
     259adding a new term in the filter (last term in \autoref{eq:TD_RA} compared to \autoref{eq:TD_asselin}). 
    233260Since the filtering of the forcing was the source of non-conservation in the classical LF-RA scheme, 
    234261the modified formulation becomes conservative \citep{leclair.madec_OM09}. 
    235 Second, the LF-RA becomes a truly quasi -second order scheme. 
     262Second, the LF-RA becomes a truly quasi-second order scheme. 
    236263Indeed, \autoref{eq:TD_forcing} used in combination with a careful treatment of static instability 
    237264(\autoref{subsec:ZDF_evd}) and of the TKE physics (\autoref{subsec:ZDF_tke_ene}) 
     
    245272even if separated by only $\rdt$ since the time filter is no longer applied to the forcing term. 
    246273 
    247 \begin{figure}[!t] 
     274\begin{figure} 
    248275  \centering 
    249   \includegraphics[width=0.66\textwidth]{Fig_MLF_forcing} 
     276  \includegraphics[width=0.66\textwidth]{TD_MLF_forcing} 
    250277  \caption[Forcing integration methods for modified leapfrog (top and bottom)]{ 
    251278    Illustration of forcing integration methods. 
     
    271298\end{listing} 
    272299 
    273 The first time step of this three level scheme when starting from initial conditions is a forward step 
    274 (Euler time integration): 
     300The first time step of this three level scheme when starting from initial conditions is 
     301a forward step (Euler time integration): 
    275302\[ 
    276303  % \label{eq:TD_DOM_euler} 
    277304  x^1 = x^0 + \rdt \ \text{RHS}^0 
    278305\] 
    279 This is done simply by keeping the leapfrog environment (\ie\ the \autoref{eq:TD} three level time stepping) but 
     306This is done simply by keeping the leapfrog environment 
     307(\ie\ the \autoref{eq:TD} three level time stepping) but 
    280308setting all $x^0$ (\textit{before}) and $x^1$ (\textit{now}) fields equal at the first time step and 
    281309using half the value of a leapfrog time step ($2 \rdt$). 
     
    286314running the model for $2N$ time steps in one go, 
    287315or by performing two consecutive experiments of $N$ steps with a restart. 
    288 This requires saving two time levels and many auxiliary data in the restart files in machine precision. 
     316This requires saving two time levels and many auxiliary data in 
     317the restart files in machine precision. 
    289318 
    290319Note that the time step $\rdt$, is also saved in the restart file. 
    291 When restarting, if the time step has been changed, or one of the prognostic variables at \textit{before} time step 
    292 is missing, an Euler time stepping scheme is imposed. A forward initial step can still be enforced by the user by setting 
    293 the namelist variable \np[=0]{nn_euler}{nn\_euler}. Other options to control the time integration of the model 
    294 are defined through the  \nam{run}{run} namelist variables. 
    295 \gmcomment{ 
     320When restarting, if the time step has been changed, or 
     321one of the prognostic variables at \textit{before} time step is missing, 
     322an Euler time stepping scheme is imposed. 
     323A forward initial step can still be enforced by the user by 
     324setting the namelist variable \np[=0]{nn_euler}{nn\_euler}. 
     325Other options to control the time integration of the model are defined through 
     326the \nam{run}{run} namelist variables. 
     327 
     328\cmtgm{ 
    296329add here how to force the restart to contain only one time step for operational purposes 
    297330 
    298331add also the idea of writing several restart for seasonal forecast : how is it done ? 
    299332 
    300 verify that all namelist pararmeters are truly described 
     333verify that all namelist parameters are truly described 
    301334 
    302335a word on the check of restart  ..... 
    303336} 
    304337 
    305 \gmcomment{       % add a subsection here 
     338\cmtgm{       % add a subsection here 
    306339 
    307340%% ================================================================================================= 
     
    309342\label{subsec:TD_time} 
    310343 
    311 Options are defined through the  \nam{dom}{dom} namelist variables. 
     344Options are defined through the\nam{dom}{dom} namelist variables. 
    312345 \colorbox{yellow}{add here a few word on nit000 and nitend} 
    313346 
    314347 \colorbox{yellow}{Write documentation on the calendar and the key variable adatrj} 
    315348 
    316 add a description of daymod, and the model calandar (leap-year and co) 
    317  
    318 }        %% end add 
    319  
    320 \gmcomment{       % add implicit in vvl case  and Crant-Nicholson scheme 
     349add a description of daymod, and the model calendar (leap-year and co) 
     350 
     351}     %% end add 
     352 
     353\cmtgm{       % add implicit in vvl case  and Crant-Nicholson scheme 
    321354 
    322355Implicit time stepping in case of variable volume thickness. 
     
    369402} 
    370403 
    371 \onlyinsubfile{\input{../../global/epilogue}} 
     404\subinc{\input{../../global/epilogue}} 
    372405 
    373406\end{document} 
Note: See TracChangeset for help on using the changeset viewer.