PDEs via Finite Differences

PDEadvanced~12 min

Solve the 1D heat equation step-by-step using an explicit finite-difference scheme implemented with for-loops.

Step 1 — Discretise the domain

Divide the spatial interval [0, 1] into N cells of width dx, and choose a time step dt small enough to satisfy the stability condition dt ≤ dx²/2 (Courant–Friedrichs–Lewy criterion).

N = 50;
dx = 1 / N;
dt = 0.4 * dx^2;
x = linspace(0, 1, N+1);
printf('dx=%.3f, dt=%.6f, CFL ratio=%.2f\n', dx, dt, dt/dx^2)
▶ Run in SimLab

Expected output: dx=0.020, dt=0.000160, CFL ratio=0.40

Step 2 — Set initial condition and boundary conditions

Start with a Gaussian heat pulse centred at x=0.5. Dirichlet boundary conditions hold u(0) = u(1) = 0 (fixed cold ends).

N = 50;
dx = 1 / N;
x = linspace(0, 1, N+1);
u = exp(-200*(x - 0.5).^2);
u(1) = 0; u(end) = 0;
plot(x, u);
xlabel('x'); ylabel('u');
title('Initial heat pulse')
▶ Run in SimLab

Expected output: Sharp Gaussian pulse centred at x=0.5

Step 3 — Time-step the heat equation

Apply the explicit update rule u_new(i) = u(i) + dt/dx² * (u(i+1) - 2u(i) + u(i-1)) for 2000 steps and plot the solution at several snapshots to watch the pulse diffuse.

N = 50; dx = 1/N; dt = 0.4*dx^2;
x = linspace(0, 1, N+1);
u = exp(-200*(x - 0.5).^2);
u(1) = 0; u(end) = 0;
hold on;
nSteps = 2000;
r = dt / dx^2;
for step = 1:nSteps
  u_new = u;
  for i = 2:N
    u_new(i) = u(i) + r*(u(i+1) - 2*u(i) + u(i-1));
  end
  u = u_new;
  if mod(step, 500) == 0
    plot(x, u);
  end
end
xlabel('x'); ylabel('u');
title('Heat Diffusion at t=0, 0.08, 0.16, 0.24, 0.32')
▶ Run in SimLab

Expected output: Four curves showing the pulse spreading and flattening over time

Related Tutorials

Try SimLab — Free MATLAB® Alternative

466 functions. Runs in your browser. No install.

Open SimLab

Stay Updated

Get notified about new simulations and tools. We send 1-2 emails per month.