s_examples/global_oce_latlon/climatalogical_ogcm.tex

% $Header: /u/gcmpack/manual/s_examples/global_oce_latlon/climatalogical_ogcm.tex,v 1.19 2010/08/30 23:09:20 jmc Exp $
% $Name:  $

\section[Global Ocean MITgcm Example]{Global Ocean Simulation at $4^\circ$ Resolution}
%\label{www:tutorials}
\label{sec:eg-global}
\begin{rawhtml}
<!-- CMIREDIR:eg-global: -->
\end{rawhtml}
\begin{center}
(in directory: {\it verification/tutorial\_global\_oce\_latlon/})
\end{center}

\bodytext{bgcolor="#FFFFFFFF"}

%\begin{center} 
%{\Large \bf Using MITgcm to Simulate Global Climatological Ocean Circulation
%At Four Degree Resolution with Asynchronous Time Stepping}
%
%\vspace*{4mm}
%
%\vspace*{3mm}
%{\large May 2001}
%\end{center}


This example experiment demonstrates using the MITgcm to simulate
the planetary ocean circulation. The simulation is configured
with realistic geography and bathymetry on a
$4^{\circ} \times 4^{\circ}$ spherical polar grid.
The files for this experiment are in the verification directory
under tutorial\_global\_oce\_latlon.
Twenty levels are used in the vertical, ranging in thickness
from $50\,{\rm m}$ at the surface to $815\,{\rm m}$ at depth,
giving a maximum model depth of $6\,{\rm km}$.
At this resolution, the configuration
can be integrated forward for thousands of years on a single 
processor desktop computer.
\\
\subsection{Overview}
%\label{www:tutorials}

The model is forced with climatological wind stress data and surface
flux data from DaSilva \cite{DaSilva94}. Climatological data
from Levitus \cite{Levitus94} is used to initialize the model hydrography.
Levitus seasonal climatology data is also used throughout the calculation
to provide additional air-sea fluxes.
These fluxes are combined with the DaSilva climatological estimates of
surface heat flux and fresh water, resulting in a mixed boundary
condition of the style described in Haney \cite{Haney}.
Altogether, this yields the following forcing applied
in the model surface layer.

\begin{eqnarray}
\label{eq:eg-global-global_forcing}
\label{eq:eg-global-global_forcing_fu}
{\cal F}_{u} & = & \frac{\tau_{x}}{\rho_{0} \Delta z_{s}}
\\
\label{eq:eg-global-global_forcing_fv}
{\cal F}_{v} & = & \frac{\tau_{y}}{\rho_{0} \Delta z_{s}}
\\
\label{eq:eg-global-global_forcing_ft}
{\cal F}_{\theta} & = & - \lambda_{\theta} ( \theta - \theta^{\ast} ) 
 - \frac{1}{C_{p} \rho_{0} \Delta z_{s}}{\cal Q}
\\
\label{eq:eg-global-global_forcing_fs}
{\cal F}_{s} & = & - \lambda_{s} ( S - S^{\ast} ) 
 + \frac{S_{0}}{\Delta z_{s}}({\cal E} - {\cal P} - {\cal R})
\end{eqnarray}

\noindent where ${\cal F}_{u}$, ${\cal F}_{v}$, ${\cal F}_{\theta}$,
${\cal F}_{s}$ are the forcing terms in the zonal and meridional
momentum and in the potential temperature and salinity
equations respectively.
The term $\Delta z_{s}$ represents the top ocean layer thickness in
meters.
It is used in conjunction with a reference density, $\rho_{0}$
(here set to $999.8\,{\rm kg\,m^{-3}}$), a
reference salinity, $S_{0}$ (here set to 35~ppt),
and a specific heat capacity, $C_{p}$ (here set to
$4000~{\rm J}~^{\circ}{\rm C}^{-1}~{\rm kg}^{-1}$), to convert
input dataset values into time tendencies of
potential temperature (with units of $^{\circ}{\rm C}~{\rm s}^{-1}$),
salinity (with units ${\rm ppt}~s^{-1}$) and
velocity (with units ${\rm m}~{\rm s}^{-2}$).
The externally supplied forcing fields used in this
experiment are $\tau_{x}$, $\tau_{y}$, $\theta^{\ast}$, $S^{\ast}$,
$\cal{Q}$ and $\cal{E}-\cal{P}-\cal{R}$. The wind stress fields ($\tau_x$, $\tau_y$)
have units of ${\rm N}~{\rm m}^{-2}$. The temperature forcing fields
($\theta^{\ast}$ and $Q$) have units of $^{\circ}{\rm C}$ and ${\rm W}~{\rm m}^{-2}$
respectively. The salinity forcing fields ($S^{\ast}$ and 
$\cal{E}-\cal{P}-\cal{R}$) have units of ${\rm ppt}$ and ${\rm m}~{\rm s}^{-1}$
respectively. The source files and procedures for ingesting this data into the
simulation are described in the experiment configuration discussion in section
\ref{sec:eg-global-clim_ocn_examp_exp_config}.


\subsection{Discrete Numerical Configuration}
%\label{www:tutorials}


 The model is configured in hydrostatic form.  The domain is discretised with 
a uniform grid spacing in latitude and longitude on the sphere
 $\Delta \phi=\Delta \lambda=4^{\circ}$, so 
that there are ninety grid cells in the zonal and forty in the 
meridional direction. The internal model coordinate variables
$x$ and $y$ are initialized according to
\begin{eqnarray}
x=r\cos(\phi),~\Delta x & = &r\cos(\Delta \phi) \\
y=r\lambda,~\Delta y &= &r\Delta \lambda 
\end{eqnarray}

