source: branches/2011/dev_NEMO_MERGE_2011/DOC/TexFiles/Chapters/Chap_LDF.tex @ 3267

Last change on this file since 3267 was 3267, checked in by agn, 9 years ago

dev_NEMO_MERGE_2011:(DOC) Finish griffies description, remove need for pstricks

File size: 25.6 KB
Line 
1
2% ================================================================
3% Chapter Ñ Lateral Ocean Physics (LDF)
4% ================================================================
5\chapter{Lateral Ocean Physics (LDF)}
6\label{LDF}
7\minitoc
8
9
10\newpage
11$\ $\newline    % force a new ligne
12
13
14The lateral physics terms in the momentum and tracer equations have been
15described in \S\ref{PE_zdf} and their discrete formulation in \S\ref{TRA_ldf} 
16and \S\ref{DYN_ldf}). In this section we further discuss each lateral physics option.
17Choosing one lateral physics scheme means for the user defining, (1) the space
18and time variations of the eddy coefficients ; (2) the direction along which the
19lateral diffusive fluxes are evaluated (model level, geopotential or isopycnal
20surfaces); and (3) the type of operator used (harmonic, or biharmonic operators,
21and for tracers only, eddy induced advection on tracers). These three aspects
22of the lateral diffusion are set through namelist parameters and CPP keys
23(see the \textit{nam\_traldf} and \textit{nam\_dynldf} below). Note
24that this chapter describes the default implementation of iso-neutral
25tracer mixing, and Griffies's implementation, which is used if
26\np{traldf\_grif}=true, is described in Appdx\ref{sec:triad}
27
28%-----------------------------------nam_traldf - nam_dynldf--------------------------------------------
29\namdisplay{namtra_ldf} 
30\namdisplay{namdyn_ldf} 
31%--------------------------------------------------------------------------------------------------------------
32
33
34% ================================================================
35% Lateral Mixing Coefficients
36% ================================================================
37\section [Lateral Mixing Coefficient (\textit{ldftra}, \textit{ldfdyn})]
38        {Lateral Mixing Coefficient (\mdl{ldftra}, \mdl{ldfdyn}) }
39\label{LDF_coef}
40
41
42Introducing a space variation in the lateral eddy mixing coefficients changes
43the model core memory requirement, adding up to four extra three-dimensional
44arrays for the geopotential or isopycnal second order operator applied to
45momentum. Six CPP keys control the space variation of eddy coefficients:
46three for momentum and three for tracer. The three choices allow:
47a space variation in the three space directions (\key{traldf\_c3d}\key{dynldf\_c3d}),
48in the horizontal plane (\key{traldf\_c2d}\key{dynldf\_c2d}),
49or in the vertical only (\key{traldf\_c1d}\key{dynldf\_c1d}).
50The default option is a constant value over the whole ocean on both momentum and tracers.
51   
52The number of additional arrays that have to be defined and the gridpoint
53position at which they are defined depend on both the space variation chosen
54and the type of operator used. The resulting eddy viscosity and diffusivity
55coefficients can be a function of more than one variable. Changes in the
56computer code when switching from one option to another have been
57minimized by introducing the eddy coefficients as statement functions
58(include file \hf{ldftra\_substitute} and \hf{ldfdyn\_substitute}). The functions
59are replaced by their actual meaning during the preprocessing step (CPP).
60The specification of the space variation of the coefficient is made in
61\mdl{ldftra} and \mdl{ldfdyn}, or more precisely in include files
62\textit{traldf\_cNd.h90} and \textit{dynldf\_cNd.h90}, with N=1, 2 or 3.
63The user can modify these include files as he/she wishes. The way the
64mixing coefficient are set in the reference version can be briefly described
65as follows:
66
67\subsubsection{Constant Mixing Coefficients (default option)}
68When none of the \textbf{key\_dynldf\_...} and \textbf{key\_traldf\_...} keys are
69defined, a constant value is used over the whole ocean for momentum and
70tracers, which is specified through the \np{rn\_ahm0} and \np{rn\_aht0} namelist
71parameters.
72
73\subsubsection{Vertically varying Mixing Coefficients (\key{traldf\_c1d} and \key{dynldf\_c1d})} 
74The 1D option is only available when using the $z$-coordinate with full step.
75Indeed in all the other types of vertical coordinate, the depth is a 3D function
76of (\textbf{i},\textbf{j},\textbf{k}) and therefore, introducing depth-dependent
77mixing coefficients will require 3D arrays. In the 1D option, a hyperbolic variation
78of the lateral mixing coefficient is introduced in which the surface value is
79\np{rn\_aht0} (\np{rn\_ahm0}), the bottom value is 1/4 of the surface value,
80and the transition takes place around z=300~m with a width of 300~m
81($i.e.$ both the depth and the width of the inflection point are set to 300~m).
82This profile is hard coded in file \hf{traldf\_c1d}, but can be easily modified by users.
83
84\subsubsection{Horizontally Varying Mixing Coefficients (\key{traldf\_c2d} and \key{dynldf\_c2d})}
85By default the horizontal variation of the eddy coefficient depends on the local mesh
86size and the type of operator used:
87\begin{equation} \label{Eq_title}
88  A_l = \left\{     
89   \begin{aligned}
90         & \frac{\max(e_1,e_2)}{e_{max}} A_o^l           & \text{for laplacian operator } \\
91         & \frac{\max(e_1,e_2)^{3}}{e_{max}^{3}} A_o^l          & \text{for bilaplacian operator } 
92   \end{aligned}    \right.
93\end{equation}
94where $e_{max}$ is the maximum of $e_1$ and $e_2$ taken over the whole masked
95ocean domain, and $A_o^l$ is the \np{rn\_ahm0} (momentum) or \np{rn\_aht0} (tracer)
96namelist parameter. This variation is intended to reflect the lesser need for subgrid
97scale eddy mixing where the grid size is smaller in the domain. It was introduced in
98the context of the DYNAMO modelling project \citep{Willebrand_al_PO01}.
99Note that such a grid scale dependance of mixing coefficients significantly increase
100the range of stability of model configurations presenting large changes in grid pacing
101such as global ocean models. Indeed, in such a case, a constant mixing coefficient
102can lead to a blow up of the model due to large coefficient compare to the smallest
103grid size (see \S\ref{STP_forward_imp}), especially when using a bilaplacian operator.
104
105Other formulations can be introduced by the user for a given configuration.
106For example, in the ORCA2 global ocean model (\key{orca\_r2}), the laplacian
107viscosity operator uses \np{rn\_ahm0}~= 4.10$^4$ m$^2$/s poleward of 20$^{\circ}$ 
108north and south and decreases linearly to \np{rn\_aht0}~= 2.10$^3$ m$^2$/s
109at the equator \citep{Madec_al_JPO96, Delecluse_Madec_Bk00}. This modification
110can be found in routine \rou{ldf\_dyn\_c2d\_orca} defined in \mdl{ldfdyn\_c2d}.
111Similar modified horizontal variations can be found with the Antarctic or Arctic
112sub-domain options of ORCA2 and ORCA05 (\key{antarctic} or \key{arctic} 
113defined, see \hf{ldfdyn\_antarctic} and \hf{ldfdyn\_arctic}).
114
115\subsubsection{Space Varying Mixing Coefficients (\key{traldf\_c3d} and \key{dynldf\_c3d})}
116
117The 3D space variation of the mixing coefficient is simply the combination of the
1181D and 2D cases, $i.e.$ a hyperbolic tangent variation with depth associated with
119a grid size dependence of the magnitude of the coefficient.
120
121\subsubsection{Space and Time Varying Mixing Coefficients}
122
123There is no default specification of space and time varying mixing coefficient.
124The only case available is specific to the ORCA2 and ORCA05 global ocean
125configurations (\key{orca\_r2} or \key{orca\_r05}). It provides only a tracer
126mixing coefficient for eddy induced velocity (ORCA2) or both iso-neutral and
127eddy induced velocity (ORCA05) that depends on the local growth rate of
128baroclinic instability. This specification is actually used when an ORCA key
129and both \key{traldf\_eiv} and \key{traldf\_c2d} are defined.
130
131$\ $\newline    % force a new ligne
132
133The following points are relevant when the eddy coefficient varies spatially:
134
135(1) the momentum diffusion operator acting along model level surfaces is
136written in terms of curl and divergent components of the horizontal current
137(see \S\ref{PE_ldf}). Although the eddy coefficient could be set to different values
138in these two terms, this option is not currently available.
139
140(2) with an horizontally varying viscosity, the quadratic integral constraints
141on enstrophy and on the square of the horizontal divergence for operators
142acting along model-surfaces are no longer satisfied
143(Appendix~\ref{Apdx_dynldf_properties}).
144
145(3) for isopycnal diffusion on momentum or tracers, an additional purely
146horizontal background diffusion with uniform coefficient can be added by
147setting a non zero value of \np{rn\_ahmb0} or \np{rn\_ahtb0}, a background horizontal
148eddy viscosity or diffusivity coefficient (namelist parameters whose default
149values are $0$). However, the technique used to compute the isopycnal
150slopes is intended to get rid of such a background diffusion, since it introduces
151spurious diapycnal diffusion (see {\S\ref{LDF_slp}).
152
153(4) when an eddy induced advection term is used (\key{traldf\_eiv}), $A^{eiv}$,
154the eddy induced coefficient has to be defined. Its space variations are controlled
155by the same CPP variable as for the eddy diffusivity coefficient ($i.e.$ 
156\textbf{key\_traldf\_cNd}).
157
158(5) the eddy coefficient associated with a biharmonic operator must be set to a \emph{negative} value.
159
160(6) it is possible to use both the laplacian and biharmonic operators concurrently.
161
162(7) it is possible to run without explicit lateral diffusion on momentum (\np{ln\_dynldf\_lap} =
163\np{ln\_dynldf\_bilap} = false). This is recommended when using the UBS advection
164scheme on momentum (\np{ln\_dynadv\_ubs} = true, see \ref{DYN_adv_ubs})
165and can be useful for testing purposes.
166
167% ================================================================
168% Direction of lateral Mixing
169% ================================================================
170\section  [Direction of Lateral Mixing (\textit{ldfslp})]
171      {Direction of Lateral Mixing (\mdl{ldfslp})}
172\label{LDF_slp}
173
174%%%
175\gmcomment{  we should emphasize here that the implementation is a rather old one.
176Better work can be achieved by using \citet{Griffies_al_JPO98, Griffies_Bk04} iso-neutral scheme. }
177
178A direction for lateral mixing has to be defined when the desired operator does
179not act along the model levels. This occurs when $(a)$ horizontal mixing is
180required on tracer or momentum (\np{ln\_traldf\_hor} or \np{ln\_dynldf\_hor})
181in $s$- or mixed $s$-$z$- coordinates, and $(b)$ isoneutral mixing is required
182whatever the vertical coordinate is. This direction of mixing is defined by its
183slopes in the \textbf{i}- and \textbf{j}-directions at the face of the cell of the
184quantity to be diffused. For a tracer, this leads to the following four slopes :
185$r_{1u}$, $r_{1w}$, $r_{2v}$, $r_{2w}$ (see \eqref{Eq_tra_ldf_iso}), while
186for momentum the slopes are  $r_{1t}$, $r_{1uw}$, $r_{2f}$, $r_{2uw}$ for
187$u$ and  $r_{1f}$, $r_{1vw}$, $r_{2t}$, $r_{2vw}$ for $v$.
188
189%gm% add here afigure of the slope in i-direction
190
191\subsection{slopes for tracer geopotential mixing in the $s$-coordinate}
192
193In $s$-coordinates, geopotential mixing ($i.e.$ horizontal mixing) $r_1$ and
194$r_2$ are the slopes between the geopotential and computational surfaces.
195Their discrete formulation is found by locally solving \eqref{Eq_tra_ldf_iso} 
196when the diffusive fluxes in the three directions are set to zero and $T$ is
197assumed to be horizontally uniform, $i.e.$ a linear function of $z_T$, the
198depth of a $T$-point.
199%gm { Steven : My version is obviously wrong since I'm left with an arbitrary constant which is the local vertical temperature gradient}
200
201\begin{equation} \label{Eq_ldfslp_geo}
202\begin{aligned}
203 r_{1u} &= \frac{e_{3u}}{ \left( e_{1u}\;\overline{\overline{e_{3w}}}^{\,i+1/2,\,k} \right)}
204           \;\delta_{i+1/2}[z_t]
205      &\approx \frac{1}{e_{1u}}\; \delta_{i+1/2}[z_t]
206\\
207 r_{2v} &= \frac{e_{3v}}{\left( e_{2v}\;\overline{\overline{e_{3w}}}^{\,j+1/2,\,k} \right)} 
208           \;\delta_{j+1/2} [z_t]
209      &\approx \frac{1}{e_{2v}}\; \delta_{j+1/2}[z_t]
210\\
211 r_{1w} &= \frac{1}{e_{1w}}\;\overline{\overline{\delta_{i+1/2}[z_t]}}^{\,i,\,k+1/2}
212      &\approx \frac{1}{e_{1w}}\; \delta_{i+1/2}[z_{uw}]
213 \\
214 r_{2w} &= \frac{1}{e_{2w}}\;\overline{\overline{\delta_{j+1/2}[z_t]}}^{\,j,\,k+1/2}
215      &\approx \frac{1}{e_{2w}}\; \delta_{j+1/2}[z_{vw}]
216 \\
217\end{aligned}
218\end{equation}
219
220%gm%  caution I'm not sure the simplification was a good idea!
221
222These slopes are computed once in \rou{ldfslp\_init} when \np{ln\_sco}=True,
223and either \np{ln\_traldf\_hor}=True or \np{ln\_dynldf\_hor}=True.
224
225\subsection{Slopes for tracer iso-neutral mixing}\label{LDF_slp_iso}
226In iso-neutral mixing  $r_1$ and $r_2$ are the slopes between the iso-neutral
227and computational surfaces. Their formulation does not depend on the vertical
228coordinate used. Their discrete formulation is found using the fact that the
229diffusive fluxes of locally referenced potential density ($i.e.$ $in situ$ density)
230vanish. So, substituting $T$ by $\rho$ in \eqref{Eq_tra_ldf_iso} and setting the
231diffusive fluxes in the three directions to zero leads to the following definition for
232the neutral slopes:
233
234\begin{equation} \label{Eq_ldfslp_iso}
235\begin{split}
236 r_{1u} &= \frac{e_{3u}}{e_{1u}}\; \frac{\delta_{i+1/2}[\rho]}
237                        {\overline{\overline{\delta_{k+1/2}[\rho]}}^{\,i+1/2,\,k}}
238\\
239 r_{2v} &= \frac{e_{3v}}{e_{2v}}\; \frac{\delta_{j+1/2}\left[\rho \right]}
240                        {\overline{\overline{\delta_{k+1/2}[\rho]}}^{\,j+1/2,\,k}}
241\\
242 r_{1w} &= \frac{e_{3w}}{e_{1w}}\; 
243         \frac{\overline{\overline{\delta_{i+1/2}[\rho]}}^{\,i,\,k+1/2}}
244             {\delta_{k+1/2}[\rho]}
245\\
246 r_{2w} &= \frac{e_{3w}}{e_{2w}}\; 
247         \frac{\overline{\overline{\delta_{j+1/2}[\rho]}}^{\,j,\,k+1/2}}
248             {\delta_{k+1/2}[\rho]}
249\\
250\end{split}
251\end{equation}
252
253%gm% rewrite this as the explanation is not very clear !!!
254%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.
255
256%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).
257
258%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.
259
260As the mixing is performed along neutral surfaces, the gradient of $\rho$ in
261\eqref{Eq_ldfslp_iso} has to be evaluated at the same local pressure (which,
262in decibars, is approximated by the depth in meters in the model). Therefore
263\eqref{Eq_ldfslp_iso} cannot be used as such, but further transformation is
264needed depending on the vertical coordinate used:
265
266\begin{description}
267
268\item[$z$-coordinate with full step : ] in \eqref{Eq_ldfslp_iso} the densities
269appearing in the $i$ and $j$ derivatives  are taken at the same depth, thus
270the $in situ$ density can be used. This is not the case for the vertical
271derivatives: $\delta_{k+1/2}[\rho]$ is replaced by $-\rho N^2/g$, where $N^2$ 
272is the local Brunt-Vais\"{a}l\"{a} frequency evaluated following
273\citet{McDougall1987} (see \S\ref{TRA_bn2}).
274
275\item[$z$-coordinate with partial step : ] this case is identical to the full step
276case except that at partial step level, the \emph{horizontal} density gradient
277is evaluated as described in \S\ref{TRA_zpshde}.
278
279\item[$s$- or hybrid $s$-$z$- coordinate : ] in the current release of \NEMO,
280iso-neutral mixing is only employed for $s$-coordinates if the
281Griffies scheme is used (\np{traldf\_grif}=true; see Appdx \ref{sec:triad}).
282In other words, iso-neutral mixing will only be accurately represented with a
283linear equation of state (\np{nn\_eos}=1 or 2). In the case of a "true" equation
284of state, the evaluation of $i$ and $j$ derivatives in \eqref{Eq_ldfslp_iso} 
285will include a pressure dependent part, leading to the wrong evaluation of
286the neutral slopes.
287
288%gm%
289Note: The solution for $s$-coordinate passes trough the use of different
290(and better) expression for the constraint on iso-neutral fluxes. Following
291\citet{Griffies_Bk04}, instead of specifying directly that there is a zero neutral
292diffusive flux of locally referenced potential density, we stay in the $T$-$S$ 
293plane and consider the balance between the neutral direction diffusive fluxes
294of potential temperature and salinity:
295\begin{equation}
296\alpha \ \textbf{F}(T) = \beta \ \textbf{F}(S)
297\end{equation}
298%gm{  where vector F is ....}
299
300This constraint leads to the following definition for the slopes:
301
302\begin{equation} \label{Eq_ldfslp_iso2}
303\begin{split}
304 r_{1u} &= \frac{e_{3u}}{e_{1u}}\; \frac
305      {\alpha_u \;\delta_{i+1/2}[T] - \beta_u \;\delta_{i+1/2}[S]}
306      {\alpha_u \;\overline{\overline{\delta_{k+1/2}[T]}}^{\,i+1/2,\,k}
307       -\beta_\;\overline{\overline{\delta_{k+1/2}[S]}}^{\,i+1/2,\,k} }
308\\
309 r_{2v} &= \frac{e_{3v}}{e_{2v}}\; \frac
310      {\alpha_v \;\delta_{j+1/2}[T] - \beta_v \;\delta_{j+1/2}[S]}
311      {\alpha_v \;\overline{\overline{\delta_{k+1/2}[T]}}^{\,j+1/2,\,k}
312       -\beta_\;\overline{\overline{\delta_{k+1/2}[S]}}^{\,j+1/2,\,k} }
313\\
314 r_{1w} &= \frac{e_{3w}}{e_{1w}}\; \frac
315      {\alpha_w \;\overline{\overline{\delta_{i+1/2}[T]}}^{\,i,\,k+1/2}
316       -\beta_\;\overline{\overline{\delta_{i+1/2}[S]}}^{\,i,\,k+1/2} }
317      {\alpha_w \;\delta_{k+1/2}[T] - \beta_w \;\delta_{k+1/2}[S]}
318\\
319 r_{2w} &= \frac{e_{3w}}{e_{2w}}\; \frac
320      {\alpha_w \;\overline{\overline{\delta_{j+1/2}[T]}}^{\,j,\,k+1/2}
321       -\beta_\;\overline{\overline{\delta_{j+1/2}[S]}}^{\,j,\,k+1/2} }
322      {\alpha_w \;\delta_{k+1/2}[T] - \beta_w \;\delta_{k+1/2}[S]}
323\\
324\end{split}
325\end{equation}
326where $\alpha$ and $\beta$, the thermal expansion and saline contraction
327coefficients introduced in \S\ref{TRA_bn2}, have to be evaluated at the three
328velocity points. In order to save computation time, they should be approximated
329by the mean of their values at $T$-points (for example in the case of $\alpha$
330$\alpha_u=\overline{\alpha_T}^{i+1/2}$$\alpha_v=\overline{\alpha_T}^{j+1/2}$ 
331and $\alpha_w=\overline{\alpha_T}^{k+1/2}$).
332
333Note that such a formulation could be also used in the $z$-coordinate and
334$z$-coordinate with partial steps cases.
335
336\end{description}
337
338This implementation is a rather old one. It is similar to the one
339proposed by Cox [1987], except for the background horizontal
340diffusion. Indeed, the Cox implementation of isopycnal diffusion in
341GFDL-type models requires a minimum background horizontal diffusion
342for numerical stability reasons.  To overcome this problem, several
343techniques have been proposed in which the numerical schemes of the
344ocean model are modified \citep{Weaver_Eby_JPO97,
345  Griffies_al_JPO98}. Griffies's scheme is now available in \NEMO if
346\np{traldf\_grif\_iso} is set true; see Appdx \ref{sec:triad}. Here,
347another strategy is presented \citep{Lazar_PhD97}: a local
348filtering of the iso-neutral slopes (made on 9 grid-points) prevents
349the development of grid point noise generated by the iso-neutral
350diffusion operator (Fig.~\ref{Fig_LDF_ZDF1}). This allows an
351iso-neutral diffusion scheme without additional background horizontal
352mixing. This technique can be viewed as a diffusion operator that acts
353along large-scale (2~$\Delta$x) \gmcomment{2deltax doesnt seem very
354  large scale} iso-neutral surfaces. The diapycnal diffusion required
355for numerical stability is thus minimized and its net effect on the
356flow is quite small when compared to the effect of an horizontal
357background mixing.
358
359Nevertheless, this iso-neutral operator does not ensure that variance cannot increase,
360contrary to the \citet{Griffies_al_JPO98} operator which has that property.
361
362%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
363\begin{figure}[!ht]      \begin{center}
364\includegraphics[width=0.70\textwidth]{./TexFiles/Figures/Fig_LDF_ZDF1.pdf}
365\caption {    \label{Fig_LDF_ZDF1}
366averaging procedure for isopycnal slope computation.}
367\end{center}    \end{figure}
368%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
369
370%There are three additional questions about the slope calculation.
371%First the expression for the rotation tensor has been obtain assuming the "small slope" approximation, so a bound has to be imposed on slopes.
372%Second, numerical stability issues also require a bound on slopes.
373%Third, the question of boundary condition specified on slopes...
374
375%from griffies: chapter 13.1....
376
377
378
379% In addition and also for numerical stability reasons \citep{Cox1987, Griffies_Bk04},
380% the slopes are bounded by $1/100$ everywhere. This limit is decreasing linearly
381% to zero fom $70$ meters depth and the surface (the fact that the eddies "feel" the
382% surface motivates this flattening of isopycnals near the surface).
383
384For numerical stability reasons \citep{Cox1987, Griffies_Bk04}, the slopes must also
385be bounded by $1/100$ everywhere. This constraint is applied in a piecewise linear
386fashion, increasing from zero at the surface to $1/100$ at $70$ metres and thereafter
387decreasing to zero at the bottom of the ocean. (the fact that the eddies "feel" the
388surface motivates this flattening of isopycnals near the surface).
389
390%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
391\begin{figure}[!ht]     \begin{center}
392\includegraphics[width=0.70\textwidth]{./TexFiles/Figures/Fig_eiv_slp.pdf}
393\caption {     \label{Fig_eiv_slp}
394Vertical profile of the slope used for lateral mixing in the mixed layer :
395\textit{(a)} in the real ocean the slope is the iso-neutral slope in the ocean interior,
396which has to be adjusted at the surface boundary (i.e. it must tend to zero at the
397surface since there is no mixing across the air-sea interface: wall boundary
398condition). Nevertheless, the profile between the surface zero value and the interior
399iso-neutral one is unknown, and especially the value at the base of the mixed layer ;
400\textit{(b)} profile of slope using a linear tapering of the slope near the surface and
401imposing a maximum slope of 1/100 ; \textit{(c)} profile of slope actually used in
402\NEMO: a linear decrease of the slope from zero at the surface to its ocean interior
403value computed just below the mixed layer. Note the huge change in the slope at the
404base of the mixed layer between  \textit{(b)}  and \textit{(c)}.}
405\end{center}   \end{figure}
406%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
407
408\colorbox{yellow}{add here a discussion about the flattening of the slopes, vs  tapering the coefficient.}
409
410\subsection{slopes for momentum iso-neutral mixing}
411
412The iso-neutral diffusion operator on momentum is the same as the one used on
413tracers but applied to each component of the velocity separately (see
414\eqref{Eq_dyn_ldf_iso} in section~\ref{DYN_ldf_iso}). The slopes between the
415surface along which the diffusion operator acts and the surface of computation
416($z$- or $s$-surfaces) are defined at $T$-, $f$-, and \textit{uw}- points for the
417$u$-component, and $T$-, $f$- and \textit{vw}- points for the $v$-component.
418They are computed from the slopes used for tracer diffusion, $i.e.$ 
419\eqref{Eq_ldfslp_geo} and \eqref{Eq_ldfslp_iso} :
420
421\begin{equation} \label{Eq_ldfslp_dyn}
422\begin{aligned}
423&r_{1t}\ \ = \overline{r_{1u}}^{\,i}       &&&    r_{1f}\ \ &= \overline{r_{1u}}^{\,i+1/2} \\
424&r_{2f} \ \ = \overline{r_{2v}}^{\,j+1/2} &&&   r_{2t}\ &= \overline{r_{2v}}^{\,j} \\
425&r_{1uw}  = \overline{r_{1w}}^{\,i+1/2} &&\ \ \text{and} \ \ &   r_{1vw}&= \overline{r_{1w}}^{\,j+1/2} \\
426&r_{2uw}= \overline{r_{2w}}^{\,j+1/2} &&&         r_{2vw}&= \overline{r_{2w}}^{\,j+1/2}\\
427\end{aligned}
428\end{equation}
429
430The major issue remaining is in the specification of the boundary conditions.
431The same boundary conditions are chosen as those used for lateral
432diffusion along model level surfaces, i.e. using the shear computed along
433the model levels and with no additional friction at the ocean bottom (see
434{\S\ref{LBC_coast}).
435
436
437% ================================================================
438% Eddy Induced Mixing
439% ================================================================
440\section  [Eddy Induced Velocity (\textit{traadv\_eiv}, \textit{ldfeiv})]
441      {Eddy Induced Velocity (\mdl{traadv\_eiv}, \mdl{ldfeiv})}
442\label{LDF_eiv}
443
444When Gent and McWilliams [1990] diffusion is used (\key{traldf\_eiv} defined),
445an eddy induced tracer advection term is added, the formulation of which
446depends on the slopes of iso-neutral surfaces. Contrary to the case of iso-neutral
447mixing, the slopes used here are referenced to the geopotential surfaces, $i.e.$ 
448\eqref{Eq_ldfslp_geo} is used in $z$-coordinates, and the sum \eqref{Eq_ldfslp_geo} 
449+ \eqref{Eq_ldfslp_iso} in $s$-coordinates. The eddy induced velocity is given by:
450\begin{equation} \label{Eq_ldfeiv}
451\begin{split}
452 u^* & = \frac{1}{e_{2u}e_{3u}}\; \delta_k \left[e_{2u} \, A_{uw}^{eiv} \; \overline{r_{1w}}^{\,i+1/2} \right]\\
453v^* & = \frac{1}{e_{1u}e_{3v}}\; \delta_k \left[e_{1v} \, A_{vw}^{eiv} \; \overline{r_{2w}}^{\,j+1/2} \right]\\
454w^* & = \frac{1}{e_{1w}e_{2w}}\; \left\{ \delta_i \left[e_{2u} \, A_{uw}^{eiv} \; \overline{r_{1w}}^{\,i+1/2} \right] + \delta_j \left[e_{1v} \, A_{vw}^{eiv} \; \overline{r_{2w}}^{\,j+1/2} \right] \right\} \\
455\end{split}
456\end{equation}
457where $A^{eiv}$ is the eddy induced velocity coefficient whose value is set
458through \np{rn\_aeiv}, a \textit{nam\_traldf} namelist parameter.
459The three components of the eddy induced velocity are computed and add
460to the eulerian velocity in \mdl{traadv\_eiv}. This has been preferred to a
461separate computation of the advective trends associated with the eiv velocity,
462since it allows us to take advantage of all the advection schemes offered for
463the tracers (see \S\ref{TRA_adv}) and not just the $2^{nd}$ order advection
464scheme as in previous releases of OPA \citep{Madec1998}. This is particularly
465useful for passive tracers where \emph{positivity} of the advection scheme is
466of paramount importance.
467
468At the surface, lateral and bottom boundaries, the eddy induced velocity,
469and thus the advective eddy fluxes of heat and salt, are set to zero.
470
471
472
473
474
Note: See TracBrowser for help on using the repository browser.