/[MITgcm]/manual/s_algorithm/text/nonlin_frsurf.tex
ViewVC logotype

Annotation of /manual/s_algorithm/text/nonlin_frsurf.tex

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


Revision 1.7 - (hide annotations) (download) (as text)
Tue Nov 13 20:27:54 2001 UTC (23 years, 8 months ago) by adcroft
Branch: MAIN
Changes since 1.6: +2 -2 lines
File MIME type: application/x-tex
More x-refs

1 adcroft 1.7 % $Header: /u/gcmpack/mitgcmdoc/part2/nonlin_frsurf.tex,v 1.6 2001/11/13 20:13:54 adcroft Exp $
2 jmc 1.1 % $Name: $
3    
4 adcroft 1.4
5 jmc 1.1
6     \subsection{Non-linear free surface}
7 adcroft 1.6 \label{sect:nonlinear-freesurface}
8 jmc 1.1
9 adcroft 1.4 Recently, two options have been added to the model (and have not yet
10     been extensively tested) that concern the free surface formulation.
11    
12 jmc 1.1
13     \subsubsection{Non-uniform linear-relation for the surface potential}
14    
15 adcroft 1.4 The linear relation between surface pressure/geo-potential
16     ($\Phi_{surf}$) and surface displacement ($\eta$) could be considered
17     to be a constant ($b_s=$ constant)
18 jmc 1.3 \marginpar{add a reference to part.1 here}
19     but is in fact dependent on the position ($x,y,r$)
20 jmc 1.2 since we linearize:
21     $$\Phi_{surf}=\int_{R_o}^{R_o+\eta} b dr \simeq b_s \eta
22     ~\mathrm{with}~ b_s = b(\theta,S,r)_{r=R_o}
23     \simeq b_s(\theta_{ref}(R_o),S_{ref}(R_o),R_o)$$
24 adcroft 1.4 Note that, for convenience, the effect on $b_s$ of the local surface
25     $\theta,S$ are not considered here, but are incorporated in to
26     $\Phi'_{hyd}$.
27    
28     For the ocean, $b_s = g \rho_{surf} / \rho_o \simeq g$ is a very good
29     approximation since the relative difference in surface density are
30     usually small and only due to local $\theta,S$ gradients (because the
31     upper surface, $R_o = 0$, is essentially flat). Therefore, they can
32     easily be incorporated in $\Phi'_{hyd}$.
33    
34     For the atmosphere, however, because of topographic effects, the
35 cnh 1.5 reference surface pressure $R_o$ has large spatial variations that
36 adcroft 1.4 are responsible for significant $b_s$ variations (from 0.8 to 1.2
37     $[m^3/kg]$). For this reason, we use a non-uniform linear coefficient
38     $b_s$.
39    
40     In practice, in an oceanic configuration or when the default value
41     (TRUE) of the parameter {\bf uniformLin\_PhiSurf} is used, then $b_s$
42     is simply set to $g$ for the ocean and $1.$ for the atmosphere.
43     Turning {\bf uniformLin\_PhiSurf} to "FALSE", tells the code to
44     evaluate $b_s$ from the reference vertical profile $\theta_{ref}$
45     ({\it S/R INI\_LINEAR\_PHISURF}) according to the reference surface
46     pressure $P_o$ ($=R_o$): $b_s = c_p \kappa (P_o / Pc)^{(\kappa - 1)}
47     \theta_{ref}(P_o)$
48    
49 jmc 1.1
50     \subsubsection{Free surface effect on column total thickness
51     (Non-linear free surface)}
52    
53 adcroft 1.4 The total thickness of the fluid column is $r_{surf} - R_{fixed} =
54     \eta + R_o - R_{fixed}$ In the linear free surface approximation
55     (detailed before), only the fixed part of it ($R_o - R_{fixed})$ is
56     considered when we integrate the continuity equation or compute tracer
57     and momentum advection term.
58    
59     This approximation is dropped when using the non-linear free surface
60     formulation. Here we discuss sections the barotropic part. In
61     sections \ref{sect:freesurf-tracer-advection} and
62     \ref{sect:freesurf-momentum-advection} we consider the baroclinic
63     component.
64    
65    
66     The continuous form of the model equations remains unchanged, except
67 adcroft 1.7 for the 2D continuity equation (\ref{eq:discrete-time-backward-free-surface}) which is now
68 adcroft 1.4 integrated from $R_{fixed}(x,y)$ up to $r_{surf}=R_o+\eta$ :
69 jmc 1.1
70     \begin{displaymath}
71     \epsilon_{fs} \partial_t \eta =
72     \left. \dot{r} \right|_{r=r_{surf}} + \epsilon_{fw} (P-E) =
73 jmc 1.3 - {\bf \nabla}_h \cdot \int_{R_{fixed}}^{R_o+\eta} \vec{\bf v} dr
74 jmc 1.1 + \epsilon_{fw} (P-E)
75     \end{displaymath}
76    
77 adcroft 1.4 Since $\eta$ has a direct effect on the horizontal velocity (through
78     $\nabla_h \Phi_{surf}$), this adds a non-linear term to the free
79     surface equation. Several options for the time discretization of this
80     non-linear part have been tested.
81    
82     If the column thickness is evaluated at time step $n$, and with
83     implicit treatment of the surface potential gradient, equations
84     (\ref{eq-solve2D} and \ref{eq-solve2D_rhs}) becomes:
85 jmc 1.1 \begin{eqnarray*}
86     \epsilon_{fs} {\eta}^{n+1} -
87 jmc 1.3 {\bf \nabla}_h \cdot \Delta t^2 (\eta^{n}+R_o-R_{fixed})
88 jmc 1.2 {\bf \nabla}_h b_s {\eta}^{n+1}
89 jmc 1.1 = {\eta}^*
90     \end{eqnarray*}
91     where
92     \begin{eqnarray*}
93     {\eta}^* = \epsilon_{fs} \: {\eta}^{n} -
94 jmc 1.3 \Delta t {\bf \nabla}_h \cdot \int_{R_{fixed}}^{R_o+\eta^n} \vec{\bf v}^* dr
95 jmc 1.1 \: + \: \epsilon_{fw} \Delta_t (P-E)^{n}
96     \end{eqnarray*}
97 adcroft 1.4 This method requires us to update the solver matrix at each time step.
98 jmc 1.1
99     Alternatively, the non-linear contribution can be evaluated fully
100     explicitly:
101     \begin{eqnarray*}
102     \epsilon_{fs} {\eta}^{n+1} -
103 jmc 1.3 {\bf \nabla}_h \cdot \Delta t^2 (R_o-R_{fixed})
104 jmc 1.2 {\bf \nabla}_h b_s {\eta}^{n+1}
105 jmc 1.1 = {\eta}^*
106 jmc 1.2 +{\bf \nabla}_h \cdot \Delta t^2 (\eta^{n})
107     {\bf \nabla}_h b_s {\eta}^{n}
108 jmc 1.1 \end{eqnarray*}
109 adcroft 1.4 This formulation allows one to keep the initial solver matrix
110     unchanged though throughout the integration, since the non-linear free
111     surface only affects the RHS.
112    
113     Finally, another option is a "linearized" formulation where the total
114     column thickness appears only in the integral term of the RHS
115     (\ref{eq-solve2D_rhs}) but not directly in the equation
116     (\ref{eq-solve2D}).
117 jmc 1.1
118    
119     \subsubsection{Free surface effect on the surface level thickness
120     (Non-linear free surface): Tracer advection}
121 adcroft 1.4 \label{sect:freesurf-tracer-advection}
122 jmc 1.1
123 adcroft 1.4 To ensure global tracer conservation (i.e., the total amount) as well
124     as local conservation, the change in the surface level thickness must
125     be consistent with the way the continuity equation is integrated, both
126     in the barotropic part (to find $\eta$) and baroclinic part (to find
127     $w = \dot{r}$).
128 jmc 1.1
129 adcroft 1.4 To illustrate this, consider the shallow water model, with uniform
130     Cartesian horizontal grid:
131 jmc 1.1 $$
132     \partial_t h + \nabla \cdot h \vec{\bf v} = 0
133     $$
134     where $h$ is the total thickness of the water column.
135     To conserve the tracer $\theta$ we have to discretize:
136     $$
137     \partial_t (h \theta) + \nabla \cdot ( h \theta \vec{\bf v})= 0
138     $$
139 adcroft 1.4 Using the implicit (non-linear) free surface described above (section
140 adcroft 1.6 \ref{sect:pressure-method-linear-backward}) we have:
141 jmc 1.1 \begin{eqnarray*}
142 jmc 1.2 h^{n+1} = h^{n} - \Delta_t \nabla \cdot (h^n \, \vec{\bf v}^{n+1} ) \\
143 jmc 1.1 \end{eqnarray*}
144 adcroft 1.4 The discretized form of the tracer equation must adopt the same
145     ``form'' in the computation of tracer fluxes, that is, the same value
146     of $h$, as used in the continuity equation:
147 jmc 1.1 \begin{eqnarray*}
148     h^{n+1} \, \theta^{n+1} = h^n \, \theta^n
149 jmc 1.2 - \Delta_t \nabla \cdot (h^n \, \theta^n \, \vec{\bf v}^{n+1})
150 jmc 1.1 \end{eqnarray*}
151    
152 adcroft 1.4 For Adams-Bashforth time-stepping, we implement this scheme slightly
153     differently from the linear free-surface method, using two steps: the
154     variation of the water column thickness (from $h^n$ to $h^{n+1}$) is
155     not incorporated directly into the tracer equation. Instead, the
156     model uses the $G_\theta$ terms (first step) as in the linear free
157     surface formulation (with the "{\it surface correction}" turned "on",
158     see tracer section):
159 jmc 1.1 $$
160 jmc 1.2 G_\theta^n = \left(- \nabla \cdot (h^n \, \theta^n \, \vec{\bf v}^{n+1})
161     - \dot{r}_{surf}^{n+1} \theta^n \right) / h^n
162 jmc 1.1 $$
163 adcroft 1.4 Then, in a second step, the thickness variation (expansion/reduction)
164     is taken into account:
165 jmc 1.1 $$
166 jmc 1.2 \theta^{n+1} = \theta^n + \Delta_t \frac{h^n}{h^{n+1}} G_\theta^{(n+1/2)}
167 jmc 1.1 $$
168     Note that with a simple forward time step (no Adams-Bashforth),
169     since
170     $
171 jmc 1.2 \dot{r}_{surf}^{n+1}
172     = - \nabla \cdot (h^n \, \vec{\bf v}^{n+1} ) = (h^{n+1} - h^{n})
173 jmc 1.1 / \Delta_t
174     $
175 adcroft 1.4 these two formulations are equivalent.
176    
177     Implementation in the MITgcm is as follows. The model ``geometry''
178     (here the {\bf hFacC,W,S}) is updated just before entering {\it
179     SOLVE\_FOR\_PRESSURE}, using the current $\eta$ field. Then, at the
180     end of the time step, the variables are advanced in time, so that
181     $\eta^n$ replaces $\eta^{n-1}$. At the next time step, the tracer
182     tendency ($G_\theta$) is computed using the same geometry, which is
183     now consistent with $\eta^{n-1}$. Finally, in S/R {\it
184     TIMESTEP\_TRACER}, the expansion/reduction ratio is applied to the
185     surface level to compute the new tracer field.
186 jmc 1.1
187    
188     \subsubsection{Free surface effect on the surface level thickness
189     (Non-linear free surface): Momentum advection}
190 adcroft 1.4 \label{sect:freesurf-momentum-advection}
191 jmc 1.1
192     Regarding momentum advection,
193     the vector invariant formulation is similar to the
194     advective form (as opposed to the flux form) and therefore
195     does not need specific modification to include the
196     free surface effect on the surface level thickness.
197     Updating the {\bf hFacC,W,S} and the {\bf recip\_hFac}(s)
198     at one given place (like describe before) is sufficient.
199    
200     With the flux form formulation, advection of momentum
201 jmc 1.2 can be treated exactly as the tracer advection is.
202 jmc 1.1 Here the expansion/reduction factors ($hFacW^{n+1}/hFacW^n$ for $u$
203     and $hFacS^{n+1}/hFacS^n$ for $v$) are simply applied in the
204     subroutine {\it TIMESTEP}.
205    

  ViewVC Help
Powered by ViewVC 1.1.22