Kronecker sum of discrete Laplacians

From HandWiki

In mathematics, the Kronecker sum of discrete Laplacians, named after Leopold Kronecker, is a discrete version of the separation of variables for the continuous Laplacian in a rectangular cuboid domain.

General form of the Kronecker sum of discrete Laplacians

In a general situation of the separation of variables in the discrete case, the multidimensional discrete Laplacian is a Kronecker sum of 1D discrete Laplacians.

Example: 2D discrete Laplacian on a regular grid with the homogeneous Dirichlet boundary condition

Mathematically, using the Kronecker sum:

[math]\displaystyle{ L = \mathbf{D_{xx}}\oplus\mathbf{D_{yy}}=\mathbf{D_{xx}}\otimes\mathbf{I}+\mathbf{I}\otimes\mathbf{D_{yy}}, \, }[/math]

where [math]\displaystyle{ \mathbf{D_{xx}} }[/math] and [math]\displaystyle{ \mathbf{D_{yy}} }[/math] are 1D discrete Laplacians in the x- and y-directions, correspondingly, and [math]\displaystyle{ \mathbf{I} }[/math] are the identities of appropriate sizes. Both [math]\displaystyle{ \mathbf{D_{xx}} }[/math] and [math]\displaystyle{ \mathbf{D_{yy}} }[/math] must correspond to the case of the homogeneous Dirichlet boundary condition at end points of the x- and y-intervals, in order to generate the 2D discrete Laplacian L corresponding to the homogeneous Dirichlet boundary condition everywhere on the boundary of the rectangular domain.

Here is a sample OCTAVE/MATLAB code to compute L on the regular 10×15 2D grid:

nx = 10; % number of grid points in the x-direction;
ny = 15; % number of grid points in the y-direction;
ex = ones(nx,1);
Dxx = spdiags([ex -2*ex ex], [-1 0 1], nx, nx); %1D discrete Laplacian in the x-direction ;
ey = ones(ny,1);
Dyy = spdiags([ey, -2*ey ey], [-1 0 1], ny, ny); %1D discrete Laplacian in the y-direction ;
L = kron(Dyy, speye(nx)) + kron(speye(ny), Dxx) ;

Eigenvalues and eigenvectors of multidimensional discrete Laplacian on a regular grid

Knowing all eigenvalues and eigenvectors of the factors, all eigenvalues and eigenvectors of the Kronecker product can be explicitly calculated. Based on this, eigenvalues and eigenvectors of the Kronecker sum can also be explicitly calculated.

The eigenvalues and eigenvectors of the standard central difference approximation of the second derivative on an interval for traditional combinations of boundary conditions at the interval end points are well known. Combining these expressions with the formulas of eigenvalues and eigenvectors for the Kronecker sum, one can easily obtain the required answer.

Example: 3D discrete Laplacian on a regular grid with the homogeneous Dirichlet boundary condition

[math]\displaystyle{ L = \mathbf{D_{xx}}\otimes\mathbf{I}\otimes\mathbf{I}+\mathbf{I}\otimes\mathbf{D_{yy}}\otimes\mathbf{I}+\mathbf{I}\otimes\mathbf{I}\otimes\mathbf{D_{zz}}, \, }[/math]

where [math]\displaystyle{ \mathbf{D_{xx}},\,\mathbf{D_{yy}} }[/math] and [math]\displaystyle{ \mathbf{D_{zz}} }[/math] are 1D discrete Laplacians in every of the 3 directions, and [math]\displaystyle{ \mathbf{I} }[/math] are the identities of appropriate sizes. Each 1D discrete Laplacian must correspond to the case of the homogeneous Dirichlet boundary condition, in order to generate the 3D discrete Laplacian L corresponding to the homogeneous Dirichlet boundary condition everywhere on the boundary. The eigenvalues are

[math]\displaystyle{ \lambda_{jx,jy,jz} = -\frac{4}{h_x^2} \sin\left(\frac{\pi j_x}{2(n_x + 1)}\right)^2 -\frac{4}{h_y^2} \sin\left(\frac{\pi j_y}{2(n_y + 1)}\right)^2 -\frac{4}{h_z^2} \sin\left(\frac{\pi j_z}{2(n_z + 1)}\right)^2 }[/math]

where [math]\displaystyle{ j_x=1,\ldots,n_x,\,j_y=1,\ldots,n_y,\,j_z=1,\ldots,n_z,\, }[/math], and the corresponding eigenvectors are

[math]\displaystyle{ v_{ix,iy,iz,jx,jy,jz} = \sqrt{\frac{2}{n_x+1}} \sin\left(\frac{i_x j_x \pi}{n_x+1}\right) \sqrt{\frac{2}{n_y+1}} \sin\left(\frac{i_y j_y \pi}{n_y+1}\right) \sqrt{\frac{2}{n_z+1}} \sin\left(\frac{i_z j_z \pi}{n_z+1}\right) }[/math]

where the multi-index [math]\displaystyle{ {jx,jy,jz} }[/math] pairs the eigenvalues and the eigenvectors, while the multi-index [math]\displaystyle{ {ix,iy,iz} }[/math] determines the location of the value of every eigenvector at the regular grid. The boundary points, where the homogeneous Dirichlet boundary condition are imposed, are just outside the grid.

Available software

An OCTAVE/MATLAB code http://www.mathworks.com/matlabcentral/fileexchange/27279-laplacian-in-1d-2d-or-3d is available under a BSD License, which computes the sparse matrix of the 1, 2D, and 3D negative Laplacians on a rectangular grid for combinations of Dirichlet, Neumann, and Periodic boundary conditions using Kronecker sums of discrete 1D Laplacians. The code also provides the exact eigenvalues and eigenvectors using the explicit formulas given above.