--- manual/s_examples/baroclinic_gyre/fourlayer.tex 2001/10/24 23:14:44 1.5 +++ manual/s_examples/baroclinic_gyre/fourlayer.tex 2002/05/16 15:54:37 1.13 @@ -1,8 +1,9 @@ -% $Header: /home/ubuntu/mnt/e9_copy/manual/s_examples/baroclinic_gyre/fourlayer.tex,v 1.5 2001/10/24 23:14:44 cnh Exp $ +% $Header: /home/ubuntu/mnt/e9_copy/manual/s_examples/baroclinic_gyre/fourlayer.tex,v 1.13 2002/05/16 15:54:37 adcroft Exp $ % $Name: $ -\section{Example: Four layer Baroclinic Ocean Gyre In Spherical Coordinates} -\label{sec:eg-fourlayer} +\section{Four Layer Baroclinic Ocean Gyre In Spherical Coordinates} +\label{www:tutorials} +\label{sect:eg-fourlayer} \bodytext{bgcolor="#FFFFFFFF"} @@ -19,17 +20,18 @@ This document describes an example experiment using MITgcm to simulate a baroclinic ocean gyre in spherical polar coordinates. The barotropic -example experiment in section \ref{sec:eg-baro} -ilustrated how to configure the code for a single layer -simulation in a cartesian grid. In this example a similar physical problem +example experiment in section \ref{sect:eg-baro} +illustrated how to configure the code for a single layer +simulation in a Cartesian grid. In this example a similar physical problem is simulated, but the code is now configured for four layers and in a spherical polar coordinate system. \subsection{Overview} +\label{www:tutorials} This example experiment demonstrates using the MITgcm to simulate a baroclinic, wind-forced, ocean gyre circulation. The experiment -is a numerical rendition of the gyre circulation problem simliar +is a numerical rendition of the gyre circulation problem similar to the problems described analytically by Stommel in 1966 \cite{Stommel66} and numerically in Holland et. al \cite{Holland75}. \\ @@ -43,7 +45,7 @@ according to latitude, $\varphi$ \begin{equation} -\label{EQ:fcori} +\label{EQ:eg-fourlayer-fcori} f(\varphi) = 2 \Omega \sin( \varphi ) \end{equation} @@ -61,14 +63,14 @@ $\tau_0$ is set to $0.1N m^{-2}$. \\ -Figure \ref{FIG:simulation_config} -summarises the configuration simulated. -In contrast to the example in section \ref{sec:eg-baro}, the +Figure \ref{FIG:eg-fourlayer-simulation_config} +summarizes the configuration simulated. +In contrast to the example in section \ref{sect:eg-baro}, the current experiment simulates a spherical polar domain. As indicated by the axes in the lower left of the figure the model code works internally -in a locally orthoganal coordinate $(x,y,z)$. For this experiment description -of this document the local orthogonal model coordinate $(x,y,z)$ is synonomous -with the spherical polar coordinate shown in figure +in a locally orthogonal coordinate $(x,y,z)$. For this experiment description +the local orthogonal model coordinate $(x,y,z)$ is synonymous +with the coordinates $(\lambda,\varphi,r)$ shown in figure \ref{fig:spherical-polar-coord} \\ @@ -82,14 +84,14 @@ linear \begin{equation} -\label{EQ:linear1_eos} +\label{EQ:eg-fourlayer-linear1_eos} \rho = \rho_{0} ( 1 - \alpha_{\theta}\theta^{'} ) \end{equation} \noindent which is implemented in the model as a density anomaly equation \begin{equation} -\label{EQ:linear1_eos_pert} +\label{EQ:eg-fourlayer-linear1_eos_pert} \rho^{'} = -\rho_{0}\alpha_{\theta}\theta^{'} \end{equation} @@ -114,26 +116,27 @@ imposed by setting the potential temperature, $\theta$, in each layer. The vertical spacing, $\Delta z$, is constant and equal to $500$m. } -\label{FIG:simulation_config} +\label{FIG:eg-fourlayer-simulation_config} \end{figure} \subsection{Equations solved} - -The implicit free surface {\bf HPE} form of the -equations described in Marshall et. al \cite{Marshall97a} is +\label{www:tutorials} +For this problem +the implicit free surface, {\bf HPE} (see section \ref{sect:hydrostatic_and_quasi-hydrostatic_forms}) form of the +equations described in Marshall et. al \cite{marshall:97a} are employed. The flow is three-dimensional with just temperature, $\theta$, as an active tracer. The equation of state is linear. -A horizontal laplacian operator $\nabla_{h}^2$ provides viscous +A horizontal Laplacian operator $\nabla_{h}^2$ provides viscous dissipation and provides a diffusive sub-grid scale closure for the temperature equation. A wind-stress momentum forcing is added to the momentum equation for the zonal flow, $u$. Other terms in the model -are explicitly switched off for this experiement configuration (see section +are explicitly switched off for this experiment configuration (see section \ref{SEC:eg_fourl_code_config} ). This yields an active set of equations solved in this configuration, written in spherical polar coordinates as follows \begin{eqnarray} -\label{EQ:model_equations} +\label{EQ:eg-fourlayer-model_equations} \frac{Du}{Dt} - fv + \frac{1}{\rho}\frac{\partial p^{\prime}}{\partial \lambda} - A_{h}\nabla_{h}^2u - A_{z}\frac{\partial^{2}u}{\partial z^{2}} @@ -185,7 +188,7 @@ part due to variations in density, $\rho^{\prime}$, integrated through the water column. -The suffices ${s},{i}$ indicate surface and interior of the domain. +The suffices ${s},{i}$ indicate surface layer and the interior of the domain. The windstress forcing, ${\cal F}_{\lambda}$, is applied in the surface layer by a source term in the zonal momentum equation. In the ocean interior this term is zero. @@ -202,6 +205,7 @@ \subsection{Discrete Numerical Configuration} +\label{www:tutorials} The domain is discretised with a uniform grid spacing in latitude and longitude @@ -210,7 +214,7 @@ Vertically the model is configured with four layers with constant depth, $\Delta z$, of $500$~m. The internal, locally orthogonal, model coordinate -variables $x$ and $y$ are initialised from the values of +variables $x$ and $y$ are initialized from the values of $\lambda$, $\varphi$, $\Delta \lambda$ and $\Delta \varphi$ in radians according to @@ -221,7 +225,7 @@ The procedure for generating a set of internal grid variables from a spherical polar grid specification is discussed in section -\ref{sec:spatial_discrete_horizontal_grid}. +\ref{sect:spatial_discrete_horizontal_grid}. \noindent\fbox{ \begin{minipage}{5.5in} {\em S/R INI\_SPHERICAL\_POLAR\_GRID} ({\em @@ -242,52 +246,60 @@ -As described in \ref{sec:tracer_equations}, the time evolution of potential +As described in \ref{sect:tracer_equations}, the time evolution of potential temperature, $\theta$, (equation \ref{eq:eg_fourl_theta}) is evaluated prognostically. The centered second-order scheme with Adams-Bashforth time stepping described in section -\ref{sec:tracer_equations_abII} is used to step forward the temperature -equation. The pressure forces that drive the fluid motions, ( +\ref{sect:tracer_equations_abII} is used to step forward the temperature +equation. Prognostic terms in +the momentum equations are solved using flux form as +described in section \ref{sect:flux-form_momentum_eqautions}. +The pressure forces that drive the fluid motions, ( $\frac{\partial p^{'}}{\partial \lambda}$ and $\frac{\partial p^{'}}{\partial \varphi}$), are found by summing pressure due to surface elevation $\eta$ and the hydrostatic pressure. The hydrostatic part of the -pressure is evaluated explicitly by integrating density. The sea-surface -height, $\eta$, is solved for implicitly as described in section -\ref{sect:pressure-method-linear-backward}. +pressure is diagnosed explicitly by integrating density. The sea-surface +height, $\eta$, is diagnosed using an implicit scheme. The pressure +field solution method is described in sections +\ref{sect:pressure-method-linear-backward} and +\ref{sect:finding_the_pressure_field}. \subsubsection{Numerical Stability Criteria} +\label{www:tutorials} -The laplacian dissipation coefficient, $A_{h}$, is set to $400 m s^{-1}$. -This value is chosen to yield a Munk layer width \cite{Adcroft_thesis}, +The Laplacian viscosity coefficient, $A_{h}$, is set to $400 m s^{-1}$. +This value is chosen to yield a Munk layer width, \begin{eqnarray} -\label{EQ:munk_layer} +\label{EQ:eg-fourlayer-munk_layer} M_{w} = \pi ( \frac { A_{h} }{ \beta } )^{\frac{1}{3}} \end{eqnarray} \noindent of $\approx 100$km. This is greater than the model -resolution in mid-latitudes $\Delta x$, ensuring that the frictional +resolution in mid-latitudes +$\Delta x=r \cos(\varphi) \Delta \lambda \approx 80~{\rm km}$ at +$\varphi=45^{\circ}$, ensuring that the frictional boundary layer is well resolved. \\ \noindent The model is stepped forward with a time step $\delta t=1200$secs. With this time step the stability -parameter to the horizontal laplacian friction \cite{Adcroft_thesis} +parameter to the horizontal Laplacian friction \begin{eqnarray} -\label{EQ:laplacian_stability} +\label{EQ:eg-fourlayer-laplacian_stability} S_{l} = 4 \frac{A_{h} \delta t}{{\Delta x}^2} \end{eqnarray} \noindent evaluates to 0.012, which is well below the 0.3 upper limit -for stability. +for stability for this term under ABII time-stepping. \\ \noindent The vertical dissipation coefficient, $A_{z}$, is set to $1\times10^{-2} {\rm m}^2{\rm s}^{-1}$. The associated stability limit \begin{eqnarray} -\label{EQ:laplacian_stability_z} +\label{EQ:eg-fourlayer-laplacian_stability_z} S_{l} = 4 \frac{A_{z} \delta t}{{\Delta z}^2} \end{eqnarray} @@ -298,10 +310,9 @@ \\ \noindent The numerical stability for inertial oscillations -\cite{Adcroft_thesis} \begin{eqnarray} -\label{EQ:inertial_stability} +\label{EQ:eg-fourlayer-inertial_stability} S_{i} = f^{2} {\delta t}^2 \end{eqnarray} @@ -309,35 +320,36 @@ limit for stability. \\ -\noindent The advective CFL \cite{Adcroft_thesis} for a extreme maximum +\noindent The advective CFL for a extreme maximum horizontal flow speed of $ | \vec{u} | = 2 ms^{-1}$ \begin{eqnarray} -\label{EQ:cfl_stability} -S_{a} = \frac{| \vec{u} | \delta t}{ \Delta x} +\label{EQ:eg-fourlayer-cfl_stability} +C_{a} = \frac{| \vec{u} | \delta t}{ \Delta x} \end{eqnarray} \noindent evaluates to $5 \times 10^{-2}$. This is well below the stability limit of 0.5. \\ -\noindent The stability parameter for internal gravity waves -\cite{Adcroft_thesis} +\noindent The stability parameter for internal gravity waves +propagating at $2~{\rm m}~{\rm s}^{-1}$ \begin{eqnarray} -\label{EQ:igw_stability} +\label{EQ:eg-fourlayer-igw_stability} S_{c} = \frac{c_{g} \delta t}{ \Delta x} \end{eqnarray} -\noindent evaluates to $5 \times 10^{-2}$. This is well below the linear +\noindent evaluates to $\approx 5 \times 10^{-2}$. This is well below the linear stability limit of 0.25. \subsection{Code Configuration} +\label{www:tutorials} \label{SEC:eg_fourl_code_config} The model configuration for this experiment resides under the -directory {\it verification/exp1/}. The experiment files +directory {\it verification/exp2/}. The experiment files \begin{itemize} \item {\it input/data} \item {\it input/data.pkg} @@ -349,10 +361,11 @@ \item {\it code/SIZE.h}. \end{itemize} contain the code customisations and parameter settings for this -experiements. Below we describe the customisations +experiments. Below we describe the customisations to these files associated with this experiment. \subsubsection{File {\it input/data}} +\label{www:tutorials} This file, reproduced completely below, specifies the main parameters for the experiment. The parameters that are significant for this configuration @@ -366,7 +379,7 @@ the initial and reference values of potential temperature at each model level in units of $^{\circ}$C. The entries are ordered from surface to depth. For each -depth level the inital and reference profiles will be uniform in +depth level the initial and reference profiles will be uniform in $x$ and $y$. The values specified here are read into the variable {\bf @@ -412,7 +425,7 @@ \item Line 6, \begin{verbatim} viscAz=1.E-2, \end{verbatim} -this line sets the vertical laplacian dissipation coefficient to +this line sets the vertical Laplacian dissipation coefficient to $1 \times 10^{-2} {\rm m^{2}s^{-1}}$. Boundary conditions for this operator are specified later. The variable @@ -432,7 +445,9 @@ \begin{rawhtml} \end{rawhtml} viscAr \begin{rawhtml} \end{rawhtml} -}. +}. At each time step, the viscous term contribution to the momentum equations +is calculated in routine +{\it S/R CALC\_DIFFUSIVITY}. \fbox{ \begin{minipage}{5.0in} @@ -463,7 +478,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +} and applied in routines {\it CALC\_MOM\_RHS} and {\it CALC\_GW}. \fbox{ \begin{minipage}{5.0in} @@ -506,7 +521,8 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +} and the boundary condition is evaluated in routine +{\it S/R CALC\_MOM\_RHS}. \fbox{ @@ -538,7 +554,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +} and is applied in the routine {\it S/R CALC\_MOM\_RHS}. \fbox{ \begin{minipage}{5.0in} @@ -570,7 +586,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +} and used in routine {\it S/R CALC\_GT}. \fbox{ \begin{minipage}{5.0in} {\it S/R CALC\_GT}({\it calc\_gt.F}) @@ -606,7 +622,7 @@ \begin{rawhtml} \end{rawhtml} diffKrT \begin{rawhtml} \end{rawhtml} -}. +} which is used in routine {\it S/R CALC\_DIFFUSIVITY}. \fbox{ \begin{minipage}{5.0in} {\it S/R CALC\_DIFFUSIVITY}({\it calc\_diffusivity.F}) @@ -637,7 +653,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +}. The routine {\it S/R FIND\_RHO} makes use of {\bf tAlpha}. \fbox{ \begin{minipage}{5.0in} @@ -666,7 +682,8 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +}. The values of {\bf eosType} sets which formula in routine +{\it FIND\_RHO} is used to calculate density. \fbox{ \begin{minipage}{5.0in} @@ -687,8 +704,8 @@ \end{verbatim} This line requests that the simulation be performed in a spherical polar coordinate system. It affects the interpretation of -grid inoput parameters, for exampl {\bf delX} and {\bf delY} and -causes the grid generation routines to initialise an internal grid based +grid input parameters, for example {\bf delX} and {\bf delY} and +causes the grid generation routines to initialize an internal grid based on spherical polar geometry. The variable {\bf @@ -701,7 +718,9 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +}. When set to {\bf .TRUE.} the settings of {\bf delX} and {\bf delY} are +taken to be in degrees. These values are used in the +routine {\it INI\_SPEHRICAL\_POLAR\_GRID}. \fbox{ \begin{minipage}{5.0in} @@ -721,7 +740,7 @@ This line sets the southern boundary of the modeled domain to $0^{\circ}$ latitude. This value affects both the generation of the locally orthogonal grid that the model -uses internally and affects the initialisation of the coriolis force. +uses internally and affects the initialization of the coriolis force. Note - it is not required to set a longitude boundary, since the absolute longitude does not alter the kernel equation discretisation. @@ -736,7 +755,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +} and is used in routine {\it INI\_SPEHRICAL\_POLAR\_GRID}. \fbox{ \begin{minipage}{5.0in} @@ -766,7 +785,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +} and is used in routine {\it INI\_SPEHRICAL\_POLAR\_GRID}. \fbox{ \begin{minipage}{5.0in} @@ -796,7 +815,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +} and is used in routine {\it INI\_SPEHRICAL\_POLAR\_GRID}. \fbox{ \begin{minipage}{5.0in} @@ -834,7 +853,7 @@ \begin{rawhtml} \end{rawhtml} delR \begin{rawhtml} \end{rawhtml} -}. +} which is used in routine {\it INI\_VERTICAL\_GRID}. \fbox{ \begin{minipage}{5.0in} @@ -873,7 +892,7 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +}. The bathymetry file is read in the routine {\it INI\_DEPTHS}. \fbox{ \begin{minipage}{5.0in} @@ -892,7 +911,7 @@ zonalWindFile='windx.sin_y' \end{verbatim} This line specifies the name of the file from which the x-direction -surface wind stress is read. This file is also a two-dimensional +(zonal) surface wind stress is read. This file is also a two-dimensional ($x,y$) map and is enumerated and formatted in the same manner as the bathymetry file. The matlab program {\it input/gendata.m} includes example code to generate a valid @@ -909,7 +928,8 @@ \begin{rawhtml} \end{rawhtml} INI\_PARMS \begin{rawhtml} \end{rawhtml} -}. +}. The wind-stress file is read in the routine +{\it EXTERNAL\_FIELDS\_LOAD}. \fbox{ \begin{minipage}{5.0in} @@ -924,9 +944,7 @@ \end{itemize} -\noindent other lines in the file {\it input/data} are standard values -that are described in the MITgcm Getting Started and MITgcm Parameters -notes. +\noindent other lines in the file {\it input/data} are standard values. \begin{rawhtml}
\end{rawhtml}
 \begin{small}
@@ -935,37 +953,45 @@
 \begin{rawhtml}
\end{rawhtml} \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 +map of wind stress ,$\tau_{x}$, values. The units used are $Nm^{-2}$ (the +default for MITgcm). +Although $\tau_{x}$ is only a function of latitude, $y$, +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 +in MITgcm (routine {\it EXTERNAL\_FIELDS\_LOAD}. +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 +$0~{\rm m}$ 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 @@ -992,17 +1018,20 @@ \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} @@ -1015,22 +1044,26 @@ \end{itemize} \subsection{Running The Example} +\label{www:tutorials} \label{SEC:running_the_example} \subsubsection{Code Download} +\label{www:tutorials} In order to run the examples you must first download the code distribution. -Instructions for downloading the code can be found in the Getting Started -Guide \cite{MITgcm_Getting_Started}. +Instructions for downloading the code can be found in section +\ref{sect:obtainingCode}. \subsubsection{Experiment Location} +\label{www:tutorials} This example experiments is located under the release sub-directory \vspace{5mm} -{\it verification/exp1/ } +{\it verification/exp2/ } \subsubsection{Running the Experiment} +\label{www:tutorials} To run the experiment @@ -1047,9 +1080,9 @@ % pwd \end{verbatim} - You shold see a response on the screen ending in + You should see a response on the screen ending in -{\it verification/exp1/input } +{\it verification/exp2/input } \item Run the genmake script to create the experiment {\it Makefile}