Arctic polar regions are not
included in this experiment. Meridionally the model extends from
$80^{\circ}{\rm S}$ to $80^{\circ}{\rm N}$.
Vertically the model is configured with twenty layers with the 
following thicknesses
$\Delta z_{1} = 50\,{\rm m},\,
 \Delta z_{2} = 50\,{\rm m},\,
 \Delta z_{3} = 55\,{\rm m},\,
 \Delta z_{4} = 60\,{\rm m},\,
 \Delta z_{5} = 65\,{\rm m},\,
$
$
 \Delta z_{6}~=~70\,{\rm m},\,
 \Delta z_{7}~=~80\,{\rm m},\,
 \Delta z_{8}~=95\,{\rm m},\,
 \Delta z_{9}=120\,{\rm m},\,
 \Delta z_{10}=155\,{\rm m},\,
$
$
 \Delta z_{11}=200\,{\rm m},\,
 \Delta z_{12}=260\,{\rm m},\,
 \Delta z_{13}=320\,{\rm m},\,
 \Delta z_{14}=400\,{\rm m},\,
 \Delta z_{15}=480\,{\rm m},\,
$
$
 \Delta z_{16}=570\,{\rm m},\,
 \Delta z_{17}=655\,{\rm m},\,
 \Delta z_{18}=725\,{\rm m},\,
 \Delta z_{19}=775\,{\rm m},\,
 \Delta z_{20}=815\,{\rm m}
$ (here the numeric subscript indicates the model level index number, ${\tt k}$) to
give a total depth, $H$, of $-5450{\rm m}$.
The implicit free surface form of the pressure equation described in Marshall et. al 
\cite{marshall:97a} is employed. A Laplacian operator, $\nabla^2$, provides viscous
dissipation. Thermal and haline diffusion is also represented by a Laplacian operator.

Wind-stress forcing is added to the momentum equations in (\ref{eq:eg-global-model_equations}) 
for both the zonal flow, $u$ and the meridional flow $v$, according to equations 
(\ref{eq:eg-global-global_forcing_fu}) and (\ref{eq:eg-global-global_forcing_fv}).
Thermodynamic forcing inputs are added to the equations 
in (\ref{eq:eg-global-model_equations}) for
potential temperature, $\theta$, and salinity, $S$, according to equations 
(\ref{eq:eg-global-global_forcing_ft}) and (\ref{eq:eg-global-global_forcing_fs}).
This produces a set of equations solved in this configuration as follows:

