# Discretization matrix for 3D Poisson equation

It is known that the 2D Poisson equation defined on a domain $$\Omega$$ (let's say $$\Omega := (0,1)^2$$) with Dirichlet boundary conditions $$u(x,y)_{|\partial \Omega}=g(x,y)$$, $$u_{xx} + u_{yy}=f$$ can be discretized, using finite differences, to obtain a system of linear equations of the form

$$A\vec{u} = \vec{f}+\vec{b}$$ where $$A$$ is a coefficient matrix, $$\vec{u}$$ the discrete solution to solve for, and $$\vec{b}$$ contains terms from the boundary conditions.

The matrix $$A$$ for the 2D Poisson equation has the block form

$$\begin{bmatrix} D & I & \dots &0 \\ I & \ddots & \ddots & \vdots\\ \vdots & \ddots & \ddots & I\\ 0 & \dots & I & D \end{bmatrix}$$

where

$$D=\begin{bmatrix} -4 & 1 & 1 &0 & \dots& 0 \\ 1 & -4 & 1 & 1 & \dots & 0 \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots\\ \vdots & \ddots & \vdots & \dots & -4 & 1\\ 0 & \dots & \dots & 1 & 1 & -4 \end{bmatrix}$$

I want to discretize the 3D Poisson equation

$$u_{xx} + u_{yy}+u_{zz}=f$$

on $$\Omega := (0,1)^3$$.

So, I suppose, that the discretized set $$\Omega_{\Delta x}$$ must be a 3D mesh. The resulting discretization equations will be

$$\Delta_h u_{i,j,k} = \frac{u_{i+1,j,k}+u_{i,j+1,k}+u_{i,j,k+1}-6u_{i,j,k}+u_{i-1,j,k}+u_{i,j-1,k}+u_{i,j,k-1}}{(\Delta x)^2}=f_{i,j,k}+b_{i,j,k}$$

So that $$A$$ should now look like this:

$$A=\begin{bmatrix} -6&1 & 1 & 1 &0 & 0 & 0&\dots& 0 \\ 1 & 1 & -6 & 1 & 1 & 1 &0&\dots & 0 \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots &\vdots&\ddots&0\\ \vdots & \vdots & \vdots & \vdots & \dots & \dots & \dots&-6 & 1\\ 0 & \dots & 0 & 0 & 0 & 1 & 1 & 1 & -6 \end{bmatrix}$$

I'm not quite sure, however, what the form of $$D$$ exactly is. I.e. is $$A$$ still block-tridiagonal, with two block diagonals comprised of $$I$$'s, and the main diagonal comprised of $$D$$'s, but with higher dimensions?

Would appreciate some insight.

I like a lot the Kronecker product representation. If $$A_1$$ is the one-dimensional second-derivative approximation, i.e. $$A_1 = \begin{pmatrix} 2 & -1 & \\ -1 & 2 & \ddots \\ & \ddots & \ddots \end{pmatrix},$$ then you can write the $$3d$$ Laplacian as $$A_3 = I \otimes I \otimes A_1 + I \otimes A_1 \otimes I + A_1 \otimes I \otimes I.$$ For the $$2d$$ case it is the same, i.e., $$A_2 = I \otimes A_1 + A_1 \otimes I$$, and this generalizes to higher dimensions. You can have a look into a small instance of this matrix, and its properties (e.g., sparsity pattern) will be clear.
• The Kronecker delta is a simple function: $\delta_{ij} = 0$ if $i \neq j$, and $\delta_{ii} = 1$. The Kronecker product is instead a matrix product, whose definition can be found here (en.wikipedia.org/wiki/Kronecker_product) and which is extremely helpful in multivariate and tensor analysis. I strongly suggest you to have a look into it. Commented Feb 21, 2019 at 7:46