Solving the Diffusion Equation Explicitly

This post is part of a series of articles about the Finite Difference Method, which includes other topics such as:

• Approximation of derivatives using the Finite Difference Method
• Explicit method for solving diffusion equations
• Implicit Crank-Nicolson method
• The solver for a tridiagonal matrix equation using Thomas's algorithm

In the first article of this series, it has been shown that the derivatives of continuous functions can be approximated on a discontinuous domain.The next step is to apply these derivatives to parabolic partial differential equations (PDEs), with the heat equation being a standard example, and we will discretize this equation.

Heat equation and discontinuity

The one-dimensional heat equation has the form:

∂u/∂t ​= D(∂2u/∂x2)​

Where u(x,t) is the temperature at position x and time t, D is the diffusivity coefficient.

We will approximate the derivative:

Use Forward Difference for the time derivative:

∂u/∂t ​≈ (uin+1​−uin)​​/Δt

Use the second-order Centered Difference for the spatial derivative:

∂x2∂2u​≈(Δx)2ui+1n​−2uin​+ui−1n​​

Substituting the values into the equation, we get:

Δtuin+1​−uin​​=D⋅(Δx)2ui+1n​−2uin​+ui−1n​​

Rearrange to find uin+1:

uin+1​=uin​+(Δx)2DΔt(ui+1n−2uin+ui−1n)

Let α=(Δx)2DΔt will yield:

uin+1​=uin​+α(ui+1n​−2uin​+ui−1n​)

Computation of Solution

This equation can be used to directly compute uin+1 from the known values of uin, ui−1n, and ui+1n, hence it is called an explicit equation.

If the initial conditions are set, for example:

u0=[0,0,1,0,0,1,0,0,1,0,0]

With "heat wave peaks" at 3 points and selecting the parameters:

Δx=1,Δt=0.001,D=1⇒α=1⋅0.00112=0.001\Delta x = 1,\quadDelta t = 0.001,\quad D = 1 \Rightarrow \alpha = \frac{1 \cdot 0.001}{1^2} = 0.001Δx=1,Δt=0.001,D=1⇒α=121⋅0.001​=0.001

Can use the formula:

Δx=1,Δt=0.001,D=1⇒α=121⋅0.001​=0.001

To calculate the time-marching continuous values

Stability Issue

If you choose Δt=0.01, you will get α=0.01, which causes the answer to oscillate violently:

• The wave pattern obtained will have unnatural fluctuations.
• Not consistent with the actual answer according to the principles of physics, which show that heat should gradually spread.

Stability Criterion

To ensure the stability of this explicit method, the condition must be met:

α=(Δx)2DΔt​≤21​

If this is not followed, the answer will be unstable and may result in distorted waves.

Example of stable calculation

Choose:

Δx=1,Δt=0.001,D=1⇒α=0.001≤0.5

Formula:

uin+1​=uin​+0.001(ui+1n​−2uin​+ui−1n​)

The result will be:

No negative values.

Accurately reflects the characteristics of heat propagation.

In summary

Although this explicit method is easy to compute, the main drawback is the limitation of the time step:

• If Δx is very small (fine mesh), Δt must also be very small accordingly.
• Resulting in the calculation taking a long time.

The next article will introduce implicit methods, such as Crank-Nicolson, to avoid this limitation.

Reference: Solving the Diffusion Equation Explicitly

From https://www.quantstart.com/articles/Solving-the-Diffusion-Equation-Explicitly/