Lesson 1: Explicit Heat-Equation Schemes

1 Modelling context

The heat equation is the canonical parabolic PDE. It describes diffusion of heat, of concentration, and of any other quantity whose flux is proportional to the negative of its gradient. It is the simplest non-trivial PDE with a time-dependent solution, and it serves as the test problem for numerical schemes.

This lesson develops the forward-time central-space (FTCS) scheme for the one-dimensional heat equation. The scheme is explicit: the value at each grid point at the new time level is computed directly from the values at the previous time level, without solving a system of equations. The scheme is easy to implement, and its behaviour is easy to analyse.

The scheme is not unconditionally stable. The stability restriction on the time step is the subject of Lesson 2, but even at this stage it is worth noting: the FTCS scheme is only usable when the diffusion number r is less than one-half.

The simplest useful PDE solver

Understanding the FTCS scheme thoroughly — its derivation, its implementation, its stability, and its limitations — provides the foundation for every other PDE method you will encounter. Higher-order and implicit schemes are extensions, not replacements.

2 Learning objectives

By the end of this lesson you will be able to

Derive the FTCS scheme for the one-dimensional heat equation.
Implement the scheme in a spreadsheet with documented inputs.
Interpret the diffusion number r and its role in the scheme.
Verify the scheme against an analytical solution.
Identify the dominant error source in the numerical solution.

3 Prerequisite review

Heat equation: ∂u/∂t = α ∂²u/∂x², where α is the thermal diffusivity, with units of length² per time.
Forward difference in time: (u^{n+1} − u^n) / Δt, a first-order approximation to ∂u/∂t at time level n.
Central second difference in space: (u_{i+1} − 2 u_i + u_{i−1}) / Δx², a second-order approximation to ∂²u/∂x² at grid point i.
Explicit scheme: A scheme in which the new value at each grid point is expressed entirely in terms of values at previous time levels, and can therefore be computed directly.

4 Mathematical formulation

Derivation

Consider the one-dimensional heat equation on the domain [0, L] × [0, T]:

∂u/∂t = α ∂²u/∂x²

Replace the time derivative with a forward difference and the spatial second derivative with a central second difference:

( u_i^{n+1} − u_i^n ) / Δt = α ( u_{i+1}^n − 2 u_i^n + u_{i−1}^n ) / Δx²

Solving for u_i^{n+1}:

u_i^{n+1} = u_i^n + ( α Δt / Δx² ) ( u_{i+1}^n − 2 u_i^n + u_{i−1}^n )

Define the diffusion number

r = α Δt / Δx²

Then the scheme becomes

u_i^{n+1} = r u_{i−1}^n + ( 1 − 2r ) u_i^n + r u_{i+1}^n

This is the forward-time central-space (FTCS) scheme. The new value at grid point i is a weighted average of the old value at i and the old values at its two neighbours. The weights sum to r + (1 − 2r) + r = 1, which is consistent with the interpretation of the scheme as a weighted average.

Interpretation of the diffusion number

The diffusion number r is dimensionless. It combines the physical diffusivity α, the time step, and the spatial step into a single quantity that governs the behaviour of the scheme.

Small r. The time step is small relative to the spatial step. The scheme is very stable but takes many steps to advance the solution.
r = 1/6. The scheme is exact for certain polynomial initial conditions, though this is a special property, not a general rule.
r = 1/2. The stability boundary. At this value, the scheme is marginally stable.
r > 1/2. The scheme is unstable. The numerical solution develops growing oscillations.

Implementation structure

The scheme is naturally implemented as a two-dimensional array of grid values, with rows indexed by time and columns by space. Starting from the initial condition at n = 0, each row is computed from the previous row using the update rule. At each step, boundary values are imposed on the first and last entries of the row.

Figure 3.6. Time-marching structure of the FTCS scheme. The solution advances one time level at a time. The boundary columns are fixed by the boundary conditions; the interior values at each new row are computed from the row below.

5 Worked example

FTCS scheme on a sine initial condition

Consider the heat equation on [0, 1] × [0, ∞) with

α = 0.1
u(x, 0) = sin(π x)
u(0, t) = u(1, t) = 0

The exact solution is u(x, t) = e^{−α π² t} sin(π x).

Grid and scheme parameters.

