/[MITgcm]/manual/s_phys_pkgs/text/gridalt.tex
ViewVC logotype

Contents of /manual/s_phys_pkgs/text/gridalt.tex

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


Revision 1.8 - (show annotations) (download) (as text)
Wed Apr 5 02:27:33 2006 UTC (19 years, 3 months ago) by edhill
Branch: MAIN
Changes since 1.7: +2 -2 lines
File MIME type: application/x-tex
refer to the model uniformly as "MITgcm" and treat it as a proper noun

1 \subsection{Gridalt - Alternate Grid Package}
2 \label{sec:pkg:gridalt}
3 \begin{rawhtml}
4 <!-- CMIREDIR:package_gridalt: -->
5 \end{rawhtml}
6
7 \subsubsection {Introduction}
8
9 The gridalt package is designed to allow different components of MITgcm to
10 be run using horizontal and/or vertical grids which are different from the main
11 model grid. The gridalt routines handle the definition of the all the various
12 alternative grid(s) and the mappings between them and the MITgcm grid.
13 The implementation of the gridalt package which allows the high end atmospheric
14 physics (fizhi) to be run on a high resolution and quasi terrain-following vertical
15 grid is documented here. The package has also (with some user modifications) been used
16 for other calculations within the GCM.
17
18 The rationale for implementing the atmospheric physics on a high resolution vertical
19 grid involves the fact that the MITgcm $p^*$ (or any pressure-type) coordinate cannot
20 maintain the vertical resolution near the surface as the bottom topography rises above
21 sea level. The vertical length scales near the ground are small and can vary
22 on small time scales, and the vertical grid must be adequate to resolve them.
23 Many studies with both regional and global atmospheric models have demonstrated the
24 improvements in the simulations when the vertical resolution near the surface is
25 increased (\cite{bm:99,Inn:01,wo:98,breth:99}). Some of the benefit of increased resolution
26 near the surface is realized by employing the higher resolution for the computation of the
27 forcing due to turbulent and convective processes in the atmosphere.
28
29 The parameterizations of atmospheric subgrid scale processes are all essentially
30 one-dimensional in nature, and the computation of the terms in the equations of
31 motion due to these processes can be performed for the air column over one grid point
32 at a time. The vertical grid on which these computations take place can therefore be
33 entirely independant of the grid on which the equations of motion are integrated, and
34 the 'tendency' terms can be interpolated to the vertical grid on which the equations
35 of motion are integrated. A modified $p^*$ coordinate, which adjusts to the local
36 terrain and adds additional levels between the lower levels of the existing $p^*$ grid
37 (and perhaps between the levels near the tropopause as well), is implemented. The
38 vertical discretization is different for each grid point, although it consist of the
39 same number of levels. Additional 'sponge' levels aloft are added when needed. The levels
40 of the physics grid are constrained to fit exactly into the existing $p^*$ grid, simplifying
41 the mapping between the two vertical coordinates. This is illustrated as follows:
42
43 \begin{figure}[htbp]
44 \vspace*{-0.4in}
45 \begin{center}
46 \includegraphics[height=2.4in]{part6/vertical.eps}
47 \caption{Vertical discretization for MITgcm (dark grey lines) and for the
48 atmospheric physics (light grey lines). In this implementation, all MITgcm level
49 interfaces must coincide with atmospheric physics level interfaces.}
50 \end{center}
51 \end{figure}
52
53 The algorithm presented here retains the state variables on the high resolution 'physics'
54 grid as well as on the coarser resolution 'dynamics` grid, and ensures that the two
55 estimates of the state 'agree' on the coarse resolution grid. It would have been possible
56 to implement a technique in which the tendencies due to atmospheric physics are computed
57 on the high resolution grid and the state variables are retained at low resolution only.
58 This, however, for the case of the turbulence parameterization, would mean that the
59 turbulent kinetic energy source terms, and all the turbulence terms that are written
60 in terms of gradients of the mean flow, cannot really be computed making use of the fine
61 structure in the vertical.
62
63 \subsubsection{Equations on Both Grids}
64
65 In addition to computing the physical forcing terms of the momentum, thermodynamic and humidity
66 equations on the modified (higher resolution) grid, the higher resolution structure of the
67 atmosphere (the boundary layer) is retained between physics calculations. This neccessitates
68 a second set of evolution equations for the atmospheric state variables on the modified grid.
69 If the equation for the evolution of $U$ on $p^*$ can be expressed as:
70 \[
71 \left . {\partial U \over {\partial t}} \right |_{p^*}^{total} =
72 \left . {\partial U \over {\partial t}} \right |_{p^*}^{dynamics} +
73 \left . {\partial U \over {\partial t}} \right |_{p^*}^{physics}
74 \]
75 where the physics forcing terms on $p^*$ have been mapped from the modified grid, then an additional
76 equation to govern the evolution of $U$ (for example) on the modified grid is written:
77 \[
78 \left . {\partial U \over {\partial t}} \right |_{p^{*m}}^{total} =
79 \left . {\partial U \over {\partial t}} \right |_{p^{*m}}^{dynamics} +
80 \left . {\partial U \over {\partial t}} \right |_{p^{*m}}^{physics} +
81 \gamma ({\left . U \right |_{p^*}} - {\left . U \right |_{p^{*m}}})
82 \]
83 where $p^{*m}$ refers to the modified higher resolution grid, and the dynamics forcing terms have
84 been mapped from $p^*$ space. The last term on the RHS is a relaxation term, meant to constrain
85 the state variables on the modified vertical grid to `track' the state variables on the $p^*$ grid
86 on some time scale, governed by $\gamma$. In the present implementation, $\gamma = 1$, requiring
87 an immediate agreement between the two 'states'.
88
89 \subsubsection{Time stepping Sequence}
90 If we write $T_{phys}$ as the temperature (or any other state variable) on the high
91 resolution physics grid, and $T_{dyn}$ as the temperature on the coarse vertical resolution
92 dynamics grid, then:
93
94 \begin{enumerate}
95 %\itemsep{-0.05in}
96
97 \item{Compute the tendency due to physics processes.}
98
99 \item{Advance the physics state: ${{T^{n+1}}^{**}}_{phys}(l) = {T^n}_{phys}(l) + \delta T_{phys}$.}
100
101 \item{Interpolate the physics tendency to the dynamics grid, and advance the dynamics
102 state by physics and dynamics tendencies:
103 ${T^{n+1}}_{dyn}(L) = {T^n}_{dyn}(L) + \delta T_{dyn}(L) + [\delta T _{phys}(l)](L)$.}
104
105 \item{Interpolate the dynamics tendency to the physics grid, and update the physics
106 grid due to dynamics tendencies:
107 ${{T^{n+1}}^*}_{phys}(l)$ = ${{T^{n+1}}^{**}}_{phys}(l) + {\delta T_{dyn}(L)}(l)$.}
108
109 \item{Apply correction term to physics state to account for divergence from dynamics state:
110 ${T^{n+1}}_{phys}(l)$ = ${{T^{n+1}}^*}_{phys}(l) + \gamma \{ T_{dyn}(L) - [T_{phys}(l)](L) \}(l)$.} \\
111 Where $\gamma=1$ here.
112
113 \end{enumerate}
114
115 \subsubsection{Interpolation}
116 In order to minimize the correction terms for the state variables on the alternative,
117 higher resolution grid, the vertical interpolation scheme must be constructed so that
118 a dynamics-to-physics interpolation can be exactly reversed with a physics-to-dynamics mapping.
119 The simple scheme employed to achieve this is:\\
120
121 Coarse to fine:\
122 For all physics layers l in dynamics layer L, $ T_{phys}(l) = \{T_{dyn}(L)\} = T_{dyn}(L) $.
123
124 Fine to coarse:\
125 For all physics layers l in dynamics layer L, $T_{dyn}(L) = [T_{phys}(l)] = \int{T_{phys} dp } $.\\
126
127 Where $\{\}$ is defined as the dynamics-to-physics operator and $[ ]$ is the physics-to-dynamics operator, $T$ stands for any state variable, and the subscripts $phys$ and $dyn$ stand for variables on
128 the physics and dynamics grids, respectively.
129
130 \subsubsection {Key subroutines, parameters and files }
131
132 \noindent
133 One of the central elements of the gridalt package is the routine which
134 is called from subroutine gridalt\_initialise to define the grid to be
135 used for the high end physics calculations. Routine make\_phys\_grid
136 passes back the parameters which define the grid, ultimately stored
137 in the common block gridalt\_mapping.
138
139 \begin{verbatim}
140 subroutine make_phys_grid(drF,hfacC,im1,im2,jm1,jm2,Nr,
141 . Nsx,Nsy,i1,i2,j1,j2,bi,bj,Nrphys,Lbot,dpphys,numlevphys,nlperdyn)
142 c***********************************************************************
143 c Purpose: Define the grid that the will be used to run the high-end
144 c atmospheric physics.
145 c
146 c Algorithm: Fit additional levels of some (~) known thickness in
147 c between existing levels of the grid used for the dynamics
148 c
149 c Need: Information about the dynamics grid vertical spacing
150 c
151 c Input: drF - delta r (p*) edge-to-edge
152 c hfacC - fraction of grid box above topography
153 c im1, im2 - beginning and ending i - dimensions
154 c jm1, jm2 - beginning and ending j - dimensions
155 c Nr - number of levels in dynamics grid
156 c Nsx,Nsy - number of processes in x and y direction
157 c i1, i2 - beginning and ending i - index to fill
158 c j1, j2 - beginning and ending j - index to fill
159 c bi, bj - x-dir and y-dir index of process
160 c Nrphys - number of levels in physics grid
161 c
162 c Output: dpphys - delta r (p*) edge-to-edge of physics grid
163 c numlevphys - number of levels used in the physics
164 c nlperdyn - physics level number atop each dynamics layer
165 c
166 c NOTES: 1) Pressure levs are built up from bottom, using p0, ps and dp:
167 c p(i,j,k)=p(i,j,k-1) + dp(k)*ps(i,j)/p0(i,j)
168 c 2) Output dp's are aligned to fit EXACTLY between existing
169 c levels of the dynamics vertical grid
170 c 3) IMPORTANT! This routine assumes the levels are numbered
171 c from the bottom up, ie, level 1 is the surface.
172 c IT WILL NOT WORK OTHERWISE!!!
173 c 4) This routine does NOT work for surface pressures less
174 c (ie, above in the atmosphere) than about 350 mb
175 c***********************************************************************
176 \end{verbatim}
177
178 \noindent In the case of the grid used to compute the atmospheric physical
179 forcing (fizhi package), the locations of the grid points move in time with
180 the MITgcm $p^*$ coordinate, and subroutine gridalt\_update is called during
181 the run to update the locations of the grid points:
182
183 \begin{verbatim}
184 subroutine gridalt_update(myThid)
185 c***********************************************************************
186 c Purpose: Update the pressure thicknesses of the layers of the
187 c alternative vertical grid (used now for atmospheric physics).
188 c
189 c Calculate: dpphys - new delta r (p*) edge-to-edge of physics grid
190 c using dpphys0 (initial value) and rstarfacC
191 c***********************************************************************
192 \end{verbatim}
193
194 \noindent The gridalt package also supplies utility routines which perform
195 the mappings from one grid to the other. These routines are called from the
196 code which computes the fields on the alternative (fizhi) grid.
197
198 \begin{verbatim}
199 subroutine dyn2phys(qdyn,pedyn,im1,im2,jm1,jm2,lmdyn,Nsx,Nsy,
200 . idim1,idim2,jdim1,jdim2,bi,bj,windphy,pephy,Lbot,lmphy,nlperdyn,
201 . flg,qphy)
202 C***********************************************************************
203 C Purpose:
204 C To interpolate an arbitrary quantity from the 'dynamics' eta (pstar)
205 C grid to the higher resolution physics grid
206 C Algorithm:
207 C Routine works one layer (edge to edge pressure) at a time.
208 C Dynamics -> Physics retains the dynamics layer mean value,
209 C weights the field either with the profile of the physics grid
210 C wind speed (for U and V fields), or uniformly (T and Q)
211 C
212 C Input:
213 C qdyn..... [im,jm,lmdyn] Arbitrary Quantity on Input Grid
214 C pedyn.... [im,jm,lmdyn+1] Pressures at bottom edges of input levels
215 C im1,2 ... Limits for Longitude Dimension of Input
216 C jm1,2 ... Limits for Latitude Dimension of Input
217 C lmdyn.... Vertical Dimension of Input
218 C Nsx...... Number of processes in x-direction
219 C Nsy...... Number of processes in y-direction
220 C idim1,2.. Beginning and ending i-values to calculate
221 C jdim1,2.. Beginning and ending j-values to calculate
222 C bi....... Index of process number in x-direction
223 C bj....... Index of process number in x-direction
224 C windphy.. [im,jm,lmphy] Magnitude of the wind on the output levels
225 C pephy.... [im,jm,lmphy+1] Pressures at bottom edges of output levels
226 C lmphy.... Vertical Dimension of Output
227 C nlperdyn. [im,jm,lmdyn] Highest Physics level in each dynamics level
228 C flg...... Flag to indicate field type (0 for T or Q, 1 for U or V)
229 C
230 C Output:
231 C qphy..... [im,jm,lmphy] Quantity at output grid (physics grid)
232 C
233 C Notes:
234 C 1) This algorithm assumes that the output (physics) grid levels
235 C fit exactly into the input (dynamics) grid levels
236 C***********************************************************************
237 \end{verbatim}
238
239 \noindent And similarly, gridalt contains subroutine phys2dyn.
240
241 \subsubsection {Dos and donts}
242
243 \subsubsection {Gridalt Reference}

  ViewVC Help
Powered by ViewVC 1.1.22