/[MITgcm]/manual/s_phys_pkgs/text/gmredi.tex
ViewVC logotype

Diff of /manual/s_phys_pkgs/text/gmredi.tex

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.1 by adcroft, Thu Sep 27 19:43:36 2001 UTC revision 1.4 by adcroft, Tue Nov 13 14:54:50 2001 UTC
# Line 167  In the instance that $\kappa_{GM} = \kap Line 167  In the instance that $\kappa_{GM} = \kap
167  \end{array}  \end{array}
168  \right)  \right)
169  \end{equation}  \end{equation}
170  which differs from the variable laplacian diffusion tensor by only  which differs from the variable Laplacian diffusion tensor by only
171  two non-zero elements in the $z$-row.  two non-zero elements in the $z$-row.
172    
173    \fbox{ \begin{minipage}{4.75in}
174    {\em S/R GMREDI\_CALC\_TENSOR} ({\em pkg/gmredi/gmredi\_calc\_tensor.F})
175    
176    $\sigma_x$: {\bf SlopeX} (argument on entry)
177    
178    $\sigma_y$: {\bf SlopeY} (argument on entry)
179    
180    $\sigma_z$: {\bf SlopeY} (argument)
181    
182    $S_x$: {\bf SlopeX} (argument on exit)
183    
184    $S_y$: {\bf SlopeY} (argument on exit)
185    
186    \end{minipage} }
187    
188    
189    
190  \subsection{Variable $\kappa_{GM}$}  \subsection{Variable $\kappa_{GM}$}
191    
192  Visbeck et al., 1996, suggest making the eddy coefficient,  Visbeck et al., 1996, suggest making the eddy coefficient,
# Line 201  Substituting into the formula for $\kapp Line 218  Substituting into the formula for $\kapp
218  Experience with the GFDL model showed that the GM scheme has to be  Experience with the GFDL model showed that the GM scheme has to be
219  matched to the convective parameterization. This was originally  matched to the convective parameterization. This was originally
220  expressed in connection with the introduction of the KPP boundary  expressed in connection with the introduction of the KPP boundary
221  layer scheme (Large et al., 97) but infact, as subsequent experience  layer scheme (Large et al., 97) but in fact, as subsequent experience
222  with the MIT model has found, is necessary for any convective  with the MIT model has found, is necessary for any convective
223  parameterization.  parameterization.
224    
# Line 219  $z_\sigma^{*}$: {\bf dRdSigmaLtd} (argum Line 236  $z_\sigma^{*}$: {\bf dRdSigmaLtd} (argum
236    
237  \end{minipage} }  \end{minipage} }
238    
239    \begin{figure}
240    \begin{center}
241    \resizebox{5.0in}{3.0in}{\includegraphics{part6/tapers.eps}}
242    \end{center}
243    \caption{Taper functions used in GKW91 and DM95.}
244    \label{fig:tapers}
245    \end{figure}
246    
247    \begin{figure}
248    \begin{center}
249    \resizebox{5.0in}{3.0in}{\includegraphics{part6/effective_slopes.eps}}
250    \end{center}
251    \caption{Effective slope as a function of ``true'' slope using Cox
252    slope clipping, GKW91 limiting and DM95 limiting.}
253    \label{fig:effective_slopes}
254    \end{figure}
255    
256    
257  \subsubsection{Slope clipping}  \subsubsection{Slope clipping}
258    
# Line 227  homogenized, unstable or nearly unstable Line 261  homogenized, unstable or nearly unstable
261  such regions can be either infinite, very large with a sign reversal  such regions can be either infinite, very large with a sign reversal
262  or simply very large. From a numerical point of view, large slopes  or simply very large. From a numerical point of view, large slopes
263  lead to large variations in the tensor elements (implying large bolus  lead to large variations in the tensor elements (implying large bolus
264  flow) and can be numerically unstable. This was first reognized by  flow) and can be numerically unstable. This was first recognized by
265  Cox, 1987, who implemented ``slope clipping'' in the isopycnal mixing  Cox, 1987, who implemented ``slope clipping'' in the isopycnal mixing
266  tensor. Here, the slope magnitude is simply restricted by an upper  tensor. Here, the slope magnitude is simply restricted by an upper
267  limit:  limit:
# Line 262  a) using the GM scheme with clipping and Line 296  a) using the GM scheme with clipping and
296  diffusion). The classic result of dramatically reduced mixed layers is  diffusion). The classic result of dramatically reduced mixed layers is
297  evident. Indeed, the deep convection sites to just one or two points  evident. Indeed, the deep convection sites to just one or two points
298  each and are much shallower than we might prefer. This, it turns out,  each and are much shallower than we might prefer. This, it turns out,
299  is due to the over zealous restratification due to the bolus transport  is due to the over zealous re-stratification due to the bolus transport
300  parameterization. Limiting the slopes also breaks the adiabatic nature  parameterization. Limiting the slopes also breaks the adiabatic nature
301  of the GM/Redi parameterization, re-introducing diabatic fluxes in  of the GM/Redi parameterization, re-introducing diabatic fluxes in
302  regions where the limiting is in effect.  regions where the limiting is in effect.
# Line 326  GM\_tap\-er\_scheme = 'ldd97'} in {\em d Line 360  GM\_tap\-er\_scheme = 'ldd97'} in {\em d
360    
361    
362    
363    
364  \begin{figure}  \begin{figure}
365    \begin{center}
366  %\includegraphics{mixedlayer-cox.eps}  %\includegraphics{mixedlayer-cox.eps}
367  %\includegraphics{mixedlayer-diff.eps}  %\includegraphics{mixedlayer-diff.eps}
368    Figure missing.
369    \end{center}
370  \caption{Mixed layer depth using GM parameterization with a) Cox slope  \caption{Mixed layer depth using GM parameterization with a) Cox slope
371  clipping and for comparison b) using horizontal constant diffusion.}  clipping and for comparison b) using horizontal constant diffusion.}
372  \ref{fig-mixedlayer}  \label{fig-mixedlayer}
 \end{figure}  
   
 \begin{figure}  
 %\includegraphics{slopelimits.eps}  
 \caption{Effective slope as a function of ``true'' slope using a) Cox  
 slope clipping, b) GKW91 limiting, c) DM95 limiting and d) LDD97  
 limiting.}  
