2.44 Generate sparse matrix for the Laplacian differential operator \(\nabla ^{2}u\) for 2D grid

\(\nabla ^{2}u=f\) in 2D is defined as \(\frac {\partial ^{2}u}{\partial x^{2}}+\frac {\partial ^{2}u}{\partial y^{2}}=f\) and the Laplacian operator using second order standard differences results in \(\left ( \frac {u_{i-1,j}+u_{i+1,j}+u_{i,j-1}+u_{i,j+1}-4u_{i,j}}{h^{2}}\right ) =f_{i,j}\) where \(h\) is the grid size. The above is solved for \(x\) in the form \(Ax=f\) by generating the \(A\) matrix and taking into consideration the boundary conditions. The follows show how to generate the sparse representation of \(A\). Assuming the number of unknowns \(n=3\) in one direction, there are 9 unknowns to solve for and the \(A\) matrix will be \(9\) by \(9\).

Matlab internalPoints=3; e = ones(internalPoints,1); spe = spdiags([e -2*e e], -1:1,... internalPoints,internalPoints); Iz = speye(internalPoints); sp = kron(Iz,spe)+kron(spe,Iz); full(sp)

ans = -4 1 0 1 0 0 0 0 0 1 -4 1 0 1 0 0 0 0 0 1 -4 0 0 1 0 0 0 1 0 0 -4 1 0 1 0 0 0 1 0 1 -4 1 0 1 0 0 0 1 0 1 -4 0 0 1 0 0 0 1 0 0 -4 1 0 0 0 0 0 1 0 1 -4 1 0 0 0 0 0 1 0 1 -4