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

enter image description here

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). $$

  1. 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!

  • 1
    $\begingroup$ I'd start by looking for a coordinate transformation (using the eigenvalues and eigenvectors of $D$) which transforms your anisotropic diffusion into isotropic diffusion. Then you can transform the classical heat kernel solution back into the original anisotropic space. But I haven't worked out the details. $\endgroup$
    – user3883
    Jul 20, 2018 at 10:46
  • $\begingroup$ This may be a stupid comment, but have you thought about just manufacturing an analytical solution using the method of manufactured solutions? This is generally what I do for anisotropic equations... $\endgroup$
    – BlaB
    Aug 22, 2018 at 14:07

1 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}}$$

  • $\begingroup$ There is one step missing here in getting kernel in $x$ space. The solution involves a convolution so you have to transform the $d$-dimensional measure $d \tilde{x}$ to $dx$. In $x$-space the kernel becomes $\exp(-x^\top D^{-1}x/(4t))/(4\pi t |D|)^{d/2}$. $\endgroup$
    – cfdlab
    Aug 22, 2018 at 11:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.