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

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

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

revision 1.3 by adcroft, Tue Sep 25 20:13:42 2001 UTC revision 1.7 by cnh, Wed Oct 24 23:14:44 2001 UTC
# Line 2  Line 2 
2  % $Name$  % $Name$
3    
4  \section{Tracer equations}  \section{Tracer equations}
5    \label{sec:tracer_equations}
6    
7  The basic discretization used for the tracer equations is the second  The basic discretization used for the tracer equations is the second
8  order piece-wise constant finite volume form of the forced  order piece-wise constant finite volume form of the forced
# Line 14  part of the tracer equations and the var Line 15  part of the tracer equations and the var
15  described here.  described here.
16    
17  \subsection{Time-stepping of tracers: ABII}  \subsection{Time-stepping of tracers: ABII}
18    \label{sec:tracer_equations_abII}
19    
20  The default advection scheme is the centered second order method which  The default advection scheme is the centered second order method which
21  requires a second order or quasi-second order time-stepping scheme to  requires a second order or quasi-second order time-stepping scheme to
# Line 66  $F_r$: {\bf fVerT} (argument) Line 68  $F_r$: {\bf fVerT} (argument)
68    
69  \end{minipage} }  \end{minipage} }
70    
   
71  The space and time discretizations are treated seperately (method of  The space and time discretizations are treated seperately (method of
72  lines). The Adams-Bashforth time discretization reads:  lines). Tendancies are calculated at time levels $n$ and $n-1$ and
73    extrapolated to $n+1/2$ using the Adams-Bashforth method:
74  \marginpar{$\epsilon$: {\bf AB\_eps}}  \marginpar{$\epsilon$: {\bf AB\_eps}}
 \marginpar{$\Delta t$: {\bf deltaTtracer}}  
75  \begin{equation}  \begin{equation}
76  \tau^{(n+1)} = \tau^{(n)} + \Delta t \left(  G^{(n+1/2)} =
77  (\frac{3}{2} + \epsilon) G^{(n)} - (\frac{1}{2} + \epsilon) G^{(n-1)}  (\frac{3}{2} + \epsilon) G^{(n)} - (\frac{1}{2} + \epsilon) G^{(n-1)}
 \right)  
78  \end{equation}  \end{equation}
79  where $G^{(n)} = G_{adv}^\tau + G_{diff}^\tau + G_{src}^\tau$ at time  where $G^{(n)} = G_{adv}^\tau + G_{diff}^\tau + G_{src}^\tau$ at time
80  step $n$.  step $n$. The tendancy at $n-1$ is not re-calculated but rather the
81    tendancy at $n$ is stored in a global array for later re-use.
82    
83  Strictly speaking the ABII scheme should be applied only to the  \fbox{ \begin{minipage}{4.75in}
84  advection terms. However, this scheme is only used in conjuction with  {\em S/R ADAMS\_BASHFORTH2} ({\em model/src/adams\_bashforth2.F})
85  the standard second, third and fourth order advection  
86  schemes. Selection of any other advection scheme disables  $G^{(n+1/2)}$: {\bf gTracer} (argument on exit)
87  Adams-Bashforth for tracers so that explicit diffusion and forcing use  
88  the forward method.  $G^{(n)}$: {\bf gTracer} (argument on entry)
89    
90    $G^{(n-1)}$: {\bf gTrNm1} (argument)
91    
92    $\epsilon$: {\bf ABeps} (PARAMS.h)
93    
94    \end{minipage} }
95    
96    The tracers are stepped forward in time using the extrapolated tendancy:
97    \begin{equation}
98    \tau^{(n+1)} = \tau^{(n)} + \Delta t G^{(n+1/2)}
99    \end{equation}
100    \marginpar{$\Delta t$: {\bf deltaTtracer}}
101    
102  \fbox{ \begin{minipage}{4.75in}  \fbox{ \begin{minipage}{4.75in}
103  {\em S/R TIMESTEP\_TRACER} ({\em model/src/timestep\_tracer.F})  {\em S/R TIMESTEP\_TRACER} ({\em model/src/timestep\_tracer.F})
104    
105  $\tau$: {\bf tracer} (argument)  $\tau^{(n+1)}$: {\bf gTracer} (argument on exit)
106    
107  $G^{(n)}$: {\bf gTracer} (argument)  $\tau^{(n)}$: {\bf tracer} (argument on entry)
108    
109  $G^{(n-1)}$: {\bf gTrNm1} (argument)  $G^{(n+1/2)}$: {\bf gTracer} (argument)
110    
111  $\Delta t$: {\bf deltaTtracer} (PARAMS.h)  $\Delta t$: {\bf deltaTtracer} (PARAMS.h)
112    
113  \end{minipage} }  \end{minipage} }
114    
115    Strictly speaking the ABII scheme should be applied only to the
116    advection terms. However, this scheme is only used in conjuction with
117    the standard second, third and fourth order advection
118    schemes. Selection of any other advection scheme disables
119    Adams-Bashforth for tracers so that explicit diffusion and forcing use
120    the forward method.
121    
122    
123    
124    
125    \section{Linear advection schemes}
126    
127  \begin{figure}  \begin{figure}
128  \resizebox{5.5in}{!}{\includegraphics{part2/advect-1d-lo.eps}}  \resizebox{5.5in}{!}{\includegraphics{part2/advect-1d-lo.eps}}
129  \caption{  \caption{
# Line 136  $\mu=c/(1-c)$. Line 161  $\mu=c/(1-c)$.
161  }  }
162  \end{figure}  \end{figure}
163    
 \section{Linear advection schemes}  
   
