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 11512 for NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/doc/latex/NEMO/subfiles/chap_TRA.tex – NEMO

Ignore:
Timestamp:
2019-09-09T12:05:20+02:00 (5 years ago)
Author:
smasson
Message:

dev_r10984_HPC-13 : merge with trunk@11511, see #2285

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/doc/latex/NEMO/subfiles/chap_TRA.tex

    r11263 r11512  
    88\label{chap:TRA} 
    99 
    10 \minitoc 
    11  
    12 % missing/update  
     10\chaptertoc 
     11 
     12% missing/update 
    1313% traqsr: need to coordinate with SBC module 
    1414 
    15 %STEVEN :  is the use of the word "positive" to describe a scheme enough, or should it be "positive definite"? I added a comment to this effect on some instances of this below 
     15%STEVEN :  is the use of the word "positive" to describe a scheme enough, or should it be "positive definite"? 
     16%I added a comment to this effect on some instances of this below 
    1617 
    1718Using the representation described in \autoref{chap:DOM}, several semi -discrete space forms of 
     
    3536The terms QSR, BBC, BBL and DMP are optional. 
    3637The external forcings and parameterisations require complex inputs and complex calculations 
    37 (\eg bulk formulae, estimation of mixing coefficients) that are carried out in the SBC, 
     38(\eg\ bulk formulae, estimation of mixing coefficients) that are carried out in the SBC, 
    3839LDF and ZDF modules and described in \autoref{chap:SBC}, \autoref{chap:LDF} and 
    3940\autoref{chap:ZDF}, respectively. 
     
    4748associated modules \mdl{eosbn2} and \mdl{phycst}). 
    4849 
    49 The different options available to the user are managed by namelist logicals or CPP keys. 
     50The different options available to the user are managed by namelist logicals. 
    5051For each equation term \textit{TTT}, the namelist logicals are \textit{ln\_traTTT\_xxx}, 
    5152where \textit{xxx} is a 3 or 4 letter acronym corresponding to each optional scheme. 
    52 The CPP key (when it exists) is \key{traTTT}. 
    5353The equivalent code can be found in the \textit{traTTT} or \textit{traTTT\_xxx} module, 
    5454in the \path{./src/OCE/TRA} directory. 
     
    6868%------------------------------------------------------------------------------------------------------------- 
    6969 
    70 When considered (\ie when \np{ln\_traadv\_NONE} is not set to \forcode{.true.}), 
     70When considered (\ie\ when \np{ln\_traadv\_OFF} is not set to \forcode{.true.}), 
    7171the advection tendency of a tracer is expressed in flux form, 
    72 \ie as the divergence of the advective fluxes. 
     72\ie\ as the divergence of the advective fluxes. 
    7373Its discrete expression is given by : 
    7474\begin{equation} 
     
    8282Indeed, it is obtained by using the following equality: $\nabla \cdot (\vect U \, T) = \vect U \cdot \nabla T$ which 
    8383results from the use of the continuity equation, $\partial_t e_3 + e_3 \; \nabla \cdot \vect U = 0$ 
    84 (which reduces to $\nabla \cdot \vect U = 0$ in linear free surface, \ie \np{ln\_linssh}\forcode{ = .true.}). 
     84(which reduces to $\nabla \cdot \vect U = 0$ in linear free surface, \ie\ \np{ln\_linssh}\forcode{ = .true.}). 
    8585Therefore it is of paramount importance to design the discrete analogue of the advection tendency so that 
    8686it is consistent with the continuity equation in order to enforce the conservation properties of 
     
    110110%>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    111111 
    112 The key difference between the advection schemes available in \NEMO is the choice made in space and 
     112The key difference between the advection schemes available in \NEMO\ is the choice made in space and 
    113113time interpolation to define the value of the tracer at the velocity points 
    114114(\autoref{fig:adv_scheme}). 
     
    125125  the moving surface $z = \eta$. 
    126126  There is a non-zero advective flux which is set for all advection schemes as 
    127   $\tau_w|_{k = 1/2} = T_{k = 1}$, \ie the product of surface velocity (at $z = 0$) by 
     127  $\tau_w|_{k = 1/2} = T_{k = 1}$, \ie\ the product of surface velocity (at $z = 0$) by 
    128128  the first level tracer value. 
    129129\item[non-linear free surface:] 
     
    139139two quantities that are not correlated \citep{roullet.madec_JGR00, griffies.pacanowski.ea_MWR01, campin.adcroft.ea_OM04}. 
    140140 
    141 The velocity field that appears in (\autoref{eq:tra_adv} and \autoref{eq:tra_adv_zco?}) is 
    142 the centred (\textit{now}) \textit{effective} ocean velocity, \ie the \textit{eulerian} velocity 
     141The velocity field that appears in (\autoref{eq:tra_adv} is  
     142the centred (\textit{now}) \textit{effective} ocean velocity, \ie\ the \textit{eulerian} velocity 
    143143(see \autoref{chap:DYN}) plus the eddy induced velocity (\textit{eiv}) and/or 
    144144the mixed layer eddy induced velocity (\textit{eiv}) when those parameterisations are used 
     
    149149Conservative Laws scheme (MUSCL), a $3^{rd}$ Upstream Biased Scheme (UBS, also often called UP3), 
    150150and a Quadratic Upstream Interpolation for Convective Kinematics with Estimated Streaming Terms scheme (QUICKEST). 
    151 The choice is made in the \ngn{namtra\_adv} namelist, by setting to \forcode{.true.} one of 
     151The choice is made in the \nam{tra\_adv} namelist, by setting to \forcode{.true.} one of 
    152152the logicals \textit{ln\_traadv\_xxx}. 
    153153The corresponding code can be found in the \textit{traadv\_xxx.F90} module, where 
    154154\textit{xxx} is a 3 or 4 letter acronym corresponding to each scheme. 
    155 By default (\ie in the reference namelist, \textit{namelist\_ref}), all the logicals are set to \forcode{.false.}. 
     155By default (\ie\ in the reference namelist, \textit{namelist\_ref}), all the logicals are set to \forcode{.false.}. 
    156156If the user does not select an advection scheme in the configuration namelist (\textit{namelist\_cfg}), 
    157157the tracers will \textit{not} be advected! 
     
    188188\label{subsec:TRA_adv_cen} 
    189189 
    190 %        2nd order centred scheme   
     190%        2nd order centred scheme 
    191191 
    192192The centred advection scheme (CEN) is used when \np{ln\_traadv\_cen}\forcode{ = .true.}. 
     
    203203\end{equation} 
    204204 
    205 CEN2 is non diffusive (\ie it conserves the tracer variance, $\tau^2$) but dispersive 
    206 (\ie it may create false extrema). 
     205CEN2 is non diffusive (\ie\ it conserves the tracer variance, $\tau^2$) but dispersive 
     206(\ie\ it may create false extrema). 
    207207It is therefore notoriously noisy and must be used in conjunction with an explicit diffusion operator to 
    208208produce a sensible solution. 
     
    213213both (\autoref{eq:tra_adv}) and (\autoref{eq:tra_adv_cen2}) have this order of accuracy. 
    214214 
    215 %        4nd order centred scheme   
     215%        4nd order centred scheme 
    216216 
    217217In the $4^{th}$ order formulation (CEN4), tracer values are evaluated at u- and v-points as 
     
    225225a $4^{th}$ COMPACT interpolation has been prefered \citep{demange_phd14}. 
    226226In the COMPACT scheme, both the field and its derivative are interpolated, which leads, after a matrix inversion, 
    227 spectral characteristics similar to schemes of higher order \citep{lele_JCP92}.  
     227spectral characteristics similar to schemes of higher order \citep{lele_JCP92}. 
    228228 
    229229Strictly speaking, the CEN4 scheme is not a $4^{th}$ order advection scheme but 
     
    237237 
    238238A direct consequence of the pseudo-fourth order nature of the scheme is that it is not non-diffusive, 
    239 \ie the global variance of a tracer is not preserved using CEN4. 
     239\ie\ the global variance of a tracer is not preserved using CEN4. 
    240240Furthermore, it must be used in conjunction with an explicit diffusion operator to produce a sensible solution. 
    241241As in CEN2 case, the time-stepping is performed using a leapfrog scheme in conjunction with an Asselin time-filter, 
     
    250250 
    251251% ------------------------------------------------------------------------------------------------------------- 
    252 %        FCT scheme   
     252%        FCT scheme 
    253253% ------------------------------------------------------------------------------------------------------------- 
    254254\subsection[FCT: Flux Corrected Transport scheme (\forcode{ln_traadv_fct = .true.})] 
     
    278278where $c_u$ is a flux limiter function taking values between 0 and 1. 
    279279The FCT order is the one of the centred scheme used 
    280 (\ie it depends on the setting of \np{nn\_fct\_h} and \np{nn\_fct\_v}). 
     280(\ie\ it depends on the setting of \np{nn\_fct\_h} and \np{nn\_fct\_v}). 
    281281There exist many ways to define $c_u$, each corresponding to a different FCT scheme. 
    282 The one chosen in \NEMO is described in \citet{zalesak_JCP79}. 
     282The one chosen in \NEMO\ is described in \citet{zalesak_JCP79}. 
    283283$c_u$ only departs from $1$ when the advective term produces a local extremum in the tracer field. 
    284284The resulting scheme is quite expensive but \textit{positive}. 
     
    286286A comparison of FCT-2 with MUSCL and a MPDATA scheme can be found in \citet{levy.estublier.ea_GRL01}. 
    287287 
    288 An additional option has been added controlled by \np{nn\_fct\_zts}. 
    289 By setting this integer to a value larger than zero, 
    290 a $2^{nd}$ order FCT scheme is used on both horizontal and vertical direction, but on the latter, 
    291 a split-explicit time stepping is used, with a number of sub-timestep equals to \np{nn\_fct\_zts}. 
    292 This option can be useful when the size of the timestep is limited by vertical advection \citep{lemarie.debreu.ea_OM15}. 
    293 Note that in this case, a similar split-explicit time stepping should be used on vertical advection of momentum to 
    294 insure a better stability (see \autoref{subsec:DYN_zad}). 
    295288 
    296289For stability reasons (see \autoref{chap:STP}), 
     
    301294 
    302295% ------------------------------------------------------------------------------------------------------------- 
    303 %        MUSCL scheme   
     296%        MUSCL scheme 
    304297% ------------------------------------------------------------------------------------------------------------- 
    305298\subsection[MUSCL: Monotone Upstream Scheme for Conservative Laws (\forcode{ln_traadv_mus = .true.})] 
     
    310303MUSCL implementation can be found in the \mdl{traadv\_mus} module. 
    311304 
    312 MUSCL has been first implemented in \NEMO by \citet{levy.estublier.ea_GRL01}. 
     305MUSCL has been first implemented in \NEMO\ by \citet{levy.estublier.ea_GRL01}. 
    313306In its formulation, the tracer at velocity points is evaluated assuming a linear tracer variation between 
    314307two $T$-points (\autoref{fig:adv_scheme}). 
     
    338331 
    339332% ------------------------------------------------------------------------------------------------------------- 
    340 %        UBS scheme   
     333%        UBS scheme 
    341334% ------------------------------------------------------------------------------------------------------------- 
    342335\subsection[UBS a.k.a. UP3: Upstream-Biased Scheme (\forcode{ln_traadv_ubs = .true.})] 
     
    374367\citep{shchepetkin.mcwilliams_OM05, demange_phd14}. 
    375368Therefore the vertical flux is evaluated using either a $2^nd$ order FCT scheme or a $4^th$ order COMPACT scheme 
    376 (\np{nn\_cen\_v}\forcode{ = 2 or 4}). 
     369(\np{nn\_ubs\_v}\forcode{ = 2 or 4}). 
    377370 
    378371For stability reasons (see \autoref{chap:STP}), the first term  in \autoref{eq:tra_adv_ubs} 
     
    408401Thirdly, the diffusion term is in fact a biharmonic operator with an eddy coefficient which 
    409402is simply proportional to the velocity: $A_u^{lm} = \frac{1}{12} \, {e_{1u}}^3 \, |u|$. 
    410 Note the current version of NEMO uses the computationally more efficient formulation \autoref{eq:tra_adv_ubs}. 
    411  
    412 % ------------------------------------------------------------------------------------------------------------- 
    413 %        QCK scheme   
     403Note the current version of \NEMO\ uses the computationally more efficient formulation \autoref{eq:tra_adv_ubs}. 
     404 
     405% ------------------------------------------------------------------------------------------------------------- 
     406%        QCK scheme 
    414407% ------------------------------------------------------------------------------------------------------------- 
    415408\subsection[QCK: QuiCKest scheme (\forcode{ln_traadv_qck = .true.})] 
     
    423416QUICKEST is the third order Godunov scheme which is associated with the ULTIMATE QUICKEST limiter 
    424417\citep{leonard_CMAME91}. 
    425 It has been implemented in NEMO by G. Reffray (MERCATOR-ocean) and can be found in the \mdl{traadv\_qck} module. 
     418It has been implemented in \NEMO\ by G. Reffray (MERCATOR-ocean) and can be found in the \mdl{traadv\_qck} module. 
    426419The resulting scheme is quite expensive but \textit{positive}. 
    427420It can be used on both active and passive tracers. 
     
    444437\nlst{namtra_ldf} 
    445438%------------------------------------------------------------------------------------------------------------- 
    446   
    447 Options are defined through the \ngn{namtra\_ldf} namelist variables. 
    448 They are regrouped in four items, allowing to specify  
     439 
     440Options are defined through the \nam{tra\_ldf} namelist variables. 
     441They are regrouped in four items, allowing to specify 
    449442$(i)$   the type of operator used (none, laplacian, bilaplacian), 
    450443$(ii)$  the direction along which the operator acts (iso-level, horizontal, iso-neutral), 
    451 $(iii)$ some specific options related to the rotated operators (\ie non-iso-level operator), and 
     444$(iii)$ some specific options related to the rotated operators (\ie\ non-iso-level operator), and 
    452445$(iv)$  the specification of eddy diffusivity coefficient (either constant or variable in space and time). 
    453446Item $(iv)$ will be described in \autoref{chap:LDF}. 
     
    457450 
    458451The lateral diffusion of tracers is evaluated using a forward scheme, 
    459 \ie the tracers appearing in its expression are the \textit{before} tracers in time, 
     452\ie\ the tracers appearing in its expression are the \textit{before} tracers in time, 
    460453except for the pure vertical component that appears when a rotation tensor is used. 
    461454This latter component is solved implicitly together with the vertical diffusion term (see \autoref{chap:STP}). 
     
    466459%        Type of operator 
    467460% ------------------------------------------------------------------------------------------------------------- 
    468 \subsection[Type of operator (\texttt{ln\_traldf}\{\texttt{\_NONE,\_lap,\_blp}\})] 
    469 {Type of operator (\protect\np{ln\_traldf\_NONE}, \protect\np{ln\_traldf\_lap}, or \protect\np{ln\_traldf\_blp}) }  
     461\subsection[Type of operator (\texttt{ln\_traldf}\{\texttt{\_OFF,\_lap,\_blp}\})] 
     462{Type of operator (\protect\np{ln\_traldf\_OFF}, \protect\np{ln\_traldf\_lap}, or \protect\np{ln\_traldf\_blp}) } 
    470463\label{subsec:TRA_ldf_op} 
    471464 
     
    473466 
    474467\begin{description} 
    475 \item[\np{ln\_traldf\_NONE}\forcode{ = .true.}:] 
     468\item[\np{ln\_traldf\_OFF}\forcode{ = .true.}:] 
    476469  no operator selected, the lateral diffusive tendency will not be applied to the tracer equation. 
    477470  This option can be used when the selected advection scheme is diffusive enough (MUSCL scheme for example). 
     
    494487minimizing the impact on the larger scale features. 
    495488The main difference between the two operators is the scale selectiveness. 
    496 The bilaplacian damping time (\ie its spin down time) scales like $\lambda^{-4}$ for 
     489The bilaplacian damping time (\ie\ its spin down time) scales like $\lambda^{-4}$ for 
    497490disturbances of wavelength $\lambda$ (so that short waves damped more rapidelly than long ones), 
    498491whereas the laplacian damping time scales only like $\lambda^{-2}$. 
     
    502495% ------------------------------------------------------------------------------------------------------------- 
    503496\subsection[Action direction (\texttt{ln\_traldf}\{\texttt{\_lev,\_hor,\_iso,\_triad}\})] 
    504 {Direction of action (\protect\np{ln\_traldf\_lev}, \protect\np{ln\_traldf\_hor}, \protect\np{ln\_traldf\_iso}, or \protect\np{ln\_traldf\_triad}) }  
     497{Direction of action (\protect\np{ln\_traldf\_lev}, \protect\np{ln\_traldf\_hor}, \protect\np{ln\_traldf\_iso}, or \protect\np{ln\_traldf\_triad}) } 
    505498\label{subsec:TRA_ldf_dir} 
    506499 
     
    508501The operator is a simple (re-entrant) laplacian acting in the (\textbf{i},\textbf{j}) plane when 
    509502iso-level option is used (\np{ln\_traldf\_lev}\forcode{ = .true.}) or 
    510 when a horizontal (\ie geopotential) operator is demanded in \textit{z}-coordinate 
     503when a horizontal (\ie\ geopotential) operator is demanded in \textit{z}-coordinate 
    511504(\np{ln\_traldf\_hor} and \np{ln\_zco} equal \forcode{.true.}). 
    512505The associated code can be found in the \mdl{traldf\_lap\_blp} module. 
     
    516509(\np{ln\_traldf\_iso} or \np{ln\_traldf\_triad} equals \forcode{.true.}, 
    517510see \mdl{traldf\_iso} or \mdl{traldf\_triad} module, resp.), or 
    518 when a horizontal (\ie geopotential) operator is demanded in \textit{s}-coordinate 
     511when a horizontal (\ie\ geopotential) operator is demanded in \textit{s}-coordinate 
    519512(\np{ln\_traldf\_hor} and \np{ln\_sco} equal \forcode{.true.}) 
    520513\footnote{In this case, the standard iso-neutral operator will be automatically selected}. 
     
    532525\label{subsec:TRA_ldf_lev} 
    533526 
    534 The laplacian diffusion operator acting along the model (\textit{i,j})-surfaces is given by:  
     527The laplacian diffusion operator acting along the model (\textit{i,j})-surfaces is given by: 
    535528\begin{equation} 
    536529  \label{eq:tra_ldf_lap} 
     
    541534where zero diffusive fluxes is assumed across solid boundaries, 
    542535first (and third in bilaplacian case) horizontal tracer derivative are masked. 
    543 It is implemented in the \rou{traldf\_lap} subroutine found in the \mdl{traldf\_lap} module. 
    544 The module also contains \rou{traldf\_blp}, the subroutine calling twice \rou{traldf\_lap} in order to 
     536It is implemented in the \rou{tra\_ldf\_lap} subroutine found in the \mdl{traldf\_lap\_blp}} module. 
     537The module also contains \rou{tra\_ldf\_blp}, the subroutine calling twice \rou{tra\_ldf\_lap} in order to 
    545538compute the iso-level bilaplacian operator. 
    546539 
     
    584577where $b_t = e_{1t} \, e_{2t} \, e_{3t}$  is the volume of $T$-cells, 
    585578$r_1$ and $r_2$ are the slopes between the surface of computation ($z$- or $s$-surfaces) and 
    586 the surface along which the diffusion operator acts (\ie horizontal or iso-neutral surfaces). 
     579the surface along which the diffusion operator acts (\ie\ horizontal or iso-neutral surfaces). 
    587580It is thus used when, in addition to \np{ln\_traldf\_lap}\forcode{ = .true.}, 
    588581we have \np{ln\_traldf\_iso}\forcode{ = .true.}, 
     
    613606\label{subsec:TRA_ldf_triad} 
    614607 
    615 If the Griffies triad scheme is employed (\np{ln\_traldf\_triad}\forcode{ = .true.}; see \autoref{apdx:triad}) 
    616  
    617608An alternative scheme developed by \cite{griffies.gnanadesikan.ea_JPO98} which ensures tracer variance decreases 
    618 is also available in \NEMO (\np{ln\_traldf\_grif}\forcode{ = .true.}). 
     609is also available in \NEMO\ (\np{ln\_traldf\_triad}\forcode{ = .true.}). 
    619610A complete description of the algorithm is given in \autoref{apdx:triad}. 
    620611 
     
    637628\item \np{rn\_slpmax} = slope limit (both operators) 
    638629\item \np{ln\_triad\_iso} = pure horizontal mixing in ML (triad only) 
    639 \item \np{rn\_sw\_triad} $= 1$ switching triad; $= 0$ all 4 triads used (triad only)  
     630\item \np{rn\_sw\_triad} $= 1$ switching triad; $= 0$ all 4 triads used (triad only) 
    640631\item \np{ln\_botmix\_triad} = lateral mixing on bottom (triad only) 
    641632\end{itemize} 
     
    652643%-------------------------------------------------------------------------------------------------------------- 
    653644 
    654 Options are defined through the \ngn{namzdf} namelist variables. 
     645Options are defined through the \nam{zdf} namelist variables. 
    655646The formulation of the vertical subgrid scale tracer physics is the same for all the vertical coordinates, 
    656647and is based on a laplacian operator. 
     
    664655respectively. 
    665656Generally, $A_w^{vT} = A_w^{vS}$ except when double diffusive mixing is parameterised 
    666 (\ie \key{zdfddm} is defined). 
     657(\ie\ \np{ln\_zdfddm} equals \forcode{.true.},). 
    667658The way these coefficients are evaluated is given in \autoref{chap:ZDF} (ZDF). 
    668659Furthermore, when iso-neutral mixing is used, both mixing coefficients are increased by 
     
    676667 
    677668The large eddy coefficient found in the mixed layer together with high vertical resolution implies that 
    678 in the case of explicit time stepping (\np{ln\_zdfexp}\forcode{ = .true.}) 
    679 there would be too restrictive a constraint on the time step. 
    680 Therefore, the default implicit time stepping is preferred for the vertical diffusion since 
     669there would be too restrictive constraint on the time step if we use explicit time stepping. 
     670Therefore an implicit time stepping is preferred for the vertical diffusion since 
    681671it overcomes the stability constraint. 
    682 A forward time differencing scheme (\np{ln\_zdfexp}\forcode{ = .true.}) using 
    683 a time splitting technique (\np{nn\_zdfexp} $> 1$) is provided as an alternative. 
    684 Namelist variables \np{ln\_zdfexp} and \np{nn\_zdfexp} apply to both tracers and dynamics. 
    685672 
    686673% ================================================================ 
     
    704691 
    705692Due to interactions and mass exchange of water ($F_{mass}$) with other Earth system components 
    706 (\ie atmosphere, sea-ice, land), the change in the heat and salt content of the surface layer of the ocean is due 
     693(\ie\ atmosphere, sea-ice, land), the change in the heat and salt content of the surface layer of the ocean is due 
    707694both to the heat and salt fluxes crossing the sea surface (not linked with $F_{mass}$) and 
    708695to the heat and salt content of the mass exchange. 
     
    716703\item 
    717704  $Q_{ns}$, the non-solar part of the net surface heat flux that crosses the sea surface 
    718   (\ie the difference between the total surface heat flux and the fraction of the short wave flux that 
     705  (\ie\ the difference between the total surface heat flux and the fraction of the short wave flux that 
    719706  penetrates into the water column, see \autoref{subsec:TRA_qsr}) 
    720707  plus the heat content associated with of the mass exchange with the atmosphere and lands. 
     
    758745          &\overline{(\textit{sfx} -        \textit{emp} \lt. S \rt|_{k = 1})}^t 
    759746  \end{alignedat} 
    760 \end{equation}  
     747\end{equation} 
    761748Note that an exact conservation of heat and salt content is only achieved with non-linear free surface. 
    762749In the linear free surface case, there is a small imbalance. 
     
    765752 
    766753% ------------------------------------------------------------------------------------------------------------- 
    767 %        Solar Radiation Penetration  
     754%        Solar Radiation Penetration 
    768755% ------------------------------------------------------------------------------------------------------------- 
    769756\subsection[Solar radiation penetration (\textit{traqsr.F90})] 
     
    775762%-------------------------------------------------------------------------------------------------------------- 
    776763 
    777 Options are defined through the \ngn{namtra\_qsr} namelist variables. 
    778 When the penetrative solar radiation option is used (\np{ln\_flxqsr}\forcode{ = .true.}), 
     764Options are defined through the \nam{tra\_qsr} namelist variables. 
     765When the penetrative solar radiation option is used (\np{ln\_traqsr}\forcode{ = .true.}), 
    779766the solar radiation penetrates the top few tens of meters of the ocean. 
    780 If it is not used (\np{ln\_flxqsr}\forcode{ = .false.}) all the heat flux is absorbed in the first ocean level. 
     767If it is not used (\np{ln\_traqsr}\forcode{ = .false.}) all the heat flux is absorbed in the first ocean level. 
    781768Thus, in the former case a term is added to the time evolution equation of temperature \autoref{eq:PE_tra_T} and 
    782 the surface boundary condition is modified to take into account only the non-penetrative part of the surface  
     769the surface boundary condition is modified to take into account only the non-penetrative part of the surface 
    783770heat flux: 
    784771\begin{equation} 
     
    789776  \end{gathered} 
    790777\end{equation} 
    791 where $Q_{sr}$ is the penetrative part of the surface heat flux (\ie the shortwave radiation) and 
     778where $Q_{sr}$ is the penetrative part of the surface heat flux (\ie\ the shortwave radiation) and 
    792779$I$ is the downward irradiance ($\lt. I \rt|_{z = \eta} = Q_{sr}$). 
    793780The additional term in \autoref{eq:PE_qsr} is discretized as follows: 
     
    803790(specified through namelist parameter \np{rn\_abs}). 
    804791It is assumed to penetrate the ocean with a decreasing exponential profile, with an e-folding depth scale, $\xi_0$, 
    805 of a few tens of centimetres (typically $\xi_0 = 0.35~m$ set as \np{rn\_si0} in the \ngn{namtra\_qsr} namelist). 
     792of a few tens of centimetres (typically $\xi_0 = 0.35~m$ set as \np{rn\_si0} in the \nam{tra\_qsr} namelist). 
    806793For shorter wavelengths (400-700~nm), the ocean is more transparent, and solar energy propagates to 
    807794larger depths where it contributes to local heating. 
     
    836823 
    837824The RGB formulation is used when \np{ln\_qsr\_rgb}\forcode{ = .true.}. 
    838 The RGB attenuation coefficients (\ie the inverses of the extinction length scales) are tabulated over 
     825The RGB attenuation coefficients (\ie\ the inverses of the extinction length scales) are tabulated over 
    83982661 nonuniform chlorophyll classes ranging from 0.01 to 10 g.Chl/L 
    840827(see the routine \rou{trc\_oce\_rgb} in \mdl{trc\_oce} module). 
     
    842829 
    843830\begin{description} 
    844 \item[\np{nn\_chdta}\forcode{ = 0}] 
    845   a constant 0.05 g.Chl/L value everywhere ;  
    846 \item[\np{nn\_chdta}\forcode{ = 1}] 
     831\item[\np{nn\_chldta}\forcode{ = 0}] 
     832  a constant 0.05 g.Chl/L value everywhere ; 
     833\item[\np{nn\_chldta}\forcode{ = 1}] 
    847834  an observed time varying chlorophyll deduced from satellite surface ocean color measurement spread uniformly in 
    848835  the vertical direction; 
    849 \item[\np{nn\_chdta}\forcode{ = 2}] 
     836\item[\np{nn\_chldta}\forcode{ = 2}] 
    850837  same as previous case except that a vertical profile of chlorophyl is used. 
    851838  Following \cite{morel.berthon_LO89}, the profile is computed from the local surface chlorophyll value; 
     
    853840  simulated time varying chlorophyll by TOP biogeochemical model. 
    854841  In this case, the RGB formulation is used to calculate both the phytoplankton light limitation in 
    855   PISCES or LOBSTER and the oceanic heating rate. 
    856 \end{description}  
     842  PISCES and the oceanic heating rate. 
     843\end{description} 
    857844 
    858845The trend in \autoref{eq:tra_qsr} associated with the penetration of the solar radiation is added to 
     
    862849the depth of $w-$levels does not significantly vary with location. 
    863850The level at which the light has been totally absorbed 
    864 (\ie it is less than the computer precision) is computed once, 
     851(\ie\ it is less than the computer precision) is computed once, 
    865852and the trend associated with the penetration of the solar radiation is only added down to that level. 
    866853Finally, note that when the ocean is shallow ($<$ 200~m), part of the solar radiation can reach the ocean floor. 
    867854In this case, we have chosen that all remaining radiation is absorbed in the last ocean level 
    868 (\ie $I$ is masked). 
     855(\ie\ $I$ is masked). 
    869856 
    870857%>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     
    889876%        Bottom Boundary Condition 
    890877% ------------------------------------------------------------------------------------------------------------- 
    891 \subsection[Bottom boundary condition (\textit{trabbc.F90})] 
     878\subsection[Bottom boundary condition (\textit{trabbc.F90}) - \forcode{ln_trabbc = .true.})] 
    892879{Bottom boundary condition (\protect\mdl{trabbc})} 
    893880\label{subsec:TRA_bbc} 
     
    910897 
    911898Usually it is assumed that there is no exchange of heat or salt through the ocean bottom, 
    912 \ie a no flux boundary condition is applied on active tracers at the bottom. 
     899\ie\ a no flux boundary condition is applied on active tracers at the bottom. 
    913900This is the default option in \NEMO, and it is implemented using the masking technique. 
    914901However, there is a non-zero heat flux across the seafloor that is associated with solid earth cooling. 
     
    916903but it warms systematically the ocean and acts on the densest water masses. 
    917904Taking this flux into account in a global ocean model increases the deepest overturning cell 
    918 (\ie the one associated with the Antarctic Bottom Water) by a few Sverdrups \citep{emile-geay.madec_OS09}. 
    919  
    920 Options are defined through the  \ngn{namtra\_bbc} namelist variables. 
     905(\ie\ the one associated with the Antarctic Bottom Water) by a few Sverdrups \citep{emile-geay.madec_OS09}. 
     906 
     907Options are defined through the \nam{bbc} namelist variables. 
    921908The presence of geothermal heating is controlled by setting the namelist parameter \np{ln\_trabbc} to true. 
    922909Then, when \np{nn\_geoflx} is set to 1, a constant geothermal heating is introduced whose value is given by 
    923 the \np{nn\_geoflx\_cst}, which is also a namelist parameter. 
     910the \np{rn\_geoflx\_cst}, which is also a namelist parameter. 
    924911When \np{nn\_geoflx} is set to 2, a spatially varying geothermal heat flux is introduced which is provided in 
    925912the \ifile{geothermal\_heating} NetCDF file (\autoref{fig:geothermal}) \citep{emile-geay.madec_OS09}. 
     
    928915% Bottom Boundary Layer 
    929916% ================================================================ 
    930 \section[Bottom boundary layer (\textit{trabbl.F90} - \texttt{\textbf{key\_trabbl}})] 
    931 {Bottom boundary layer (\protect\mdl{trabbl} - \protect\key{trabbl})} 
     917\section[Bottom boundary layer (\textit{trabbl.F90} - \forcode{ln_trabbl = .true.})] 
     918{Bottom boundary layer (\protect\mdl{trabbl} - \protect\np{ln\_trabbl}\forcode{ = .true.})} 
    932919\label{sec:TRA_bbl} 
    933920%--------------------------------------------nambbl--------------------------------------------------------- 
     
    936923%-------------------------------------------------------------------------------------------------------------- 
    937924 
    938 Options are defined through the \ngn{nambbl} namelist variables. 
     925Options are defined through the \nam{bbl} namelist variables. 
    939926In a $z$-coordinate configuration, the bottom topography is represented by a series of discrete steps. 
    940927This is not adequate to represent gravity driven downslope flows. 
     
    965952\label{subsec:TRA_bbl_diff} 
    966953 
    967 When applying sigma-diffusion (\key{trabbl} defined and \np{nn\_bbl\_ldf} set to 1), 
    968 the diffusive flux between two adjacent cells at the ocean floor is given by  
     954When applying sigma-diffusion (\np{ln\_trabbl}\forcode{ = .true.} and \np{nn\_bbl\_ldf} set to 1), 
     955the diffusive flux between two adjacent cells at the ocean floor is given by 
    969956\[ 
    970957  % \label{eq:tra_bbl_diff} 
     
    974961$A_l^\sigma$ the lateral diffusivity in the BBL. 
    975962Following \citet{beckmann.doscher_JPO97}, the latter is prescribed with a spatial dependence, 
    976 \ie in the conditional form 
     963\ie\ in the conditional form 
    977964\begin{equation} 
    978965  \label{eq:tra_bbl_coef} 
     
    990977(see green arrow in \autoref{fig:bbl}). 
    991978In practice, this constraint is applied separately in the two horizontal directions, 
    992 and the density gradient in \autoref{eq:tra_bbl_coef} is evaluated with the log gradient formulation:  
     979and the density gradient in \autoref{eq:tra_bbl_coef} is evaluated with the log gradient formulation: 
    993980\[ 
    994981  % \label{eq:tra_bbl_Drho} 
     
    10411028(see black arrow in \autoref{fig:bbl}) \citep{beckmann.doscher_JPO97}. 
    10421029It is a \textit{conditional advection}, that is, advection is allowed only 
    1043 if dense water overlies less dense water on the slope (\ie $\nabla_\sigma \rho \cdot \nabla H < 0$) and 
    1044 if the velocity is directed towards greater depth (\ie $\vect U \cdot \nabla H > 0$). 
     1030if dense water overlies less dense water on the slope (\ie\ $\nabla_\sigma \rho \cdot \nabla H < 0$) and 
     1031if the velocity is directed towards greater depth (\ie\ $\vect U \cdot \nabla H > 0$). 
    10451032 
    10461033\np{nn\_bbl\_adv}\forcode{ = 2}: 
     
    10481035the density difference between the higher cell and lower cell densities \citep{campin.goosse_T99}. 
    10491036The advection is allowed only  if dense water overlies less dense water on the slope 
    1050 (\ie $\nabla_\sigma \rho \cdot \nabla H < 0$). 
     1037(\ie\ $\nabla_\sigma \rho \cdot \nabla H < 0$). 
    10511038For example, the resulting transport of the downslope flow, here in the $i$-direction (\autoref{fig:bbl}), 
    10521039is simply given by the following expression: 
     
    10701057The advective BBL scheme modifies the tracer time tendency of the ocean cells near the topographic step by 
    10711058the downslope flow \autoref{eq:bbl_dw}, the horizontal \autoref{eq:bbl_hor} and 
    1072 the upward \autoref{eq:bbl_up} return flows as follows:  
     1059the upward \autoref{eq:bbl_up} return flows as follows: 
    10731060\begin{alignat}{3} 
    10741061  \label{eq:bbl_dw} 
     
    11081095    \pd[S]{t} = \cdots - \gamma (S - S_o) 
    11091096  \end{gathered} 
    1110 \end{equation}  
     1097\end{equation} 
    11111098where $\gamma$ is the inverse of a time scale, and $T_o$ and $S_o$ are given temperature and salinity fields 
    11121099(usually a climatology). 
    1113 Options are defined through the  \ngn{namtra\_dmp} namelist variables. 
     1100Options are defined through the  \nam{tra\_dmp} namelist variables. 
    11141101The restoring term is added when the namelist parameter \np{ln\_tradmp} is set to true. 
    1115 It also requires that both \np{ln\_tsd\_init} and \np{ln\_tsd\_tradmp} are set to true in 
    1116 \ngn{namtsd} namelist as well as \np{sn\_tem} and \np{sn\_sal} structures are correctly set 
    1117 (\ie that $T_o$ and $S_o$ are provided in input files and read using \mdl{fldread}, 
     1102It also requires that both \np{ln\_tsd\_init} and \np{ln\_tsd\_dmp} are set to true in 
     1103\nam{tsd} namelist as well as \np{sn\_tem} and \np{sn\_sal} structures are correctly set 
     1104(\ie\ that $T_o$ and $S_o$ are provided in input files and read using \mdl{fldread}, 
    11181105see \autoref{subsec:SBC_fldread}). 
    11191106The restoring coefficient $\gamma$ is a three-dimensional array read in during the \rou{tra\_dmp\_init} routine. 
     
    11571144%-------------------------------------------------------------------------------------------------------------- 
    11581145 
    1159 Options are defined through the \ngn{namdom} namelist variables. 
     1146Options are defined through the \nam{dom} namelist variables. 
    11601147The general framework for tracer time stepping is a modified leap-frog scheme \citep{leclair.madec_OM09}, 
    1161 \ie a three level centred time scheme associated with a Asselin time filter (cf. \autoref{sec:STP_mLF}): 
     1148\ie\ a three level centred time scheme associated with a Asselin time filter (cf. \autoref{sec:STP_mLF}): 
    11621149\begin{equation} 
    11631150  \label{eq:tra_nxt} 
     
    11651152    &(e_{3t}T)^{t + \rdt} &&= (e_{3t}T)_f^{t - \rdt} &&+ 2 \, \rdt \,e_{3t}^t \ \text{RHS}^t \\ 
    11661153    &(e_{3t}T)_f^t        &&= (e_{3t}T)^t            &&+ \, \gamma \, \lt[ (e_{3t}T)_f^{t - \rdt} - 2(e_{3t}T)^t + (e_{3t}T)^{t + \rdt} \rt] \\ 
    1167     &                     &&                         &&- \, \gamma \, \rdt \, \lt[ Q^{t + \rdt/2} - Q^{t - \rdt/2} \rt]   
     1154    &                     &&                         &&- \, \gamma \, \rdt \, \lt[ Q^{t + \rdt/2} - Q^{t - \rdt/2} \rt] 
    11681155  \end{alignedat} 
    1169 \end{equation}  
     1156\end{equation} 
    11701157where RHS is the right hand side of the temperature equation, the subscript $f$ denotes filtered values, 
    11711158$\gamma$ is the Asselin coefficient, and $S$ is the total forcing applied on $T$ 
    1172 (\ie fluxes plus content in mass exchanges). 
     1159(\ie\ fluxes plus content in mass exchanges). 
    11731160$\gamma$ is initialized as \np{rn\_atfp} (\textbf{namelist} parameter). 
    11741161Its default value is \np{rn\_atfp}\forcode{ = 10.e-3}. 
    11751162Note that the forcing correction term in the filter is not applied in linear free surface 
    1176 (\jp{lk\_vvl}\forcode{ = .false.}) (see \autoref{subsec:TRA_sbc}). 
     1163(\jp{ln\_linssh}\forcode{ = .true.}) (see \autoref{subsec:TRA_sbc}). 
    11771164Not also that in constant volume case, the time stepping is performed on $T$, not on its content, $e_{3t}T$. 
    11781165 
     
    11851172 
    11861173% ================================================================ 
    1187 % Equation of State (eosbn2)  
     1174% Equation of State (eosbn2) 
    11881175% ================================================================ 
    11891176\section[Equation of state (\textit{eosbn2.F90})] 
     
    11981185%        Equation of State 
    11991186% ------------------------------------------------------------------------------------------------------------- 
    1200 \subsection[Equation of seawater (\forcode{nn_eos = {-1,1}})] 
    1201 {Equation of seawater (\protect\np{nn\_eos}\forcode{ = {-1,1}})} 
     1187\subsection[Equation of seawater (\texttt{ln}\{\texttt{\_teso10,\_eos80,\_seos}\})] 
     1188{Equation of seawater (\protect\np{ln\_teos10}, \protect\np{ln\_teos80}, or \protect\np{ln\_seos}) } 
    12021189\label{subsec:TRA_eos} 
     1190 
    12031191 
    12041192The Equation Of Seawater (EOS) is an empirical nonlinear thermodynamic relationship linking seawater density, 
     
    12171205\textit{(ii)}  it is more accurate, being based on an updated database of laboratory measurements, and 
    12181206\textit{(iii)} it uses Conservative Temperature and Absolute Salinity (instead of potential temperature and 
    1219 practical salinity for EOS-980, both variables being more suitable for use as model variables 
     1207practical salinity for EOS-80, both variables being more suitable for use as model variables 
    12201208\citep{ioc.iapso_bk10, graham.mcdougall_JPO13}. 
    1221 EOS-80 is an obsolescent feature of the NEMO system, kept only for backward compatibility. 
     1209EOS-80 is an obsolescent feature of the \NEMO\ system, kept only for backward compatibility. 
    12221210For process studies, it is often convenient to use an approximation of the EOS. 
    12231211To that purposed, a simplified EOS (S-EOS) inspired by \citet{vallis_bk06} is also available. 
     
    12291217density in the World Ocean varies by no more than 2$\%$ from that value \citep{gill_bk82}. 
    12301218 
    1231 Options are defined through the \ngn{nameos} namelist variables, and in particular \np{nn\_eos} which 
    1232 controls the EOS used (\forcode{= -1} for TEOS10 ; \forcode{= 0} for EOS-80 ; \forcode{= 1} for S-EOS). 
     1219Options which control the EOS used are defined through the \ngn{nameos} namelist variables. 
    12331220 
    12341221\begin{description} 
    1235 \item[\np{nn\_eos}\forcode{ = -1}] 
     1222\item[\np{ln\_teos10}\forcode{ = .true.}] 
    12361223  the polyTEOS10-bsq equation of seawater \citep{roquet.madec.ea_OM15} is used. 
    12371224  The accuracy of this approximation is comparable to the TEOS-10 rational function approximation, 
     
    12491236  In particular, the initial state deined by the user have to be given as \textit{Conservative} Temperature and 
    12501237  \textit{Absolute} Salinity. 
    1251   In addition, setting \np{ln\_useCT} to \forcode{.true.} convert the Conservative SST to potential SST prior to 
     1238  In addition, when using TEOS10, the Conservative SST is converted to potential SST prior to 
    12521239  either computing the air-sea and ice-sea fluxes (forced mode) or 
    12531240  sending the SST field to the atmosphere (coupled mode). 
    1254 \item[\np{nn\_eos}\forcode{ = 0}] 
     1241\item[\np{ln\_eos80}\forcode{ = .true.}] 
    12551242  the polyEOS80-bsq equation of seawater is used. 
    12561243  It takes the same polynomial form as the polyTEOS10, but the coefficients have been optimized to 
     
    12641251  Nevertheless, a severe assumption is made in order to have a heat content ($C_p T_p$) which 
    12651252  is conserved by the model: $C_p$ is set to a constant value, the TEOS10 value. 
    1266 \item[\np{nn\_eos}\forcode{ = 1}] 
     1253\item[\np{ln\_seos}\forcode{ = .true.}] 
    12671254  a simplified EOS (S-EOS) inspired by \citet{vallis_bk06} is chosen, 
    12681255  the coefficients of which has been optimized to fit the behavior of TEOS10 
     
    12741261  as well as between \textit{absolute} and \textit{practical} salinity. 
    12751262  S-EOS takes the following expression: 
     1263 
    12761264  \begin{gather*} 
    12771265    % \label{eq:tra_S-EOS} 
    12781266    \begin{alignedat}{2} 
    12791267    &d_a(T,S,z) = \frac{1}{\rho_o} \big[ &- a_0 \; ( 1 + 0.5 \; \lambda_1 \; T_a + \mu_1 \; z ) * &T_a \big. \\ 
    1280     &                                    &+ b_0 \; ( 1 - 0.5 \; \lambda_2 \; S_a - \mu_2 \; z ) * &S_a       \\   
     1268    &                                    &+ b_0 \; ( 1 - 0.5 \; \lambda_2 \; S_a - \mu_2 \; z ) * &S_a       \\ 
    12811269    &                              \big. &- \nu \;                           T_a                  &S_a \big] \\ 
    12821270    \end{alignedat} 
     
    13261314%        Brunt-V\"{a}is\"{a}l\"{a} Frequency 
    13271315% ------------------------------------------------------------------------------------------------------------- 
    1328 \subsection[Brunt-V\"{a}is\"{a}l\"{a} frequency (\forcode{nn_eos = [0-2]})] 
    1329 {Brunt-V\"{a}is\"{a}l\"{a} frequency (\protect\np{nn\_eos}\forcode{ = [0-2]})} 
     1316\subsection[Brunt-V\"{a}is\"{a}l\"{a} frequency] 
     1317{Brunt-V\"{a}is\"{a}l\"{a} frequency} 
    13301318\label{subsec:TRA_bn2} 
    13311319 
     
    13361324In particular, $N^2$ has to be computed at the local pressure 
    13371325(pressure in decibar being approximated by the depth in meters). 
    1338 The expression for $N^2$  is given by:  
     1326The expression for $N^2$  is given by: 
    13391327\[ 
    13401328  % \label{eq:tra_bn2} 
     
    13581346  \begin{split} 
    13591347    &T_f (S,p) = \lt( a + b \, \sqrt{S} + c \, S \rt) \, S + d \, p \\ 
    1360     &\text{where~} a = -0.0575, \, b = 1.710523~10^{-3}, \, c = -2.154996~10^{-4} \\  
     1348    &\text{where~} a = -0.0575, \, b = 1.710523~10^{-3}, \, c = -2.154996~10^{-4} \\ 
    13611349    &\text{and~} d = -7.53~10^{-3} 
    13621350    \end{split} 
     
    13641352 
    13651353\autoref{eq:tra_eos_fzp} is only used to compute the potential freezing point of sea water 
    1366 (\ie referenced to the surface $p = 0$), 
     1354(\ie\ referenced to the surface $p = 0$), 
    13671355thus the pressure dependent terms in \autoref{eq:tra_eos_fzp} (last term) have been dropped. 
    13681356The freezing point is computed through \textit{eos\_fzp}, 
     
    13701358 
    13711359% ------------------------------------------------------------------------------------------------------------- 
    1372 %        Potential Energy      
     1360%        Potential Energy 
    13731361% ------------------------------------------------------------------------------------------------------------- 
    13741362%\subsection{Potential Energy anomalies} 
     
    13791367 
    13801368% ================================================================ 
    1381 % Horizontal Derivative in zps-coordinate  
     1369% Horizontal Derivative in zps-coordinate 
    13821370% ================================================================ 
    13831371\section[Horizontal derivative in \textit{zps}-coordinate (\textit{zpshde.F90})] 
     
    13851373\label{sec:TRA_zpshde} 
    13861374 
    1387 \gmcomment{STEVEN: to be consistent with earlier discussion of differencing and averaging operators,  
     1375\gmcomment{STEVEN: to be consistent with earlier discussion of differencing and averaging operators, 
    13881376I've changed "derivative" to "difference" and "mean" to "average"} 
    13891377 
     
    14271415  \rt. 
    14281416\] 
    1429 and the resulting forms for the horizontal difference and the horizontal average value of $T$ at a $U$-point are:  
     1417and the resulting forms for the horizontal difference and the horizontal average value of $T$ at a $U$-point are: 
    14301418\begin{equation} 
    14311419  \label{eq:zps_hde} 
     
    14531441Instead of forming a linear approximation of density, we compute $\widetilde \rho$ from the interpolated values of 
    14541442$T$ and $S$, and the pressure at a $u$-point 
    1455 (in the equation of state pressure is approximated by depth, see \autoref{subsec:TRA_eos}):  
     1443(in the equation of state pressure is approximated by depth, see \autoref{subsec:TRA_eos}): 
    14561444\[ 
    14571445  % \label{eq:zps_hde_rho} 
Note: See TracChangeset for help on using the changeset viewer.