The Heat Equation on ℝⁿ

This chapter solves the heat equation on all of space. The solution is a convolution with one function, the heat kernel

and almost everything about the equation can be read off from it. It is positive, it has total mass $1$, it is a Gaussian whose width grows like $\sqrt t$, and it is smooth. So solutions are averages of their initial data with Gaussian weights: they are instantly smooth, they spread at infinite speed, and their maxima can only decrease.

The heat kernel is also a probability law: the distribution of a Brownian particle. That is how Einstein used it in 1905 and how Perrin confirmed it, and it is the sense in which Perelman's and Bamler's heat kernels on Ricci flows are probability measures. And the Gaussian is the function Perelman puts at the centre of his entropy: this chapter's rehearsal computes it.

By the end of this chapter you will be able to:

• derive the heat kernel by self-similarity and by the Fourier transform, and list its properties;
• solve the Cauchy problem by convolution and prove that solutions are smooth for $t > 0$;
• prove a maximum principle and uniqueness on bounded domains, and explain why uniqueness on $\mathbb{R}^n$ needs a growth condition;
• use energy methods for uniqueness and decay;
• interpret the heat kernel as the law of Brownian motion, and estimate diffusion from physical data.

3.1Brownian motion

3.2The heat kernel

There are two derivations, both already half done.

By scaling (6A.1 What a PDE Is). A solution of the form $t^{-n/2}F(|x|/\sqrt t)$ conserves $\int u\,dx$, and substituting into the heat equation forces $F(r) = Ce^{-r^2/4}$. Choosing $C$ to make the total integral $1$ gives $\Phi$ (Exercise 3.5).

By the Fourier transform (4A.5 The Fourier Transform). With $\hat u(\xi, t) = \int u(x, t)e^{-2\pi ix\cdot\xi}\,dx$, the heat equation becomes $\partial_t\hat u = -4\pi^2|\xi|^2\hat u$, so $\hat u(\xi, t) = e^{-4\pi^2t|\xi|^2}\hat g(\xi)$. A product of transforms is the transform of a convolution, and the inverse transform of the Gaussian $e^{-4\pi^2t|\xi|^2}$ is $\Phi(\cdot, t)$ (computed by contour shifting in the last exercise of 5A.3 Residues and Fourier Transforms). So $u = \Phi(\cdot, t) * g$.

The semigroup property says that running the heat equation for time $s$ and then for time $t$ is the same as running it for time $t + s$. The approximate-identity property says that as $t \to 0$ the kernel concentrates all its mass at the origin, so $\Phi(\cdot, t) \to \delta_0$: the heat kernel is the temperature after a unit of heat is released at a point at time $0$ (Figure 3.1).

3.3Solving the Cauchy problem

Three features distinguish the heat equation from the wave equation and from ODEs.

Instant smoothing. The initial data $g$ need only be bounded and continuous (in fact bounded and measurable will do, with convergence almost everywhere), yet $u$ is $C^\infty$ for every $t > 0$. Moreover, differentiating the kernel gives explicit bounds:

Each derivative costs a factor $t^{-1/2}$, exactly as parabolic scaling predicts (Exercise 3.7). These smoothing estimates are the model for Shi's derivative estimates for the Ricci flow (11A.3 Short-Time Existence and Uniqueness).

Infinite speed of propagation. If $g \geq 0$ is not identically zero, then $u(x, t) > 0$ for every $x$ and every $t > 0$, because $\Phi > 0$ everywhere. A disturbance is felt everywhere immediately, though only with weight $e^{-|x|^2/4t}$ at distance $|x|$. The heat equation is a model, accurate over the scales where the random-walk picture applies, not a law that violates relativity.

Inhomogeneous equations. The solution of $u_t - \Delta u = f$ with $u(\cdot, 0) = 0$ is given by Duhamel's principle: superpose the solutions started at each earlier time $s$ with data $f(\cdot, s)$,

the same idea as variation of constants for linear ODE (2B.10 Ordinary Differential Equations). It is the starting point of the fixed-point arguments for nonlinear equations in 6A.7 Nonlinear Parabolic Equations.

