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

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

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

-revision 1.1 by jmc,
Thu Aug  9 16:22:09 2001 UTC
+revision 1.2 by jmc,
Fri Aug 17 18:38:10 2001 UTC
 Line 15 
 that concern the free surface formulatio
  The linear relation between
  surface pressure / geo- potential ($\Phi_{surf}$)
  and surface displacement ($\eta$)
- has been considered as uniform ($b_{surf} = Constant$)
+ has been considered as uniform ($b_s =$ Constant)
  but is in fact
- dependent on the position
+ dependent on the position ($x,y,r$)
- since we have
+ since we linearize:
- $\Phi_{surf}=\int_{R_o}^{R_o+\eta} b dr \simeq b_{surf} \eta$
+ $$\Phi_{surf}=\int_{R_o}^{R_o+\eta} b dr \simeq b_s \eta
+ ~\mathrm{with}~ b_s = b(\theta,S,r)_{r=R_o}
+ \simeq b_s(\theta_{ref}(R_o),S_{ref}(R_o),R_o)$$
+ Note that, for convinience, the effect of the local
+ surface $\theta,S$ on $b_s$
+ are not considered here, but incorporated in $\Phi'_{hyd}$.
- For the ocean, $b_{surf} = g \rho_{surf} / \rho_o \simeq g$
+ For the ocean, $b_s = g \rho_{surf} / \rho_o \simeq g$
  is a fairly good approximation since the relative difference
- in surface density are usually small.
+ in surface density are usually small and only due to
+ local $\theta,S$ gradient (because $R_o = 0$);
+ Therefore, they can easely be incorporated in $\Phi'_{hyd}$.
  For the atmosphere, because of topographic effects,
  the reference surface pressure $R_o$ has large spacial differences
- that are responsible for significant $b_surf$ variations
+ that are responsible for significant $b_s$ variations
- (from 0.8 to 1.2 $[]$).
+ (from 0.8 to 1.2 $[m^3/kg]$). For this reason,
+ we use a non-uniform linear coefficient $b_s$.
  Practically, in an oceanic configuration or
- when the default value (TRUE) of the parameter {\bf uniformLin\_PhiSurf}
+ when the default value (TRUE) of the parameter
- is used )
+ {\bf uniformLin\_PhiSurf} is used
- then $b_{surf}$ is simply set to $g$ for the ocean
+ then $b_s$ is simply set to $g$ for the ocean
  and $1.$ for the atmosphere.\\
  Turning {\bf uniformLin\_PhiSurf} to "FALSE", allows to
- evaluate $b_{surf}$ from the reference Theta,Q vertical profile
+ evaluate $b_s$ from the reference vertical profile $\theta_{ref}$
  ({\it S/R INI\_LINEAR\_PHISURF})
  according to the reference surface pressure $P_o$ ($=R_o$):
