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 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/DOC/TexFiles/Chapters/Chap_LDF.tex – NEMO

Ignore:
Timestamp:
2016-01-08T10:35:19+01:00 (8 years ago)
Author:
jamesharle
Message:

Update MPP_BDY_UPDATE branch to be consistent with head of trunk

File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/DOC/TexFiles/Chapters/Chap_LDF.tex

    r4147 r6225  
    11 
    22% ================================================================ 
    3 % Chapter � Lateral Ocean Physics (LDF) 
     3% Chapter ——— Lateral Ocean Physics (LDF) 
    44% ================================================================ 
    55\chapter{Lateral Ocean Physics (LDF)} 
    … …  
    1515described in \S\ref{PE_zdf} and their discrete formulation in \S\ref{TRA_ldf}  
    1616and \S\ref{DYN_ldf}). In this section we further discuss each lateral physics option.  
    17 Choosing one lateral physics scheme means for the user defining, (1) the space  
    18 and time variations of the eddy coefficients ; (2) the direction along which the  
    19 lateral diffusive fluxes are evaluated (model level, geopotential or isopycnal  
    20 surfaces); and (3) the type of operator used (harmonic, or biharmonic operators,  
    21 and for tracers only, eddy induced advection on tracers). These three aspects  
    22 of the lateral diffusion are set through namelist parameters and CPP keys  
    23 (see the \textit{\ngn{nam\_traldf}} and \textit{\ngn{nam\_dynldf}} below). Note 
    24 that this chapter describes the default implementation of iso-neutral 
     17Choosing one lateral physics scheme means for the user defining,  
     18(1) the type of operator used (laplacian or bilaplacian operators, or no lateral mixing term) ;  
     19(2) the direction along which the lateral diffusive fluxes are evaluated (model level, geopotential or isopycnal surfaces) ; and  
     20(3) the space and time variations of the eddy coefficients.  
     21These three aspects of the lateral diffusion are set through namelist parameters  
     22(see the \textit{\ngn{nam\_traldf}} and \textit{\ngn{nam\_dynldf}} below).  
     23Note that this chapter describes the standard implementation of iso-neutral 
    2524tracer mixing, and Griffies's implementation, which is used if 
    2625\np{traldf\_grif}=true, is described in Appdx\ref{sec:triad} 
    … …  
    3332 
    3433% ================================================================ 
     34% Direction of lateral Mixing 
     35% ================================================================ 
     36\section  [Direction of Lateral Mixing (\textit{ldfslp})] 
     37      {Direction of Lateral Mixing (\mdl{ldfslp})} 
     38\label{LDF_slp} 
     39 
     40%%% 
     41\gmcomment{  we should emphasize here that the implementation is a rather old one.  
     42Better work can be achieved by using \citet{Griffies_al_JPO98, Griffies_Bk04} iso-neutral scheme. } 
     43 
     44A direction for lateral mixing has to be defined when the desired operator does  
     45not act along the model levels. This occurs when $(a)$ horizontal mixing is  
     46required on tracer or momentum (\np{ln\_traldf\_hor} or \np{ln\_dynldf\_hor})  
     47in $s$- or mixed $s$-$z$- coordinates, and $(b)$ isoneutral mixing is required  
     48whatever the vertical coordinate is. This direction of mixing is defined by its  
     49slopes in the \textbf{i}- and \textbf{j}-directions at the face of the cell of the  
     50quantity to be diffused. For a tracer, this leads to the following four slopes :  
     51$r_{1u}$, $r_{1w}$, $r_{2v}$, $r_{2w}$ (see \eqref{Eq_tra_ldf_iso}), while  
     52for momentum the slopes are  $r_{1t}$, $r_{1uw}$, $r_{2f}$, $r_{2uw}$ for  
     53$u$ and  $r_{1f}$, $r_{1vw}$, $r_{2t}$, $r_{2vw}$ for $v$.  
     54 
     55%gm% add here afigure of the slope in i-direction 
     56 
     57\subsection{slopes for tracer geopotential mixing in the $s$-coordinate} 
     58 
     59In $s$-coordinates, geopotential mixing ($i.e.$ horizontal mixing) $r_1$ and  
     60$r_2$ are the slopes between the geopotential and computational surfaces.  
     61Their discrete formulation is found by locally solving \eqref{Eq_tra_ldf_iso}  
     62when the diffusive fluxes in the three directions are set to zero and $T$ is  
     63assumed to be horizontally uniform, $i.e.$ a linear function of $z_T$, the  
     64depth of a $T$-point.  
     65%gm { Steven : My version is obviously wrong since I'm left with an arbitrary constant which is the local vertical temperature gradient} 
     66 
     67\begin{equation} \label{Eq_ldfslp_geo} 
     68\begin{aligned} 
     69 r_{1u} &= \frac{e_{3u}}{ \left( e_{1u}\;\overline{\overline{e_{3w}}}^{\,i+1/2,\,k} \right)} 
     70           \;\delta_{i+1/2}[z_t]  
     71      &\approx \frac{1}{e_{1u}}\; \delta_{i+1/2}[z_t]  
     72\\ 
     73 r_{2v} &= \frac{e_{3v}}{\left( e_{2v}\;\overline{\overline{e_{3w}}}^{\,j+1/2,\,k} \right)}  
     74           \;\delta_{j+1/2} [z_t]  
     75      &\approx \frac{1}{e_{2v}}\; \delta_{j+1/2}[z_t]  
     76\\ 
     77 r_{1w} &= \frac{1}{e_{1w}}\;\overline{\overline{\delta_{i+1/2}[z_t]}}^{\,i,\,k+1/2} 
     78      &\approx \frac{1}{e_{1w}}\; \delta_{i+1/2}[z_{uw}]  
     79 \\ 
     80 r_{2w} &= \frac{1}{e_{2w}}\;\overline{\overline{\delta_{j+1/2}[z_t]}}^{\,j,\,k+1/2} 
     81      &\approx \frac{1}{e_{2w}}\; \delta_{j+1/2}[z_{vw}]  
     82 \\ 
     83\end{aligned} 
     84\end{equation} 
     85 
     86%gm%  caution I'm not sure the simplification was a good idea!  
     87 
     88These slopes are computed once in \rou{ldfslp\_init} when \np{ln\_sco}=True,  
     89and either \np{ln\_traldf\_hor}=True or \np{ln\_dynldf\_hor}=True.  
     90 
     91\subsection{Slopes for tracer iso-neutral mixing}\label{LDF_slp_iso} 
     92In iso-neutral mixing  $r_1$ and $r_2$ are the slopes between the iso-neutral  
     93and computational surfaces. Their formulation does not depend on the vertical  
     94coordinate used. Their discrete formulation is found using the fact that the  
     95diffusive fluxes of locally referenced potential density ($i.e.$ $in situ$ density)  
     96vanish. So, substituting $T$ by $\rho$ in \eqref{Eq_tra_ldf_iso} and setting the  
     97diffusive fluxes in the three directions to zero leads to the following definition for  
     98the neutral slopes: 
     99 
     100\begin{equation} \label{Eq_ldfslp_iso} 
     101\begin{split} 
     102 r_{1u} &= \frac{e_{3u}}{e_{1u}}\; \frac{\delta_{i+1/2}[\rho]} 
     103                        {\overline{\overline{\delta_{k+1/2}[\rho]}}^{\,i+1/2,\,k}} 
     104\\ 
     105 r_{2v} &= \frac{e_{3v}}{e_{2v}}\; \frac{\delta_{j+1/2}\left[\rho \right]} 
     106                        {\overline{\overline{\delta_{k+1/2}[\rho]}}^{\,j+1/2,\,k}} 
     107\\ 
     108 r_{1w} &= \frac{e_{3w}}{e_{1w}}\;  
     109         \frac{\overline{\overline{\delta_{i+1/2}[\rho]}}^{\,i,\,k+1/2}} 
     110             {\delta_{k+1/2}[\rho]} 
     111\\ 
     112 r_{2w} &= \frac{e_{3w}}{e_{2w}}\;  
     113         \frac{\overline{\overline{\delta_{j+1/2}[\rho]}}^{\,j,\,k+1/2}} 
     114             {\delta_{k+1/2}[\rho]} 
     115\\ 
     116\end{split} 
     117\end{equation} 
     118 
     119%gm% rewrite this as the explanation is not very clear !!! 
     120%In practice, \eqref{Eq_ldfslp_iso} is of little help in evaluating the neutral surface slopes. Indeed, for an unsimplified equation of state, the density has a strong dependancy on pressure (here approximated as the depth), therefore applying \eqref{Eq_ldfslp_iso} using the $in situ$ density, $\rho$, computed at T-points leads to a flattening of slopes as the depth increases. This is due to the strong increase of the $in situ$ density with depth.  
     121 
     122%By definition, neutral surfaces are tangent to the local $in situ$ density \citep{McDougall1987}, therefore in \eqref{Eq_ldfslp_iso}, all the derivatives have to be evaluated at the same local pressure (which in decibars is approximated by the depth in meters). 
     123 
     124%In the $z$-coordinate, the derivative of the  \eqref{Eq_ldfslp_iso} numerator is evaluated at the same depth \nocite{as what?} ($T$-level, which is the same as the $u$- and $v$-levels), so  the $in situ$ density can be used for its evaluation.  
     125 
     126As the mixing is performed along neutral surfaces, the gradient of $\rho$ in  
     127\eqref{Eq_ldfslp_iso} has to be evaluated at the same local pressure (which,  
     128in decibars, is approximated by the depth in meters in the model). Therefore  
     129\eqref{Eq_ldfslp_iso} cannot be used as such, but further transformation is  
     130needed depending on the vertical coordinate used: 
     131 
     132\begin{description} 
     133 
     134\item[$z$-coordinate with full step : ] in \eqref{Eq_ldfslp_iso} the densities  
     135appearing in the $i$ and $j$ derivatives  are taken at the same depth, thus  
     136the $in situ$ density can be used. This is not the case for the vertical  
     137derivatives: $\delta_{k+1/2}[\rho]$ is replaced by $-\rho N^2/g$, where $N^2$  
     138is the local Brunt-Vais\"{a}l\"{a} frequency evaluated following  
     139\citet{McDougall1987} (see \S\ref{TRA_bn2}).  
     140 
     141\item[$z$-coordinate with partial step : ] this case is identical to the full step  
     142case except that at partial step level, the \emph{horizontal} density gradient  
     143is evaluated as described in \S\ref{TRA_zpshde}. 
     144 
     145\item[$s$- or hybrid $s$-$z$- coordinate : ] in the current release of \NEMO,  
     146iso-neutral mixing is only employed for $s$-coordinates if the 
     147Griffies scheme is used (\np{traldf\_grif}=true; see Appdx \ref{sec:triad}).  
     148In other words, iso-neutral mixing will only be accurately represented with a  
     149linear equation of state (\np{nn\_eos}=1 or 2). In the case of a "true" equation  
     150of state, the evaluation of $i$ and $j$ derivatives in \eqref{Eq_ldfslp_iso}  
     151will include a pressure dependent part, leading to the wrong evaluation of  
     152the neutral slopes. 
     153 
     154%gm%  
     155Note: The solution for $s$-coordinate passes trough the use of different  
     156(and better) expression for the constraint on iso-neutral fluxes. Following  
     157\citet{Griffies_Bk04}, instead of specifying directly that there is a zero neutral  
     158diffusive flux of locally referenced potential density, we stay in the $T$-$S$  
     159plane and consider the balance between the neutral direction diffusive fluxes  
     160of potential temperature and salinity: 
     161\begin{equation} 
     162\alpha \ \textbf{F}(T) = \beta \ \textbf{F}(S) 
     163\end{equation} 
     164%gm{  where vector F is ....} 
     165 
     166This constraint leads to the following definition for the slopes: 
     167 
     168\begin{equation} \label{Eq_ldfslp_iso2} 
     169\begin{split} 
     170 r_{1u} &= \frac{e_{3u}}{e_{1u}}\; \frac 
     171      {\alpha_u \;\delta_{i+1/2}[T] - \beta_u \;\delta_{i+1/2}[S]} 
     172      {\alpha_u \;\overline{\overline{\delta_{k+1/2}[T]}}^{\,i+1/2,\,k} 
     173       -\beta_u  \;\overline{\overline{\delta_{k+1/2}[S]}}^{\,i+1/2,\,k} } 
     174\\ 
     175 r_{2v} &= \frac{e_{3v}}{e_{2v}}\; \frac 
     176      {\alpha_v \;\delta_{j+1/2}[T] - \beta_v \;\delta_{j+1/2}[S]} 
     177      {\alpha_v \;\overline{\overline{\delta_{k+1/2}[T]}}^{\,j+1/2,\,k} 
     178       -\beta_v  \;\overline{\overline{\delta_{k+1/2}[S]}}^{\,j+1/2,\,k} } 
     179\\ 
     180 r_{1w} &= \frac{e_{3w}}{e_{1w}}\; \frac 
     181      {\alpha_w \;\overline{\overline{\delta_{i+1/2}[T]}}^{\,i,\,k+1/2} 
     182       -\beta_w  \;\overline{\overline{\delta_{i+1/2}[S]}}^{\,i,\,k+1/2} } 
     183      {\alpha_w \;\delta_{k+1/2}[T] - \beta_w \;\delta_{k+1/2}[S]} 
     184\\ 
     185 r_{2w} &= \frac{e_{3w}}{e_{2w}}\; \frac 
     186      {\alpha_w \;\overline{\overline{\delta_{j+1/2}[T]}}^{\,j,\,k+1/2} 
     187       -\beta_w  \;\overline{\overline{\delta_{j+1/2}[S]}}^{\,j,\,k+1/2} } 
     188      {\alpha_w \;\delta_{k+1/2}[T] - \beta_w \;\delta_{k+1/2}[S]} 
     189\\ 
     190\end{split} 
     191\end{equation} 
     192where $\alpha$ and $\beta$, the thermal expansion and saline contraction  
     193coefficients introduced in \S\ref{TRA_bn2}, have to be evaluated at the three  
     194velocity points. In order to save computation time, they should be approximated  
     195by the mean of their values at $T$-points (for example in the case of $\alpha$:   
     196$\alpha_u=\overline{\alpha_T}^{i+1/2}$,  $\alpha_v=\overline{\alpha_T}^{j+1/2}$  
     197and $\alpha_w=\overline{\alpha_T}^{k+1/2}$). 
     198 
     199Note that such a formulation could be also used in the $z$-coordinate and  
     200$z$-coordinate with partial steps cases. 
     201 
     202\end{description} 
     203 
     204This implementation is a rather old one. It is similar to the one 
     205proposed by Cox [1987], except for the background horizontal 
     206diffusion. Indeed, the Cox implementation of isopycnal diffusion in 
     207GFDL-type models requires a minimum background horizontal diffusion 
     208for numerical stability reasons.  To overcome this problem, several 
     209techniques have been proposed in which the numerical schemes of the 
     210ocean model are modified \citep{Weaver_Eby_JPO97, 
     211  Griffies_al_JPO98}. Griffies's scheme is now available in \NEMO if 
     212\np{traldf\_grif\_iso} is set true; see Appdx \ref{sec:triad}. Here, 
     213another strategy is presented \citep{Lazar_PhD97}: a local 
     214filtering of the iso-neutral slopes (made on 9 grid-points) prevents 
     215the development of grid point noise generated by the iso-neutral 
     216diffusion operator (Fig.~\ref{Fig_LDF_ZDF1}). This allows an 
     217iso-neutral diffusion scheme without additional background horizontal 
     218mixing. This technique can be viewed as a diffusion operator that acts 
     219along large-scale (2~$\Delta$x) \gmcomment{2deltax doesnt seem very 
     220  large scale} iso-neutral surfaces. The diapycnal diffusion required 
     221for numerical stability is thus minimized and its net effect on the 
     222flow is quite small when compared to the effect of an horizontal 
     223background mixing. 
     224 
     225Nevertheless, this iso-neutral operator does not ensure that variance cannot increase,  
     226contrary to the \citet{Griffies_al_JPO98} operator which has that property.  
     227 
     228%>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     229\begin{figure}[!ht]      \begin{center} 
     230\includegraphics[width=0.70\textwidth]{./TexFiles/Figures/Fig_LDF_ZDF1.pdf} 
     231\caption {    \label{Fig_LDF_ZDF1} 
     232averaging procedure for isopycnal slope computation.} 
     233\end{center}    \end{figure} 
     234%>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     235 
     236%There are three additional questions about the slope calculation.  
     237%First the expression for the rotation tensor has been obtain assuming the "small slope" approximation, so a bound has to be imposed on slopes.  
     238%Second, numerical stability issues also require a bound on slopes.  
     239%Third, the question of boundary condition specified on slopes... 
     240 
     241%from griffies: chapter 13.1.... 
     242 
     243 
     244 
     245% In addition and also for numerical stability reasons \citep{Cox1987, Griffies_Bk04},  
     246% the slopes are bounded by $1/100$ everywhere. This limit is decreasing linearly  
     247% to zero fom $70$ meters depth and the surface (the fact that the eddies "feel" the  
     248% surface motivates this flattening of isopycnals near the surface). 
     249 
     250For numerical stability reasons \citep{Cox1987, Griffies_Bk04}, the slopes must also  
     251be bounded by $1/100$ everywhere. This constraint is applied in a piecewise linear  
     252fashion, increasing from zero at the surface to $1/100$ at $70$ metres and thereafter  
     253decreasing to zero at the bottom of the ocean. (the fact that the eddies "feel" the  
     254surface motivates this flattening of isopycnals near the surface). 
     255 
     256%>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     257\begin{figure}[!ht]     \begin{center} 
     258\includegraphics[width=0.70\textwidth]{./TexFiles/Figures/Fig_eiv_slp.pdf} 
     259\caption {     \label{Fig_eiv_slp} 
     260Vertical profile of the slope used for lateral mixing in the mixed layer :  
     261\textit{(a)} in the real ocean the slope is the iso-neutral slope in the ocean interior,  
     262which has to be adjusted at the surface boundary (i.e. it must tend to zero at the  
     263surface since there is no mixing across the air-sea interface: wall boundary  
     264condition). Nevertheless, the profile between the surface zero value and the interior  
     265iso-neutral one is unknown, and especially the value at the base of the mixed layer ;  
     266\textit{(b)} profile of slope using a linear tapering of the slope near the surface and  
     267imposing a maximum slope of 1/100 ; \textit{(c)} profile of slope actually used in  
     268\NEMO: a linear decrease of the slope from zero at the surface to its ocean interior  
     269value computed just below the mixed layer. Note the huge change in the slope at the  
     270base of the mixed layer between  \textit{(b)}  and \textit{(c)}.} 
     271\end{center}   \end{figure} 
     272%>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     273 
     274\colorbox{yellow}{add here a discussion about the flattening of the slopes, vs  tapering the coefficient.} 
     275 
     276\subsection{slopes for momentum iso-neutral mixing} 
     277 
     278The iso-neutral diffusion operator on momentum is the same as the one used on  
     279tracers but applied to each component of the velocity separately (see  
     280\eqref{Eq_dyn_ldf_iso} in section~\ref{DYN_ldf_iso}). The slopes between the  
     281surface along which the diffusion operator acts and the surface of computation  
     282($z$- or $s$-surfaces) are defined at $T$-, $f$-, and \textit{uw}- points for the  
     283$u$-component, and $T$-, $f$- and \textit{vw}- points for the $v$-component.  
     284They are computed from the slopes used for tracer diffusion, $i.e.$  
     285\eqref{Eq_ldfslp_geo} and \eqref{Eq_ldfslp_iso} : 
     286 
     287\begin{equation} \label{Eq_ldfslp_dyn} 
     288\begin{aligned} 
     289&r_{1t}\ \ = \overline{r_{1u}}^{\,i}       &&&    r_{1f}\ \ &= \overline{r_{1u}}^{\,i+1/2} \\ 
     290&r_{2f} \ \ = \overline{r_{2v}}^{\,j+1/2} &&&   r_{2t}\ &= \overline{r_{2v}}^{\,j} \\ 
     291&r_{1uw}  = \overline{r_{1w}}^{\,i+1/2} &&\ \ \text{and} \ \ &   r_{1vw}&= \overline{r_{1w}}^{\,j+1/2} \\ 
     292&r_{2uw}= \overline{r_{2w}}^{\,j+1/2} &&&         r_{2vw}&= \overline{r_{2w}}^{\,j+1/2}\\ 
     293\end{aligned} 
     294\end{equation} 
     295 
     296The major issue remaining is in the specification of the boundary conditions.  
     297The same boundary conditions are chosen as those used for lateral  
     298diffusion along model level surfaces, i.e. using the shear computed along  
     299the model levels and with no additional friction at the ocean bottom (see  
     300{\S\ref{LBC_coast}). 
     301 
     302 
     303% ================================================================ 
     304% Lateral Mixing Operator 
     305% ================================================================ 
     306\section [Lateral Mixing Operators (\textit{ldftra}, \textit{ldfdyn})]  
     307        {Lateral Mixing Operators (\mdl{traldf}, \mdl{traldf}) } 
     308\label{LDF_op} 
     309 
     310 
     311    
     312% ================================================================ 
    35313% Lateral Mixing Coefficients 
    36314% ================================================================ 
    … …  
    38316        {Lateral Mixing Coefficient (\mdl{ldftra}, \mdl{ldfdyn}) } 
    39317\label{LDF_coef} 
    40  
    41318 
    42319Introducing a space variation in the lateral eddy mixing coefficients changes  
    … …  
    165442 
    166443% ================================================================ 
    167 % Direction of lateral Mixing 
    168 % ================================================================ 
    169 \section  [Direction of Lateral Mixing (\textit{ldfslp})] 
    170       {Direction of Lateral Mixing (\mdl{ldfslp})} 
    171 \label{LDF_slp} 
    172  
    173 %%% 
    174 \gmcomment{  we should emphasize here that the implementation is a rather old one.  
    175 Better work can be achieved by using \citet{Griffies_al_JPO98, Griffies_Bk04} iso-neutral scheme. } 
    176  
    177 A direction for lateral mixing has to be defined when the desired operator does  
    178 not act along the model levels. This occurs when $(a)$ horizontal mixing is  
    179 required on tracer or momentum (\np{ln\_traldf\_hor} or \np{ln\_dynldf\_hor})  
    180 in $s$- or mixed $s$-$z$- coordinates, and $(b)$ isoneutral mixing is required  
    181 whatever the vertical coordinate is. This direction of mixing is defined by its  
    182 slopes in the \textbf{i}- and \textbf{j}-directions at the face of the cell of the  
    183 quantity to be diffused. For a tracer, this leads to the following four slopes :  
    184 $r_{1u}$, $r_{1w}$, $r_{2v}$, $r_{2w}$ (see \eqref{Eq_tra_ldf_iso}), while  
    185 for momentum the slopes are  $r_{1t}$, $r_{1uw}$, $r_{2f}$, $r_{2uw}$ for  
    186 $u$ and  $r_{1f}$, $r_{1vw}$, $r_{2t}$, $r_{2vw}$ for $v$.  
    187  
    188 %gm% add here afigure of the slope in i-direction 
    189  
    190 \subsection{slopes for tracer geopotential mixing in the $s$-coordinate} 
    191  
    192 In $s$-coordinates, geopotential mixing ($i.e.$ horizontal mixing) $r_1$ and  
    193 $r_2$ are the slopes between the geopotential and computational surfaces.  
    194 Their discrete formulation is found by locally solving \eqref{Eq_tra_ldf_iso}  
    195 when the diffusive fluxes in the three directions are set to zero and $T$ is  
    196 assumed to be horizontally uniform, $i.e.$ a linear function of $z_T$, the  
    197 depth of a $T$-point.  
    198 %gm { Steven : My version is obviously wrong since I'm left with an arbitrary constant which is the local vertical temperature gradient} 
    199  
    200 \begin{equation} \label{Eq_ldfslp_geo} 
    201 \begin{aligned} 
    202  r_{1u} &= \frac{e_{3u}}{ \left( e_{1u}\;\overline{\overline{e_{3w}}}^{\,i+1/2,\,k} \right)} 
    203            \;\delta_{i+1/2}[z_t]  
    204       &\approx \frac{1}{e_{1u}}\; \delta_{i+1/2}[z_t]  
    205 \\ 
    206  r_{2v} &= \frac{e_{3v}}{\left( e_{2v}\;\overline{\overline{e_{3w}}}^{\,j+1/2,\,k} \right)}  
    207            \;\delta_{j+1/2} [z_t]  
    208       &\approx \frac{1}{e_{2v}}\; \delta_{j+1/2}[z_t]  
    209 \\ 
    210  r_{1w} &= \frac{1}{e_{1w}}\;\overline{\overline{\delta_{i+1/2}[z_t]}}^{\,i,\,k+1/2} 
    211       &\approx \frac{1}{e_{1w}}\; \delta_{i+1/2}[z_{uw}]  
    212  \\ 
    213  r_{2w} &= \frac{1}{e_{2w}}\;\overline{\overline{\delta_{j+1/2}[z_t]}}^{\,j,\,k+1/2} 
    214       &\approx \frac{1}{e_{2w}}\; \delta_{j+1/2}[z_{vw}]  
    215  \\ 
    216 \end{aligned} 
    217 \end{equation} 
    218  
    219 %gm%  caution I'm not sure the simplification was a good idea!  
    220  
    221 These slopes are computed once in \rou{ldfslp\_init} when \np{ln\_sco}=True,  
    222 and either \np{ln\_traldf\_hor}=True or \np{ln\_dynldf\_hor}=True.  
    223  
    224 \subsection{Slopes for tracer iso-neutral mixing}\label{LDF_slp_iso} 
    225 In iso-neutral mixing  $r_1$ and $r_2$ are the slopes between the iso-neutral  
    226 and computational surfaces. Their formulation does not depend on the vertical  
    227 coordinate used. Their discrete formulation is found using the fact that the  
    228 diffusive fluxes of locally referenced potential density ($i.e.$ $in situ$ density)  
    229 vanish. So, substituting $T$ by $\rho$ in \eqref{Eq_tra_ldf_iso} and setting the  
    230 diffusive fluxes in the three directions to zero leads to the following definition for  
    231 the neutral slopes: 
    232  
    233 \begin{equation} \label{Eq_ldfslp_iso} 
    234 \begin{split} 
    235  r_{1u} &= \frac{e_{3u}}{e_{1u}}\; \frac{\delta_{i+1/2}[\rho]} 
    236                         {\overline{\overline{\delta_{k+1/2}[\rho]}}^{\,i+1/2,\,k}} 
    237 \\ 
    238  r_{2v} &= \frac{e_{3v}}{e_{2v}}\; \frac{\delta_{j+1/2}\left[\rho \right]} 
    239                         {\overline{\overline{\delta_{k+1/2}[\rho]}}^{\,j+1/2,\,k}} 
    240 \\ 
    241  r_{1w} &= \frac{e_{3w}}{e_{1w}}\;  
    242          \frac{\overline{\overline{\delta_{i+1/2}[\rho]}}^{\,i,\,k+1/2}} 
    243              {\delta_{k+1/2}[\rho]} 
    244 \\ 
    245  r_{2w} &= \frac{e_{3w}}{e_{2w}}\;  
    246          \frac{\overline{\overline{\delta_{j+1/2}[\rho]}}^{\,j,\,k+1/2}} 
    247              {\delta_{k+1/2}[\rho]} 
    248 \\ 
    249 \end{split} 
    250 \end{equation} 
    251  
    252 %gm% rewrite this as the explanation is not very clear !!! 
    253 %In practice, \eqref{Eq_ldfslp_iso} is of little help in evaluating the neutral surface slopes. Indeed, for an unsimplified equation of state, the density has a strong dependancy on pressure (here approximated as the depth), therefore applying \eqref{Eq_ldfslp_iso} using the $in situ$ density, $\rho$, computed at T-points leads to a flattening of slopes as the depth increases. This is due to the strong increase of the $in situ$ density with depth.  
    254  
    255 %By definition, neutral surfaces are tangent to the local $in situ$ density \citep{McDougall1987}, therefore in \eqref{Eq_ldfslp_iso}, all the derivatives have to be evaluated at the same local pressure (which in decibars is approximated by the depth in meters). 
    256  
    257 %In the $z$-coordinate, the derivative of the  \eqref{Eq_ldfslp_iso} numerator is evaluated at the same depth \nocite{as what?} ($T$-level, which is the same as the $u$- and $v$-levels), so  the $in situ$ density can be used for its evaluation.  
    258  
    259 As the mixing is performed along neutral surfaces, the gradient of $\rho$ in  
    260 \eqref{Eq_ldfslp_iso} has to be evaluated at the same local pressure (which,  
    261 in decibars, is approximated by the depth in meters in the model). Therefore  
    262 \eqref{Eq_ldfslp_iso} cannot be used as such, but further transformation is  
    263 needed depending on the vertical coordinate used: 
    264  
    265 \begin{description} 
    266  
    267 \item[$z$-coordinate with full step : ] in \eqref{Eq_ldfslp_iso} the densities  
    268 appearing in the $i$ and $j$ derivatives  are taken at the same depth, thus  
    269 the $in situ$ density can be used. This is not the case for the vertical  
    270 derivatives: $\delta_{k+1/2}[\rho]$ is replaced by $-\rho N^2/g$, where $N^2$  
    271 is the local Brunt-Vais\"{a}l\"{a} frequency evaluated following  
    272 \citet{McDougall1987} (see \S\ref{TRA_bn2}).  
    273  
    274 \item[$z$-coordinate with partial step : ] this case is identical to the full step  
    275 case except that at partial step level, the \emph{horizontal} density gradient  
    276 is evaluated as described in \S\ref{TRA_zpshde}. 
    277  
    278 \item[$s$- or hybrid $s$-$z$- coordinate : ] in the current release of \NEMO,  
    279 iso-neutral mixing is only employed for $s$-coordinates if the 
    280 Griffies scheme is used (\np{traldf\_grif}=true; see Appdx \ref{sec:triad}).  
    281 In other words, iso-neutral mixing will only be accurately represented with a  
    282 linear equation of state (\np{nn\_eos}=1 or 2). In the case of a "true" equation  
    283 of state, the evaluation of $i$ and $j$ derivatives in \eqref{Eq_ldfslp_iso}  
    284 will include a pressure dependent part, leading to the wrong evaluation of  
    285 the neutral slopes. 
    286  
    287 %gm%  
    288 Note: The solution for $s$-coordinate passes trough the use of different  
    289 (and better) expression for the constraint on iso-neutral fluxes. Following  
    290 \citet{Griffies_Bk04}, instead of specifying directly that there is a zero neutral  
    291 diffusive flux of locally referenced potential density, we stay in the $T$-$S$  
    292 plane and consider the balance between the neutral direction diffusive fluxes  
    293 of potential temperature and salinity: 
    294 \begin{equation} 
    295 \alpha \ \textbf{F}(T) = \beta \ \textbf{F}(S) 
    296 \end{equation} 
    297 %gm{  where vector F is ....} 
    298  
    299 This constraint leads to the following definition for the slopes: 
    300  
    301 \begin{equation} \label{Eq_ldfslp_iso2} 
    302 \begin{split} 
    303  r_{1u} &= \frac{e_{3u}}{e_{1u}}\; \frac 
    304       {\alpha_u \;\delta_{i+1/2}[T] - \beta_u \;\delta_{i+1/2}[S]} 
    305       {\alpha_u \;\overline{\overline{\delta_{k+1/2}[T]}}^{\,i+1/2,\,k} 
    306        -\beta_u  \;\overline{\overline{\delta_{k+1/2}[S]}}^{\,i+1/2,\,k} } 
    307 \\ 
    308  r_{2v} &= \frac{e_{3v}}{e_{2v}}\; \frac 
    309       {\alpha_v \;\delta_{j+1/2}[T] - \beta_v \;\delta_{j+1/2}[S]} 
    310       {\alpha_v \;\overline{\overline{\delta_{k+1/2}[T]}}^{\,j+1/2,\,k} 
    311        -\beta_v  \;\overline{\overline{\delta_{k+1/2}[S]}}^{\,j+1/2,\,k} } 
    312 \\ 
    313  r_{1w} &= \frac{e_{3w}}{e_{1w}}\; \frac 
    314       {\alpha_w \;\overline{\overline{\delta_{i+1/2}[T]}}^{\,i,\,k+1/2} 
    315        -\beta_w  \;\overline{\overline{\delta_{i+1/2}[S]}}^{\,i,\,k+1/2} } 
    316       {\alpha_w \;\delta_{k+1/2}[T] - \beta_w \;\delta_{k+1/2}[S]} 
    317 \\ 
    318  r_{2w} &= \frac{e_{3w}}{e_{2w}}\; \frac 
    319       {\alpha_w \;\overline{\overline{\delta_{j+1/2}[T]}}^{\,j,\,k+1/2} 
    320        -\beta_w  \;\overline{\overline{\delta_{j+1/2}[S]}}^{\,j,\,k+1/2} } 
    321       {\alpha_w \;\delta_{k+1/2}[T] - \beta_w \;\delta_{k+1/2}[S]} 
    322 \\ 
    323 \end{split} 
    324 \end{equation} 
    325 where $\alpha$ and $\beta$, the thermal expansion and saline contraction  
    326 coefficients introduced in \S\ref{TRA_bn2}, have to be evaluated at the three  
    327 velocity points. In order to save computation time, they should be approximated  
    328 by the mean of their values at $T$-points (for example in the case of $\alpha$:   
    329 $\alpha_u=\overline{\alpha_T}^{i+1/2}$,  $\alpha_v=\overline{\alpha_T}^{j+1/2}$  
    330 and $\alpha_w=\overline{\alpha_T}^{k+1/2}$). 
    331  
    332 Note that such a formulation could be also used in the $z$-coordinate and  
    333 $z$-coordinate with partial steps cases. 
    334  
    335 \end{description} 
    336  
    337 This implementation is a rather old one. It is similar to the one 
    338 proposed by Cox [1987], except for the background horizontal 
    339 diffusion. Indeed, the Cox implementation of isopycnal diffusion in 
    340 GFDL-type models requires a minimum background horizontal diffusion 
    341 for numerical stability reasons.  To overcome this problem, several 
    342 techniques have been proposed in which the numerical schemes of the 
    343 ocean model are modified \citep{Weaver_Eby_JPO97, 
    344   Griffies_al_JPO98}. Griffies's scheme is now available in \NEMO if 
    345 \np{traldf\_grif\_iso} is set true; see Appdx \ref{sec:triad}. Here, 
    346 another strategy is presented \citep{Lazar_PhD97}: a local 
    347 filtering of the iso-neutral slopes (made on 9 grid-points) prevents 
    348 the development of grid point noise generated by the iso-neutral 
    349 diffusion operator (Fig.~\ref{Fig_LDF_ZDF1}). This allows an 
    350 iso-neutral diffusion scheme without additional background horizontal 
    351 mixing. This technique can be viewed as a diffusion operator that acts 
    352 along large-scale (2~$\Delta$x) \gmcomment{2deltax doesnt seem very 
    353   large scale} iso-neutral surfaces. The diapycnal diffusion required 
    354 for numerical stability is thus minimized and its net effect on the 
    355 flow is quite small when compared to the effect of an horizontal 
    356 background mixing. 
    357  
    358 Nevertheless, this iso-neutral operator does not ensure that variance cannot increase,  
    359 contrary to the \citet{Griffies_al_JPO98} operator which has that property.  
    360  
    361 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    362 \begin{figure}[!ht]      \begin{center} 
    363 \includegraphics[width=0.70\textwidth]{./TexFiles/Figures/Fig_LDF_ZDF1.pdf} 
    364 \caption {    \label{Fig_LDF_ZDF1} 
    365 averaging procedure for isopycnal slope computation.} 
    366 \end{center}    \end{figure} 
    367 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    368  
    369 %There are three additional questions about the slope calculation.  
    370 %First the expression for the rotation tensor has been obtain assuming the "small slope" approximation, so a bound has to be imposed on slopes.  
    371 %Second, numerical stability issues also require a bound on slopes.  
    372 %Third, the question of boundary condition specified on slopes... 
    373  
    374 %from griffies: chapter 13.1.... 
    375  
    376  
    377  
    378 % In addition and also for numerical stability reasons \citep{Cox1987, Griffies_Bk04},  
    379 % the slopes are bounded by $1/100$ everywhere. This limit is decreasing linearly  
    380 % to zero fom $70$ meters depth and the surface (the fact that the eddies "feel" the  
    381 % surface motivates this flattening of isopycnals near the surface). 
    382  
    383 For numerical stability reasons \citep{Cox1987, Griffies_Bk04}, the slopes must also  
    384 be bounded by $1/100$ everywhere. This constraint is applied in a piecewise linear  
    385 fashion, increasing from zero at the surface to $1/100$ at $70$ metres and thereafter  
    386 decreasing to zero at the bottom of the ocean. (the fact that the eddies "feel" the  
    387 surface motivates this flattening of isopycnals near the surface). 
    388  
    389 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    390 \begin{figure}[!ht]     \begin{center} 
    391 \includegraphics[width=0.70\textwidth]{./TexFiles/Figures/Fig_eiv_slp.pdf} 
    392 \caption {     \label{Fig_eiv_slp} 
    393 Vertical profile of the slope used for lateral mixing in the mixed layer :  
    394 \textit{(a)} in the real ocean the slope is the iso-neutral slope in the ocean interior,  
    395 which has to be adjusted at the surface boundary (i.e. it must tend to zero at the  
    396 surface since there is no mixing across the air-sea interface: wall boundary  
    397 condition). Nevertheless, the profile between the surface zero value and the interior  
    398 iso-neutral one is unknown, and especially the value at the base of the mixed layer ;  
    399 \textit{(b)} profile of slope using a linear tapering of the slope near the surface and  
    400 imposing a maximum slope of 1/100 ; \textit{(c)} profile of slope actually used in  
    401 \NEMO: a linear decrease of the slope from zero at the surface to its ocean interior  
    402 value computed just below the mixed layer. Note the huge change in the slope at the  
    403 base of the mixed layer between  \textit{(b)}  and \textit{(c)}.} 
    404 \end{center}   \end{figure} 
    405 %>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    406  
    407 \colorbox{yellow}{add here a discussion about the flattening of the slopes, vs  tapering the coefficient.} 
    408  
    409 \subsection{slopes for momentum iso-neutral mixing} 
    410  
    411 The iso-neutral diffusion operator on momentum is the same as the one used on  
    412 tracers but applied to each component of the velocity separately (see  
    413 \eqref{Eq_dyn_ldf_iso} in section~\ref{DYN_ldf_iso}). The slopes between the  
    414 surface along which the diffusion operator acts and the surface of computation  
    415 ($z$- or $s$-surfaces) are defined at $T$-, $f$-, and \textit{uw}- points for the  
    416 $u$-component, and $T$-, $f$- and \textit{vw}- points for the $v$-component.  
    417 They are computed from the slopes used for tracer diffusion, $i.e.$  
    418 \eqref{Eq_ldfslp_geo} and \eqref{Eq_ldfslp_iso} : 
    419  
    420 \begin{equation} \label{Eq_ldfslp_dyn} 
    421 \begin{aligned} 
    422 &r_{1t}\ \ = \overline{r_{1u}}^{\,i}       &&&    r_{1f}\ \ &= \overline{r_{1u}}^{\,i+1/2} \\ 
    423 &r_{2f} \ \ = \overline{r_{2v}}^{\,j+1/2} &&&   r_{2t}\ &= \overline{r_{2v}}^{\,j} \\ 
    424 &r_{1uw}  = \overline{r_{1w}}^{\,i+1/2} &&\ \ \text{and} \ \ &   r_{1vw}&= \overline{r_{1w}}^{\,j+1/2} \\ 
    425 &r_{2uw}= \overline{r_{2w}}^{\,j+1/2} &&&         r_{2vw}&= \overline{r_{2w}}^{\,j+1/2}\\ 
    426 \end{aligned} 
    427 \end{equation} 
    428  
    429 The major issue remaining is in the specification of the boundary conditions.  
    430 The same boundary conditions are chosen as those used for lateral  
    431 diffusion along model level surfaces, i.e. using the shear computed along  
    432 the model levels and with no additional friction at the ocean bottom (see  
    433 {\S\ref{LBC_coast}). 
    434  
    435  
    436 % ================================================================ 
    437444% Eddy Induced Mixing 
    438445% ================================================================ 
Note: See TracChangeset for help on using the changeset viewer.