Further to my posts on solving the Laplace Equation and Poisson Equation in $R \subset \mathbb{R}^{2}$ this post is describing mathematical extensions to these methods allowing us to consider numerical solvers for regions $R \subset \mathbb{R}^{n}$, for arbitrary $n \in \mathbb{N}$.

We begin by once again considering the multidimensional Poisson equation in $\mathbb{R}^{n}$, where $\psi,f : \mathbb{R}^{n} \to \mathbb{R}$:

$$\nabla^{2} \psi = \sum_{k=1}^{n}\frac{\partial^{2} \psi}{\partial x_{k}^{2}}=f(\vec{r})$$

Recalling the central difference method from before: $f^{\prime \prime}(x) \approx \frac{f(x+h) – 2f(x) + f(x-h)}{h^{2}}$, we can write:

$$\nabla^{2} \psi \approx \sum_{k=1}^{n}\frac{\psi(\dots,x_{k}+\delta x_{k},\dots) – 2\psi(\dots,x_{k},\dots) + \psi(\dots,x_{k}-\delta x_{k},\dots)}{(\delta x_{k})^{2}}=f(x_{1},\dots,x_{n})$$

Or in a more readable vector notation, where $\vec{r} = \left(\begin{smallmatrix}x_{1} & \dots & x_{n}\end{smallmatrix}\right)^{T}$ and $\delta \vec{x}_{k} = \delta x_{k} \boldsymbol{\hat{e}}_{k}$, and $\boldsymbol{\hat{e}}_{k}$ is the $k$th unit vector:

$$\nabla^{2} \psi \approx \sum_{k=1}^{n}\frac{\psi(\vec{r}+\delta \vec{x}_{k}) – 2\psi(\vec{r}) + \psi(\vec{r} – \delta \vec{x}_{k})}{\|\delta \vec{x}_{k}\|^{2}}=f(\vec{r})$$

Rearranging for $\psi(\vec{r})$, we have:

$$\psi(\vec{r}) = \frac{1}{2\xi} \left(\sum_{k=1}^{n}\frac{\psi(\vec{r}+\delta \vec{x}_{k}) + \psi(\vec{r} – \delta \vec{x}_{k})}{\|\delta \vec{x}_{k}\|^{2}} -f(\vec{r})\right), \qquad \xi \triangleq \sum_{k=1}^{n}\frac{1}{\|\delta \vec{x}_{k}\|^{2}}$$

As before, if we take a uniform grid, we have $\|\delta \vec{x}_{i}\| = \|\delta \vec{x}_{j}\| = \delta$, $\ \forall i,j\in\{1,\dots,n\}$, then this simplifies as before to give the overrelaxation formula:

$$\psi^{m+1}(\vec{r}) = \psi^{m}(\vec{r}) + \frac{\alpha}{2n}\left(\sum_{k=1}^{n}\left(\psi^{m+1}(\vec{r}+\delta \vec{x}_{k}) + \psi^{m}(\vec{r} – \delta \vec{x}_{k})\right) – 2n\psi^{m}(\vec{r}) – \delta^{2} f(\vec{r})\right)$$

So the only real difference between the 2-dimensional case analyzed before and the $n$-dimensional generalization here is the introduction of the $\frac{1}{n}$ as a prefactor for the residual and the $2n$ factor for the $\psi^{m}(\vec{r})$ term in the residual.

As we have seen, successive over-relaxation is a versatile and easily-programmable method for the solution of the multidimensional Poisson equation; however, it should be noted that there are more efficient methods such as multigrid methods. I will address these concerns in a future post, in which I will detail a C++ implementation for the numeric solution of the multidimensional Poisson equation by successive overrelaxation using a parallelized red-black method.