164  The advection schemes known as centered second order, centered fourth  The advection schemes known as centered second order, centered fourth
165  order, first order upwind and upwind biased third order are known as  order, first order upwind and upwind biased third order are known as
166  linear advection schemes because the coefficient for interpolation of  linear advection schemes because the coefficient for interpolation of
# Line 228  At boundaries, $\delta_{\hat{n}} \tau$ i Line 251  At boundaries, $\delta_{\hat{n}} \tau$ i
251  $\delta_{nn}$ to be evaluated. We are currently examing the accuracy  $\delta_{nn}$ to be evaluated. We are currently examing the accuracy
252  of this boundary condition and the effect on the solution.  of this boundary condition and the effect on the solution.
253    
254    \fbox{ \begin{minipage}{4.75in}
255    {\em S/R GAD\_U3\_ADV\_X} ({\em gad\_u3\_adv\_x.F})
256    
257    $F_x$: {\bf uT} (argument)
258    
259    $U$: {\bf uTrans} (argument)
260    
261    $\tau$: {\bf tracer} (argument)
262    
263    {\em S/R GAD\_U3\_ADV\_Y} ({\em gad\_u3\_adv\_y.F})
264    
265    $F_y$: {\bf vT} (argument)
266    
267    $V$: {\bf vTrans} (argument)
268    
269    $\tau$: {\bf tracer} (argument)
270    
271    {\em S/R GAD\_U3\_ADV\_R} ({\em gad\_u3\_adv\_r.F})
272    
273    $F_r$: {\bf wT} (argument)
274    
275    $W$: {\bf rTrans} (argument)
276    
277    $\tau$: {\bf tracer} (argument)
278    
279    \end{minipage} }
280    
281  \subsection{Centered fourth order advection}  \subsection{Centered fourth order advection}
282    
# Line 250  F_r & = & W \overline{\tau - \frac{1}{6} Line 299  F_r & = & W \overline{\tau - \frac{1}{6}
299  As for the third order scheme, the best discretization near boundaries  As for the third order scheme, the best discretization near boundaries
300  is under investigation but currenlty $\delta_i \tau=0$ on a boundary.  is under investigation but currenlty $\delta_i \tau=0$ on a boundary.
301    
302    \fbox{ \begin{minipage}{4.75in}
303    {\em S/R GAD\_C4\_ADV\_X} ({\em gad\_c4\_adv\_x.F})
304    
305    $F_x$: {\bf uT} (argument)
306    
307    $U$: {\bf uTrans} (argument)
308    
309    $\tau$: {\bf tracer} (argument)
310    
311    {\em S/R GAD\_C4\_ADV\_Y} ({\em gad\_c4\_adv\_y.F})
312    
313    $F_y$: {\bf vT} (argument)
314    
315    $V$: {\bf vTrans} (argument)
316    
317    $\tau$: {\bf tracer} (argument)
318    
319    {\em S/R GAD\_C4\_ADV\_R} ({\em gad\_c4\_adv\_r.F})
320    
321    $F_r$: {\bf wT} (argument)
322    
323    $W$: {\bf rTrans} (argument)
324    
325    $\tau$: {\bf tracer} (argument)
326    
327    \end{minipage} }
328    
329    
330  \subsection{First order upwind advection}  \subsection{First order upwind advection}
331    
332  Although the upwind scheme is the underlying scheme for the robust or  Although the upwind scheme is the underlying scheme for the robust or
# Line 315  only provide the Superbee limiter \cite{ Line 392  only provide the Superbee limiter \cite{
392  \psi(r) = \max[0,\min[1,2r],\min[2,r]]  \psi(r) = \max[0,\min[1,2r],\min[2,r]]
393  \end{equation}  \end{equation}
394    
395    \fbox{ \begin{minipage}{4.75in}
396    {\em S/R GAD\_FLUXLIMIT\_ADV\_X} ({\em gad\_fluxlimit\_adv\_x.F})
397    
398    $F_x$: {\bf uT} (argument)
399    
400    $U$: {\bf uTrans} (argument)
401    
402    $\tau$: {\bf tracer} (argument)
403    
404    {\em S/R GAD\_FLUXLIMIT\_ADV\_Y} ({\em gad\_fluxlimit\_adv\_y.F})
405    
406    $F_y$: {\bf vT} (argument)
407    
408    $V$: {\bf vTrans} (argument)
409    
410    $\tau$: {\bf tracer} (argument)
411    
412    {\em S/R GAD\_FLUXLIMIT\_ADV\_R} ({\em gad\_fluxlimit\_adv\_r.F})
413    
414    $F_r$: {\bf wT} (argument)
415    
416    $W$: {\bf rTrans} (argument)
417    
418    $\tau$: {\bf tracer} (argument)
419    
420    \end{minipage} }
421    
422    
423  \subsection{Third order direct space time}  \subsection{Third order direct space time}
424    
# Line 353  large Courant number, where the linear u Line 457  large Courant number, where the linear u
457  unstable, the scheme is extremely accurate  unstable, the scheme is extremely accurate
458  (Fig.~\ref{fig:advect-1d-hi}) with only minor overshoots.  (Fig.~\ref{fig:advect-1d-hi}) with only minor overshoots.
459    
460    \fbox{ \begin{minipage}{4.75in}
461    {\em S/R GAD\_DST3\_ADV\_X} ({\em gad\_dst3\_adv\_x.F})
462    
463    $F_x$: {\bf uT} (argument)
464    
465    $U$: {\bf uTrans} (argument)
466    
467    $\tau$: {\bf tracer} (argument)
468    
469    {\em S/R GAD\_DST3\_ADV\_Y} ({\em gad\_dst3\_adv\_y.F})
470    
471    $F_y$: {\bf vT} (argument)
472    
473    $V$: {\bf vTrans} (argument)
474    
475    $\tau$: {\bf tracer} (argument)
476    
477    {\em S/R GAD\_DST3\_ADV\_R} ({\em gad\_dst3\_adv\_r.F})
478    
479    $F_r$: {\bf wT} (argument)
480    
481    $W$: {\bf rTrans} (argument)
482    
483    $\tau$: {\bf tracer} (argument)
484    
485    \end{minipage} }
486    
487    
488  \subsection{Third order direct space time with flux limiting}  \subsection{Third order direct space time with flux limiting}
489    
490  The overshoots in the DST3 method can be controlled with a flux limiter.  The overshoots in the DST3 method can be controlled with a flux limiter.
# Line 373  and the limiter is the Sweby limiter: Line 505  and the limiter is the Sweby limiter:
505  \psi(r) = \max[0, \min[\min(1,d_0+d_1r],\frac{1-c}{c}r ]]  \psi(r) = \max[0, \min[\min(1,d_0+d_1r],\frac{1-c}{c}r ]]
506  \end{equation}  \end{equation}
507    
508    \fbox{ \begin{minipage}{4.75in}
509    {\em S/R GAD\_DST3FL\_ADV\_X} ({\em gad\_dst3\_adv\_x.F})
510    
511    $F_x$: {\bf uT} (argument)
512    
513    $U$: {\bf uTrans} (argument)
514    
515    $\tau$: {\bf tracer} (argument)
516    
517    {\em S/R GAD\_DST3FL\_ADV\_Y} ({\em gad\_dst3\_adv\_y.F})
518    
519    $F_y$: {\bf vT} (argument)
520    
521    $V$: {\bf vTrans} (argument)
522    
523    $\tau$: {\bf tracer} (argument)
524    
525    {\em S/R GAD\_DST3FL\_ADV\_R} ({\em gad\_dst3\_adv\_r.F})
526    
527    $F_r$: {\bf wT} (argument)
528    
529    $W$: {\bf rTrans} (argument)
530    
531    $\tau$: {\bf tracer} (argument)
532    
533    \end{minipage} }
534    
535    
536  \subsection{Multi-dimensional advection}  \subsection{Multi-dimensional advection}
537    
538    \begin{figure}
539    \resizebox{5.5in}{!}{\includegraphics{part2/advect-2d-lo-diag.eps}}
540    \caption{
541    Comparison of advection schemes in two dimensions; diagonal advection
542    of a resolved Guassian feature. Courant number is 0.01 with
543    30$\times$30 points and solutions are shown for T=1/2. White lines
544    indicate zero crossing (ie. the presence of false minima).  The left
545    column shows the second order schemes; top) centered second order with
546    Adams-Bashforth, middle) Lax-Wendroff and bottom) Superbee flux
547    limited. The middle column shows the third order schemes; top) upwind
548    biased third order with Adams-Bashforth, middle) third order direct
549    space-time method and bottom) the same with flux limiting. The top
550    right panel shows the centered fourth order scheme with
551    Adams-Bashforth and right middle panel shows a fourth order variant on
552    the DST method. Bottom right panel shows the Superbee flux limiter
553    (second order) applied independantly in each direction (method of
554    lines).
555    \label{fig:advect-2d-lo-diag}
556    }
557    \end{figure}
558    
559    \begin{figure}
560    \resizebox{5.5in}{!}{\includegraphics{part2/advect-2d-mid-diag.eps}}
561    \caption{
562    Comparison of advection schemes in two dimensions; diagonal advection
563    of a resolved Guassian feature. Courant number is 0.27 with
564    30$\times$30 points and solutions are shown for T=1/2. White lines
565    indicate zero crossing (ie. the presence of false minima).  The left
566    column shows the second order schemes; top) centered second order with
567    Adams-Bashforth, middle) Lax-Wendroff and bottom) Superbee flux
568    limited. The middle column shows the third order schemes; top) upwind
569    biased third order with Adams-Bashforth, middle) third order direct
570    space-time method and bottom) the same with flux limiting. The top
571    right panel shows the centered fourth order scheme with
572    Adams-Bashforth and right middle panel shows a fourth order variant on
573    the DST method. Bottom right panel shows the Superbee flux limiter
574    (second order) applied independantly in each direction (method of
575    lines).
576    \label{fig:advect-2d-mid-diag}
577    }
578    \end{figure}
579    
580    \begin{figure}
581    \resizebox{5.5in}{!}{\includegraphics{part2/advect-2d-hi-diag.eps}}
582    \caption{
583    Comparison of advection schemes in two dimensions; diagonal advection
584    of a resolved Guassian feature. Courant number is 0.47 with
585    30$\times$30 points and solutions are shown for T=1/2. White lines
586    indicate zero crossings and initial maximum values (ie. the presence
587    of false extrema).  The left column shows the second order schemes;
588    top) centered second order with Adams-Bashforth, middle) Lax-Wendroff
589    and bottom) Superbee flux limited. The middle column shows the third
590    order schemes; top) upwind biased third order with Adams-Bashforth,
591    middle) third order direct space-time method and bottom) the same with
592    flux limiting. The top right panel shows the centered fourth order
593    scheme with Adams-Bashforth and right middle panel shows a fourth
594    order variant on the DST method. Bottom right panel shows the Superbee
595    flux limiter (second order) applied independantly in each direction
596    (method of lines).
597    \label{fig:advect-2d-hi-diag}
598    }
599    \end{figure}
600    
601    
602    
603  In many of the aforementioned advection schemes the behaviour in  In many of the aforementioned advection schemes the behaviour in
604  multiple dimensions is not necessarily as good as the one dimensional  multiple dimensions is not necessarily as good as the one dimensional
605  behaviour. For instance, a shape preserving monotonic scheme in one  behaviour. For instance, a shape preserving monotonic scheme in one
# Line 408  So that the over all time-stepping looks Line 633  So that the over all time-stepping looks
633  \begin{equation}  \begin{equation}
634  \tau^{n+1} = \tau^{n} + \Delta t \left( G^{n+1/2}_{adv} + G_{diff}(\tau^{n}) + G^{n}_{forcing} \right)  \tau^{n+1} = \tau^{n} + \Delta t \left( G^{n+1/2}_{adv} + G_{diff}(\tau^{n}) + G^{n}_{forcing} \right)
635  \end{equation}  \end{equation}
636    
637    \fbox{ \begin{minipage}{4.75in}
638    {\em S/R GAD\_ADVECTION} ({\em gad\_advection.F})
639    
640    $\tau$: {\bf Tracer} (argument)
641    
642    $G^{n+1/2}_{adv}$: {\bf Gtracer} (argument)
643    
644    $F_x, F_y, F_r$: {\bf af} (local)
645    
646    $U$: {\bf uTrans} (local)
647    
648    $V$: {\bf vTrans} (local)
649    
650    $W$: {\bf rTrans} (local)
651    
652    \end{minipage} }
653    
654    
655    \section{Comparison of advection schemes}
656    
657    Figs.~\ref{fig:advect-2d-lo-diag}, \ref{fig:advect-2d-mid-diag} and
658    \ref{fig:advect-2d-hi-diag} show solutions to a simple diagonal
659    advection problem using a selection of schemes for low, moderate and
660    high Courant numbers, respectively. The top row shows the linear
661    schemes, integrated with the Adams-Bashforth method. Theses schemes
662    are clearly unstable for the high Courant number and weakly unstable
663    for the moderate Courant number. The presence of false extrema is very
664    apparent for all Courant numbers. The middle row shows solutions
665    obtained with the unlimited but multi-dimensional schemes. These
666    solutions also exhibit false extrema though the pattern now shows
667    symmetry due to the multi-dimensional scheme. Also, the schemes are
668    stable at high Courant number where the linear schemes weren't. The
669    bottom row (left and middle) shows the limited schemes and most
670    obvious is the absence of false extrema. The accuracy and stability of
671    the unlimited non-linear schemes is retained at high Courant number
672    but at low Courant number the tendancy is to loose amplitude in sharp
673    peaks due to diffusion. The one dimensional tests shown in
674    Figs.~\ref{fig:advect-1d-lo} and \ref{fig:advect-1d-hi} showed this
675    phenomenum.
676    
677    Finally, the bottom left and right panels use the same advection
678    scheme but the right does not use the mutli-dimensional method. At low
679    Courant number this appears to not matter but for moderate Courant
680    number severe distortion of the feature is apparent. Moreoever, the
681    stability of the multi-dimensional scheme is determined by the maximum
682    Courant number applied of each dimension while the stability of the
683    method of lines is determined by the sum. Hence, in the high Courant
684    number plot, the scheme is unstable.
685    
686    With many advection schemes implemented in the code two questions
687    arise: ``Which scheme is best?'' and ``Why don't you just offer the
688    best advection scheme?''. Unfortunately, no one advection scheme is
689    ``the best'' for all particular applications and for new applications
690    it is often a matter of trial to determine which is most
691    suitable. Here are some guidelines but these are not the rule;
692    \begin{itemize}
693    \item If you have a coarsely resolved model, using a
694    positive or upwind biased scheme will introduce significant diffusion
695    to the solution and using a centered higher order scheme will
696    introduce more noise. In this case, simplest may be best.
697    \item If you have a high resolution model, using a higher order
698    scheme will give a more accurate solution but scale-selective
699    diffusion might need to be employed. The flux limited methods
700    offer similar accuracy in this regime.
701    \item If your solution has shocks or propagatin fronts then a
702    flux limited scheme is almost essential.
703    \item If your time-step is limited by advection, the multi-dimensional
704    non-linear schemes have the most stablility (upto Courant number 1).
705    \item If you need to know how much diffusion/dissipation has occured you
706    will have a lot of trouble figuring it out with a non-linear method.
707    \item The presence of false extrema is unphysical and this alone is the
708    strongest argument for using a positive scheme.
709    \end{itemize}

Legend:
Removed from v.1.3  
changed lines
  Added in v.1.7

  ViewVC Help
Powered by ViewVC 1.1.22