Δx = 0.1, so N_x = 10 and the grid has 11 points.
Choose r = 0.4 (below the stability limit of 0.5).
Then Δt = r Δx² / α = 0.4 · 0.01 / 0.1 = 0.04.
Run to t = 0.4, so N_t = 10.

Initial condition.

u_i^0 = sin(π x_i) for i = 0, …, 10

At i = 0 and i = 10, the value is zero, consistent with the boundary conditions.

First time step.

Apply the update rule at the interior points:

u_i^1 = 0.4 u_{i−1}^0 + 0.2 u_i^0 + 0.4 u_{i+1}^0

For example, at i = 5 (x = 0.5, the peak of the sine):

u_5^0 = sin(π/2) = 1.000
u_4^0 = sin(0.4π) ≈ 0.951
u_6^0 = sin(0.6π) ≈ 0.951
u_5^1 = 0.4(0.951) + 0.2(1.000) + 0.4(0.951) = 0.961

The exact solution at t = 0.04 gives u(0.5, 0.04) = e^{−0.1·π²·0.04} ≈ 0.961. The numerical and exact values agree to three significant figures at this point.

Over the full interval.

Running to t = 0.4 and comparing to the exact solution:

Comparison of FTCS solution with exact solution at `x = 0.5`
t	Numerical	Exact	Error
0.04	0.961	0.961	< 10⁻³
0.12	0.891	0.891	< 10⁻³
0.20	0.827	0.827	< 10⁻³
0.40	0.711	0.713	≈ 2 × 10⁻³

Behaviour.

The numerical solution tracks the exact solution closely. The error grows slowly over time, which is characteristic of a first-order scheme in time. Reducing Δt while holding Δx fixed would reduce this temporal error at the cost of more time steps.

6 Verification check

Checks to perform

Confirm boundary values remain zero (or their specified values) at every time step. A common error is to accidentally recompute the boundary using the interior stencil.
Check the stability condition before running. Compute r = α Δt / Δx² and confirm it is below 0.5.
Compare against the analytical solution where one exists. For the sine initial condition, the exact solution is available. Compute the maximum error over the grid at a fixed time.
Verify convergence under grid refinement. Halve Δx and Δt simultaneously and confirm the error decreases. The observed order should be consistent with the scheme's truncation error.
Check for symmetry. The sine initial condition is symmetric about x = 0.5. The numerical solution should preserve this symmetry to within floating-point tolerance.
Test the steady-state limit. For the given problem, the solution decays to zero as t → ∞. Run for a long time and confirm the numerical solution approaches zero.

7 Interpretation

The FTCS scheme has a clear physical interpretation. The weighted average r u_{i−1} + (1 − 2r) u_i + r u_{i+1} smooths the solution locally. This is exactly what the heat equation does: it diffuses sharp features and drives the solution toward a uniform state.

The diffusion number r controls the strength of the smoothing. Larger r means a stronger smoothing per step; smaller r means a weaker smoothing but more steps to reach the same time. The stability limit r ≤ 1/2 ensures that the smoothing never overshoots and produces values outside the range of the previous step.

Physical interpretation of stability

When r ≤ 1/2, the weight 1 − 2r on the central point is non-negative. The new value is a convex combination of the three previous values, so it lies between their minimum and maximum. When r > 1/2, the central weight becomes negative, and the new value can lie outside that range. This is the origin of the stability restriction.

8 Uncertainty and limitations

The FTCS scheme has two sources of truncation error, one temporal and one spatial.

Truncation error of the FTCS scheme
Component	Approximation	Error order
Time derivative	Forward difference	`O(Δt)`
Space derivative	Central second difference	`O(Δx²)`
Overall scheme	Forward-time central-space	`O(Δt) + O(Δx²)`

The scheme is first-order accurate in time and second-order accurate in space. Refining Δx without refining Δt will not improve the overall accuracy beyond the temporal limit. Both must be refined together to obtain a higher-order overall result.

Stability limits the time step

The CFL condition r ≤ 1/2 imposes an upper bound on the time step: Δt ≤ Δx² / (2α). Refining Δx reduces this bound quadratically, so grid refinement dramatically increases the required number of time steps. This is the principal cost of the FTCS scheme, and the motivation for considering implicit methods in more advanced treatments.

9 Common implementation errors

Watch for these

