1 |
\section{Gent/McWiliams/Redi SGS Eddy parameterization} |
\subsection{GMREDI: Gent-McWilliams/Redi SGS Eddy Parameterization} |
2 |
|
\label{sec:pkg:gmredi} |
3 |
|
\begin{rawhtml} |
4 |
|
<!-- CMIREDIR:gmredi: --> |
5 |
|
\end{rawhtml} |
6 |
|
|
7 |
There are two parts to the Redi/GM parameterization of geostrophic |
There are two parts to the Redi/GM parameterization of geostrophic |
8 |
eddies. The first aims to mix tracer properties along isentropes |
eddies. The first, the Redi scheme \citep{re82}, aims to mix tracer properties |
9 |
(neutral surfaces) by means of a diffusion operator oriented along the |
along isentropes (neutral surfaces) by means of a diffusion operator oriented |
10 |
local isentropic surface (Redi). The second part, adiabatically |
along the local isentropic surface. |
11 |
|
The second part, GM \citep{gen-mcw:90,gen-eta:95}, adiabatically |
12 |
re-arranges tracers through an advective flux where the advecting flow |
re-arranges tracers through an advective flux where the advecting flow |
13 |
is a function of slope of the isentropic surfaces (GM). |
is a function of slope of the isentropic surfaces. |
14 |
|
|
15 |
The first GCM implementation of the Redi scheme was by Cox 1987 in the |
The first GCM implementation of the Redi scheme was by \cite{Cox87} in the |
16 |
GFDL ocean circulation model. The original approach failed to |
GFDL ocean circulation model. The original approach failed to |
17 |
distinguish between isopycnals and surfaces of locally referenced |
distinguish between isopycnals and surfaces of locally referenced |
18 |
potential density (now called neutral surfaces) which are proper |
potential density (now called neutral surfaces) which are proper |
27 |
to the boundary condition of zero value on upper and lower |
to the boundary condition of zero value on upper and lower |
28 |
boundaries. The horizontal bolus velocities are then the vertical |
boundaries. The horizontal bolus velocities are then the vertical |
29 |
derivative of these functions. Here in lies a problem highlighted by |
derivative of these functions. Here in lies a problem highlighted by |
30 |
Griffies et al., 1997: the bolus velocities involve multiple |
\cite{gretal:98}: the bolus velocities involve multiple |
31 |
derivatives on the potential density field, which can consequently |
derivatives on the potential density field, which can consequently |
32 |
give rise to noise. Griffies et al. point out that the GM bolus fluxes |
give rise to noise. Griffies et al. point out that the GM bolus fluxes |
33 |
can be identically written as a skew flux which involves fewer |
can be identically written as a skew flux which involves fewer |
36 |
that the horizontal fluxes are unmodified from the lateral diffusion |
that the horizontal fluxes are unmodified from the lateral diffusion |
37 |
parameterization. |
parameterization. |
38 |
|
|
39 |
\subsection{Redi scheme: Isopycnal diffusion} |
\subsubsection{Redi scheme: Isopycnal diffusion} |
40 |
|
|
41 |
The Redi scheme diffuses tracers along isopycnals and introduces a |
The Redi scheme diffuses tracers along isopycnals and introduces a |
42 |
term in the tendency (rhs) of such a tracer (here $\tau$) of the form: |
term in the tendency (rhs) of such a tracer (here $\tau$) of the form: |
47 |
$\bf{K}_{Redi}$ is a rank 2 tensor that projects the gradient of |
$\bf{K}_{Redi}$ is a rank 2 tensor that projects the gradient of |
48 |
$\tau$ onto the isopycnal surface. The unapproximated projection tensor is: |
$\tau$ onto the isopycnal surface. The unapproximated projection tensor is: |
49 |
\begin{equation} |
\begin{equation} |
50 |
\bf{K}_{Redi} = \left( |
\bf{K}_{Redi} = \frac{1}{1 + |S|^2} \left( |
51 |
\begin{array}{ccc} |
\begin{array}{ccc} |
52 |
1 + S_x& S_x S_y & S_x \\ |
1 + S_y^2& -S_x S_y & S_x \\ |
53 |
S_x S_y & 1 + S_y & S_y \\ |
-S_x S_y & 1 + S_x^2 & S_y \\ |
54 |
S_x & S_y & |S|^2 \\ |
S_x & S_y & |S|^2 \\ |
55 |
\end{array} |
\end{array} |
56 |
\right) |
\right) |
76 |
\end{equation} |
\end{equation} |
77 |
|
|
78 |
|
|
79 |
\subsection{GM parameterization} |
\subsubsection{GM parameterization} |
80 |
|
|
81 |
The GM parameterization aims to parameterise the ``advective'' or |
The GM parameterization aims to represent the ``advective'' or |
82 |
``transport'' effect of geostrophic eddies by means of a ``bolus'' |
``transport'' effect of geostrophic eddies by means of a ``bolus'' |
83 |
velocity, $\bf{u}^*$. The divergence of this advective flux is added |
velocity, $\bf{u}^\star$. The divergence of this advective flux is added |
84 |
to the tracer tendency equation (on the rhs): |
to the tracer tendency equation (on the rhs): |
85 |
\begin{equation} |
\begin{equation} |
86 |
- \bf{\nabla} \cdot \tau \bf{u}^* |
- \bf{\nabla} \cdot \tau \bf{u}^\star |
87 |
\end{equation} |
\end{equation} |
88 |
|
|
89 |
The bolus velocity is defined as: |
The bolus velocity $\bf{u}^\star$ is defined as the rotational of a |
90 |
\begin{eqnarray} |
streamfunction $\bf{F}^\star$=$(F_x^\star,F_y^\star,0)$: |
91 |
u^* & = & - \partial_z F_x \\ |
\begin{equation} |
92 |
v^* & = & - \partial_z F_y \\ |
\bf{u}^\star = \nabla \times \bf{F}^\star = |
93 |
w^* & = & \partial_x F_x + \partial_y F_y |
\left( \begin{array}{c} |
94 |
\end{eqnarray} |
- \partial_z F_y^\star \\ |
95 |
where $F_x$ and $F_y$ are stream-functions with boundary conditions |
+ \partial_z F_x^\star \\ |
96 |
$F_x=F_y=0$ on upper and lower boundaries. The virtue of casting the |
\partial_x F_y^\star - \partial_y F_x^\star |
97 |
bolus velocity in terms of these stream-functions is that they are |
\end{array} \right), |
98 |
automatically non-divergent ($\partial_x u^* + \partial_y v^* = - |
\end{equation} |
99 |
\partial_{xz} F_x - \partial_{yz} F_y = - \partial_z w^*$). $F_x$ and |
and thus is automatically non-divergent. In the GM parameterization, the streamfunction is |
100 |
$F_y$ are specified in terms of the isoneutral slopes $S_x$ and $S_y$: |
specified in terms of the isoneutral slopes $S_x$ and $S_y$: |
101 |
\begin{eqnarray} |
\begin{eqnarray} |
102 |
F_x & = & \kappa_{GM} S_x \\ |
F_x^\star & = & -\kappa_{GM} S_y \\ |
103 |
F_y & = & \kappa_{GM} S_y |
F_y^\star & = & \kappa_{GM} S_x |
104 |
\end{eqnarray} |
\end{eqnarray} |
105 |
|
with boundary conditions $F_x^\star=F_y^\star=0$ on upper and lower boundaries. |
106 |
|
In the end, the bolus transport in the GM parameterization is given by: |
107 |
|
\begin{equation} |
108 |
|
\bf{u}^\star = \left( |
109 |
|
\begin{array}{c} |
110 |
|
u^\star \\ |
111 |
|
v^\star \\ |
112 |
|
w^\star |
113 |
|
\end{array} |
114 |
|
\right) = \left( |
115 |
|
\begin{array}{c} |
116 |
|
- \partial_z (\kappa_{GM} S_x) \\ |
117 |
|
- \partial_z (\kappa_{GM} S_y) \\ |
118 |
|
\partial_x (\kappa_{GM} S_x) + \partial_y (\kappa_{GM} S_y) |
119 |
|
\end{array} |
120 |
|
\right) |
121 |
|
\end{equation} |
122 |
|
|
123 |
This is the form of the GM parameterization as applied by Donabasaglu, |
This is the form of the GM parameterization as applied by Donabasaglu, |
124 |
1997, in MOM versions 1 and 2. |
1997, in MOM versions 1 and 2. |
125 |
|
|
126 |
\subsection{Griffies Skew Flux} |
Note that in the MITgcm, the variables containing the GM bolus streamfunction are: |
127 |
|
\begin{equation} |
128 |
|
\left( |
129 |
|
\begin{array}{c} |
130 |
|
GM\_PsiX \\ |
131 |
|
GM\_PsiY |
132 |
|
\end{array} |
133 |
|
\right) = \left( |
134 |
|
\begin{array}{c} |
135 |
|
\kappa_{GM} S_x \\ |
136 |
|
\kappa_{GM} S_y |
137 |
|
\end{array} |
138 |
|
\right)= \left( |
139 |
|
\begin{array}{c} |
140 |
|
F_y^\star \\ |
141 |
|
-F_x^\star |
142 |
|
\end{array} |
143 |
|
\right). |
144 |
|
\end{equation} |
145 |
|
|
146 |
|
\subsubsection{Griffies Skew Flux} |
147 |
|
|
148 |
Griffies notes that the discretisation of bolus velocities involves |
\cite{gr:98} notes that the discretisation of bolus velocities involves |
149 |
multiple layers of differencing and interpolation that potentially |
multiple layers of differencing and interpolation that potentially |
150 |
lead to noisy fields and computational modes. He pointed out that the |
lead to noisy fields and computational modes. He pointed out that the |
151 |
bolus flux can be re-written in terms of a non-divergent flux and a |
bolus flux can be re-written in terms of a non-divergent flux and a |
152 |
skew-flux: |
skew-flux: |
153 |
\begin{eqnarray*} |
\begin{eqnarray*} |
154 |
\bf{u}^* \tau |
\bf{u}^\star \tau |
155 |
& = & |
& = & |
156 |
\left( \begin{array}{c} |
\left( \begin{array}{c} |
157 |
- \partial_z ( \kappa_{GM} S_x ) \tau \\ |
- \partial_z ( \kappa_{GM} S_x ) \tau \\ |
163 |
\left( \begin{array}{c} |
\left( \begin{array}{c} |
164 |
- \partial_z ( \kappa_{GM} S_x \tau) \\ |
- \partial_z ( \kappa_{GM} S_x \tau) \\ |
165 |
- \partial_z ( \kappa_{GM} S_y \tau) \\ |
- \partial_z ( \kappa_{GM} S_y \tau) \\ |
166 |
\partial_x ( \kappa_{GM} S_x \tau) + \partial_y ( \kappa_{GM} S_y) \tau) |
\partial_x ( \kappa_{GM} S_x \tau) + \partial_y ( \kappa_{GM} S_y \tau) |
167 |
\end{array} \right) |
\end{array} \right) |
168 |
+ \left( \begin{array}{c} |
+ \left( \begin{array}{c} |
169 |
\kappa_{GM} S_x \partial_z \tau \\ |
\kappa_{GM} S_x \partial_z \tau \\ |
170 |
\kappa_{GM} S_y \partial_z \tau \\ |
\kappa_{GM} S_y \partial_z \tau \\ |
171 |
- \kappa_{GM} S_x \partial_x \tau - \kappa_{GM} S_y) \partial_y \tau |
- \kappa_{GM} S_x \partial_x \tau - \kappa_{GM} S_y \partial_y \tau |
172 |
\end{array} \right) |
\end{array} \right) |
173 |
\end{eqnarray*} |
\end{eqnarray*} |
174 |
The first vector is non-divergent and thus has no effect on the tracer |
The first vector is non-divergent and thus has no effect on the tracer |
175 |
field and can be dropped. The remaining flux can be written: |
field and can be dropped. The remaining flux can be written: |
176 |
\begin{equation} |
\begin{equation} |
177 |
\bf{u}^* \tau = - \kappa_{GM} \bf{K}_{GM} \bf{\nabla} \tau |
\bf{u}^\star \tau = - \kappa_{GM} \bf{K}_{GM} \bf{\nabla} \tau |
178 |
\end{equation} |
\end{equation} |
179 |
where |
where |
180 |
\begin{equation} |
\begin{equation} |
196 |
with the Redi isoneutral mixing scheme: |
with the Redi isoneutral mixing scheme: |
197 |
\begin{equation} |
\begin{equation} |
198 |
\kappa_\rho \bf{K}_{Redi} \bf{\nabla} \tau |
\kappa_\rho \bf{K}_{Redi} \bf{\nabla} \tau |
199 |
- u^* \tau = |
- u^\star \tau = |
200 |
( \kappa_\rho \bf{K}_{Redi} + \kappa_{GM} \bf{K}_{GM} ) \bf{\nabla} \tau |
( \kappa_\rho \bf{K}_{Redi} + \kappa_{GM} \bf{K}_{GM} ) \bf{\nabla} \tau |
201 |
\end{equation} |
\end{equation} |
202 |
In the instance that $\kappa_{GM} = \kappa_{\rho}$ then |
In the instance that $\kappa_{GM} = \kappa_{\rho}$ then |
210 |
\end{array} |
\end{array} |
211 |
\right) |
\right) |
212 |
\end{equation} |
\end{equation} |
213 |
which differs from the variable laplacian diffusion tensor by only |
which differs from the variable Laplacian diffusion tensor by only |
214 |
two non-zero elements in the $z$-row. |
two non-zero elements in the $z$-row. |
215 |
|
|
216 |
\subsection{Variable $\kappa_{GM}$} |
\fbox{ \begin{minipage}{4.75in} |
217 |
|
{\em S/R GMREDI\_CALC\_TENSOR} ({\em pkg/gmredi/gmredi\_calc\_tensor.F}) |
218 |
|
|
219 |
|
$\sigma_x$: {\bf SlopeX} (argument on entry) |
220 |
|
|
221 |
|
$\sigma_y$: {\bf SlopeY} (argument on entry) |
222 |
|
|
223 |
|
$\sigma_z$: {\bf SlopeY} (argument) |
224 |
|
|
225 |
|
$S_x$: {\bf SlopeX} (argument on exit) |
226 |
|
|
227 |
|
$S_y$: {\bf SlopeY} (argument on exit) |
228 |
|
|
229 |
|
\end{minipage} } |
230 |
|
|
231 |
|
|
232 |
Visbeck et al., 1996, suggest making the eddy coefficient, |
|
233 |
|
\subsubsection{Variable $\kappa_{GM}$} |
234 |
|
|
235 |
|
\cite{visbeck:97} suggest making the eddy coefficient, |
236 |
$\kappa_{GM}$, a function of the Eady growth rate, |
$\kappa_{GM}$, a function of the Eady growth rate, |
237 |
$|f|/\sqrt{Ri}$. The formula involves a non-dimensional constant, |
$|f|/\sqrt{Ri}$. The formula involves a non-dimensional constant, |
238 |
$\alpha$, and a length-scale $L$: |
$\alpha$, and a length-scale $L$: |
256 |
\end{displaymath} |
\end{displaymath} |
257 |
|
|
258 |
|
|
259 |
\subsection{Tapering and stability} |
\subsubsection{Tapering and stability} |
260 |
|
|
261 |
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 |
262 |
matched to the convective parameterization. This was originally |
matched to the convective parameterization. This was originally |
263 |
expressed in connection with the introduction of the KPP boundary |
expressed in connection with the introduction of the KPP boundary |
264 |
layer scheme (Large et al., 97) but infact, as subsequent experience |
layer scheme \citep{lar-eta:94} but in fact, as subsequent experience |
265 |
with the MIT model has found, is necessary for any convective |
with the MIT model has found, is necessary for any convective |
266 |
parameterization. |
parameterization. |
267 |
|
|
279 |
|
|
280 |
\end{minipage} } |
\end{minipage} } |
281 |
|
|
282 |
|
\begin{figure} |
283 |
|
\begin{center} |
284 |
|
\resizebox{5.0in}{3.0in}{\includegraphics{s_phys_pkgs/figs/tapers.eps}} |
285 |
|
\end{center} |
286 |
|
\caption{Taper functions used in GKW91 and DM95.} |
287 |
|
\label{fig:tapers} |
288 |
|
\end{figure} |
289 |
|
|
290 |
|
\begin{figure} |
291 |
|
\begin{center} |
292 |
|
\resizebox{5.0in}{3.0in}{\includegraphics{s_phys_pkgs/figs/effective_slopes.eps}} |
293 |
|
\end{center} |
294 |
|
\caption{Effective slope as a function of ``true'' slope using Cox |
295 |
|
slope clipping, GKW91 limiting and DM95 limiting.} |
296 |
|
\label{fig:effective_slopes} |
297 |
|
\end{figure} |
298 |
|
|
299 |
\subsubsection{Slope clipping} |
|
300 |
|
\subsubsection*{Slope clipping} |
301 |
|
|
302 |
Deep convection sites and the mixed layer are indicated by |
Deep convection sites and the mixed layer are indicated by |
303 |
homogenized, unstable or nearly unstable stratification. The slopes in |
homogenized, unstable or nearly unstable stratification. The slopes in |
304 |
such regions can be either infinite, very large with a sign reversal |
such regions can be either infinite, very large with a sign reversal |
305 |
or simply very large. From a numerical point of view, large slopes |
or simply very large. From a numerical point of view, large slopes |
306 |
lead to large variations in the tensor elements (implying large bolus |
lead to large variations in the tensor elements (implying large bolus |
307 |
flow) and can be numerically unstable. This was first reognized by |
flow) and can be numerically unstable. This was first recognized by |
308 |
Cox, 1987, who implemented ``slope clipping'' in the isopycnal mixing |
\cite{Cox87} who implemented ``slope clipping'' in the isopycnal mixing |
309 |
tensor. Here, the slope magnitude is simply restricted by an upper |
tensor. Here, the slope magnitude is simply restricted by an upper |
310 |
limit: |
limit: |
311 |
\begin{eqnarray} |
\begin{eqnarray} |
339 |
diffusion). The classic result of dramatically reduced mixed layers is |
diffusion). The classic result of dramatically reduced mixed layers is |
340 |
evident. Indeed, the deep convection sites to just one or two points |
evident. Indeed, the deep convection sites to just one or two points |
341 |
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, |
342 |
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 |
343 |
parameterization. Limiting the slopes also breaks the adiabatic nature |
parameterization. Limiting the slopes also breaks the adiabatic nature |
344 |
of the GM/Redi parameterization, re-introducing diabatic fluxes in |
of the GM/Redi parameterization, re-introducing diabatic fluxes in |
345 |
regions where the limiting is in effect. |
regions where the limiting is in effect. |
346 |
|
|
347 |
\subsubsection{Tapering: Gerdes, Koberle and Willebrand, Clim. Dyn. 1991} |
\subsubsection*{Tapering: Gerdes, Koberle and Willebrand, Clim. Dyn. 1991} |
348 |
|
|
349 |
The tapering scheme used in Gerdes et al., 1991, (\cite{gkw91}) |
The tapering scheme used in \cite{gkw:91} |
350 |
addressed two issues with the clipping method: the introduction of |
addressed two issues with the clipping method: the introduction of |
351 |
large vertical fluxes in addition to convective adjustment fluxes is |
large vertical fluxes in addition to convective adjustment fluxes is |
352 |
avoided by tapering the GM/Redi slopes back to zero in |
avoided by tapering the GM/Redi slopes back to zero in |
365 |
that the effective vertical diffusivity term $\kappa f_1(S) |S|^2 = |
that the effective vertical diffusivity term $\kappa f_1(S) |S|^2 = |
366 |
\kappa S_{max}^2$. |
\kappa S_{max}^2$. |
367 |
|
|
368 |
The GKW tapering scheme is activated in the model by setting {\bf |
The GKW91 tapering scheme is activated in the model by setting {\bf |
369 |
GM\_tap\-er\_scheme = 'gkw91'} in {\em data.gmredi}. |
GM\_tap\-er\_scheme = 'gkw91'} in {\em data.gmredi}. |
370 |
|
|
371 |
\subsection{Tapering: Danabasoglu and McWilliams, J. Clim. 1995} |
\subsubsection*{Tapering: Danabasoglu and McWilliams, J. Clim. 1995} |
372 |
|
|
373 |
The tapering scheme used by Danabasoglu and McWilliams, 1995, |
The tapering scheme used by \cite{dm:95} followed a similar procedure but used a different |
|
\cite{DM95}, followed a similar procedure but used a different |
|
374 |
tapering function, $f_1(S)$: |
tapering function, $f_1(S)$: |
375 |
\begin{equation} |
\begin{equation} |
376 |
f_1(S) = \frac{1}{2} \left( 1+\tanh \left[ \frac{S_c - |S|}{S_d} \right] \right) |
f_1(S) = \frac{1}{2} \left( 1+\tanh \left[ \frac{S_c - |S|}{S_d} \right] \right) |
381 |
cut-off, turning off the GM/Redi SGS parameterization for weaker |
cut-off, turning off the GM/Redi SGS parameterization for weaker |
382 |
slopes. |
slopes. |
383 |
|
|
384 |
The DM tapering scheme is activated in the model by setting {\bf |
The DM95 tapering scheme is activated in the model by setting {\bf |
385 |
GM\_tap\-er\_scheme = 'dm95'} in {\em data.gmredi}. |
GM\_tap\-er\_scheme = 'dm95'} in {\em data.gmredi}. |
386 |
|
|
387 |
\subsection{Tapering: Large, Danabasoglu and Doney, JPO 1997} |
\subsubsection*{Tapering: Large, Danabasoglu and Doney, JPO 1997} |
388 |
|
|
389 |
The tapering used in Large et al., 1997, \cite{ldd97}, is based on the |
The tapering used in \cite{lar-eta:97} is based on the |
390 |
DM95 tapering scheme, but also tapers the scheme with an additional |
DM95 tapering scheme, but also tapers the scheme with an additional |
391 |
function of height, $f_2(z)$, so that the GM/Redi SGS fluxes are |
function of height, $f_2(z)$, so that the GM/Redi SGS fluxes are |
392 |
reduced near the surface: |
reduced near the surface: |
393 |
\begin{equation} |
\begin{equation} |
394 |
f_2(S) = \frac{1}{2} \left( 1 + \sin(\pi \frac{z}{D} - \pi/2)\right) |
f_2(z) = \frac{1}{2} \left( 1 + \sin(\pi \frac{z}{D} - \frac{\pi}{2})\right) |
395 |
\end{equation} |
\end{equation} |
396 |
where $D = L_\rho |S|$ is a depth-scale and $L_\rho=c/f$ with |
where $D = L_\rho |S|$ is a depth-scale and $L_\rho=c/f$ with |
397 |
$c=2$~m~s$^{-1}$. This tapering with height was introduced to fix |
$c=2$~m~s$^{-1}$. This tapering with height was introduced to fix |
398 |
some spurious interaction with the mixed-layer KPP parameterization. |
some spurious interaction with the mixed-layer KPP parameterization. |
399 |
|
|
400 |
The LDD tapering scheme is activated in the model by setting {\bf |
The LDD97 tapering scheme is activated in the model by setting {\bf |
401 |
GM\_tap\-er\_scheme = 'ldd97'} in {\em data.gmredi}. |
GM\_tap\-er\_scheme = 'ldd97'} in {\em data.gmredi}. |
402 |
|
|
403 |
|
|
404 |
|
|
405 |
|
|
406 |
\begin{figure} |
\begin{figure} |
407 |
|
\begin{center} |
408 |
%\includegraphics{mixedlayer-cox.eps} |
%\includegraphics{mixedlayer-cox.eps} |
409 |
%\includegraphics{mixedlayer-diff.eps} |
%\includegraphics{mixedlayer-diff.eps} |
410 |
|
Figure missing. |
411 |
|
\end{center} |
412 |
\caption{Mixed layer depth using GM parameterization with a) Cox slope |
\caption{Mixed layer depth using GM parameterization with a) Cox slope |
413 |
clipping and for comparison b) using horizontal constant diffusion.} |
clipping and for comparison b) using horizontal constant diffusion.} |
414 |
\ref{fig-mixedlayer} |
\label{fig-mixedlayer} |
415 |
\end{figure} |
\end{figure} |
416 |
|
|
417 |
\begin{figure} |
\subsubsection{Package Reference} |
418 |
%\includegraphics{slopelimits.eps} |
\label{sec:pkg:gmredi:diagnostics} |
|
\caption{Effective slope as a function of ``true'' slope using a) Cox |
|
|
slope clipping, b) GKW91 limiting, c) DM95 limiting and d) LDD97 |
|
|
limiting.} |
|
|
\end{figure} |
|
|
|
|
|
|
|
|
%\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()} |
|
419 |
|
|
420 |
{\small |
{\footnotesize |
421 |
\begin{verbatim} |
\begin{verbatim} |
422 |
subroutine RediTensor(Temp,Salt,Kredigm,K31,K32,K33, nIter,DumpFlag) |
------------------------------------------------------------------------ |
423 |
|---in--| |-------out-------| |
<-Name->|Levs|<-parsing code->|<-- Units -->|<- Tile (max=80c) |
424 |
! Input |
------------------------------------------------------------------------ |
425 |
real Temp(Nx,Ny,Nz) ! Potential temperature |
GM_VisbK| 1 |SM P M1 |m^2/s |Mixing coefficient from Visbeck etal parameterization |
426 |
real Salt(Nx,Ny,Nz) ! Salinity |
GM_Kux | 15 |UU P 177MR |m^2/s |K_11 element (U.point, X.dir) of GM-Redi tensor |
427 |
! Output |
GM_Kvy | 15 |VV P 176MR |m^2/s |K_22 element (V.point, Y.dir) of GM-Redi tensor |
428 |
real Kredigm(Nx,Ny,Nz) ! Redi/GM eddy coefficient |
GM_Kuz | 15 |UU 179MR |m^2/s |K_13 element (U.point, Z.dir) of GM-Redi tensor |
429 |
real K31(Nx,Ny,Nz) ! Redi/GM (3,1) tensor component |
GM_Kvz | 15 |VV 178MR |m^2/s |K_23 element (V.point, Z.dir) of GM-Redi tensor |
430 |
real K32(Nx,Ny,Nz) ! Redi/GM (3,2) tensor component |
GM_Kwx | 15 |UM 181LR |m^2/s |K_31 element (W.point, X.dir) of GM-Redi tensor |
431 |
real K33(Nx,Ny,Nz) ! Redi/GM (3,3) tensor component |
GM_Kwy | 15 |VM 180LR |m^2/s |K_32 element (W.point, Y.dir) of GM-Redi tensor |
432 |
! Auxiliary input |
GM_Kwz | 15 |WM P LR |m^2/s |K_33 element (W.point, Z.dir) of GM-Redi tensor |
433 |
integer nIter ! interation/time-step number |
GM_PsiX | 15 |UU 184LR |m^2/s |GM Bolus transport stream-function : X component |
434 |
logical DumpFlag ! flag to indicate routine should ``dump'' |
GM_PsiY | 15 |VV 183LR |m^2/s |GM Bolus transport stream-function : Y component |
435 |
|
GM_KuzTz| 15 |UU 186MR |degC.m^3/s |Redi Off-diagonal Tempetature flux: X component |
436 |
|
GM_KvzTz| 15 |VV 185MR |degC.m^3/s |Redi Off-diagonal Tempetature flux: Y component |
437 |
\end{verbatim} |
\end{verbatim} |
438 |
} |
} |
439 |
|
|
440 |
The subroutine {\bf RediTensor()} is called from {\bf model()} with |
\subsubsection{Experiments and tutorials that use gmredi} |
441 |
input arguments $T$ and $S$. It returns the 3D-arrays {\tt Kredigm}, |
\label{sec:pkg:gmredi:experiments} |
|
{\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()} |
|
442 |
|
|
443 |
{\bf inc\-tracer()} is called from {\bf model()} and has {\em four |
\begin{itemize} |
444 |
new} arguments: |
\item{Global Ocean tutorial, in tutorial\_global\_oce\_latlon verification directory, |
445 |
\begin{verbatim} |
described in section \ref{sec:eg-global} } |
446 |
subroutine inc_tracer( ...,Kredigm,K31,K32,K33, ... ) |
\item{ Front Relax experiment, in front\_relax verification directory.} |
447 |
real Kredigm(Nx,Ny,Nz) ! Eddy coefficient |
\item{ Ideal 2D Ocean experiment, in ideal\_2D\_oce verification directory.} |
448 |
real K31(Nx,Ny,Nz) ! (3,1) tensor coefficient |
\end{itemize} |
|
real K32(Nx,Ny,Nz) ! (3,2) tensor coefficient |
|
|
real K33(Nx,Ny,Nz) ! (3,3) tensor coefficient |
|
|
\end{verbatim} |
|
449 |
|
|
450 |
Within the routine, the lateral fluxes, {\tt fluxWest} and {\tt |
% DO NOT EDIT HERE |
|
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*} |
|
451 |
|
|
|
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. |
|
452 |
|
|
453 |
|
|