/[MITgcm]/manual/s_examples/tracer_adjsens/co2sens.tex
ViewVC logotype

Contents of /manual/s_examples/tracer_adjsens/co2sens.tex

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


Revision 1.4 - (show annotations) (download) (as text)
Tue Nov 13 18:22:24 2001 UTC (23 years, 8 months ago) by adcroft
Branch: MAIN
Changes since 1.3: +2 -2 lines
File MIME type: application/x-tex
Change cite to ref

1 % $Header: /u/gcmpack/mitgcmdoc/part3/case_studies/carbon_outgassing_sensitivity/co2sens.tex,v 1.3 2001/11/13 18:19:18 adcroft Exp $
2 % $Name: $
3
4 \section{Example: Centennial Time Scale Sensitivities}
5
6 \bodytext{bgcolor="#FFFFFFFF"}
7
8 %\begin{center}
9 %{\Large \bf Using MITgcm to Look at Centennial Time Scale
10 %Sensitivities}
11 %
12 %\vspace*{4mm}
13 %
14 %\vspace*{3mm}
15 %{\large May 2001}
16 %\end{center}
17
18 \subsection{Introduction}
19
20 This document describes the fourth example MITgcm experiment.
21 This example illustrates the use of
22 the MITgcm to perform sensitivity analysis in a
23 large scale ocean circulation simulation.
24
25 \subsection{Overview}
26
27 This example experiment demonstrates using the MITgcm to simulate
28 the planetary ocean circulation. The simulation is configured
29 with realistic geography and bathymetry on a
30 $4^{\circ} \times 4^{\circ}$ spherical polar grid.
31 Twenty vertical layers are used in the vertical, ranging in thickness
32 from $50\,{\rm m}$ at the surface to $815\,{\rm m}$ at depth,
33 giving a maximum model depth of $6\,{\rm km}$.
34 At this resolution, the configuration
35 can be integrated forward for thousands of years on a single
36 processor desktop computer.
37 \\
38
39 The model is forced with climatological wind stress data and surface
40 flux data from Da Silva \cite{DaSilva94}. Climatological data
41 from Levitus \cite{Levitus94} is used to initialize the model hydrography.
42 Levitus data is also used throughout the calculation
43 to derive air-sea fluxes of heat at the ocean surface.
44 These fluxes are combined with climatological estimates of
45 surface heat flux and fresh water, resulting in a mixed boundary
46 condition of the style described in Haney \cite{Haney}.
47 Altogether, this yields the following forcing applied
48 in the model surface layer.
49
50 \begin{eqnarray}
51 \label{EQ:global_forcing}
52 \label{EQ:global_forcing_fu}
53 {\cal F}_{u} & = & \frac{\tau_{x}}{\rho_{0} \Delta z_{s}}
54 \\
55 \label{EQ:global_forcing_fv}
56 {\cal F}_{v} & = & \frac{\tau_{y}}{\rho_{0} \Delta z_{s}}
57 \\
58 \label{EQ:global_forcing_ft}
59 {\cal F}_{\theta} & = & - \lambda_{\theta} ( \theta - \theta^{\ast} )
60 - \frac{1}{C_{p} \rho_{0} \Delta z_{s}}{\cal Q}
61 \\
62 \label{EQ:global_forcing_fs}
63 {\cal F}_{s} & = & - \lambda_{s} ( S - S^{\ast} )
64 + \frac{S_{0}}{\Delta z_{s}}({\cal E} - {\cal P} - {\cal R})
65 \end{eqnarray}
66
67 \noindent where ${\cal F}_{u}$, ${\cal F}_{v}$, ${\cal F}_{\theta}$,
68 ${\cal F}_{s}$ are the forcing terms in the zonal and meridional
69 momentum and in the potential temperature and salinity
70 equations respectively.
71 The term $\Delta z_{s}$ represents the top ocean layer thickness.
72 It is used in conjunction with the reference density, $\rho_{0}$
73 (here set to $999.8\,{\rm kg\,m^{-3}}$), the
74 reference salinity, $S_{0}$ (here set to 35ppt),
75 and a specific heat capacity $C_{p}$ to convert
76 wind-stress fluxes given in ${\rm N}\,m^{-2}$,
77 \\
78
79
80 The configuration is illustrated in figure \ref{simulation_config}.
81
82
83 \subsection{Discrete Numerical Configuration}
84
85
86 The model is configured in hydrostatic form. The domain is discretised with
87 a uniform grid spacing in latitude and longitude of
88 $\Delta x=\Delta y=4^{\circ}$, so
89 that there are ninety grid cells in the $x$ and forty in the
90 $y$ direction (Arctic polar regions are not
91 included in this experiment). Vertically the
92 model is configured with twenty layers with the following thicknesses
93 $\Delta z_{1} = 50\,{\rm m},\,
94 \Delta z_{2} = 50\,{\rm m},\,
95 \Delta z_{3} = 55\,{\rm m},\,
96 \Delta z_{4} = 60\,{\rm m},\,
97 \Delta z_{5} = 65\,{\rm m},\,
98 $
99 $
100 \Delta z_{6}~=~70\,{\rm m},\,
101 \Delta z_{7}~=~80\,{\rm m},\,
102 \Delta z_{8}~=95\,{\rm m},\,
103 \Delta z_{9}=120\,{\rm m},\,
104 \Delta z_{10}=155\,{\rm m},\,
105 $
106 $
107 \Delta z_{11}=200\,{\rm m},\,
108 \Delta z_{12}=260\,{\rm m},\,
109 \Delta z_{13}=320\,{\rm m},\,
110 \Delta z_{14}=400\,{\rm m},\,
111 \Delta z_{15}=480\,{\rm m},\,
112 $
113 $
114 \Delta z_{16}=570\,{\rm m},\,
115 \Delta z_{17}=655\,{\rm m},\,
116 \Delta z_{18}=725\,{\rm m},\,
117 \Delta z_{19}=775\,{\rm m},\,
118 \Delta z_{20}=815\,{\rm m}
119 $ (here the numeric subscript indicates the model level index number, ${\tt k}$).
120 The implicit free surface form of the pressure equation described in Marshall et. al
121 \cite{Marshall97a} is employed. A Laplacian operator, $\nabla^2$, provides viscous
122 dissipation. Thermal and haline diffusion is also represented by a Laplacian operator.
123 \\
124
125 Wind-stress momentum inputs are added to the momentum equations for both
126 the zonal flow, $u$ and the meridional flow $v$, according to equations
127 (\ref{EQ:global_forcing_fu}) and (\ref{EQ:global_forcing_fv}).
128 Thermodynamic forcing inputs are added to the equations for
129 potential temperature, $\theta$, and salinity, $S$, according to equations
130 (\ref{EQ:global_forcing_ft}) and (\ref{EQ:global_forcing_fs}).
131 This produces a set of equations solved in this configuration as follows:
132 % {\fracktur}
133
134
135 \begin{eqnarray}
136 \label{EQ:model_equations}
137 \frac{Du}{Dt} - fv +
138 \frac{1}{\rho}\frac{\partial p^{'}}{\partial x} -
139 A_{h}\nabla_{h}^2u - A_{z}\frac{\partial^{2}u}{\partial z^{2}}
140 & = &
141 {\cal F}_{u}
142 \\
143 \frac{Dv}{Dt} + fu +
144 \frac{1}{\rho}\frac{\partial p^{'}}{\partial y} -
145 A_{h}\nabla_{h}^2v - A_{z}\frac{\partial^{2}v}{\partial z^{2}}
146 & = &
147 {\cal F}_{v}
148 \\
149 \frac{\partial \eta}{\partial t} + \nabla_{h}\cdot \vec{u}
150 &=&
151 0
152 \\
153 \frac{D\theta}{Dt} -
154 K_{h}\nabla_{h}^2\theta - \Gamma(K_{z})\frac{\partial^{2}\theta}{\partial z^{2}}
155 & = &
156 {\cal F}_{\theta}
157 \\
158 \frac{D s}{Dt} -
159 K_{h}\nabla_{h}^2 s - \Gamma(K_{z})\frac{\partial^{2} s}{\partial z^{2}}
160 & = &
161 {\cal F}_{s}
162 \\
163 g\rho_{0} \eta + \int^{0}_{-z}\rho^{'} dz & = & p^{'}
164 \\
165 \end{eqnarray}
166
167 \noindent where $u$ and $v$ are the $x$ and $y$ components of the
168 flow vector $\vec{u}$. The suffices ${s},{i}$ indicate surface and
169 interior model levels respectively. As described in
170 MITgcm Numerical Solution Procedure \ref{chap:discretization}, the time
171 evolution of potential temperature, $\theta$, equation is solved prognostically.
172 The total pressure, $p$, is diagnosed by summing pressure due to surface
173 elevation $\eta$ and the hydrostatic pressure.
174 \\
175
176 \subsubsection{Numerical Stability Criteria}
177
178 The Laplacian dissipation coefficient, $A_{h}$, is set to $400 m s^{-1}$.
179 This value is chosen to yield a Munk layer width \cite{adcroft:95},
180
181 \begin{eqnarray}
182 \label{EQ:munk_layer}
183 M_{w} = \pi ( \frac { A_{h} }{ \beta } )^{\frac{1}{3}}
184 \end{eqnarray}
185
186 \noindent of $\approx 100$km. This is greater than the model
187 resolution in mid-latitudes $\Delta x$, ensuring that the frictional
188 boundary layer is well resolved.
189 \\
190
191 \noindent The model is stepped forward with a
192 time step $\delta t=1200$secs. With this time step the stability
193 parameter to the horizontal Laplacian friction \cite{adcroft:95}
194
195 \begin{eqnarray}
196 \label{EQ:laplacian_stability}
197 S_{l} = 4 \frac{A_{h} \delta t}{{\Delta x}^2}
198 \end{eqnarray}
199
200 \noindent evaluates to 0.012, which is well below the 0.3 upper limit
201 for stability.
202 \\
203
204 \noindent The vertical dissipation coefficient, $A_{z}$, is set to
205 $1\times10^{-2} {\rm m}^2{\rm s}^{-1}$. The associated stability limit
206
207 \begin{eqnarray}
208 \label{EQ:laplacian_stability_z}
209 S_{l} = 4 \frac{A_{z} \delta t}{{\Delta z}^2}
210 \end{eqnarray}
211
212 \noindent evaluates to $4.8 \times 10^{-5}$ which is again well below
213 the upper limit.
214 The values of $A_{h}$ and $A_{z}$ are also used for the horizontal ($K_{h}$)
215 and vertical ($K_{z}$) diffusion coefficients for temperature respectively.
216 \\
217
218 \noindent The numerical stability for inertial oscillations
219 \cite{adcroft:95}
220
221 \begin{eqnarray}
222 \label{EQ:inertial_stability}
223 S_{i} = f^{2} {\delta t}^2
224 \end{eqnarray}
225
226 \noindent evaluates to $0.0144$, which is well below the $0.5$ upper
227 limit for stability.
228 \\
229
230 \noindent The advective CFL \cite{adcroft:95} for a extreme maximum
231 horizontal flow
232 speed of $ | \vec{u} | = 2 ms^{-1}$
233
234 \begin{eqnarray}
235 \label{EQ:cfl_stability}
236 S_{a} = \frac{| \vec{u} | \delta t}{ \Delta x}
237 \end{eqnarray}
238
239 \noindent evaluates to $5 \times 10^{-2}$. This is well below the stability
240 limit of 0.5.
241 \\
242
243 \noindent The stability parameter for internal gravity waves
244 \cite{adcroft:95}
245
246 \begin{eqnarray}
247 \label{EQ:cfl_stability}
248 S_{c} = \frac{c_{g} \delta t}{ \Delta x}
249 \end{eqnarray}
250
251 \noindent evaluates to $5 \times 10^{-2}$. This is well below the linear
252 stability limit of 0.25.
253
254 \subsection{Code Configuration}
255 \label{SEC:code_config}
256
257 The model configuration for this experiment resides under the
258 directory {\it verification/exp1/}. The experiment files
259 \begin{itemize}
260 \item {\it input/data}
261 \item {\it input/data.pkg}
262 \item {\it input/eedata},
263 \item {\it input/windx.sin\_y},
264 \item {\it input/topog.box},
265 \item {\it code/CPP\_EEOPTIONS.h}
266 \item {\it code/CPP\_OPTIONS.h},
267 \item {\it code/SIZE.h}.
268 \end{itemize}
269 contain the code customizations and parameter settings for this
270 experiments. Below we describe the customizations
271 to these files associated with this experiment.
272
273 \subsubsection{File {\it input/data}}
274
275 This file, reproduced completely below, specifies the main parameters
276 for the experiment. The parameters that are significant for this configuration
277 are
278
279 \begin{itemize}
280
281 \item Line 4,
282 \begin{verbatim} tRef=20.,10.,8.,6., \end{verbatim}
283 this line sets
284 the initial and reference values of potential temperature at each model
285 level in units of $^{\circ}$C.
286 The entries are ordered from surface to depth. For each
287 depth level the initial and reference profiles will be uniform in
288 $x$ and $y$.
289
290 \fbox{
291 \begin{minipage}{5.0in}
292 {\it S/R INI\_THETA}({\it ini\_theta.F})
293 \end{minipage}
294 }
295
296
297 \item Line 6,
298 \begin{verbatim} viscAz=1.E-2, \end{verbatim}
299 this line sets the vertical Laplacian dissipation coefficient to
300 $1 \times 10^{-2} {\rm m^{2}s^{-1}}$. Boundary conditions
301 for this operator are specified later. This variable is copied into
302 model general vertical coordinate variable {\bf viscAr}.
303
304 \fbox{
305 \begin{minipage}{5.0in}
306 {\it S/R CALC\_DIFFUSIVITY}({\it calc\_diffusivity.F})
307 \end{minipage}
308 }
309
310 \item Line 7,
311 \begin{verbatim}
312 viscAh=4.E2,
313 \end{verbatim}
314 this line sets the horizontal Laplacian frictional dissipation coefficient to
315 $1 \times 10^{-2} {\rm m^{2}s^{-1}}$. Boundary conditions
316 for this operator are specified later.
317
318 \item Lines 8,
319 \begin{verbatim}
320 no_slip_sides=.FALSE.
321 \end{verbatim}
322 this line selects a free-slip lateral boundary condition for
323 the horizontal Laplacian friction operator
324 e.g. $\frac{\partial u}{\partial y}$=0 along boundaries in $y$ and
325 $\frac{\partial v}{\partial x}$=0 along boundaries in $x$.
326
327 \item Lines 9,
328 \begin{verbatim}
329 no_slip_bottom=.TRUE.
330 \end{verbatim}
331 this line selects a no-slip boundary condition for bottom
332 boundary condition in the vertical Laplacian friction operator
333 e.g. $u=v=0$ at $z=-H$, where $H$ is the local depth of the domain.
334
335 \item Line 10,
336 \begin{verbatim}
337 diffKhT=4.E2,
338 \end{verbatim}
339 this line sets the horizontal diffusion coefficient for temperature
340 to $400\,{\rm m^{2}s^{-1}}$. The boundary condition on this
341 operator is $\frac{\partial}{\partial x}=\frac{\partial}{\partial y}=0$ at
342 all boundaries.
343
344 \item Line 11,
345 \begin{verbatim}
346 diffKzT=1.E-2,
347 \end{verbatim}
348 this line sets the vertical diffusion coefficient for temperature
349 to $10^{-2}\,{\rm m^{2}s^{-1}}$. The boundary condition on this
350 operator is $\frac{\partial}{\partial z}$ = 0 on all boundaries.
351
352 \item Line 13,
353 \begin{verbatim}
354 tAlpha=2.E-4,
355 \end{verbatim}
356 This line sets the thermal expansion coefficient for the fluid
357 to $2 \times 10^{-4}\,{\rm degrees}^{-1}$
358
359 \fbox{
360 \begin{minipage}{5.0in}
361 {\it S/R FIND\_RHO}({\it find\_rho.F})
362 \end{minipage}
363 }
364
365 \item Line 18,
366 \begin{verbatim}
367 eosType='LINEAR'
368 \end{verbatim}
369 This line selects the linear form of the equation of state.
370
371 \fbox{
372 \begin{minipage}{5.0in}
373 {\it S/R FIND\_RHO}({\it find\_rho.F})
374 \end{minipage}
375 }
376
377
378
379 \item Line 40,
380 \begin{verbatim}
381 usingSphericalPolarGrid=.TRUE.,
382 \end{verbatim}
383 This line requests that the simulation be performed in a
384 spherical polar coordinate system. It affects the interpretation of
385 grid input parameters, for example {\bf delX} and {\bf delY} and
386 causes the grid generation routines to initialize an internal grid based
387 on spherical polar geometry.
388
389 \fbox{
390 \begin{minipage}{5.0in}
391 {\it S/R INI\_SPEHRICAL\_POLAR\_GRID}({\it ini\_spherical\_polar\_grid.F})
392 \end{minipage}
393 }
394
395 \item Line 41,
396 \begin{verbatim}
397 phiMin=0.,
398 \end{verbatim}
399 This line sets the southern boundary of the modeled
400 domain to $0^{\circ}$ latitude. This value affects both the
401 generation of the locally orthogonal grid that the model
402 uses internally and affects the initialization of the coriolis force.
403 Note - it is not required to set
404 a longitude boundary, since the absolute longitude does
405 not alter the kernel equation discretisation.
406
407 \item Line 42,
408 \begin{verbatim}
409 delX=60*1.,
410 \end{verbatim}
411 This line sets the horizontal grid spacing between each y-coordinate line
412 in the discrete grid to $1^{\circ}$ in longitude.
413
414 \item Line 43,
415 \begin{verbatim}
416 delY=60*1.,
417 \end{verbatim}
418 This line sets the horizontal grid spacing between each y-coordinate line
419 in the discrete grid to $1^{\circ}$ in latitude.
420
421 \item Line 44,
422 \begin{verbatim}
423 delZ=500.,500.,500.,500.,
424 \end{verbatim}
425 This line sets the vertical grid spacing between each z-coordinate line
426 in the discrete grid to $500\,{\rm m}$, so that the total model depth
427 is $2\,{\rm km}$. The variable {\bf delZ} is copied into the internal
428 model coordinate variable {\bf delR}
429
430 \fbox{
431 \begin{minipage}{5.0in}
432 {\it S/R INI\_VERTICAL\_GRID}({\it ini\_vertical\_grid.F})
433 \end{minipage}
434 }
435
436 \item Line 47,
437 \begin{verbatim}
438 bathyFile='topog.box'
439 \end{verbatim}
440 This line specifies the name of the file from which the domain
441 bathymetry is read. This file is a two-dimensional ($x,y$) map of
442 depths. This file is assumed to contain 64-bit binary numbers
443 giving the depth of the model at each grid cell, ordered with the x
444 coordinate varying fastest. The points are ordered from low coordinate
445 to high coordinate for both axes. The units and orientation of the
446 depths in this file are the same as used in the MITgcm code. In this
447 experiment, a depth of $0m$ indicates a solid wall and a depth
448 of $-2000m$ indicates open ocean. The matlab program
449 {\it input/gendata.m} shows an example of how to generate a
450 bathymetry file.
451
452
453 \item Line 50,
454 \begin{verbatim}
455 zonalWindFile='windx.sin_y'
456 \end{verbatim}
457 This line specifies the name of the file from which the x-direction
458 surface wind stress is read. This file is also a two-dimensional
459 ($x,y$) map and is enumerated and formatted in the same manner as the
460 bathymetry file. The matlab program {\it input/gendata.m} includes example
461 code to generate a valid
462 {\bf zonalWindFile}
463 file.
464
465 \end{itemize}
466
467 \noindent other lines in the file {\it input/data} are standard values
468 that are described in the MITgcm Getting Started and MITgcm Parameters
469 notes.
470
471 \begin{small}
472 % \input{part3/case_studies/carbon_outgassing_sensitivity/input/data}
473 \end{small}
474
475 \subsubsection{File {\it input/data.pkg}}
476
477 This file uses standard default values and does not contain
478 customizations for this experiment.
479
480 \subsubsection{File {\it input/eedata}}
481
482 This file uses standard default values and does not contain
483 customizations for this experiment.
484
485 \subsubsection{File {\it input/windx.sin\_y}}
486
487 The {\it input/windx.sin\_y} file specifies a two-dimensional ($x,y$)
488 map of wind stress ,$\tau_{x}$, values. The units used are $Nm^{-2}$.
489 Although $\tau_{x}$ is only a function of $y$n in this experiment
490 this file must still define a complete two-dimensional map in order
491 to be compatible with the standard code for loading forcing fields
492 in MITgcm. The included matlab program {\it input/gendata.m} gives a complete
493 code for creating the {\it input/windx.sin\_y} file.
494
495 \subsubsection{File {\it input/topog.box}}
496
497
498 The {\it input/topog.box} file specifies a two-dimensional ($x,y$)
499 map of depth values. For this experiment values are either
500 $0m$ or $-2000\,{\rm m}$, corresponding respectively to a wall or to deep
501 ocean. The file contains a raw binary stream of data that is enumerated
502 in the same way as standard MITgcm two-dimensional, horizontal arrays.
503 The included matlab program {\it input/gendata.m} gives a complete
504 code for creating the {\it input/topog.box} file.
505
506 \subsubsection{File {\it code/SIZE.h}}
507
508 Two lines are customized in this file for the current experiment
509
510 \begin{itemize}
511
512 \item Line 39,
513 \begin{verbatim} sNx=60, \end{verbatim} this line sets
514 the lateral domain extent in grid points for the
515 axis aligned with the x-coordinate.
516
517 \item Line 40,
518 \begin{verbatim} sNy=60, \end{verbatim} this line sets
519 the lateral domain extent in grid points for the
520 axis aligned with the y-coordinate.
521
522 \item Line 49,
523 \begin{verbatim} Nr=4, \end{verbatim} this line sets
524 the vertical domain extent in grid points.
525
526 \end{itemize}
527
528 \begin{small}
529 % \include{code/SIZE.h}
530 \end{small}
531
532 \subsubsection{File {\it code/CPP\_OPTIONS.h}}
533
534 This file uses standard default values and does not contain
535 customizations for this experiment.
536
537
538 \subsubsection{File {\it code/CPP\_EEOPTIONS.h}}
539
540 This file uses standard default values and does not contain
541 customizations for this experiment.
542
543 \subsubsection{Other Files }
544
545 Other files relevant to this experiment are
546 \begin{itemize}
547 \item {\it model/src/ini\_cori.F}. This file initializes the model
548 coriolis variables {\bf fCorU}.
549 \item {\it model/src/ini\_spherical\_polar\_grid.F}
550 \item {\it model/src/ini\_parms.F},
551 \item {\it input/windx.sin\_y},
552 \end{itemize}
553 contain the code customizations and parameter settings for this
554 experiments. Below we describe the customizations
555 to these files associated with this experiment.

  ViewVC Help
Powered by ViewVC 1.1.22