- $b_{surf} = c_p \kappa (P_o / Pc)^{(\kappa - 1)} \theta_{ref}(P_o)$
+ $b_s = c_p \kappa (P_o / Pc)^{(\kappa - 1)} \theta_{ref}(P_o)$
  %------------------------------------------
  \subsubsection{Free surface effect on column total thickness
-Line 52 
 it ($R_o - R_{min})$ is considered when
+Line 61 
 it ($R_o - R_{min})$ is considered when
  continuity equation or compute tracer and momentum advection term.
  This approximation is dropped when using
- the non linear free surface formulation.
+ the non-linear free surface formulation.
  Details follow here after for the barotropic part
- and in the 2 next subsection for the baroclinic
+ and in the 2 next subsections for the baroclinic
  part.
  %------------------------------------------
- % Non Linear Barotropic part
+ % Non-Linear Barotropic part
  The continuous form of the model equations remains
  unchanged, except for the 2D continuity equation
- (\ref{eq-cont-2D}) that is now integrated
+ (\ref{eq-tCsC-eta}) that is now integrated
  from $R_{min}(x,y)$ up to $r_{surf}=R_o+\eta$ :
  \begin{displaymath}
  \epsilon_{fs} \partial_t \eta =
  \left. \dot{r} \right|_{r=r_{surf}} + \epsilon_{fw} (P-E) =
- - {\bf \nabla}_r \cdot \int_{R_{min}}^{R_o+\eta} \vec{\bf v} dr
+ - {\bf \nabla}_h \cdot \int_{R_{min}}^{R_o+\eta} \vec{\bf v} dr
  + \epsilon_{fw} (P-E)
  \end{displaymath}
  Since $\eta$ has a direct effect on the horizontal
- velocity (through $\nabla_r \Phi_{surf}$), this
+ velocity (through $\nabla_h \Phi_{surf}$), this
- add a non linear term to the free surface equation.
+ adds a non-linear term to the free surface equation.
- Regarding the time discretization of this non linear part,
+ Regarding the time discretization of this non-linear part,
- several option are now being tested:
+ several options are now being tested:
  With the column thickness evaluated at time step $n$,
  and the surface potential gradient still implicit,
- equation (\ref{solve_2d} \& \ref{solve_2d_rhs})
+ equation (\ref{eq-solve2D} \& \ref{eq-solve2D_rhs})
  become:
  \begin{eqnarray*}
  \epsilon_{fs} {\eta}^{n+1} -
- {\bf \nabla}_r \cdot \Delta t^2 (\eta^{n}+R_o-R_{min})
+ {\bf \nabla}_h \cdot \Delta t^2 (\eta^{n}+R_o-R_{min})
- {\bf \nabla}_r B_o {\eta}^{n+1}
+ {\bf \nabla}_h b_s {\eta}^{n+1}
  = {\eta}^*
  %\label{solve_2d}
  \end{eqnarray*}
  where
  \begin{eqnarray*}
  {\eta}^* = \epsilon_{fs} \: {\eta}^{n} -
- \Delta t {\bf \nabla}_r \cdot \int_{R_{min}}^{R_o+\eta} \vec{\bf v}^* dr
+ \Delta t {\bf \nabla}_h \cdot \int_{R_{min}}^{R_o+\eta^n} \vec{\bf v}^* dr
  \: + \: \epsilon_{fw} \Delta_t (P-E)^{n}
  %\label{solve_2d_rhs}
  \end{eqnarray*}
-Line 103 
 Alternatively, the non-linear contributi
+Line 112 
 Alternatively, the non-linear contributi
  explicitly:
  \begin{eqnarray*}
  \epsilon_{fs} {\eta}^{n+1} -
- {\bf \nabla}_r \cdot \Delta t^2 (R_o-R_{min})
+ {\bf \nabla}_h \cdot \Delta t^2 (R_o-R_{min})
- {\bf \nabla}_r B_o {\eta}^{n+1}
+ {\bf \nabla}_h b_s {\eta}^{n+1}
  = {\eta}^*
- +{\bf \nabla}_r \cdot \Delta t^2 (\eta^{n})
+ +{\bf \nabla}_h \cdot \Delta t^2 (\eta^{n})
- {\bf \nabla}_r B_o {\eta}^{n}
+ {\bf \nabla}_h b_s {\eta}^{n}
  \end{eqnarray*}
  This formulation allows to keep the initial solver matrix
  since the non-linear free surface only affects the RHS.
  An other option is a "linearized" formulation where the
  total column thickness appears only in the integral term of
- the RHS (\ref{solve_2d_rhs}) but not directly in
+ the RHS (\ref{eq-solve2D_rhs}) but not directly in
- the equation (\ref{solve_2d}).
+ the equation (\ref{eq-solve2D}).
  %------------------------------------------
  \subsubsection{Free surface effect on the surface level thickness
-Line 126 
 as well as the local one (see tracer sec
+Line 135 
 as well as the local one (see tracer sec
  the change in the surface level thickness must be consistent with
  the way the continuity equation is integrated, both in
  in the barotropic part (to find $\eta$) and baroclinic part
- (to find $\dot{r}$).
+ (to find $w = \dot{r}$).
  To illustrate this, let's consider the shallow water model,
  with uniform cartesian horizontal grid:
-Line 140 
 $$
+Line 149 
 $$
  $$
  Using the implicit (non-linear) free surface described before, we have:
  \begin{eqnarray*}
- h^{n+1} = h^{n} - \Delta_t \nabla (h^n \, \vec{\bf v}^{n+1} ) \\
+ h^{n+1} = h^{n} - \Delta_t \nabla \cdot (h^n \, \vec{\bf v}^{n+1} ) \\
  \end{eqnarray*}
  The discretized form of the tracer equation must use the same
  "geometry" to compute the tracer fluxes, that is, the same value of
  h, as the continuity equation did:
  \begin{eqnarray*}
  h^{n+1} \, \theta^{n+1} = h^n \, \theta^n
-         - \Delta_t \nabla (h^n \, \theta^n \, \vec{\bf v}^{n+1})
+         - \Delta_t \nabla \cdot (h^n \, \theta^n \, \vec{\bf v}^{n+1})
  \end{eqnarray*}
  In order to deal with the Adams-Bashforth time stepping,
- we implement this scheme slightly differently:
+ we implement this scheme slightly differently, in two step:
  the variation of the water column thickness (from
  $h^n$ to $h^{n+1}$)
  is not incorporated directly to the tracer equation.
  Instead,
- the model continues to evaluate the $G_\theta$ term
+ the model continues to evaluate the $G_\theta$ term (first step)
- exactly as it did with the Linear free surface formulation,
+ as it use to do with the Linear free surface formulation
- with the "{\it surface correction}" turned "on" (see tracer section):
+ (with the "{\it surface correction}" turned "on", see tracer section):
  $$
- G_\theta = \left(- \nabla (h^n \, \theta^n \, \vec{\bf v}^{n+1})
+ G_\theta^n = \left(- \nabla \cdot (h^n \, \theta^n \, \vec{\bf v}^{n+1})
-          + w_{surf}^{n+1} \theta^n \right) / h^n
+          - \dot{r}_{surf}^{n+1} \theta^n \right) / h^n
  $$
- Then after,
+ Then in a second step,
  thickness variation (expansion/reduction) is taken into account :
  $$
- \theta^{n+1} = (\theta^n + \Delta_t G_\theta^{(n+1)}) \frac{h^{n+1}}{h^n}
+ \theta^{n+1} = \theta^n + \Delta_t \frac{h^n}{h^{n+1}} G_\theta^{(n+1/2)}
  $$
  Note that with a simple forward time step (no Adams-Bashforth),
  since
  $
- w_{surf} = - \nabla (h^n \, \vec{\bf v}^{n+1} ) = (h^{n+1} - h^{n})
+ \dot{r}_{surf}^{n+1}
+ = - \nabla \cdot (h^n \, \vec{\bf v}^{n+1} ) = (h^{n+1} - h^{n})
  / \Delta_t
  $
  those 2 formulations are equivalent.
  The implementation in the MITgcm follows this scheme.
  The model "geometry" (here the {\bf hFacC,W,S}) is updated
- at the beginning of {\it SOLVE\_FOR\_PRESSURE}
+ just before entering {\it SOLVE\_FOR\_PRESSURE},
  according to the current $\eta$ field.
  Then, at the end of the time step, the variables are
  advanced in time, so that $\eta^n$ becomes $\eta^{n-1}$.
-Line 201 
 Updating the {\bf hFacC,W,S} and the {\b
+Line 211 
 Updating the {\bf hFacC,W,S} and the {\b
  at one given place (like describe before) is sufficient.
  With the flux form formulation, advection of momentum
- can be treated exactly as the tracer advection is,
+ can be treated exactly as the tracer advection is.
  Here the expansion/reduction factors ($hFacW^{n+1}/hFacW^n$ for $u$
  and $hFacS^{n+1}/hFacS^n$ for $v$) are simply applied in the
  subroutine {\it TIMESTEP}.

 Legend:



Removed from v.1.1
 


changed lines


 
Added in v.1.2
 Legend:



Removed from v.1.1
 


changed lines


 
Added in v.1.2
-Removed from v.1.1
+Added in v.1.2

	ViewVC Help
Powered by ViewVC 1.1.22