Applying the update rule to boundary points. The interior stencil requires both neighbours. Boundary points must be assigned from the boundary condition, not computed.
Using the new time level's values when updating. The FTCS scheme is explicit: the new value depends on the old values at all three points. Using the newly computed values in the update would convert the scheme to an implicit one and produce incorrect results.
Reversing the indices. The stencil is symmetric, so reversing i−1 and i+1 does not change the result. But reversing the sign in the (1 − 2r) term does.
Forgetting to compute the diffusion number. The stability check should be performed before every run. A change in α, Δx, or Δt changes r.
Using the wrong formula for r. The correct form is r = α Δt / Δx². Using Δx instead of Δx² produces a scheme that is dimensionally inconsistent and numerically incorrect.
Imposing a boundary condition only at the initial time. Boundary values must be imposed at every time step, not just at n = 0.
Ignoring the shape of the boundary column in the spreadsheet. A common layout error is to put boundary values in a separate row or column that is not updated at each time step.

10 Independent modelling task

Your turn

Choose a one-dimensional heat-equation problem with a known analytical solution. The sine initial condition used in the worked example is a good default.
Choose a spatial grid Δx and a diffusion number r = 0.4. Compute the corresponding time step Δt.
Implement the FTCS scheme in a spreadsheet, with inputs, grid parameters, and computed values in clearly separated regions.
Run the simulation to a final time T. Compute the maximum absolute error at the final time by comparing against the analytical solution.
Repeat the simulation at Δx/2 and Δx/4, adjusting Δt to keep r fixed. Tabulate the errors and compute the observed order of convergence.
Tabulate the number of grid points in space, the number of time steps, and the total number of operations for each grid.
Write a short paragraph describing the tradeoff between accuracy and computational cost for this problem.

This exercise forms the basis of Laboratory 6 and the PDE simulation component of Assignment 3.

11 Reflection prompt

Reflect before continuing

The FTCS scheme represents the heat equation as a weighted average of neighbouring values. Is this always a valid interpretation of the scheme?

Two or three sentences. Consider what happens when r > 1/2, and what the negative weight on the central point means.

12 Lesson summary

Key takeaways

The FTCS scheme combines a forward time difference with a central space difference. It is the simplest explicit scheme for the heat equation.
The update rule is u_i^{n+1} = r u_{i−1}^n + (1 − 2r) u_i^n + r u_{i+1}^n. The diffusion number r = α Δt / Δx² governs the scheme's behaviour.
The scheme is first-order in time and second-order in space. Both Δt and Δx must be refined for the overall error to decrease.
The stability limit is r ≤ 1/2. Above this limit, the numerical solution becomes unstable and grows without bound.
Boundary values are imposed at every time step. They are not computed by the interior stencil.

Coming next

Lesson 2 develops the stability analysis of the FTCS scheme in detail, deriving the CFL condition from the von Neumann method and demonstrating the consequences of violating it.

§ Self-quiz, practice, and homework

Practice exercises

For the sine initial condition on [0, 1] with α = 0.1, compute one FTCS step by hand at i = 5 with Δx = 0.1, r = 0.4.
Verify the exact solution at (0.5, Δt) for the same problem and compare to the numerical value.
For a Gaussian initial condition u(x, 0) = exp(−100(x − 0.5)²), compute the initial FTCS update at i = 5 and i = 4.
Explain in two sentences why the FTCS scheme is described as "smoothing".

Homework problems

FTCS implementation. Implement the FTCS scheme for ∂u/∂t = α ∂²u/∂x² on [0, 1] × [0, 0.4] with u(x, 0) = sin(πx), u(0, t) = u(1, t) = 0, α = 0.1. Run at r = 0.4. Compare the numerical and exact solutions at x = 0.5 for t = 0.1, 0.2, 0.3, 0.4. Produce a table and a plot.
Convergence. Refine the grid at Δx/2 and Δx/4 (holding r fixed). Tabulate the errors and compute the observed order.
Cost scaling. Compute the total number of grid updates at each resolution. Show that the cost grows approximately as 1/Δx³.

Mini project

Implement the FTCS scheme for a heat-equation problem of your choice, with a nontrivial initial condition (e.g. a triangular profile or two overlapping Gaussians). Verify against the analytical solution if one is available, or against a very fine numerical solution otherwise. Produce a technical report with the model, the grid, the scheme, the verification evidence, and the interpretation. This is the basis for Laboratory 6 and Assignment 3.