\begin{eqnarray}
\label{eq:eg-global-model_equations}
\frac{Du}{Dt} - fv + 
  \frac{1}{\rho}\frac{\partial p^{'}}{\partial x} - 
  \nabla_{h}\cdot A_{h}\nabla_{h}u - 
  \frac{\partial}{\partial z}A_{z}\frac{\partial u}{\partial z} 
 & = &
\begin{cases}
{\cal F}_u & \text{(surface)} \\
0 & \text{(interior)}
\end{cases}
\\
\frac{Dv}{Dt} + fu + 
  \frac{1}{\rho}\frac{\partial p^{'}}{\partial y} - 
  \nabla_{h}\cdot A_{h}\nabla_{h}v - 
  \frac{\partial}{\partial z}A_{z}\frac{\partial v}{\partial z} 
& = &
\begin{cases}
{\cal F}_v & \text{(surface)} \\
0 & \text{(interior)}
\end{cases}
\\
\frac{\partial \eta}{\partial t} + \nabla_{h}\cdot \vec{u}
&=&
0
\\
\frac{D\theta}{Dt} -
 \nabla_{h}\cdot K_{h}\nabla_{h}\theta
 - \frac{\partial}{\partial z}\Gamma(K_{z})\frac{\partial\theta}{\partial z} 
& = &
\begin{cases}
{\cal F}_\theta & \text{(surface)} \\
0 & \text{(interior)}
\end{cases}
\\
\frac{D s}{Dt} -
 \nabla_{h}\cdot K_{h}\nabla_{h}s
 - \frac{\partial}{\partial z}\Gamma(K_{z})\frac{\partial s}{\partial z} 
& = &
\begin{cases}
{\cal F}_s & \text{(surface)} \\
0 & \text{(interior)}
\end{cases}
\\
g\rho_{0} \eta + \int^{0}_{-z}\rho^{'} dz & = & p^{'}
\end{eqnarray}

\noindent where $u=\frac{Dx}{Dt}=r \cos(\phi)\frac{D \lambda}{Dt}$ and 
$v=\frac{Dy}{Dt}=r \frac{D \phi}{Dt}$ 
are the zonal and meridional components of the
flow vector, $\vec{u}$, on the sphere. As described in
MITgcm Numerical Solution Procedure \ref{chap:discretization}, the time 
evolution of potential temperature, $\theta$, equation is solved prognostically.
The total pressure, $p$, is diagnosed by summing pressure due to surface 
elevation $\eta$ and the hydrostatic pressure.
\\

\subsubsection{Numerical Stability Criteria}
%\label{www:tutorials}

The Laplacian dissipation coefficient, $A_{h}$, is set to $5 \times 10^5 m s^{-1}$.
This value is chosen to yield a Munk layer width \cite{adcroft:95},
\begin{eqnarray}
\label{eq:eg-global-munk_layer}
&& M_{w} = \pi ( \frac { A_{h} }{ \beta } )^{\frac{1}{3}}
\end{eqnarray}

\noindent  of $\approx 600$km. This is greater than the model
resolution in low-latitudes, $\Delta x \approx 400{\rm km}$, ensuring that the frictional 
boundary layer is adequately resolved.
\\

\noindent The model is stepped forward with a 
time step $\delta t_{\theta}=30~{\rm hours}$ for thermodynamic variables and
$\delta t_{v}=40~{\rm minutes}$ for momentum terms. With this time step, the stability 
parameter to the horizontal Laplacian friction \cite{adcroft:95}
\begin{eqnarray}
\label{eq:eg-global-laplacian_stability}
&& S_{l} = 4 \frac{A_{h} \delta t_{v}}{{\Delta x}^2}
\end{eqnarray}

\noindent evaluates to 0.16 at a latitude of $\phi=80^{\circ}$, which is below the 
0.3 upper limit for stability. The zonal grid spacing $\Delta x$ is smallest at
$\phi=80^{\circ}$ where $\Delta x=r\cos(\phi)\Delta \phi\approx 77{\rm km}$.
\\

\noindent The vertical dissipation coefficient, $A_{z}$, is set to 
$1\times10^{-3} {\rm m}^2{\rm s}^{-1}$. The associated stability limit
\begin{eqnarray}
\label{eq:eg-global-laplacian_stability_z}
S_{l} = 4 \frac{A_{z} \delta t_{v}}{{\Delta z}^2}
\end{eqnarray}

\noindent evaluates to $0.015$ for the smallest model
level spacing ($\Delta z_{1}=50{\rm m}$) which is again well below
the upper stability limit.
\\

The values of the horizontal ($K_{h}$) and vertical ($K_{z}$) diffusion coefficients 
for both temperature and salinity are set to $1 \times 10^{3}~{\rm m}^{2}{\rm s}^{-1}$ 
and $3 \times 10^{-5}~{\rm m}^{2}{\rm s}^{-1}$ respectively. The stability limit 
related to $K_{h}$ will be at $\phi=80^{\circ}$ where $\Delta x \approx 77 {\rm km}$. 
Here the stability parameter 
\begin{eqnarray} 
\label{eq:eg-global-laplacian_stability_xtheta}
S_{l} = \frac{4 K_{h} \delta t_{\theta}}{{\Delta x}^2} 
\end{eqnarray}
evaluates to $0.07$, well below the stability limit of $S_{l} \approx 0.5$. The 
stability parameter related to $K_{z}$
\begin{eqnarray} 
\label{eq:eg-global-laplacian_stability_ztheta}
S_{l} = \frac{4 K_{z} \delta t_{\theta}}{{\Delta z}^2} 
\end{eqnarray}
evaluates to $0.005$ for $\min(\Delta z)=50{\rm m}$, well below the stability limit 
of $S_{l} \approx 0.5$.
\\

\noindent The numerical stability for inertial oscillations
\cite{adcroft:95} 

\begin{eqnarray}
\label{eq:eg-global-inertial_stability}
S_{i} = f^{2} {\delta t_v}^2
\end{eqnarray}

\noindent evaluates to $0.24$ for $f=2\omega\sin(80^{\circ})=1.43\times10^{-4}~{\rm s}^{-1}$, which is close to 
the $S_{i} < 1$ upper limit for stability.
\\

\noindent The advective CFL \cite{adcroft:95} for a extreme maximum 
horizontal flow
speed of $ | \vec{u} | = 2 ms^{-1}$

\begin{eqnarray}
\label{eq:eg-global-cfl_stability}
S_{a} = \frac{| \vec{u} | \delta t_{v}}{ \Delta x}
\end{eqnarray}

\noindent evaluates to $6 \times 10^{-2}$. This is well below the stability 
limit of 0.5.
\\

\noindent The stability parameter for internal gravity waves propagating
with a maximum speed of $c_{g}=10~{\rm ms}^{-1}$
\cite{adcroft:95}

\begin{eqnarray}
\label{eq:eg-global-gfl_stability}
S_{c} = \frac{c_{g} \delta t_{v}}{ \Delta x}
\end{eqnarray}

\noindent evaluates to $3 \times 10^{-1}$. This is close to the linear
stability limit of 0.5.
  
\subsection{Experiment Configuration}
%\label{www:tutorials}
\label{sec:eg-global-clim_ocn_examp_exp_config}

The model configuration for this experiment resides under the 
directory {\it tutorial\_examples/global\_ocean\_circulation/}.  
The experiment files 

\begin{itemize}
\item {\it input/data}
\item {\it input/data.pkg}
\item {\it input/eedata},
\item {\it input/windx.bin},
\item {\it input/windy.bin},
\item {\it input/salt.bin},
\item {\it input/theta.bin},
\item {\it input/SSS.bin},
\item {\it input/SST.bin},
\item {\it input/topog.bin},
\item {\it code/CPP\_EEOPTIONS.h}
\item {\it code/CPP\_OPTIONS.h},
\item {\it code/SIZE.h}. 
\end{itemize}
contain the code customizations and parameter settings for these
experiments. Below we describe the customizations
to these files associated with this experiment.

\subsubsection{Driving Datasets}
%\label{www:tutorials}

Figures ({\it --- missing figures ---}) 
%(\ref{fig:sim_config_tclim}-\ref{fig:sim_config_empmr})
show the relaxation temperature ($\theta^{\ast}$) and salinity ($S^{\ast}$) 
fields, the wind stress components ($\tau_x$ and $\tau_y$), the heat flux ($Q$)
and the net fresh water flux (${\cal E} - {\cal P} - {\cal R}$) used
in equations 
(\ref{eq:eg-global-global_forcing_fu}-\ref{eq:eg-global-global_forcing_fs}).
The figures also indicate the lateral extent and coastline used in the 
experiment. Figure ({\it --- missing figure --- }) %ref{fig:model_bathymetry}) 
shows the depth contours of the model domain.

\subsubsection{File {\it input/data}}
%\label{www:tutorials}

\input{s_examples/global_oce_latlon/inp_data}

\subsubsection{File {\it input/data.pkg}}
%\label{www:tutorials}

This file uses standard default values and does not contain
customisations for this experiment.

\subsubsection{File {\it input/eedata}}
%\label{www:tutorials}

This file uses standard default values and does not contain
customisations for this experiment.

\subsubsection{File {\it input/windx.sin\_y}}
%\label{www:tutorials}

The {\it input/windx.sin\_y} file specifies a two-dimensional ($x,y$) 
map of wind stress ,$\tau_{x}$, values. The units used are $Nm^{-2}$.
Although $\tau_{x}$ is only a function of $y$n in this experiment
this file must still define a complete two-dimensional map in order
to be compatible with the standard code for loading forcing fields 
in MITgcm. The included matlab program {\it input/gendata.m} gives a complete
code for creating the {\it input/windx.sin\_y} file.

\subsubsection{File {\it input/topog.box}}
%\label{www:tutorials}


The {\it input/topog.box} file specifies a two-dimensional ($x,y$) 
map of depth values. For this experiment values are either
$0m$ or $-2000\,{\rm m}$, corresponding respectively to a wall or to deep
ocean. The file contains a raw binary stream of data that is enumerated
in the same way as standard MITgcm two-dimensional, horizontal arrays.
The included matlab program {\it input/gendata.m} gives a complete
code for creating the {\it input/topog.box} file.

\subsubsection{File {\it code/SIZE.h}}
%\label{www:tutorials}

Two lines are customized in this file for the current experiment

\begin{itemize}

\item Line 39, 
\begin{verbatim} sNx=60, \end{verbatim} this line sets
the lateral domain extent in grid points for the
axis aligned with the x-coordinate.

\item Line 40, 
\begin{verbatim} sNy=60, \end{verbatim} this line sets
the lateral domain extent in grid points for the
axis aligned with the y-coordinate.

\item Line 49, 
\begin{verbatim} Nr=4,   \end{verbatim} this line sets
the vertical domain extent in grid points.

\end{itemize}

\begin{small}
\input{s_examples/global_oce_latlon/code/SIZE.h}
\end{small}

\subsubsection{File {\it code/CPP\_OPTIONS.h}}
%\label{www:tutorials}

This file uses standard default values and does not contain
customisations for this experiment.


\subsubsection{File {\it code/CPP\_EEOPTIONS.h}}
%\label{www:tutorials}

This file uses standard default values and does not contain
customisations for this experiment.

\subsubsection{Other Files }
%\label{www:tutorials}

Other files relevant to this experiment are
\begin{itemize}
\item {\it model/src/ini\_cori.F}. This file initializes the model
coriolis variables {\bf fCorU}.
\item {\it model/src/ini\_spherical\_polar\_grid.F}
\item {\it model/src/ini\_parms.F},
\item {\it input/windx.sin\_y},
\end{itemize}
contain the code customisations and parameter settings for this 
experiments. Below we describe the customisations
to these files associated with this experiment.
1	jmc	1.20	% $Header: /u/gcmpack/manual/s_examples/global_oce_latlon/climatalogical_ogcm.tex,v 1.19 2010/08/30 23:09:20 jmc Exp $
2	cnh	1.2	% $Name: $
3	adcroft	1.1
4	jmc	1.17	\section[Global Ocean MITgcm Example]{Global Ocean Simulation at $4^\circ$ Resolution}
5	jmc	1.19	%\label{www:tutorials}
6			\label{sec:eg-global}
7	edhill	1.12	\begin{rawhtml}
8			<!-- CMIREDIR:eg-global: -->
9			\end{rawhtml}
10	jmc	1.16	\begin{center}
11			(in directory: {\it verification/tutorial\_global\_oce\_latlon/})
12			\end{center}
13	adcroft	1.1
14			\bodytext{bgcolor="#FFFFFFFF"}
15
16			%\begin{center}
17	cnh	1.3	%{\Large \bf Using MITgcm to Simulate Global Climatological Ocean Circulation
18	adcroft	1.1	%At Four Degree Resolution with Asynchronous Time Stepping}
19			%
20			%\vspace*{4mm}
21			%
22			%\vspace*{3mm}
23			%{\large May 2001}
24			%\end{center}
25
26
27			This example experiment demonstrates using the MITgcm to simulate
28			the planetary ocean circulation. The simulation is configured
29			with realistic geography and bathymetry on a
30			$4^{\circ} \times 4^{\circ}$ spherical polar grid.
31	molod	1.14	The files for this experiment are in the verification directory
32			under tutorial\_global\_oce\_latlon.
33	adcroft	1.1	Twenty levels are used in the vertical, ranging in thickness
34			from $50\,{\rm m}$ at the surface to $815\,{\rm m}$ at depth,
35			giving a maximum model depth of $6\,{\rm km}$.
36			At this resolution, the configuration
37			can be integrated forward for thousands of years on a single
38			processor desktop computer.
39			\\
40	cnh	1.8	\subsection{Overview}
41	jmc	1.19	%\label{www:tutorials}
42	adcroft	1.1
43	cnh	1.3	The model is forced with climatological wind stress data and surface
44			flux data from DaSilva \cite{DaSilva94}. Climatological data
45			from Levitus \cite{Levitus94} is used to initialize the model hydrography.
46			Levitus seasonal climatology data is also used throughout the calculation
47	adcroft	1.1	to provide additional air-sea fluxes.
48	cnh	1.3	These fluxes are combined with the DaSilva climatological estimates of
49	adcroft	1.1	surface heat flux and fresh water, resulting in a mixed boundary
50	cnh	1.3	condition of the style described in Haney \cite{Haney}.
51	adcroft	1.1	Altogether, this yields the following forcing applied
52			in the model surface layer.
53
54			\begin{eqnarray}
55	jmc	1.19	\label{eq:eg-global-global_forcing}
56			\label{eq:eg-global-global_forcing_fu}
57	adcroft	1.1	{\cal F}_{u} & = & \frac{\tau_{x}}{\rho_{0} \Delta z_{s}}
58			\\
59	jmc	1.19	\label{eq:eg-global-global_forcing_fv}
60	adcroft	1.1	{\cal F}_{v} & = & \frac{\tau_{y}}{\rho_{0} \Delta z_{s}}
61			\\
62	jmc	1.19	\label{eq:eg-global-global_forcing_ft}
63	adcroft	1.1	{\cal F}_{\theta} & = & - \lambda_{\theta} ( \theta - \theta^{\ast} )
64			- \frac{1}{C_{p} \rho_{0} \Delta z_{s}}{\cal Q}
65			\\
66	jmc	1.19	\label{eq:eg-global-global_forcing_fs}
67	adcroft	1.1	{\cal F}_{s} & = & - \lambda_{s} ( S - S^{\ast} )
68			+ \frac{S_{0}}{\Delta z_{s}}({\cal E} - {\cal P} - {\cal R})
69			\end{eqnarray}
70
71			\noindent where ${\cal F}_{u}$, ${\cal F}_{v}$, ${\cal F}_{\theta}$,
72			${\cal F}_{s}$ are the forcing terms in the zonal and meridional
73			momentum and in the potential temperature and salinity
74			equations respectively.
75			The term $\Delta z_{s}$ represents the top ocean layer thickness in
76			meters.
77			It is used in conjunction with a reference density, $\rho_{0}$
78			(here set to $999.8\,{\rm kg\,m^{-3}}$), a
79			reference salinity, $S_{0}$ (here set to 35~ppt),
80			and a specific heat capacity, $C_{p}$ (here set to
81			$4000~{\rm J}~^{\circ}{\rm C}^{-1}~{\rm kg}^{-1}$), to convert
82			input dataset values into time tendencies of
83			potential temperature (with units of $^{\circ}{\rm C}~{\rm s}^{-1}$),
84			salinity (with units ${\rm ppt}~s^{-1}$) and
85			velocity (with units ${\rm m}~{\rm s}^{-2}$).
86			The externally supplied forcing fields used in this
87			experiment are $\tau_{x}$, $\tau_{y}$, $\theta^{\ast}$, $S^{\ast}$,
88			$\cal{Q}$ and $\cal{E}-\cal{P}-\cal{R}$. The wind stress fields ($\tau_x$, $\tau_y$)
89			have units of ${\rm N}~{\rm m}^{-2}$. The temperature forcing fields
90			($\theta^{\ast}$ and $Q$) have units of $^{\circ}{\rm C}$ and ${\rm W}~{\rm m}^{-2}$
91			respectively. The salinity forcing fields ($S^{\ast}$ and
92			$\cal{E}-\cal{P}-\cal{R}$) have units of ${\rm ppt}$ and ${\rm m}~{\rm s}^{-1}$
93	cnh	1.8	respectively. The source files and procedures for ingesting this data into the
94			simulation are described in the experiment configuration discussion in section
95	jmc	1.19	\ref{sec:eg-global-clim_ocn_examp_exp_config}.
96	adcroft	1.1
97
98			\subsection{Discrete Numerical Configuration}
99	jmc	1.19	%\label{www:tutorials}
100	adcroft	1.1
101
102			The model is configured in hydrostatic form. The domain is discretised with
103			a uniform grid spacing in latitude and longitude on the sphere
104			$\Delta \phi=\Delta \lambda=4^{\circ}$, so
105			that there are ninety grid cells in the zonal and forty in the
106			meridional direction. The internal model coordinate variables
107	cnh	1.3	$x$ and $y$ are initialized according to
108	adcroft	1.1	\begin{eqnarray}
109			x=r\cos(\phi),~\Delta x & = &r\cos(\Delta \phi) \\
110	cnh	1.8	y=r\lambda,~\Delta y &= &r\Delta \lambda
111	adcroft	1.1	\end{eqnarray}
112
113			Arctic polar regions are not
114			included in this experiment. Meridionally the model extends from
115			$80^{\circ}{\rm S}$ to $80^{\circ}{\rm N}$.
116			Vertically the model is configured with twenty layers with the
117			following thicknesses
118			$\Delta z_{1} = 50\,{\rm m},\,
119			\Delta z_{2} = 50\,{\rm m},\,
120			\Delta z_{3} = 55\,{\rm m},\,
121			\Delta z_{4} = 60\,{\rm m},\,
122			\Delta z_{5} = 65\,{\rm m},\,
123			$
124			$
125			\Delta z_{6}~=~70\,{\rm m},\,
126			\Delta z_{7}~=~80\,{\rm m},\,
127			\Delta z_{8}~=95\,{\rm m},\,
128			\Delta z_{9}=120\,{\rm m},\,
129			\Delta z_{10}=155\,{\rm m},\,
130			$
131			$
132			\Delta z_{11}=200\,{\rm m},\,
133			\Delta z_{12}=260\,{\rm m},\,
134			\Delta z_{13}=320\,{\rm m},\,
135			\Delta z_{14}=400\,{\rm m},\,
136			\Delta z_{15}=480\,{\rm m},\,
137			$
138			$
139			\Delta z_{16}=570\,{\rm m},\,
140			\Delta z_{17}=655\,{\rm m},\,
141			\Delta z_{18}=725\,{\rm m},\,
142			\Delta z_{19}=775\,{\rm m},\,
143			\Delta z_{20}=815\,{\rm m}
144	cnh	1.8	$ (here the numeric subscript indicates the model level index number, ${\tt k}$) to
145			give a total depth, $H$, of $-5450{\rm m}$.
146	adcroft	1.1	The implicit free surface form of the pressure equation described in Marshall et. al
147	adcroft	1.6	\cite{marshall:97a} is employed. A Laplacian operator, $\nabla^2$, provides viscous
148	cnh	1.3	dissipation. Thermal and haline diffusion is also represented by a Laplacian operator.
149	adcroft	1.1
150	jmc	1.19	Wind-stress forcing is added to the momentum equations in (\ref{eq:eg-global-model_equations})
151	cnh	1.8	for both the zonal flow, $u$ and the meridional flow $v$, according to equations
152	jmc	1.19	(\ref{eq:eg-global-global_forcing_fu}) and (\ref{eq:eg-global-global_forcing_fv}).
153	cnh	1.8	Thermodynamic forcing inputs are added to the equations
154	jmc	1.19	in (\ref{eq:eg-global-model_equations}) for
155	adcroft	1.1	potential temperature, $\theta$, and salinity, $S$, according to equations
156	jmc	1.19	(\ref{eq:eg-global-global_forcing_ft}) and (\ref{eq:eg-global-global_forcing_fs}).
157	adcroft	1.1	This produces a set of equations solved in this configuration as follows:
158
159			\begin{eqnarray}
160	jmc	1.19	\label{eq:eg-global-model_equations}
161	adcroft	1.1	\frac{Du}{Dt} - fv +
162			\frac{1}{\rho}\frac{\partial p^{'}}{\partial x} -
163			\nabla_{h}\cdot A_{h}\nabla_{h}u -
164			\frac{\partial}{\partial z}A_{z}\frac{\partial u}{\partial z}
165			& = &
166			\begin{cases}
167			{\cal F}_u & \text{(surface)} \\
168			0 & \text{(interior)}
169			\end{cases}
170			\\
171			\frac{Dv}{Dt} + fu +
172			\frac{1}{\rho}\frac{\partial p^{'}}{\partial y} -
173			\nabla_{h}\cdot A_{h}\nabla_{h}v -
174			\frac{\partial}{\partial z}A_{z}\frac{\partial v}{\partial z}
175			& = &
176			\begin{cases}
177			{\cal F}_v & \text{(surface)} \\
178			0 & \text{(interior)}
179			\end{cases}
180			\\
181			\frac{\partial \eta}{\partial t} + \nabla_{h}\cdot \vec{u}
182			&=&
183			0
184			\\
185			\frac{D\theta}{Dt} -
186			\nabla_{h}\cdot K_{h}\nabla_{h}\theta
187			- \frac{\partial}{\partial z}\Gamma(K_{z})\frac{\partial\theta}{\partial z}
188			& = &
189			\begin{cases}
190			{\cal F}_\theta & \text{(surface)} \\
191			0 & \text{(interior)}
192			\end{cases}
193			\\
194			\frac{D s}{Dt} -
195			\nabla_{h}\cdot K_{h}\nabla_{h}s
196			- \frac{\partial}{\partial z}\Gamma(K_{z})\frac{\partial s}{\partial z}
197			& = &
198			\begin{cases}
199			{\cal F}_s & \text{(surface)} \\
200			0 & \text{(interior)}
201			\end{cases}
202			\\
203			g\rho_{0} \eta + \int^{0}_{-z}\rho^{'} dz & = & p^{'}
204			\end{eqnarray}
205
206			\noindent where $u=\frac{Dx}{Dt}=r \cos(\phi)\frac{D \lambda}{Dt}$ and
207			$v=\frac{Dy}{Dt}=r \frac{D \phi}{Dt}$
208			are the zonal and meridional components of the
209			flow vector, $\vec{u}$, on the sphere. As described in
210	adcroft	1.5	MITgcm Numerical Solution Procedure \ref{chap:discretization}, the time
211	adcroft	1.1	evolution of potential temperature, $\theta$, equation is solved prognostically.
212			The total pressure, $p$, is diagnosed by summing pressure due to surface
213			elevation $\eta$ and the hydrostatic pressure.
214			\\
215
216			\subsubsection{Numerical Stability Criteria}
217	jmc	1.19	%\label{www:tutorials}
218	adcroft	1.1
219	cnh	1.3	The Laplacian dissipation coefficient, $A_{h}$, is set to $5 \times 10^5 m s^{-1}$.
220	adcroft	1.4	This value is chosen to yield a Munk layer width \cite{adcroft:95},
221	adcroft	1.1	\begin{eqnarray}
222	jmc	1.19	\label{eq:eg-global-munk_layer}
223	adcroft	1.10	&& M_{w} = \pi ( \frac { A_{h} }{ \beta } )^{\frac{1}{3}}
224	adcroft	1.1	\end{eqnarray}
225
226			\noindent of $\approx 600$km. This is greater than the model
227			resolution in low-latitudes, $\Delta x \approx 400{\rm km}$, ensuring that the frictional
228			boundary layer is adequately resolved.
229			\\
230
231			\noindent The model is stepped forward with a
232			time step $\delta t_{\theta}=30~{\rm hours}$ for thermodynamic variables and
233			$\delta t_{v}=40~{\rm minutes}$ for momentum terms. With this time step, the stability
234	adcroft	1.4	parameter to the horizontal Laplacian friction \cite{adcroft:95}
235	adcroft	1.1	\begin{eqnarray}
236	jmc	1.19	\label{eq:eg-global-laplacian_stability}
237	adcroft	1.10	&& S_{l} = 4 \frac{A_{h} \delta t_{v}}{{\Delta x}^2}
238	adcroft	1.1	\end{eqnarray}
239
240			\noindent evaluates to 0.16 at a latitude of $\phi=80^{\circ}$, which is below the
241			0.3 upper limit for stability. The zonal grid spacing $\Delta x$ is smallest at
242			$\phi=80^{\circ}$ where $\Delta x=r\cos(\phi)\Delta \phi\approx 77{\rm km}$.
243			\\
244
245			\noindent The vertical dissipation coefficient, $A_{z}$, is set to
246			$1\times10^{-3} {\rm m}^2{\rm s}^{-1}$. The associated stability limit
247			\begin{eqnarray}
248	jmc	1.19	\label{eq:eg-global-laplacian_stability_z}
249	adcroft	1.1	S_{l} = 4 \frac{A_{z} \delta t_{v}}{{\Delta z}^2}
250			\end{eqnarray}
251
252			\noindent evaluates to $0.015$ for the smallest model
253	cnh	1.3	level spacing ($\Delta z_{1}=50{\rm m}$) which is again well below
254	adcroft	1.1	the upper stability limit.
255			\\
256
257			The values of the horizontal ($K_{h}$) and vertical ($K_{z}$) diffusion coefficients
258			for both temperature and salinity are set to $1 \times 10^{3}~{\rm m}^{2}{\rm s}^{-1}$
259			and $3 \times 10^{-5}~{\rm m}^{2}{\rm s}^{-1}$ respectively. The stability limit
260			related to $K_{h}$ will be at $\phi=80^{\circ}$ where $\Delta x \approx 77 {\rm km}$.
261			Here the stability parameter
262			\begin{eqnarray}
263	jmc	1.19	\label{eq:eg-global-laplacian_stability_xtheta}
264	adcroft	1.1	S_{l} = \frac{4 K_{h} \delta t_{\theta}}{{\Delta x}^2}
265			\end{eqnarray}
266	cnh	1.3	evaluates to $0.07$, well below the stability limit of $S_{l} \approx 0.5$. The
267	adcroft	1.1	stability parameter related to $K_{z}$
268			\begin{eqnarray}
269	jmc	1.19	\label{eq:eg-global-laplacian_stability_ztheta}
270	adcroft	1.1	S_{l} = \frac{4 K_{z} \delta t_{\theta}}{{\Delta z}^2}
271			\end{eqnarray}
272			evaluates to $0.005$ for $\min(\Delta z)=50{\rm m}$, well below the stability limit
273			of $S_{l} \approx 0.5$.
274			\\
275
276			\noindent The numerical stability for inertial oscillations
277	adcroft	1.4	\cite{adcroft:95}
278	adcroft	1.1
279			\begin{eqnarray}
280	jmc	1.19	\label{eq:eg-global-inertial_stability}
281	adcroft	1.1	S_{i} = f^{2} {\delta t_v}^2
282			\end{eqnarray}
283
284			\noindent evaluates to $0.24$ for $f=2\omega\sin(80^{\circ})=1.43\times10^{-4}~{\rm s}^{-1}$, which is close to
285			the $S_{i} < 1$ upper limit for stability.
286			\\
287
288	adcroft	1.4	\noindent The advective CFL \cite{adcroft:95} for a extreme maximum
289	adcroft	1.1	horizontal flow
290			speed of $ \| \vec{u} \| = 2 ms^{-1}$
291
292			\begin{eqnarray}
293	jmc	1.19	\label{eq:eg-global-cfl_stability}
294	adcroft	1.1	S_{a} = \frac{\| \vec{u} \| \delta t_{v}}{ \Delta x}
295			\end{eqnarray}
296
297			\noindent evaluates to $6 \times 10^{-2}$. This is well below the stability
298			limit of 0.5.
299			\\
300
301	cnh	1.3	\noindent The stability parameter for internal gravity waves propagating
302	adcroft	1.1	with a maximum speed of $c_{g}=10~{\rm ms}^{-1}$
303	adcroft	1.4	\cite{adcroft:95}
304	adcroft	1.1
305			\begin{eqnarray}
306	jmc	1.19	\label{eq:eg-global-gfl_stability}
307	adcroft	1.1	S_{c} = \frac{c_{g} \delta t_{v}}{ \Delta x}
308			\end{eqnarray}
309
310			\noindent evaluates to $3 \times 10^{-1}$. This is close to the linear
311			stability limit of 0.5.
312
313			\subsection{Experiment Configuration}
314	jmc	1.19	%\label{www:tutorials}
315			\label{sec:eg-global-clim_ocn_examp_exp_config}
316	adcroft	1.1
317			The model configuration for this experiment resides under the
318	cnh	1.8	directory {\it tutorial\_examples/global\_ocean\_circulation/}.
319			The experiment files
320
321	adcroft	1.1	\begin{itemize}
322			\item {\it input/data}
323			\item {\it input/data.pkg}
324			\item {\it input/eedata},
325			\item {\it input/windx.bin},
326			\item {\it input/windy.bin},
327			\item {\it input/salt.bin},
328			\item {\it input/theta.bin},
329			\item {\it input/SSS.bin},
330			\item {\it input/SST.bin},
331			\item {\it input/topog.bin},
332			\item {\it code/CPP\_EEOPTIONS.h}
333			\item {\it code/CPP\_OPTIONS.h},
334			\item {\it code/SIZE.h}.
335			\end{itemize}
336	cnh	1.3	contain the code customizations and parameter settings for these
337			experiments. Below we describe the customizations
338	adcroft	1.1	to these files associated with this experiment.
339	cnh	1.8
340			\subsubsection{Driving Datasets}
341	jmc	1.19	%\label{www:tutorials}
342	cnh	1.8
343	jmc	1.19	Figures ({\it --- missing figures ---})
344			%(\ref{fig:sim_config_tclim}-\ref{fig:sim_config_empmr})
345			show the relaxation temperature ($\theta^{\ast}$) and salinity ($S^{\ast}$)
346			fields, the wind stress components ($\tau_x$ and $\tau_y$), the heat flux ($Q$)
347	cnh	1.8	and the net fresh water flux (${\cal E} - {\cal P} - {\cal R}$) used
348	jmc	1.19	in equations
349			(\ref{eq:eg-global-global_forcing_fu}-\ref{eq:eg-global-global_forcing_fs}).
350			The figures also indicate the lateral extent and coastline used in the
351			experiment. Figure ({\it --- missing figure --- }) %ref{fig:model_bathymetry})
352			shows the depth contours of the model domain.
353	adcroft	1.1
354			\subsubsection{File {\it input/data}}
355	jmc	1.19	%\label{www:tutorials}
356	adcroft	1.1
357	jmc	1.20	\input{s_examples/global_oce_latlon/inp_data}
358	adcroft	1.1
359			\subsubsection{File {\it input/data.pkg}}
360	jmc	1.19	%\label{www:tutorials}
361	adcroft	1.1
362			This file uses standard default values and does not contain
363			customisations for this experiment.
364
365			\subsubsection{File {\it input/eedata}}
366	jmc	1.19	%\label{www:tutorials}
367	adcroft	1.1
368			This file uses standard default values and does not contain
369			customisations for this experiment.
370
371			\subsubsection{File {\it input/windx.sin\_y}}
372	jmc	1.19	%\label{www:tutorials}
373	adcroft	1.1
374			The {\it input/windx.sin\_y} file specifies a two-dimensional ($x,y$)
375			map of wind stress ,$\tau_{x}$, values. The units used are $Nm^{-2}$.
376			Although $\tau_{x}$ is only a function of $y$n in this experiment
377			this file must still define a complete two-dimensional map in order
378			to be compatible with the standard code for loading forcing fields
379			in MITgcm. The included matlab program {\it input/gendata.m} gives a complete
380			code for creating the {\it input/windx.sin\_y} file.
381
382			\subsubsection{File {\it input/topog.box}}
383	jmc	1.19	%\label{www:tutorials}
384	adcroft	1.1
385
386			The {\it input/topog.box} file specifies a two-dimensional ($x,y$)
387			map of depth values. For this experiment values are either
388			$0m$ or $-2000\,{\rm m}$, corresponding respectively to a wall or to deep
389			ocean. The file contains a raw binary stream of data that is enumerated
390			in the same way as standard MITgcm two-dimensional, horizontal arrays.
391			The included matlab program {\it input/gendata.m} gives a complete
392			code for creating the {\it input/topog.box} file.
393
394			\subsubsection{File {\it code/SIZE.h}}
395	jmc	1.19	%\label{www:tutorials}
396	adcroft	1.1
397			Two lines are customized in this file for the current experiment
398
399			\begin{itemize}
400
401			\item Line 39,
402			\begin{verbatim} sNx=60, \end{verbatim} this line sets
403			the lateral domain extent in grid points for the
404			axis aligned with the x-coordinate.
405
406			\item Line 40,
407			\begin{verbatim} sNy=60, \end{verbatim} this line sets
408			the lateral domain extent in grid points for the
409			axis aligned with the y-coordinate.
410
411			\item Line 49,
412			\begin{verbatim} Nr=4, \end{verbatim} this line sets
413			the vertical domain extent in grid points.
414
415			\end{itemize}
416
417			\begin{small}
418	jmc	1.18	\input{s_examples/global_oce_latlon/code/SIZE.h}
419	adcroft	1.1	\end{small}
420
421			\subsubsection{File {\it code/CPP\_OPTIONS.h}}
422	jmc	1.19	%\label{www:tutorials}
423	adcroft	1.1
424			This file uses standard default values and does not contain
425			customisations for this experiment.
426
427
428			\subsubsection{File {\it code/CPP\_EEOPTIONS.h}}
429	jmc	1.19	%\label{www:tutorials}
430	adcroft	1.1
431			This file uses standard default values and does not contain
432			customisations for this experiment.
433
434			\subsubsection{Other Files }
435	jmc	1.19	%\label{www:tutorials}
436	adcroft	1.1
437			Other files relevant to this experiment are
438			\begin{itemize}
439			\item {\it model/src/ini\_cori.F}. This file initializes the model
440			coriolis variables {\bf fCorU}.
441			\item {\it model/src/ini\_spherical\_polar\_grid.F}
442			\item {\it model/src/ini\_parms.F},
443			\item {\it input/windx.sin\_y},
444			\end{itemize}
445			contain the code customisations and parameter settings for this
446	cnh	1.3	experiments. Below we describe the customisations
447	adcroft	1.1	to these files associated with this experiment.