/[MITgcm]/manual/s_overview/text/manual.tex
ViewVC logotype

Annotation of /manual/s_overview/text/manual.tex

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


Revision 1.22 - (hide annotations) (download) (as text)
Sun Oct 17 04:15:13 2004 UTC (20 years, 8 months ago) by jmc
Branch: MAIN
Changes since 1.21: +8 -7 lines
File MIME type: application/x-tex
add a labbel for reference in part2.

1 jmc 1.22 % $Header: /u/gcmpack/manual/part1/manual.tex,v 1.21 2004/10/16 03:40:12 edhill Exp $
2 cnh 1.2 % $Name: $
3 cnh 1.1
4 adcroft 1.4 %tci%\documentclass[12pt]{book}
5     %tci%\usepackage{amsmath}
6     %tci%\usepackage{html}
7     %tci%\usepackage{epsfig}
8     %tci%\usepackage{graphics,subfigure}
9     %tci%\usepackage{array}
10     %tci%\usepackage{multirow}
11     %tci%\usepackage{fancyhdr}
12     %tci%\usepackage{psfrag}
13    
14     %tci%%TCIDATA{OutputFilter=Latex.dll}
15     %tci%%TCIDATA{LastRevised=Thursday, October 04, 2001 14:41:22}
16     %tci%%TCIDATA{<META NAME="GraphicsSave" CONTENT="32">}
17     %tci%%TCIDATA{Language=American English}
18    
19     %tci%\fancyhead{}
20     %tci%\fancyhead[LO]{\slshape \rightmark}
21     %tci%\fancyhead[RE]{\slshape \leftmark}
22     %tci%\fancyhead[RO,LE]{\thepage}
23     %tci%\fancyfoot[CO,CE]{\today}
24     %tci%\fancyfoot[RO,LE]{ }
25     %tci%\renewcommand{\headrulewidth}{0.4pt}
26     %tci%\renewcommand{\footrulewidth}{0.4pt}
27     %tci%\setcounter{secnumdepth}{3}
28     %tci%\input{tcilatex}
29    
30     %tci%\begin{document}
31    
32     %tci%\tableofcontents
33    
34    
35 cnh 1.1 % Section: Overview
36    
37 jmc 1.22 % $Header: /u/gcmpack/manual/part1/manual.tex,v 1.21 2004/10/16 03:40:12 edhill Exp $
38 cnh 1.1 % $Name: $
39    
40 cnh 1.16 This document provides the reader with the information necessary to
41 cnh 1.1 carry out numerical experiments using MITgcm. It gives a comprehensive
42     description of the continuous equations on which the model is based, the
43     numerical algorithms the model employs and a description of the associated
44     program code. Along with the hydrodynamical kernel, physical and
45     biogeochemical parameterizations of key atmospheric and oceanic processes
46     are available. A number of examples illustrating the use of the model in
47     both process and general circulation studies of the atmosphere and ocean are
48     also presented.
49    
50 cnh 1.16 \section{Introduction}
51 afe 1.18 \begin{rawhtml}
52 afe 1.19 <!-- CMIREDIR:innovations: -->
53 afe 1.18 \end{rawhtml}
54    
55 cnh 1.16
56 cnh 1.1 MITgcm has a number of novel aspects:
57    
58     \begin{itemize}
59     \item it can be used to study both atmospheric and oceanic phenomena; one
60     hydrodynamical kernel is used to drive forward both atmospheric and oceanic
61 cnh 1.7 models - see fig \ref{fig:onemodel}
62 cnh 1.1
63 cnh 1.3 %% CNHbegin
64     \input{part1/one_model_figure}
65     %% CNHend
66    
67 cnh 1.1 \item it has a non-hydrostatic capability and so can be used to study both
68 cnh 1.7 small-scale and large scale processes - see fig \ref{fig:all-scales}
69 cnh 1.1
70 cnh 1.3 %% CNHbegin
71     \input{part1/all_scales_figure}
72     %% CNHend
73    
74 cnh 1.1 \item finite volume techniques are employed yielding an intuitive
75     discretization and support for the treatment of irregular geometries using
76 cnh 1.7 orthogonal curvilinear grids and shaved cells - see fig \ref{fig:finite-volumes}
77 cnh 1.3
78     %% CNHbegin
79     \input{part1/fvol_figure}
80     %% CNHend
81 cnh 1.1
82     \item tangent linear and adjoint counterparts are automatically maintained
83     along with the forward model, permitting sensitivity and optimization
84     studies.
85    
86     \item the model is developed to perform efficiently on a wide variety of
87     computational platforms.
88     \end{itemize}
89    
90 cnh 1.16 Key publications reporting on and charting the development of the model are
91     \cite{hill:95,marshall:97a,marshall:97b,adcroft:97,marshall:98,adcroft:99,hill:99,maro-eta:99}:
92 cnh 1.12
93     \begin{verbatim}
94     Hill, C. and J. Marshall, (1995)
95     Application of a Parallel Navier-Stokes Model to Ocean Circulation in
96     Parallel Computational Fluid Dynamics
97     In Proceedings of Parallel Computational Fluid Dynamics: Implementations
98     and Results Using Parallel Computers, 545-552.
99     Elsevier Science B.V.: New York
100    
101     Marshall, J., C. Hill, L. Perelman, and A. Adcroft, (1997)
102 cnh 1.16 Hydrostatic, quasi-hydrostatic, and nonhydrostatic ocean modeling
103 cnh 1.12 J. Geophysical Res., 102(C3), 5733-5752.
104    
105     Marshall, J., A. Adcroft, C. Hill, L. Perelman, and C. Heisey, (1997)
106     A finite-volume, incompressible Navier Stokes model for studies of the ocean
107     on parallel computers,
108     J. Geophysical Res., 102(C3), 5753-5766.
109    
110     Adcroft, A.J., Hill, C.N. and J. Marshall, (1997)
111     Representation of topography by shaved cells in a height coordinate ocean
112     model
113     Mon Wea Rev, vol 125, 2293-2315
114    
115     Marshall, J., Jones, H. and C. Hill, (1998)
116     Efficient ocean modeling using non-hydrostatic algorithms
117     Journal of Marine Systems, 18, 115-134
118    
119     Adcroft, A., Hill C. and J. Marshall: (1999)
120     A new treatment of the Coriolis terms in C-grid models at both high and low
121     resolutions,
122     Mon. Wea. Rev. Vol 127, pages 1928-1936
123    
124     Hill, C, Adcroft,A., Jamous,D., and J. Marshall, (1999)
125     A Strategy for Terascale Climate Modeling.
126 cnh 1.14 In Proceedings of the Eighth ECMWF Workshop on the Use of Parallel Processors
127     in Meteorology, pages 406-425
128     World Scientific Publishing Co: UK
129 cnh 1.12
130     Marotzke, J, Giering,R., Zhang, K.Q., Stammer,D., Hill,C., and T.Lee, (1999)
131     Construction of the adjoint MIT ocean general circulation model and
132     application to Atlantic heat transport variability
133     J. Geophysical Res., 104(C12), 29,529-29,547.
134    
135     \end{verbatim}
136 cnh 1.1
137     We begin by briefly showing some of the results of the model in action to
138     give a feel for the wide range of problems that can be addressed using it.
139    
140 jmc 1.22 % $Header: /u/gcmpack/manual/part1/manual.tex,v 1.21 2004/10/16 03:40:12 edhill Exp $
141 cnh 1.1 % $Name: $
142    
143     \section{Illustrations of the model in action}
144    
145     The MITgcm has been designed and used to model a wide range of phenomena,
146     from convection on the scale of meters in the ocean to the global pattern of
147 cnh 1.7 atmospheric winds - see figure \ref{fig:all-scales}. To give a flavor of the
148 cnh 1.1 kinds of problems the model has been used to study, we briefly describe some
149     of them here. A more detailed description of the underlying formulation,
150     numerical algorithm and implementation that lie behind these calculations is
151 cnh 1.2 given later. Indeed many of the illustrative examples shown below can be
152     easily reproduced: simply download the model (the minimum you need is a PC
153 cnh 1.10 running Linux, together with a FORTRAN\ 77 compiler) and follow the examples
154 cnh 1.2 described in detail in the documentation.
155 cnh 1.1
156     \subsection{Global atmosphere: `Held-Suarez' benchmark}
157 afe 1.18 \begin{rawhtml}
158 afe 1.19 <!-- CMIREDIR:atmospheric_example: -->
159 afe 1.18 \end{rawhtml}
160    
161    
162 cnh 1.1
163 cnh 1.7 A novel feature of MITgcm is its ability to simulate, using one basic algorithm,
164     both atmospheric and oceanographic flows at both small and large scales.
165 cnh 1.2
166 cnh 1.7 Figure \ref{fig:eddy_cs} shows an instantaneous plot of the 500$mb$
167 cnh 1.2 temperature field obtained using the atmospheric isomorph of MITgcm run at
168     2.8$^{\circ }$ resolution on the cubed sphere. We see cold air over the pole
169     (blue) and warm air along an equatorial band (red). Fully developed
170     baroclinic eddies spawned in the northern hemisphere storm track are
171     evident. There are no mountains or land-sea contrast in this calculation,
172     but you can easily put them in. The model is driven by relaxation to a
173     radiative-convective equilibrium profile, following the description set out
174     in Held and Suarez; 1994 designed to test atmospheric hydrodynamical cores -
175     there are no mountains or land-sea contrast.
176    
177 cnh 1.3 %% CNHbegin
178     \input{part1/cubic_eddies_figure}
179     %% CNHend
180    
181 cnh 1.2 As described in Adcroft (2001), a `cubed sphere' is used to discretize the
182 cnh 1.10 globe permitting a uniform griding and obviated the need to Fourier filter.
183 cnh 1.2 The `vector-invariant' form of MITgcm supports any orthogonal curvilinear
184     grid, of which the cubed sphere is just one of many choices.
185 cnh 1.1
186 cnh 1.7 Figure \ref{fig:hs_zave_u} shows the 5-year mean, zonally averaged zonal
187     wind from a 20-level configuration of
188 cnh 1.2 the model. It compares favorable with more conventional spatial
189 cnh 1.7 discretization approaches. The two plots show the field calculated using the
190     cube-sphere grid and the flow calculated using a regular, spherical polar
191     latitude-longitude grid. Both grids are supported within the model.
192 cnh 1.1
193 cnh 1.3 %% CNHbegin
194     \input{part1/hs_zave_u_figure}
195     %% CNHend
196    
197 cnh 1.2 \subsection{Ocean gyres}
198 afe 1.18 \begin{rawhtml}
199 afe 1.19 <!-- CMIREDIR:oceanic_example: -->
200 afe 1.18 \end{rawhtml}
201     \begin{rawhtml}
202 afe 1.19 <!-- CMIREDIR:ocean_gyres: -->
203 afe 1.18 \end{rawhtml}
204 cnh 1.1
205 cnh 1.2 Baroclinic instability is a ubiquitous process in the ocean, as well as the
206     atmosphere. Ocean eddies play an important role in modifying the
207     hydrographic structure and current systems of the oceans. Coarse resolution
208     models of the oceans cannot resolve the eddy field and yield rather broad,
209     diffusive patterns of ocean currents. But if the resolution of our models is
210     increased until the baroclinic instability process is resolved, numerical
211     solutions of a different and much more realistic kind, can be obtained.
212    
213 cnh 1.7 Figure \ref{fig:ocean-gyres} shows the surface temperature and velocity
214     field obtained from MITgcm run at $\frac{1}{6}^{\circ }$ horizontal
215     resolution on a $lat-lon$
216 cnh 1.2 grid in which the pole has been rotated by 90$^{\circ }$ on to the equator
217     (to avoid the converging of meridian in northern latitudes). 21 vertical
218     levels are used in the vertical with a `lopped cell' representation of
219     topography. The development and propagation of anomalously warm and cold
220 cnh 1.7 eddies can be clearly seen in the Gulf Stream region. The transport of
221 cnh 1.2 warm water northward by the mean flow of the Gulf Stream is also clearly
222     visible.
223 cnh 1.1
224 cnh 1.3 %% CNHbegin
225 cnh 1.11 \input{part1/atl6_figure}
226 cnh 1.3 %% CNHend
227    
228    
229 cnh 1.1 \subsection{Global ocean circulation}
230 afe 1.18 \begin{rawhtml}
231 afe 1.19 <!-- CMIREDIR:global_ocean_circulation: -->
232 afe 1.18 \end{rawhtml}
233 cnh 1.1
234 cnh 1.7 Figure \ref{fig:large-scale-circ} (top) shows the pattern of ocean currents at
235     the surface of a 4$^{\circ }$
236 cnh 1.2 global ocean model run with 15 vertical levels. Lopped cells are used to
237     represent topography on a regular $lat-lon$ grid extending from 70$^{\circ
238     }N $ to 70$^{\circ }S$. The model is driven using monthly-mean winds with
239     mixed boundary conditions on temperature and salinity at the surface. The
240     transfer properties of ocean eddies, convection and mixing is parameterized
241     in this model.
242    
243 cnh 1.7 Figure \ref{fig:large-scale-circ} (bottom) shows the meridional overturning
244     circulation of the global ocean in Sverdrups.
245 cnh 1.2
246 cnh 1.3 %%CNHbegin
247     \input{part1/global_circ_figure}
248     %%CNHend
249    
250 cnh 1.2 \subsection{Convection and mixing over topography}
251 afe 1.18 \begin{rawhtml}
252 afe 1.19 <!-- CMIREDIR:mixing_over_topography: -->
253 afe 1.18 \end{rawhtml}
254    
255 cnh 1.2
256     Dense plumes generated by localized cooling on the continental shelf of the
257     ocean may be influenced by rotation when the deformation radius is smaller
258     than the width of the cooling region. Rather than gravity plumes, the
259     mechanism for moving dense fluid down the shelf is then through geostrophic
260 adcroft 1.9 eddies. The simulation shown in the figure \ref{fig:convect-and-topo}
261 cnh 1.7 (blue is cold dense fluid, red is
262 cnh 1.2 warmer, lighter fluid) employs the non-hydrostatic capability of MITgcm to
263     trigger convection by surface cooling. The cold, dense water falls down the
264     slope but is deflected along the slope by rotation. It is found that
265     entrainment in the vertical plane is reduced when rotational control is
266     strong, and replaced by lateral entrainment due to the baroclinic
267     instability of the along-slope current.
268 cnh 1.1
269 cnh 1.3 %%CNHbegin
270     \input{part1/convect_and_topo}
271     %%CNHend
272    
273 cnh 1.1 \subsection{Boundary forced internal waves}
274 afe 1.18 \begin{rawhtml}
275 afe 1.19 <!-- CMIREDIR:boundary_forced_internal_waves: -->
276 afe 1.18 \end{rawhtml}
277 cnh 1.1
278 cnh 1.2 The unique ability of MITgcm to treat non-hydrostatic dynamics in the
279     presence of complex geometry makes it an ideal tool to study internal wave
280     dynamics and mixing in oceanic canyons and ridges driven by large amplitude
281     barotropic tidal currents imposed through open boundary conditions.
282    
283 cnh 1.7 Fig. \ref{fig:boundary-forced-wave} shows the influence of cross-slope
284     topographic variations on
285 cnh 1.2 internal wave breaking - the cross-slope velocity is in color, the density
286     contoured. The internal waves are excited by application of open boundary
287 cnh 1.7 conditions on the left. They propagate to the sloping boundary (represented
288 cnh 1.2 using MITgcm's finite volume spatial discretization) where they break under
289     nonhydrostatic dynamics.
290    
291 cnh 1.3 %%CNHbegin
292     \input{part1/boundary_forced_waves}
293     %%CNHend
294    
295 cnh 1.2 \subsection{Parameter sensitivity using the adjoint of MITgcm}
296 afe 1.18 \begin{rawhtml}
297 afe 1.19 <!-- CMIREDIR:parameter_sensitivity: -->
298 afe 1.18 \end{rawhtml}
299 cnh 1.2
300     Forward and tangent linear counterparts of MITgcm are supported using an
301     `automatic adjoint compiler'. These can be used in parameter sensitivity and
302     data assimilation studies.
303    
304 cnh 1.7 As one example of application of the MITgcm adjoint, Figure \ref{fig:hf-sensitivity}
305     maps the gradient $\frac{\partial J}{\partial \mathcal{H}}$where $J$ is the magnitude
306 cnh 1.10 of the overturning stream-function shown in figure \ref{fig:large-scale-circ}
307 cnh 1.7 at 60$^{\circ }$N and $
308     \mathcal{H}(\lambda,\varphi)$ is the mean, local air-sea heat flux over
309     a 100 year period. We see that $J$ is
310 cnh 1.2 sensitive to heat fluxes over the Labrador Sea, one of the important sources
311     of deep water for the thermohaline circulations. This calculation also
312     yields sensitivities to all other model parameters.
313    
314 cnh 1.3 %%CNHbegin
315     \input{part1/adj_hf_ocean_figure}
316     %%CNHend
317    
318 cnh 1.2 \subsection{Global state estimation of the ocean}
319 afe 1.18 \begin{rawhtml}
320 afe 1.19 <!-- CMIREDIR:global_state_estimation: -->
321 afe 1.18 \end{rawhtml}
322    
323 cnh 1.2
324     An important application of MITgcm is in state estimation of the global
325     ocean circulation. An appropriately defined `cost function', which measures
326     the departure of the model from observations (both remotely sensed and
327 cnh 1.10 in-situ) over an interval of time, is minimized by adjusting `control
328 cnh 1.2 parameters' such as air-sea fluxes, the wind field, the initial conditions
329 cnh 1.15 etc. Figure \ref{fig:assimilated-globes} shows the large scale planetary
330     circulation and a Hopf-Muller plot of Equatorial sea-surface height.
331     Both are obtained from assimilation bringing the model in to
332 cnh 1.7 consistency with altimetric and in-situ observations over the period
333 cnh 1.15 1992-1997.
334 cnh 1.2
335 cnh 1.3 %% CNHbegin
336 cnh 1.13 \input{part1/assim_figure}
337 cnh 1.3 %% CNHend
338    
339 cnh 1.2 \subsection{Ocean biogeochemical cycles}
340 afe 1.18 \begin{rawhtml}
341 afe 1.19 <!-- CMIREDIR:ocean_biogeo_cycles: -->
342 afe 1.18 \end{rawhtml}
343 cnh 1.2
344     MITgcm is being used to study global biogeochemical cycles in the ocean. For
345     example one can study the effects of interannual changes in meteorological
346     forcing and upper ocean circulation on the fluxes of carbon dioxide and
347 cnh 1.7 oxygen between the ocean and atmosphere. Figure \ref{fig:biogeo} shows
348     the annual air-sea flux of oxygen and its relation to density outcrops in
349     the southern oceans from a single year of a global, interannually varying
350     simulation. The simulation is run at $1^{\circ}\times1^{\circ}$ resolution
351     telescoping to $\frac{1}{3}^{\circ}\times\frac{1}{3}^{\circ}$ in the tropics (not shown).
352 cnh 1.2
353 cnh 1.3 %%CNHbegin
354     \input{part1/biogeo_figure}
355     %%CNHend
356 cnh 1.2
357     \subsection{Simulations of laboratory experiments}
358 afe 1.18 \begin{rawhtml}
359 afe 1.19 <!-- CMIREDIR:classroom_exp: -->
360 afe 1.18 \end{rawhtml}
361 cnh 1.2
362 cnh 1.7 Figure \ref{fig:lab-simulation} shows MITgcm being used to simulate a
363 edhill 1.17 laboratory experiment inquiring into the dynamics of the Antarctic Circumpolar Current (ACC). An
364 cnh 1.2 initially homogeneous tank of water ($1m$ in diameter) is driven from its
365     free surface by a rotating heated disk. The combined action of mechanical
366     and thermal forcing creates a lens of fluid which becomes baroclinically
367     unstable. The stratification and depth of penetration of the lens is
368 cnh 1.7 arrested by its instability in a process analogous to that which sets the
369 cnh 1.2 stratification of the ACC.
370 cnh 1.1
371 cnh 1.3 %%CNHbegin
372     \input{part1/lab_figure}
373     %%CNHend
374    
375 jmc 1.22 % $Header: /u/gcmpack/manual/part1/manual.tex,v 1.21 2004/10/16 03:40:12 edhill Exp $
376 cnh 1.1 % $Name: $
377    
378     \section{Continuous equations in `r' coordinates}
379 afe 1.18 \begin{rawhtml}
380 afe 1.19 <!-- CMIREDIR:z-p_isomorphism: -->
381 afe 1.18 \end{rawhtml}
382 cnh 1.1
383     To render atmosphere and ocean models from one dynamical core we exploit
384     `isomorphisms' between equation sets that govern the evolution of the
385 cnh 1.7 respective fluids - see figure \ref{fig:isomorphic-equations}.
386     One system of hydrodynamical equations is written down
387 cnh 1.1 and encoded. The model variables have different interpretations depending on
388     whether the atmosphere or ocean is being studied. Thus, for example, the
389     vertical coordinate `$r$' is interpreted as pressure, $p$, if we are
390 edhill 1.17 modeling the atmosphere (right hand side of figure \ref{fig:isomorphic-equations})
391     and height, $z$, if we are modeling the ocean (left hand side of figure
392 cnh 1.7 \ref{fig:isomorphic-equations}).
393 cnh 1.1
394 cnh 1.3 %%CNHbegin
395     \input{part1/zandpcoord_figure.tex}
396     %%CNHend
397    
398 cnh 1.1 The state of the fluid at any time is characterized by the distribution of
399     velocity $\vec{\mathbf{v}}$, active tracers $\theta $ and $S$, a
400     `geopotential' $\phi $ and density $\rho =\rho (\theta ,S,p)$ which may
401     depend on $\theta $, $S$, and $p$. The equations that govern the evolution
402     of these fields, obtained by applying the laws of classical mechanics and
403     thermodynamics to a Boussinesq, Navier-Stokes fluid are, written in terms of
404 cnh 1.7 a generic vertical coordinate, $r$, so that the appropriate
405     kinematic boundary conditions can be applied isomorphically
406     see figure \ref{fig:zandp-vert-coord}.
407 cnh 1.1
408 cnh 1.3 %%CNHbegin
409     \input{part1/vertcoord_figure.tex}
410     %%CNHend
411    
412 jmc 1.20 \begin{equation}
413 adcroft 1.4 \frac{D\vec{\mathbf{v}_{h}}}{Dt}+\left( 2\vec{\Omega}\times \vec{\mathbf{v}}
414     \right) _{h}+\mathbf{\nabla }_{h}\phi =\mathcal{F}_{\vec{\mathbf{v}_{h}}}
415 cnh 1.8 \text{ horizontal mtm} \label{eq:horizontal_mtm}
416 jmc 1.20 \end{equation}
417 cnh 1.1
418 cnh 1.8 \begin{equation}
419 adcroft 1.4 \frac{D\dot{r}}{Dt}+\widehat{k}\cdot \left( 2\vec{\Omega}\times \vec{\mathbf{
420 cnh 1.1 v}}\right) +\frac{\partial \phi }{\partial r}+b=\mathcal{F}_{\dot{r}}\text{
421 cnh 1.8 vertical mtm} \label{eq:vertical_mtm}
422     \end{equation}
423 cnh 1.1
424     \begin{equation}
425 adcroft 1.4 \mathbf{\nabla }_{h}\cdot \vec{\mathbf{v}}_{h}+\frac{\partial \dot{r}}{
426 cnh 1.8 \partial r}=0\text{ continuity} \label{eq:continuity}
427 cnh 1.1 \end{equation}
428    
429 cnh 1.8 \begin{equation}
430     b=b(\theta ,S,r)\text{ equation of state} \label{eq:equation_of_state}
431     \end{equation}
432 cnh 1.1
433 cnh 1.8 \begin{equation}
434 cnh 1.2 \frac{D\theta }{Dt}=\mathcal{Q}_{\theta }\text{ potential temperature}
435 cnh 1.8 \label{eq:potential_temperature}
436     \end{equation}
437 cnh 1.1
438 cnh 1.8 \begin{equation}
439 cnh 1.2 \frac{DS}{Dt}=\mathcal{Q}_{S}\text{ humidity/salinity}
440 adcroft 1.9 \label{eq:humidity_salt}
441 cnh 1.8 \end{equation}
442 cnh 1.1
443     Here:
444    
445     \begin{equation*}
446 cnh 1.2 r\text{ is the vertical coordinate}
447 cnh 1.1 \end{equation*}
448    
449     \begin{equation*}
450     \frac{D}{Dt}=\frac{\partial }{\partial t}+\vec{\mathbf{v}}\cdot \nabla \text{
451 cnh 1.2 is the total derivative}
452 cnh 1.1 \end{equation*}
453    
454     \begin{equation*}
455 adcroft 1.4 \mathbf{\nabla }=\mathbf{\nabla }_{h}+\widehat{k}\frac{\partial }{\partial r}
456 cnh 1.2 \text{ is the `grad' operator}
457 cnh 1.1 \end{equation*}
458 adcroft 1.4 with $\mathbf{\nabla }_{h}$ operating in the horizontal and $\widehat{k}
459 cnh 1.1 \frac{\partial }{\partial r}$ operating in the vertical, where $\widehat{k}$
460     is a unit vector in the vertical
461    
462     \begin{equation*}
463 cnh 1.2 t\text{ is time}
464 cnh 1.1 \end{equation*}
465    
466     \begin{equation*}
467     \vec{\mathbf{v}}=(u,v,\dot{r})=(\vec{\mathbf{v}}_{h},\dot{r})\text{ is the
468 cnh 1.2 velocity}
469 cnh 1.1 \end{equation*}
470    
471     \begin{equation*}
472 cnh 1.2 \phi \text{ is the `pressure'/`geopotential'}
473 cnh 1.1 \end{equation*}
474    
475     \begin{equation*}
476 cnh 1.2 \vec{\Omega}\text{ is the Earth's rotation}
477 cnh 1.1 \end{equation*}
478    
479     \begin{equation*}
480 cnh 1.2 b\text{ is the `buoyancy'}
481 cnh 1.1 \end{equation*}
482    
483     \begin{equation*}
484 cnh 1.2 \theta \text{ is potential temperature}
485 cnh 1.1 \end{equation*}
486    
487     \begin{equation*}
488 cnh 1.2 S\text{ is specific humidity in the atmosphere; salinity in the ocean}
489 cnh 1.1 \end{equation*}
490    
491     \begin{equation*}
492 adcroft 1.4 \mathcal{F}_{\vec{\mathbf{v}}}\text{ are forcing and dissipation of }\vec{
493 cnh 1.1 \mathbf{v}}
494     \end{equation*}
495    
496     \begin{equation*}
497 cnh 1.2 \mathcal{Q}_{\theta }\mathcal{\ }\text{are forcing and dissipation of }\theta
498 cnh 1.1 \end{equation*}
499    
500     \begin{equation*}
501     \mathcal{Q}_{S}\mathcal{\ }\text{are forcing and dissipation of }S
502     \end{equation*}
503    
504     The $\mathcal{F}^{\prime }s$ and $\mathcal{Q}^{\prime }s$ are provided by
505 cnh 1.7 `physics' and forcing packages for atmosphere and ocean. These are described
506     in later chapters.
507 cnh 1.1
508     \subsection{Kinematic Boundary conditions}
509    
510     \subsubsection{vertical}
511    
512 cnh 1.7 at fixed and moving $r$ surfaces we set (see figure \ref{fig:zandp-vert-coord}):
513 cnh 1.1
514     \begin{equation}
515 edhill 1.17 \dot{r}=0 \text{\ at\ } r=R_{fixed}(x,y)\text{ (ocean bottom, top of the atmosphere)}
516 cnh 1.1 \label{eq:fixedbc}
517     \end{equation}
518    
519     \begin{equation}
520 edhill 1.17 \dot{r}=\frac{Dr}{Dt} \text{\ at\ } r=R_{moving}\text{ \
521 cnh 1.10 (ocean surface,bottom of the atmosphere)} \label{eq:movingbc}
522 cnh 1.1 \end{equation}
523    
524     Here
525    
526     \begin{equation*}
527 cnh 1.2 R_{moving}=R_{o}+\eta
528 cnh 1.1 \end{equation*}
529     where $R_{o}(x,y)$ is the `$r-$value' (height or pressure, depending on
530     whether we are in the atmosphere or ocean) of the `moving surface' in the
531     resting fluid and $\eta $ is the departure from $R_{o}(x,y)$ in the presence
532     of motion.
533    
534     \subsubsection{horizontal}
535    
536     \begin{equation}
537     \vec{\mathbf{v}}\cdot \vec{\mathbf{n}}=0 \label{eq:noflow}
538 adcroft 1.4 \end{equation}
539 cnh 1.1 where $\vec{\mathbf{n}}$ is the normal to a solid boundary.
540    
541     \subsection{Atmosphere}
542    
543 cnh 1.7 In the atmosphere, (see figure \ref{fig:zandp-vert-coord}), we interpret:
544 cnh 1.1
545     \begin{equation}
546     r=p\text{ is the pressure} \label{eq:atmos-r}
547     \end{equation}
548    
549     \begin{equation}
550     \dot{r}=\frac{Dp}{Dt}=\omega \text{ is the vertical velocity in }p\text{
551     coordinates} \label{eq:atmos-omega}
552     \end{equation}
553    
554     \begin{equation}
555     \phi =g\,z\text{ is the geopotential height} \label{eq:atmos-phi}
556     \end{equation}
557    
558     \begin{equation}
559     b=\frac{\partial \Pi }{\partial p}\theta \text{ is the buoyancy}
560     \label{eq:atmos-b}
561     \end{equation}
562    
563     \begin{equation}
564     \theta =T(\frac{p_{c}}{p})^{\kappa }\text{ is potential temperature}
565     \label{eq:atmos-theta}
566     \end{equation}
567    
568     \begin{equation}
569     S=q,\text{ is the specific humidity} \label{eq:atmos-s}
570     \end{equation}
571     where
572    
573     \begin{equation*}
574     T\text{ is absolute temperature}
575 adcroft 1.4 \end{equation*}
576 cnh 1.1 \begin{equation*}
577     p\text{ is the pressure}
578 adcroft 1.4 \end{equation*}
579 cnh 1.1 \begin{eqnarray*}
580     &&z\text{ is the height of the pressure surface} \\
581     &&g\text{ is the acceleration due to gravity}
582     \end{eqnarray*}
583    
584     In the above the ideal gas law, $p=\rho RT$, has been expressed in terms of
585     the Exner function $\Pi (p)$ given by (see Appendix Atmosphere)
586     \begin{equation}
587     \Pi (p)=c_{p}(\frac{p}{p_{c}})^{\kappa } \label{eq:exner}
588 adcroft 1.4 \end{equation}
589 cnh 1.1 where $p_{c}$ is a reference pressure and $\kappa =R/c_{p}$ with $R$ the gas
590     constant and $c_{p}$ the specific heat of air at constant pressure.
591    
592     At the top of the atmosphere (which is `fixed' in our $r$ coordinate):
593    
594     \begin{equation*}
595 cnh 1.2 R_{fixed}=p_{top}=0
596 cnh 1.1 \end{equation*}
597     In a resting atmosphere the elevation of the mountains at the bottom is
598     given by
599     \begin{equation*}
600 cnh 1.2 R_{moving}=R_{o}(x,y)=p_{o}(x,y)
601 cnh 1.1 \end{equation*}
602     i.e. the (hydrostatic) pressure at the top of the mountains in a resting
603     atmosphere.
604    
605     The boundary conditions at top and bottom are given by:
606    
607     \begin{eqnarray}
608     &&\omega =0~\text{at }r=R_{fixed} \text{ (top of the atmosphere)}
609     \label{eq:fixed-bc-atmos} \\
610     \omega &=&\frac{Dp_{s}}{Dt}\text{; at }r=R_{moving}\text{ (bottom of the
611     atmosphere)} \label{eq:moving-bc-atmos}
612     \end{eqnarray}
613    
614 edhill 1.21 Then the (hydrostatic form of) equations
615     (\ref{eq:horizontal_mtm}-\ref{eq:humidity_salt}) yields a consistent
616     set of atmospheric equations which, for convenience, are written out
617     in $p$ coordinates in Appendix Atmosphere - see
618     eqs(\ref{eq:atmos-prime}).
619 cnh 1.1
620     \subsection{Ocean}
621    
622     In the ocean we interpret:
623     \begin{eqnarray}
624     r &=&z\text{ is the height} \label{eq:ocean-z} \\
625     \dot{r} &=&\frac{Dz}{Dt}=w\text{ is the vertical velocity}
626     \label{eq:ocean-w} \\
627     \phi &=&\frac{p}{\rho _{c}}\text{ is the pressure} \label{eq:ocean-p} \\
628     b(\theta ,S,r) &=&\frac{g}{\rho _{c}}\left( \rho (\theta ,S,r)-\rho
629     _{c}\right) \text{ is the buoyancy} \label{eq:ocean-b}
630     \end{eqnarray}
631     where $\rho _{c}$ is a fixed reference density of water and $g$ is the
632     acceleration due to gravity.\noindent
633    
634     In the above
635    
636     At the bottom of the ocean: $R_{fixed}(x,y)=-H(x,y)$.
637    
638     The surface of the ocean is given by: $R_{moving}=\eta $
639    
640 adcroft 1.4 The position of the resting free surface of the ocean is given by $
641 cnh 1.1 R_{o}=Z_{o}=0$.
642    
643     Boundary conditions are:
644    
645     \begin{eqnarray}
646     w &=&0~\text{at }r=R_{fixed}\text{ (ocean bottom)} \label{eq:fixed-bc-ocean}
647     \\
648 adcroft 1.4 w &=&\frac{D\eta }{Dt}\text{ at }r=R_{moving}=\eta \text{ (ocean surface)
649 cnh 1.1 \label{eq:moving-bc-ocean}}
650     \end{eqnarray}
651     where $\eta $ is the elevation of the free surface.
652    
653 adcroft 1.9 Then equations (\ref{eq:horizontal_mtm}-\ref{eq:humidity_salt}) yield a consistent set
654 cnh 1.8 of oceanic equations
655 cnh 1.1 which, for convenience, are written out in $z$ coordinates in Appendix Ocean
656     - see eqs(\ref{eq:ocean-mom}) to (\ref{eq:ocean-salt}).
657    
658     \subsection{Hydrostatic, Quasi-hydrostatic, Quasi-nonhydrostatic and
659     Non-hydrostatic forms}
660 afe 1.18 \begin{rawhtml}
661 afe 1.19 <!-- CMIREDIR:non_hydrostatic: -->
662 afe 1.18 \end{rawhtml}
663    
664 cnh 1.1
665     Let us separate $\phi $ in to surface, hydrostatic and non-hydrostatic terms:
666    
667     \begin{equation}
668     \phi (x,y,r)=\phi _{s}(x,y)+\phi _{hyd}(x,y,r)+\phi _{nh}(x,y,r)
669     \label{eq:phi-split}
670 adcroft 1.4 \end{equation}
671 jmc 1.20 %and write eq(\ref{eq:incompressible}) in the form:
672     % ^- this eq is missing (jmc) ; replaced with:
673     and write eq( \ref{eq:horizontal_mtm}) in the form:
674 cnh 1.1
675     \begin{equation}
676     \frac{\partial \vec{\mathbf{v}_{h}}}{\partial t}+\mathbf{\nabla }_{h}\phi
677     _{s}+\mathbf{\nabla }_{h}\phi _{hyd}+\epsilon _{nh}\mathbf{\nabla }_{h}\phi
678     _{nh}=\vec{\mathbf{G}}_{\vec{v}_{h}} \label{eq:mom-h}
679     \end{equation}
680    
681     \begin{equation}
682     \frac{\partial \phi _{hyd}}{\partial r}=-b \label{eq:hydrostatic}
683     \end{equation}
684    
685     \begin{equation}
686 adcroft 1.4 \epsilon _{nh}\frac{\partial \dot{r}}{\partial t}+\frac{\partial \phi _{nh}}{
687 cnh 1.1 \partial r}=G_{\dot{r}} \label{eq:mom-w}
688     \end{equation}
689     Here $\epsilon _{nh}$ is a non-hydrostatic parameter.
690    
691 adcroft 1.4 The $\left( \vec{\mathbf{G}}_{\vec{v}},G_{\dot{r}}\right) $ in eq(\ref
692 cnh 1.1 {eq:mom-h}) and (\ref{eq:mom-w}) represent advective, metric and Coriolis
693 adcroft 1.4 terms in the momentum equations. In spherical coordinates they take the form
694     \footnote{
695 cnh 1.1 In the hydrostatic primitive equations (\textbf{HPE}) all underlined terms
696 adcroft 1.4 in (\ref{eq:gu-speherical}), (\ref{eq:gv-spherical}) and (\ref
697 cnh 1.1 {eq:gw-spherical}) are omitted; the singly-underlined terms are included in
698 adcroft 1.4 the quasi-hydrostatic model (\textbf{QH}). The fully non-hydrostatic model (
699 cnh 1.1 \textbf{NH}) includes all terms.} - see Marshall et al 1997a for a full
700     discussion:
701    
702     \begin{equation}
703     \left.
704     \begin{tabular}{l}
705     $G_{u}=-\vec{\mathbf{v}}.\nabla u$ \\
706 cnh 1.6 $-\left\{ \underline{\frac{u\dot{r}}{{r}}}-\frac{uv\tan \varphi}{{r}}\right\} $
707 cnh 1.1 \\
708 cnh 1.6 $-\left\{ -2\Omega v\sin \varphi+\underline{2\Omega \dot{r}\cos \varphi}\right\} $
709 cnh 1.1 \\
710 adcroft 1.4 $+\mathcal{F}_{u}$
711     \end{tabular}
712 cnh 1.1 \ \right\} \left\{
713     \begin{tabular}{l}
714     \textit{advection} \\
715     \textit{metric} \\
716     \textit{Coriolis} \\
717 adcroft 1.4 \textit{\ Forcing/Dissipation}
718     \end{tabular}
719 cnh 1.2 \ \right. \qquad \label{eq:gu-speherical}
720 cnh 1.1 \end{equation}
721    
722     \begin{equation}
723     \left.
724     \begin{tabular}{l}
725     $G_{v}=-\vec{\mathbf{v}}.\nabla v$ \\
726 cnh 1.6 $-\left\{ \underline{\frac{v\dot{r}}{{r}}}-\frac{u^{2}\tan \varphi}{{r}}\right\}
727 cnh 1.1 $ \\
728 cnh 1.6 $-\left\{ -2\Omega u\sin \varphi \right\} $ \\
729 adcroft 1.4 $+\mathcal{F}_{v}$
730     \end{tabular}
731 cnh 1.1 \ \right\} \left\{
732     \begin{tabular}{l}
733     \textit{advection} \\
734     \textit{metric} \\
735     \textit{Coriolis} \\
736 adcroft 1.4 \textit{\ Forcing/Dissipation}
737     \end{tabular}
738 cnh 1.2 \ \right. \qquad \label{eq:gv-spherical}
739 adcroft 1.4 \end{equation}
740 cnh 1.2 \qquad \qquad \qquad \qquad \qquad
741 cnh 1.1
742     \begin{equation}
743     \left.
744     \begin{tabular}{l}
745     $G_{\dot{r}}=-\underline{\underline{\vec{\mathbf{v}}.\nabla \dot{r}}}$ \\
746     $+\left\{ \underline{\frac{u^{_{^{2}}}+v^{2}}{{r}}}\right\} $ \\
747 cnh 1.6 ${+}\underline{{2\Omega u\cos \varphi}}$ \\
748 adcroft 1.4 $\underline{\underline{\mathcal{F}_{\dot{r}}}}$
749     \end{tabular}
750 cnh 1.1 \ \right\} \left\{
751     \begin{tabular}{l}
752     \textit{advection} \\
753     \textit{metric} \\
754     \textit{Coriolis} \\
755 adcroft 1.4 \textit{\ Forcing/Dissipation}
756     \end{tabular}
757 cnh 1.2 \ \right. \label{eq:gw-spherical}
758 adcroft 1.4 \end{equation}
759 cnh 1.2 \qquad \qquad \qquad \qquad \qquad
760 cnh 1.1
761 cnh 1.6 In the above `${r}$' is the distance from the center of the earth and `$\varphi$
762 cnh 1.1 ' is latitude.
763    
764     Grad and div operators in spherical coordinates are defined in appendix
765 adcroft 1.4 OPERATORS.
766 cnh 1.1
767 cnh 1.3 %%CNHbegin
768     \input{part1/sphere_coord_figure.tex}
769     %%CNHend
770    
771 cnh 1.1 \subsubsection{Shallow atmosphere approximation}
772    
773     Most models are based on the `hydrostatic primitive equations' (HPE's) in
774     which the vertical momentum equation is reduced to a statement of
775     hydrostatic balance and the `traditional approximation' is made in which the
776     Coriolis force is treated approximately and the shallow atmosphere
777     approximation is made.\ The MITgcm need not make the `traditional
778     approximation'. To be able to support consistent non-hydrostatic forms the
779 adcroft 1.4 shallow atmosphere approximation can be relaxed - when dividing through by $
780 cnh 1.2 r $ in, for example, (\ref{eq:gu-speherical}), we do not replace $r$ by $a$,
781 cnh 1.1 the radius of the earth.
782    
783     \subsubsection{Hydrostatic and quasi-hydrostatic forms}
784 cnh 1.7 \label{sec:hydrostatic_and_quasi-hydrostatic_forms}
785 cnh 1.1
786     These are discussed at length in Marshall et al (1997a).
787    
788     In the `hydrostatic primitive equations' (\textbf{HPE)} all the underlined
789     terms in Eqs. (\ref{eq:gu-speherical} $\rightarrow $\ \ref{eq:gw-spherical})
790     are neglected and `${r}$' is replaced by `$a$', the mean radius of the
791     earth. Once the pressure is found at one level - e.g. by inverting a 2-d
792     Elliptic equation for $\phi _{s}$ at $r=R_{moving}$ - the pressure can be
793 adcroft 1.4 computed at all other levels by integration of the hydrostatic relation, eq(
794 cnh 1.1 \ref{eq:hydrostatic}).
795    
796     In the `quasi-hydrostatic' equations (\textbf{QH)} strict balance between
797     gravity and vertical pressure gradients is not imposed. The $2\Omega u\cos
798 cnh 1.6 \varphi $ Coriolis term are not neglected and are balanced by a non-hydrostatic
799 adcroft 1.4 contribution to the pressure field: only the terms underlined twice in Eqs. (
800 cnh 1.1 \ref{eq:gu-speherical}$\rightarrow $\ \ref{eq:gw-spherical}) are set to zero
801     and, simultaneously, the shallow atmosphere approximation is relaxed. In
802     \textbf{QH}\ \textit{all} the metric terms are retained and the full
803     variation of the radial position of a particle monitored. The \textbf{QH}\
804     vertical momentum equation (\ref{eq:mom-w}) becomes:
805    
806     \begin{equation*}
807 cnh 1.6 \frac{\partial \phi _{nh}}{\partial r}=2\Omega u\cos \varphi
808 cnh 1.1 \end{equation*}
809     making a small correction to the hydrostatic pressure.
810    
811     \textbf{QH} has good energetic credentials - they are the same as for
812     \textbf{HPE}. Importantly, however, it has the same angular momentum
813     principle as the full non-hydrostatic model (\textbf{NH)} - see Marshall
814     et.al., 1997a. As in \textbf{HPE }only a 2-d elliptic problem need be solved.
815    
816     \subsubsection{Non-hydrostatic and quasi-nonhydrostatic forms}
817    
818     The MIT model presently supports a full non-hydrostatic ocean isomorph, but
819     only a quasi-non-hydrostatic atmospheric isomorph.
820    
821     \paragraph{Non-hydrostatic Ocean}
822    
823 adcroft 1.4 In the non-hydrostatic ocean model all terms in equations Eqs.(\ref
824 cnh 1.1 {eq:gu-speherical} $\rightarrow $\ \ref{eq:gw-spherical}) are retained. A
825     three dimensional elliptic equation must be solved subject to Neumann
826     boundary conditions (see below). It is important to note that use of the
827     full \textbf{NH} does not admit any new `fast' waves in to the system - the
828 cnh 1.8 incompressible condition eq(\ref{eq:continuity}) has already filtered out
829 cnh 1.1 acoustic modes. It does, however, ensure that the gravity waves are treated
830     accurately with an exact dispersion relation. The \textbf{NH} set has a
831     complete angular momentum principle and consistent energetics - see White
832     and Bromley, 1995; Marshall et.al.\ 1997a.
833    
834     \paragraph{Quasi-nonhydrostatic Atmosphere}
835    
836 adcroft 1.4 In the non-hydrostatic version of our atmospheric model we approximate $\dot{
837 cnh 1.1 r}$ in the vertical momentum eqs(\ref{eq:mom-w}) and (\ref{eq:gv-spherical})
838     (but only here) by:
839    
840     \begin{equation}
841     \dot{r}=\frac{Dp}{Dt}=\frac{1}{g}\frac{D\phi }{Dt} \label{eq:quasi-nh-w}
842 adcroft 1.4 \end{equation}
843 cnh 1.1 where $p_{hy}$ is the hydrostatic pressure.
844    
845     \subsubsection{Summary of equation sets supported by model}
846    
847     \paragraph{Atmosphere}
848    
849     Hydrostatic, and quasi-hydrostatic and quasi non-hydrostatic forms of the
850     compressible non-Boussinesq equations in $p-$coordinates are supported.
851    
852     \subparagraph{Hydrostatic and quasi-hydrostatic}
853    
854     The hydrostatic set is written out in $p-$coordinates in appendix Atmosphere
855     - see eq(\ref{eq:atmos-prime}).
856    
857     \subparagraph{Quasi-nonhydrostatic}
858    
859     A quasi-nonhydrostatic form is also supported.
860    
861     \paragraph{Ocean}
862    
863     \subparagraph{Hydrostatic and quasi-hydrostatic}
864    
865     Hydrostatic, and quasi-hydrostatic forms of the incompressible Boussinesq
866     equations in $z-$coordinates are supported.
867    
868     \subparagraph{Non-hydrostatic}
869    
870 adcroft 1.4 Non-hydrostatic forms of the incompressible Boussinesq equations in $z-$
871     coordinates are supported - see eqs(\ref{eq:ocean-mom}) to (\ref
872 cnh 1.1 {eq:ocean-salt}).
873    
874     \subsection{Solution strategy}
875    
876 adcroft 1.4 The method of solution employed in the \textbf{HPE}, \textbf{QH} and \textbf{
877 cnh 1.8 NH} models is summarized in Figure \ref{fig:solution-strategy}.
878     Under all dynamics, a 2-d elliptic equation is
879 cnh 1.1 first solved to find the surface pressure and the hydrostatic pressure at
880     any level computed from the weight of fluid above. Under \textbf{HPE} and
881     \textbf{QH} dynamics, the horizontal momentum equations are then stepped
882     forward and $\dot{r}$ found from continuity. Under \textbf{NH} dynamics a
883     3-d elliptic equation must be solved for the non-hydrostatic pressure before
884     stepping forward the horizontal momentum equations; $\dot{r}$ is found by
885     stepping forward the vertical momentum equation.
886    
887 cnh 1.3 %%CNHbegin
888     \input{part1/solution_strategy_figure.tex}
889     %%CNHend
890    
891 cnh 1.1 There is no penalty in implementing \textbf{QH} over \textbf{HPE} except, of
892 cnh 1.6 course, some complication that goes with the inclusion of $\cos \varphi \ $
893 cnh 1.1 Coriolis terms and the relaxation of the shallow atmosphere approximation.
894     But this leads to negligible increase in computation. In \textbf{NH}, in
895     contrast, one additional elliptic equation - a three-dimensional one - must
896     be inverted for $p_{nh}$. However the `overhead' of the \textbf{NH} model is
897     essentially negligible in the hydrostatic limit (see detailed discussion in
898     Marshall et al, 1997) resulting in a non-hydrostatic algorithm that, in the
899     hydrostatic limit, is as computationally economic as the \textbf{HPEs}.
900    
901     \subsection{Finding the pressure field}
902 cnh 1.7 \label{sec:finding_the_pressure_field}
903 cnh 1.1
904     Unlike the prognostic variables $u$, $v$, $w$, $\theta $ and $S$, the
905     pressure field must be obtained diagnostically. We proceed, as before, by
906     dividing the total (pressure/geo) potential in to three parts, a surface
907     part, $\phi _{s}(x,y)$, a hydrostatic part $\phi _{hyd}(x,y,r)$ and a
908     non-hydrostatic part $\phi _{nh}(x,y,r)$, as in (\ref{eq:phi-split}), and
909     writing the momentum equation as in (\ref{eq:mom-h}).
910    
911     \subsubsection{Hydrostatic pressure}
912    
913     Hydrostatic pressure is obtained by integrating (\ref{eq:hydrostatic})
914     vertically from $r=R_{o}$ where $\phi _{hyd}(r=R_{o})=0$, to yield:
915    
916     \begin{equation*}
917 adcroft 1.4 \int_{r}^{R_{o}}\frac{\partial \phi _{hyd}}{\partial r}dr=\left[ \phi _{hyd}
918 cnh 1.2 \right] _{r}^{R_{o}}=\int_{r}^{R_{o}}-bdr
919 cnh 1.1 \end{equation*}
920     and so
921    
922     \begin{equation}
923     \phi _{hyd}(x,y,r)=\int_{r}^{R_{o}}bdr \label{eq:hydro-phi}
924     \end{equation}
925    
926     The model can be easily modified to accommodate a loading term (e.g
927     atmospheric pressure pushing down on the ocean's surface) by setting:
928    
929     \begin{equation}
930     \phi _{hyd}(r=R_{o})=loading \label{eq:loading}
931     \end{equation}
932    
933     \subsubsection{Surface pressure}
934    
935 cnh 1.8 The surface pressure equation can be obtained by integrating continuity,
936     (\ref{eq:continuity}), vertically from $r=R_{fixed}$ to $r=R_{moving}$
937 cnh 1.1
938     \begin{equation*}
939 adcroft 1.4 \int_{R_{fixed}}^{R_{moving}}\left( \mathbf{\nabla }_{h}\cdot \vec{\mathbf{v}
940 cnh 1.2 }_{h}+\partial _{r}\dot{r}\right) dr=0
941 cnh 1.1 \end{equation*}
942    
943     Thus:
944    
945     \begin{equation*}
946     \frac{\partial \eta }{\partial t}+\vec{\mathbf{v}}.\nabla \eta
947 adcroft 1.4 +\int_{R_{fixed}}^{R_{moving}}\mathbf{\nabla }_{h}\cdot \vec{\mathbf{v}}
948 cnh 1.2 _{h}dr=0
949 cnh 1.1 \end{equation*}
950 adcroft 1.4 where $\eta =R_{moving}-R_{o}$ is the free-surface $r$-anomaly in units of $
951 cnh 1.1 r $. The above can be rearranged to yield, using Leibnitz's theorem:
952    
953     \begin{equation}
954     \frac{\partial \eta }{\partial t}+\mathbf{\nabla }_{h}\cdot
955     \int_{R_{fixed}}^{R_{moving}}\vec{\mathbf{v}}_{h}dr=\text{source}
956     \label{eq:free-surface}
957 adcroft 1.4 \end{equation}
958 cnh 1.1 where we have incorporated a source term.
959    
960     Whether $\phi $ is pressure (ocean model, $p/\rho _{c}$) or geopotential
961 cnh 1.8 (atmospheric model), in (\ref{eq:mom-h}), the horizontal gradient term can
962 cnh 1.1 be written
963     \begin{equation}
964 cnh 1.2 \mathbf{\nabla }_{h}\phi _{s}=\mathbf{\nabla }_{h}\left( b_{s}\eta \right)
965 cnh 1.1 \label{eq:phi-surf}
966 adcroft 1.4 \end{equation}
967 cnh 1.1 where $b_{s}$ is the buoyancy at the surface.
968    
969 cnh 1.8 In the hydrostatic limit ($\epsilon _{nh}=0$), equations (\ref{eq:mom-h}), (\ref
970 cnh 1.1 {eq:free-surface}) and (\ref{eq:phi-surf}) can be solved by inverting a 2-d
971     elliptic equation for $\phi _{s}$ as described in Chapter 2. Both `free
972     surface' and `rigid lid' approaches are available.
973    
974     \subsubsection{Non-hydrostatic pressure}
975    
976 cnh 1.8 Taking the horizontal divergence of (\ref{eq:mom-h}) and adding
977     $\frac{\partial }{\partial r}$ of (\ref{eq:mom-w}), invoking the continuity equation
978     (\ref{eq:continuity}), we deduce that:
979 cnh 1.1
980     \begin{equation}
981 adcroft 1.4 \nabla _{3}^{2}\phi _{nh}=\nabla .\vec{\mathbf{G}}_{\vec{v}}-\left( \mathbf{
982     \nabla }_{h}^{2}\phi _{s}+\mathbf{\nabla }^{2}\phi _{hyd}\right) =\nabla .
983 cnh 1.1 \vec{\mathbf{F}} \label{eq:3d-invert}
984     \end{equation}
985    
986     For a given rhs this 3-d elliptic equation must be inverted for $\phi _{nh}$
987     subject to appropriate choice of boundary conditions. This method is usually
988     called \textit{The Pressure Method} [Harlow and Welch, 1965; Williams, 1969;
989     Potter, 1976]. In the hydrostatic primitive equations case (\textbf{HPE}),
990     the 3-d problem does not need to be solved.
991    
992     \paragraph{Boundary Conditions}
993    
994     We apply the condition of no normal flow through all solid boundaries - the
995     coasts (in the ocean) and the bottom:
996    
997     \begin{equation}
998     \vec{\mathbf{v}}.\widehat{n}=0 \label{nonormalflow}
999     \end{equation}
1000     where $\widehat{n}$ is a vector of unit length normal to the boundary. The
1001     kinematic condition (\ref{nonormalflow}) is also applied to the vertical
1002 adcroft 1.4 velocity at $r=R_{moving}$. No-slip $\left( v_{T}=0\right) \ $or slip $
1003 cnh 1.1 \left( \partial v_{T}/\partial n=0\right) \ $conditions are employed on the
1004     tangential component of velocity, $v_{T}$, at all solid boundaries,
1005     depending on the form chosen for the dissipative terms in the momentum
1006     equations - see below.
1007    
1008 cnh 1.8 Eq.(\ref{nonormalflow}) implies, making use of (\ref{eq:mom-h}), that:
1009 cnh 1.1
1010     \begin{equation}
1011     \widehat{n}.\nabla \phi _{nh}=\widehat{n}.\vec{\mathbf{F}}
1012     \label{eq:inhom-neumann-nh}
1013     \end{equation}
1014     where
1015    
1016     \begin{equation*}
1017     \vec{\mathbf{F}}=\vec{\mathbf{G}}_{\vec{v}}-\left( \mathbf{\nabla }_{h}\phi
1018     _{s}+\mathbf{\nabla }\phi _{hyd}\right)
1019 adcroft 1.4 \end{equation*}
1020 cnh 1.1 presenting inhomogeneous Neumann boundary conditions to the Elliptic problem
1021     (\ref{eq:3d-invert}). As shown, for example, by Williams (1969), one can
1022     exploit classical 3D potential theory and, by introducing an appropriately
1023 cnh 1.2 chosen $\delta $-function sheet of `source-charge', replace the
1024     inhomogeneous boundary condition on pressure by a homogeneous one. The
1025 adcroft 1.4 source term $rhs$ in (\ref{eq:3d-invert}) is the divergence of the vector $
1026     \vec{\mathbf{F}}.$ By simultaneously setting $
1027 cnh 1.1 \begin{array}{l}
1028 adcroft 1.4 \widehat{n}.\vec{\mathbf{F}}
1029     \end{array}
1030 cnh 1.1 =0$\ and $\widehat{n}.\nabla \phi _{nh}=0\ $on the boundary the following
1031 cnh 1.2 self-consistent but simpler homogenized Elliptic problem is obtained:
1032 cnh 1.1
1033     \begin{equation*}
1034 cnh 1.2 \nabla ^{2}\phi _{nh}=\nabla .\widetilde{\vec{\mathbf{F}}}\qquad
1035 adcroft 1.4 \end{equation*}
1036 cnh 1.1 where $\widetilde{\vec{\mathbf{F}}}$ is a modified $\vec{\mathbf{F}}$ such
1037 adcroft 1.4 that $\widetilde{\vec{\mathbf{F}}}.\widehat{n}=0$. As is implied by (\ref
1038 cnh 1.1 {eq:inhom-neumann-nh}) the modified boundary condition becomes:
1039    
1040     \begin{equation}
1041     \widehat{n}.\nabla \phi _{nh}=0 \label{eq:hom-neumann-nh}
1042     \end{equation}
1043    
1044     If the flow is `close' to hydrostatic balance then the 3-d inversion
1045     converges rapidly because $\phi _{nh}\ $is then only a small correction to
1046     the hydrostatic pressure field (see the discussion in Marshall et al, a,b).
1047    
1048 cnh 1.8 The solution $\phi _{nh}\ $to (\ref{eq:3d-invert}) and (\ref{eq:inhom-neumann-nh})
1049 cnh 1.1 does not vanish at $r=R_{moving}$, and so refines the pressure there.
1050    
1051     \subsection{Forcing/dissipation}
1052    
1053     \subsubsection{Forcing}
1054    
1055     The forcing terms $\mathcal{F}$ on the rhs of the equations are provided by
1056 cnh 1.8 `physics packages' and forcing packages. These are described later on.
1057 cnh 1.1
1058     \subsubsection{Dissipation}
1059    
1060     \paragraph{Momentum}
1061    
1062     Many forms of momentum dissipation are available in the model. Laplacian and
1063     biharmonic frictions are commonly used:
1064    
1065     \begin{equation}
1066 adcroft 1.4 D_{V}=A_{h}\nabla _{h}^{2}v+A_{v}\frac{\partial ^{2}v}{\partial z^{2}}
1067 cnh 1.1 +A_{4}\nabla _{h}^{4}v \label{eq:dissipation}
1068     \end{equation}
1069     where $A_{h}$ and $A_{v}\ $are (constant) horizontal and vertical viscosity
1070     coefficients and $A_{4}\ $is the horizontal coefficient for biharmonic
1071     friction. These coefficients are the same for all velocity components.
1072    
1073     \paragraph{Tracers}
1074    
1075     The mixing terms for the temperature and salinity equations have a similar
1076     form to that of momentum except that the diffusion tensor can be
1077 adcroft 1.4 non-diagonal and have varying coefficients. $\qquad $
1078 cnh 1.1 \begin{equation}
1079     D_{T,S}=\nabla .[\underline{\underline{K}}\nabla (T,S)]+K_{4}\nabla
1080     _{h}^{4}(T,S) \label{eq:diffusion}
1081     \end{equation}
1082 adcroft 1.4 where $\underline{\underline{K}}\ $is the diffusion tensor and the $K_{4}\ $
1083 cnh 1.1 horizontal coefficient for biharmonic diffusion. In the simplest case where
1084     the subgrid-scale fluxes of heat and salt are parameterized with constant
1085     horizontal and vertical diffusion coefficients, $\underline{\underline{K}}$,
1086     reduces to a diagonal matrix with constant coefficients:
1087    
1088     \begin{equation}
1089     \qquad \qquad \qquad \qquad K=\left(
1090     \begin{array}{ccc}
1091     K_{h} & 0 & 0 \\
1092     0 & K_{h} & 0 \\
1093 adcroft 1.4 0 & 0 & K_{v}
1094 cnh 1.1 \end{array}
1095     \right) \qquad \qquad \qquad \label{eq:diagonal-diffusion-tensor}
1096     \end{equation}
1097     where $K_{h}\ $and $K_{v}\ $are the horizontal and vertical diffusion
1098     coefficients. These coefficients are the same for all tracers (temperature,
1099     salinity ... ).
1100    
1101     \subsection{Vector invariant form}
1102    
1103 edhill 1.21 For some purposes it is advantageous to write momentum advection in
1104     eq(\ref {eq:horizontal_mtm}) and (\ref{eq:vertical_mtm}) in the
1105     (so-called) `vector invariant' form:
1106 cnh 1.1
1107     \begin{equation}
1108 adcroft 1.4 \frac{D\vec{\mathbf{v}}}{Dt}=\frac{\partial \vec{\mathbf{v}}}{\partial t}
1109     +\left( \nabla \times \vec{\mathbf{v}}\right) \times \vec{\mathbf{v}}+\nabla
1110 cnh 1.2 \left[ \frac{1}{2}(\vec{\mathbf{v}}\cdot \vec{\mathbf{v}})\right]
1111 cnh 1.1 \label{eq:vi-identity}
1112 adcroft 1.4 \end{equation}
1113 cnh 1.1 This permits alternative numerical treatments of the non-linear terms based
1114     on their representation as a vorticity flux. Because gradients of coordinate
1115     vectors no longer appear on the rhs of (\ref{eq:vi-identity}), explicit
1116 adcroft 1.4 representation of the metric terms in (\ref{eq:gu-speherical}), (\ref
1117 cnh 1.1 {eq:gv-spherical}) and (\ref{eq:gw-spherical}), can be avoided: information
1118     about the geometry is contained in the areas and lengths of the volumes used
1119     to discretize the model.
1120    
1121     \subsection{Adjoint}
1122    
1123 cnh 1.8 Tangent linear and adjoint counterparts of the forward model are described
1124 cnh 1.2 in Chapter 5.
1125 cnh 1.1
1126 jmc 1.22 % $Header: /u/gcmpack/manual/part1/manual.tex,v 1.21 2004/10/16 03:40:12 edhill Exp $
1127 cnh 1.1 % $Name: $
1128    
1129     \section{Appendix ATMOSPHERE}
1130    
1131     \subsection{Hydrostatic Primitive Equations for the Atmosphere in pressure
1132     coordinates}
1133    
1134     \label{sect-hpe-p}
1135    
1136     The hydrostatic primitive equations (HPEs) in p-coordinates are:
1137     \begin{eqnarray}
1138 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1139 cnh 1.2 _{h}+\mathbf{\nabla }_{p}\phi &=&\vec{\mathbf{\mathcal{F}}}
1140 cnh 1.1 \label{eq:atmos-mom} \\
1141 cnh 1.2 \frac{\partial \phi }{\partial p}+\alpha &=&0 \label{eq-p-hydro-start} \\
1142 adcroft 1.4 \mathbf{\nabla }_{p}\cdot \vec{\mathbf{v}}_{h}+\frac{\partial \omega }{
1143 cnh 1.1 \partial p} &=&0 \label{eq:atmos-cont} \\
1144 cnh 1.2 p\alpha &=&RT \label{eq:atmos-eos} \\
1145 cnh 1.1 c_{v}\frac{DT}{Dt}+p\frac{D\alpha }{Dt} &=&\mathcal{Q} \label{eq:atmos-heat}
1146 adcroft 1.4 \end{eqnarray}
1147 cnh 1.1 where $\vec{\mathbf{v}}_{h}=(u,v,0)$ is the `horizontal' (on pressure
1148     surfaces) component of velocity,$\frac{D}{Dt}=\vec{\mathbf{v}}_{h}\cdot
1149     \mathbf{\nabla }_{p}+\omega \frac{\partial }{\partial p}$ is the total
1150 cnh 1.6 derivative, $f=2\Omega \sin \varphi$ is the Coriolis parameter, $\phi =gz$ is
1151 adcroft 1.4 the geopotential, $\alpha =1/\rho $ is the specific volume, $\omega =\frac{Dp
1152     }{Dt}$ is the vertical velocity in the $p-$coordinate. Equation(\ref
1153     {eq:atmos-heat}) is the first law of thermodynamics where internal energy $
1154     e=c_{v}T$, $T$ is temperature, $Q$ is the rate of heating per unit mass and $
1155 cnh 1.1 p\frac{D\alpha }{Dt}$ is the work done by the fluid in compressing.
1156    
1157     It is convenient to cast the heat equation in terms of potential temperature
1158     $\theta $ so that it looks more like a generic conservation law.
1159     Differentiating (\ref{eq:atmos-eos}) we get:
1160     \begin{equation*}
1161     p\frac{D\alpha }{Dt}+\alpha \frac{Dp}{Dt}=R\frac{DT}{Dt}
1162 adcroft 1.4 \end{equation*}
1163     which, when added to the heat equation (\ref{eq:atmos-heat}) and using $
1164 cnh 1.1 c_{p}=c_{v}+R$, gives:
1165     \begin{equation}
1166     c_{p}\frac{DT}{Dt}-\alpha \frac{Dp}{Dt}=\mathcal{Q}
1167     \label{eq-p-heat-interim}
1168 adcroft 1.4 \end{equation}
1169 cnh 1.1 Potential temperature is defined:
1170     \begin{equation}
1171     \theta =T(\frac{p_{c}}{p})^{\kappa } \label{eq:potential-temp}
1172 adcroft 1.4 \end{equation}
1173 cnh 1.1 where $p_{c}$ is a reference pressure and $\kappa =R/c_{p}$. For convenience
1174     we will make use of the Exner function $\Pi (p)$ which defined by:
1175     \begin{equation}
1176     \Pi (p)=c_{p}(\frac{p}{p_{c}})^{\kappa } \label{Exner}
1177 adcroft 1.4 \end{equation}
1178 cnh 1.1 The following relations will be useful and are easily expressed in terms of
1179     the Exner function:
1180     \begin{equation*}
1181     c_{p}T=\Pi \theta \;\;;\;\;\frac{\partial \Pi }{\partial p}=\frac{\kappa \Pi
1182 adcroft 1.4 }{p}\;\;;\;\;\alpha =\frac{\kappa \Pi \theta }{p}=\frac{\partial \ \Pi }{
1183     \partial p}\theta \;\;;\;\;\frac{D\Pi }{Dt}=\frac{\partial \Pi }{\partial p}
1184 cnh 1.1 \frac{Dp}{Dt}
1185 adcroft 1.4 \end{equation*}
1186 cnh 1.1 where $b=\frac{\partial \ \Pi }{\partial p}\theta $ is the buoyancy.
1187    
1188     The heat equation is obtained by noting that
1189     \begin{equation*}
1190     c_{p}\frac{DT}{Dt}=\frac{D(\Pi \theta )}{Dt}=\Pi \frac{D\theta }{Dt}+\theta
1191 cnh 1.2 \frac{D\Pi }{Dt}=\Pi \frac{D\theta }{Dt}+\alpha \frac{Dp}{Dt}
1192 cnh 1.1 \end{equation*}
1193     and on substituting into (\ref{eq-p-heat-interim}) gives:
1194     \begin{equation}
1195     \Pi \frac{D\theta }{Dt}=\mathcal{Q}
1196     \label{eq:potential-temperature-equation}
1197     \end{equation}
1198     which is in conservative form.
1199    
1200 adcroft 1.4 For convenience in the model we prefer to step forward (\ref
1201 cnh 1.1 {eq:potential-temperature-equation}) rather than (\ref{eq:atmos-heat}).
1202    
1203     \subsubsection{Boundary conditions}
1204    
1205     The upper and lower boundary conditions are :
1206     \begin{eqnarray}
1207     \mbox{at the top:}\;\;p=0 &&\text{, }\omega =\frac{Dp}{Dt}=0 \\
1208     \mbox{at the surface:}\;\;p=p_{s} &&\text{, }\phi =\phi _{topo}=g~Z_{topo}
1209     \label{eq:boundary-condition-atmosphere}
1210     \end{eqnarray}
1211     In $p$-coordinates, the upper boundary acts like a solid boundary ($\omega
1212     =0 $); in $z$-coordinates and the lower boundary is analogous to a free
1213     surface ($\phi $ is imposed and $\omega \neq 0$).
1214    
1215     \subsubsection{Splitting the geo-potential}
1216 jmc 1.22 \label{sec:hpe-p-geo-potential-split}
1217 cnh 1.1
1218     For the purposes of initialization and reducing round-off errors, the model
1219     deals with perturbations from reference (or ``standard'') profiles. For
1220     example, the hydrostatic geopotential associated with the resting atmosphere
1221     is not dynamically relevant and can therefore be subtracted from the
1222     equations. The equations written in terms of perturbations are obtained by
1223     substituting the following definitions into the previous model equations:
1224     \begin{eqnarray}
1225     \theta &=&\theta _{o}+\theta ^{\prime } \label{eq:atmos-ref-prof-theta} \\
1226     \alpha &=&\alpha _{o}+\alpha ^{\prime } \label{eq:atmos-ref-prof-alpha} \\
1227     \phi &=&\phi _{o}+\phi ^{\prime } \label{eq:atmos-ref-prof-phi}
1228     \end{eqnarray}
1229     The reference state (indicated by subscript ``0'') corresponds to
1230     horizontally homogeneous atmosphere at rest ($\theta _{o},\alpha _{o},\phi
1231     _{o}$) with surface pressure $p_{o}(x,y)$ that satisfies $\phi
1232     _{o}(p_{o})=g~Z_{topo}$, defined:
1233     \begin{eqnarray*}
1234     \theta _{o}(p) &=&f^{n}(p) \\
1235     \alpha _{o}(p) &=&\Pi _{p}\theta _{o} \\
1236     \phi _{o}(p) &=&\phi _{topo}-\int_{p_{0}}^{p}\alpha _{o}dp
1237     \end{eqnarray*}
1238     %\begin{eqnarray*}
1239     %\phi'_\alpha & = & \int^p_{p_o} (\alpha_o -\alpha) dp \\
1240     %\phi'_s(x,y,t) & = & \int_{p_o}^{p_s} \alpha dp
1241     %\end{eqnarray*}
1242    
1243     The final form of the HPE's in p coordinates is then:
1244     \begin{eqnarray}
1245 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1246 edhill 1.21 _{h}+\mathbf{\nabla }_{p}\phi ^{\prime } &=&\vec{\mathbf{\mathcal{F}}}
1247     \label{eq:atmos-prime} \\
1248 cnh 1.1 \frac{\partial \phi ^{\prime }}{\partial p}+\alpha ^{\prime } &=&0 \\
1249 adcroft 1.4 \mathbf{\nabla }_{p}\cdot \vec{\mathbf{v}}_{h}+\frac{\partial \omega }{
1250 cnh 1.1 \partial p} &=&0 \\
1251     \frac{\partial \Pi }{\partial p}\theta ^{\prime } &=&\alpha ^{\prime } \\
1252 cnh 1.8 \frac{D\theta }{Dt} &=&\frac{\mathcal{Q}}{\Pi }
1253 cnh 1.1 \end{eqnarray}
1254    
1255 jmc 1.22 % $Header: /u/gcmpack/manual/part1/manual.tex,v 1.21 2004/10/16 03:40:12 edhill Exp $
1256 cnh 1.1 % $Name: $
1257    
1258     \section{Appendix OCEAN}
1259    
1260     \subsection{Equations of motion for the ocean}
1261    
1262     We review here the method by which the standard (Boussinesq, incompressible)
1263     HPE's for the ocean written in z-coordinates are obtained. The
1264     non-Boussinesq equations for oceanic motion are:
1265     \begin{eqnarray}
1266 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1267 cnh 1.1 _{h}+\frac{1}{\rho }\mathbf{\nabla }_{z}p &=&\vec{\mathbf{\mathcal{F}}} \\
1268     \epsilon _{nh}\frac{Dw}{Dt}+g+\frac{1}{\rho }\frac{\partial p}{\partial z}
1269     &=&\epsilon _{nh}\mathcal{F}_{w} \\
1270 adcroft 1.4 \frac{1}{\rho }\frac{D\rho }{Dt}+\mathbf{\nabla }_{z}\cdot \vec{\mathbf{v}}
1271 cnh 1.8 _{h}+\frac{\partial w}{\partial z} &=&0 \label{eq-zns-cont}\\
1272     \rho &=&\rho (\theta ,S,p) \label{eq-zns-eos}\\
1273     \frac{D\theta }{Dt} &=&\mathcal{Q}_{\theta } \label{eq-zns-heat}\\
1274     \frac{DS}{Dt} &=&\mathcal{Q}_{s} \label{eq-zns-salt}
1275     \label{eq:non-boussinesq}
1276 adcroft 1.4 \end{eqnarray}
1277 cnh 1.1 These equations permit acoustics modes, inertia-gravity waves,
1278 cnh 1.10 non-hydrostatic motions, a geostrophic (Rossby) mode and a thermohaline
1279 cnh 1.1 mode. As written, they cannot be integrated forward consistently - if we
1280     step $\rho $ forward in (\ref{eq-zns-cont}), the answer will not be
1281 adcroft 1.4 consistent with that obtained by stepping (\ref{eq-zns-heat}) and (\ref
1282 cnh 1.1 {eq-zns-salt}) and then using (\ref{eq-zns-eos}) to yield $\rho $. It is
1283     therefore necessary to manipulate the system as follows. Differentiating the
1284     EOS (equation of state) gives:
1285    
1286     \begin{equation}
1287     \frac{D\rho }{Dt}=\left. \frac{\partial \rho }{\partial \theta }\right|
1288     _{S,p}\frac{D\theta }{Dt}+\left. \frac{\partial \rho }{\partial S}\right|
1289     _{\theta ,p}\frac{DS}{Dt}+\left. \frac{\partial \rho }{\partial p}\right|
1290     _{\theta ,S}\frac{Dp}{Dt} \label{EOSexpansion}
1291     \end{equation}
1292    
1293 edhill 1.21 Note that $\frac{\partial \rho }{\partial p}=\frac{1}{c_{s}^{2}}$ is
1294     the reciprocal of the sound speed ($c_{s}$) squared. Substituting into
1295     \ref{eq-zns-cont} gives:
1296 cnh 1.1 \begin{equation}
1297 adcroft 1.4 \frac{1}{\rho c_{s}^{2}}\frac{Dp}{Dt}+\mathbf{\nabla }_{z}\cdot \vec{\mathbf{
1298 cnh 1.1 v}}+\partial _{z}w\approx 0 \label{eq-zns-pressure}
1299     \end{equation}
1300     where we have used an approximation sign to indicate that we have assumed
1301     adiabatic motion, dropping the $\frac{D\theta }{Dt}$ and $\frac{DS}{Dt}$.
1302     Replacing \ref{eq-zns-cont} with \ref{eq-zns-pressure} yields a system that
1303     can be explicitly integrated forward:
1304     \begin{eqnarray}
1305 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1306 cnh 1.1 _{h}+\frac{1}{\rho }\mathbf{\nabla }_{z}p &=&\vec{\mathbf{\mathcal{F}}}
1307     \label{eq-cns-hmom} \\
1308     \epsilon _{nh}\frac{Dw}{Dt}+g+\frac{1}{\rho }\frac{\partial p}{\partial z}
1309     &=&\epsilon _{nh}\mathcal{F}_{w} \label{eq-cns-hydro} \\
1310 adcroft 1.4 \frac{1}{\rho c_{s}^{2}}\frac{Dp}{Dt}+\mathbf{\nabla }_{z}\cdot \vec{\mathbf{
1311 cnh 1.1 v}}_{h}+\frac{\partial w}{\partial z} &=&0 \label{eq-cns-cont} \\
1312     \rho &=&\rho (\theta ,S,p) \label{eq-cns-eos} \\
1313     \frac{D\theta }{Dt} &=&\mathcal{Q}_{\theta } \label{eq-cns-heat} \\
1314     \frac{DS}{Dt} &=&\mathcal{Q}_{s} \label{eq-cns-salt}
1315     \end{eqnarray}
1316    
1317     \subsubsection{Compressible z-coordinate equations}
1318    
1319     Here we linearize the acoustic modes by replacing $\rho $ with $\rho _{o}(z)$
1320     wherever it appears in a product (ie. non-linear term) - this is the
1321     `Boussinesq assumption'. The only term that then retains the full variation
1322     in $\rho $ is the gravitational acceleration:
1323     \begin{eqnarray}
1324 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1325 cnh 1.1 _{h}+\frac{1}{\rho _{o}}\mathbf{\nabla }_{z}p &=&\vec{\mathbf{\mathcal{F}}}
1326     \label{eq-zcb-hmom} \\
1327 adcroft 1.4 \epsilon _{nh}\frac{Dw}{Dt}+\frac{g\rho }{\rho _{o}}+\frac{1}{\rho _{o}}
1328 cnh 1.1 \frac{\partial p}{\partial z} &=&\epsilon _{nh}\mathcal{F}_{w}
1329     \label{eq-zcb-hydro} \\
1330 adcroft 1.4 \frac{1}{\rho _{o}c_{s}^{2}}\frac{Dp}{Dt}+\mathbf{\nabla }_{z}\cdot \vec{
1331 cnh 1.1 \mathbf{v}}_{h}+\frac{\partial w}{\partial z} &=&0 \label{eq-zcb-cont} \\
1332     \rho &=&\rho (\theta ,S,p) \label{eq-zcb-eos} \\
1333     \frac{D\theta }{Dt} &=&\mathcal{Q}_{\theta } \label{eq-zcb-heat} \\
1334     \frac{DS}{Dt} &=&\mathcal{Q}_{s} \label{eq-zcb-salt}
1335     \end{eqnarray}
1336     These equations still retain acoustic modes. But, because the
1337 adcroft 1.4 ``compressible'' terms are linearized, the pressure equation \ref
1338 cnh 1.1 {eq-zcb-cont} can be integrated implicitly with ease (the time-dependent
1339     term appears as a Helmholtz term in the non-hydrostatic pressure equation).
1340     These are the \emph{truly} compressible Boussinesq equations. Note that the
1341     EOS must have the same pressure dependency as the linearized pressure term,
1342 adcroft 1.4 ie. $\left. \frac{\partial \rho }{\partial p}\right| _{\theta ,S}=\frac{1}{
1343 cnh 1.1 c_{s}^{2}}$, for consistency.
1344    
1345     \subsubsection{`Anelastic' z-coordinate equations}
1346    
1347     The anelastic approximation filters the acoustic mode by removing the
1348 adcroft 1.4 time-dependency in the continuity (now pressure-) equation (\ref{eq-zcb-cont}
1349     ). This could be done simply by noting that $\frac{Dp}{Dt}\approx -g\rho _{o}
1350 cnh 1.1 \frac{Dz}{Dt}=-g\rho _{o}w$, but this leads to an inconsistency between
1351     continuity and EOS. A better solution is to change the dependency on
1352     pressure in the EOS by splitting the pressure into a reference function of
1353     height and a perturbation:
1354     \begin{equation*}
1355 cnh 1.2 \rho =\rho (\theta ,S,p_{o}(z)+\epsilon _{s}p^{\prime })
1356 cnh 1.1 \end{equation*}
1357     Remembering that the term $\frac{Dp}{Dt}$ in continuity comes from
1358     differentiating the EOS, the continuity equation then becomes:
1359     \begin{equation*}
1360 adcroft 1.4 \frac{1}{\rho _{o}c_{s}^{2}}\left( \frac{Dp_{o}}{Dt}+\epsilon _{s}\frac{
1361     Dp^{\prime }}{Dt}\right) +\mathbf{\nabla }_{z}\cdot \vec{\mathbf{v}}_{h}+
1362 cnh 1.2 \frac{\partial w}{\partial z}=0
1363 cnh 1.1 \end{equation*}
1364     If the time- and space-scales of the motions of interest are longer than
1365 adcroft 1.4 those of acoustic modes, then $\frac{Dp^{\prime }}{Dt}<<(\frac{Dp_{o}}{Dt},
1366 cnh 1.1 \mathbf{\nabla }\cdot \vec{\mathbf{v}}_{h})$ in the continuity equations and
1367 adcroft 1.4 $\left. \frac{\partial \rho }{\partial p}\right| _{\theta ,S}\frac{
1368 cnh 1.1 Dp^{\prime }}{Dt}<<\left. \frac{\partial \rho }{\partial p}\right| _{\theta
1369     ,S}\frac{Dp_{o}}{Dt}$ in the EOS (\ref{EOSexpansion}). Thus we set $\epsilon
1370     _{s}=0$, removing the dependency on $p^{\prime }$ in the continuity equation
1371     and EOS. Expanding $\frac{Dp_{o}(z)}{Dt}=-g\rho _{o}w$ then leads to the
1372     anelastic continuity equation:
1373     \begin{equation}
1374 adcroft 1.4 \mathbf{\nabla }_{z}\cdot \vec{\mathbf{v}}_{h}+\frac{\partial w}{\partial z}-
1375 cnh 1.1 \frac{g}{c_{s}^{2}}w=0 \label{eq-za-cont1}
1376     \end{equation}
1377     A slightly different route leads to the quasi-Boussinesq continuity equation
1378 adcroft 1.4 where we use the scaling $\frac{\partial \rho ^{\prime }}{\partial t}+
1379     \mathbf{\nabla }_{3}\cdot \rho ^{\prime }\vec{\mathbf{v}}<<\mathbf{\nabla }
1380 cnh 1.1 _{3}\cdot \rho _{o}\vec{\mathbf{v}}$ yielding:
1381     \begin{equation}
1382 adcroft 1.4 \mathbf{\nabla }_{z}\cdot \vec{\mathbf{v}}_{h}+\frac{1}{\rho _{o}}\frac{
1383 cnh 1.1 \partial \left( \rho _{o}w\right) }{\partial z}=0 \label{eq-za-cont2}
1384     \end{equation}
1385     Equations \ref{eq-za-cont1} and \ref{eq-za-cont2} are in fact the same
1386     equation if:
1387     \begin{equation}
1388     \frac{1}{\rho _{o}}\frac{\partial \rho _{o}}{\partial z}=\frac{-g}{c_{s}^{2}}
1389     \end{equation}
1390     Again, note that if $\rho _{o}$ is evaluated from prescribed $\theta _{o}$
1391 adcroft 1.4 and $S_{o}$ profiles, then the EOS dependency on $p_{o}$ and the term $\frac{
1392 cnh 1.1 g}{c_{s}^{2}}$ in continuity should be referred to those same profiles. The
1393     full set of `quasi-Boussinesq' or `anelastic' equations for the ocean are
1394     then:
1395     \begin{eqnarray}
1396 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1397 cnh 1.1 _{h}+\frac{1}{\rho _{o}}\mathbf{\nabla }_{z}p &=&\vec{\mathbf{\mathcal{F}}}
1398     \label{eq-zab-hmom} \\
1399 adcroft 1.4 \epsilon _{nh}\frac{Dw}{Dt}+\frac{g\rho }{\rho _{o}}+\frac{1}{\rho _{o}}
1400 cnh 1.1 \frac{\partial p}{\partial z} &=&\epsilon _{nh}\mathcal{F}_{w}
1401     \label{eq-zab-hydro} \\
1402 adcroft 1.4 \mathbf{\nabla }_{z}\cdot \vec{\mathbf{v}}_{h}+\frac{1}{\rho _{o}}\frac{
1403 cnh 1.1 \partial \left( \rho _{o}w\right) }{\partial z} &=&0 \label{eq-zab-cont} \\
1404     \rho &=&\rho (\theta ,S,p_{o}(z)) \label{eq-zab-eos} \\
1405     \frac{D\theta }{Dt} &=&\mathcal{Q}_{\theta } \label{eq-zab-heat} \\
1406     \frac{DS}{Dt} &=&\mathcal{Q}_{s} \label{eq-zab-salt}
1407     \end{eqnarray}
1408    
1409     \subsubsection{Incompressible z-coordinate equations}
1410    
1411     Here, the objective is to drop the depth dependence of $\rho _{o}$ and so,
1412     technically, to also remove the dependence of $\rho $ on $p_{o}$. This would
1413     yield the ``truly'' incompressible Boussinesq equations:
1414     \begin{eqnarray}
1415 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1416 cnh 1.1 _{h}+\frac{1}{\rho _{c}}\mathbf{\nabla }_{z}p &=&\vec{\mathbf{\mathcal{F}}}
1417     \label{eq-ztb-hmom} \\
1418 adcroft 1.4 \epsilon _{nh}\frac{Dw}{Dt}+\frac{g\rho }{\rho _{c}}+\frac{1}{\rho _{c}}
1419 cnh 1.1 \frac{\partial p}{\partial z} &=&\epsilon _{nh}\mathcal{F}_{w}
1420     \label{eq-ztb-hydro} \\
1421     \mathbf{\nabla }_{z}\cdot \vec{\mathbf{v}}_{h}+\frac{\partial w}{\partial z}
1422     &=&0 \label{eq-ztb-cont} \\
1423     \rho &=&\rho (\theta ,S) \label{eq-ztb-eos} \\
1424     \frac{D\theta }{Dt} &=&\mathcal{Q}_{\theta } \label{eq-ztb-heat} \\
1425     \frac{DS}{Dt} &=&\mathcal{Q}_{s} \label{eq-ztb-salt}
1426     \end{eqnarray}
1427     where $\rho _{c}$ is a constant reference density of water.
1428    
1429     \subsubsection{Compressible non-divergent equations}
1430    
1431     The above ``incompressible'' equations are incompressible in both the flow
1432     and the density. In many oceanic applications, however, it is important to
1433     retain compressibility effects in the density. To do this we must split the
1434     density thus:
1435     \begin{equation*}
1436     \rho =\rho _{o}+\rho ^{\prime }
1437 adcroft 1.4 \end{equation*}
1438 cnh 1.1 We then assert that variations with depth of $\rho _{o}$ are unimportant
1439     while the compressible effects in $\rho ^{\prime }$ are:
1440     \begin{equation*}
1441     \rho _{o}=\rho _{c}
1442 adcroft 1.4 \end{equation*}
1443 cnh 1.1 \begin{equation*}
1444     \rho ^{\prime }=\rho (\theta ,S,p_{o}(z))-\rho _{o}
1445 adcroft 1.4 \end{equation*}
1446 cnh 1.1 This then yields what we can call the semi-compressible Boussinesq
1447     equations:
1448     \begin{eqnarray}
1449 adcroft 1.4 \frac{D\vec{\mathbf{v}}_{h}}{Dt}+f\hat{\mathbf{k}}\times \vec{\mathbf{v}}
1450     _{h}+\frac{1}{\rho _{c}}\mathbf{\nabla }_{z}p^{\prime } &=&\vec{\mathbf{
1451 cnh 1.1 \mathcal{F}}} \label{eq:ocean-mom} \\
1452     \epsilon _{nh}\frac{Dw}{Dt}+\frac{g\rho ^{\prime }}{\rho _{c}}+\frac{1}{\rho
1453     _{c}}\frac{\partial p^{\prime }}{\partial z} &=&\epsilon _{nh}\mathcal{F}_{w}
1454     \label{eq:ocean-wmom} \\
1455     \mathbf{\nabla }_{z}\cdot \vec{\mathbf{v}}_{h}+\frac{\partial w}{\partial z}
1456     &=&0 \label{eq:ocean-cont} \\
1457     \rho ^{\prime } &=&\rho (\theta ,S,p_{o}(z))-\rho _{c} \label{eq:ocean-eos}
1458     \\
1459     \frac{D\theta }{Dt} &=&\mathcal{Q}_{\theta } \label{eq:ocean-theta} \\
1460     \frac{DS}{Dt} &=&\mathcal{Q}_{s} \label{eq:ocean-salt}
1461 adcroft 1.4 \end{eqnarray}
1462 cnh 1.1 Note that the hydrostatic pressure of the resting fluid, including that
1463     associated with $\rho _{c}$, is subtracted out since it has no effect on the
1464     dynamics.
1465    
1466     Though necessary, the assumptions that go into these equations are messy
1467     since we essentially assume a different EOS for the reference density and
1468     the perturbation density. Nevertheless, it is the hydrostatic ($\epsilon
1469     _{nh}=0$ form of these equations that are used throughout the ocean modeling
1470     community and referred to as the primitive equations (HPE).
1471    
1472 jmc 1.22 % $Header: /u/gcmpack/manual/part1/manual.tex,v 1.21 2004/10/16 03:40:12 edhill Exp $
1473 cnh 1.1 % $Name: $
1474    
1475     \section{Appendix:OPERATORS}
1476    
1477     \subsection{Coordinate systems}
1478    
1479     \subsubsection{Spherical coordinates}
1480    
1481     In spherical coordinates, the velocity components in the zonal, meridional
1482     and vertical direction respectively, are given by (see Fig.2) :
1483    
1484     \begin{equation*}
1485 cnh 1.6 u=r\cos \varphi \frac{D\lambda }{Dt}
1486 cnh 1.1 \end{equation*}
1487    
1488     \begin{equation*}
1489 cnh 1.6 v=r\frac{D\varphi }{Dt}\qquad
1490 cnh 1.1 \end{equation*}
1491     $\qquad \qquad \qquad \qquad $
1492    
1493     \begin{equation*}
1494 cnh 1.2 \dot{r}=\frac{Dr}{Dt}
1495 cnh 1.1 \end{equation*}
1496    
1497 cnh 1.6 Here $\varphi $ is the latitude, $\lambda $ the longitude, $r$ the radial
1498 cnh 1.1 distance of the particle from the center of the earth, $\Omega $ is the
1499     angular speed of rotation of the Earth and $D/Dt$ is the total derivative.
1500    
1501     The `grad' ($\nabla $) and `div' ($\nabla $.) operators are defined by, in
1502     spherical coordinates:
1503    
1504     \begin{equation*}
1505 cnh 1.6 \nabla \equiv \left( \frac{1}{r\cos \varphi }\frac{\partial }{\partial \lambda }
1506     ,\frac{1}{r}\frac{\partial }{\partial \varphi },\frac{\partial }{\partial r}
1507 cnh 1.2 \right)
1508 cnh 1.1 \end{equation*}
1509    
1510     \begin{equation*}
1511 cnh 1.6 \nabla .v\equiv \frac{1}{r\cos \varphi }\left\{ \frac{\partial u}{\partial
1512     \lambda }+\frac{\partial }{\partial \varphi }\left( v\cos \varphi \right) \right\}
1513 cnh 1.2 +\frac{1}{r^{2}}\frac{\partial \left( r^{2}\dot{r}\right) }{\partial r}
1514 cnh 1.1 \end{equation*}
1515    
1516 adcroft 1.4 %tci%\end{document}

  ViewVC Help
Powered by ViewVC 1.1.22