Changeset 6391

Timestamp:

2016-03-15T16:14:04+01:00 (8 years ago)

Author:

acc

Message:

Branch 2016/dev_r6325_SIMPLIF_1. Documentation changes for the Smagorinsky reimplementation (#1689. 2016WP-SIMPLIF-5). Note this documentation has been purposely submitted to the 2016WP-SIMPLIF-1 development branch at gm's request to avoid complications or loss when Chap_LDF.tex is rewritten on that branch. This submission affects Chapters/Chap_LDF.tex and Namelist/namdyn_ldf only.

Location:

branches/2016/dev_r6325_SIMPLIF_1/DOC/TexFiles

Files:

• Chapters/Chap_LDF.tex (modified) (2 diffs)
• Namelist/namdyn_ldf (modified) (1 diff)

branches/2016/dev_r6325_SIMPLIF_1/DOC/TexFiles/Chapters/Chap_LDF.tex

@@ -404 +404 @@
 and both \key{traldf\_eiv} and \key{traldf\_c2d} are defined.
 
-\subsubsection{Smagorinsky viscosity (\key{dynldf\_c3d} and \key{dynldf\_smag})}
-
-The \key{dynldf\_smag} key activates a 3D, time-varying viscosity that depends on the
-resolved motions. Following \citep{Smagorinsky_93} the viscosity coefficient is set
-proportional to a local deformation rate based on the horizontal shear and tension,
-namely:
-
-\begin{equation}
-A_{m_{Smag}} = \left(\frac{{\sf CM_{Smag}}}{\pi}\right)^2L^2\vert{D}\vert
-\end{equation}
-
-\noindent where the deformation rate $\vert{D}\vert$ is given by 
-
-\begin{equation}
-\vert{D}\vert=\sqrt{\left({\frac{\partial{u}} {\partial{x}}}
-                         -{\frac{\partial{v}} {\partial{y}}}\right)^2
-                 +  \left({\frac{\partial{u}} {\partial{y}}}
-                         +{\frac{\partial{v}} {\partial{x}}}\right)^2} 
-\end{equation}
-
-\noindent and $L$ is the local gridscale given by:
-
-\begin{equation}
-L^2 = \frac{2{e_1}^2 {e_2}^2}{\left ( {e_1}^2 + {e_2}^2 \right )}
-\end{equation}
+\subsubsection{Smagorinsky viscosity (\np{nn\_ahm\_ijk\_t}=32)}
+
+Setting \np{nn\_ahm\_ijk\_t}=32 activates a 3D, time-varying viscosity that depends on the
+resolved motions. Three control parameters are described in the following sections which
+may be set in \textit{\ngn{namdyn\_ldf}}. Following \citep{Smagorinsky_93} and
+\citep{Griffies_Hallberg_MWR00} the viscosity coefficient is set proportional to a local
+deformation rate based on the horizontal shear and tension, namely:
+
+\begin{equation}
+A_{hm} = \left(\frac{{\sf C_{Smag}}}{\pi}\right)^2L^2\vert{D}\vert
+\end{equation}
+where $D$ is the deformation rate and $L$ is the local gridscale given by:
+\begin{equation}
+L = \frac{2{e_1} {e_2}}{\left ( {e_1} + {e_2} \right )}
+\end{equation}
+In orthogonal curvilinear coordinate, we have the tension and shearing strains ($D_T$ and $D_S$) given by
+\begin{equation}       
+\begin{split}   \label{EVP_strain}
+D_T &= \frac{e_2}{e_1} \,\frac{\partial}{\partial i} \left( \frac{u} {e_2} \right) 
+     - \frac{e_1}{e_2} \,\frac{\partial}{\partial j} \left( \frac{v} {e_1} \right) \\
+D_S &= \frac{e_1}{e_2} \,\frac{\partial}{\partial j} \left( \frac{u} {e_1} \right) 
+     + \frac{e_2}{e_1} \,\frac{\partial}{\partial i} \left( \frac{v} {e_2} \right) \\
+\end{split}   \end{equation}
+and the deformation rate, $D$, is given by:
+\begin{equation}       
+D = \sqrt{  {D_T}^2 + {D_S}^2 }
+\end{equation}
+On an Arakawa C-grid, $D_T$ is naturally defined at $t$-points and $D_S$ at $f$-points.
+Their discretizations are:
+\begin{equation} \label{VP_strain_discrete} \begin{split}
+%
+{D_T} &\equiv \frac{1}{e_{1t} e_{2t}} \left(  e_{2t}^2 \,\delta_i 
+                                          \left[ \frac{u}{e_{2u}} \right]  
+                    - e_{1t}^2 \,\delta_j \left[ \frac{v}{e_{1v}} \right]   \right)  \\
+%
+{D_S} &\equiv \frac{1}{e_{1f} e_{2f}} \left(  e_{1f}^2 \,\delta_{j+1/2} 
+                                         \left[ \frac{u}{e_{1u}} \right]  
+                    - e_{2f}^2 \,\delta_{i+1/2} \left[ \frac{v}{e_{2v}} \right] \right)  
+\end{split} \end{equation}
+The deformation rate, $D$ needs to be defined at both $t$-points and $f$-points:
+\begin{equation} \begin{split}
+D_t &= \sqrt{  {D_T}^2  +   \overline{ \overline{ {D_S}^2 }}^{\,i,\,j} \,  }   \\
+D_f &= \sqrt{               \overline{ \overline{ {D_T}^2 }}^{\,i,\,j} + {D_S}^2   \, }
+\end{split} \end{equation}
 
 \citep{Griffies_Hallberg_MWR00} suggest values in the range 2.2 to 4.0 of the coefficient
-$\sf CM_{Smag}$ for oceanic flows. This value is set via the \np{rn\_cmsmag\_1} namelist
-parameter. An additional parameter: \np{rn\_cmsh} is included in NEMO for experimenting
-with the contribution of the shear term. A value of 1.0 (the default) calculates the
-deformation rate as above; a value of 0.0 will discard the shear term entirely.
-
-For numerical stability, the calculated viscosity is bounded according to the following:
-
-\begin{equation}
-{\rm MIN}\left ({ L^2\over {8\Delta{t}}}, rn\_ahm\_m\_lap\right ) \geq A_{m_{Smag}} 
-                                                                  \geq rn\_ahm\_0\_lap
-\end{equation}
-
-\noindent with both parameters for the upper and lower bounds being provided via the
-indicated namelist parameters.
-
-\bigskip When $ln\_dynldf\_bilap = .true.$, a biharmonic version of the Smagorinsky
-viscosity is also available which sets a coefficient for the biharmonic viscosity as:
-
-\begin{equation}
-B_{m_{Smag}} = - \left(\frac{{\sf CM_{bSmag}}}{\pi}\right)^2 {L^4\over 8}\vert{D}\vert
-\end{equation}
-
-\noindent which is bounded according to:
-
-\begin{equation}
-{\rm MAX}\left (-{ L^4\over {64\Delta{t}}}, rn\_ahm\_m\_blp\right ) \leq B_{m_{Smag}} 
-                                                                    \leq rn\_ahm\_0\_blp
-\end{equation}
-
-\noindent Note the reversal of the inequalities here because NEMO requires the biharmonic
-coefficients as negative numbers. $\sf CM_{bSmag}$ is set via the \np{rn\_cmsmag\_2}
-namelist parameter and the bounding values have corresponding entries in the namelist too.
-
-\bigskip The current implementation in NEMO also allows for 3D, time-varying diffusivities
-to be set using the Smagorinsky approach. Users should note that this option is not
-recommended for many applications since diffusivities will tend to be largest near
-boundaries (where shears are greatest) leading to spurious upwellings
-(\citep{Griffies_Bk04}, chapter 18.3.4). Nevertheless the option is there for those
-wishing to experiment. This choice requires both \key{traldf\_c3d} and \key{traldf\_smag}
-and uses the \np{rn\_chsmag} (${\sf CH_{Smag}}$), \np{rn\_smsh} and \np{rn\_aht\_m}
-namelist parameters in an analogous way to \np{rn\_cmsmag\_1}, \np{rn\_cmsh} and
-\np{rn\_ahm\_m\_lap} (see above) to set the diffusion coefficient:
-
-\begin{equation}
-A_{h_{Smag}} = \left(\frac{{\sf CH_{Smag}}}{\pi}\right)^2L^2\vert{D}\vert
-\end{equation}
-
- 
-For numerical stability, the calculated diffusivity is bounded according to the following:
-
-\begin{equation}
-{\rm MIN}\left ({ L^2\over {8\Delta{t}}}, rn\_aht\_m\right ) \geq A_{h_{Smag}} 
-                                                             \geq rn\_aht\_0
-\end{equation}
-
+$\sf C_{Smag}$ for oceanic flows. This value is set via the \np{rn\_csmc} namelist
+parameter.  The same reference also derives the numerical stability criteria that suggest
+the calculated viscosity should be bounded according to the following:
+
+\begin{equation}
+{ \vert U \vert L\over 2} \leq A_{hm} \leq {L^2 \over 4 (2\Delta t)}
+\end{equation}
+which differ from \citep{Griffies_Hallberg_MWR00} only in the explicit acknowledgement of
+the time interval as $2\Delta t$ to reflect the LeapFrog scheme used in NEMO.  These
+bounds are imposed on the coefficients calculated at T and F points with the speed
+calculated using the appropriate average of the squares of the surrounding velocities.
+
+Some flexibility in these bounds may be appropriate in certain situations. For example,
+when modelling the equatorial region the presence of a strong equatorial current may
+impose a higher lower bound than is desirable or strictly necessary. Flexability is
+provided via two extra namelist parameters which can be used to adjust either or both
+bounds by constant factors.  These parameters are: \np{rn\_minfac} and \np{rn\_maxfac}.
+Their default values are 1.0 which provide the exact bounds specified above but values
+other than 1.0 will adjust the bounds according to:
+
+\begin{equation}
+\np{rn\_minfac} \times { \vert U \vert L\over 2} \leq A_{hm} \leq 
+                                   {L^2 \over 8 \Delta t} \times \np{rn\_maxfac}
+\end{equation}
+
+\bigskip 
+
+When $ln\_dynldf\_bilap = .true.$, a biharmonic version of the Smagorinsky viscosity is
+also available which sets a coefficient for the biharmonic viscosity ( following
+\citep{Griffies_Hallberg_MWR00}) as:
+
+\begin{equation}\begin{split}
+B_{hm} &= \left(\frac{{\sf C_{Smag}}}{\pi}\right)^2 {L^4\over 8}\vert{D}\vert \\
+             &= A_{hm} {L^2\over 8}
+\end{split}\end{equation}
+where $A_{hm}$ is the laplacian form of the Smagorinsky viscosity.  This coefficient is
+computed from the laplacian values with bounds already applied so the effective bounds on
+$B_{hm}$ would be:
+\begin{equation}
+\np{rn\_minfac} \times { \vert U \vert L^3\over 16} \leq B_{hm} \leq 
+                                   {L^4 \over 32 \Delta t} \times \np{rn\_maxfac}
+\end{equation}
+but this lower limit does not reflect the behaviour of the 4th order advection schemes
+available in NEMO that have a diffusive component with a biharmonic operator whose eddy
+coefficient is equivalent to:
+\begin{equation}
+{ \vert U \vert L^3\over 12} 
+\end{equation}
+(see equation \eqref{Eq_traadv_ubs2} and surrounding discussion). A safer lower bound for
+the biharmonic Smagorinsky eddy coefficient in NEMO is therefore implemented by
+introducing a hard-coded multiplier of $4/3$ to the lower bound such that the bounds on
+the biharmonic coefficient are:
+\begin{equation}
+\np{rn\_minfac} \times { \vert U \vert L^3\over 12} \leq B_{hm} \leq 
+                                   {L^4 \over 32 \Delta t} \times \np{rn\_maxfac}
+\end{equation} 
+Note that the bilaplacian operator is implemented as a re-entrant laplacian
+which means the actual values of the coefficient returned by \np{ldf\_dyn} in the
+biharmonic case are $\sqrt{B_{hm}}$.
 
 $\ $\newline    % force a new ligne
@@ -516 +539 @@
 (5) the eddy coefficient associated with a biharmonic operator must be set to a \emph{negative} value.
 
-(6) it is possible to use both the laplacian and biharmonic operators concurrently.
+(6) it is not possible to use both the laplacian and biharmonic operators concurrently.
 
 (7) it is possible to run without explicit lateral diffusion on momentum (\np{ln\_dynldf\_lap} = 

branches/2016/dev_r6325_SIMPLIF_1/DOC/TexFiles/Namelist/namdyn_ldf

@@ -19 +19 @@
    !                                !  = 30  F(i,j,k)=c2d*c1d
    !                                !  = 31  F(i,j,k)=F(grid spacing and local velocity)
+   !                                !  = 32  F(i,j,k)=F(local gridscale and deformation rate)
+   ! Caution in 20 and 30 cases the coefficient have to be given for a 1 degree grid (~111km)
    rn_ahm_0      =  40000.     !  horizontal laplacian eddy viscosity   [m2/s]
    rn_ahm_b      =      0.     !  background eddy viscosity for ldf_iso [m2/s]
    rn_bhm_0      = 1.e+12      !  horizontal bilaplacian eddy viscosity [m4/s]
-   !
-   ! Caution in 20 and 30 cases the coefficient have to be given for a 1 degree grid (~111km)
+   !                       !  Smagorinsky settings (nn_ahm_ijk_t  = 32) :
+   rn_csmc       = 3.5         !  Smagorinsky constant of proportionality
+   rn_minfac     = 1.0         !  multiplier of theorectical lower limit
+   rn_maxfac     = 1.0         !  multiplier of theorectical upper limit
 /