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

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

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

revision 1.14 by adcroft, Wed Nov 14 21:07:13 2001 UTC revision 1.18 by jmc, Thu Oct 14 22:22:30 2004 UTC
# Line 1  Line 1 
1  % $Header$  % $Header$
2  % $Name$  % $Name$
3    
4    This chapter lays out the numerical schemes that are
5    employed in the core MITgcm algorithm. Whenever possible
6    links are made to actual program code in the MITgcm implementation.
7    The chapter begins with a discussion of the temporal discretization
8    used in MITgcm. This discussion is followed by sections that
9    describe the spatial discretization. The schemes employed for momentum
10    terms are described first, afterwards the schemes that apply to
11    passive and dynamically active tracers are described.
12    
13    
14    \section{Time-stepping}
15  The equations of motion integrated by the model involve four  The equations of motion integrated by the model involve four
16  prognostic equations for flow, $u$ and $v$, temperature, $\theta$, and  prognostic equations for flow, $u$ and $v$, temperature, $\theta$, and
17  salt/moisture, $S$, and three diagnostic equations for vertical flow,  salt/moisture, $S$, and three diagnostic equations for vertical flow,
# Line 70  temporary.} Line 81  temporary.}
81  \begin{figure}  \begin{figure}
82  \begin{center} \fbox{ \begin{minipage}{4.5in} \begin{tabbing}  \begin{center} \fbox{ \begin{minipage}{4.5in} \begin{tabbing}
83  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill
84  FORWARD\_STEP \\  \filelink{FORWARD\_STEP}{model-src-forward_step.F} \\
85  \> DYNAMICS \\  \> DYNAMICS \\
86  \>\> TIMESTEP \` $u^*$,$v^*$ (\ref{eq:ustar-rigid-lid},\ref{eq:vstar-rigid-lid}) \\  \>\> TIMESTEP \` $u^*$,$v^*$ (\ref{eq:ustar-rigid-lid},\ref{eq:vstar-rigid-lid}) \\
87  \> SOLVE\_FOR\_PRESSURE \\  \> SOLVE\_FOR\_PRESSURE \\
88  \>\> CALC\_DIV\_GHAT \` $H\widehat{u^*}$,$H\widehat{v^*}$ (\ref{eq:elliptic}) \\  \>\> CALC\_DIV\_GHAT \` $H\widehat{u^*}$,$H\widehat{v^*}$ (\ref{eq:elliptic}) \\
89  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic}) \\  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic}) \\
90  \> THE\_CORRECTION\_STEP  \\  \> MOMENTUM\_CORRECTION\_STEP  \\
91  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\
92  \>\> CORRECTION\_STEP \` $u^{n+1}$,$v^{n+1}$ (\ref{eq:un+1-rigid-lid},\ref{eq:vn+1-rigid-lid})  \>\> CORRECTION\_STEP \` $u^{n+1}$,$v^{n+1}$ (\ref{eq:un+1-rigid-lid},\ref{eq:vn+1-rigid-lid})
93  \end{tabbing} \end{minipage} } \end{center}  \end{tabbing} \end{minipage} } \end{center}
94  \caption{Calling tree for the pressure method algorihtm}  \caption{Calling tree for the pressure method algorithm
95      (\filelink{FORWARD\_STEP}{model-src-forward_step.F})}
96  \label{fig:call-tree-pressure-method}  \label{fig:call-tree-pressure-method}
97  \end{figure}  \end{figure}
98    
# Line 162  The correspondence to the code is as fol Line 174  The correspondence to the code is as fol
174  \item  \item
175  the prognostic phase, equations \ref{eq:ustar-rigid-lid} and \ref{eq:vstar-rigid-lid},  the prognostic phase, equations \ref{eq:ustar-rigid-lid} and \ref{eq:vstar-rigid-lid},
176  stepping forward $u^n$ and $v^n$ to $u^{*}$ and $v^{*}$ is coded in  stepping forward $u^n$ and $v^n$ to $u^{*}$ and $v^{*}$ is coded in
177  {\em TIMESTEP.F}  \filelink{TIMESTEP()}{model-src-timestep.F}
178  \item  \item
179  the vertical integration, $H \widehat{u^*}$ and $H  the vertical integration, $H \widehat{u^*}$ and $H
180  \widehat{v^*}$, divergence and inversion of the elliptic operator in  \widehat{v^*}$, divergence and inversion of the elliptic operator in
181  equation \ref{eq:elliptic} is coded in {\em  equation \ref{eq:elliptic} is coded in
182  SOLVE\_FOR\_PRESSURE.F}  \filelink{SOLVE\_FOR\_PRESSURE()}{model-src-solve_for_pressure.F}
183  \item  \item
184  finally, the new flow field at time level $n+1$ given by equations  finally, the new flow field at time level $n+1$ given by equations
185  \ref{eq:un+1-rigid-lid} and \ref{eq:vn+1-rigid-lid} is calculated in {\em CORRECTION\_STEP.F}.  \ref{eq:un+1-rigid-lid} and \ref{eq:vn+1-rigid-lid} is calculated in
186    \filelink{CORRECTION\_STEP()}{model-src-correction_step.F}.
187  \end{itemize}  \end{itemize}
188  The calling tree for these routines is given in  The calling tree for these routines is given in
189  Fig.~\ref{fig:call-tree-pressure-method}.  Fig.~\ref{fig:call-tree-pressure-method}.
# Line 281  FORWARD\_STEP \\ Line 294  FORWARD\_STEP \\
294  \> THERMODYNAMICS \\  \> THERMODYNAMICS \\
295  \>\> CALC\_GT \\  \>\> CALC\_GT \\
296  \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$ \\  \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$ \\
297  \>\>\> EXTERNAL\_FORCING \` $G_\theta^n = G_\theta^n + {\cal Q}$ \\  \>either\>\> EXTERNAL\_FORCING \` $G_\theta^n = G_\theta^n + {\cal Q}$ \\
298  \>\>\> ADAMS\_BASHFORTH2 \` $G_\theta^{(n+1/2)}$ (\ref{eq:adams-bashforth2}) \\  \>\>\> ADAMS\_BASHFORTH2 \` $G_\theta^{(n+1/2)}$ (\ref{eq:adams-bashforth2}) \\
299    \>or\>\> EXTERNAL\_FORCING \` $G_\theta^{(n+1/2)} = G_\theta^{(n+1/2)} + {\cal Q}$ \\
300  \>\> TIMESTEP\_TRACER \` $\tau^*$ (\ref{eq:taustar}) \\  \>\> TIMESTEP\_TRACER \` $\tau^*$ (\ref{eq:taustar}) \\
301  \>\> IMPLDIFF \` $\tau^{(n+1)}$ (\ref{eq:tau-n+1-implicit})  \>\> IMPLDIFF \` $\tau^{(n+1)}$ (\ref{eq:tau-n+1-implicit})
302  \end{tabbing} \end{minipage} } \end{center}  \end{tabbing} \end{minipage} } \end{center}
# Line 319  simpler to include these terms and this Line 333  simpler to include these terms and this
333  and forcing evolves smoothly. Problems can, and do, arise when forcing  and forcing evolves smoothly. Problems can, and do, arise when forcing
334  or motions are high frequency and this corresponds to a reduced  or motions are high frequency and this corresponds to a reduced
335  stability compared to a simple forward time-stepping of such terms.  stability compared to a simple forward time-stepping of such terms.
336    The model offers the possibility to leave the forcing term outside the
337    Adams-Bashforth extrapolation, by turning off the logical flag
338    {\bf forcing\_In\_AB } (parameter file {\em data}, namelist {\em PARM01},
339    default value = True).
340    
341  A stability analysis for an oscillation equation should be given at this point.  A stability analysis for an oscillation equation should be given at this point.
342  \marginpar{AJA needs to find his notes on this...}  \marginpar{AJA needs to find his notes on this...}
# Line 330  A stability analysis for a relaxation eq Line 348  A stability analysis for a relaxation eq
348  \section{Implicit time-stepping: backward method}  \section{Implicit time-stepping: backward method}
349    
350  Vertical diffusion and viscosity can be treated implicitly in time  Vertical diffusion and viscosity can be treated implicitly in time
351  using the backward method which is an intrinsic scheme. For tracers,  using the backward method which is an intrinsic scheme.
352    Recently, the option to treat the vertical advection
353    implicitly has been added, but not yet tested; therefore,
354    the description hereafter is limited to diffusion and viscosity.
355    For tracers,
356  the time discretized equation is:  the time discretized equation is:
357  \begin{equation}  \begin{equation}
358  \tau^{n+1} - \Delta t \partial_r \kappa_v \partial_r \tau^{n+1} =  \tau^{n+1} - \Delta t \partial_r \kappa_v \partial_r \tau^{n+1} =
# Line 361  Fig.~\ref{fig:call-tree-adams-bashforth} Line 383  Fig.~\ref{fig:call-tree-adams-bashforth}
383  stepping forward a tracer variable such as temperature.  stepping forward a tracer variable such as temperature.
384    
385  In order to fit within the pressure method, the implicit viscosity  In order to fit within the pressure method, the implicit viscosity
386  must not alter the barotropic flow. In other words, it can on ly  must not alter the barotropic flow. In other words, it can only
387  redistribute momentum in the vertical. The upshot of this is that  redistribute momentum in the vertical. The upshot of this is that
388  although vertical viscosity may be backward implicit and  although vertical viscosity may be backward implicit and
389  unconditionally stable, no-slip boundary conditions may not be made  unconditionally stable, no-slip boundary conditions may not be made
# Line 389  is solved to yield the state variables a Line 411  is solved to yield the state variables a
411  \end{figure}  \end{figure}
412    
413  \begin{figure}  \begin{figure}
414  \begin{center} \fbox{ \begin{minipage}{4.5in} \begin{tabbing}  \begin{center} \fbox{ \begin{minipage}{4.7in} \begin{tabbing}
415  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill
416  FORWARD\_STEP \\  FORWARD\_STEP \\
417    \>\> EXTERNAL\_FIELDS\_LOAD\\
418    \>\> DO\_ATMOSPHERIC\_PHYS \\
419    \>\> DO\_OCEANIC\_PHYS \\
420  \> THERMODYNAMICS \\  \> THERMODYNAMICS \\
421  \>\> CALC\_GT \\  \>\> CALC\_GT \\
422  \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$ (\ref{eq:Gt-n-sync})\\  \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$ (\ref{eq:Gt-n-sync})\\
# Line 404  FORWARD\_STEP \\ Line 429  FORWARD\_STEP \\
429  \>\> MOM\_FLUXFORM or MOM\_VECINV \` $G_{\vec{\bf v}}^n$ (\ref{eq:Gv-n-sync})\\  \>\> MOM\_FLUXFORM or MOM\_VECINV \` $G_{\vec{\bf v}}^n$ (\ref{eq:Gv-n-sync})\\
430  \>\> TIMESTEP \` $\vec{\bf v}^*$ (\ref{eq:Gv-n+5-sync}, \ref{eq:vstar-sync}) \\  \>\> TIMESTEP \` $\vec{\bf v}^*$ (\ref{eq:Gv-n+5-sync}, \ref{eq:vstar-sync}) \\
431  \>\> IMPLDIFF \` $\vec{\bf v}^{**}$ (\ref{eq:vstarstar-sync}) \\  \>\> IMPLDIFF \` $\vec{\bf v}^{**}$ (\ref{eq:vstarstar-sync}) \\
432    \> UPDATE\_R\_STAR or UPDATE\_SURF\_DR \` (NonLin-FS only)\\
433  \> SOLVE\_FOR\_PRESSURE \\  \> SOLVE\_FOR\_PRESSURE \\
434  \>\> CALC\_DIV\_GHAT \` $\eta^*$ (\ref{eq:nstar-sync}) \\  \>\> CALC\_DIV\_GHAT \` $\eta^*$ (\ref{eq:nstar-sync}) \\
435  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic-sync}) \\  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic-sync}) \\
436  \> THE\_CORRECTION\_STEP  \\  \> MOMENTUM\_CORRECTION\_STEP  \\
437  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\
438  \>\> CORRECTION\_STEP \` $u^{n+1}$,$v^{n+1}$ (\ref{eq:v-n+1-sync})  \>\> CORRECTION\_STEP \` $u^{n+1}$,$v^{n+1}$ (\ref{eq:v-n+1-sync})\\
439    \> TRACERS\_CORRECTION\_STEP  \\
440    \>\> CYCLE\_TRACER \` $\theta^{n+1}$ \\
441    \>\> FILTER \` Shapiro Filter, Zonal Filter (FFT) \\
442    \>\> CONVECTIVE\_ADJUSTMENT \` \\
443  \end{tabbing} \end{minipage} } \end{center}  \end{tabbing} \end{minipage} } \end{center}
444  \caption{  \caption{
445  Calling tree for the overall synchronous algorithm using  Calling tree for the overall synchronous algorithm using
446  Adams-Bashforth time-stepping.}  Adams-Bashforth time-stepping.
447    The place where the model geometry
448    ({\em hFac} factors) is updated is added here but is only relevant
449    for the non-linear free-surface algorithm.
450    For completeness, the external forcing,
451    ocean and atmospheric physics have been added, although they are mainly
452    optional}
453  \label{fig:call-tree-adams-bashforth-sync}  \label{fig:call-tree-adams-bashforth-sync}
454  \end{figure}  \end{figure}
455    
# Line 442  G_{\theta,S}^{(n+1/2)} & = & (3/2+\epsil Line 478  G_{\theta,S}^{(n+1/2)} & = & (3/2+\epsil
478  \label{eq:vstar-sync} \\  \label{eq:vstar-sync} \\
479  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )
480  \label{eq:vstarstar-sync} \\  \label{eq:vstarstar-sync} \\
481  \eta^* & = & \epsilon_{fs} \left( \eta^{n} +P-E+R \right)- \Delta t  \eta^* & = & \epsilon_{fs} \left( \eta^{n} + \Delta t (P-E) \right)- \Delta t
482    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }
483  \label{eq:nstar-sync} \\  \label{eq:nstar-sync} \\
484  \nabla \cdot g H \nabla \eta^{n+1} - \frac{\epsilon_{fs} \eta^{n+1}}{\Delta t^2}  \nabla \cdot g H \nabla \eta^{n+1} & - & \frac{\epsilon_{fs} \eta^{n+1}}{\Delta t^2}
485  & = & - \frac{\eta^*}{\Delta t^2}  ~ = ~ - \frac{\eta^*}{\Delta t^2}
486  \label{eq:elliptic-sync} \\  \label{eq:elliptic-sync} \\
487  \vec{\bf v}^{n+1} & = & \vec{\bf v}^{*} - \Delta t g \nabla \eta^{n+1}  \vec{\bf v}^{n+1} & = & \vec{\bf v}^{*} - \Delta t g \nabla \eta^{n+1}
488  \label{eq:v-n+1-sync}  \label{eq:v-n+1-sync}
# Line 465  accelerations, stepping forward and solv Line 501  accelerations, stepping forward and solv
501  surface pressure gradient terms, corresponding to equations  surface pressure gradient terms, corresponding to equations
502  \ref{eq:Gv-n-sync} to \ref{eq:v-n+1-sync}.  \ref{eq:Gv-n-sync} to \ref{eq:v-n+1-sync}.
503  These operations are carried out in subroutines {\em DYNAMCIS}, {\em  These operations are carried out in subroutines {\em DYNAMCIS}, {\em
504  SOLVE\_FOR\_PRESSURE} and {\em THE\_CORRECTION\_STEP}. This, then,  SOLVE\_FOR\_PRESSURE} and {\em MOMENTUM\_CORRECTION\_STEP}. This, then,
505  represents an entire algorithm for stepping forward the model one  represents an entire algorithm for stepping forward the model one
506  time-step. The corresponding calling tree is given in  time-step. The corresponding calling tree is given in
507  \ref{fig:call-tree-adams-bashforth-sync}.  \ref{fig:call-tree-adams-bashforth-sync}.
# Line 480  time-step. The corresponding calling tre Line 516  time-step. The corresponding calling tre
516  \caption{  \caption{
517  A schematic of the explicit Adams-Bashforth and implicit time-stepping  A schematic of the explicit Adams-Bashforth and implicit time-stepping
518  phases of the algorithm but with staggering in time of thermodynamic  phases of the algorithm but with staggering in time of thermodynamic
519  variables with the flow. Explicit thermodynamics tendencies are  variables with the flow.
520  evaluated at time level $n-1/2$ as a function of the thermodynamics  Explicit momentum tendencies are evaluated at time level $n-1/2$ as a
521  state at that time level $n$ and flow at time $n$ (dotted arrow). The  function of the flow field at that time level $n-1/2$.
522  explicit tendency from the previous time level, $n-3/2$, is used to  The explicit tendency from the previous time level, $n-3/2$, is used to
523  extrapolate tendencies to $n$ (dashed arrow). This extrapolated  extrapolate tendencies to $n$ (dashed arrow).
524  tendency allows thermo-dynamics variables to be stably integrated  The hydrostatic pressure/geo-potential $\phi_{hyd}$ is evaluated directly
525  forward-in-time to render an estimate ($*$-variables) at the $n+1/2$  at time level $n$ (vertical arrows) and used with the extrapolated tendencies
526  time level (solid arc-arrow). The implicit-in-time operator ${\cal  to step forward the flow variables from $n-1/2$ to $n+1/2$ (solid arc-arrow).
527  L_{\theta,S}}$ is solved to yield the thermodynamic variables at time  The implicit-in-time operator ${\cal L_{u,v}}$ (vertical arrows) is
528  level $n+1/2$. These are then used to calculate the hydrostatic  then applied to the previous estimation of the the flow field ($*$-variables)
529  pressure/geo-potential, $\phi_{hyd}$ (vertical arrows). The  and yields to the two velocity components $u,v$ at time level $n+1/2$.
530  hydrostatic pressure gradient is evaluated directly an time level  These are then used to calculate the advection term (dashed arc-arrow)
531  $n+1/2$ in stepping forward the flow variables from $n$ to $n+1$  of the thermo-dynamics tendencies at time step $n$.
532  (solid arc-arrow). }  The extrapolated thermodynamics tendency, from time level $n-1$ and $n$
533    to $n+1/2$, allows thermodynamic variables to be stably integrated
534    forward-in-time (solid arc-arrow) up to time level $n+1$.
535    }
536  \label{fig:adams-bashforth-staggered}  \label{fig:adams-bashforth-staggered}
537  \end{figure}  \end{figure}
538    
# Line 503  circumstance, it is more efficient to st Line 542  circumstance, it is more efficient to st
542  thermodynamic variables with the flow  thermodynamic variables with the flow
543  variables. Fig.~\ref{fig:adams-bashforth-staggered} illustrates the  variables. Fig.~\ref{fig:adams-bashforth-staggered} illustrates the
544  staggering and algorithm. The key difference between this and  staggering and algorithm. The key difference between this and
545  Fig.~\ref{fig:adams-bashforth-sync} is that the new thermodynamics  Fig.~\ref{fig:adams-bashforth-sync} is that the thermodynamic variables
546  fields are used to compute the hydrostatic pressure at time level  are solved after the dynamics, using the recently updated flow field.
547  $n+1/2$. The essentially allows the gravity wave terms to leap-frog in  This essentially allows the gravity wave terms to leap-frog in
548  time giving second order accuracy and more stability.  time giving second order accuracy and more stability.
549    
550  The essential change in the staggered algorithm is the calculation of  The essential change in the staggered algorithm is that the
551  hydrostatic pressure which, in the context of the synchronous  thermodynamics solver is delayed from half a time step,
552  algorithm involves replacing equation \ref{eq:phi-hyd-sync} with  allowing the use of the most recent velocities to compute
553  \begin{displaymath}  the advection terms. Once the thermodynamics fields are
554  \phi_{hyd}^n = \int b(\theta^{n+1},S^{n+1}) dr  updated, the hydrostatic pressure is computed
555  \end{displaymath}  to step frowrad the dynamics
556  but the pressure gradient must also be taken out of the  Note that the pressure gradient must also be taken out of the
557  Adams-Bashforth extrapolation. Also, retaining the integer time-levels,  Adams-Bashforth extrapolation. Also, retaining the integer time-levels,
558  $n$ and $n+1$, does not give a user the sense of where variables are  $n$ and $n+1$, does not give a user the sense of where variables are
559  located in time.  Instead, we re-write the entire algorithm,  located in time.  Instead, we re-write the entire algorithm,
560  \ref{eq:Gt-n-sync} to \ref{eq:v-n+1-sync}, annotating the  \ref{eq:Gt-n-sync} to \ref{eq:v-n+1-sync}, annotating the
561  position in time of variables appropriately:  position in time of variables appropriately:
562  \begin{eqnarray}  \begin{eqnarray}
563  G_{\theta,S}^{n-1/2} & = & G_{\theta,S} ( u^n, \theta^{n-1/2}, S^{n-1/2} )  \vec{\bf G}_{\vec{\bf v}}^{n-1/2} & = & \vec{\bf G}_{\vec{\bf v}} ( \vec{\bf v}^{n-1/2} )
 \label{eq:Gt-n-staggered} \\  
 G_{\theta,S}^{(n)} & = & (3/2+\epsilon_{AB}) G_{\theta,S}^{n-1/2}-(1/2+\epsilon_{AB}) G_{\theta,S}^{n-3/2}  
 \label{eq:Gt-n+5-staggered} \\  
 (\theta^*,S^*) & = & (\theta^{n},S^{n}) + \Delta t G_{\theta,S}^{(n)}  
 \label{eq:tstar-staggered} \\  
 (\theta^{n+1/2},S^{n+1/2}) & = & {\cal L}^{-1}_{\theta,S} (\theta^*,S^*)  
 \label{eq:t-n+1-staggered} \\  
 \phi^{n+1/2}_{hyd} & = & \int b(\theta^{n+1/2},S^{n+1/2}) dr  
 \label{eq:phi-hyd-staggered} \\  
 \vec{\bf G}_{\vec{\bf v}}^{n} & = & \vec{\bf G}_{\vec{\bf v}} ( \vec{\bf v}^n )  
564  \label{eq:Gv-n-staggered} \\  \label{eq:Gv-n-staggered} \\
565  \vec{\bf G}_{\vec{\bf v}}^{(n+1/2)} & = & (3/2 + \epsilon_{AB} ) \vec{\bf G}_{\vec{\bf v}}^{n} - (1/2 + \epsilon_{AB} ) \vec{\bf G}_{\vec{\bf v}}^{n-1}  \vec{\bf G}_{\vec{\bf v}}^{(n)} & = & (3/2 + \epsilon_{AB} ) \vec{\bf G}_{\vec{\bf v}}^{n-1/2} - (1/2 + \epsilon_{AB} ) \vec{\bf G}_{\vec{\bf v}}^{n-3/2}
566  \label{eq:Gv-n+5-staggered} \\  \label{eq:Gv-n+5-staggered} \\
567  \vec{\bf v}^{*} & = & \vec{\bf v}^{n} + \Delta t \left( \vec{\bf G}_{\vec{\bf v}}^{(n+1/2)} - \nabla \phi_{hyd}^{n+1/2} \right)  \phi^{n}_{hyd} & = & \int b(\theta^{n},S^{n}) dr
568    \label{eq:phi-hyd-staggered} \\
569    \vec{\bf v}^{*} & = & \vec{\bf v}^{n-1/2} + \Delta t \left( \vec{\bf G}_{\vec{\bf v}}^{(n)} - \nabla \phi_{hyd}^{n} \right)
570  \label{eq:vstar-staggered} \\  \label{eq:vstar-staggered} \\
571  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )
572  \label{eq:vstarstar-staggered} \\  \label{eq:vstarstar-staggered} \\
573  \eta^* & = & \epsilon_{fs} \left( \eta^{n} +P-E+R \right)- \Delta t  \eta^* & = & \epsilon_{fs} \left( \eta^{n-1/2} + \Delta t (P-E)^n \right)- \Delta t
574    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }
575  \label{eq:nstar-staggered} \\  \label{eq:nstar-staggered} \\
576  \nabla \cdot g H \nabla \eta^{n+1} - \frac{\epsilon_{fs} \eta^{n+1}}{\Delta t^2}  \nabla \cdot g H \nabla \eta^{n+1/2} & - & \frac{\epsilon_{fs} \eta^{n+1/2}}{\Delta t^2}
577  & = & - \frac{\eta^*}{\Delta t^2}  ~ = ~ - \frac{\eta^*}{\Delta t^2}
578  \label{eq:elliptic-staggered} \\  \label{eq:elliptic-staggered} \\
579  \vec{\bf v}^{n+1} & = & \vec{\bf v}^{*} - \Delta t g \nabla \eta^{n+1}  \vec{\bf v}^{n+1/2} & = & \vec{\bf v}^{*} - \Delta t g \nabla \eta^{n+1/2}
580  \label{eq:v-n+1-staggered}  \label{eq:v-n+1-staggered} \\
581    G_{\theta,S}^{n} & = & G_{\theta,S} ( u^{n+1/2}, \theta^{n}, S^{n} )
582    \label{eq:Gt-n-staggered} \\
583    G_{\theta,S}^{(n+1/2)} & = & (3/2+\epsilon_{AB}) G_{\theta,S}^{n}-(1/2+\epsilon_{AB}) G_{\theta,S}^{n-1}
584    \label{eq:Gt-n+5-staggered} \\
585    (\theta^*,S^*) & = & (\theta^{n},S^{n}) + \Delta t G_{\theta,S}^{(n+1/2)}
586    \label{eq:tstar-staggered} \\
587    (\theta^{n+1},S^{n+1}) & = & {\cal L}^{-1}_{\theta,S} (\theta^*,S^*)
588    \label{eq:t-n+1-staggered} \\
589  \end{eqnarray}  \end{eqnarray}
590  The calling sequence is unchanged from  The corresponding calling tree is given in
591  Fig.~\ref{fig:call-tree-adams-bashforth-sync}. The staggered algorithm  \ref{fig:call-tree-adams-bashforth-staggered}.
592  is activated with the run-time flag {\bf staggerTimeStep=.TRUE.} in  The staggered algorithm is activated with the run-time flag
593  {\em PARM01} of {\em data}.  {\bf staggerTimeStep=.TRUE.} in parameter file {\em data},
594    namelist {\em PARM01}.
595    
596    \begin{figure}
597    \begin{center} \fbox{ \begin{minipage}{4.7in} \begin{tabbing}
598    aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill
599    FORWARD\_STEP \\
600    \>\> EXTERNAL\_FIELDS\_LOAD\\
601    \>\> DO\_ATMOSPHERIC\_PHYS \\
602    \>\> DO\_OCEANIC\_PHYS \\
603    \> DYNAMICS \\
604    \>\> CALC\_PHI\_HYD \` $\phi_{hyd}^n$ (\ref{eq:phi-hyd-staggered}) \\
605    \>\> MOM\_FLUXFORM or MOM\_VECINV \` $G_{\vec{\bf v}}^{n-1/2}$
606        (\ref{eq:Gv-n-staggered})\\
607    \>\> TIMESTEP \` $\vec{\bf v}^*$ (\ref{eq:Gv-n+5-staggered},
608                                      \ref{eq:vstar-staggered}) \\
609    \>\> IMPLDIFF \` $\vec{\bf v}^{**}$ (\ref{eq:vstarstar-staggered}) \\
610    \> UPDATE\_R\_STAR or UPDATE\_SURF\_DR \` (NonLin-FS only)\\
611    \> SOLVE\_FOR\_PRESSURE \\
612    \>\> CALC\_DIV\_GHAT \` $\eta^*$ (\ref{eq:nstar-staggered}) \\
613    \>\> CG2D \` $\eta^{n+1/2}$ (\ref{eq:elliptic-staggered}) \\
614    \> MOMENTUM\_CORRECTION\_STEP  \\
615    \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1/2}$ \\
616    \>\> CORRECTION\_STEP \` $u^{n+1/2}$,$v^{n+1/2}$ (\ref{eq:v-n+1-staggered})\\
617    \> THERMODYNAMICS \\
618    \>\> CALC\_GT \\
619    \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$
620         (\ref{eq:Gt-n-staggered})\\
621    \>\>\> EXTERNAL\_FORCING \` $G_\theta^n = G_\theta^n + {\cal Q}$ \\
622    \>\>\> ADAMS\_BASHFORTH2 \` $G_\theta^{(n+1/2)}$ (\ref{eq:Gt-n+5-staggered}) \\
623    \>\> TIMESTEP\_TRACER \` $\tau^*$ (\ref{eq:tstar-staggered}) \\
624    \>\> IMPLDIFF \` $\tau^{(n+1)}$ (\ref{eq:t-n+1-staggered}) \\
625    \> TRACERS\_CORRECTION\_STEP  \\
626    \>\> CYCLE\_TRACER \` $\theta^{n+1}$ \\
627    \>\> FILTER \` Shapiro Filter, Zonal Filter (FFT) \\
628    \>\> CONVECTIVE\_ADJUSTMENT \` \\
629    \end{tabbing} \end{minipage} } \end{center}
630    \caption{
631    Calling tree for the overall staggered algorithm using
632    Adams-Bashforth time-stepping.
633    The place where the model geometry
634    ({\em hFac} factors) is updated is added here but is only relevant
635    for the non-linear free-surface algorithm.
636    }
637    \label{fig:call-tree-adams-bashforth-staggered}
638    \end{figure}
639    
640  The only difficulty with this approach is apparent in equation  The only difficulty with this approach is apparent in equation
641  \ref{eq:Gt-n-staggered} and illustrated by the dotted arrow  \ref{eq:Gt-n-staggered} and illustrated by the dotted arrow
642  connecting $u,v^n$ with $G_\theta^{n-1/2}$. The flow used to advect  connecting $u,v^{n+1/2}$ with $G_\theta^{n}$. The flow used to advect
643  tracers around is not naturally located in time. This could be avoided  tracers around is not naturally located in time. This could be avoided
644  by applying the Adams-Bashforth extrapolation to the tracer field  by applying the Adams-Bashforth extrapolation to the tracer field
645  itself and advecting that around but this approach is not yet  itself and advecting that around but this approach is not yet
# Line 638  the following equations: Line 722  the following equations:
722  u^{*} & = & u^{n} + \Delta t G_u^{(n+1/2)} \label{eq:ustar-nh} \\  u^{*} & = & u^{n} + \Delta t G_u^{(n+1/2)} \label{eq:ustar-nh} \\
723  v^{*} & = & v^{n} + \Delta t G_v^{(n+1/2)} \label{eq:vstar-nh} \\  v^{*} & = & v^{n} + \Delta t G_v^{(n+1/2)} \label{eq:vstar-nh} \\
724  w^{*} & = & w^{n} + \Delta t G_w^{(n+1/2)} \label{eq:wstar-nh} \\  w^{*} & = & w^{n} + \Delta t G_w^{(n+1/2)} \label{eq:wstar-nh} \\
725  \eta^* & = & \epsilon_{fs} \left( \eta^{n} +P-E+R \right)- \Delta t  \eta^* ~ = ~ \epsilon_{fs} \left( \eta^{n} + \Delta t (P-E) \right)
726    & - & \Delta t
727    \partial_x H \widehat{u^{*}}    \partial_x H \widehat{u^{*}}
728  + \partial_y H \widehat{v^{*}}  + \partial_y H \widehat{v^{*}}
729  \\  \\
730    \partial_x g H \partial_x \eta^{n+1}    \partial_x g H \partial_x \eta^{n+1}
731  + \partial_y g H \partial_y \eta^{n+1}  + \partial_y g H \partial_y \eta^{n+1}
732  - \frac{\epsilon_{fs} \eta^{n+1}}{\Delta t^2}  & - & \frac{\epsilon_{fs} \eta^{n+1}}{\Delta t^2}
733  & = &  ~ = ~
734  - \frac{\eta^*}{\Delta t^2}  - \frac{\eta^*}{\Delta t^2}
735  \label{eq:elliptic-nh}  \label{eq:elliptic-nh}
736  \\  \\

Legend:
Removed from v.1.14  
changed lines
  Added in v.1.18

  ViewVC Help
Powered by ViewVC 1.1.22