373  \end{figure}  \end{figure}
374    
375    
 %\begin{figure}  
 %\includegraphics{coxslope.eps}  
 %\includegraphics{gkw91slope.eps}  
 %\includegraphics{dm95slope.eps}  
 %\includegraphics{ldd97slope.eps}  
 %\caption{Effective slope magnitude at 100~m depth evaluated using a)  
 %Cox slope clipping, b) GKW91 limiting, c) DM95 limiting and d) LDD97  
 %limiting.}  
 %\end{figure}  
   
 \section{Discretisation and code}  
   
 This is the old documentation.....has to be brought upto date with MITgcm.  
   
   
 The Gent-McWilliams-Redi parameterization is implemented through the  
 package ``gmredi''. There are two necessary calls to ``gmredi''  
 routines other than initialization; 1) to calculate the slope tensor  
 as a function of the current model state ({\bf gmredi\_calc\_tensor})  
 and 2) evaluation of the lateral and vertical fluxes due to gradients  
 along isopycnals or bolus transport ({\bf gmredi\_xtransport}, {\bf  
 gmredi\_ytransport} and {\bf gmredi-rtransport}).  
   
 Each element of the tensor is discretised to be adiabatic and so that  
 there would be no flux if the gmredi operator is applied to buoyancy.  
 To acheive this we have to consider both these constraints for each  
 row of the tensor, each row corresponding to a 'u', 'v' or 'w' point  
 on the model grid.  
   
 The code that implements the Redi/GM/Griffies schemes involves an  
 original core routine {\bf inc\_tracer()} that is used to calculate  
 the tendency in the tracers (namely, salt and potential temperature)  
 and a new routine {\bf RediTensor()} that calculates the tensor  
 components and $\kappa_{GM}$.  
   
 \subsection{subroutine RediTensor()}  
   
 {\small  
 \begin{verbatim}  
 subroutine RediTensor(Temp,Salt,Kredigm,K31,K32,K33, nIter,DumpFlag)  
                       |---in--| |-------out-------|  
 ! Input  
 real Temp(Nx,Ny,Nz)     ! Potential temperature  
 real Salt(Nx,Ny,Nz)     ! Salinity  
 ! Output  
 real Kredigm(Nx,Ny,Nz)  ! Redi/GM eddy coefficient  
 real K31(Nx,Ny,Nz)      ! Redi/GM (3,1) tensor component  
 real K32(Nx,Ny,Nz)      ! Redi/GM (3,2) tensor component  
 real K33(Nx,Ny,Nz)      ! Redi/GM (3,3) tensor component  
 ! Auxiliary input  
 integer nIter           ! interation/time-step number  
 logical DumpFlag        ! flag to indicate routine should ``dump''  
 \end{verbatim}  
 }  
   
 The subroutine {\bf RediTensor()} is called from {\bf model()} with  
 input arguments $T$ and $S$. It returns the 3D-arrays {\tt Kredigm},  
 {\t K31}, {\tt K32} and {\tt K33} which represent $\kappa_{GM}$ (at  
 $T/S$ points) and the three components of the bottom row in the  
 Redi/GM tensor; $2 S_x$, $2 S_y$ and $|S|^2$ respectively, all at $W$  
 points.  
   
 The discretisations and algorithm within {\bf RediTensor()} are as  
 follows.  The routine first calculates the locally reference potential  
 density $\sigma_\theta$ from $T$ and $S$ and calculates the potential  
 density gradients in subroutine {\bf gradSigma()}:  
   
 \centerline{\begin{tabular}{ccl}  
 & & \\  
 Array & Grid-point & Definition \\  
 {\tt SigX} & U &  
 $\sigma_x = \frac{1}{\Delta x} \delta_x \sigma|_{z(k)}$  
 \\  
 {\tt SigY} & V &  
 $\sigma_y = \frac{1}{\Delta y} \delta_y \sigma|_{z(k)}$  
 \\  
 {\tt SigZ} & W &  
 $\sigma_z = \frac{1}{\Delta z}  
 [ \sigma|_{z(k)}(k-1/2) - \sigma|_{z(k)}(k+1/2) ]$  
 \\  
 \\  
 \end{tabular}}  
   
 Note that $\sigma_z$ is the static stability because the potential  
 densities are referenced to the same reference level ($W$-level).  
   
 The next step calculates the three tensor components {\tt K13}, {\tt  
 K23} and {\tt K33} in subroutine {\bf KtensorWface()}. First, the  
 lateral gradients $\sigma_x$ and $\sigma_y$ are interpolated to the  
 $W$ points and stored in intermediate variables:  
 \begin{eqnarray*}  
 \mbox{\tt Sx} & = & \overline{ \overline{ \sigma_x }^x }^z \\  
 \mbox{\tt Sy} & = & \overline{ \overline{ \sigma_y }^y }^z  
 \end{eqnarray*}  
 Next, the magnitude of ${\bf \nabla}_z \sigma$ is stored in an intermediate  
 variable:  
 \begin{displaymath}  
 \mbox{\tt Sxy2} = \sqrt{ {\tt Sx}^2 + {\tt Sy}^2 }  
 \end{displaymath}  
 The stratification ($\sigma_z$) is ``checked'' such that the slope  
 vector has magnitude less than or equal to {\tt Smax} and stored in  
 an intermediate variable:  
 \begin{displaymath}  
 \mbox{\tt Sz} = \max ( \sigma_z , - \mbox{\tt Sxy2/Smax} )  
 \end{displaymath}  
 This guarantees stability and at the same time retains the lateral  
 orientation of the slope vector. The tensor components are then calculated:  
 \begin{eqnarray*}  
 \mbox{\tt K13} & = & -2 {\tt Sx/Sz} \\  
 \mbox{\tt K23} & = & -2 {\tt Sx/Sz} \\  
 \mbox{\tt K33} & = & ({\tt Sx/Sz})^2 + ({\tt Sy/Sz})^2  
 \end{eqnarray*}  
   
 Finally, {\tt Kredigm} ($\kappa_{GM}$) is calculated in subroutine  
 {\bf GMRediCoefficient()}. First, all the gradients are interpolated  
 to the $T/S$ points and stored in intermediate variables:  
 \begin{eqnarray*}  
 \mbox{\tt Sx} & = & \overline{ \sigma_x }^x \\  
 \mbox{\tt Sy} & = & \overline{ \sigma_y }^y \\  
 \mbox{\tt Sz} & = & \overline{ \sigma_z }^z  
 \end{eqnarray*}  
 Again, a nominal stratification is found by ``check'' the magnitude of  
 the slope vector but here is converted to a Brunt-Vasala frequency:  
 \begin{eqnarray*}  
 {\tt M2} & = & \sqrt{ {\tt Sx}^2 + {\tt Sy}^2} \\  
 {\tt N2} & = & - \frac{g}{\rho_o} \max ( {\tt Sz} , -{\tt M2 / Smax}  
 \end{eqnarray*}  
 The magnitude of the slope is then $|S| = {\tt M2}/{\tt N2}$. The Eady  
 growth rate is defined as $|f|/\sqrt(Ri) = |S| N$ and is calculated  
 as:  
 \begin{displaymath}  
 {\tt FrRi} = \frac{\tt M2}{\tt N2} ( - \frac{g}{\rho} {\tt Sz} )  
 \end{displaymath}  
 The Eady growth rate is then averaged over the upper layers (about  
 1100m) and $\kappa_{GM}$ specified from this 2D-variable:  
 \begin{displaymath}  
 {\tt Kredigm} = 0.02 * (200d3 **2) * {\tt FrRi}  
 \end{displaymath}  
   
 \subsection{subroutine inc\_tracer()}  
   
 {\bf inc\-tracer()} is called from {\bf model()} and has {\em four  
 new} arguments:  
 \begin{verbatim}  
 subroutine inc_tracer( ...,Kredigm,K31,K32,K33, ... )  
 real Kredigm(Nx,Ny,Nz)  ! Eddy coefficient  
 real K31(Nx,Ny,Nz)      ! (3,1) tensor coefficient  
 real K32(Nx,Ny,Nz)      ! (3,2) tensor coefficient  
 real K33(Nx,Ny,Nz)      ! (3,3) tensor coefficient  
 \end{verbatim}  
   
 Within the routine, the lateral fluxes, {\tt fluxWest} and {\tt  
 fluxSouth}, in the Redi/GM/Griffies scheme are very similar to the  
 conventional horizontal diffusion terms except that the diffusion  
 coefficient is a function of space and must be interpolated from the  
 $T/S$ points:  
 \begin{eqnarray*}  
 {\tt fluxWest}(\tau) & = & \ldots +  
 \overline{\tt Kredigm}^x \partial_x \tau \\  
 {\tt fluxSouth}(\tau) & = & \ldots +  
 \overline{\tt Kredigm}^y \partial_y \tau  
 \end{eqnarray*}  
   
 The Redi/GM/Griffies scheme adds three terms to the vertical flux  
 ({\tt fluxUpper}) in the tracer equation. It is discretise simply:  
 \begin{displaymath}  
 {\tt fluxUpper}(\tau) = \ldots + \overline{\tt Kredigm}^z  
 \left(  
 {\tt K13} \overline{\partial_x \tau}^{xz} +  
 {\tt K23} \overline{\partial_y \tau}^{yz} +  
 {\tt K33} \partial_z \tau  
 \right)  
 \end{displaymath}  
 On boundaries, {\tt fluxUpper} is set to zero.  
376    
377    

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.4

  ViewVC Help
Powered by ViewVC 1.1.22