previous up next SYLLABUS  Previous: 2.5 Implicit Crank-Nicholson  Up: 2.5 Implicit Crank-Nicholson  Next: 2.5.2 Schrödinger equation

2.5.1 Advection-diffusion equation

With Crank-Nicholson, the centered differences in space are evaluated in the same manner and with equal weights for the present and the future values of the function:

Figure 2.5#fig.5: Implicit Crank-Nicholson.

$\displaystyle \frac{f_j^{t+\Delta t} -f_j^t}{\Delta t}$ $\displaystyle +$ $\displaystyle \frac{u}{2}\left\{
\frac{f_{j+1}^{t+\Delta t}-f_{j-1}^{t+\Delta t}}{2\Delta x}
+\frac{f_{j+1}^{t }-f_{j-1}^{t }}{2\Delta x}
  $\displaystyle -$ $\displaystyle \frac{D}{2}\left\{
\frac{f_{j+1}^{t+\Delta t}-2f_j^{t+\Delta t}+f...
...ta t}}{\Delta x^2}
+\frac{f_{j+1}^t -2f_j^t +f_{j-1}^t}{\Delta x^2}
\right\} =0$  

The scheme can conveniently be recast into a linear system with $ j=1,N$ equations

$\displaystyle \left(\begin{array}{c} -\alpha/2 -\beta/4 \\ 1+\alpha \\ -\alpha/...
...\cdot \left(\begin{array}{c} f_{j-1}^t \\ f_j^t \\ f_{j+1}^t \end{array}\right)$ (2.5.1#eq.11)

where the band matrix on the left has three diagonals and the right hand side can be evaluated explicitly from already known quantities. Performing a von Neumann stability analysis for a pure diffusive process ( $ u=\alpha=0$ ) yields the amplification factor

$\displaystyle G=\frac{1-2\alpha\sin^2\left(\frac{k\Delta x}{2}\right)} {1+2\alpha\sin^2\left(\frac{k\Delta x}{2}\right)}$ (2.5.1#eq.12)

This shows that the Crank-Nicholson scheme is unconditionally stable ( $ \forall \Delta x$ , $ \forall\Delta t$ ) with phase errors affecting mainly short wavelengths $ k\Delta x\sim 1$ . A similar conclusion can also be reached for the advective part. The Crank-Nicholson scheme (2.5.1#eq.12) has been implemented in JBONE with
      BandMatrix a = new BandMatrix(3, f.length);
      BandMatrix b = new BandMatrix(3, f.length);
      double[]   c = new double[f.length];
      for (int j=0; j<=n; j++) {
        a.setL(j,-0.25*beta -0.5*alpha);             //Matrix elements
        a.setD(j,            1. +alpha);
        a.setR(j, 0.25*beta -0.5*alpha);
        b.setL(j, 0.25*beta +0.5*alpha);             //Right hand side
        b.setD(j,            1. -alpha);
        b.setR(j,-0.25*beta +0.5*alpha);
      };                                    //Right hand side
      c[0]=c[0]+b.getL(0)*f[n];                      // with periodicity

      fp=a.solve3(c);                                //Solve linear problem
The BandMatrix.solve3() method solves the linear system efficiently in $ \mathcal{O}(N)$ operations, using an LU-factorization that will be discussed later in sect.3.

The favorable stability property (2.5.1#eq.13) can nicely be exploited in diffusion dominated problems dealing with the evolution of large scale features $ \lambda\gg\Delta x$ . Starting from a relatively smooth Gaussian pulse that is subject to both advection and diffusion $ u=D=1$ , the applet below shows that a reasonably accurate solution (12 % for the valley to peak ratio when the time reaches 100) can be computed using extremely large time steps $ \Delta t=5$ , where $ \alpha=\beta=5$ .

JBONE applet:  press Start/Stop to simulate the advection of a box function and to test how the implicit 2 levels Crank-Nicholson scheme affects the dispersion and the damping of superposed short and long wavelengths.

Numerical experiments: Crank-Nicholson for advection-diffusion
  1. Switch to an initial Box function; explain what happens after one time step.
  2. Complete the Implicit 2-level calculation with a Box and terminate with a few tiny TimeStep=0.1 using the Explicit 2-level scheme (2.1#eq.1). Compare the result with the best solution you have obtained so far.

SYLLABUS  Previous: 2.5 Implicit Crank-Nicholson  Up: 2.5 Implicit Crank-Nicholson  Next: 2.5.2 Schrödinger equation

back up next contents bibliography Copyright © Lifelong-learners at 03:13:01, October 21st, 2017