/[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.15 by cnh, Thu Feb 28 19:32:19 2002 UTC revision 1.19 by edhill, Sat Oct 16 03:40:12 2004 UTC
# Line 12  passive and dynamically active tracers a Line 12  passive and dynamically active tracers a
12    
13    
14  \section{Time-stepping}  \section{Time-stepping}
15    \begin{rawhtml}
16    <!-- CMIREDIR:time-stepping: -->
17    \end{rawhtml}
18    
19  The equations of motion integrated by the model involve four  The equations of motion integrated by the model involve four
20  prognostic equations for flow, $u$ and $v$, temperature, $\theta$, and  prognostic equations for flow, $u$ and $v$, temperature, $\theta$, and
21  salt/moisture, $S$, and three diagnostic equations for vertical flow,  salt/moisture, $S$, and three diagnostic equations for vertical flow,
# Line 61  treated more exactly, including non-line Line 65  treated more exactly, including non-line
65  described in section \ref{sect:nonlinear-freesurface}.  described in section \ref{sect:nonlinear-freesurface}.
66    
67    
68  \section{Pressure method with rigid-lid} \label{sect:pressure-method-rigid-lid}  \section{Pressure method with rigid-lid}
69    \label{sect:pressure-method-rigid-lid}
70    \begin{rawhtml}
71    <!-- CMIREDIR:pressure_method_rigid_lid: -->
72    \end{rawhtml}
73    
74  \begin{figure}  \begin{figure}
75  \begin{center}  \begin{center}
# Line 81  temporary.} Line 89  temporary.}
89  \begin{figure}  \begin{figure}
90  \begin{center} \fbox{ \begin{minipage}{4.5in} \begin{tabbing}  \begin{center} \fbox{ \begin{minipage}{4.5in} \begin{tabbing}
91  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill
92  \proclink{FORWARD\_STEP}{../code/._model_src_forward_step.F} \\  \filelink{FORWARD\_STEP}{model-src-forward_step.F} \\
93  \> DYNAMICS \\  \> DYNAMICS \\
94  \>\> 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}) \\
95  \> SOLVE\_FOR\_PRESSURE \\  \> SOLVE\_FOR\_PRESSURE \\
96  \>\> CALC\_DIV\_GHAT \` $H\widehat{u^*}$,$H\widehat{v^*}$ (\ref{eq:elliptic}) \\  \>\> CALC\_DIV\_GHAT \` $H\widehat{u^*}$,$H\widehat{v^*}$ (\ref{eq:elliptic}) \\
97  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic}) \\  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic}) \\
98  \> THE\_CORRECTION\_STEP  \\  \> MOMENTUM\_CORRECTION\_STEP  \\
99  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\
100  \>\> 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})
101  \end{tabbing} \end{minipage} } \end{center}  \end{tabbing} \end{minipage} } \end{center}
102  \caption{Calling tree for the pressure method algorihtm}  \caption{Calling tree for the pressure method algorithm
103      (\filelink{FORWARD\_STEP}{model-src-forward_step.F})}
104  \label{fig:call-tree-pressure-method}  \label{fig:call-tree-pressure-method}
105  \end{figure}  \end{figure}
106    
# Line 173  The correspondence to the code is as fol Line 182  The correspondence to the code is as fol
182  \item  \item
183  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},
184  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
185  \proclink{TIMESTEP}{../code/._model_src_timestep.F}  \filelink{TIMESTEP()}{model-src-timestep.F}
186  \item  \item
187  the vertical integration, $H \widehat{u^*}$ and $H  the vertical integration, $H \widehat{u^*}$ and $H
188  \widehat{v^*}$, divergence and inversion of the elliptic operator in  \widehat{v^*}$, divergence and inversion of the elliptic operator in
189  equation \ref{eq:elliptic} is coded in  equation \ref{eq:elliptic} is coded in
190  \proclink{SOLVE\_FOR\_PRESSURE}{../code/._model_src_solve_for_pressure.F}  \filelink{SOLVE\_FOR\_PRESSURE()}{model-src-solve_for_pressure.F}
191  \item  \item
192  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
193  \ref{eq:un+1-rigid-lid} and \ref{eq:vn+1-rigid-lid} is calculated in  \ref{eq:un+1-rigid-lid} and \ref{eq:vn+1-rigid-lid} is calculated in
194  \proclink{CORRECTION\_STEP}{../code/._model_src_correction_step.F}.  \filelink{CORRECTION\_STEP()}{model-src-correction_step.F}.
195  \end{itemize}  \end{itemize}
196  The calling tree for these routines is given in  The calling tree for these routines is given in
197  Fig.~\ref{fig:call-tree-pressure-method}.  Fig.~\ref{fig:call-tree-pressure-method}.
# Line 202  G_u^{(n+1/2)} Line 211  G_u^{(n+1/2)}
211    
212  \section{Pressure method with implicit linear free-surface}  \section{Pressure method with implicit linear free-surface}
213  \label{sect:pressure-method-linear-backward}  \label{sect:pressure-method-linear-backward}
214    \begin{rawhtml}
215    <!-- CMIREDIR:pressure_method_linear_backward: -->
216    \end{rawhtml}
217    
218  The rigid-lid approximation filters out external gravity waves  The rigid-lid approximation filters out external gravity waves
219  subsequently modifying the dispersion relation of barotropic Rossby  subsequently modifying the dispersion relation of barotropic Rossby
# Line 278  pressure-method. Line 290  pressure-method.
290    
291  \section{Explicit time-stepping: Adams-Bashforth}  \section{Explicit time-stepping: Adams-Bashforth}
292  \label{sect:adams-bashforth}  \label{sect:adams-bashforth}
293    \begin{rawhtml}
294    <!-- CMIREDIR:adams_bashforth: -->
295    \end{rawhtml}
296    
297  In describing the the pressure method above we deferred describing the  In describing the the pressure method above we deferred describing the
298  time discretization of the explicit terms. We have historically used  time discretization of the explicit terms. We have historically used
# Line 293  FORWARD\_STEP \\ Line 308  FORWARD\_STEP \\
308  \> THERMODYNAMICS \\  \> THERMODYNAMICS \\
309  \>\> CALC\_GT \\  \>\> CALC\_GT \\
310  \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$ \\  \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$ \\
311  \>\>\> EXTERNAL\_FORCING \` $G_\theta^n = G_\theta^n + {\cal Q}$ \\  \>either\>\> EXTERNAL\_FORCING \` $G_\theta^n = G_\theta^n + {\cal Q}$ \\
312  \>\>\> ADAMS\_BASHFORTH2 \` $G_\theta^{(n+1/2)}$ (\ref{eq:adams-bashforth2}) \\  \>\>\> ADAMS\_BASHFORTH2 \` $G_\theta^{(n+1/2)}$ (\ref{eq:adams-bashforth2}) \\
313    \>or\>\> EXTERNAL\_FORCING \` $G_\theta^{(n+1/2)} = G_\theta^{(n+1/2)} + {\cal Q}$ \\
314  \>\> TIMESTEP\_TRACER \` $\tau^*$ (\ref{eq:taustar}) \\  \>\> TIMESTEP\_TRACER \` $\tau^*$ (\ref{eq:taustar}) \\
315  \>\> IMPLDIFF \` $\tau^{(n+1)}$ (\ref{eq:tau-n+1-implicit})  \>\> IMPLDIFF \` $\tau^{(n+1)}$ (\ref{eq:tau-n+1-implicit})
316  \end{tabbing} \end{minipage} } \end{center}  \end{tabbing} \end{minipage} } \end{center}
# Line 331  simpler to include these terms and this Line 347  simpler to include these terms and this
347  and forcing evolves smoothly. Problems can, and do, arise when forcing  and forcing evolves smoothly. Problems can, and do, arise when forcing
348  or motions are high frequency and this corresponds to a reduced  or motions are high frequency and this corresponds to a reduced
349  stability compared to a simple forward time-stepping of such terms.  stability compared to a simple forward time-stepping of such terms.
350    The model offers the possibility to leave the forcing term outside the
351    Adams-Bashforth extrapolation, by turning off the logical flag
352    {\bf forcing\_In\_AB } (parameter file {\em data}, namelist {\em PARM01},
353    default value = True).
354    
355  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.
356  \marginpar{AJA needs to find his notes on this...}  \marginpar{AJA needs to find his notes on this...}
# Line 340  A stability analysis for a relaxation eq Line 360  A stability analysis for a relaxation eq
360    
361    
362  \section{Implicit time-stepping: backward method}  \section{Implicit time-stepping: backward method}
363    \begin{rawhtml}
364    <!-- CMIREDIR:implicit_time-stepping_backward: -->
365    \end{rawhtml}
366    
367  Vertical diffusion and viscosity can be treated implicitly in time  Vertical diffusion and viscosity can be treated implicitly in time
368  using the backward method which is an intrinsic scheme. For tracers,  using the backward method which is an intrinsic scheme.
369    Recently, the option to treat the vertical advection
370    implicitly has been added, but not yet tested; therefore,
371    the description hereafter is limited to diffusion and viscosity.
372    For tracers,
373  the time discretized equation is:  the time discretized equation is:
374  \begin{equation}  \begin{equation}
375  \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 373  Fig.~\ref{fig:call-tree-adams-bashforth} Line 400  Fig.~\ref{fig:call-tree-adams-bashforth}
400  stepping forward a tracer variable such as temperature.  stepping forward a tracer variable such as temperature.
401    
402  In order to fit within the pressure method, the implicit viscosity  In order to fit within the pressure method, the implicit viscosity
403  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
404  redistribute momentum in the vertical. The upshot of this is that  redistribute momentum in the vertical. The upshot of this is that
405  although vertical viscosity may be backward implicit and  although vertical viscosity may be backward implicit and
406  unconditionally stable, no-slip boundary conditions may not be made  unconditionally stable, no-slip boundary conditions may not be made
# Line 381  implicit and are thus cast as a an expli Line 408  implicit and are thus cast as a an expli
408    
409  \section{Synchronous time-stepping: variables co-located in time}  \section{Synchronous time-stepping: variables co-located in time}
410  \label{sect:adams-bashforth-sync}  \label{sect:adams-bashforth-sync}
411    \begin{rawhtml}
412    <!-- CMIREDIR:adams_bashforth_sync: -->
413    \end{rawhtml}
414    
415  \begin{figure}  \begin{figure}
416  \begin{center}  \begin{center}
# Line 401  is solved to yield the state variables a Line 431  is solved to yield the state variables a
431  \end{figure}  \end{figure}
432    
433  \begin{figure}  \begin{figure}
434  \begin{center} \fbox{ \begin{minipage}{4.5in} \begin{tabbing}  \begin{center} \fbox{ \begin{minipage}{4.7in} \begin{tabbing}
435  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill  aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill
436  FORWARD\_STEP \\  FORWARD\_STEP \\
437    \>\> EXTERNAL\_FIELDS\_LOAD\\
438    \>\> DO\_ATMOSPHERIC\_PHYS \\
439    \>\> DO\_OCEANIC\_PHYS \\
440  \> THERMODYNAMICS \\  \> THERMODYNAMICS \\
441  \>\> CALC\_GT \\  \>\> CALC\_GT \\
442  \>\>\> 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 416  FORWARD\_STEP \\ Line 449  FORWARD\_STEP \\
449  \>\> 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})\\
450  \>\> 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}) \\
451  \>\> IMPLDIFF \` $\vec{\bf v}^{**}$ (\ref{eq:vstarstar-sync}) \\  \>\> IMPLDIFF \` $\vec{\bf v}^{**}$ (\ref{eq:vstarstar-sync}) \\
452    \> UPDATE\_R\_STAR or UPDATE\_SURF\_DR \` (NonLin-FS only)\\
453  \> SOLVE\_FOR\_PRESSURE \\  \> SOLVE\_FOR\_PRESSURE \\
454  \>\> CALC\_DIV\_GHAT \` $\eta^*$ (\ref{eq:nstar-sync}) \\  \>\> CALC\_DIV\_GHAT \` $\eta^*$ (\ref{eq:nstar-sync}) \\
455  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic-sync}) \\  \>\> CG2D \` $\eta^{n+1}$ (\ref{eq:elliptic-sync}) \\
456  \> THE\_CORRECTION\_STEP  \\  \> MOMENTUM\_CORRECTION\_STEP  \\
457  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\  \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1}$ \\
458  \>\> 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})\\
459    \> TRACERS\_CORRECTION\_STEP  \\
460    \>\> CYCLE\_TRACER \` $\theta^{n+1}$ \\
461    \>\> FILTER \` Shapiro Filter, Zonal Filter (FFT) \\
462    \>\> CONVECTIVE\_ADJUSTMENT \` \\
463  \end{tabbing} \end{minipage} } \end{center}  \end{tabbing} \end{minipage} } \end{center}
464  \caption{  \caption{
465  Calling tree for the overall synchronous algorithm using  Calling tree for the overall synchronous algorithm using
466  Adams-Bashforth time-stepping.}  Adams-Bashforth time-stepping.
467    The place where the model geometry
468    ({\em hFac} factors) is updated is added here but is only relevant
469    for the non-linear free-surface algorithm.
470    For completeness, the external forcing,
471    ocean and atmospheric physics have been added, although they are mainly
472    optional}
473  \label{fig:call-tree-adams-bashforth-sync}  \label{fig:call-tree-adams-bashforth-sync}
474  \end{figure}  \end{figure}
475    
# Line 454  G_{\theta,S}^{(n+1/2)} & = & (3/2+\epsil Line 498  G_{\theta,S}^{(n+1/2)} & = & (3/2+\epsil
498  \label{eq:vstar-sync} \\  \label{eq:vstar-sync} \\
499  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )
500  \label{eq:vstarstar-sync} \\  \label{eq:vstarstar-sync} \\
501  \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
502    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }
503  \label{eq:nstar-sync} \\  \label{eq:nstar-sync} \\
504  \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}
505  & = & - \frac{\eta^*}{\Delta t^2}  ~ = ~ - \frac{\eta^*}{\Delta t^2}
506  \label{eq:elliptic-sync} \\  \label{eq:elliptic-sync} \\
507  \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}
508  \label{eq:v-n+1-sync}  \label{eq:v-n+1-sync}
# Line 477  accelerations, stepping forward and solv Line 521  accelerations, stepping forward and solv
521  surface pressure gradient terms, corresponding to equations  surface pressure gradient terms, corresponding to equations
522  \ref{eq:Gv-n-sync} to \ref{eq:v-n+1-sync}.  \ref{eq:Gv-n-sync} to \ref{eq:v-n+1-sync}.
523  These operations are carried out in subroutines {\em DYNAMCIS}, {\em  These operations are carried out in subroutines {\em DYNAMCIS}, {\em
524  SOLVE\_FOR\_PRESSURE} and {\em THE\_CORRECTION\_STEP}. This, then,  SOLVE\_FOR\_PRESSURE} and {\em MOMENTUM\_CORRECTION\_STEP}. This, then,
525  represents an entire algorithm for stepping forward the model one  represents an entire algorithm for stepping forward the model one
526  time-step. The corresponding calling tree is given in  time-step. The corresponding calling tree is given in
527  \ref{fig:call-tree-adams-bashforth-sync}.  \ref{fig:call-tree-adams-bashforth-sync}.
528    
529  \section{Staggered baroclinic time-stepping}  \section{Staggered baroclinic time-stepping}
530  \label{sect:adams-bashforth-staggered}  \label{sect:adams-bashforth-staggered}
531    \begin{rawhtml}
532    <!-- CMIREDIR:adams_bashforth_staggered: -->
533    \end{rawhtml}
534    
535  \begin{figure}  \begin{figure}
536  \begin{center}  \begin{center}
# Line 492  time-step. The corresponding calling tre Line 539  time-step. The corresponding calling tre
539  \caption{  \caption{
540  A schematic of the explicit Adams-Bashforth and implicit time-stepping  A schematic of the explicit Adams-Bashforth and implicit time-stepping
541  phases of the algorithm but with staggering in time of thermodynamic  phases of the algorithm but with staggering in time of thermodynamic
542  variables with the flow. Explicit thermodynamics tendencies are  variables with the flow.
543  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
544  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$.
545  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
546  extrapolate tendencies to $n$ (dashed arrow). This extrapolated  extrapolate tendencies to $n$ (dashed arrow).
547  tendency allows thermo-dynamics variables to be stably integrated  The hydrostatic pressure/geo-potential $\phi_{hyd}$ is evaluated directly
548  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
549  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).
550  L_{\theta,S}}$ is solved to yield the thermodynamic variables at time  The implicit-in-time operator ${\cal L_{u,v}}$ (vertical arrows) is
551  level $n+1/2$. These are then used to calculate the hydrostatic  then applied to the previous estimation of the the flow field ($*$-variables)
552  pressure/geo-potential, $\phi_{hyd}$ (vertical arrows). The  and yields to the two velocity components $u,v$ at time level $n+1/2$.
553  hydrostatic pressure gradient is evaluated directly an time level  These are then used to calculate the advection term (dashed arc-arrow)
554  $n+1/2$ in stepping forward the flow variables from $n$ to $n+1$  of the thermo-dynamics tendencies at time step $n$.
555  (solid arc-arrow). }  The extrapolated thermodynamics tendency, from time level $n-1$ and $n$
556    to $n+1/2$, allows thermodynamic variables to be stably integrated
557    forward-in-time (solid arc-arrow) up to time level $n+1$.
558    }
559  \label{fig:adams-bashforth-staggered}  \label{fig:adams-bashforth-staggered}
560  \end{figure}  \end{figure}
561    
# Line 515  circumstance, it is more efficient to st Line 565  circumstance, it is more efficient to st
565  thermodynamic variables with the flow  thermodynamic variables with the flow
566  variables. Fig.~\ref{fig:adams-bashforth-staggered} illustrates the  variables. Fig.~\ref{fig:adams-bashforth-staggered} illustrates the
567  staggering and algorithm. The key difference between this and  staggering and algorithm. The key difference between this and
568  Fig.~\ref{fig:adams-bashforth-sync} is that the new thermodynamics  Fig.~\ref{fig:adams-bashforth-sync} is that the thermodynamic variables
569  fields are used to compute the hydrostatic pressure at time level  are solved after the dynamics, using the recently updated flow field.
570  $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
571  time giving second order accuracy and more stability.  time giving second order accuracy and more stability.
572    
573  The essential change in the staggered algorithm is the calculation of  The essential change in the staggered algorithm is that the
574  hydrostatic pressure which, in the context of the synchronous  thermodynamics solver is delayed from half a time step,
575  algorithm involves replacing equation \ref{eq:phi-hyd-sync} with  allowing the use of the most recent velocities to compute
576  \begin{displaymath}  the advection terms. Once the thermodynamics fields are
577  \phi_{hyd}^n = \int b(\theta^{n+1},S^{n+1}) dr  updated, the hydrostatic pressure is computed
578  \end{displaymath}  to step frowrad the dynamics
579  but the pressure gradient must also be taken out of the  Note that the pressure gradient must also be taken out of the
580  Adams-Bashforth extrapolation. Also, retaining the integer time-levels,  Adams-Bashforth extrapolation. Also, retaining the integer time-levels,
581  $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
582  located in time.  Instead, we re-write the entire algorithm,  located in time.  Instead, we re-write the entire algorithm,
583  \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
584  position in time of variables appropriately:  position in time of variables appropriately:
585  \begin{eqnarray}  \begin{eqnarray}
586  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 )  
587  \label{eq:Gv-n-staggered} \\  \label{eq:Gv-n-staggered} \\
588  \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}
589  \label{eq:Gv-n+5-staggered} \\  \label{eq:Gv-n+5-staggered} \\
590  \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
591    \label{eq:phi-hyd-staggered} \\
592    \vec{\bf v}^{*} & = & \vec{\bf v}^{n-1/2} + \Delta t \left( \vec{\bf G}_{\vec{\bf v}}^{(n)} - \nabla \phi_{hyd}^{n} \right)
593  \label{eq:vstar-staggered} \\  \label{eq:vstar-staggered} \\
594  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )  \vec{\bf v}^{**} & = & {\cal L}_{\vec{\bf v}}^{-1} ( \vec{\bf v}^* )
595  \label{eq:vstarstar-staggered} \\  \label{eq:vstarstar-staggered} \\
596  \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
597    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }    \nabla \cdot H \widehat{ \vec{\bf v}^{**} }
598  \label{eq:nstar-staggered} \\  \label{eq:nstar-staggered} \\
599  \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}
600  & = & - \frac{\eta^*}{\Delta t^2}  ~ = ~ - \frac{\eta^*}{\Delta t^2}
601  \label{eq:elliptic-staggered} \\  \label{eq:elliptic-staggered} \\
602  \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}
603  \label{eq:v-n+1-staggered}  \label{eq:v-n+1-staggered} \\
604    G_{\theta,S}^{n} & = & G_{\theta,S} ( u^{n+1/2}, \theta^{n}, S^{n} )
605    \label{eq:Gt-n-staggered} \\
606    G_{\theta,S}^{(n+1/2)} & = & (3/2+\epsilon_{AB}) G_{\theta,S}^{n}-(1/2+\epsilon_{AB}) G_{\theta,S}^{n-1}
607    \label{eq:Gt-n+5-staggered} \\
608    (\theta^*,S^*) & = & (\theta^{n},S^{n}) + \Delta t G_{\theta,S}^{(n+1/2)}
609    \label{eq:tstar-staggered} \\
610    (\theta^{n+1},S^{n+1}) & = & {\cal L}^{-1}_{\theta,S} (\theta^*,S^*)
611    \label{eq:t-n+1-staggered} \\
612  \end{eqnarray}  \end{eqnarray}
613  The calling sequence is unchanged from  The corresponding calling tree is given in
614  Fig.~\ref{fig:call-tree-adams-bashforth-sync}. The staggered algorithm  \ref{fig:call-tree-adams-bashforth-staggered}.
615  is activated with the run-time flag {\bf staggerTimeStep=.TRUE.} in  The staggered algorithm is activated with the run-time flag
616  {\em PARM01} of {\em data}.  {\bf staggerTimeStep=.TRUE.} in parameter file {\em data},
617    namelist {\em PARM01}.
618    
619    \begin{figure}
620    \begin{center} \fbox{ \begin{minipage}{4.7in} \begin{tabbing}
621    aaa \= aaa \= aaa \= aaa \= aaa \= aaa \kill
622    FORWARD\_STEP \\
623    \>\> EXTERNAL\_FIELDS\_LOAD\\
624    \>\> DO\_ATMOSPHERIC\_PHYS \\
625    \>\> DO\_OCEANIC\_PHYS \\
626    \> DYNAMICS \\
627    \>\> CALC\_PHI\_HYD \` $\phi_{hyd}^n$ (\ref{eq:phi-hyd-staggered}) \\
628    \>\> MOM\_FLUXFORM or MOM\_VECINV \` $G_{\vec{\bf v}}^{n-1/2}$
629        (\ref{eq:Gv-n-staggered})\\
630    \>\> TIMESTEP \` $\vec{\bf v}^*$ (\ref{eq:Gv-n+5-staggered},
631                                      \ref{eq:vstar-staggered}) \\
632    \>\> IMPLDIFF \` $\vec{\bf v}^{**}$ (\ref{eq:vstarstar-staggered}) \\
633    \> UPDATE\_R\_STAR or UPDATE\_SURF\_DR \` (NonLin-FS only)\\
634    \> SOLVE\_FOR\_PRESSURE \\
635    \>\> CALC\_DIV\_GHAT \` $\eta^*$ (\ref{eq:nstar-staggered}) \\
636    \>\> CG2D \` $\eta^{n+1/2}$ (\ref{eq:elliptic-staggered}) \\
637    \> MOMENTUM\_CORRECTION\_STEP  \\
638    \>\> CALC\_GRAD\_PHI\_SURF \` $\nabla \eta^{n+1/2}$ \\
639    \>\> CORRECTION\_STEP \` $u^{n+1/2}$,$v^{n+1/2}$ (\ref{eq:v-n+1-staggered})\\
640    \> THERMODYNAMICS \\
641    \>\> CALC\_GT \\
642    \>\>\> GAD\_CALC\_RHS \` $G_\theta^n = G_\theta( u, \theta^n )$
643         (\ref{eq:Gt-n-staggered})\\
644    \>\>\> EXTERNAL\_FORCING \` $G_\theta^n = G_\theta^n + {\cal Q}$ \\
645    \>\>\> ADAMS\_BASHFORTH2 \` $G_\theta^{(n+1/2)}$ (\ref{eq:Gt-n+5-staggered}) \\
646    \>\> TIMESTEP\_TRACER \` $\tau^*$ (\ref{eq:tstar-staggered}) \\
647    \>\> IMPLDIFF \` $\tau^{(n+1)}$ (\ref{eq:t-n+1-staggered}) \\
648    \> TRACERS\_CORRECTION\_STEP  \\
649    \>\> CYCLE\_TRACER \` $\theta^{n+1}$ \\
650    \>\> FILTER \` Shapiro Filter, Zonal Filter (FFT) \\
651    \>\> CONVECTIVE\_ADJUSTMENT \` \\
652    \end{tabbing} \end{minipage} } \end{center}
653    \caption{
654    Calling tree for the overall staggered algorithm using
655    Adams-Bashforth time-stepping.
656    The place where the model geometry
657    ({\em hFac} factors) is updated is added here but is only relevant
658    for the non-linear free-surface algorithm.
659    }
660    \label{fig:call-tree-adams-bashforth-staggered}
661    \end{figure}
662    
663  The only difficulty with this approach is apparent in equation  The only difficulty with this approach is apparent in equation
664  \ref{eq:Gt-n-staggered} and illustrated by the dotted arrow  \ref{eq:Gt-n-staggered} and illustrated by the dotted arrow
665  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
666  tracers around is not naturally located in time. This could be avoided  tracers around is not naturally located in time. This could be avoided
667  by applying the Adams-Bashforth extrapolation to the tracer field  by applying the Adams-Bashforth extrapolation to the tracer field
668  itself and advecting that around but this approach is not yet  itself and advecting that around but this approach is not yet
# Line 578  time-level variables and terms correspon Line 673  time-level variables and terms correspon
673    
674  \section{Non-hydrostatic formulation}  \section{Non-hydrostatic formulation}
675  \label{sect:non-hydrostatic}  \label{sect:non-hydrostatic}
676    \begin{rawhtml}
677    <!-- CMIREDIR:non-hydrostatic_formulation: -->
678    \end{rawhtml}
679    
680  The non-hydrostatic formulation re-introduces the full vertical  The non-hydrostatic formulation re-introduces the full vertical
681  momentum equation and requires the solution of a 3-D elliptic  momentum equation and requires the solution of a 3-D elliptic
# Line 650  the following equations: Line 748  the following equations:
748  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} \\
749  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} \\
750  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} \\
751  \eta^* & = & \epsilon_{fs} \left( \eta^{n} +P-E+R \right)- \Delta t  \eta^* ~ = ~ \epsilon_{fs} \left( \eta^{n} + \Delta t (P-E) \right)
752    & - & \Delta t
753    \partial_x H \widehat{u^{*}}    \partial_x H \widehat{u^{*}}
754  + \partial_y H \widehat{v^{*}}  + \partial_y H \widehat{v^{*}}
755  \\  \\
756    \partial_x g H \partial_x \eta^{n+1}    \partial_x g H \partial_x \eta^{n+1}
757  + \partial_y g H \partial_y \eta^{n+1}  + \partial_y g H \partial_y \eta^{n+1}
758  - \frac{\epsilon_{fs} \eta^{n+1}}{\Delta t^2}  & - & \frac{\epsilon_{fs} \eta^{n+1}}{\Delta t^2}
759  & = &  ~ = ~
760  - \frac{\eta^*}{\Delta t^2}  - \frac{\eta^*}{\Delta t^2}
761  \label{eq:elliptic-nh}  \label{eq:elliptic-nh}
762  \\  \\

Legend:
Removed from v.1.15  
changed lines
  Added in v.1.19

  ViewVC Help
Powered by ViewVC 1.1.22