Solutions to the heat equation converge to the similarity solution

The heat equation $$ \begin{align} u_t &= u_{xx}, \quad x \in \mathbb{R}, \quad t > 0, \\ u(x, 0) &= u_0(x) \end{align} $$ admits a self-similar solution of the form $$ K(x, t) = \frac{1}{(4 \pi t)^{1/2}} \exp\left(\frac{-x^2}{4t}\right) $$ when the initial condition $u_0(x) = \delta(x)$ is a function localised at $x = 0$ with "magnitude" $\int_{-\infty}^\infty u_0(x)\,dx = 1$. The form of the solution is motivated by scaling/dimensional arguments, but why do we pay so much attention to this specific solution? Especially since the delta function initial condition seems unphysical?

The similarity solution is important for two main reasons. The first is that it is the "Green's function" or "kernel" for the diffusion operator; the linearity of the heat equation allows us to write down the general solution for arbitrary $u_0$ as the convolution of $u_0$ with $K$: $$ u(x, t) = \int_{-\infty}^\infty u_0(\xi) K(x-\xi, t)\, d\xi. $$ This property is a powerful technique for solving linear equations, although it does not generalise to nonlinear diffusion equations.

The second, more physically relevant reason is that any such solution asymptotically approaches the similarity solution, in the following sense. The convolution theorem for Fourier transforms implies that $$ u(x, t) = \frac{1}{2\pi} \int_{-\infty}^\infty \hat{u_0}(k) \exp\left( ikx - k^2 t \right) \,dk $$ where $\hat{u_0}(k) = \int_{-\infty}^\infty u_0(x) \exp(-ikx) \,dx$ is the Fourier transform of $u_0$. A little algebra shows that $$ u(x, t) = \frac{1}{2\pi}\, \exp\left(\frac{-x^2}{4t}\right) \int_{-\infty}^\infty \hat{u_0}(k) \exp\left[ -t \left(k - \frac{ix}{2t}\right)^2\right] \,dk $$

Take the limit $t\rightarrow\infty$ and $|x|\rightarrow\infty$ with $\eta \propto t^{-1/2} x$ fixed, so that the factor outside the integral is fixed. Make the substitution¹ $$ k = \frac{ix}{2t} + \frac{\lambda}{t^{1/2}}, \quad dk = t^{-1/2} \,d\lambda. $$ We then have $$ u(x, t) = \frac{1}{2\pi t^{1/2}}\, \exp\left(\frac{-x^2}{4t}\right) \int_{-\infty}^\infty \hat{u_0}\left(\frac{\lambda}{t^{1/2}} + \frac{ix}{2t}\right) \exp(-\lambda^2) \,d\lambda. $$ In the limit $t \rightarrow \infty$, this is approximately²

$$ \begin{align} u(x, t) &\sim \frac{1}{2\pi t^{1/2}}\, \exp\left(\frac{-x^2}{4t}\right) \int_{-\infty}^\infty \hat{u_0}(0) \exp(-\lambda^2) \,d\lambda \\ &= \frac{\hat{u_0}(0)}{(4\pi t)^{1/2}} \, \exp\left(\frac{-x^2}{4t}\right) \\ &= \hat{u_0}(0) K(x, t) \end{align} $$

as required, where $\hat{u_0}(0) = \int_{-\infty}^\infty u(x) \, dx$ is the "magnitude" of the initial condition – scaled appropriately.

Thus, to leading order, the solution $u(x, t)$ is approximately proportional to the similarity solution obtained from an initial condition of the same magnitude that is perfectly localised at $x=0$.

Thus, diffusion tends to "destroy information" about the initial conditions over time: regardless of the choice of $u_0$, such details are "smeared out" and the solution asymptotically approaches the same profile as an initially localised distribution.³

This behaviour is typical of diffusion processes and generalises to nonlinear diffusion equations, whose solutions also tend to approach self-similar forms as $t\rightarrow\infty$; although the proof is more difficult since in the nonlinear case it is not possible to express the solution in terms of a Fourier transform as we have done above.

1. Introducing a complex-valued substitution like this is a little shady, but can be rigorously justified by considering a contour integral. ↩
2. Some care is needed here. With $x^2/t$ held fixed, the $x/t$ term converges to 0. However, the numerator $\lambda$ in the $\lambda / t^{1/2}$ term ranges across all real values, including very large values of $\lambda$, and so $\lambda/t^{1/2}$ is not uniformly small. The resolution is to argue that for large values of $\lambda$, the $\exp(-\lambda^2)$ factor becomes asymptotically small and so the contributions to the integral from these tails is small. The argument can be made rigorous by splitting the integral into regions and putting appropriate bounds on each region. ↩
3. The first- and second-order corrections to this approximation account for the possibility that the initial condition might be localised about some $x = a \neq 0$, and has some spread. ↩