The Explicit Method for the 1-D Diffusion Equation

We are considering the model problem for u(x, t):

\begin{matrix}u_{t} & = u_{xx}  in 0 < x < 1, t > 0 \\ with u(0, t) & = 0 , u(1, t) = 1 , u(x, 0) = f(x) \end{matrix} \} .  (3.1)

We define a regular, rectangular grid (x_{n}, t^{j}) for 0 ⩽ n ⩽ N, 0 ⩽ j ⩽ J of lengths (h, k), so that Nh = 1, x_{n} = nh, t^{j} = jk. We shall seek an approximation U_{n}^{j} to the exact solution evaluated on the grid points, u_{n}^{j} ≡ u(nh, jk). The boundary conditions require u_{0}^{j} = 0, u_{N}^{j} = 1 and u_{n}^{0} = f_{n} ≡ f(x_{n}). Recalling the results from §1, we have

\frac{u_{n}^{j+1} − u_{n}^{j}}{k} = \frac{∂u_{n}^{j}}{∂t} + \frac{1}{2}k \frac{∂^{2}u_{n}^{j}}{∂t^{2}} + O(k^{2}) ,

and

\frac{u_{n+1}^{j} − 2u_{n}^{j} + u_{n−1}^{j}}{h^{2}} = \frac{∂^{2}u_{n}^{j}}{∂x^{2}} + \frac{1}{12}h^{2} \frac{∂^{4}u_{n}^{j}}{∂x^{4}} + O(h^{4}) .

Using the above in the equation u_{t} = u_{xx}, we have

\frac{u_{n}^{j+1} − u_{n}^{j}}{k} = \frac{u_{n+1}^{j} − 2u_{n}^{j} + u_{n−1}^{j}}{h^{2}} + R_{n}^{j} ,  (3.2)

where the Truncation Error, R_{n}^{j} is

R_{n}^{j} = \frac{1}{12}h^{2} \frac{∂^{4}u_{n}^{j}}{∂x^{4}} − \frac{1}{2}k \frac{∂^{2}u_{n}^{j}}{∂t^{2}} + O(k^{2}, h^{4}) .  (3.3)

The simplest explicit method neglects terms of O(k, h^{2}). If we write

r = k/h^{2} ,

neglect R_{n}^{j} and replace u by U in (3.2), we obtain the scheme

U_{n}^{j+1} = rU_{n+1}^{j} + (1 − 2r)U_{n}^{j} + rU_{n−1}^{j} ,  (3.4)

with the boundary conditions

U_{0}^{j} = 0 ,  U_{n}^{j} = 1  and U_{n}^{0} = f_{n} ≡ f(x_{n}) .

With a little thought, we see that these boundary conditions and repeated use of (3.4) enable us to calculate U_{n}^{j} everywhere. We must now consider how accurate the approximation is. We shall perform a simple “Maximum Principle Analysis” to show that under some conditions U_{n}^{j} can be made as close as we choose to the real solution u_{n}^{j}, for all n and j.

Truncation Error and Solution Error

In general, suppose a PDE Lu = f, where L is some differential operator, is approximated on a suitable grid (nh, jk), at which points u and f take the values u_{n}^{j} and f_{n}^{j}, by the FDM MU_{n}^{j} = f_{n}^{j}, where M is a difference operator. The solution error, z_{n}^{j} and the truncation error, R_{n}^{j}, are defined by

z_{n}^{j} ≡ u_{n}^{j} − U_{n}^{j} ,  and R_{n}^{j} ≡ (Lu)_{n}^{j} − Mu_{n}^{j} .  (3.5)

Returning to our specific equation, subtracting (3.4) from (3.2), and using (3.5), we obtain

z_{n}^{j+1} = rz_{n+1}^{j} + (1 − 2r)z_{n}^{j} + rz_{n−1}^{j} + kR_{n}^{j}  and z_{0}^{j} = z_{N}^{j} = z_{n}^{0} = 0 .  (3.6)

Now

|z_{n}^{j+1}| ⩽ |r| |z_{n+1}^{j}| + |(1 − 2r)| |z_{n}^{j}| + |r| |z_{n−1}^{j}| + |kR_{n}^{j}|.

Since R_{n}^{j} = O(k, h^{2}), over the finite interval 0 < t < T ≡ Jk we can find a positive constant A such that |R_{n}^{j}| ⩽ A(|k| + h^{2}). We shall also define the norm

||z^{j}|| ≡ max_{n=0...N} |z_{n}^{j}| ,

which is the maximum error over all the points at a fixed time-level j. Then |z_{n}^{j}| ⩽ ||z^{j}|| for all n, and so we have

|z_{n}^{j+1}| ⩽ (|r| + |1 − 2r| + |r|)||z^{j}|| + A(k^{2} + |k|h^{2})

As this is true for all values of n, it is true for that value which maximises its LHS, and so

||z^{j+1}|| ⩽ (|r| + |1 − 2r| + |r|)||z^{j}|| + A(k^{2} + |k|h^{2})

We now assume that 0 < r ⩽ \frac{1}{2}, so that the quantities inside modulus signs are positive. Then the maximum possible error at the time-level (j + 1) is related to the maximum at time j by

||z^{j+1}|| ⩽ ||z^{j}|| + A(k^{2} + |k|h^{2})

Now the initial error, ||z^{0}||, is zero because we know the solution exactly at t = 0. Applying the above repeatedly therefore implies

||z^{1}|| ⩽ |k|A(|k| + h^{2}) ,  ||z^{2}|| ⩽ 2|k|A(|k| + h^{2})  and ||z^{j}|| ⩽ j|k|A(|k| + h^{2}) .

We have therefore shown that provided 0 < r ⩽ \frac{1}{2}, the maximum possible error at time t^{j} ≡ jk is O(k, h^{2})t^{j}, which can be made as small as we choose by choosing small enough steplengths, k and h. Note that r < 0 would correspond to k < 0, which involves stepping backwards in time, which we saw in §2 leads to an ill-posed problem for this Parabolic equation. We have proved nothing yet about the case r > \frac{1}{2}, but we will find that for such values of r the FDM (3.2) is numerically unstable.