3.4Maximum principle and uniqueness

On a bounded region the heat equation is well posed with initial values and boundary values. Let $U \subset \mathbb{R}^n$ be bounded and open, and let $U_T = U\times(0, T]$ be the space-time cylinder. Its parabolic boundary $\Gamma_T$ is the bottom and the sides, $\big(U\times\{0\}\big)\cup\big(\partial U\times[0, T]\big)$: everything except the top, the part of the boundary where data are prescribed (Figure 3.2).

Applied to the difference of two solutions and its negative, the maximum principle gives uniqueness and stability: two solutions with the same initial and boundary values agree, and if the data differ by at most $\varepsilon$ the solutions differ by at most $\varepsilon$.

On all of $\mathbb{R}^n$ there is no lateral boundary, and something must replace it. Uniqueness holds among solutions that do not grow too fast: if $u$ solves the heat equation on $\mathbb{R}^n\times(0, T]$, is continuous up to $t = 0$ with $u(\cdot, 0) = 0$, and satisfies $|u(x, t)| \leq Ae^{a|x|^2}$ for some constants $A$, $a$, then $u \equiv 0$ (Evans, section 2.3.3). The growth condition cannot be dropped. Andrey Tychonoff constructed in 1935 a non-zero solution on $\mathbb{R}\times(0, \infty)$ with zero initial values,

which grows faster than any $e^{a|x|^2}$ as $|x| \to \infty$ (Exercise 3.10): heat that arrives "from infinity" in no time.

3.5Energy methods

A second route to uniqueness, which needs no maximum principle and works for systems: measure the size of a solution by an integral. For a solution on the bounded region $U$ with $u = 0$ on $\partial U$, let $e(t) = \int_Uu^2\,dx$. Then, integrating by parts (1A.10 Divergence, Curl and the Integral Theorems),

So $e$ never increases; if two solutions have the same data, the energy of their difference is $0$ at the start and so always. More is true. By the Poincaré inequality $\int_U|\nabla u|^2 \geq \lambda_1\int_Uu^2$ (4A.9 Sobolev Spaces), $e' \leq -2\lambda_1e$, so $e(t) \leq e^{-2\lambda_1t}e(0)$: the solution decays exponentially at a rate set by the first Dirichlet eigenvalue of $U$ (4A.7 Compact Operators and Spectra). Showing that $\log e(t)$ is a convex function of $t$, by one more differentiation, proves backward uniqueness: two solutions that agree at time $T$ agreed at all earlier times (Evans, section 2.3.4), even though the backward problem is ill-posed (6A.1 What a PDE Is).

3.6The heat kernel is a probability law

The factor $2$ is a convention: probabilists normalise Brownian motion so that its generator is $\frac12\Delta$, while the heat equation here has $\Delta$. In words: the temperature at $x$ at time $t$ is the average of the initial temperature over the endpoints of random paths started at $x$. The density $p(\cdot, t)$ of a diffusing particle evolves by the heat equation, which in this role is called the Fokker–Planck equation of Brownian motion (Figure 3.3). Harmonic functions have the same interpretation, with stopping at the boundary in place of a fixed time: the solution of the Dirichlet problem at $x$ is the expected boundary value where a Brownian path from $x$ first exits (the continuous version of the gambler's ruin in 6A.2 Harmonic Functions's exercise on graphs).

3.7History

Joseph Fourier solved the heat equation on the line with what is in effect the Gaussian kernel in his 1822 Théorie analytique de la chaleur. Kelvin's estimate of the age of the Earth appeared in 1862 in "On the secular cooling of the Earth". Einstein's paper on Brownian motion was one of his papers of 1905; Marian Smoluchowski reached similar results independently in 1906; Perrin's experiments date from 1908–09. Andrey Tychonoff's non-uniqueness example appeared in 1935. The Black–Scholes and Merton papers were published in 1973.

3.8Exercises

Solution

