Laplace Equation 2D Relaxation Method
In two-dimensional Cartesian coordinates, the Laplace equation \( \nabla^2 V = 0 \) is
$$ \frac{\partial^2 V}{\partial x^2} + \frac{\partial^2 V}{\partial y^2} = 0 $$The relaxation method gives a way to approximate the solution to this boundary value problem on a discrete grid. The starting point is the finite-difference approximation to the second derivative:
$$ \frac{\partial^2 V}{\partial x^2} = \frac{V(x - h) - 2V(x) + V(x + h)}{h^2} + O(h^2) $$Applying the same idea in both coordinate directions gives the following approximation to the two-dimensional Laplacian:
$$ \nabla^2 V = \frac{ V(x - h, y) + V(x + h, y) + V(x, y - h) + V(x, y + h) - 4V(x, y) }{h^2} + O(h^2) $$Since Laplace's equation says \( \nabla^2 V = 0 \), the value of \( V \) at an interior grid point should be approximately equal to the average of its four nearest neighbours:
$$ V(x,y) = \frac{1}{4} \left[ V(x - h, y) + V(x + h, y) + V(x, y - h) + V(x, y + h) \right] $$Thus harmonic functions are locally determined by neighbouring values.
This suggests an iterative algorithm. We begin with an initial guess \( V^{(0)}(x,y) \), chosen so that it already satisfies the boundary conditions. The interior values can initially be chosen arbitrarily, often as zero. Then, while keeping the boundary values fixed, we repeatedly update each interior point by replacing it with the average of its four neighbouring points:
$$ V^{(n+1)}(x,y) = \frac{1}{4} \left[ V^{(n)}(x - h, y) + V^{(n)}(x + h, y) + V^{(n)}(x, y - h) + V^{(n)}(x, y + h) \right] $$As the iteration continues, the grid values converge toward a discrete approximation of the solution to Laplace's equation. In this sense, the relaxation method works by gradually smoothing out the interior values until each interior point is consistent with the average of its neighbours.
In the notebook below, we implement this iteration on the rectangle \( [0,2]\times[0,1] \). The vertical sides are held at \( V=1 \), while the horizontal sides are held at \( V = -1 \):
$$ V(0,y) = V(2,y) = 1 $$ $$ V(x,0) = V(x,1) = -1 $$These fixed boundary values act as sources that determine the harmonic potential inside the rectangle.
Starting from an initial guess that satisfies these boundary values, the code repeatedly replaces each interior grid point by the average of its four nearest neighbours. After many iterations, the grid approaches a discrete approximation to the harmonic function satisfying the boundary conditions. We then visualize the result as a colour map and also examine one-dimensional slices through the middle of the rectangle.
The colour map below shows the relaxed solution on the full rectangular domain. The influence of the positive vertical boundaries and negative horizontal boundaries is visible in the way the potential smoothly interpolates across the interior.
To examine the solution more closely, we can take one-dimensional slices through the centre of the rectangle. The first slice fixes \( y = 0.5 \), showing how the solution varies horizontally from one vertical boundary to the other.
The second slice fixes \( x = 1 \), showing how the solution varies vertically between the two horizontal boundaries. This makes the transition between the lower and upper boundary values especially clear.
To see the code for these plots, click the link below:
Open Notebook Environment