2 |
% $Name$ |
% $Name$ |
3 |
|
|
4 |
\section{Example: Four layer Baroclinic Ocean Gyre In Spherical Coordinates} |
\section{Example: Four layer Baroclinic Ocean Gyre In Spherical Coordinates} |
5 |
|
\label{sec:eg-fourlayer} |
6 |
|
|
7 |
\bodytext{bgcolor="#FFFFFFFF"} |
\bodytext{bgcolor="#FFFFFFFF"} |
8 |
|
|
16 |
%{\large May 2001} |
%{\large May 2001} |
17 |
%\end{center} |
%\end{center} |
18 |
|
|
19 |
\subsection{Introduction} |
This document describes an example experiment using MITgcm |
20 |
|
to simulate a baroclinic ocean gyre in spherical |
21 |
This document describes the second example MITgcm experiment. The first |
polar coordinates. The barotropic |
22 |
example experiment ilustrated how to configure the code for a single layer |
example experiment in section \ref{sec:eg-baro} |
23 |
|
ilustrated how to configure the code for a single layer |
24 |
simulation in a cartesian grid. In this example a similar physical problem |
simulation in a cartesian grid. In this example a similar physical problem |
25 |
is simulated, but the code is now configured |
is simulated, but the code is now configured |
26 |
for four layers and in a spherical polar coordinate system. |
for four layers and in a spherical polar coordinate system. |
37 |
In this experiment the model is configured to represent a mid-latitude |
In this experiment the model is configured to represent a mid-latitude |
38 |
enclosed sector of fluid on a sphere, $60^{\circ} \times 60^{\circ}$ in |
enclosed sector of fluid on a sphere, $60^{\circ} \times 60^{\circ}$ in |
39 |
lateral extent. The fluid is $2$~km deep and is forced |
lateral extent. The fluid is $2$~km deep and is forced |
40 |
by a constant in time zonal wind stress, $\tau_x$, that varies sinusoidally |
by a constant in time zonal wind stress, $\tau_{\lambda}$, that varies |
41 |
in the north-south direction. Topologically the simulated |
sinusoidally in the north-south direction. Topologically the simulated |
42 |
domain is a sector on a sphere and the coriolis parameter, $f$, is defined |
domain is a sector on a sphere and the coriolis parameter, $f$, is defined |
43 |
according to latitude, $\phi$ |
according to latitude, $\varphi$ |
44 |
|
|
45 |
\begin{equation} |
\begin{equation} |
46 |
\label{EQ:fcori} |
\label{EQ:fcori} |
47 |
f(\phi) = 2 \Omega \sin( \phi ) |
f(\varphi) = 2 \Omega \sin( \varphi ) |
48 |
\end{equation} |
\end{equation} |
49 |
|
|
50 |
\noindent with the rotation rate, $\Omega$ set to $\frac{2 \pi}{86400s}$. |
\noindent with the rotation rate, $\Omega$ set to $\frac{2 \pi}{86400s}$. |
54 |
|
|
55 |
\begin{equation} |
\begin{equation} |
56 |
\label{EQ:taux} |
\label{EQ:taux} |
57 |
\tau_x(\phi) = \tau_{0}\sin(\pi \frac{\phi}{L_{\phi}}) |
\tau_{\lambda}(\varphi) = \tau_{0}\sin(\pi \frac{\varphi}{L_{\varphi}}) |
58 |
\end{equation} |
\end{equation} |
59 |
|
|
60 |
\noindent where $L_{\phi}$ is the lateral domain extent ($60^{\circ}$) and |
\noindent where $L_{\varphi}$ is the lateral domain extent ($60^{\circ}$) and |
61 |
$\tau_0$ is set to $0.1N m^{-2}$. |
$\tau_0$ is set to $0.1N m^{-2}$. |
62 |
\\ |
\\ |
63 |
|
|
64 |
Figure \ref{FIG:simulation_config} |
Figure \ref{FIG:simulation_config} |
65 |
summarises the configuration simulated. |
summarises the configuration simulated. |
66 |
In contrast to example (1) \cite{baro_gyre_case_study}, the current |
In contrast to the example in section \ref{sec:eg-baro}, the |
67 |
experiment simulates a spherical polar domain. However, as indicated |
current experiment simulates a spherical polar domain. As indicated |
68 |
by the axes in the lower left of the figure the model code works internally |
by the axes in the lower left of the figure the model code works internally |
69 |
in a locally orthoganal coordinate $(x,y,z)$. In the remainder of this |
in a locally orthoganal coordinate $(x,y,z)$. For this experiment description |
70 |
document the local coordinate $(x,y,z)$ will be adopted. |
of this document the local orthogonal model coordinate $(x,y,z)$ is synonomous |
71 |
|
with the spherical polar coordinate shown in figure |
72 |
|
\ref{fig:spherical-polar-coord} |
73 |
\\ |
\\ |
74 |
|
|
75 |
The experiment has four levels in the vertical, each of equal thickness, |
The experiment has four levels in the vertical, each of equal thickness, |
95 |
|
|
96 |
\noindent with $\rho_{0}=999.8\,{\rm kg\,m}^{-3}$ and |
\noindent with $\rho_{0}=999.8\,{\rm kg\,m}^{-3}$ and |
97 |
$\alpha_{\theta}=2\times10^{-4}\,{\rm degrees}^{-1} $. Integrated forward in |
$\alpha_{\theta}=2\times10^{-4}\,{\rm degrees}^{-1} $. Integrated forward in |
98 |
this configuration the model state variable {\bf theta} is synonomous with |
this configuration the model state variable {\bf theta} is equivalent to |
99 |
either in-situ temperature, $T$, or potential temperature, $\theta$. For |
either in-situ temperature, $T$, or potential temperature, $\theta$. For |
100 |
consistency with later examples, in which the equation of state is |
consistency with later examples, in which the equation of state is |
101 |
non-linear, we use $\theta$ to represent temperature here. This is |
non-linear, we use $\theta$ to represent temperature here. This is |
102 |
the quantity that is carried in the model core equations. |
the quantity that is carried in the model core equations. |
103 |
|
|
104 |
\begin{figure} |
\begin{figure} |
105 |
\centerline{ |
\begin{center} |
106 |
\resizebox{7.5in}{5.5in}{ |
\resizebox{7.5in}{5.5in}{ |
107 |
\includegraphics*[0.2in,0.7in][10.5in,10.5in] |
\includegraphics*[0.2in,0.7in][10.5in,10.5in] |
108 |
{part3/case_studies/fourlayer_gyre/simulation_config.eps} } |
{part3/case_studies/fourlayer_gyre/simulation_config.eps} } |
109 |
} |
\end{center} |
110 |
\caption{Schematic of simulation domain and wind-stress forcing function |
\caption{Schematic of simulation domain and wind-stress forcing function |
111 |
for the four-layer gyre numerical experiment. The domain is enclosed by solid |
for the four-layer gyre numerical experiment. The domain is enclosed by solid |
112 |
walls at $0^{\circ}$~E, $60^{\circ}$~E, $0^{\circ}$~N and $60^{\circ}$~N. |
walls at $0^{\circ}$~E, $60^{\circ}$~E, $0^{\circ}$~N and $60^{\circ}$~N. |
113 |
In the four-layer case an initial temperature stratification is |
An initial stratification is |
114 |
imposed by setting the potential temperature, $\theta$, in each layer. |
imposed by setting the potential temperature, $\theta$, in each layer. |
115 |
The vertical spacing, $\Delta z$, is constant and equal to $500$m. |
The vertical spacing, $\Delta z$, is constant and equal to $500$m. |
116 |
} |
} |
117 |
\label{FIG:simulation_config} |
\label{FIG:simulation_config} |
118 |
\end{figure} |
\end{figure} |
119 |
|
|
120 |
\subsection{Discrete Numerical Configuration} |
\subsection{Equations solved} |
121 |
|
|
122 |
The model is configured in hydrostatic form. The domain is discretised with |
The implicit free surface {\bf HPE} form of the |
123 |
a uniform grid spacing in latitude and longitude |
equations described in Marshall et. al \cite{Marshall97a} is |
124 |
$\Delta x=\Delta y=1^{\circ}$, so |
employed. The flow is three-dimensional with just temperature, $\theta$, as |
125 |
that there are sixty grid cells in the $x$ and $y$ directions. Vertically the |
an active tracer. The equation of state is linear. |
|
model is configured with a four layers with constant depth, |
|
|
$\Delta z$, of $500$~m. |
|
|
The implicit free surface form of the |
|
|
pressure equation described in Marshall et. al \cite{Marshall97a} is |
|
|
employed. |
|
126 |
A horizontal laplacian operator $\nabla_{h}^2$ provides viscous |
A horizontal laplacian operator $\nabla_{h}^2$ provides viscous |
127 |
dissipation. The wind-stress momentum input is added to the momentum equation |
dissipation and provides a diffusive sub-grid scale closure for the |
128 |
for the ``zonal flow'', $u$. Other terms in the model |
temperature equation. A wind-stress momentum forcing is added to the momentum |
129 |
|
equation for the zonal flow, $u$. Other terms in the model |
130 |
are explicitly switched off for this experiement configuration (see section |
are explicitly switched off for this experiement configuration (see section |
131 |
\ref{SEC:code_config} ), yielding an active set of equations solved in this |
\ref{SEC:eg_fourl_code_config} ). This yields an active set of equations |
132 |
configuration as follows |
solved in this configuration, written in spherical polar coordinates as |
133 |
|
follows |
134 |
|
|
135 |
\begin{eqnarray} |
\begin{eqnarray} |
136 |
\label{EQ:model_equations} |
\label{EQ:model_equations} |
137 |
\frac{Du}{Dt} - fv + |
\frac{Du}{Dt} - fv + |
138 |
\frac{1}{\rho}\frac{\partial p^{'}}{\partial x} - |
\frac{1}{\rho}\frac{\partial p^{\prime}}{\partial \lambda} - |
139 |
A_{h}\nabla_{h}^2u - A_{z}\frac{\partial^{2}u}{\partial z^{2}} |
A_{h}\nabla_{h}^2u - A_{z}\frac{\partial^{2}u}{\partial z^{2}} |
140 |
& = & |
& = & |
141 |
\cal{F} |
\cal{F}_{\lambda} |
142 |
\\ |
\\ |
143 |
\frac{Dv}{Dt} + fu + |
\frac{Dv}{Dt} + fu + |
144 |
\frac{1}{\rho}\frac{\partial p^{'}}{\partial y} - |
\frac{1}{\rho}\frac{\partial p^{\prime}}{\partial \varphi} - |
145 |
A_{h}\nabla_{h}^2v - A_{z}\frac{\partial^{2}v}{\partial z^{2}} |
A_{h}\nabla_{h}^2v - A_{z}\frac{\partial^{2}v}{\partial z^{2}} |
146 |
& = & |
& = & |
147 |
0 |
0 |
148 |
\\ |
\\ |
149 |
\frac{\partial \eta}{\partial t} + \nabla_{h}\cdot \vec{u} |
\frac{\partial \eta}{\partial t} + \frac{\partial H \widehat{u}}{\partial \lambda} + |
150 |
|
\frac{\partial H \widehat{v}}{\partial \varphi} |
151 |
&=& |
&=& |
152 |
0 |
0 |
153 |
|
\label{eq:fourl_example_continuity} |
154 |
\\ |
\\ |
155 |
\frac{D\theta}{Dt} - |
\frac{D\theta}{Dt} - |
156 |
K_{h}\nabla_{h}^2\theta - K_{z}\frac{\partial^{2}\theta}{\partial z^{2}} |
K_{h}\nabla_{h}^2\theta - K_{z}\frac{\partial^{2}\theta}{\partial z^{2}} |
157 |
& = & |
& = & |
158 |
0 |
0 |
159 |
|
\label{eq:eg_fourl_theta} |
160 |
\\ |
\\ |
161 |
g\rho_{0} \eta + \int^{0}_{-z}\rho^{'} dz & = & p^{'} |
p^{\prime} & = & |
162 |
|
g\rho_{0} \eta + \int^{0}_{-z}\rho^{\prime} dz |
163 |
\\ |
\\ |
164 |
{\cal F} |_{s} & = & \frac{\tau_{x}}{\rho_{0}\Delta z_{s}} |
\rho^{\prime} & = &- \alpha_{\theta}\rho_{0}\theta^{\prime} |
165 |
\\ |
\\ |
166 |
{\cal F} |_{i} & = & 0 |
{\cal F}_{\lambda} |_{s} & = & \frac{\tau_{\lambda}}{\rho_{0}\Delta z_{s}} |
167 |
|
\\ |
168 |
|
{\cal F}_{\lambda} |_{i} & = & 0 |
169 |
\end{eqnarray} |
\end{eqnarray} |
170 |
|
|
171 |
\noindent where $u$ and $v$ are the $x$ and $y$ components of the |
\noindent where $u$ and $v$ are the components of the horizontal |
172 |
flow vector $\vec{u}$. The suffices ${s},{i}$ indicate surface and |
flow vector $\vec{u}$ on the sphere ($u=\dot{\lambda},v=\dot{\varphi}$). |
173 |
interior model levels respectively. As described in |
The terms $H\widehat{u}$ and $H\widehat{v}$ are the components of the vertical |
174 |
MITgcm Numerical Solution Procedure \cite{MITgcm_Numerical_Scheme}, the time |
integral term given in equation \ref{eq:free-surface} and |
175 |
evolution of potential temperature, $\theta$, equation is solved prognostically. |
explained in more detail in section \ref{sect:pressure-method-linear-backward}. |
176 |
The total pressure, $p$, is diagnosed by summing pressure due to surface |
However, for the problem presented here, the continuity relation (equation |
177 |
elevation $\eta$ and the hydrostatic pressure. |
\ref{eq:fourl_example_continuity}) differs from the general form given |
178 |
\\ |
in section \ref{sect:pressure-method-linear-backward}, |
179 |
|
equation \ref{eq:linear-free-surface=P-E+R}, because the source terms |
180 |
|
${\cal P}-{\cal E}+{\cal R}$ |
181 |
|
are all $0$. |
182 |
|
|
183 |
|
The pressure field, $p^{\prime}$, is separated into a barotropic part |
184 |
|
due to variations in sea-surface height, $\eta$, and a hydrostatic |
185 |
|
part due to variations in density, $\rho^{\prime}$, integrated |
186 |
|
through the water column. |
187 |
|
|
188 |
|
The suffices ${s},{i}$ indicate surface and interior of the domain. |
189 |
|
The windstress forcing, ${\cal F}_{\lambda}$, is applied in the surface layer |
190 |
|
by a source term in the zonal momentum equation. In the ocean interior |
191 |
|
this term is zero. |
192 |
|
|
193 |
|
In the momentum equations |
194 |
|
lateral and vertical boundary conditions for the $\nabla_{h}^{2}$ |
195 |
|
and $\frac{\partial^{2}}{\partial z^{2}}$ operators are specified |
196 |
|
when the numerical simulation is run - see section |
197 |
|
\ref{SEC:eg_fourl_code_config}. For temperature |
198 |
|
the boundary condition is ``zero-flux'' |
199 |
|
e.g. $\frac{\partial \theta}{\partial \varphi}= |
200 |
|
\frac{\partial \theta}{\partial \lambda}=\frac{\partial \theta}{\partial z}=0$. |
201 |
|
|
202 |
|
|
203 |
|
|
204 |
|
\subsection{Discrete Numerical Configuration} |
205 |
|
|
206 |
|
The domain is discretised with |
207 |
|
a uniform grid spacing in latitude and longitude |
208 |
|
$\Delta \lambda=\Delta \varphi=1^{\circ}$, so |
209 |
|
that there are sixty grid cells in the zonal and meridional directions. |
210 |
|
Vertically the |
211 |
|
model is configured with four layers with constant depth, |
212 |
|
$\Delta z$, of $500$~m. The internal, locally orthogonal, model coordinate |
213 |
|
variables $x$ and $y$ are initialised from the values of |
214 |
|
$\lambda$, $\varphi$, $\Delta \lambda$ and $\Delta \varphi$ in |
215 |
|
radians according to |
216 |
|
|
217 |
|
\begin{eqnarray} |
218 |
|
x=r\cos(\varphi)\lambda,~\Delta x & = &r\cos(\varphi)\Delta \lambda \\ |
219 |
|
y=r\varphi,~\Delta y &= &r\Delta \varphi |
220 |
|
\end{eqnarray} |
221 |
|
|
222 |
|
The procedure for generating a set of internal grid variables from a |
223 |
|
spherical polar grid specification is discussed in section |
224 |
|
\ref{sec:spatial_discrete_horizontal_grid}. |
225 |
|
|
226 |
|
\noindent\fbox{ \begin{minipage}{5.5in} |
227 |
|
{\em S/R INI\_SPHERICAL\_POLAR\_GRID} ({\em |
228 |
|
model/src/ini\_spherical\_polar\_grid.F}) |
229 |
|
|
230 |
|
$A_c$, $A_\zeta$, $A_w$, $A_s$: {\bf rAc}, {\bf rAz}, {\bf rAw}, {\bf rAs} |
231 |
|
({\em GRID.h}) |
232 |
|
|
233 |
|
$\Delta x_g$, $\Delta y_g$: {\bf DXg}, {\bf DYg} ({\em GRID.h}) |
234 |
|
|
235 |
|
$\Delta x_c$, $\Delta y_c$: {\bf DXc}, {\bf DYc} ({\em GRID.h}) |
236 |
|
|
237 |
|
$\Delta x_f$, $\Delta y_f$: {\bf DXf}, {\bf DYf} ({\em GRID.h}) |
238 |
|
|
239 |
|
$\Delta x_v$, $\Delta y_u$: {\bf DXv}, {\bf DYu} ({\em GRID.h}) |
240 |
|
|
241 |
|
\end{minipage} }\\ |
242 |
|
|
243 |
|
|
244 |
|
|
245 |
|
As described in \ref{sec:tracer_equations}, the time evolution of potential |
246 |
|
temperature, |
247 |
|
$\theta$, (equation \ref{eq:eg_fourl_theta}) |
248 |
|
is evaluated prognostically. The centered second-order scheme with |
249 |
|
Adams-Bashforth time stepping described in section |
250 |
|
\ref{sec:tracer_equations_abII} is used to step forward the temperature |
251 |
|
equation. The pressure forces that drive the fluid motions, ( |
252 |
|
$\frac{\partial p^{'}}{\partial \lambda}$ and $\frac{\partial p^{'}}{\partial \varphi}$), are found by summing pressure due to surface |
253 |
|
elevation $\eta$ and the hydrostatic pressure. The hydrostatic part of the |
254 |
|
pressure is evaluated explicitly by integrating density. The sea-surface |
255 |
|
height, $\eta$, is solved for implicitly as described in section |
256 |
|
\ref{sect:pressure-method-linear-backward}. |
257 |
|
|
258 |
\subsubsection{Numerical Stability Criteria} |
\subsubsection{Numerical Stability Criteria} |
259 |
|
|
268 |
\noindent of $\approx 100$km. This is greater than the model |
\noindent of $\approx 100$km. This is greater than the model |
269 |
resolution in mid-latitudes $\Delta x$, ensuring that the frictional |
resolution in mid-latitudes $\Delta x$, ensuring that the frictional |
270 |
boundary layer is well resolved. |
boundary layer is well resolved. |
271 |
|
|
272 |
\\ |
\\ |
273 |
|
|
274 |
\noindent The model is stepped forward with a |
\noindent The model is stepped forward with a |
335 |
stability limit of 0.25. |
stability limit of 0.25. |
336 |
|
|
337 |
\subsection{Code Configuration} |
\subsection{Code Configuration} |
338 |
\label{SEC:code_config} |
\label{SEC:eg_fourl_code_config} |
339 |
|
|
340 |
The model configuration for this experiment resides under the |
The model configuration for this experiment resides under the |
341 |
directory {\it verification/exp1/}. The experiment files |
directory {\it verification/exp1/}. The experiment files |