(a) The integral factors into $n$ one-dimensional ones, each $\int e^{-s^2/4t}ds = \sqrt{4\pi t}$. (b) In one dimension, $\int s^2e^{-s^2/4t}ds = 2t\sqrt{4\pi t}$, by integrating by parts or differentiating $\int e^{-as^2}ds = \sqrt{\pi/a}$ in $a$ at $a = 1/4t$: $\int s^2e^{-as^2} = \frac{\sqrt\pi}{2}a^{-3/2}$, and $\frac{\sqrt\pi}{2}(4t)^{3/2}/\sqrt{4\pi t} = 2t$. Summing over coordinates gives $2nt$.

Solution

$\partial_t\Phi = \big(-\frac{n}{2t} + \frac{|x|^2}{4t^2}\big)\Phi$. $\Delta\Phi = \operatorname{div}\big(-\frac{x}{2t}\Phi\big) = -\frac{n}{2t}\Phi + \frac{|x|^2}{4t^2}\Phi$.

Solution

$\int|\partial_x\Phi| = 2\int_0^\infty\frac{x}{2t}\Phi\,dx = \frac1t\cdot\frac{1}{\sqrt{4\pi t}}\int_0^\infty xe^{-x^2/4t}dx = \frac1t\cdot\frac{2t}{\sqrt{4\pi t}} = \frac{1}{\sqrt{\pi t}}$. Then $|u_x| = |\int\partial_x\Phi(x - y, t)g(y)\,dy| \leq \sup|g|\int|\partial_x\Phi|$. The bound must be invariant under the scaling $u(\lambda x, \lambda^2t)$, which multiplies $u_x$ by $\lambda$ and $t$ by $\lambda^{-2}$, so it must be $Ct^{-1/2}\sup|g|$.

Solution

$k_BT = 4.05 \times 10^{-21}$ J; $6\pi\eta a = 9.42 \times 10^{-9}$ kg/s; $D \approx 4.3 \times 10^{-13}$ m²/s; $\sqrt{2D\cdot30} \approx 5.1 \times 10^{-6}$ m. From $\langle x^2\rangle = 2Dt$, $D = 3.0 \times 10^{-11}/60 = 5.0 \times 10^{-13}$, and $N_A = \frac{RT}{6\pi\eta aD} = \frac{8.314 \times 293}{9.42 \times 10^{-9}\times5.0 \times 10^{-13}} \approx 5.2 \times 10^{23}$. Small errors in $a$ and $\eta$ change the answer proportionally, which is one reason Perrin's careful selection of uniform particles mattered.

Solution

(a) $\partial_x\operatorname{erf}(x/2\sqrt{\kappa t}) = \frac{2}{\sqrt\pi}\cdot\frac{1}{2\sqrt{\kappa t}}$ at $x = 0$. (b) $G = 1/15.24 \approx 0.0656$ °F/m; $t = 98 \times 10^6 \times 3.156 \times 10^7 \approx 3.09 \times 10^{15}$ s; $\kappa = V^2/(\pi tG^2) \approx 4.9 \times 10^7/(3.1416 \times 3.09 \times 10^{15}\times0.0043) \approx 1.2 \times 10^{-6}$ m²/s. (c) Halved; multiplied by $(10/7)^2 \approx 2$, about $200$ million years, which is the figure Kelvin also gives for that melting temperature.

Solution

(a) $\partial_t = -\partial_\tau$, and $\partial_\tau\big((4\pi\tau)^{-n/2}e^{-|x|^2/4\tau}\big) = \Delta(\cdots)$ because it is the heat kernel in the variable $\tau$; so $\partial_tu = -\Delta u$. (b) $\nabla f = \frac{x}{2\tau}$, $|\nabla f|^2 = \frac{|x|^2}{4\tau^2} = \frac f\tau$, $\nabla^2f = \frac{1}{2\tau}I$, $\Delta f = \frac{n}{2\tau}$. Then $2\Delta f - |\nabla f|^2 + \frac{f - n}{\tau} = \frac n\tau - \frac f\tau + \frac f\tau - \frac n\tau = 0$.