Analytical testcase for 2D/3D anisotropic Diffusion (Heat Kernel)

Question

I want to verify and compare different Discretizations of the anisotropic diffusion equation in 2D / 3D.

In order to both test the timestepping and the spatial discretisations I had a look at using the heat kernels as an analytical solution to diffusion equations. The one dimensional heat kernel (wiki-link) looks like this: $$ u(x,t) = \frac{1}{(4\pi t)^{d/2}} e^{-x^2/4t}\ $$ and solves following diffusion equation of dimension $d$ analytically: $$ \frac{\partial u}{\partial t} = -\Delta u $$

The idea is that for an analytical solution I can verify my convergence both spatially and temporal. It is important in my case, that the diffusion is anisotropic, and might not be aligned with the underlying grids.

My anisotropic diffusion equation reads:

$$ \frac{\partial u}{\partial t} = \nabla (D \nabla u) $$

In order to expand the analytical solution to my case, where the diffusion is not scalar and unity, but a full (symmetric positive definite) tensor D, I had a look at multivariate gaußian distributions (wiki-link).

My question is the following:

1. Is there somewhere in the literature a heat kernel for anisotropic (but homogeneous) diffusion which I can use to verify my numerics?

By analogy I would come up with something like:

$$ u(\vec{x},\vec{\mu},t) = \frac{1}{(4\pi~ |\mathbf{D}| t)^{d/2}} \exp \large(\frac{-(\vec{x}-\vec{\mu})^{T} \mathbf{D}^{-1}(\vec{x}-\vec{\mu})}{4~|\mathbf{D}|~t}\large). $$

2. Is the above equation a proper analytical solution of the anisotropic diffusion equation?

I am not sure if the exponent $d/2$ is correct for the $|D|$. I am aware that the question might touch both physics, and numerics, but I thought this might be useful to other people trying to verify their codes.

Any hints or help is appreciated!

Answer 1

Since Laplacian is an elliptic operator you are looking for the Cholesky decomposition of the assumed constant diffusion matrix $D$: $$D = LL^T$$

Therefore the parabolic equation may be written as: $$\partial_t u=-\nabla^T(D\nabla u)=-(L^T\nabla)^T(L^T\nabla u)=-\tilde{\nabla}^T(\tilde{\nabla} u)=-\tilde{\nabla}^2u$$

Naming $\tilde{\nabla}=L^T\nabla$ and therefore $\tilde{x}=L^{-1}x$. This last PDE is totally decoupled (using Fourier transform and separation of variables).

The solution for the unnormalised heat kernel will be (using $\tilde{x}=L^{-1}x$): $$\Phi=\prod_{n=1}^d\frac{\exp{(-\tilde{x}_n^2/(4t))}}{\sqrt{4\pi t}}=\frac{\exp{(-\tilde{x}^T\tilde{x}/(4t))}}{(4\pi t)^{d/2}}=\frac{\exp{(-x^TD^{-1}x/(4t))}}{(4\pi t)^{d/2}}$$