Calculus of several variables

Differentiation, extrema and integration in several variables.

In this chapter 7

Tick a unit to add it to your PDF basket.

Differentiation and extrema
Integration

Differentiation and extrema

Method

Partial derivatives and the chain rule

Adapted

Add to PDF

Defines partial derivatives and the total derivative of a map between Euclidean spaces, identifies the total derivative with the Jacobian matrix, and gives the chain rule as a matrix product together with the $C^2$ condition under which mixed partials commute.

Definition(Partial derivative)AdaptedAdapted from Lee (2013), App. C, p. 644; Riley, Hobson & Bence (2006), §5.1, Eqs. (5.1)–(5.2) and p. 152; do Carmo (2016), Appendix, p. 127

Let $U \subseteq \mathbb{R}^n$ be open, let $f : U \to \mathbb{R}$ , let $a \in U$ and let $e_j$ be the $j$ th standard basis vector. The $j$ th partial derivative of $f$ at $a$ is

\frac{\partial f}{\partial x_j}(a) \;=\; \lim_{h \to 0} \frac{f(a + h e_j) - f(a)}{h} \in \mathbb{R},

whenever that limit exists, $h$ ranging over the real numbers with $a + he_j \in U$ . It is also written $\partial_j f(a)$ or $f_{x_j}(a)$ . For a map $F = (F_1,\dots,F_m) : U \to \mathbb{R}^m$ , the partial derivatives of $F$ are those of its components, $\partial F_i/\partial x_j$ .

Higher partial derivatives are the partial derivatives of partial derivatives, with

\frac{\partial^2 f}{\partial x_j \partial x_k}(a) = \frac{\partial}{\partial x_j}\!\left(\frac{\partial f}{\partial x_k}\right)\!(a) :

the rightmost index is differentiated first. We read the subscript notation the same way round, $f_{x_j x_k} = \partial^2 f/\partial x_j \partial x_k$ , which is the convention of Riley, Hobson and Bence (p. 152) and of do Carmo (p. 127): both write $f_{xy} = \partial^2 f/\partial x\,\partial y = \partial_x(\partial_y f)$ . Other texts read the subscripts in the order of differentiation instead, so that $f_{xy}$ means $\partial_y(\partial_x f) = \partial^2 f/\partial y\, \partial x$ . The two readings agree whenever $f$ is $C^2$ , which is why the difference is usually passed over; for a function that is not $C^2$ they can disagree (see the remark).

$F$ is of class $C^k$ on $U$ ( $k \ge 1$ ) if every partial derivative of $F$ of order at most $k$ exists at every point of $U$ and is continuous on $U$ ; $C^0$ means continuous, and $C^\infty$ (smooth) means $C^k$ for every $k$ .

Definition(Differentiability and the total derivative)AdaptedAdapted from Lee (2013), App. C, Eqs. (C.1)–(C.2) and Exercise C.1, p. 642

With $U \subseteq \mathbb{R}^n$ open and $a \in U$ , a map $F : U \to \mathbb{R}^m$ is differentiable at $a$ if there is a linear map $L : \mathbb{R}^n \to \mathbb{R}^m$ with

\lim_{v \to 0} \frac{\lvert F(a + v) - F(a) - Lv \rvert}{\lvert v \rvert} = 0 ,

the norms being Euclidean and $v$ ranging over $\mathbb{R}^n \setminus \{0\}$ with $a + v \in U$ . Such an $L$ is unique; it is written $DF(a)$ and called the total derivative of $F$ at $a$ . Equivalently $F(a+v) = F(a) + DF(a)v + R(v)$ with $\lvert R(v)\rvert / \lvert v \rvert \to 0$ as $v \to 0$ , so $DF(a)$ is the linear map that approximates the increment of $F$ to first order.

Definition(Jacobian matrix, total differential, directional derivative)AdaptedAdapted from Lee (2013), App. C, p. 644 and Eq. (C.6), p. 647; Riley, Hobson & Bence (2006), §5.2, Eqs. (5.5)–(5.6), p. 154 and §6.4.3, p. 205

Let $U \subseteq \mathbb{R}^n$ be open and $F : U \to \mathbb{R}^m$ .

The Jacobian matrix of $F$ at $a \in U$ is the $m \times n$ matrix $J_F(a) = \bigl(\partial F_i/\partial x_j (a)\bigr)_{i \le m,\ j \le n}$ , whenever these partial derivatives exist. If $m = n$ , its determinant is the Jacobian determinant, also written $\partial(F_1,\dots,F_n)/\partial(x_1,\dots,x_n)$ .
For real-valued $f$ differentiable at $a$ , the linear form $Df(a)$ written in the coordinate differentials $dx_j$ (the linear forms $v \mapsto v_j$ ) is the total differential $df = \sum_{j=1}^{n} (\partial f/\partial x_j)\, dx_j$ .
For $v \in \mathbb{R}^n$ , the directional derivative of $f$ at $a$ along $v$ is $D_v f(a) = \frac{d}{dt} f(a + tv)\big|_{t=0}$ , when that one-variable derivative exists. No normalisation of $v$ is assumed.

Method(Partial derivatives, continuity and the total derivative)AdaptedAdapted from Lee (2013), App. C, Exercise C.1 (p. 642), Exercise C.2(a) (p. 643), Prop. C.8 and Exercise C.9 (p. 646), Prop. C.10 and the display after Eq. (C.6) (p. 647)

Let $U \subseteq \mathbb{R}^n$ be open, $F : U \to \mathbb{R}^m$ and $a \in U$ .

If $F$ is differentiable at $a$ , then $DF(a)$ is unique, $F$ is continuous at $a$ , every partial derivative $\partial F_i/\partial x_j(a)$ exists, and the matrix of $DF(a)$ in the standard bases is the Jacobian matrix: $DF(a)\,v = J_F(a)\,v$ for all $v \in \mathbb{R}^n$ . In particular $D_v f(a) = \sum_{j} v_j\, \partial f/\partial x_j (a)$ for real-valued $f$ differentiable at $a$ .
If $F$ is of class $C^1$ on $U$ — every $\partial F_i/\partial x_j$ exists at every point of $U$ and is continuous on $U$ — then $F$ is differentiable at every point of $U$ .
$F$ is differentiable at $a$ if and only if each component $F_1, \dots, F_m$ is differentiable at $a$ , and then the rows of $J_F(a)$ are the $J_{F_i}(a)$ .

Neither of the following converses holds: existence of all partial derivatives at $a$ does not give differentiability at $a$ , and differentiability at every point of $U$ does not give $C^1$ (see the remark).

Idea of the proofAdaptedAdapted from Lee (2013), App. C, Prop. C.8, p. 646

Everything in part 1 comes from testing the definition of differentiability along a single direction: put $v = te_j$ , divide by $t$ , and the remainder term dies because it is small compared with $\lvert v\rvert = \lvert t\rvert$ . That one substitution gives uniqueness, the existence of the partials and the identification of the matrix at once. Part 2 goes the other way and is the only part that needs work: the increment is split into $n$ steps, each changing one coordinate, so that the one-variable mean value theorem applies to each; continuity of the partials at $a$ is then exactly what makes the $n$ intermediate values collapse onto their values at $a$ .

ProofAdaptedAdapted from Lee (2013), App. C, Prop. C.8, p. 646

Part 1. Let $F$ be differentiable at $a$ with $F(a+v) = F(a) + Lv + R(v)$ , $\lvert R(v)\rvert/\lvert v\rvert \to 0$ .

Uniqueness. Suppose $L$ and $L'$ both satisfy the definition. Fix $u \in \mathbb{R}^n$ with $\lvert u \rvert = 1$ and put $v = tu$ with $t \to 0$ . Subtracting the two expansions, $\lvert (L - L')(tu)\rvert/\lvert tu\rvert \to 0$ ; but by linearity that quotient equals $\lvert (L-L')u \rvert$ , a number independent of $t$ . So $(L - L')u = 0$ for every unit $u$ , hence $L = L'$ .

Continuity. A linear map between Euclidean spaces satisfies $\lvert Lv\rvert \le C\lvert v\rvert$ for some constant $C$ (take $C$ to be the square root of the sum of the squares of the matrix entries, by the Cauchy-Schwarz inequality). Then $\lvert F(a+v) - F(a)\rvert \le C\lvert v\rvert + \lvert R(v)\rvert \to 0$ as $v \to 0$ , which is continuity at $a$ .

The partial derivatives and the matrix. Put $v = he_j$ with $h$ real and small, so that $\lvert v\rvert = \lvert h\rvert$ . Then

\frac{F(a + he_j) - F(a)}{h} = Le_j + \frac{R(he_j)}{h}, \qquad \left\lvert \frac{R(he_j)}{h}\right\rvert = \frac{\lvert R(he_j)\rvert}{\lvert he_j\rvert} \longrightarrow 0 .

So the difference quotient converges, which means each component limit exists: $\partial F_i/\partial x_j(a)$ exists and equals the $i$ th component of $Le_j$ . Since $Le_j$ is the $j$ th column of the matrix of $L$ in the standard bases, that matrix is $\bigl(\partial F_i/\partial x_j(a)\bigr) = J_F(a)$ .

The directional derivative. For real-valued $f$ differentiable at $a$ and $v \ne 0$ , the same substitution with $v$ in place of $e_j$ gives $\bigl(f(a+tv) - f(a)\bigr)/t = Df(a)v + R(tv)/t$ , and $\lvert R(tv)/t\rvert = \lvert v\rvert\,\lvert R(tv)\rvert/\lvert tv\rvert \to 0$ . Hence $D_vf(a)$ exists and equals $Df(a)v = \sum_j v_j\,\partial f/\partial x_j(a)$ , the last step by the matrix identification. For $v = 0$ both sides are $0$ .

Part 2. By part 3, proved below, it is enough to treat a real-valued $f$ of class $C^1$ . Fix $a \in U$ and $r > 0$ with the closed ball $\lvert v\rvert \le r$ contained in $U$ . For $\lvert v \rvert < r$ set $v^{(0)} = 0$ and $v^{(k)} = (v_1,\dots,v_k,0,\dots,0)$ , and telescope the increment one coordinate at a time:

f(a+v) - f(a) = \sum_{k=1}^{n}\Bigl[ f\bigl(a + v^{(k)}\bigr) - f\bigl(a + v^{(k-1)}\bigr) \Bigr].

The $k$ th bracket changes only the $k$ th coordinate, so it is the increment of the one-variable function $t \mapsto f\bigl(a + v^{(k-1)} + te_k\bigr)$ over $0 \le t \le v_k$ (or $v_k \le t \le 0$ ). That function is differentiable, because $\partial f/\partial x_k$ exists at every point of $U$ — this is where the existence half of the $C^1$ hypothesis is used — so the mean value theorem gives a point $c_k$ on the segment from $a + v^{(k-1)}$ to $a + v^{(k)}$ with

f\bigl(a+v^{(k)}\bigr) - f\bigl(a+v^{(k-1)}\bigr) = v_k\,\frac{\partial f}{\partial x_k}(c_k) .

Every point of that segment has the form $a + (v_1,\dots,v_{k-1},s,0,\dots,0)$ with $\lvert s\rvert \le \lvert v_k\rvert$ , so $\lvert c_k - a\rvert \le \lvert v\rvert$ and $c_k$ lies in the ball. Subtracting the candidate linear map,

\Bigl\lvert f(a+v) - f(a) - \sum_{k=1}^{n} v_k \frac{\partial f}{\partial x_k}(a)\Bigr\rvert = \Bigl\lvert \sum_{k=1}^{n} v_k\Bigl(\frac{\partial f}{\partial x_k}(c_k) - \frac{\partial f}{\partial x_k}(a)\Bigr)\Bigr\rvert \le \lvert v\rvert \sum_{k=1}^{n}\Bigl\lvert \frac{\partial f}{\partial x_k}(c_k) - \frac{\partial f}{\partial x_k}(a)\Bigr\rvert ,

using $\lvert v_k\rvert \le \lvert v\rvert$ . Divide by $\lvert v\rvert$ and let $v \to 0$ : each $c_k \to a$ because $\lvert c_k - a\rvert \le \lvert v\rvert$ , and each $\partial f/\partial x_k$ is continuous at $a$ — this is where continuity, as opposed to mere existence, is used, and it is the only place — so the sum tends to $0$ . Hence $f$ is differentiable at $a$ with $Df(a)v = \sum_k v_k\,\partial f/\partial x_k(a)$ .

Part 3. The Euclidean norm of a vector dominates the modulus of each of its components and is dominated by their sum. So if $F$ is differentiable at $a$ with remainder $R$ , then each component remainder satisfies $\lvert R_i(v)\rvert \le \lvert R(v)\rvert = o(\lvert v\rvert)$ , and $F_i$ is differentiable with $DF_i(a)$ the $i$ th row. Conversely, if each $F_i$ is differentiable at $a$ with linear map $L_i$ and remainder $R_i$ , let $L$ be the linear map whose $i$ th row is the matrix of $L_i$ ; then $\lvert R(v)\rvert \le \sum_i \lvert R_i(v)\rvert$ , a sum of $m$ terms each $o(\lvert v\rvert)$ , so $F$ is differentiable at $a$ with $DF(a) = L$ . In either direction the rows of $J_F(a)$ are the $J_{F_i}(a)$ , by the matrix identification in part 1. $\square$

Method(Chain rule)AdaptedAdapted from Lee (2013), App. C, Prop. C.3, p. 643; Lee (2013), App. C, Cor. C.11, p. 647

Let $V \subseteq \mathbb{R}^n$ and $W \subseteq \mathbb{R}^m$ be open, $F : V \to W$ and $G : W \to \mathbb{R}^p$ .

Total-derivative form. If $F$ is differentiable at $a \in V$ and $G$ is differentiable at $F(a)$ , then $G \circ F$ is differentiable at $a$ and $D(G \circ F)(a) = DG(F(a)) \circ DF(a)$ , which in matrices is the product of a $p \times m$ with an $m \times n$ matrix, $J_{G \circ F}(a) = J_G\bigl(F(a)\bigr)\, J_F(a).$
Partial-derivative form. If $F$ is of class $C^1$ on $V$ and $G$ is of class $C^1$ on $W$ , then $G \circ F$ is of class $C^1$ on $V$ and, writing $y = F(x)$ for the intermediate variables, $\frac{\partial (G_i \circ F)}{\partial x_j}(x) = \sum_{k=1}^{m} \frac{\partial G_i}{\partial y_k}\bigl(F(x)\bigr)\, \frac{\partial F_k}{\partial x_j}(x) \qquad (1 \le i \le p,\ 1 \le j \le n).$ If $F$ and $G$ are both of class $C^k$ for some $k \ge 1$ , or both smooth, so is $G \circ F$ . (Lee's Cor. C.11 states the two ends of this, parts (a) and (b): the cases $k = 1$ and smooth. The paragraph after it, p. 647, goes no further than to observe that, where there is differentiability enough to support it, induction on the chain rule yields whatever second- and higher-order formulas one needs; it does not assert that a composite of $C^k$ maps is again $C^k$ . The intermediate $k$ is therefore ours, and the proof below establishes it.)

Idea of the proofAdaptedAdapted from Lee (2013), App. C, Prop. C.3, pp. 643-644

Composing the two first-order expansions is almost the whole proof: substituting the increment of $F$ into the expansion of $G$ produces the wanted linear term $DG(F(a))\,DF(a)\,v$ plus two errors. One is the linear map $DG(F(a))$ applied to the error of $F$ , which is small because a linear map cannot inflate a vector by more than a fixed factor. The other is the error of $G$ evaluated at an increment that is itself at most a constant times $\lvert v\rvert$ , so it too is small compared with $\lvert v\rvert$ . The partial-derivative form is then the matrix product written out entry by entry.

ProofAdaptedAdapted from Lee (2013), App. C, Prop. C.3, pp. 643-644; Lee (2013), App. C, Cor. C.11, p. 647

Part 1. Put $b = F(a)$ , $A = DF(a)$ and $B = DG(b)$ , and name the two remainders:

F(a+v) = b + Av + r(v), \qquad G(b+w) = G(b) + Bw + s(w),

with $\lvert r(v)\rvert/\lvert v\rvert \to 0$ and $\lvert s(w)\rvert/\lvert w\rvert \to 0$ ; set $s(0) = 0$ , consistent with the definition. Let $w(v) = F(a+v) - b = Av + r(v)$ , which is the increment at which $G$ is evaluated. Substituting,

G\bigl(F(a+v)\bigr) - G\bigl(F(a)\bigr) - BAv \;=\; B\,r(v) \; + \; s\bigl(w(v)\bigr) .

The first error. As in part 1 of the previous statement, $\lvert By\rvert \le C'\lvert y\rvert$ for a constant $C'$ , so $\lvert B\,r(v)\rvert \le C'\lvert r(v)\rvert$ , which is $o(\lvert v\rvert)$ because $F$ is differentiable at $a$ — the first hypothesis.

The second error. Likewise $\lvert Av\rvert \le C\lvert v\rvert$ , and $\lvert r(v)\rvert \le \lvert v\rvert$ once $\lvert v\rvert$ is small enough, so $\lvert w(v)\rvert \le (C+1)\lvert v\rvert$ for small $\lvert v\rvert$ . Let $\varepsilon > 0$ . Because $G$ is differentiable at $b$ — the second hypothesis — there is $\delta > 0$ with $\lvert s(w)\rvert \le \varepsilon\lvert w\rvert$ whenever $\lvert w\rvert \le \delta$ . For $\lvert v\rvert$ small enough that $\lvert w(v)\rvert \le \delta$ we then have $\lvert s(w(v))\rvert \le \varepsilon(C+1)\lvert v\rvert$ . Since $\varepsilon$ was arbitrary, $\lvert s(w(v))\rvert/\lvert v\rvert \to 0$ .

Conclusion. The right-hand side of the display is therefore $o(\lvert v\rvert)$ , so $G \circ F$ is differentiable at $a$ with total derivative the composite linear map $BA = DG(F(a)) \circ DF(a)$ . Passing to matrices in the standard bases, where composition of linear maps is matrix multiplication, gives $J_{G\circ F}(a) = J_G(F(a))\,J_F(a)$ , a $p \times m$ matrix times an $m \times n$ one.

Part 2. If $F$ and $G$ are $C^1$ , then by part 2 of the previous statement each is differentiable at every point, so part 1 applies at every $x \in V$ . Reading off the $(i,j)$ entry of the matrix product, and using part 1 of the previous statement to identify entries of a Jacobian with partial derivatives,

\frac{\partial (G_i\circ F)}{\partial x_j}(x) = \sum_{k=1}^{m}\frac{\partial G_i}{\partial y_k}\bigl(F(x)\bigr)\,\frac{\partial F_k}{\partial x_j}(x) .

Each summand is a product of $\partial G_i/\partial y_k \circ F$ , continuous as a composite of continuous maps, with $\partial F_k/\partial x_j$ , continuous by hypothesis; a finite sum of products of continuous functions is continuous, so every first partial derivative of $G \circ F$ is continuous and $G\circ F$ is $C^1$ .

For higher orders, argue by induction on the order $k$ , writing the summation index over the intermediate variables as $\ell$ to keep the two uses of the letter apart. The displayed formula writes each first partial derivative of $G\circ F$ as a sum of products of (i) first partial derivatives $\partial G_i/\partial y_\ell$ of $G$ composed with $F$ and (ii) first partial derivatives $\partial F_\ell/\partial x_j$ of $F$ . If $F$ and $G$ are $C^k$ with $k \ge 2$ , those factors are $C^{k-1}$ maps — the composite $(\partial G_i/\partial y_\ell)\circ F$ by the induction hypothesis, applied to the $C^{k-1}$ maps $\partial G_i/\partial y_\ell$ and $F$ — and sums and products of $C^{k-1}$ functions are $C^{k-1}$ . So the first partial derivatives of $G\circ F$ are $C^{k-1}$ , which says that $G\circ F$ is $C^k$ . The smooth case is the case of every $k$ . $\square$

Method(Equality of mixed partial derivatives)AdaptedAdapted from Lee (2013), App. C, Prop. C.6 and Cor. C.7, p. 646

Let $U \subseteq \mathbb{R}^n$ be open and let $F : U \to \mathbb{R}^m$ be of class $C^2$ . Then

\frac{\partial^2 F_i}{\partial x_j \, \partial x_k}(x) = \frac{\partial^2 F_i}{\partial x_k \, \partial x_j}(x) \qquad (i \le m,\ j, k \le n,\ x \in U).

Lee's Cor. C.7 draws the consequence for smooth $F$ . By the same argument a function of class $C^r$ has all its mixed partial derivatives of order at most $r$ independent of the order of differentiation, and a smooth function has this at every order; the extension to finite $r$ is ours, and the proof below gives it. For $n = 2$ this leaves three independent second derivatives: $f_{xx}$ , $f_{yy}$ and the common value $f_{xy} = f_{yx}$ .

Idea of the proofOriginal

The two mixed partials are two ways of evaluating one object: the second difference of $F$ over a small coordinate rectangle, $\Delta = F(a + se_j + te_k) - F(a+se_j) - F(a+te_k) + F(a)$ . Applying the one-variable mean value theorem first in one direction and then in the other turns $\Delta/st$ into a second partial derivative at some interior point of the rectangle, and doing the two directions in the opposite order turns it into the other second partial derivative at another interior point. Shrinking the rectangle sends both interior points to $a$ , and continuity of the second partials — the $C^2$ hypothesis — forces the two limits to agree.

ProofOriginal

It is enough to prove the identity for a real-valued function of two variables. The assertion involves only the component $F_i$ and only the two coordinates $x_j$ , $x_k$ , so fix $i$ and fix the remaining $n-2$ coordinates at their values at the point $x$ ; write $f$ for the resulting function of the two free coordinates, which is $C^2$ on an open set because $F$ is, and write $\partial_1$ , $\partial_2$ for the two partial derivatives, with $\partial_2\partial_1 f$ meaning $\partial_2(\partial_1 f)$ , so that no convention about the order of the indices is needed. The claim is $\partial_1\partial_2 f = \partial_2\partial_1 f$ at the point in question, which we take to be $a$ .

Fix $r > 0$ with the square $\lvert s\rvert, \lvert t\rvert \le r$ about $a$ inside the open set, and for $0 < \lvert s\rvert, \lvert t\rvert \le r$ define the second difference

\Delta(s,t) = f(a + se_1 + te_2) - f(a + se_1) - f(a + te_2) + f(a) .

First order of differentiation. Put $\varphi(u) = f(a + ue_1 + te_2) - f(a + ue_1)$ , so that $\Delta(s,t) = \varphi(s) - \varphi(0)$ . Since $\partial_1 f$ exists on the square, $\varphi$ is differentiable with $\varphi'(u) = \partial_1 f(a+ue_1+te_2) - \partial_1 f(a+ue_1)$ , and the mean value theorem gives a $\sigma$ between $0$ and $s$ with $\Delta(s,t) = s\,\varphi'(\sigma)$ . Now apply the mean value theorem again, to the function $\tau \mapsto \partial_1 f(a + \sigma e_1 + \tau e_2)$ , which is differentiable because $\partial_2\partial_1 f$ exists on the square: there is a $\tau$ between $0$ and $t$ with $\varphi'(\sigma) = t\,\partial_2\partial_1 f(a + \sigma e_1 + \tau e_2)$ . Hence

\frac{\Delta(s,t)}{st} = \partial_2\partial_1 f(a + \sigma e_1 + \tau e_2), \qquad \lvert \sigma\rvert \le \lvert s\rvert,\ \lvert\tau\rvert \le \lvert t\rvert .

The other order. $\Delta(s,t)$ is symmetric under exchanging the roles of the two coordinates together with $s$ and $t$ : grouping its four terms the other way, $\Delta(s,t) = \psi(t) - \psi(0)$ with $\psi(u) = f(a + se_1 + ue_2) - f(a + ue_2)$ . The same two applications of the mean value theorem, in the reverse order, give points $\sigma'$ between $0$ and $s$ and $\tau'$ between $0$ and $t$ with

\frac{\Delta(s,t)}{st} = \partial_1\partial_2 f(a + \sigma' e_1 + \tau' e_2) .

Let the rectangle shrink. Equating the two expressions and letting $(s,t) \to (0,0)$ : both arguments converge to $a$ , since $\lvert\sigma\rvert,\lvert\sigma'\rvert \le \lvert s\rvert$ and $\lvert\tau\rvert,\lvert\tau'\rvert \le \lvert t\rvert$ , and both $\partial_2\partial_1 f$ and $\partial_1\partial_2 f$ are continuous at $a$ by the $C^2$ hypothesis. Therefore $\partial_2\partial_1 f(a) = \partial_1\partial_2 f(a)$ . This is the only step that uses continuity, and the remark below shows that the conclusion fails without it.

Higher orders for a $C^r$ function. Let $F$ be $C^r$ and consider a mixed partial derivative of order $m \le r$ , taken in some order of the indices. Any permutation of $m$ symbols is a product of transpositions of adjacent symbols, so it suffices to show that exchanging two consecutive differentiations leaves the result unchanged. Those two act on a partial derivative of $F$ of order $m - 2 \le r-2$ , which is a $C^2$ function because $F$ is $C^r$ ; the case already proved applies to it, and the exchange is legitimate. The smooth case is the case of every $r$ . $\square$

NoteOriginal

What it means. Differentiability in several variables is the existence of a linear approximation, not the existence of the partial derivatives, which record only what $F$ does along the $n$ coordinate directions. Class $C^1$ is the practical bridge between the two, and it is the hypothesis almost every later theorem uses. Once $F$ is differentiable, everything first-order about it sits in one matrix, and the chain rule is matrix multiplication in the order "outer map first, evaluated at the image point".

Method(Differentiating a composite map)AdaptedAdapted from Lee (2013), App. C, Cor. C.11, p. 647; Riley, Hobson & Bence (2006), §5.2, Eq. (5.8), p. 154

To differentiate $x \mapsto G(F(x))$ :

Name the stages: $y = F(x)$ with $x \in V \subseteq \mathbb{R}^n$ , $y \in W \subseteq \mathbb{R}^m$ , and $z = G(y) \in \mathbb{R}^p$ , recording which letters are independent at each stage.
Check that $F$ is $C^1$ near $x$ and $G$ is $C^1$ near $F(x)$ , or at least that each is differentiable at the relevant point.
Multiply in the right order and at the right point: $J_{G\circ F}(x) = J_G(F(x))\,J_F(x)$ , the outer factor evaluated at $y = F(x)$ . The $(i,j)$ entry is the sum over intermediate variables of part 2 of the chain rule, one term per path from $x_j$ to $z_i$ .
Specialise if the chain is a curve. With $F = \gamma$ a differentiable curve, $\frac{d}{dt}G(\gamma(t)) = J_G(\gamma(t))\,\gamma'(t)$ , and $\gamma(t) = a + tv$ recovers $D_v G(a)$ .
Iterate for higher derivatives, applying step 3 to each first-order expression together with the product rule; once both maps are $C^2$ the order of differentiation does not matter.

When it fails. Step 3 needs differentiability, not just the existence of partial derivatives (see the remark). It also fails silently if the outer Jacobian is evaluated at the wrong point, or if a variable is held fixed at one stage and varied at another.

RemarkOriginal

Where the hypotheses are used, and what fails without them.

Partial derivatives alone do not give continuity. Put $f(x,y) = xy/(x^2+y^2)$ for $(x,y) \ne (0,0)$ and $f(0,0) = 0$ . Then $f(x,0) = f(0,y) = 0$ , so $f_x(0,0) = f_y(0,0) = 0$ exist; but $f = \tfrac12$ on the line $y = x$ away from the origin, so $f$ is not continuous there, hence not differentiable, and any first-order formula built from the two partials alone — the chain rule included — is false for it at the origin.
Continuity of the partials in part 2 is sufficient, not necessary. Let $f(x,y) = (x^2+y^2)\sin(1/(x^2+y^2))$ with $f(0,0) = 0$ . Since $\lvert f(v)\rvert \le \lvert v\rvert^2$ , $f$ is differentiable at the origin with $Df(0,0) = 0$ ; but $f_x$ contains $-2x\cos(1/(x^2+y^2))/(x^2+y^2)$ , unbounded near the origin, so $f$ is differentiable everywhere without being $C^1$ .
$C^2$ cannot be weakened to mere existence of the mixed partials. Let $f(x,y) = xy(x^2-y^2)/(x^2+y^2)$ away from the origin, $f(0,0) = 0$ . From the definition $f_x(0,y) = -y$ and $f_y(x,0) = x$ , so both mixed second partials exist at the origin yet $\partial^2 f/\partial y\,\partial x(0,0) = -1 \ne 1 = \partial^2 f/\partial x\,\partial y(0,0)$ ; what fails is their continuity. For such an $f$ the two readings of $f_{xy}$ named in the first definition also disagree.
A weaker hypothesis in one source. Riley, Hobson and Bence (p. 153) assert the equality of mixed partials under a weaker condition than the $C^2$ hypothesis of Lee's Prop. C.6 used above: they ask only for continuity of the second partials at the point concerned. The statement here follows Lee.
Notation that differs between fields. Physics and thermodynamics texts write $(\partial f/\partial x)_y$ , listing the held variables as subscripts (Riley, Hobson and Bence, §5.1, p. 152), because one symbol often denotes different functions of different variable sets. Those books call $df = \sum_j (\partial f/\partial x_j)dx_j$ the total differential and $df/dx_1$ along a curve the total derivative; here these are $Df(a)$ in coordinate differentials and the case $n=1$ of the chain rule. The symbol $\partial(x,y,z)/\partial(u,v,w)$ is standard in both fields and returns in the change-of-variables formula. A reader turning to Riley, Hobson and Bence (p. 205) should note that they write that determinant with the array transposed relative to the Jacobian matrix defined above: their rows are indexed by $u$ , $v$ , $w$ and their columns by $x$ , $y$ , $z$ , so their first row is $\partial x/\partial u$ , $\partial y/\partial u$ , $\partial z/\partial u$ , whereas the first row of $J_F$ here is $\partial x/\partial u$ , $\partial x/\partial v$ , $\partial x/\partial w$ . Transposition leaves the determinant unchanged but not the matrix.
Smoothness assumed by the source. Lee's Eq. (C.6) defines the directional derivative for smooth $f$ , and the identity $D_v f(a) = \sum_j v_j \partial f/\partial x_j(a) = Df(a)v$ is the display just after it (p. 647). Differentiability of $f$ at the single point $a$ already implies it, by the chain rule applied to $t \mapsto f(a+tv)$ , which is the form stated in part 1; the implication does not reverse, since all directional derivatives at $a$ can exist without $f$ being differentiable there.

Example(A chain rule computed both ways)Original

Take

F : \mathbb{R}^2 \to \mathbb{R}^3, \quad F(u,v) = \bigl(u + v^2,\; uv,\; e^{\,u-v}\bigr), \qquad G : \mathbb{R}^3 \to \mathbb{R}^2, \quad G(x,y,z) = \bigl(x^2y,\; y + xz\bigr),

and compute $J_{G\circ F}$ at $(u,v) = (1,2)$ , following the method.

Stages. $(u,v) \mapsto (x,y,z) = F(u,v) \mapsto G(x,y,z)$ , with $n = 2$ , $m = 3$ , $p = 2$ . At the point, $F(1,2) = (1 + 4,\; 2,\; e^{-1}) = (5,\, 2,\, e^{-1})$ , so the outer factor must be evaluated at $(x,y,z) = (5, 2, e^{-1})$ and nowhere else.
Hypotheses. Every component of $F$ and of $G$ is a polynomial in the coordinates and the exponential, so all partial derivatives of all orders exist and are continuous: both maps are smooth, hence $C^1$ , hence differentiable, and part 2 of the chain rule applies.
The two Jacobians.
$J_F(u,v) = \begin{pmatrix} 1 & 2v \\ v & u \\ e^{\,u-v} & -e^{\,u-v}\end{pmatrix} \ \xrightarrow{(1,2)}\ \begin{pmatrix} 1 & 4 \\ 2 & 1 \\ e^{-1} & -e^{-1}\end{pmatrix}, \qquad J_G(x,y,z) = \begin{pmatrix} 2xy & x^2 & 0 \\ z & 1 & x \end{pmatrix} \ \xrightarrow{(5,2,e^{-1})}\ \begin{pmatrix} 20 & 25 & 0 \\ e^{-1} & 1 & 5\end{pmatrix}.$
Multiply, outer factor first.
$J_{G\circ F}(1,2) = \begin{pmatrix} 20 & 25 & 0 \\ e^{-1} & 1 & 5\end{pmatrix} \begin{pmatrix} 1 & 4 \\ 2 & 1 \\ e^{-1} & -e^{-1}\end{pmatrix} = \begin{pmatrix} 70 & 105 \\ 2 + 6e^{-1} & 1 - e^{-1} \end{pmatrix} = \begin{pmatrix} 70 & 105 \\ 4.207277 & 0.632121\end{pmatrix}.$
For instance the $(1,1)$ entry is $20\cdot1 + 25\cdot2 + 0\cdot e^{-1} = 70$ : one term for each intermediate variable, and the term through $z$ drops out because $\partial G_1/\partial z = 0$ .
Check by composing first. Writing $H = G \circ F$ out,
$H_1(u,v) = (u+v^2)^2\,uv, \qquad H_2(u,v) = uv + (u+v^2)\,e^{\,u-v},$
and differentiating directly: $\partial H_1/\partial u = 2(u+v^2)uv + (u+v^2)^2v$ , which at $(1,2)$ is $2\cdot5\cdot2 + 25\cdot2 = 70$ ; $\partial H_1/\partial v = 4v(u+v^2)uv + (u+v^2)^2u = 4\cdot2\cdot5\cdot2 + 25 = 105$ ; $\partial H_2/\partial u = v + (1 + u + v^2)e^{u-v} = 2 + 6e^{-1}$ ; and $\partial H_2/\partial v = u + (2v - u - v^2)e^{u-v} = 1 + (4-5)e^{-1} = 1 - e^{-1}$ . All four entries agree.
A directional derivative. By part 1 of the first statement, the derivative of $H_1$ at $(1,2)$ along $w = (1,-1)$ is $J_{H_1}(1,2)\,w = 70 - 105 = -35$ . Equivalently $\frac{d}{dt}H_1(1+t, 2-t)\big|_{t=0} = -35$ , which is step 4 of the method with the curve $\gamma(t) = (1,2) + t(1,-1)$ .
Mixed partials. $H_2$ is smooth, so the equality of mixed partials applies. Differentiating the expressions in step 5 once more, $\partial^2H_2/\partial u\,\partial v = \partial_u\bigl(u + (2v-u-v^2)e^{u-v}\bigr) = 1 + (2v - u - v^2 - 1)e^{u-v}$ , which is $1 - 2e^{-1} = 0.264241$ at $(1,2)$ ; and $\partial^2H_2/\partial v\,\partial u = \partial_v\bigl(v + (1+u+v^2)e^{u-v}\bigr) = 1 + (2v - 1 - u - v^2)e^{u-v}$ , the same expression and the same value $1 - 2e^{-1}$ . For a function that is not $C^2$ this agreement fails, as the third bullet of the remark shows.

ApplicationOriginal

Maxwell relations in thermodynamics. The equality of mixed partials is the whole content of the Maxwell relations. The Helmholtz free energy $A(T,V)$ satisfies $dA = -S\,dT - p\,dV$ , so $\partial A/\partial T = -S$ and $\partial A/\partial V = -p$ as functions of $(T,V)$ . If $A$ is $C^2$ on the region of the phase diagram concerned, then $\partial^2A/\partial T\,\partial V = \partial^2A/\partial V\,\partial T$ , that is

\left(\frac{\partial S}{\partial V}\right)_T = \left(\frac{\partial p}{\partial T}\right)_V ,

which converts a quantity that cannot be measured directly (the volume dependence of the entropy) into one that can (the pressure rise at fixed volume). Each thermodynamic potential gives one such relation the same way. The $C^2$ hypothesis is not idle: it fails at a phase boundary, where second derivatives of the potential jump or diverge, and the classification of transitions into first and second order is precisely a statement about which derivative of the potential first misbehaves.

Changing coordinates. Every curvilinear calculation is the chain rule with $F$ the coordinate map. For plane polars, $F(r,\theta) = (r\cos\theta, r\sin\theta)$ has $J_F = \begin{pmatrix} \cos\theta & -r\sin\theta\\ \sin\theta & r\cos\theta\end{pmatrix}$ with determinant $r$ ; the chain rule then expresses $\partial/\partial r$ and $\partial/\partial\theta$ in terms of $\partial/\partial x$ and $\partial/\partial y$ , and inverting that relation (possible exactly where the Jacobian determinant is non-zero, so not at the origin) is how the Laplacian is transformed. The determinant $r$ reappears as the area element, which is the subject of the change-of-variables formula.

Material derivatives in continuum mechanics. A quantity $\phi$ carried by a moving fluid is observed along a particle path $t \mapsto x(t)$ , so the rate of change seen by the particle is the chain rule applied to $t \mapsto \phi(x(t), t)$ : $\frac{D\phi}{Dt} = \partial\phi/\partial t + \sum_j \dot x_j\,\partial\phi/\partial x_j = \phi_t + v\cdot\nabla\phi$ . Step 4 of the method, the case of a curve, is exactly this computation, and the warning in the method about which variables are held fixed at which stage is the reason the convective term is so often dropped by mistake: $\partial\phi/\partial t$ at fixed position and $D\phi/Dt$ following the fluid are different functions written with similar symbols.

Theorem

Taylor expansion and critical points in several variables

Adapted

Add to PDF

Gives the second-order Taylor expansion of a $C^2$ function of $n$ variables in Peano and mean-value form, and classifies a critical point as a strict local minimum, a strict local maximum or a saddle according as the Hessian there is positive definite, negative definite or indefinite, with no conclusion when it is singular and semidefinite.

Prerequisites: Partial derivatives and the chain rule, Taylor's theorem, Quadratic forms and the Rayleigh quotient

Definition(Gradient, Hessian, critical point, definiteness)AdaptedAdapted from Riley, Hobson & Bence (2006), §5.8, pp. 162–163 and pp. 166–167; Lee (2013), App. C, Prop. C.6, p. 646

Let $U \subseteq \mathbb{R}^n$ be open and $f : U \to \mathbb{R}$ .

If $f$ is differentiable at $a \in U$ , its gradient there is the column vector $\nabla f(a) = \bigl(\partial f/\partial x_1(a), \dots, \partial f/\partial x_n(a)\bigr)^{\mathsf T}$ , so that $Df(a)h = \nabla f(a) \cdot h$ .
If $f$ is of class $C^2$ on $U$ , its Hessian at $a$ is the $n \times n$ matrix $H_f(a) = \left( \frac{\partial^2 f}{\partial x_i \, \partial x_j}(a) \right)_{i,j \le n},$ the entry in row $i$ and column $j$ being $\partial_i\bigl(\partial_j f\bigr)(a)$ — the rightmost index differentiated first, the convention fixed in the unit on partial derivatives. That unit's equality of mixed partials makes $H_f(a)$ symmetric because $f$ is $C^2$ , so which index is read first makes no difference here.
$a \in U$ is a critical point (or stationary point) of $f$ if $f$ is differentiable at $a$ and $\nabla f(a) = 0$ .
$a$ is a local minimum of $f$ if $f(x) \ge f(a)$ for all $x$ in some ball $B_r(a) \subseteq U$ , and a strict local minimum if $f(x) > f(a)$ for all $x \in B_r(a)$ with $x \neq a$ ; local maxima are defined with the inequalities reversed. A saddle point is a critical point that is neither a local minimum nor a local maximum.
A real symmetric $n \times n$ matrix $H$ is positive definite if $h^{\mathsf T} H h > 0$ for every $h \ne 0$ , positive semidefinite if $h^{\mathsf T} H h \ge 0$ for every $h$ , negative (semi)definite if $-H$ is positive (semi)definite, and indefinite if $u^{\mathsf T} H u > 0$ and $w^{\mathsf T} H w < 0$ for some $u, w \in \mathbb{R}^n$ . Its eigenvalues are real, and the eigenvalue test of quadratic forms makes $H$ positive definite exactly when all of them are positive and positive semidefinite exactly when none of them is negative. Applied to $-H$ , whose eigenvalues are $-\lambda$ for each eigenvalue $\lambda$ of $H$ , the same test makes $H$ negative definite exactly when all are negative and negative semidefinite exactly when none is positive. The indefinite case is not among the ones that unit states, but it follows from those two: being indefinite is exactly failing both semidefiniteness conditions, so $H$ is indefinite exactly when some eigenvalue is positive and some is negative.

Theorem(Second-order Taylor expansion)AdaptedAdapted from Riley, Hobson & Bence (2006), §5.7, Eqs. (5.18) and (5.20), pp. 160–162; Lee (2013), App. C, Thm C.15, p. 648

Let $U \subseteq \mathbb{R}^n$ be open, let $f : U \to \mathbb{R}$ be of class $C^2$ on $U$ and let $a \in U$ . Then, as $h \to 0$ ,

f(a + h) = f(a) + \sum_{j=1}^{n} \frac{\partial f}{\partial x_j}(a)\, h_j + \frac{1}{2} \sum_{j=1}^{n} \sum_{k=1}^{n} \frac{\partial^2 f}{\partial x_j \, \partial x_k}(a)\, h_j h_k + o\bigl(\lvert h \rvert^2\bigr),

that is, the difference between $f(a+h)$ and the quadratic polynomial written here, divided by $\lvert h \rvert^2$ , tends to $0$ . The single sum is $\nabla f(a)\cdot h$ and the double sum is $h^{\mathsf T} H_f(a)\, h$ .

If $f$ is of class $C^{k+1}$ on $U$ and $W \subseteq U$ is convex with $a \in W$ , the expansion continues to order $k$ for every $x \in W$ : $f(x)$ equals its $k$ th Taylor polynomial at $a$ plus a remainder given by an explicit integral of the $(k+1)$ st partial derivatives along the segment from $a$ to $x$ .

Neither cited source states the $o(\lvert h\rvert^2)$ remainder; what each of them does state, and how the display follows, is set out in the remark on provenance below.

Theorem(Mean-value form of the second-order expansion)Original

Let $U \subseteq \mathbb{R}^n$ be open, let $f : U \to \mathbb{R}$ be of class $C^2$ on $U$ , let $a \in U$ and let $h \in \mathbb{R}^n$ be such that the closed segment $\{a + th : 0 \le t \le 1\}$ lies in $U$ . Then there is some $\theta \in (0,1)$ , depending on $a$ and $h$ , with

f(a + h) = f(a) + \nabla f(a) \cdot h + \tfrac{1}{2}\, h^{\mathsf T} H_f(a + \theta h)\, h .

Neither cited source states this form; see the remark on provenance below.

Idea of the proof(Both expansions, from one variable)Original

Restrict $f$ to the segment through $a$ in the direction $h$ . The function $g(t) = f(a+th)$ of the single variable $t$ has $g'(t) = \nabla f(a+th)\cdot h$ and $g''(t) = h^{\mathsf T} H_f(a+th)\,h$ , by the chain rule, so one-variable Taylor's theorem with the Lagrange remainder applied to $g$ on $[0,1]$ is the mean-value form. The $o(\lvert h\rvert^2)$ form follows from it by replacing $H_f(a+\theta h)$ with $H_f(a)$ : what that costs is at most $\lvert h \rvert^2$ times the amount by which the second partial derivatives vary on the ball of radius $\lvert h\rvert$ about $a$ , and that amount tends to $0$ as $h \to 0$ precisely because $f$ is $C^2$ .

Proof(Mean-value form)AdaptedAdapted from Lee (2013), App. C, Thm C.15 and its proof, p. 648

Let $U$ , $f$ , $a$ and $h$ be as in the statement. If $h = 0$ every $\theta$ works, so assume $h \neq 0$ .

The one-variable function is defined on an open set containing $[0,1]$ . The set $T = \{t \in \mathbb{R} : a + th \in U\}$ is the preimage of the open set $U$ under the continuous map $t \mapsto a+th$ , hence open, and it contains $[0,1]$ because the closed segment lies in $U$ — the hypothesis on $h$ . Every $t \in [0,1]$ , the end points included, therefore has an open interval about it inside $T$ , so derivatives taken below at $t = 0$ and $t = 1$ are two-sided ones, exactly like those at interior points. No other property of $T$ is used.
Its derivatives are the directional derivatives of $f$ . On $T$ put $g(t) = f(a+th)$ . The map $t \mapsto a+th$ is affine, so the chain rule applies and gives $g'(t) = \sum_{j=1}^{n} \frac{\partial f}{\partial x_j}(a+th)\, h_j = \nabla f(a+th)\cdot h ,$ and applying it once more to each term, $g''(t) = \sum_{j=1}^{n}\sum_{k=1}^{n} \frac{\partial^2 f}{\partial x_k \,\partial x_j}(a+th)\, h_j h_k = h^{\mathsf T} H_f(a+th)\, h .$ Because $f$ is $C^2$ on $U$ , all the first and second partial derivatives appearing here are continuous, so $g''$ is continuous on $T$ and in particular on $[0,1]$ .
Apply the one-variable theorem. Taylor's theorem with $N = 1$ , expansion point $0$ and evaluation point $1$ requires exactly that $g''$ exist and be continuous on $[0,1]$ , which the previous step supplies. Its Lagrange remainder gives a point $\theta$ between $0$ and $1$ with $g(1) = g(0) + g'(0) + \tfrac{1}{2} g''(\theta) .$ As that unit states the conclusion, $\theta$ lies in the closed interval $[0,1]$ , and that is all the rest of this unit uses: every later appeal to the mean-value form needs only $\lvert \theta h\rvert \le \lvert h\rvert$ . The strict inclusion $\theta \in (0,1)$ asserted in the statement comes from the sharper version recorded in that unit's remark on its hypotheses, where it is enough that $g'$ be continuous on $[0,1]$ and $g''$ exist at every interior point, and $\theta$ may then be taken strictly between the end points. Both conditions hold here, since $g$ is $C^2$ on the open set $T \supseteq [0,1]$ .
Read the terms back. $g(1) = f(a+h)$ , $g(0) = f(a)$ , $g'(0) = \nabla f(a)\cdot h$ and $g''(\theta) = h^{\mathsf T}H_f(a+\theta h)h$ , which is the display of the statement. $\square$

The reduction to the single variable $t$ along the segment is the route of Lee's Thm C.15; he carries it to order $k$ and keeps an integral remainder, whereas the Lagrange remainder of our one-variable unit is what produces the intermediate point $\theta$ .

Proof(The

o(\lvert h \rvert^2)

form)Original

Write $H = H_f(a)$ and let $P(h) = f(a) + \nabla f(a)\cdot h + \tfrac12 h^{\mathsf T}Hh$ be the quadratic polynomial of the statement, whose two sums are the ones displayed there.

A ball on which the mean-value form is available. $U$ is open, so fix $r > 0$ with $B_r(a) \subseteq U$ . For $0 < \lvert h\rvert < r$ the closed segment from $a$ to $a+h$ lies in the convex set $B_r(a) \subseteq U$ , so the mean-value form applies and yields a $\theta = \theta(h) \in (0,1)$ with $f(a+h) - P(h) = \tfrac12\, h^{\mathsf T}\bigl[ H_f(a+\theta h) - H \bigr] h .$
A crude bound for a quadratic form. For any real $n \times n$ matrix $A$ , $\bigl\lvert h^{\mathsf T} A h \bigr\rvert \le \sum_{j,k} \lvert A_{jk}\rvert \lvert h_j\rvert \lvert h_k \rvert \le \Bigl(\max_{j,k} \lvert A_{jk}\rvert\Bigr) \Bigl( \sum_{j} \lvert h_j \rvert \Bigr)^{2} \le n \Bigl(\max_{j,k}\lvert A_{jk}\rvert\Bigr) \lvert h \rvert^2 ,$ the last step by the Cauchy–Schwarz inequality $\bigl(\sum_j \lvert h_j\rvert\bigr)^2 \le n \sum_j h_j^2$ .
The matrix in question is small. For $0 \le s < r$ put $\omega(s) = \max_{j,k} \; \sup_{\lvert u \rvert \le s} \left\lvert \frac{\partial^2 f}{\partial x_j\,\partial x_k}(a+u) - \frac{\partial^2 f}{\partial x_j \,\partial x_k}(a) \right\rvert .$ Each second partial derivative is continuous at $a$ , because $f$ is $C^2$ ; there are finitely many of them, so $\omega(s) \to 0$ as $s \to 0^+$ . Since $\lvert \theta h\rvert \le \lvert h \rvert$ , the entries of $H_f(a+\theta h) - H$ are bounded in modulus by $\omega(\lvert h\rvert)$ .
Conclude. Combining the three previous steps, $\lvert f(a+h) - P(h)\rvert \le \tfrac12 n\, \omega(\lvert h\rvert)\, \lvert h\rvert^2$ for $0 < \lvert h \rvert < r$ , so $\frac{\lvert f(a+h) - P(h)\rvert}{\lvert h\rvert^{2}} \le \tfrac12 n\, \omega(\lvert h\rvert) \xrightarrow[h \to 0]{} 0 ,$ which is the assertion $f(a+h) = P(h) + o(\lvert h\rvert^2)$ . $\square$

Only a ball about $a$ was used, so this form needs no convexity hypothesis on $U$ , as the statement records.

RemarkAdaptedAdapted from Riley, Hobson & Bence (2006), §5.7, Eqs. (5.18) and (5.20), pp. 160–162, and §5.8, pp. 164 and 167; Lee (2013), App. C, Thm C.15 and Eqs. (C.7)–(C.9), p. 648

What the sources state, and what is ours.

The $o(\lvert h\rvert^2)$ remainder. Riley, Hobson and Bence write Eqs. (5.18) and (5.20) as a formal series ending in " $+\cdots$ ", with no remainder estimate at all. Lee's Thm C.15 does give an exact remainder, but it is not the one displayed above: taking $k = 1$ there — the case whose hypothesis is exactly $f \in C^2$ — produces the first-order Taylor polynomial, with $f(a+h) = f(a) + \nabla f(a)\cdot h + R_1(h), \qquad R_1(h) = \sum_{j=1}^{n}\sum_{k=1}^{n} h_j h_k \int_0^1 (1-t)\, \frac{\partial^2 f}{\partial x_j \, \partial x_k}(a+th)\, dt ,$ valid whenever the segment from $a$ to $a+h$ lies in a convex subset of $U$ , so for all small $\lvert h\rvert$ . This $R_1(h)$ is not $o(\lvert h\rvert^2)$ : it carries the quadratic term itself. (For $n = 1$ , $f(x) = x^2$ and $a = 0$ the first-order polynomial is $0$ while $R_1(h) = h^2$ .)
How the display follows from it, which is ours. $R_1$ splits. Since $\int_0^1 (1-t)\,dt = \tfrac12$ , subtracting and adding the value of each second partial derivative at $a$ gives $R_1(h) = \tfrac12\, h^{\mathsf T} H_f(a)\, h \;+\; \sum_{j=1}^{n}\sum_{k=1}^{n} h_j h_k \int_0^1 (1-t)\left[ \frac{\partial^2 f}{\partial x_j\,\partial x_k}(a+th) - \frac{\partial^2 f}{\partial x_j\,\partial x_k}(a) \right] dt ,$ and it is the second piece alone that is $o(\lvert h\rvert^2)$ , because the second partial derivatives are continuous at $a$ and so the brackets tend to $0$ uniformly in $t \in [0,1]$ as $h \to 0$ . The displayed expansion is therefore a consequence of Lee's case $k = 1$ and not that case itself. The splitting is ours, and so is the route actually taken in the two proofs above, which reaches the same $o(\lvert h\rvert^2)$ statement from the mean-value form instead of from the integral remainder.
The mean-value form. Neither source states it. It is Taylor's theorem with the Lagrange remainder, applied to the one-variable function $t \mapsto f(a+th)$ on $[0,1]$ , and the reduction to that one variable is the only thing it shares with Lee's Thm C.15.
The word "strict" in the classification. Riley, Hobson and Bence give the two-variable criteria on p. 164 and the eigenvalue criteria on p. 167, and both are phrased with "minima" and "maxima"; neither says strict. The strengthening to strict local extrema in parts 1 and 2 of the next statement goes beyond those pages and is ours. What supplies it is the remainder estimate in the proof below: it bounds $\lvert R(h)\rvert$ by a fraction of $\lambda_1\lvert h\rvert^2$ , which gives a sign of $f(a+h) - f(a)$ at every $h \ne 0$ in a punctured ball, not merely a weak inequality.

Theorem(Classification of a critical point by the Hessian)AdaptedAdapted from Riley, Hobson & Bence (2006), §5.8, pp. 164 and 167

Let $U \subseteq \mathbb{R}^n$ be open, let $f : U \to \mathbb{R}$ be of class $C^2$ on $U$ , let $a \in U$ satisfy $\nabla f(a) = 0$ , and put $H = H_f(a)$ .

If $H$ is positive definite, $a$ is a strict local minimum of $f$ .
If $H$ is negative definite, $a$ is a strict local maximum of $f$ .
If $H$ is indefinite, $a$ is a saddle point; more precisely, every ball $B_r(a) \subseteq U$ contains points where $f > f(a)$ and points where $f < f(a)$ .

In two variables, with $D = \det H = f_{xx}f_{yy} - f_{xy}^2$ evaluated at $a$ :

$D > 0$ and $f_{xx} > 0$ : strict local minimum (then $f_{yy} > 0$ automatically);
$D > 0$ and $f_{xx} < 0$ : strict local maximum (then $f_{yy} < 0$ automatically);
$D < 0$ : saddle point;
$D = 0$ : the test does not classify the point. (It is not quite silent: the next statement says what is still excluded and what is not.)

Idea of the proof(Why the Hessian decides)AdaptedAdapted from Riley, Hobson & Bence (2006), §5.8, pp. 166–167

At a critical point the linear term of the expansion is absent, so $f(a+h) - f(a)$ is the quadratic form $\tfrac12 h^{\mathsf T}Hh$ plus an error small compared with $\lvert h\rvert^2$ . In an orthonormal eigenbasis of $H$ the form becomes $\tfrac12\sum_r \lambda_r a_r^2$ , so the signs of the eigenvalues say which signs the form takes. A definite form does more: on the unit sphere it is bounded away from zero, by the extreme eigenvalue, and that margin is a fixed multiple of $\lvert h \rvert^2$ , which the error cannot overcome. When the form takes both signs one need only walk along the two directions that witness this.

Proof(Classification of a critical point)AdaptedAdapted from Riley, Hobson & Bence (2006), §5.8, pp. 164, 166–167

Throughout, $\nabla f(a) = 0$ and $H = H_f(a)$ , which is symmetric because $f$ is $C^2$ . The second-order expansion proved above then reads

f(a+h) - f(a) = \tfrac12\, h^{\mathsf T}Hh + R(h), \qquad \frac{R(h)}{\lvert h\rvert^2} \to 0 \text{ as } h \to 0 .

The five numbered items below are the steps of the argument, not the parts of the statement: items 1 and 2 together prove part 1, item 3 proves part 2, item 4 proves part 3, and item 5 is the two-variable specialisation. Each item names the part it serves.

Part 1, the quadratic form has a margin. The eigenvalues $\lambda_1 \le \cdots \le \lambda_n$ of $H$ are real and, by the eigenvalue test for definiteness, all positive, so $\lambda_1 > 0$ . The Rayleigh bound of quadratic forms gives $h^{\mathsf T}Hh \ge \lambda_1 \lvert h\rvert^2$ for every $h$ .
Part 1 concluded: positive definite gives a strict local minimum. Because $R(h) = o(\lvert h\rvert^2)$ and $U$ is open, there is $r > 0$ with $B_r(a) \subseteq U$ and $\lvert R(h)\rvert \le \tfrac14 \lambda_1 \lvert h\rvert^2$ whenever $\lvert h \rvert < r$ — this is the definition of the limit in the display, with $\tfrac14\lambda_1$ as the tolerance, and it is here that positive definiteness, not mere semidefiniteness, is used: the tolerance has to be a positive multiple of $\lvert h\rvert^2$ . Hence for $0 < \lvert h \rvert < r$ $f(a+h) - f(a) \ge \tfrac12\lambda_1\lvert h\rvert^2 - \tfrac14\lambda_1\lvert h\rvert^2 = \tfrac14 \lambda_1 \lvert h \rvert^2 > 0 ,$ so $a$ is a strict local minimum.
Part 2: negative definite. Apply part 1 to $-f$ , which is $C^2$ , has $\nabla(-f)(a) = 0$ and $H_{-f}(a) = -H$ positive definite; a strict local minimum of $-f$ at $a$ is a strict local maximum of $f$ .
Part 3: indefinite gives a saddle. Choose $u$ with $u^{\mathsf T}Hu > 0$ ; then $u \neq 0$ . For $t \neq 0$ with $a + tu \in U$ , $f(a+tu) - f(a) = t^2\left( \tfrac12 u^{\mathsf T}Hu + \frac{R(tu)}{t^2} \right), \qquad \frac{R(tu)}{t^{2}} = \lvert u\rvert^{2}\,\frac{R(tu)}{\lvert tu\rvert^{2}} \xrightarrow[t \to 0]{} 0 ,$ so the bracket is positive for all small enough $\lvert t\rvert \neq 0$ and $f(a+tu) > f(a)$ there. The same argument with $w$ , for which $w^{\mathsf T}Hw < 0$ , gives points with $f < f(a)$ . Every ball $B_\rho(a) \subseteq U$ contains $a + tu$ and $a+tw$ for all small $\lvert t \rvert \ne 0$ , so each contains points of both kinds; in particular $a$ is neither a local minimum nor a local maximum, that is, a saddle point.
The two-variable criteria, from parts 1–3. Here $H$ $H$ has eigenvalues $\lambda_1, \lambda_2$ $λ_{1}, λ_{2}$ with $\lambda_1\lambda_2 = \det H = D$ $λ_{1} λ_{2} = det H = D$ and $\lambda_1 + \lambda_2 = \operatorname{tr} H = f_{xx} + f_{yy}$ $λ_{1} + λ_{2} = tr H = f_{xx} + f_{y y}$ , and $f_{xx} = e_1^{\mathsf T}He_1$ $f_{xx} = e_{1}^{T} H e_{1}$ with $e_1 = (1,0)^{\mathsf T}$ $e_{1} = (1, 0)^{T}$ .
- If $D > 0$ the two eigenvalues share a sign. Were both negative, $H$ would be negative definite and $f_{xx} = e_1^{\mathsf T}He_1$ would be negative; so if also $f_{xx} > 0$ , both are positive, $H$ is positive definite and part 1 applies. Moreover $f_{yy} = (D + f_{xy}^2)/f_{xx} > 0$ , which is the parenthetical claim.
- If $D > 0$ and $f_{xx} < 0$ , apply the previous case to $-f$ , whose determinant is again $D$ and whose $(1,1)$ entry is $-f_{xx} > 0$ .
- If $D < 0$ the eigenvalues have opposite signs, so $H$ is indefinite and part 3 applies.
- If $D = 0$ then $0$ is an eigenvalue. If $\operatorname{tr}H \neq 0$ the eigenvalues are $0$ and $\operatorname{tr}H$ , of one sign only, so $H$ is semidefinite; if $\operatorname{tr}H = 0$ both eigenvalues vanish and $H = 0$ , since a symmetric matrix is orthogonally diagonalisable. Either way $H$ is singular and semidefinite, so none of parts 1–3 applies, and the next statement shows that the three-way classification cannot be completed from $H$ alone — although, when $H \neq 0$ , it still excludes one of the two kinds of extremum. $\square$

Reading the signs of the quadratic form off the eigenvalues of $H$ , after expanding the displacement in an orthonormal eigenbasis, is the route of Riley, Hobson and Bence (pp. 166–167, with the two-variable criteria on p. 164). Step 1 reaches the one consequence of that expansion it needs — the lower bound $\lambda_1\lvert h\rvert^2$ — by quoting the Rayleigh bound of quadratic forms instead of expanding; step 5 uses the eigenvalue signs directly, as they do. They argue with the quadratic approximation itself; the remainder estimates in steps 2 and 4, which are what turn a sign of the form into a sign of $f(a+h)-f(a)$ , are ours.

Theorem(The converse, and the cases the Hessian does not decide)Original

Let $U \subseteq \mathbb{R}^n$ be open, let $f : U \to \mathbb{R}$ be of class $C^2$ on $U$ , let $a \in U$ satisfy $\nabla f(a) = 0$ and put $H = H_f(a)$ , as in the previous statement.

If $a$ is a local minimum of $f$ , then $H$ is positive semidefinite; if $a$ is a local maximum, $H$ is negative semidefinite.
Suppose $H$ $H$ is semidefinite and singular, so that none of parts 1–3 applies (part 4 may still apply). Then:
- if $H$ is positive semidefinite and $H \neq 0$ , $a$ is not a local maximum, by part 4; but each of the three remaining possibilities occurs — a strict local minimum, a local minimum that is not strict, and no local extremum at all — and $H$ does not distinguish them, since all three occur among functions with the same $\nabla f(a) = 0$ and the same $H_f(a)$ ;
- if $H$ is negative semidefinite and $H \neq 0$ , symmetrically, $a$ is not a local minimum, and the three mirror-image possibilities — a strict local maximum, a non-strict local maximum, and no local extremum — all occur;
- if $H = 0$ , nothing is excluded: $a$ may be a strict local minimum, a strict local maximum, a non-strict local extremum, or no local extremum at all.

So in the singular semidefinite case the Hessian fixes at most the sign of the extremum that is still possible, never whether there is one.

Idea of the proof(The converse, and why the degenerate cases are undecidable)Original

Part 4 is the expansion read along a single ray: if $f$ does not fall below $f(a)$ nearby, then the quadratic term, which dominates for small $t$ , cannot be negative in any direction. Part 5 asks for the opposite of a theorem — functions that the Hessian cannot tell apart — so it has to be settled by constructing them, and the eigenbasis of $H$ hands them over. Split the coordinates into the directions where the form is positive and the directions where it vanishes; adding a quartic in the vanishing directions changes neither the gradient nor the Hessian at $a$ , and its sign, or its absence, decides what happens there.

Proof(The converse, and the cases the Hessian does not decide)Original

Write $H = H_f(a)$ , symmetric since $f$ is $C^2$ , and keep the expansion $f(a+h) - f(a) = \tfrac12 h^{\mathsf T}Hh + R(h)$ with $R(h) = o(\lvert h\rvert^2)$ , as in the previous proof.

Part 4. Let $a$ be a local minimum, so $f(a+h) \ge f(a)$ for $\lvert h\rvert < \rho$ , some $\rho > 0$ . Fix $h \ne 0$ . For $0 < t < \rho/\lvert h\rvert$ the point $a + th$ lies in that ball, so

0 \le f(a+th) - f(a) = t^2\left( \tfrac12 h^{\mathsf T}Hh + \lvert h\rvert^2 \frac{R(th)}{\lvert th\rvert^2} \right).

Dividing by $t^2 > 0$ and letting $t \to 0^+$ , the last term vanishes and $h^{\mathsf T}Hh \ge 0$ remains. As $h$ was arbitrary, $H$ is positive semidefinite. If $a$ is a local maximum, apply this to $-f$ to get $-H$ positive semidefinite, that is $H$ negative semidefinite.

Part 5, the constructions. Let $H$ be positive semidefinite and singular. By the principal-axes theorem (quadratic forms) choose an orthonormal basis $q_1, \dots, q_n$ of eigenvectors, with eigenvalues $\lambda_1, \dots, \lambda_n$ ; none is negative, by semidefiniteness, and at least one is $0$ , by singularity, so we may order them as $\lambda_1, \dots, \lambda_p > 0$ and $\lambda_{p+1} = \cdots = \lambda_n = 0$ with $0 \le p \le n-1$ . Put $y_j(x) = q_j \cdot (x - a)$ , so that $y = (y_1, \dots, y_n)$ are orthonormal coordinates centred at $a$ , and let

q(y) = \tfrac12\sum_{j \le p} \lambda_j y_j^{2},

which is the quadratic form of $H$ in these coordinates: $q(y(x)) = \tfrac12 (x-a)^{\mathsf T}H(x-a)$ .

Each $y_j$ is an affine function of $x$ with $\nabla y_j = q_j$ , so for a polynomial $\varphi$ in $y$ with no constant and no linear term, the function $x \mapsto f(a) + \varphi(y(x))$ has gradient $0$ at $a$ , and its Hessian at $a$ is $\sum_{j,k}\partial^2\varphi/\partial y_j\partial y_k(0)\, q_jq_k^{\mathsf T}$ . In particular each of

\Phi_1 = f(a) + q(y) + \sum_{j>p} y_j^{4}, \qquad \Phi_2 = f(a) + q(y), \qquad \Phi_3 = f(a) + q(y) - y_n^{4}

is a polynomial in $x$ of class $C^\infty$ with gradient $0$ and Hessian exactly $H$ at $a$ : the quartic terms contribute nothing to either, and $q$ contributes $\sum_{j\le p}\lambda_j q_jq_j^{\mathsf T} = H$ .

Now suppose in addition $H \neq 0$ , that is $p \ge 1$ .

$\Phi_1$ has a strict local minimum at $a$ — indeed a strict global one. For $y \neq 0$ both $q(y) \ge 0$ and $\sum_{j>p}y_j^4 \ge 0$ , and they vanish together only if $y_j = 0$ for $j \le p$ (as $\lambda_j > 0$ there) and for $j > p$ , that is only if $y = 0$ .
$\Phi_2$ has a local minimum at $a$ that is not strict: $q \ge 0$ always, while $\Phi_2 = f(a)$ on the whole subspace $\{y_1 = \cdots = y_p = 0\}$ , of dimension $n - p \ge 1$ .
$\Phi_3$ has no local extremum at $a$ . Along $x = a + s q_n$ (a direction with $\lambda_n = 0$ , available because $p < n$ ) it equals $f(a) - s^4 < f(a)$ for $s \neq 0$ , so $a$ is not a local minimum; along $x = a + sq_1$ it equals $f(a) + \tfrac12\lambda_1 s^2 > f(a)$ for $s \ne 0$ , so $a$ is not a local maximum either.

All three have the same gradient $0$ and the same Hessian $H$ at $a$ , which proves that no criterion depending on those two can separate them; and by part 4 none of them has a local maximum at $a$ , since $H$ has the positive eigenvalue $\lambda_1$ and so is not negative semidefinite. The negative semidefinite case is this one applied to $-f$ .

Part 5, the case $H = 0$ . Then $p = 0$ and $q \equiv 0$ . The four functions

f(a) + \lvert x - a\rvert^{4}, \qquad f(a) - \lvert x-a\rvert^{4}, \qquad f(a), \qquad f(a) + y_1(x)^{3}

are all $C^\infty$ with gradient $0$ and Hessian $0$ at $a$ , and have there, respectively, a strict local minimum, a strict local maximum, a local extremum that is both a non-strict minimum and a non-strict maximum, and no local extremum at all — the last because $y_1(x)^3$ takes both signs arbitrarily close to $a$ . So nothing is excluded. $\square$

NoteOriginal

What it means. Near a critical point a $C^2$ function looks, to the accuracy that matters, like the constant $f(a)$ plus the quadratic form $\tfrac12 h^{\mathsf T} H h$ : the first-order term has been switched off, so the second-order term decides the local shape. Classifying the critical point is then nothing but asking what signs that quadratic form takes, a question about the eigenvalues of a symmetric matrix. The test is silent exactly when $H$ is semidefinite and singular — the form vanishes along some direction without ever changing sign — because the quadratic approximation is then too coarse along those null directions and the higher-order terms take over. A form that does change sign still decides, however degenerate it is.

RemarkOriginal

Where the hypotheses are used, and what fails without them.

The segment must lie in $U$ . The mean-value form is the one-variable theorem applied to $g(t) = f(a+th)$ on $[0,1]$ , and $g$ is defined only where $a + th \in U$ ; on a domain that is not convex there may be no intermediate point at which to evaluate the Hessian. The $o(\lvert h\rvert^2)$ form is local and needs only a ball around $a$ .
$C^2$ is used twice. It makes $H$ symmetric — without it "definiteness" could refer only to the symmetric part, since $h^{\mathsf T} H h = h^{\mathsf T}\bigl(\tfrac12(H + H^{\mathsf T})\bigr)h$ for every $h$ — and it is what gives the second-order expansion a remainder small compared with $\lvert h \rvert^2$ , which is what turns a sign of the quadratic form into a sign of $f(a+h) - f(a)$ in parts 1–3.
Singular and semidefinite really is undecidable. Both $f(x,y) = x^2 + y^4$ and $g(x,y) = x^2 - y^4$ have a critical point at the origin with Hessian $\operatorname{diag}(2,0)$ , positive semidefinite and singular. The origin is a strict local minimum of $f$ , whereas $g(0,y) = -y^4 < 0$ for $y \ne 0$ and $g(x,0) = x^2 > 0$ for $x \neq 0$ , so it is a saddle point of $g$ . Neither can be a local maximum, which is part 5. The third possibility is realised by $p(x,y) = x^2$ , with the same gradient and the same Hessian at the origin again: there $p$ has a local minimum that is not strict, since $p = 0$ along the whole $y$ -axis. For $H = 0$ take $u(x,y) = (x^2+y^2)^2$ , with a strict local minimum at the origin, $-u$ , with a strict local maximum, and $v(x,y) = x^3 - 3xy^2$ , whose value on the line $y = 0$ is $x^3$ and so changes sign at the origin.
The test classifies, it does not locate. Points where $f$ fails to be differentiable, and boundary points of the region of interest, are not critical points and are not covered; on a compact set the maximum and minimum may well be attained there.
The same criteria in other books. Riley, Hobson and Bence (p. 164) give the minimum test as " $f_{xx}$ and $f_{yy}$ both positive and $f_{xy}^2 < f_{xx}f_{yy}$ ", equivalent to the determinant form above because $f_{xx} > 0$ with $f_{xx}f_{yy} > f_{xy}^2 \ge 0$ forces $f_{yy} > 0$ ; p. 167 gives the $n$ -variable test by the signs of the eigenvalues of the Hessian and notes that it may fail when some eigenvalues vanish and the non-zero ones share a sign, which is the case of part 5.
Constrained problems are different. Parts 1–5 describe a critical point of $f$ on an open set. If the variables are tied by equations, the first-order condition is not $\nabla f = 0$ but a multiplier condition. That is the subject of Lagrange multipliers, a later unit which takes this one as a prerequisite; the pointer is forward, and nothing in the present unit uses it.

Example(Two critical-point problems, in two and in three variables)Original

(a) $f(x,y) = x^4 + y^4 - 4xy$ on $\mathbb{R}^2$ .

Critical points. $\partial f/\partial x = 4x^3 - 4y$ and $\partial f/\partial y = 4y^3 - 4x$ , so $\nabla f = 0$ means $y = x^3$ and $x = y^3 = x^9$ , that is $x(x^8-1) = 0$ . Over the reals $x \in \{0, 1, -1\}$ , giving the three critical points $(0,0)$ , $(1,1)$ , $(-1,-1)$ .

Hessian. $f_{xx} = 12x^2$ , $f_{yy} = 12y^2$ , $f_{xy} = -4$ , so

H_f(x,y) = \begin{pmatrix} 12x^2 & -4 \\ -4 & 12y^2\end{pmatrix}.

Classification. At $(0,0)$ : $H = \begin{pmatrix} 0 & -4\\ -4 & 0\end{pmatrix}$ and $D = \det H = -16 < 0$ , so by the two-variable criteria $(0,0)$ is a saddle point; $f(0,0) = 0$ . At $(\pm1,\pm1)$ (both signs together): $H = \begin{pmatrix}12 & -4\\ -4 & 12\end{pmatrix}$ , $D = 144 - 16 = 128 > 0$ and $f_{xx} = 12 > 0$ , so both are strict local minima, with the common value

f(1,1) = f(-1,-1) = 1 + 1 - 4 = -2 .

The eigenvalues of that $H$ are $12 \mp 4$ , namely $8$ and $16$ , both positive, as part 1 requires.

Check that $-2$ is the global minimum. From $4xy \le 2(x^2+y^2)$ (which is $2(x-y)^2 \ge 0$ ),

f(x,y) \ge x^4 + y^4 - 2x^2 - 2y^2 = (x^2-1)^2 + (y^2-1)^2 - 2 \ge -2 ,

with equality only if $x^2 = y^2 = 1$ and $x = y$ , that is at the two minima just found. So $\min_{\mathbb{R}^2} f = -2$ .

(b) $f(x,y,z) = x^2+y^2+z^2-2xyz$ on $\mathbb{R}^3$ , where the two-variable determinant test has nothing to say and the eigenvalues must be used.

Critical points. $\nabla f = 2(x - yz,\; y - xz,\; z - xy)$ , so we need $x = yz$ , $y = xz$ , $z = xy$ . Multiplying the first two gives $xy = xyz^2$ , so $xy(1-z^2) = 0$ ; with $z = xy$ this reads $z(1-z^2) = 0$ , so $z \in \{0,1,-1\}$ , and the symmetric argument applies to $x$ and $y$ . If one coordinate is $0$ then all three are, by the equations; otherwise each is $\pm 1$ and $x = yz$ forces the product $xyz = 1$ . The critical points are therefore the origin and the four points $(1,1,1)$ , $(1,-1,-1)$ , $(-1,1,-1)$ , $(-1,-1,1)$ .

Hessian.

H_f(x,y,z) = 2\begin{pmatrix} 1 & -z & -y \\ -z & 1 & -x \\ -y & -x & 1 \end{pmatrix}.

At the origin $H = 2I$ , positive definite, so the origin is a strict local minimum, with $f = 0$ . (It is not a global minimum: $f(t,t,t) = 3t^2 - 2t^3 \to -\infty$ .)

At $(1,1,1)$ ,

H = \begin{pmatrix} 2 & -2 & -2 \\ -2 & 2 & -2 \\ -2 & -2 & 2\end{pmatrix} = 4I - 2J, \qquad J = \begin{pmatrix}1&1&1\\1&1&1\\1&1&1\end{pmatrix}.

$J$ has eigenvalue $3$ on $(1,1,1)^{\mathsf T}$ and eigenvalue $0$ on the plane $x+y+z=0$ , so $H$ has eigenvalues $4 - 2\cdot 3 = -2$ (once) and $4 - 2\cdot 0 = 4$ (twice). These have both signs, so $H$ is indefinite and $(1,1,1)$ is a saddle point, with $f(1,1,1) = 3 - 2 = 1$ . The other three points are the same by the symmetry $(x,y,z)\mapsto(-x,-y,z)$ and its relatives, which fix $f$ .

Note that $\det H = (-2)(4)(4) = -32 < 0$ here while $\det(2I) = 8 > 0$ at the origin: in three variables a negative determinant only says the eigenvalues' product is negative, which it also would be for a negative definite $3\times 3$ matrix. The sign of the determinant is a test only when $n = 2$ .

Application(Stable equilibria and transition states on a potential-energy surface)Original

A mechanical or molecular system with coordinates $x = (x_1,\dots,x_n)$ and potential energy $V(x)$ is in equilibrium exactly where $\nabla V(x) = 0$ , since the generalised forces are $-\partial V/\partial x_j$ . Which kind of equilibrium it is, is the question this unit answers with $H = H_V(x_*)$ :

$H$ positive definite: a strict local minimum of the energy, a stable equilibrium. Expanding $V$ to second order, $V \approx V(x_*) + \tfrac12 h^{\mathsf T}Hh$ , is exactly the harmonic approximation, and the squared frequencies of the small oscillations are the eigenvalues of the mass-weighted Hessian $M^{-1/2}HM^{-1/2}$ (normal modes); positive definiteness of $H$ is what makes them positive, so that the motion oscillates rather than runs away.
$H$ indefinite with exactly one negative eigenvalue and all the others positive: an index-1 saddle, the transition state of chemical kinetics. The system is stable against all displacements except along the eigenvector of the negative eigenvalue, which is the reaction coordinate; the corresponding "frequency" is imaginary, and that is the quantity transition-state theory removes from the partition function of the saddle.
$H$ singular and semidefinite: the harmonic approximation fails along the null directions. Genuine instances are the rigid translations and rotations of an isolated molecule of $N$ atoms, in the $n = 3N$ nuclear coordinates: they leave $V$ unchanged, so they are exactly zero eigenvalues of $H$ at every stationary point. The three translations always count; of the rotations, three are independent if the molecule is nonlinear, but only two if it is linear, because a rotation about the molecular axis moves no nucleus. The number of zero eigenvalues is therefore $6$ for a nonlinear molecule, $5$ for a linear one — every diatomic, and $\mathrm{CO_2}$ — and $3$ for a single atom, which has no rotations at all. The remaining eigenvalues are the vibrations, which is the familiar count of $3N - 6$ normal modes, or $3N - 5$ for a linear molecule. These zero modes are projected out before the vibrational analysis rather than classified.

As a two-coordinate model, take $V(x,y) = (x^2-1)^2 + 2y^2$ (energies in units of some $\epsilon$ , lengths in units of some $a$ ). Then $\nabla V = (4x(x^2-1),\, 4y)$ , so the stationary points are $(\pm1, 0)$ and $(0,0)$ , and $H_V = \operatorname{diag}(12x^2-4,\,4)$ .

At $(\pm 1, 0)$ : $H_V = \operatorname{diag}(8, 4)$ , positive definite — two stable equilibria, of equal energy $V = 0$ , with harmonic frequencies $\sqrt{8\epsilon/(ma^2)}$ and $\sqrt{4\epsilon/(ma^2)}$ for a particle of mass $m$ .
At $(0,0)$ : $H_V = \operatorname{diag}(-4, 4)$ , eigenvalues $-4$ and $4$ — an index-1 saddle, the transition state between the two wells, at energy $V = 1$ . The barrier is $\Delta V = 1 - 0 = 1$ (that is, $\epsilon$ ), the reaction coordinate is the $x$ direction, and the $y$ direction remains a stable vibration of the transition state with frequency $\sqrt{4\epsilon/(ma^2)}$ .

An Arrhenius rate $\propto e^{-\Delta V/k_BT}$ needs only the barrier, a number obtained from the critical points; the Hessian supplies the prefactor.

Theorem

Implicit and inverse function theorems

Adapted

Add to PDF

A map of class $C^k$ between open subsets of $\mathbb{R}^n$ whose Jacobian determinant is non-zero at a point is a $C^k$ diffeomorphism of some neighbourhood of that point onto an open set, and an equation $\Phi(x,y) = c$ whose partial Jacobian in $y$ is invertible at a solution determines $y$ locally as a $C^k$ function $g$ of $x$ , with $Dg(x) = -(\partial\Phi/\partial y)^{-1}(\partial\Phi/\partial x)$ evaluated at $(x, g(x))$ .

Prerequisites: Partial derivatives and the chain rule, Limits of sequences and completeness, Connectedness and simple connectivity

Definition(Diffeomorphism and regular point)AdaptedAdapted from Folland (1999), §2.6, p. 74; Lee (2013), App. C, p. 645

Let $\Omega, \Omega' \subseteq \mathbb{R}^n$ be open and let $k \in \{1, 2, \dots\} \cup \{\infty\}$ .

A map $F : \Omega \to \Omega'$ is a $C^k$ diffeomorphism if $F$ is bijective, $F$ is of class $C^k$ on $\Omega$ , and $F^{-1}$ is of class $C^k$ on $\Omega'$ . For $k = \infty$ this is simply a diffeomorphism.
A point $a \in \Omega$ is a regular point of a $C^1$ map $F : \Omega \to \mathbb{R}^n$ if $DF(a)$ is invertible, equivalently if the Jacobian determinant $\det J_F(a)$ is non-zero.
Following Folland, a $C^1$ map $G$ on an open $\Omega \subseteq \mathbb{R}^n$ is called a $C^1$ diffeomorphism when $G$ is injective on $\Omega$ and $D_xG$ is invertible at every $x \in \Omega$ ; by the second statement below this agrees with the definition above, with $\Omega' = G(\Omega)$ .

Theorem(Inverse function theorem)AdaptedAdapted from Lee (2013), App. C, Thm C.34, p. 657 and Prop. C.4, p. 645; Folland (1999), §2.6, p. 74

Let $\Omega \subseteq \mathbb{R}^n$ be open, let $F : \Omega \to \mathbb{R}^n$ be of class $C^k$ for some $k$ with $1 \le k \le \infty$ , and let $a \in \Omega$ be a regular point of $F$ , that is $\det J_F(a) \neq 0$ . Then there exist open sets $\Omega_0$ and $V_0$ with $a \in \Omega_0 \subseteq \Omega$ and $F(a) \in V_0$ such that

$F$ restricted to $\Omega_0$ is a bijection onto $V_0$ whose inverse $G = (F|_{\Omega_0})^{-1} : V_0 \to \Omega_0$ is of class $C^k$ ; that is, $F|_{\Omega_0} : \Omega_0 \to V_0$ is a $C^k$ diffeomorphism, and $\Omega_0$ may be taken connected;
for every $y \in V_0$ , writing $x = G(y)$ , $DG(y) = \bigl[DF(x)\bigr]^{-1}, \qquad\text{so}\qquad \det J_{G}(y) = \frac{1}{\det J_F(x)} .$

The map $F$ need not be injective on the whole of $\Omega$ , and $\Omega_0$ cannot in general be taken to be $\Omega$ (see the remark).

Idea of the proof(Inverse function theorem)AdaptedAdapted from Lee (2013), App. C, proof of Thm C.34, pp. 657–659

Normalise so that the point is the origin and the derivative there is the identity. Then $F$ differs from the identity by a map $\varphi$ whose derivative vanishes at the origin, so $\varphi$ moves points by at most half their separation on a small ball. Solving $F(x) = y$ is the same as finding a fixed point of $x \mapsto y + \varphi(x)$ , which that bound makes a contraction of a closed ball into itself, so the contraction mapping theorem supplies exactly one solution for each small $y$ — existence and uniqueness of the inverse at one stroke. The same bound shows that $F$ can shrink a distance by no more than a factor of two, $\lvert F(x') - F(x)\rvert \ge \tfrac12\lvert x' - x\rvert$ , which makes the inverse Lipschitz with constant $2$ ; the derivative of the inverse is then read off the difference quotient, and its smoothness class follows from $DG = (DF\circ G)^{-1}$ , matrix inversion being smooth.

Proof(Inverse function theorem)AdaptedAdapted from Lee (2013), App. C, Lemma C.35 and the proof of Thm C.34, pp. 657–660

Normalisation. Let $A = DF(a)$ , invertible by hypothesis, and set

\tilde F(v) = A^{-1}\bigl(F(a+v) - F(a)\bigr), \qquad v \in \tilde\Omega = \{v \in \mathbb{R}^n : a + v \in \Omega\} .

$\tilde\Omega$ is open and contains $0$ ; $\tilde F$ is $C^k$ on it, being $F$ composed on both sides with affine maps; $\tilde F(0) = 0$ ; and $D\tilde F(0) = A^{-1}A = I$ by the chain rule. The maps $v \mapsto a+v$ and $w \mapsto F(a) + Aw$ are $C^\infty$ bijections of $\mathbb{R}^n$ with $C^\infty$ inverses, so if $\tilde F$ restricted to an open connected $\tilde\Omega_0 \ni 0$ is a $C^k$ diffeomorphism onto an open $\tilde V_0 \ni 0$ , then $F$ restricted to $\Omega_0 = a + \tilde\Omega_0$ is a $C^k$ diffeomorphism onto $V_0 = F(a) + A\tilde V_0$ , and $\Omega_0$ , $V_0$ are open with $\Omega_0$ connected. Part 1 for $F$ therefore follows from part 1 for $\tilde F$ , and we may assume $a = 0$ , $F(0) = 0$ , $DF(0) = I$ .

Step 1: $F$ differs from the identity by a $\tfrac12$ -contraction. Put $\varphi(x) = x - F(x)$ on $\Omega$ , so that $\varphi$ is $C^k$ with $D\varphi(0) = I - I = 0$ . Because $F$ is of class $C^1$ — this is where continuity of the derivative, and not merely its existence, is used — the partial derivatives of $\varphi$ are continuous, so there is $r > 0$ with $\overline{B_r(0)} \subseteq \Omega$ and

\Bigl\lvert \frac{\partial \varphi_j}{\partial x_k}(x) \Bigr\rvert \le \frac{1}{2n} \quad \text{for all } j,k \text{ and all } \lvert x \rvert \le r .

The closed ball is convex, so for $x, x' \in \overline{B_r(0)}$ the fundamental theorem of calculus along the segment between them gives, for each component,

\varphi_j(x') - \varphi_j(x) = \int_0^1 \sum_{k=1}^n \frac{\partial \varphi_j}{\partial x_k}\bigl(x + t(x'-x)\bigr)\,(x'_k - x_k)\, dt ,

whence $\lvert \varphi_j(x') - \varphi_j(x)\rvert \le \frac{1}{2n}\sum_k \lvert x'_k - x_k\rvert \le \frac{1}{2\sqrt n}\lvert x' - x\rvert$ by the Cauchy–Schwarz inequality. Squaring and summing over the $n$ components,

\lvert \varphi(x') - \varphi(x) \rvert \le \tfrac12 \lvert x' - x\rvert \qquad (x, x' \in \overline{B_r(0)}). \tag{1}

Step 2: $F$ expands distances, so it is injective there. For $x, x' \in \overline{B_r(0)}$ write $x' - x = \bigl(F(x') - F(x)\bigr) + \bigl(\varphi(x') - \varphi(x)\bigr)$ . The triangle inequality and $(1)$ give $\lvert x'-x\rvert \le \lvert F(x')-F(x)\rvert + \tfrac12\lvert x'-x\rvert$ , that is

\lvert F(x') - F(x)\rvert \ge \tfrac12 \lvert x' - x\rvert . \tag{2}

In particular $F(x') = F(x)$ forces $x' = x$ .

Step 3: each $y$ with $\lvert y \rvert < r/2$ has exactly one preimage. Fix such a $y$ and put $T_y(x) = y + \varphi(x)$ on $K = \overline{B_r(0)}$ ; then $T_y(x) = x$ if and only if $F(x) = y$ . Taking $x = 0$ in $(1)$ and using $\varphi(0) = 0$ gives $\lvert \varphi(x)\rvert \le \tfrac12\lvert x\rvert$ , so for $x \in K$

\lvert T_y(x) \rvert \le \lvert y \rvert + \tfrac12\lvert x \rvert < \tfrac{r}{2} + \tfrac{r}{2} = r ,

and $T_y$ maps $K$ into $B_r(0) \subseteq K$ . By $(1)$ , $\lvert T_y(x') - T_y(x)\rvert = \lvert \varphi(x') - \varphi(x)\rvert \le \tfrac12 \lvert x' - x\rvert$ for $x, x' \in K$ , so $T_y$ is a contraction of constant $\tfrac12$ on the non-empty closed bounded set $K$ . That such a map has exactly one fixed point is the contraction mapping theorem; since that unit comes later in the handbook, here is the one instance needed. Set $x_0 = 0$ and $x_{m+1} = T_y(x_m)$ , all in $K$ . The contraction bound gives $\lvert x_{m+1} - x_m \rvert \le 2^{-m}\lvert x_1 - x_0\rvert$ by induction on $m$ , so $\sum_{m \ge 0}\lvert x_{m+1}-x_m\rvert$ converges and the sequence is Cauchy. It therefore converges in $\mathbb{R}^n$ , one coordinate at a time: $\lvert (x_m)_j - (x_{m'})_j\rvert \le \lvert x_m - x_{m'}\rvert$ makes each coordinate sequence Cauchy in $\mathbb{R}$ , hence convergent by the Cauchy criterion of limits of sequences, and if $x$ is the vector of the $n$ limits then $\lvert x_m - x\rvert \le \sum_{j}\lvert (x_m)_j - x_j\rvert \to 0$ . The limit lies in $K$ because $K$ is closed; boundedness of $K$ plays no part in this step. Letting $m \to \infty$ in $x_{m+1} = T_y(x_m)$ , with $T_y$ continuous by the same bound, gives $T_y(x) = x$ . If $x' \in K$ also satisfies $T_y(x') = x'$ then $\lvert x - x'\rvert \le \tfrac12\lvert x - x'\rvert$ , forcing $x = x'$ . Finally the displayed inequality shows $\lvert x \rvert = \lvert T_y(x)\rvert < r$ , so $x \in B_r(0)$ .

Step 4: the neighbourhoods, and continuity of the inverse. Put

V_0 = B_{r/2}(0), \qquad \Omega_0 = B_r(0) \cap F^{-1}(V_0),

the latter open because $F$ is continuous, and both containing $0$ because $F(0) = 0$ . By step 3, $F$ maps $\Omega_0$ onto $V_0$ (every $y \in V_0$ has a preimage in $B_r(0)$ , which then lies in $\Omega_0$ ), and by step 2 it does so injectively. Let $G : V_0 \to \Omega_0$ be the inverse. Substituting $x = G(y)$ and $x' = G(y')$ in $(2)$ ,

\lvert G(y') - G(y) \rvert \le 2 \lvert y' - y \rvert , \tag{3}

so $G$ is continuous and $F|_{\Omega_0}$ is a homeomorphism onto $V_0$ . Finally $\Omega_0$ is connected. The ball $V_0$ is: it is the union of the segments joining its centre to its points; each such segment is the image of an interval under a continuous map, so it is connected by items 1 and 2 of "Properties of connected spaces" in connectedness and simple connectivity; and they all contain the centre, so item 4 of that statement glues them into one connected set. Item 2 of the same statement — connectedness transported through a continuous map — then makes the image $\Omega_0 = G(V_0)$ connected.

Step 5: $G$ is differentiable, with $DG(y) = [DF(G(y))]^{-1}$ . First, $DF(x)$ is invertible at every $x \in \overline{B_r(0)}$ : if $DF(x)v = 0$ then $v = D\varphi(x)v$ , and the bound of step 1 with the same two inequalities as in $(1)$ gives $\lvert D\varphi(x)v\rvert \le \tfrac12 \lvert v\rvert$ , so $\lvert v \rvert \le \tfrac12\lvert v\rvert$ and $v = 0$ . Now fix $y \in V_0$ , put $x = G(y)$ and $L = DF(x)$ . For $y' \in V_0$ with $y' \neq y$ , put $x' = G(y')$ , which differs from $x$ because $G$ is injective; using $y' - y = F(x') - F(x)$ ,

\frac{G(y') - G(y) - L^{-1}(y'-y)}{\lvert y' - y\rvert} = - L^{-1}\!\left( \frac{F(x') - F(x) - L(x'-x)}{\lvert x' - x\rvert} \right) \frac{\lvert x' - x\rvert}{\lvert y' - y\rvert} .

The last factor is at most $2$ , by $(3)$ . The vector in brackets tends to $0$ as $x' \to x$ , because $F$ is differentiable at $x$ with derivative $L$ ; and $x' \to x$ as $y' \to y$ , again by $(3)$ . Since the linear map $L^{-1}$ is bounded, the whole right-hand side tends to $0$ as $y' \to y$ . So $G$ is differentiable at $y$ with $DG(y) = L^{-1}$ .

Step 6: $G$ is of class $C^k$ . By step 5, $DG = \iota \circ DF \circ G$ , where $\iota(B) = B^{-1}$ on the open set of invertible matrices. Cramer's rule expresses each entry of $B^{-1}$ as a polynomial in the entries of $B$ divided by $\det B \neq 0$ , so $\iota$ is $C^\infty$ . Now induct on the order. $G$ is continuous (step 4) and $DF$ is continuous ( $F$ is $C^1$ ), so $DG$ is continuous, that is $G$ is $C^1$ . Suppose $G$ is $C^j$ with $1 \le j \le k-1$ ; then $DF$ is $C^{k-1}$ and hence $C^j$ , so the composition $DG = \iota\circ DF\circ G$ is $C^j$ , which says that $G$ is $C^{j+1}$ . For finite $k$ this ends at $G \in C^k$ ; for $k = \infty$ it gives $G \in C^j$ for every finite $j$ , that is $G \in C^\infty$ . This proves part 1.

Part 2. On $\Omega_0$ we have $G \circ F = \operatorname{id}$ with both maps differentiable, so the chain rule gives $DG(F(x))\,DF(x) = I$ for every $x \in \Omega_0$ . Hence $DF(x)$ is invertible and $DG(y) = [DF(x)]^{-1}$ at $y = F(x)$ ; taking determinants, $\det J_G(y)\,\det J_F(x) = 1$ . $\square$

Lee's proof is the model for steps 1–6. He states the theorem for smooth maps, normalises in two stages and quotes a general Lipschitz estimate for the bound $(1)$ ; here the normalisation is a single composition, $(1)$ is derived on the spot, the one instance of his contraction lemma (Lemma C.35) that is needed is proved in place rather than quoted, and step 6 keeps track of a finite differentiability class.

Theorem(Local diffeomorphism and openness)AdaptedAdapted from Lee (2013), App. C, Cor. C.36, p. 660

Let $\Omega \subseteq \mathbb{R}^n$ be open and let $F : \Omega \to \mathbb{R}^n$ be of class $C^k$ , $1 \le k \le \infty$ , with $\det J_F(x) \neq 0$ for every $x \in \Omega$ . Then

$F$ is an open map: $F(W)$ is open in $\mathbb{R}^n$ for every open $W \subseteq \Omega$ ; in particular $F(\Omega)$ is open;
if in addition $F$ is injective on $\Omega$ , then $F : \Omega \to F(\Omega)$ is a $C^k$ diffeomorphism.

Idea of the proof(Local diffeomorphism and openness)Original

The hypothesis is the one of the inverse function theorem, holding at every point at once. So around the image of each point there is a whole open set of image points, which is openness; and smoothness of an inverse is a local property, so once injectivity makes the local inverses fit together into one map, that map inherits their smoothness.

Proof(Local diffeomorphism and openness)AdaptedAdapted from Lee (2013), App. C, Cor. C.36 and its proof, p. 660

Openness. Let $W \subseteq \Omega$ be open and let $y \in F(W)$ , say $y = F(x)$ with $x \in W$ . The restriction $F|_W$ is $C^k$ on the open set $W$ and $\det J_F(x) \ne 0$ , so the inverse function theorem applied to $F|_W$ at $x$ produces open sets $\Omega_0 \subseteq W$ containing $x$ and $V_0$ containing $y$ with $F(\Omega_0) = V_0$ . Then $y \in V_0 \subseteq F(W)$ with $V_0$ open, and as $y$ was arbitrary, $F(W)$ is open. Taking $W = \Omega$ gives the last clause.
Global inverse. Suppose $F$ is also injective on $\Omega$ . Then $F : \Omega \to F(\Omega)$ is a bijection, so a map $G : F(\Omega) \to \Omega$ inverting it exists for set-theoretic reasons, and $F(\Omega)$ is open by part 1. Fix $y = F(x) \in F(\Omega)$ and take $\Omega_0 \ni x$ , $V_0 \ni y$ from the inverse function theorem at $x$ , with $(F|_{\Omega_0})^{-1} : V_0 \to \Omega_0$ of class $C^k$ . For $z \in V_0$ , both $G(z)$ and $(F|_{\Omega_0})^{-1}(z)$ are preimages of $z$ under $F$ , and $F$ is injective on all of $\Omega$ , so the two agree. Hence $G$ coincides with a $C^k$ map on a neighbourhood of each point of $F(\Omega)$ , so $G$ is $C^k$ , and $F : \Omega \to F(\Omega)$ is a $C^k$ diffeomorphism. $\square$

Injectivity is used only in part 2, and only to identify the local inverses with each other; nothing in part 1 needs it.

Theorem(Implicit function theorem)AdaptedAdapted from Lee (2013), App. C, Thm C.40, pp. 661–662

Write points of $\mathbb{R}^n \times \mathbb{R}^p$ as $(x,y)$ with $x = (x_1, \dots, x_n)$ and $y = (y_1, \dots, y_p)$ . Let $U \subseteq \mathbb{R}^n \times \mathbb{R}^p$ be open, let $\Phi : U \to \mathbb{R}^p$ be of class $C^k$ for some $k$ with $1 \le k \le \infty$ , let $(a,b) \in U$ and put $c = \Phi(a,b)$ . Suppose the $p \times p$ partial Jacobian in the $y$ variables,

\frac{\partial \Phi}{\partial y}(a,b) = \left( \frac{\partial \Phi_i}{\partial y_j}(a,b) \right)_{i,j \le p},

is invertible. Then there exist open sets $V_0 \subseteq \mathbb{R}^n$ with $a \in V_0$ and $W_0 \subseteq \mathbb{R}^p$ with $b \in W_0$ , with $V_0 \times W_0 \subseteq U$ , and a map $g : V_0 \to W_0$ of class $C^k$ , such that

\bigl\{ (x,y) \in V_0 \times W_0 \;:\; \Phi(x,y) = c \bigr\} = \bigl\{ (x, g(x)) \;:\; x \in V_0 \bigr\} ;

that is, for $(x,y) \in V_0 \times W_0$ one has $\Phi(x,y) = c$ if and only if $y = g(x)$ . In particular $g(a) = b$ , and $g$ is the only map $V_0 \to W_0$ with this property.

Moreover $V_0$ may be taken small enough that $\partial\Phi/\partial y\bigl(x, g(x)\bigr)$ is invertible at every $x \in V_0$ — replacing $V_0$ by a smaller open set containing $a$ leaves the displayed equality of sets intact, since it merely restricts both sides. With $V_0$ so chosen, for every $x \in V_0$

\frac{\partial \Phi}{\partial x}\bigl(x, g(x)\bigr) + \frac{\partial \Phi}{\partial y}\bigl(x, g(x)\bigr)\, Dg(x) = 0, \qquad Dg(x) = - \left[ \frac{\partial \Phi}{\partial y}\bigl(x, g(x)\bigr) \right]^{-1} \frac{\partial \Phi}{\partial x}\bigl(x, g(x)\bigr),

where $\partial \Phi/\partial x$ is the $p \times n$ block of derivatives in the $x$ variables. For $n = p = 1$ this is the familiar $g'(x) = -\Phi_x/\Phi_y$ , both evaluated at $(x, g(x))$ . What is asserted is invertibility along the graph of $g$ , which is all that is used; invertibility throughout $V_0 \times W_0$ is not claimed.

Of this, Lee's Thm C.40 supplies the two neighbourhoods, the $C^k$ map $g$ and the equality of the level set with its graph. The uniqueness of $g$ , the refinement of $V_0$ just described and the derivative formula are ours, obtained from the chain rule in steps 6–8 of the proof below.

Idea of the proof(Implicit function theorem)AdaptedAdapted from Lee (2013), App. C, proof of Thm C.40, pp. 661–662

Make the equation part of a change of variables. The map $\Psi(x,y) = (x, \Phi(x,y))$ leaves the $x$ coordinates untouched and replaces $y$ by the value of $\Phi$ , so its derivative is block triangular with diagonal blocks $I_n$ and $\partial\Phi/\partial y$ — invertible exactly under the hypothesis. The inverse function theorem then makes $\Psi$ a $C^k$ diffeomorphism of a neighbourhood of $(a,b)$ . Its inverse must return $x$ in the first slot, so it has the form $(x,z) \mapsto (x, \beta(x,z))$ ; freezing the second argument at the level $c$ turns $\beta$ into the wanted function $g$ , and the level set becomes its graph.

Proof(Implicit function theorem)AdaptedAdapted from Lee (2013), App. C, proof of Thm C.40, pp. 661–662

A map with invertible derivative. Define $\Psi : U \to \mathbb{R}^n\times\mathbb{R}^p$ by $\Psi(x,y) = \bigl(x, \Phi(x,y)\bigr)$ . It is of class $C^k$ , since $\Phi$ is and the first block is the identity, and in the splitting of the coordinates its derivative is $D\Psi(x,y) = \begin{pmatrix} I_n & 0 \\ \partial\Phi/\partial x & \partial\Phi/\partial y \end{pmatrix} ,$ block lower triangular, so $\det D\Psi(x,y) = \det\bigl(\partial\Phi/\partial y\,(x,y)\bigr)$ . At $(a,b)$ that determinant is non-zero: this is the hypothesis, and it is the only place it is used.
Invert it locally, on a product. By the inverse function theorem there are open sets $\Omega_1 \subseteq U$ containing $(a,b)$ and $Y_1$ containing $\Psi(a,b) = (a,c)$ with $\Psi|_{\Omega_1} : \Omega_1 \to Y_1$ a $C^k$ diffeomorphism. Choose open $V \ni a$ in $\mathbb{R}^n$ and $W_0 \ni b$ in $\mathbb{R}^p$ with $V\times W_0 \subseteq \Omega_1$ , and set $\Omega_0 = V\times W_0$ and $Y_0 = \Psi(\Omega_0)$ . Since $\Psi|_{\Omega_1}$ is a homeomorphism onto the open set $Y_1$ , the set $Y_0$ is open, it contains $(a,c)$ , and $\Psi|_{\Omega_0} : \Omega_0 \to Y_0$ is a $C^k$ diffeomorphism. Write its inverse as $\bigl(\Psi|_{\Omega_0}\bigr)^{-1}(x,z) = \bigl(\alpha(x,z),\, \beta(x,z)\bigr), \qquad (x,z) \in Y_0,$ with $\alpha$ and $\beta$ of class $C^k$ (they are the two blocks of components of a $C^k$ map).
The inverse returns $x$ unchanged. For $(x,z) \in Y_0$ , $(x,z) = \Psi\bigl(\alpha(x,z), \beta(x,z)\bigr) = \bigl(\alpha(x,z),\, \Phi(\alpha(x,z), \beta(x,z))\bigr) .$ Comparing first components gives $\alpha(x,z) = x$ ; substituting that and comparing second components gives $\Phi\bigl(x, \beta(x,z)\bigr) = z \qquad \text{for all } (x,z) \in Y_0. \tag{$*$}$
Definition of $g$ . Put $V_0 = \{x \in V : (x,c) \in Y_0\}$ , which is open because $x \mapsto (x,c)$ is continuous and $Y_0$ is open, and which contains $a$ because $(a,c) \in Y_0$ . Define $g : V_0 \to \mathbb{R}^p$ by $g(x) = \beta(x,c)$ ; it is $C^k$ , being $\beta$ composed with an affine map. For $x \in V_0$ the point $(x, g(x))$ is $(\Psi|_{\Omega_0})^{-1}(x,c) \in \Omega_0 = V\times W_0$ , so $g(x) \in W_0$ : thus $g$ maps $V_0$ into $W_0$ , as claimed. Also $V_0\times W_0 \subseteq V\times W_0 \subseteq U$ , and $g(a) = b$ because $(\Psi|_{\Omega_0})^{-1}(a,c) = (a,b)$ .
The level set is the graph of $g$ . Let $(x,y) \in V_0\times W_0$ . If $y = g(x)$ then $\Phi(x,y) = \Phi(x,\beta(x,c)) = c$ by $(*)$ . Conversely, suppose $\Phi(x,y) = c$ . Then $\Psi(x,y) = (x,c)$ , and $(x,y) \in \Omega_0$ , so applying $(\Psi|_{\Omega_0})^{-1}$ to both sides gives $(x,y) = (x, \beta(x,c)) = (x, g(x))$ , that is $y = g(x)$ . The two sets displayed in the statement are therefore equal.
Uniqueness. Any $\tilde g : V_0 \to W_0$ whose graph is that same set has $\tilde g(x) = g(x)$ for every $x \in V_0$ , since the set meets each $\{x\}\times W_0$ in the single point $(x,g(x))$ .
Shrinking $V_0$ to make $\partial\Phi/\partial y$ invertible along the graph. Consider $d(x) = \det\bigl(\partial\Phi/\partial y(x, g(x))\bigr)$ on $V_0$ . Each entry of $\partial\Phi/\partial y$ is continuous, because $\Phi$ is of class $C^k$ with $k \ge 1$ ; $g$ is continuous by step 4; and the determinant is a polynomial in the entries. So $d$ is continuous. Its value at $x = a$ is $\det\bigl(\partial\Phi/\partial y(a,b)\bigr) \neq 0$ — the hypothesis again. Hence $V_0' = \{x \in V_0 : d(x) \neq 0\} = d^{-1}\bigl(\mathbb{R}\setminus\{0\}\bigr)$ is open and contains $a$ . Replace $V_0$ by $V_0'$ : the set equality of step 5 restricts to both sides, $g$ restricted to $V_0'$ is still $C^k$ with values in $W_0$ , and step 6 applies verbatim. This secures invertibility only along the graph; to have it at every point of a product neighbourhood one would have to shrink $W_0$ as well and then re-shrink $V_0$ to keep $g$ inside the smaller $W_0$ , and nothing below needs that.
The derivative of $g$ . On $V_0$ the identity $\Phi\bigl(x, g(x)\bigr) = c$ holds, by step 5, with $\Phi$ and $g$ both of class $C^k$ , $k \ge 1$ . The map $x \mapsto (x, g(x))$ has for derivative the $(n+p)\times n$ matrix whose top block is $I_n$ and whose bottom block is $Dg(x)$ . The composite $x \mapsto \Phi(x, g(x))$ is constant, so its derivative vanishes, and the chain rule computes that derivative as $\partial\Phi/\partial x$ times the top block plus $\partial\Phi/\partial y$ times the bottom block: $\frac{\partial \Phi}{\partial x}\bigl(x, g(x)\bigr) + \frac{\partial \Phi}{\partial y}\bigl(x, g(x)\bigr)\, Dg(x) = 0 .$ By step 7 the $p \times p$ matrix $\partial\Phi/\partial y(x,g(x))$ is invertible for every $x \in V_0$ , so multiplying on the left by its inverse and rearranging gives $Dg(x) = -\bigl[\partial\Phi/\partial y(x,g(x))\bigr]^{-1}\,\partial\Phi/\partial x(x,g(x))$ , the formula of the statement. $\square$

Steps 1–5 are Lee's route: the map $\Psi$ , its block-triangular derivative, the local inverse and the identification of the level set with a graph. Steps 6–8 are ours; Thm C.40 states none of the uniqueness of $g$ , the refinement of $V_0$ , or the derivative formula.

NoteOriginal

What it means. Both theorems say the same thing twice: where the linear approximation of a map is invertible, the map itself behaves like that linear approximation — but only nearby. The inverse function theorem reads this as " $F$ can be undone locally"; the implicit function theorem reads it as "the equation $\Phi(x,y) = c$ can be solved locally for $y$ , and the solution inherits the smoothness of $\Phi$ ". Neither gives a formula for the solution, only its existence, its uniqueness in a box, and its derivative. Each follows from the other: Lee's route to the implicit function theorem applies the inverse function theorem to $(x,y) \mapsto (x, \Phi(x,y))$ , whose Jacobian is block triangular with diagonal blocks $I_n$ and $\partial \Phi/\partial y$ , and the inverse function theorem itself is usually proved by a contraction-mapping argument.

RemarkOriginal

Where the hypotheses are used, and what fails without them.

Local, never global. Let $F(x,y) = (e^x \cos y,\; e^x \sin y)$ on $\Omega = \mathbb{R}^2$ . Its Jacobian determinant is $e^{2x}$ , non-zero everywhere, so $F$ is a local diffeomorphism at every point and an open map; but $F(x, y + 2\pi) = F(x,y)$ , so $F$ is not injective. A non-vanishing Jacobian determinant on a connected open set therefore does not give global invertibility; injectivity is an extra hypothesis, as in part 2 of the local-diffeomorphism statement.
Class $C^1$ cannot be weakened to differentiability. Let $f(0) = 0$ and $f(x) = x + 2x^2 \sin(1/x)$ for $x \neq 0$ . Then $f$ is differentiable on $\mathbb{R}$ with $f'(0) = \lim_{h\to 0}(1 + 2h\sin(1/h)) = 1 \ne 0$ , while for $x \neq 0$ $f'(x) = 1 + 4x \sin(1/x) - 2 \cos(1/x),$ which equals $-1$ at every $x = 1/(2\pi j)$ and $3$ at every $x = 1/((2j+1)\pi)$ , $j \in \mathbb{N}$ . So $f'$ takes both signs in every neighbourhood of $0$ . That destroys injectivity near $0$ , by the following argument. Given a neighbourhood of $0$ , choose $j$ so large that it contains the interval $I = [l, r]$ with $l = 1/((2j+1)\pi)$ and $r = 1/(2\pi j)$ , where indeed $l < r$ because $2j+1 > 2j$ . Then $f'(l) = 3 > 0$ and $f'(r) = -1 < 0$ , so $f$ exceeds $f(l)$ just to the right of $l$ and exceeds $f(r)$ just to the left of $r$ : neither endpoint is where $f$ attains its maximum over the compact interval $I$ . That maximum is therefore attained at an interior point $\xi$ , with $f(\xi) > f(l)$ and $f(\xi) > f(r)$ . Pick $v$ with $\max\{f(l), f(r)\} < v < f(\xi)$ . The intermediate value theorem, applied to the continuous $f$ on $[l,\xi]$ and then on $[\xi, r]$ , produces one point of each interval at which $f = v$ ; neither of them is $\xi$ , since $f(\xi) > v$ , so one lies strictly left of $\xi$ and the other strictly right of it, and the two are distinct. So $f$ is injective on no neighbourhood of $0$ : the conclusion fails although $Df(0)$ is invertible. What is missing is the continuity of $f'$ at $0$ .
The invertible block decides which variables are solved for. For $\Phi(x,y) = x^2 + y^2 - 1$ on $\mathbb{R}^2$ , $\partial \Phi/\partial y = 2y$ is invertible at every point of the unit circle except $(\pm 1, 0)$ ; near those two points the circle is not the graph of a function of $x$ , but $\partial \Phi/\partial x = 2x$ is invertible there, so it is the graph of a function of $y$ . The hypothesis is on one block of the Jacobian, not on the whole of it.
Existence is not a formula. For $\Phi(x,y) = y^5 + y - x$ one has $\partial \Phi/\partial y = 5y^4 + 1 \ge 1$ everywhere, so near each solution the theorem gives $y$ as a smooth function of $x$ , with $g'(x) = 1/(5g(x)^4+1)$ . That these local solutions are pieces of a single smooth function on the whole of $\mathbb{R}$ is a further, non-local fact, which the theorem cannot supply: it holds because $y \mapsto y^5+y$ has derivative $\ge 1$ , so is strictly increasing, and tends to $\pm\infty$ , so is a bijection of $\mathbb{R}$ onto $\mathbb{R}$ — whence every $x$ has exactly one solution $g(x)$ , and the local solutions, being solutions, agree with it. But $g$ is given by no formula in radicals, since $y^5 + y - x$ has Galois group $S_5$ over $\mathbb{Q}(x)$ and so is not solvable by radicals. The derivative formula is still fully usable, which is the normal way these theorems are applied.
Smoothness class in the sources. Lee states both theorems for smooth maps; Folland uses and states the $C^1$ version, since that is what the change-of-variables theorem needs. The versions above are for $C^k$ with $1 \le k \le \infty$ , covering both: the inverse is obtained first as a continuous map, then shown differentiable, and its derivative is a composition of $F^{-1}$ , $x \mapsto DF(x)$ and matrix inversion, so each extra degree of smoothness of $F$ gives one more for $F^{-1}$ .
Why the Jacobian determinant. The hypothesis is the invertibility of the linear map $DF(a)$ , which exists as a total derivative because $F$ is $C^1$ (partial derivatives and the chain rule). Arfken's reciprocal relation $\partial(x_1,x_2,\dots)/\partial(u_1,u_2,\dots) = [\partial(u_1,u_2,\dots)/\partial(x_1,x_2,\dots)]^{-1}$ between the Jacobian of a change of variables and that of its inverse (Eq. (4.75), p. 230) is exactly the determinant identity of part 2 of the inverse function theorem.

Example(Solving a system implicitly, and inverting a map locally)Original

(a) Two equations, two unknowns, one parameter. Take $n = 1$ , $p = 2$ and

\Phi(x, y_1, y_2) = \bigl( y_1^3 + y_2 - 2x, \;\; y_1 + y_2^2 - 2x \bigr)

on $U = \mathbb{R}^3$ , a polynomial map and so of class $C^\infty$ . At $(a,b) = \bigl(1,(1,1)\bigr)$ both components vanish: $1 + 1 - 2 = 0$ twice. So $c = \Phi(a,b) = (0,0)$ .

Check the hypothesis.

\frac{\partial\Phi}{\partial y} = \begin{pmatrix} 3y_1^2 & 1 \\ 1 & 2y_2 \end{pmatrix} \quad\Longrightarrow\quad \frac{\partial\Phi}{\partial y}(a,b) = \begin{pmatrix} 3 & 1 \\ 1 & 2\end{pmatrix}, \qquad \det = 6 - 1 = 5 \neq 0 .

The theorem therefore gives a neighbourhood $V_0$ of $x = 1$ and a $C^\infty$ map $g = (g_1,g_2)$ on it with $g(1) = (1,1)$ , such that for $(x,y)$ in a box about $\bigl(1,(1,1)\bigr)$ the system $y_1^3 + y_2 = 2x$ , $y_1 + y_2^2 = 2x$ holds exactly when $y = g(x)$ .

Its derivative. With $\partial\Phi/\partial x = (-2,-2)^{\mathsf T}$ and

\begin{pmatrix} 3 & 1\\ 1 & 2\end{pmatrix}^{-1} = \frac{1}{5}\begin{pmatrix} 2 & -1 \\ -1 & 3\end{pmatrix},

the formula $Dg = -[\partial\Phi/\partial y]^{-1}\,[\partial\Phi/\partial x]$ gives

Dg(1) = -\frac{1}{5}\begin{pmatrix} 2 & -1\\ -1 & 3\end{pmatrix}\begin{pmatrix}-2\\-2\end{pmatrix} = \frac{1}{5}\begin{pmatrix}2\\4\end{pmatrix} = \begin{pmatrix} 2/5 \\ 4/5\end{pmatrix},

so $g_1'(1) = 0.4$ and $g_2'(1) = 0.8$ .

Numerical check. Solving the system at $x = 1.01$ by Newton's method from $(1,1)$ gives $y \approx (1.0039935, 1.0079715)$ , and the difference quotients $(y_i - 1)/0.01$ are $0.39935$ and $0.79715$ — the predicted slopes $0.4$ and $0.8$ to the accuracy expected of a first-order approximation over a step of $0.01$ .

(b) A map that is locally but not globally invertible. Let

F(x,y) = \bigl(x^2 - y^2,\; xy\bigr), \qquad J_F(x,y) = \begin{pmatrix} 2x & -2y \\ y & x\end{pmatrix}, \qquad \det J_F = 2x^2 + 2y^2 .

The determinant is non-zero except at the origin, so on $\Omega = \mathbb{R}^2\setminus\{(0,0)\}$ the local-diffeomorphism statement applies: $F$ is an open map there and a $C^\infty$ diffeomorphism of any open set on which it happens to be injective. It is not injective on $\Omega$ , because $F(-x,-y) = F(x,y)$ .

At $(2,1)$ : $F(2,1) = (3,2)$ and $\det J_F = 10$ , so there are neighbourhoods $\Omega_0 \ni (2,1)$ and $V_0 \ni (3,2)$ with $F|_{\Omega_0}$ a diffeomorphism onto $V_0$ , and its inverse $G$ has

DG(3,2) = \begin{pmatrix} 4 & -2 \\ 1 & 2\end{pmatrix}^{-1} = \frac{1}{10}\begin{pmatrix}2 & 2\\ -1 & 4\end{pmatrix} = \begin{pmatrix} 0.2 & 0.2 \\ -0.1 & 0.4 \end{pmatrix},

and $\det J_G(3,2) = 1/10$ . Solving $F(x,y) = (3.1, 2)$ numerically from $(2,1)$ gives $(2.019979, 0.990109)$ : the increments divided by $0.1$ are $0.19979$ and $-0.09891$ , against the first column $(0.2, -0.1)$ of $DG(3,2)$ .

Application(Thermodynamic derivatives and the cyclic relation)Original

An equation of state ties the state variables of a fluid together by one equation, $\Phi(p, V, T) = 0$ — for one mole of a van der Waals gas,

\Phi(p,V,T) = \Bigl(p + \frac{a}{V^2}\Bigr)(V - b) - RT .

Thermodynamics then treats any two of $p, V, T$ as independent and the third as a function of them, and writes partial derivatives such as $(\partial p/\partial V)_T$ . The implicit function theorem is what licenses this: wherever $\Phi$ is $C^1$ and $\partial\Phi/\partial p \neq 0$ , the equation determines $p$ locally as a $C^1$ function of $(V,T)$ , and

\Bigl(\frac{\partial p}{\partial V}\Bigr)_T = -\frac{\Phi_V}{\Phi_p}, \qquad \Bigl(\frac{\partial V}{\partial T}\Bigr)_p = -\frac{\Phi_T}{\Phi_V}, \qquad \Bigl(\frac{\partial T}{\partial p}\Bigr)_V = -\frac{\Phi_p}{\Phi_T},

each from the one-dimensional derivative formula with the other two variables playing the role of parameters. Two consequences used constantly in the subject come out at once. Solving instead for $V$ as a function of $(p,T)$ gives $(\partial V/\partial p)_T = -\Phi_p/\Phi_V$ , the reciprocal of the first expression above: this is the reciprocal relation $(\partial p/\partial V)_T = 1/(\partial V/\partial p)_T$ , the one-dimensional case of the determinant identity of the inverse function theorem. Multiplying the three displayed expressions, every factor cancels against another while all three minus signs survive, giving the cyclic relation

\Bigl(\frac{\partial p}{\partial V}\Bigr)_T \Bigl(\frac{\partial V}{\partial T}\Bigr)_p \Bigl(\frac{\partial T}{\partial p}\Bigr)_V = -1 .

It holds for any equation of state whose three partial derivatives are non-zero at the state in question; the minus sign, which surprises people who expect the three factors to cancel to $+1$ , is simply the product of the three minus signs in the derivative formula.

For the van der Waals gas, $\Phi_p = V - b$ , $\Phi_T = -R$ and $\Phi_V = p - a/V^2 + 2ab/V^3$ , so for instance

\Bigl(\frac{\partial p}{\partial T}\Bigr)_V = -\frac{\Phi_T}{\Phi_p} = \frac{R}{V-b} ,

with no need to solve the equation of state for anything. Where the hypothesis fails, the physics changes. Since $\Phi_p = V - b \neq 0$ , the first of the two conditions that fix the critical point, $(\partial p/\partial V)_T = -\Phi_V/\Phi_p = 0$ , says exactly that $\Phi_V = 0$ ; and that is the failure of the hypothesis under which $V$ is a function of $p$ at fixed $T$ . So wherever $\Phi_V = 0$ the isotherm is flat and the isothermal compressibility $-V^{-1}(\partial V/\partial p)_T$ is infinite: the single-phase description has lost its local mechanical stability. That locus is the spinodal, and it is not where the two-phase region begins. The boundary of the two-phase region is the binodal, fixed by the Maxwell equal-area construction — a global condition on the isotherm, which no derivative of $\Phi$ and no local theorem can deliver. For every $T < T_c$ the spinodal lies strictly inside the binodal; between the two the homogeneous states are metastable rather than absent, with $(\partial p/\partial V)_T < 0$ and so $\Phi_V > 0$ , and the phases have in fact already separated there. The two loci meet only at the critical point, which is why $\Phi_V = 0$ is nevertheless exactly the right condition for the computation that follows. The second condition, $(\partial^2 p/\partial V^2)_T = 0$ , is a genuinely further one, not a restatement of the first; both are needed, and solving the pair for the van der Waals gas gives the familiar

V_c = 3b, \qquad T_c = \frac{8a}{27Rb}, \qquad p_c = \frac{a}{27b^2}, \qquad \frac{p_cV_c}{RT_c} = \frac{3}{8},

the last number being the same for every van der Waals fluid, which is why it is a test of the model.

Method

Lagrange multipliers

Adapted

Add to PDF

If a $C^1$ function has a local extremum, among the points of the common zero set of $m < n$ constraints, at a point where the $m$ constraint gradients are linearly independent, then its gradient there is a unique linear combination of them; this gives the standard recipe for constrained extrema together with the points at which the recipe can miss them.

Prerequisites: Implicit and inverse function theorems, Taylor expansion and critical points in several variables, Compactness

Definition(Constraint set, regular point, Lagrangian)AdaptedAdapted from Lee (2013), §5, p. 107 and Problem 11-11, p. 301; Riley, Hobson & Bence (2006), §5.9, pp. 168–169

Let $U \subseteq \mathbb{R}^n$ be open, let $1 \le m < n$ , and let $g = (g_1, \dots, g_m) : U \to \mathbb{R}^m$ be of class $C^1$ .

The constraint set is the level set $M = \{x \in U : g(x) = 0\}$ . (A constraint written $g_i(x) = c_i$ is brought to this form by replacing $g_i$ with $g_i - c_i$ , which changes no derivative.)
A point $x_* \in M$ is a regular point of the constraint if the $m \times n$ Jacobian matrix $J_g(x_*)$ has rank $m$ ; equivalently, if the gradients $\nabla g_1(x_*), \dots, \nabla g_m(x_*)$ are linearly independent. If every point of $M$ is regular, $M$ is a regular level set of $g$ and $g$ is a defining function for it.
The tangent space to $M$ at a regular point $x_*$ is $T_{x_*}M = \ker J_g(x_*) = \bigl\{ v \in \mathbb{R}^n : \nabla g_i(x_*) \cdot v = 0 \text{ for } i = 1, \dots, m \bigr\},$ a subspace of dimension $n - m$ ; its orthogonal complement, the normal space, is spanned by the $\nabla g_i(x_*)$ .
$f : U \to \mathbb{R}$ has a local minimum on $M$ at $x_* \in M$ if there is $r > 0$ with $f(x) \ge f(x_*)$ for every $x \in M \cap B_r(x_*)$ ; a local maximum on $M$ reverses the inequality, and a local extremum on $M$ is either. Only the values of $f$ on $M$ are compared, so $x_*$ need not be a critical point of $f$ on $U$ .
The Lagrangian of the problem is the function on $U \times \mathbb{R}^m$ $L(x, \lambda) = f(x) - \sum_{i=1}^{m} \lambda_i\, g_i(x),$ and $\lambda_1, \dots, \lambda_m$ are its Lagrange multipliers (also called undetermined multipliers).

Method(First-order necessary condition)AdaptedAdapted from Lee (2013), Problem 11-11, p. 301; Riley, Hobson & Bence (2006), §5.9, Eqs. (5.27)–(5.28), p. 168 and Eq. (5.31), p. 169; Arfken, Weber & Harris (2013), §22.3, Eq. (22.67), p. 1109

Let $U \subseteq \mathbb{R}^n$ be open, let $1 \le m < n$ , let $f : U \to \mathbb{R}$ and $g = (g_1, \dots, g_m) : U \to \mathbb{R}^m$ be of class $C^1$ , and let $M = \{x \in U : g(x) = 0\}$ .

Suppose $f$ has a local extremum on $M$ at a point $x_* \in M$ which is a regular point of the constraint, that is $\operatorname{rank} J_g(x_*) = m$ . Then there is exactly one $\lambda = (\lambda_1, \dots, \lambda_m) \in \mathbb{R}^m$ with

\nabla f(x_*) = \sum_{i=1}^{m} \lambda_i\, \nabla g_i(x_*),

equivalently $\nabla_x L(x_*, \lambda) = 0$ . Since $\nabla_\lambda L(x_*,\lambda) = -g(x_*) = 0$ holds automatically for $x_* \in M$ , the pair $(x_*, \lambda)$ is a critical point of $L$ on $U \times \mathbb{R}^m$ , and the $n + m$ equations

\frac{\partial f}{\partial x_j}(x) - \sum_{i=1}^{m} \lambda_i \frac{\partial g_i}{\partial x_j}(x) = 0 \quad (j = 1, \dots, n), \qquad g_i(x) = 0 \quad (i = 1, \dots, m)

in the $n + m$ unknowns $(x, \lambda)$ hold at $(x_*, \lambda)$ . Geometrically, $\nabla f(x_*)$ lies in the normal space of $M$ at $x_*$ , that is $\nabla f(x_*) \cdot v = 0$ for every $v \in T_{x_*}M$ : the derivative of $f$ vanishes in every direction tangent to the constraint set. The condition is necessary only: a solution of the displayed system need be neither a maximum nor a minimum of $f$ on $M$ .

Idea of the proof(First-order necessary condition)AdaptedAdapted from Arfken, Weber & Harris (2013), §22.3, Eqs. (22.63)–(22.66), p. 1108

Regularity says that $m$ of the $n$ coordinates can be eliminated. Call those $m$ coordinates $y$ and the remaining $n-m$ of them $u$ , chosen so that $\partial g/\partial y$ is invertible; the implicit function theorem then presents $M$ near $x_*$ as the graph $y = h(u)$ of a $C^1$ map, and $f$ restricted to $M$ becomes an ordinary unconstrained function $F(u) = f(u,h(u))$ of the free variables, whose gradient must vanish. The chain rule turns that into one relation between the $u$ - and the $y$ -derivatives of $f$ and of $g$ . Define $\lambda$ from the $y$ -block alone, which is where $\partial g/\partial y$ can be inverted; the relation then says precisely that the $u$ -block matches as well, and the two blocks together are $\nabla f(x_*) = \sum_i\lambda_i\nabla g_i(x_*)$ .

Proof(First-order necessary condition)AdaptedAdapted from Arfken, Weber & Harris (2013), §22.3, Eqs. (22.63)–(22.66), p. 1108; Evans (2010), §8.4.1, proof of Thm 2, pp. 460–462

Split the coordinates. $J_g(x_*)$ is $m\times n$ of rank $m$ , so $m$ of its columns are linearly independent. Relabel the coordinates of $\mathbb{R}^n$ so that these are the last $m$ , writing $x = (u,y)$ with $u\in\mathbb{R}^{n-m}$ and $y\in\mathbb{R}^m$ , and $x_* = (u_*,y_*)$ . A relabelling permutes the entries of every gradient in the same way, so a conclusion of the form $\nabla f = \sum_i\lambda_i\nabla g_i$ is unaffected by it. By construction $\partial g/\partial y\,(x_*)$ is an invertible $m\times m$ matrix.
Eliminate $y$ . Apply the implicit function theorem to $g$ at $(u_*,y_*)$ with level $c = 0$ : there are open sets $V_0\ni u_*$ and $W_0 \ni y_*$ with $V_0\times W_0\subseteq U$ , and a $C^1$ map $h : V_0\to W_0$ with $h(u_*) = y_*$ , such that for $(u,y)\in V_0\times W_0$ $g(u,y) = 0 \iff y = h(u) .$
The restricted problem is unconstrained. Put $F(u) = f\bigl(u, h(u)\bigr)$ , a $C^1$ function on $V_0$ . Let $r>0$ be the radius in the hypothesis that $x_*$ is a local extremum of $f$ on $M$ , and shrink $V_0$ so that $(u,h(u)) \in B_r(x_*)$ for every $u\in V_0$ , which is possible because $u\mapsto (u,h(u))$ is continuous and maps $u_*$ to $x_*$ . Each such point lies in $M$ by step 2, so $F(u)\ge F(u_*)$ for all $u\in V_0$ in the minimum case and $F(u)\le F(u_*)$ in the maximum case. Hence $u_*$ is an unconstrained local extremum of $F$ at an interior point of its domain, so $\nabla F(u_*) = 0$ : each partial derivative is the derivative of a one-variable differentiable function with an interior local extremum, and so vanishes.
What the chain rule makes of that. Differentiating $F = f\circ(u\mapsto(u,h(u)))$ at $u_*$ with the chain rule gives, as an identity between $1\times(n-m)$ row vectors, $0 \;=\; \frac{\partial f}{\partial u}(x_*) + \frac{\partial f}{\partial y}(x_*)\, Dh(u_*) .$
What the implicit function theorem gives for $Dh$ . Its derivative formula, applied to $g$ , reads $Dh(u_*) = -\left[\frac{\partial g}{\partial y}(x_*)\right]^{-1}\frac{\partial g}{\partial u}(x_*) .$
Define the multipliers, then check both blocks. Put $\lambda^{\mathsf T} \;=\; \frac{\partial f}{\partial y}(x_*)\left[\frac{\partial g}{\partial y}(x_*)\right]^{-1} \;\in\; \mathbb{R}^{1\times m},$ which is Arfken's $\lambda = (\partial f/\partial z)\big/(\partial g/\partial z)$ with a matrix inverse in place of the division. The $y$ -block then holds by construction, $\frac{\partial f}{\partial y}(x_*) = \lambda^{\mathsf T}\frac{\partial g}{\partial y}(x_*) ,$ and substituting step 5 into step 4 gives the $u$ -block: $\frac{\partial f}{\partial u}(x_*) = -\frac{\partial f}{\partial y}(x_*)Dh(u_*) = \frac{\partial f}{\partial y}(x_*)\left[\frac{\partial g}{\partial y}(x_*)\right]^{-1}\frac{\partial g}{\partial u}(x_*) = \lambda^{\mathsf T}\frac{\partial g}{\partial u}(x_*) .$ Together the two blocks say $Df(x_*) = \lambda^{\mathsf T}J_g(x_*)$ , that is $\nabla f(x_*) = \sum_{i=1}^m \lambda_i\nabla g_i(x_*)$ ; undoing the relabelling of step 1 states it in the original coordinates.
Uniqueness, and the other forms of the conclusion. If $\sum_i \lambda_i\nabla g_i(x_*) = \sum_i\mu_i\nabla g_i(x_*)$ then $\sum_i(\lambda_i-\mu_i)\nabla g_i(x_*) = 0$ , and the rank hypothesis makes those gradients linearly independent, so $\mu = \lambda$ : the multiplier is unique. For $v\in T_{x_*}M = \ker J_g(x_*)$ we get $\nabla f(x_*)\cdot v = \lambda^{\mathsf T}J_g(x_*)v = 0$ , the geometric form. And $\nabla_xL(x_*,\lambda) = \nabla f(x_*) - \sum_i\lambda_i\nabla g_i(x_*) = 0$ while $\nabla_\lambda L(x_*,\lambda) = -g(x_*) = 0$ because $x_* \in M$ , so $(x_*,\lambda)$ is a critical point of $L$ on $U\times\mathbb{R}^m$ . $\square$

Arfken carries out exactly this elimination for $n = 3$ , $m = 1$ , and his Eq. (22.65) is step 6 in that case; Evans's Thm 2 follows the same implicit-function route in a function space. The $m$ -constraint bookkeeping and the relabelling of step 1 are ours.

NoteOriginal

What it means. On the constraint set only $n - m$ of the $n$ coordinates can be moved freely, so the right first-order condition is not that $\nabla f$ vanish but that it have no component along the directions one is still allowed to move in — and those are exactly the directions perpendicular to all the constraint gradients. The multipliers are the coefficients that express the leftover, normal part of $\nabla f$ . The reason the regularity hypothesis is needed, and the reason the result is true, is the implicit function theorem: where the constraint gradients are independent, $m$ of the variables can be eliminated as $C^1$ functions of the other $n-m$ , turning the problem into an unconstrained one in $n-m$ free variables.

Method(Constrained extrema by multipliers)AdaptedAdapted from Riley, Hobson & Bence (2006), §5.9, pp. 168–171; Arfken, Weber & Harris (2013), §22.3, pp. 1107–1109

To find the extrema of a $C^1$ function $f$ on $M = \{x \in U : g(x) = 0\}$ , with $g$ of class $C^1$ and $m < n$ constraints:

Write the constraints as equalities with zero on the right, $g_i(x) = 0$ , and count: $n$ variables, $m$ constraints, $m < n$ .
Form the Lagrangian $L(x,\lambda) = f(x) - \sum_i \lambda_i g_i(x)$ .
Solve the $n+m$ stationarity equations $\partial L/\partial x_j = 0$ for $j \le n$ together with $g_i(x) = 0$ for $i \le m$ , in the $n+m$ unknowns $x_1, \dots, x_n, \lambda_1, \dots, \lambda_m$ . Keep every branch: never divide by a quantity that could vanish without treating that case separately — the step at which Riley, Hobson and Bence's worked sphere problem would lose solutions.
List separately every point of $M$ where $\operatorname{rank} J_g < m$ . The theorem says nothing at such points, so they are extra candidates and must be examined on their own.

How it fails. Step 3 finds only regular candidates, so skipping step 4 can miss the answer outright (see the remark). The equations locate maxima, minima and saddle points alike — Arfken makes this point on p. 1109 — so they only produce candidates; and on a non-compact $M$ , or when $f$ is unbounded on $M$ , there may be no extremum at all although step 3 has solutions. The alternative of substituting the constraint to eliminate a variable, which Arfken contrasts the method with, is equally valid whenever the elimination can actually be carried out; multipliers are preferred when it cannot be, or when the symmetry of the variables should be kept.

Method(Deciding which candidates are extrema)Original

The multiplier equations do not say which candidates are extrema. Let $C \subseteq M$ be the set of candidates produced by steps 3 and 4 above.

If $M$ is compact (closed and bounded) and $f$ is continuous on $M$ , then $f$ attains a maximum and a minimum on $M$ , and by the theorem each is attained at a point of $C$ . So the largest and the smallest value of $f$ on $C$ are the maximum and the minimum of $f$ on $M$ ; when $C$ is finite this is a comparison of finitely many numbers, and when it is not, it is the easier problem of extremising $f$ on $C$ . No second-order test is needed. The multipliers may be discarded once $x_*$ is known.
Otherwise test a candidate directly. With $f$ and $g$ of class $C^2$ , $x_*$ a regular point, $\lambda$ the multiplier belonging to $x_*$ , and $H \;=\; H_f(x_*) - \sum_{i=1}^{m} \lambda_i\, H_{g_i}(x_*)$ the Hessian of $x \mapsto L(x,\lambda)$ at $x_*$ : if $h^{\mathsf T} H h > 0$ for every non-zero $h \in T_{x_*}M$ then $x_*$ is a strict local minimum of $f$ on $M$ , and if $h^{\mathsf T} H h < 0$ for every such $h$ then a strict local maximum. $H$ need only be definite on the tangent space, not on all of $\mathbb{R}^n$ ; if it is indefinite on the tangent space, $x_*$ is neither.

Idea of the proof(Why steps 5 and 6 work)Original

Step 5 needs no calculus. On a compact constraint set a continuous $f$ attains its largest and smallest values somewhere, and wherever that is, either the theorem applies there or the rank condition fails there — so the point is already on the candidate list, and comparing values on the list is the whole of the matter. Step 6 is the unconstrained second-derivative test in disguise: parametrise $M$ near $x_*$ by the free variables, and the Hessian of $f$ pulled back to those parameters turns out to be exactly $H$ restricted to the tangent space, because the first-order condition annihilates the terms that involve the curvature of the parametrisation. That is what the Hessian of the Lagrangian, rather than of $f$ , is for.

Proof(The tests of steps 5 and 6)Original

Step 5. Let $M$ be compact and $f$ continuous on $M$ . If $M = \varnothing$ there is nothing to extremise, so assume $M \neq \varnothing$ . By the extreme value theorem (compactness), which is stated there for a non-empty compact space, $f$ attains a maximum at some $\hat x\in M$ and a minimum at some $\check x\in M$ . Now $\hat x$ is in particular a local maximum of $f$ on $M$ , so:

if $\operatorname{rank}J_g(\hat x) = m$ , the first-order condition applies at $\hat x$ and produces a $\lambda$ solving the $n+m$ equations, so $\hat x$ is among the solutions found in step 3;
otherwise $\operatorname{rank}J_g(\hat x) < m$ and $\hat x$ is listed by step 4.

Either way $\hat x \in C$ , and the same for $\check x$ . Hence $\max_M f = \max_C f$ and $\min_M f = \min_C f$ , and these are attained at points of $C$ — which is what step 5 asserts, with no second-order information used. $\square$

Step 6. Let $f, g$ be of class $C^2$ , let $x_*\in M$ be a regular point with multiplier $\lambda$ , and let $H = H_f(x_*) - \sum_i\lambda_iH_{g_i}(x_*)$ .

A $C^2$ parametrisation of $M$ near $x_*$ . Repeat steps 1 and 2 of the previous proof; since $g$ is now $C^2$ , the implicit function theorem returns an $h$ of class $C^2$ . Put $\phi(u) = (u, h(u))$ for $u\in V_0$ : a $C^2$ map whose image is $M\cap(V_0\times W_0)$ , and a homeomorphism onto that set, its inverse being the restriction of the projection $(u,y)\mapsto u$ . Write $P = D\phi(u_*)$ , an $n\times(n-m)$ matrix whose top block is $I_{n-m}$ , so $\operatorname{rank}P = n-m$ .
Its image is the tangent space. Differentiating $g(\phi(u)) = 0$ at $u_*$ gives $J_g(x_*)P = 0$ , so the image of $P$ lies in $T_{x_*}M = \ker J_g(x_*)$ ; both spaces have dimension $n-m$ , so they coincide.
The pulled-back Hessian is $H$ on the tangent space. Put $F = f\circ\phi$ and $G_i = g_i\circ\phi$ , all of class $C^2$ on $V_0$ . Two applications of the chain rule give, for $\alpha,\beta \le n-m$ , $\frac{\partial^2 F}{\partial u_\alpha \partial u_\beta}(u_*) = \bigl(P^{\mathsf T}H_f(x_*)P\bigr)_{\alpha\beta} + \sum_{j=1}^{n}\frac{\partial f}{\partial x_j}(x_*)\,\frac{\partial^2\phi_j}{\partial u_\alpha\partial u_\beta}(u_*),$ and the same formula with each $g_i$ in place of $f$ — whose left-hand side is $0$ , because $G_i$ vanishes identically on $V_0$ . Multiplying the $g_i$ identity by $\lambda_i$ , summing over $i$ and subtracting from the $f$ identity, $H_F(u_*) = P^{\mathsf T}HP + \sum_{j=1}^{n}\left(\frac{\partial f}{\partial x_j}(x_*) - \sum_{i=1}^{m}\lambda_i\frac{\partial g_i}{\partial x_j}(x_*)\right)H_{\phi_j}(u_*) = P^{\mathsf T}HP ,$ since each bracket is a component of $\nabla f(x_*) - \sum_i\lambda_i\nabla g_i(x_*)$ , which vanishes by the first-order condition. This is where the multiplier earns its keep: it is exactly what makes the second derivatives of the parametrisation drop out.
And $u_*$ is a critical point of $F$ . $\nabla F(u_*) = P^{\mathsf T}\nabla f(x_*) = 0$ , because the columns of $P$ lie in $T_{x_*}M$ and $\nabla f(x_*)$ is orthogonal to that space, by the geometric form of the first-order condition.
Conclude. Suppose $h^{\mathsf T}Hh>0$ for every non-zero $h\in T_{x_*}M$ . For $w \ne 0$ in $\mathbb{R}^{n-m}$ we have $Pw\ne0$ and $Pw\in T_{x_*}M$ , so $w^{\mathsf T}H_F(u_*)w = (Pw)^{\mathsf T}H(Pw)>0$ : the matrix $H_F(u_*)$ is positive definite. By part 1 of the classification in Taylor expansion and critical points in several variables, $u_*$ is a strict local minimum of $F$ : there is $\rho>0$ with $F(u)>F(u_*)$ for $0<\lvert u-u_*\rvert<\rho$ . Because $\phi$ is a homeomorphism onto $M\cap(V_0\times W_0)$ , every point of $M$ close enough to $x_*$ is $\phi(u)$ for exactly one $u$ with $\lvert u - u_*\rvert<\rho$ , and $x\ne x_*$ corresponds to $u\ne u_*$ ; hence $f(x)>f(x_*)$ for all such $x\ne x_*$ , which is a strict local minimum of $f$ on $M$ . The maximum case is this one applied to $-f$ , whose multiplier is $-\lambda$ and whose matrix is $-H$ . Finally, if $H$ is indefinite on $T_{x_*}M$ then $H_F(u_*)$ is indefinite, so by part 3 of the same classification every neighbourhood of $u_*$ contains points where $F$ exceeds $F(u_*)$ and points where it is smaller; transporting these through $\phi$ shows $x_*$ is neither a local minimum nor a local maximum of $f$ on $M$ . $\square$

RemarkOriginal

Where the hypotheses are used, and what fails without them.

The rank condition cannot be dropped. Minimise $f(x,y) = x$ subject to $g(x,y) = y^2 - x^3 = 0$ on $U = \mathbb{R}^2$ . The constraint set is the cuspidal curve $y^2 = x^3$ , which forces $x \ge 0$ , so the minimum of $f$ on $M$ is at the origin. But $\nabla g(0,0) = (-3x^2,2y)\big|_{(0,0)} = (0,0)$ has rank $0$ , and $\nabla f(0,0) = (1,0) \ne \lambda(0,0)$ for any $\lambda$ : the multiplier equation has no solution at the true minimiser, which step 3 misses and only step 4 finds.
Uniqueness of $\lambda$ comes from the same hypothesis: linear independence of the $\nabla g_i(x_*)$ is exactly what makes the coefficients in $\nabla f(x_*) = \sum_i \lambda_i \nabla g_i(x_*)$ unique.
Necessary, not sufficient. For $f(x,y) = y$ with $g(x,y) = y - x^3 = 0$ , $\nabla g = (-3x^2,1)$ never vanishes, so every point of $M$ is regular, and $(0,1) = \lambda(-3x^2,1)$ gives $\lambda = 1$ , $x = 0$ . But $f$ on $M$ is $x \mapsto x^3$ , with no local extremum at $x = 0$ . The system cannot tell maxima, minima and neither apart; that is what steps 5 and 6 are for.
$m < n$ is needed for the problem to be non-trivial. With $m = n$ and the rank condition, $g$ is a local diffeomorphism near $x_*$ by the inverse function theorem, so $M$ is locally the single point $x_*$ ; with $m > n$ the rank condition is unsatisfiable.
Equality constraints only. For inequality constraints $g_i(x) \le 0$ the first-order conditions carry sign restrictions on the multipliers and a complementarity condition, and are not covered here.
Sign conventions. We write $L = f - \sum_i \lambda_i g_i$ , which is Arfken's convention: his Eqs. (22.66)–(22.67), pp. 1108–1109, read $\partial f/\partial x_i - \sum_j \lambda_j\,\partial g_j/\partial x_i = 0$ . Riley, Hobson and Bence use the opposite sign — $f + \lambda g$ in Eqs. (5.27)–(5.28), $f + \lambda g + \mu h$ in (5.31) — which replaces each $\lambda_i$ by $-\lambda_i$ and changes nothing about the solutions $x_*$ . The sign matters when the multiplier itself is interpreted, as a sensitivity or as a physical quantity such as a force or a chemical potential: with our convention, and under hypotheses making the optimal point and multiplier differentiable in the constraint levels $c$ , $\lambda_i$ is $\partial/\partial c_i$ of the constrained optimal value of $f$ .
The analogous hypothesis in the variational setting. Gelfand and Fomin's isoperimetric theorem (Thm 1, §12.1, p. 43) assumes that the extremal of $J$ is not an extremal of the constraint functional $K$ before producing $\lambda$ ; their footnote names the analogy with the finite-dimensional method, and that hypothesis is the infinite-dimensional counterpart of $\operatorname{rank}J_g(x_*) = m$ . Evans's Theorem 2 (§8.4.1, p. 460) carries no such hypothesis: it assumes only that $u$ minimises $I$ over the admissible class. Non-degeneracy of the constraint enters there only inside the proof, as the case assumption (13) under which the implicit function theorem is applied, and step 3 (p. 462) treats the opposite case and finds the conclusion holds trivially, for any $\lambda$ . In our finite-dimensional statement the regularity hypothesis cannot be dispensed with that way, as the cuspidal curve above shows.
Smoothness class in the sources. Lee's Problem 11-11 states the result for a smooth function on a smooth manifold with a smooth local defining function — his name for a map whose zero set is the constraint set and whose differential is surjective there, so his regularity hypothesis contains ours. The version above is for $C^1$ data on an open subset of $\mathbb{R}^n$ , which is what the implicit function theorem route needs and the generality the physics texts use.

Example(Two constraints on a sphere, and a non-compact constraint set)Original

(a) The highest point of a tilted great circle. Maximise and minimise $f(x,y,z) = z$ on

M = \bigl\{ (x,y,z) : g_1 = x^2+y^2+z^2-1 = 0, \;\; g_2 = x+y+z = 0 \bigr\},

the great circle cut from the unit sphere by a plane through the origin. Here $n = 3$ , $m = 2$ .

Step 4 first: every point is regular. $\nabla g_1 = 2(x,y,z)$ and $\nabla g_2 = (1,1,1)$ are linearly dependent only if $(x,y,z)$ is a multiple of $(1,1,1)$ ; on $M$ the plane equation then forces $3x = 0$ , hence $(x,y,z) = (0,0,0)$ , which is not on the sphere. So $\operatorname{rank}J_g = 2$ at every point of $M$ and step 4 contributes no candidates.

Step 3. With $L = z - \lambda_1(x^2+y^2+z^2-1) - \lambda_2(x+y+z)$ the stationarity equations are

0 = 2\lambda_1 x + \lambda_2, \qquad 0 = 2\lambda_1 y + \lambda_2, \qquad 1 = 2\lambda_1 z + \lambda_2 ,

together with the two constraints. Subtracting the first two gives $2\lambda_1(x-y) = 0$ . If $\lambda_1 = 0$ the first equation gives $\lambda_2 = 0$ and the third then reads $1 = 0$ ; so $\lambda_1 \ne 0$ and $x = y$ . The plane gives $z = -2x$ , and the sphere $2x^2 + 4x^2 = 1$ , so $x = \pm 1/\sqrt6$ . The two candidates are

\Bigl(\tfrac{1}{\sqrt6}, \tfrac{1}{\sqrt6}, -\tfrac{2}{\sqrt6}\Bigr) \quad\text{and}\quad \Bigl(-\tfrac{1}{\sqrt6},-\tfrac{1}{\sqrt6}, \tfrac{2}{\sqrt6}\Bigr).

Their multipliers follow from the first and third equations: $1 = 2\lambda_1(z-x)$ , and for the second point $z - x = 3/\sqrt6$ , so $\lambda_1 = 1/\sqrt6$ and $\lambda_2 = -2\lambda_1x = 1/3$ .

Step 5. $M$ is closed (an intersection of two level sets of continuous functions) and bounded (it lies on the unit sphere), hence compact, and $f$ is continuous, so the maximum and the minimum of $f$ on $M$ are attained and lie among the two candidates. Comparing the two values of $z$ ,

\max_M z = \frac{2}{\sqrt6} = \sqrt{\tfrac23} \approx 0.8165, \qquad \min_M z = -\frac{2}{\sqrt6} \approx -0.8165 ,

with no second-order test needed. Independent check: on $M$ , $x + y = -z$ and $x^2+y^2 = 1-z^2$ , while $x^2+y^2 \ge (x+y)^2/2$ always; so $1 - z^2 \ge z^2/2$ , that is $z^2 \le 2/3$ , with equality when $x = y$ — exactly the bound just found.

(b) A non-compact constraint set, where step 6 is needed. Minimise $f(x,y) = x+y$ on $M = \{xy = 1,\; x>0\}$ , taking $U = \{(x,y) : x > 0\}$ and $g = xy - 1$ . Then $\nabla g = (y,x)$ , which is non-zero on $M$ , so every point is regular and step 4 is empty. Step 3: $(1,1) = \lambda(y,x)$ gives $y = x = 1/\lambda$ , and $xy = 1$ with $x>0$ gives $x = y = 1$ , $\lambda = 1$ ; the single candidate is $(1,1)$ with $f = 2$ . $M$ is closed but unbounded, so step 5 does not apply. Step 6:

H = H_f(1,1) - \lambda H_g(1,1) = \begin{pmatrix}0&0\\0&0\end{pmatrix} - 1\cdot\begin{pmatrix}0&1\\1&0\end{pmatrix} = \begin{pmatrix}0&-1\\-1&0\end{pmatrix},

and $T_{(1,1)}M = \ker\begin{pmatrix}1&1\end{pmatrix} = \operatorname{span}\{(1,-1)^{\mathsf T}\}$ . For $h = t(1,-1)^{\mathsf T}$ with $t \ne 0$ ,

h^{\mathsf T}Hh = t^2\,(1,-1)\begin{pmatrix}0&-1\\-1&0\end{pmatrix}\begin{pmatrix}1\\-1\end{pmatrix} = t^2\,(1,-1)\begin{pmatrix}1\\-1\end{pmatrix} = 2t^2 > 0,

so $(1,1)$ is a strict local minimum of $f$ on $M$ . Note that $H$ itself is indefinite on $\mathbb{R}^2$ (its eigenvalues are $\pm1$ ): definiteness on the tangent line is all that is required, and all that holds. Check: for $x>0$ , $x + 1/x \ge 2$ by $(\sqrt x - 1/\sqrt x)^2 \ge 0$ , with equality only at $x=1$ , so $(1,1)$ is in fact the global minimiser and $\min_M f = 2$ .

Application(The Boltzmann distribution as a constrained maximum of entropy)Original

A system with $n \ge 3$ accessible microstates of energies $E_1,\dots,E_n$ , not all equal, is described by probabilities $p = (p_1,\dots,p_n)$ . Gibbs's variational principle takes the equilibrium distribution to be the one maximising the entropy

S(p) = -k\sum_{i=1}^{n} p_i\ln p_i

on the open set $U = \{p : p_i>0 \text{ for all } i\}$ , subject to the two constraints

g_1(p) = \sum_i p_i - 1 = 0, \qquad g_2(p) = \sum_i p_iE_i - \bar U = 0 ,

normalisation and a prescribed mean energy $\bar U$ . Here $m = 2 < n$ .

Regularity. $\nabla g_1 = (1,\dots,1)$ and $\nabla g_2 = (E_1,\dots,E_n)$ are linearly independent precisely because the energies are not all equal, so every point of $M$ is a regular point and the first-order condition applies at any maximiser.

The multiplier equations. $\partial S/\partial p_i = -k(\ln p_i + 1)$ , so $\nabla S = \lambda_1\nabla g_1 + \lambda_2\nabla g_2$ reads

-k(\ln p_i + 1) = \lambda_1 + \lambda_2 E_i \qquad (i = 1,\dots,n),

that is $p_i = e^{-1-\lambda_1/k}\,e^{-\beta E_i}$ with $\beta = \lambda_2/k$ . The factor in front is the same for every $i$ , so imposing $g_1 = 0$ fixes it and gives the Boltzmann distribution

p_i = \frac{e^{-\beta E_i}}{Z}, \qquad Z = \sum_{j=1}^{n} e^{-\beta E_j} ,

with $\beta$ then determined by $g_2 = 0$ , that is by $-\partial(\ln Z)/\partial\beta = \bar U$ . Two features of the method are doing real work here. The multipliers are not nuisances to be eliminated but the answer: $\lambda_2 = k\beta$ is, by the sensitivity reading of the multiplier noted in the remark above, the derivative of the maximal entropy with respect to the mean energy, which thermodynamics calls $1/T$ — so $\beta = 1/(kT)$ , and the temperature is the multiplier conjugate to the energy constraint. And the single candidate produced is the maximiser: $S$ is strictly concave on $U$ (each $-p\ln p$ is), while $M$ is the intersection of $U$ with an affine subspace and so convex, and a critical point of a strictly concave function on a convex set is its unique maximiser.

A worked instance. Take three equally spaced levels $E_i = (i-1)\epsilon$ for $i = 1,2,3$ and $\bar U = 0.8\,\epsilon$ . With $w = e^{-\beta\epsilon}$ the constraint $g_2 = 0$ becomes

\frac{w + 2w^2}{1 + w + w^2} = \frac45 \;\Longleftrightarrow\; 6w^2 + w - 4 = 0 \;\Longrightarrow\; w = \frac{\sqrt{97}-1}{12} \approx 0.7374048 ,

taking the positive root. Then $\beta\epsilon = -\ln w \approx 0.3046183$ , so $kT \approx 3.2828\,\epsilon$ ; $Z = 1+w+w^2 \approx 2.2811707$ and

p \approx (0.4383714,\; 0.3232572,\; 0.2383714), \qquad S/k \approx 1.0683834 .

The mean energy checks: $0.3232572 + 2(0.2383714) = 0.8000000$ . And $S/k$ sits just below $\ln 3 \approx 1.0986123$ , the entropy of the uniform distribution, which is the unconstrained maximum — the energy constraint costs a little entropy, as it must.

Integration

Method

Multiple integrals and change of variables

Adapted

Add to PDF

States when a multiple integral may be computed as an iterated integral and gives the change-of-variables formula for a $C^1$ diffeomorphism, with the absolute Jacobian determinant as the local volume factor, specialised to plane polar, cylindrical and spherical coordinates.

Prerequisites: Partial derivatives and the chain rule, Implicit and inverse function theorems, The Lebesgue integral, Fubini and Tonelli theorems

Definition(Iterated integral, Jacobian factor, admissible coordinate map)AdaptedAdapted from Riley, Hobson & Bence (2006), §6.1, Eqs. (6.3)–(6.4), pp. 188–190; Folland (1999), §2.6, p. 74

An iterated integral of $f$ over a region described by $a \le x \le b$ , $y_1(x) \le y \le y_2(x)$ is $\int_a^b \left( \int_{y_1(x)}^{y_2(x)} f(x,y)\, dy \right) dx ,$ the inner integral being computed with $x$ held fixed. The limits of the inner integral may depend on the outer variable; those of the outer integral may not. The convention extends to any number of variables.
Let $\Omega \subseteq \mathbb{R}^n$ be open. A map $G : \Omega \to \mathbb{R}^n$ is a $C^1$ diffeomorphism (Folland's usage) if $G$ is of class $C^1$ , injective on $\Omega$ , and $\det D_uG \neq 0$ for every $u \in \Omega$ . Then $G(\Omega)$ is open and $G^{-1}$ is $C^1$ (inverse function theorem).
For such a $G$ the Jacobian factor at $u \in \Omega$ is $\lvert \det D_uG \rvert$ , written $\bigl\lvert \partial(x_1,\dots,x_n)/\partial(u_1,\dots,u_n) \bigr\rvert$ when $x = G(u)$ . The absolute value is taken: the sign of the determinant records only the orientation of the coordinate ordering.

Method(When a multiple integral is an iterated integral)AdaptedAdapted from Folland (1999), Thm 2.37, p. 67 and Thm 2.39, p. 68; Riley, Hobson & Bence (2006), §6.1, Eqs. (6.3)–(6.4), pp. 188–190

Box, continuous integrand. If $f$ is continuous on the closed bounded box $R = [a_1,b_1]\times\cdots\times[a_n,b_n] \subseteq \mathbb{R}^n$ , then $\int_R f\, d^n x$ equals the iterated integral in each of the $n!$ orders of the variables.
Region with variable limits. Let $y_1 \le y_2$ be continuous on $[a,b]$ , let $R = \bigl\{ (x,y) : a \le x \le b,\ y_1(x) \le y \le y_2(x) \bigr\},$ and let $f$ be continuous on $R$ . Then $\int_R f \, dA = \int_a^b \left( \int_{y_1(x)}^{y_2(x)} f(x,y)\, dy \right) dx .$ If the same $R$ can also be written as $\{(x,y) : c \le y \le d,\ x_1(y) \le x \le x_2(y)\}$ with $x_1 \le x_2$ continuous, the iterated integral in the other order has the same value. The same holds in $n$ variables, with the innermost limits depending on all the outer variables.
The general condition. Let $f$ be Lebesgue measurable on $\mathbb{R}^p \times \mathbb{R}^q$ . If $f \ge 0$ , or if $\int \lvert f \rvert < \infty$ , then both iterated integrals and the multiple integral agree (for $f \ge 0$ all three may be $+\infty$ ); the inner integrals are then defined for almost every value of the outer variable, and are measurable functions of it. This is Folland's Thm 2.39, the form of Fubini–Tonelli for the completed product measure, which is the one that applies to Lebesgue measurable integrands. If neither hypothesis holds, the two iterated integrals can both exist and be different; see Fubini and Tonelli, where these hypotheses and a counterexample belong. Parts 1 and 2 are the cases in which $f$ is bounded on a bounded region, so that $\int \lvert f \rvert < \infty$ automatically.

Idea of the proof(Parts 1 and 2 are part 3 with the hypothesis checked)Original

Part 3 asks for absolute integrability; parts 1 and 2 are the situations in which that comes free, because a continuous function on a compact set is bounded and a bounded function on a bounded region has a finite integral of its modulus. So the only work is bookkeeping: extend the integrand by zero outside the region, so that it is defined on the whole of $\mathbb{R}^n$ and part 3 can be quoted, and then split off one variable at a time.

Proof(Parts 1 and 2, from part 3)AdaptedAdapted from Folland (1999), Thm 2.37, p. 67 and Thm 2.39, p. 68

Part 1. $f$ is continuous on the compact box $R$ , hence bounded, say $\lvert f\rvert \le C$ . Let $\tilde f = f$ on $R$ and $\tilde f = 0$ elsewhere. Then $\tilde f$ is Lebesgue measurable — for $c \ge 0$ the set $\{\tilde f > c\}$ is $\{x \in R: f(x) > c\}$ , the intersection of $R$ with an open set because $f$ is continuous on $R$ , and for $c < 0$ it is that set together with the open complement of $R$ — and

\int_{\mathbb{R}^n} \lvert \tilde f\rvert \, d^nx \le C \operatorname{vol}(R) < \infty .

So part 3 applies to $\tilde f$ with the splitting $\mathbb{R}^n = \mathbb{R}^{n-1}\times\mathbb{R}$ and gives

\int_R f\,d^nx = \int_{\mathbb{R}^n}\tilde f \,d^nx = \int_{\mathbb{R}^{n-1}}\left( \int_{\mathbb{R}} \tilde f(x', x_n)\, dx_n\right) d^{n-1}x' ,

where the inner integral is $\int_{a_n}^{b_n} f(x',x_n)\,dx_n$ for $x'$ in the box $R' = [a_1,b_1]\times\cdots\times[a_{n-1},b_{n-1}]$ and $0$ for $x'$ outside it. Repeating the step on $R'$ , and $n-1$ times in all, produces the iterated integral with the variables peeled off in the order $x_n, x_{n-1}, \dots, x_1$ . Lebesgue measure on $\mathbb{R}^n$ is unchanged by a permutation of the coordinates — the case $T =$ a permutation matrix, $\lvert\det T\rvert = 1$ , of the change-of-variables statement below — so relabelling the variables first and then repeating the argument gives each of the $n!$ orders, all equal to $\int_R f$ .

Part 2. The region

R = \bigl\{(x,y) \in [a,b]\times\mathbb{R} \;:\; y_1(x) \le y \le y_2(x)\bigr\}

is closed, being cut out of the closed set $[a,b]\times\mathbb{R}$ by the two non-strict inequalities $y - y_2(x) \le 0$ and $y_1(x) - y \le 0$ whose left sides are continuous there; and it is bounded, because $y_1$ and $y_2$ are continuous on a compact interval and so bounded. Hence $R$ is compact, and the extension $\tilde f$ of $f$ by zero is measurable with $\int\lvert\tilde f\rvert < \infty$ , exactly as in part 1. Part 3 with $p = q = 1$ then gives

\int_R f\, dA = \int_{\mathbb{R}^2}\tilde f \,dA = \int_{\mathbb{R}}\left(\int_{\mathbb{R}} \tilde f(x,y)\,dy\right) dx ,

and for fixed $x \in [a,b]$ the function $\tilde f(x,\cdot)$ is $f(x,\cdot)$ on $[y_1(x), y_2(x)]$ and $0$ off it, while for $x \notin [a,b]$ it vanishes identically. So the right-hand side is $\int_a^b \bigl(\int_{y_1(x)}^{y_2(x)} f(x,y)\,dy\bigr)dx$ , which is the claim. If $R$ has the second description as well, the same computation with the two variables exchanged — part 3 gives both iterated integrals, with the same value — yields the other order. In $n$ variables the argument is the same, applied once for each variable, the innermost limits being allowed to depend on all the outer variables because nothing above required the limits to be constants. $\square$

NoteOriginal

Part 3 is the Fubini–Tonelli theorem, in the form Folland states for the completed product measure (his Thm 2.39, p. 68), which is the form that applies to Lebesgue measurable integrands. It is not proved in this unit: it belongs to Fubini and Tonelli, along with the construction of the product measure and the counterexample showing that the hypothesis cannot be dropped. Everything here that needs it quotes it.

Method(Change of variables)AdaptedAdapted from Folland (1999), §2.6, Thm 2.47(a)–(b), p. 74

Let $\Omega \subseteq \mathbb{R}^n$ be open and let $G : \Omega \to \mathbb{R}^n$ be a $C^1$ diffeomorphism, so that $G$ is $C^1$ , injective, and $\det D_uG \neq 0$ for every $u \in \Omega$ . Let $f$ be a Lebesgue measurable function on $G(\Omega)$ . If $f \ge 0$ , or if $f$ is integrable on $G(\Omega)$ , then $f \circ G$ is measurable on $\Omega$ and

\int_{G(\Omega)} f(x)\, d^{n}x \;=\; \int_{\Omega} f\bigl(G(u)\bigr)\, \bigl\lvert \det D_uG \bigr\rvert \, d^{n}u .

Part (b) of the same theorem is the volume rule: for Lebesgue measurable $E \subseteq \Omega$ , the image $G(E)$ is measurable and

\operatorname{vol}\bigl(G(E)\bigr) = \int_E \lvert \det D_uG \rvert \, d^n u

— equivalently, the displayed formula applied to $f = \chi_{G(E)}$ , which is the indicator function on $G(\Omega)$ that the hypotheses admit. The Jacobian factor is thus the local ratio by which $G$ stretches $n$ -dimensional volume.

Algebra of Jacobian determinants. For composable $C^1$ changes of variables $x = x(y)$ , $y = y(z)$ , the chain rule and the multiplicativity of determinants give

\frac{\partial(x_1,\dots,x_n)}{\partial(z_1,\dots,z_n)} = \frac{\partial(x_1,\dots,x_n)}{\partial(y_1,\dots,y_n)}\; \frac{\partial(y_1,\dots,y_n)}{\partial(z_1,\dots,z_n)}, \qquad \frac{\partial(x_1,\dots,x_n)}{\partial(y_1,\dots,y_n)} = \left[ \frac{\partial(y_1,\dots,y_n)}{\partial(x_1,\dots,x_n)} \right]^{-1},

the second with both sides evaluated at corresponding points. These help when one of the two sets of partial derivatives is much easier to compute than the other.

Idea of the proof(Why the Jacobian factor is the local volume ratio)AdaptedAdapted from Folland (1999), §2.6, proof of Thm 2.47, pp. 75–76

Near a point $u_0$ the map $G$ agrees, to higher order, with the affine map $u \mapsto G(u_0) + D_{u_0}G\,(u - u_0)$ , and an affine map multiplies the volume of every set by $\lvert\det\rvert$ of its linear part — the linear case of the theorem, which follows from Fubini–Tonelli and the one-variable scaling and translation rules. So on a cube small enough that $D_uG$ hardly varies, $G$ multiplies volume by nearly $\lvert\det D_uG\rvert$ . Subdividing a cube into many such pieces and adding up gives a sum that is a Riemann-type approximation to $\int\lvert\det D_uG\rvert$ , and the mesh can be sent to zero because $D_uG$ is continuous. That gives the volume rule, and the formula for integrals then follows by approximating $f$ by simple functions.

Proof(The volume rule and the algebra of Jacobian determinants)AdaptedAdapted from Folland (1999), §2.6, Thm 2.47(b), p. 74; Arfken, Weber & Harris (2013), §4.4, Eq. (4.75), p. 230; Riley, Hobson & Bence (2006), §6.4.4, Eqs. (6.16)–(6.17), pp. 206–207

The volume rule, given the integral formula. Let $E \subseteq \Omega$ be Lebesgue measurable; that $G(E)$ is again measurable is part of Folland's statement and is not re-proved here. The function $f = \chi_{G(E)}$ is measurable and non-negative on $G(\Omega)$ , so the hypotheses of the displayed formula are met, and its left-hand side is $\int_{G(\Omega)}\chi_{G(E)}\,d^nx = \operatorname{vol}(G(E))$ . On the right, $\chi_{G(E)}(G(u)) = 1$ exactly when $G(u) \in G(E)$ , which by the injectivity of $G$ happens exactly when $u \in E$ ; so $\chi_{G(E)}\circ G = \chi_E$ and the right-hand side is $\int_E\lvert\det D_uG\rvert\,d^nu$ . This is the volume rule, and it is why injectivity, and not merely the non-vanishing of the determinant, is needed: a map that covered part of $G(\Omega)$ twice would make $\chi_{G(E)}\circ G$ larger than $\chi_E$ .

The composition rule. Let $x = x(y)$ and $y = y(z)$ be $C^1$ with the composite defined. The chain rule gives $D_z x = (D_yx)(D_zy)$ at corresponding points, and the determinant of a product of square matrices is the product of the determinants, which is the first display.

The reciprocal rule. Apply the composition rule with $z = x$ , so that the composite is the identity map, whose derivative is $I$ with determinant $1$ :

\frac{\partial(x_1,\dots,x_n)}{\partial(y_1,\dots,y_n)}\;\frac{\partial(y_1,\dots,y_n)}{\partial(x_1,\dots,x_n)} = 1 ,

both factors evaluated at corresponding points. This is the second display, and it is part 2 of the inverse function theorem read as an identity between determinants — the relation Arfken writes as $J^{-1} = \partial(u)/\partial(x)$ . In particular neither determinant can vanish where the other is defined. $\square$

NoteOriginal

The change-of-variables formula itself — the first display of the statement — is not proved in this handbook. Folland's proof (§2.6, pp. 75–76) runs along the lines of the idea above: it settles the linear case first, then bounds the measure of the image of a cube by that of the image of its affine approximation, sums over a fine subdivision and passes to the limit, and finally extends from indicator functions to simple and then to measurable functions by monotone convergence. It needs his decomposition of an open set into cubes with disjoint interiors and the regularity of Lebesgue measure, neither of which is developed here; he also notes (p. 74) that he proves a more general result by different methods in his §11.2. Everything else in this unit is proved from the formula, not alongside it.

Method(The three standard coordinate maps)AdaptedAdapted from Riley, Hobson & Bence (2006), §6.4.3, p. 205; Lee (2013), App. C, Examples C.37–C.38, pp. 660–661; Arfken, Weber & Harris (2013), §1.10, p. 73

Each map below is a $C^\infty$ diffeomorphism from the stated open set onto its image, whose complement in $\mathbb{R}^n$ is a Lebesgue null set; so by the change-of-variables theorem the resulting formula computes the integral over all of $\mathbb{R}^n$ .

Plane polar. $G(r,\varphi) = (r\cos\varphi,\; r\sin\varphi)$ on $(0,\infty) \times (-\pi,\pi)$ , with image $\mathbb{R}^2$ minus the closed half-line $\{(x,0) : x \le 0\}$ : $\det DG = r, \qquad dA = r\, dr\, d\varphi .$
Cylindrical polar. $G(\rho,\varphi,z) = (\rho\cos\varphi,\; \rho\sin\varphi,\; z)$ on $(0,\infty) \times (-\pi,\pi) \times \mathbb{R}$ , with image $\mathbb{R}^3$ minus the closed half-plane $\{(x,0,z) : x \le 0\}$ : $\det DG = \rho, \qquad dV = \rho\, d\rho\, d\varphi\, dz .$
Spherical polar (physics convention: $\theta$ the polar angle from the positive $z$ -axis, $\varphi$ the azimuth). $G(r,\theta,\varphi) = (r\sin\theta\cos\varphi,\; r\sin\theta\sin\varphi,\; r\cos\theta)$ on $(0,\infty) \times (0,\pi) \times (-\pi,\pi)$ , with the same image as in 2: $\det DG = r^2 \sin\theta, \qquad dV = r^2 \sin\theta\, dr\, d\theta\, d\varphi .$ Since $0 < \theta < \pi$ , $\sin\theta > 0$ and the absolute value may be dropped.

Idea of the proof(The three coordinate maps)Original

Each determinant is a direct computation. The real content is set-theoretic: the map is injective on the stated box because it can be inverted by hand — the radius is the norm of the point, the polar angle comes from the last coordinate, the azimuth from the first two — and the injectivity of the azimuth is just the fact that two angles in an interval of length $2\pi$ with the same sine and cosine are equal. What the box leaves out is a piece of a hyperplane, which is null because each of its one-dimensional sections is a single point, and Tonelli turns that into a statement about $n$ -dimensional measure.

Proof(The three coordinate maps)AdaptedAdapted from Riley, Hobson & Bence (2006), §6.4.3, p. 205; Lee (2013), App. C, Examples C.37–C.38, pp. 660–661

The determinants. For plane polar, $DG(r,\varphi) = \begin{pmatrix} \cos\varphi & -r\sin\varphi \\ \sin\varphi & r\cos\varphi\end{pmatrix}, \qquad \det DG = r\cos^2\varphi + r\sin^2\varphi = r .$ For cylindrical polar the matrix is $DG(\rho,\varphi,z) = \begin{pmatrix} \cos\varphi & -\rho\sin\varphi & 0 \\ \sin\varphi & \rho\cos\varphi & 0 \\ 0 & 0 & 1\end{pmatrix},$ and expanding along the last row gives $\det DG = \rho$ , the previous determinant. For spherical polar, $DG(r,\theta,\varphi) = \begin{pmatrix} \sin\theta\cos\varphi & r\cos\theta\cos\varphi & -r\sin\theta\sin\varphi \\ \sin\theta\sin\varphi & r\cos\theta\sin\varphi & r\sin\theta\cos\varphi \\ \cos\theta & -r\sin\theta & 0\end{pmatrix} ;$ expanding along the last row, the cofactor of the $(3,1)$ entry is $r^2\cos\theta\sin\theta\cos^2\varphi + r^2\sin\theta\cos\theta\sin^2\varphi = r^2\sin\theta\cos\theta ,$ and that of the $(3,2)$ entry is $-\bigl(r\sin^2\theta\cos^2\varphi + r\sin^2\theta\sin^2\varphi\bigr) = -r\sin^2\theta$ , so $\det DG = \cos\theta\cdot r^2\sin\theta\cos\theta + (-r\sin\theta)\cdot(-r\sin^2\theta) = r^2\sin\theta\,(\cos^2\theta + \sin^2\theta) = r^2\sin\theta .$ On each box the determinant is positive, since $r, \rho > 0$ and $\sin\theta > 0$ for $0 < \theta < \pi$ .
Injectivity, and the image, for plane polar. Let $x = G(r,\varphi)$ with $(r,\varphi) \in (0,\infty)\times(-\pi,\pi)$ . Then $\lvert x\rvert = r$ , so $r$ is determined by $x$ ; and $x/r = (\cos\varphi,\sin\varphi)$ , so $\varphi$ is determined too, because if $\cos\varphi = \cos\varphi'$ and $\sin\varphi = \sin\varphi'$ then $\varphi - \varphi' \in 2\pi\mathbb{Z}$ , and the only member of $2\pi\mathbb{Z}$ in $(-2\pi,2\pi)$ is $0$ . So $G$ is injective. Its image omits the origin ( $r > 0$ ) and omits the points with $x/\lvert x\rvert = (-1,0)$ , which is the case $\varphi = \pm\pi$ , excluded from the box; everything else is attained, by taking $r = \lvert x\rvert$ and $\varphi$ the unique angle in $(-\pi,\pi)$ with $(\cos\varphi,\sin\varphi) = x/\lvert x\rvert$ . The image is therefore $\mathbb{R}^2$ minus the closed half-line $\{(x,0) : x \le 0\}$ .
Cylindrical polar. The first two coordinates are the plane polar map and the third is the identity, so injectivity and the image follow from step 2 applied in each horizontal plane: the complement of the image is the closed half-plane $\{(x,0,z) : x \le 0\}$ .
Spherical polar. Let $x = G(r,\theta,\varphi)$ on the stated box. Then $\lvert x\rvert^2 = r^2\sin^2\theta(\cos^2\varphi+\sin^2\varphi) + r^2\cos^2\theta = r^2$ , so $r = \lvert x\rvert$ is determined. Next $x_3/r = \cos\theta$ , and $\cos$ is injective on $(0,\pi)$ , so $\theta$ is determined. Then $\sin\theta > 0$ and $r\sin\theta = \sqrt{r^2 - x_3^2} = \lvert(x_1,x_2)\rvert$ , so $(x_1,x_2)/(r\sin\theta)$ is a unit vector equal to $(\cos\varphi,\sin\varphi)$ , and $\varphi$ is determined as in step 2. So $G$ is injective. Conversely, let $x$ be outside the closed half-plane $\{(x,0,z) : x\le 0\}$ ; then $(x_1,x_2) \neq (0,0)$ and $(x_1,x_2)/\lvert(x_1,x_2)\rvert \ne (-1,0)$ , so putting $r = \lvert x\rvert > 0$ , $\theta = \arccos(x_3/r) \in (0,\pi)$ and $\varphi$ the unique angle in $(-\pi,\pi)$ with $(\cos\varphi,\sin\varphi) = (x_1,x_2)/\lvert(x_1,x_2)\rvert$ recovers $x = G(r,\theta,\varphi)$ . Points of that half-plane are not attained: $r>0$ excludes the origin, $0<\theta<\pi$ excludes the rest of the $z$ -axis, and $\lvert\varphi\rvert < \pi$ excludes the open half-plane $\{(x,0,z): x<0\}$ . So the image is the same as in step 3.
The omitted sets are null. In $\mathbb{R}^2$ the complement of the image lies in the line $L = \{(x,y) : y = 0\}$ , which is closed and hence measurable, and by part 3 of the first statement, applied to the non-negative function $\chi_L$ , $\operatorname{vol}(L) = \int_{\mathbb{R}}\left(\int_{\mathbb{R}} \chi_L(x,y)\, dy\right) dx = \int_{\mathbb{R}} 0 \, dx = 0 ,$ because for each fixed $x$ the inner integrand is the indicator of the single point $y = 0$ . The same computation in $\mathbb{R}^3$ , integrating in $y$ first, gives measure $0$ for the plane $\{y = 0\}$ and hence for the closed half-plane inside it. A subset of a null set is null.
The formulas. Each $G$ is $C^\infty$ , injective, and has non-vanishing Jacobian determinant on its box, so it is a $C^1$ diffeomorphism in the sense used here and the change-of-variables theorem applies with $\lvert \det DG\rvert$ equal to $r$ , $\rho$ and $r^2\sin\theta$ respectively. Since $\mathbb{R}^n \setminus G(\Omega)$ is null, $\int_{\mathbb{R}^n} f = \int_{G(\Omega)}f$ for every $f$ that is non-negative or integrable, which is the clause at the head of the statement. $\square$

NoteOriginal

What it means. Two different reductions are in play. An iterated integral turns one $n$ -dimensional integral into $n$ one-dimensional ones; what must be checked is absolute integrability, and what is gained is that each stage is one-variable calculus. A change of variables replaces the region by an easier one; what must be checked is that the new coordinates label each point of the region once, and the price is the Jacobian factor, the local volume distortion of the coordinate map. In practice: choose coordinates that make the region a box, then iterate.

Method(Computing a multiple integral)AdaptedAdapted from Riley, Hobson & Bence (2006), §6.1, pp. 188–190 and §6.4.3, p. 205; Arfken, Weber & Harris (2013), §1.10, p. 73

Describe the region $R$ by inequalities, and sketch it if $n \le 3$ . The limits are a property of $R$ , not of the integrand.
Choose coordinates. Pick a map $G$ making $R$ the image of a box or of a region with simple limits: polar or cylindrical for symmetry about an axis, spherical for symmetry about a point, a linear map for an ellipsoid or a parallelepiped.
Check the map on the open set you will integrate over: $G$ must be $C^1$ , injective and have non-vanishing Jacobian determinant there. Excluding a null set — the origin, a half-plane cut, an axis — to achieve this does not change the integral.
Compute the Jacobian factor $\lvert \det D_uG \rvert$ ; if the inverse map is the one given to you, compute $\partial(u)/\partial(x)$ and take the reciprocal.
Substitute: $f(x) \mapsto f(G(u))$ , $d^nx \mapsto \lvert \det D_uG\rvert\,d^nu$ , and $R \mapsto G^{-1}(R)$ .
Justify iterating, then iterate. Check that the integrand is non-negative, or that the integral of its absolute value is finite — Tonelli applied to $\lvert f\rvert$ is the usual way — then iterate in whichever order makes the limits simplest. Under the same condition the order may be exchanged, which often turns variable limits into constants.

How it fails. Step 3 is where most errors occur: a map that covers part of the region twice doubles that part's contribution, and one whose Jacobian determinant vanishes on a set of positive measure is not a diffeomorphism there at all. Step 6 cannot be skipped when $f$ changes sign and the region is unbounded, or when $f$ is unbounded: both iterated integrals can then exist and differ. For a region with part of its boundary at infinity the integral is a limit over an increasing family of bounded regions, and arguments that "neglect" the corners of that family are not valid — Arfken's warning on p. 73.

RemarkOriginal

Where the hypotheses are used, and conventions that differ.

Injectivity. The polar map on $(0,\infty)\times(0,4\pi)$ has Jacobian determinant $r \ne 0$ throughout but wraps the punctured plane twice, and the formula applied to it returns exactly twice the correct answer. Non-vanishing of the determinant is local and does not give injectivity; both are hypotheses.
Null sets removed to get a diffeomorphism. The polar map is not injective on $(0,\infty)\times[-\pi,\pi]$ , since the ends of the $\varphi$ -interval give the same points; on $[0,\infty)\times[-\pi,\pi]$ the determinant also vanishes, along $r = 0$ . The spherical map degenerates likewise at $\theta \in \{0,\pi\}$ . Deleting these sets costs nothing, because their images — a half-line in $\mathbb{R}^2$ , a half-plane and an axis in $\mathbb{R}^3$ — are Lebesgue null.
The integrability hypothesis in Fubini. Both iterated integrals can exist and disagree when $\int\lvert f\rvert = \infty$ ; Folland's remarks after Thm 2.37 record this and send the reader to his Exercises 47 and 48 (pp. 68–69) for the examples, and Riley, Hobson and Bence warn about unbounded regions and discontinuous integrands for the same reason. For non-negative integrands no such hypothesis is needed.
Which matrix, and whether it is transposed. Arfken (Eq. (4.74)) and Riley, Hobson and Bence (p. 205) lay the Jacobian matrix out with $\partial x_j/\partial u_i$ in row $i$ , the transpose of the convention used here. Transposition does not change a determinant, so every formula above is the same either way; the layout matters only if the matrix itself is wanted.
Spherical angle names. We use the physics convention, $\theta$ the polar angle and $\varphi$ the azimuth, giving $r^2\sin\theta$ , as do Riley, Hobson and Bence (p. 205). Lee's Example C.38 writes the same map as $(\rho,\varphi,\theta)$ with $\varphi$ from the positive $z$ -axis and $\theta$ the azimuth, and so records the determinant as $\rho^2\sin\varphi$ . Only the names differ, but mixing them silently exchanges the two angles.
Relation to the elementary statement. The theorem above is Folland's, for Lebesgue measurable $f$ ; for continuous $f$ on a bounded region with null boundary it reduces to the Riemann-integral version of calculus texts. Riley, Hobson and Bence state the $n$ -dimensional rule without proof (§6.4.4) and derive the volume element geometrically from a parallelepiped.

Example(A region that becomes a box, and an integral over a ball)Original

(a) Hyperbolic coordinates. Let $D$ be the part of the open first quadrant between the hyperbolas $xy = 1$ and $xy = 4$ and between the lines $y = x$ and $y = 4x$ . Compute the area of $D$ and $\iint_D xy \, dA$ .

Step 1, the region. $D = \{(x,y) : x>0,\ y>0,\ 1 \le xy \le 4,\ 1 \le y/x \le 4\}$ , a curvilinear quadrilateral with no description by simple limits in $x$ and $y$ : the roles of the two hyperbolas and the two lines change as $x$ increases.

Step 2, the coordinates. The two pairs of boundary curves are level sets of $u = xy$ and $v = y/x$ , so take these as the new variables. On the open first quadrant $Q = (0,\infty)^2$ the map $\Psi(x,y) = (xy,\, y/x)$ has the explicit inverse

G(u,v) = \Bigl( \sqrt{u/v},\; \sqrt{uv} \Bigr), \qquad (u,v) \in (0,\infty)^2 ,

as one checks by composing both ways, and $D = G\bigl([1,4]\times[1,4]\bigr)$ : a box.

Step 3, the map is admissible. $G$ is $C^\infty$ on $(0,\infty)^2$ and injective, having the inverse $\Psi$ . The boundary of the box is a null set, so we may work on the open box $B = (1,4)\times(1,4)$ and integrate over $G(B)$ : it differs from $D$ by the image of that boundary, which is null by the volume rule, $\operatorname{vol}(G(N)) = \int_N\lvert\det D G\rvert = 0$ for a null $N$ .

Step 4, the Jacobian factor, by the reciprocal rule. It is easier to differentiate $\Psi$ :

\frac{\partial(u,v)}{\partial(x,y)} = \det\begin{pmatrix} y & x \\ -y/x^2 & 1/x\end{pmatrix} = \frac{y}{x} + \frac{xy}{x^2} = \frac{2y}{x} = 2v ,

so by the reciprocal rule $\det D_{(u,v)}G = 1/(2v)$ , which is positive, and the Jacobian factor is $1/(2v)$ .

Steps 5 and 6, substitute and iterate. Both integrands below are non-negative, so Tonelli licenses the iteration in either order. For the area,

\operatorname{area}(D) = \iint_B \frac{du\,dv}{2v} = \left(\int_1^4 du\right)\left(\int_1^4 \frac{dv}{2v}\right) = 3 \cdot \tfrac12\ln 4 = 3\ln 2 \approx 2.079442 .

For the second integral, $xy = u$ in the new variables, so

\iint_D xy\, dA = \iint_B \frac{u}{2v}\,du\,dv = \left(\int_1^4 u\,du\right)\left(\int_1^4\frac{dv}{2v}\right) = \frac{15}{2}\ln 2 \approx 5.198604 .

Check. Integrating numerically in the original variables, over $x$ from $1/2$ to $2$ with $y$ between $\max(1/x, x)$ and $\min(4/x, 4x)$ , gives $2.0794415$ and $5.1986039$ — the closed forms to seven figures.

(b) Spherical coordinates. Compute $I = \int_{B_R}(x^2+y^2)\,dV$ over the ball $\lvert x\rvert \le R$ . In spherical polars $x^2+y^2 = r^2\sin^2\theta$ and $dV = r^2\sin\theta\,dr\,d\theta\,d\varphi$ , and the ball minus a null set is the image of the box $(0,R)\times(0,\pi)\times(-\pi,\pi)$ , so

I = \int_{-\pi}^{\pi}\!\!\int_0^{\pi}\!\!\int_0^{R} r^4\sin^3\theta \, dr\, d\theta\, d\varphi = 2\pi \cdot \frac{R^5}{5}\cdot \int_0^\pi \sin^3\theta\,d\theta = 2\pi\cdot\frac{R^5}{5}\cdot\frac43 = \frac{8\pi R^5}{15},

using $\int_0^\pi\sin^3\theta\,d\theta = \int_{-1}^{1}(1-c^2)\,dc = 4/3$ with $c = \cos\theta$ . For $R = 1$ this is $8\pi/15 \approx 1.675516$ . The three integrals separated because the integrand factorised and the limits were constants — which is the whole point of choosing coordinates that turn the region into a box.

Application(Moments of inertia, and the Maxwell speed distribution)Original

A uniform solid ellipsoid. Let $E = \{x^2/a^2 + y^2/b^2 + z^2/c^2 \le 1\}$ carry total mass $M$ at constant density $\rho = M/\operatorname{vol}(E)$ . Its moment of inertia about the $z$ -axis is $I_z = \rho\int_E (x^2+y^2)\,dV$ . Two changes of variables do the whole computation.

First the linear map $G(u) = (au_1, bu_2, cu_3)$ , with constant $\det DG = abc$ , which carries the closed unit ball $B$ onto $E$ ; the volume rule gives $\operatorname{vol}(E) = abc\operatorname{vol}(B) = \tfrac43\pi abc$ , and

\int_E (x^2+y^2)\,dV = abc\int_B \bigl(a^2u_1^2 + b^2u_2^2\bigr)\,dV .

Then spherical polars on $B$ : $\int_B \lvert u\rvert^2 dV = \int_0^1 r^2\cdot 4\pi r^2\,dr = 4\pi/5$ , and since $\int_Bu_1^2 = \int_Bu_2^2 = \int_Bu_3^2$ by the symmetry of $B$ under permuting coordinates, each equals $\tfrac13\cdot4\pi/5 = 4\pi/15$ . Hence

I_z = \frac{3M}{4\pi abc}\cdot abc\,(a^2+b^2)\,\frac{4\pi}{15} = \frac{M\,(a^2+b^2)}{5},

which for $a=b=c=R$ is the familiar $\tfrac25 MR^2$ of a uniform sphere. Nothing but the Jacobian factors $abc$ and $r^2\sin\theta$ was needed; attempting the same integral in Cartesian coordinates means three nested square-root limits.

The Maxwell speed distribution. In the kinetic theory of gases the velocity $v \in \mathbb{R}^3$ of a molecule of mass $m$ at temperature $T$ has probability density

f(v) = \left(\frac{m}{2\pi kT}\right)^{3/2} e^{-m\lvert v\rvert^2/(2kT)} ,

which is largest at $v = 0$ . The distribution of the speed $\lvert v\rvert$ is obtained by changing to spherical coordinates in velocity space, where $d^3v = v^2\sin\theta\,dv\,d\theta\,d\varphi$ : since $f$ depends only on $v = \lvert v\rvert$ , the two angular integrations give the factor $\int_{-\pi}^{\pi}\int_0^\pi \sin\theta\,d\theta\,d\varphi = 4\pi$ and

F(v) = 4\pi v^2\left(\frac{m}{2\pi kT}\right)^{3/2} e^{-mv^2/(2kT)} .

The famous $v^2$ in front is precisely the Jacobian factor, with the surface of the sphere of radius $v$ folded into it. It is also what makes $F(0) = 0$ and puts the most probable speed at

\frac{dF}{dv} = 0 \;\Longrightarrow\; 2v = \frac{mv^2}{kT}\,v \;\Longrightarrow\; v_{\mathrm{mp}} = \sqrt{\frac{2kT}{m}} ,

strictly positive, even though the most probable velocity vector is zero. The apparent paradox is entirely a matter of the volume element: few molecules have small speed because there is little room in velocity space at small $\lvert v\rvert$ .

Theorem

Differentiation under the integral sign

Adapted

Add to PDF

If $f$ and all its first partial derivatives in the parameters $x \in U \subseteq \mathbb{R}^n$ are continuous on $U \times [a,b]$ , then $F(x) = \int_a^b f(x,t)\,dt$ is of class $C^1$ with $\partial F/\partial x_i = \int_a^b \partial f/\partial x_i\,dt$ ; with differentiable limits of integration the derivative acquires the two boundary terms of the Leibniz rule.

Prerequisites: Uniform convergence, Partial derivatives and the chain rule, Compactness

Theorem(Fixed limits of integration)AdaptedAdapted from Lee (2013), App. C, Thm C.14, p. 648

Let $U \subseteq \mathbb{R}^n$ be open, let $a < b$ be real, and let $f : U \times [a,b] \to \mathbb{R}$ be continuous. Suppose that for each $i = 1,\dots,n$ the partial derivative $\partial f/\partial x_i$ exists at every point of $U \times [a,b]$ and is continuous on $U \times [a,b]$ . Define

F(x) = \int_a^b f(x,t)\, dt \qquad (x \in U).

Then $F$ is of class $C^1$ on $U$ and, for every $x \in U$ and every $i$ ,

\frac{\partial F}{\partial x_i}(x) = \int_a^b \frac{\partial f}{\partial x_i}(x,t)\, dt .

For $n = 1$ this reads $F'(x) = \int_a^b \partial f/\partial x\,(x,t)\, dt$ .

The following corollary is ours, not part of Lee's theorem, which concludes only $C^1$ ; it is proved in Step 7 below. Iteration. If $f$ and all its partial derivatives in $x$ up to order $k$ are continuous on $U \times [a,b]$ , then $F$ is of class $C^k$ and each such derivative may be taken inside the integral.

Theorem(Leibniz rule with variable limits)AdaptedAdapted from Stakgold & Holst (2011), §1.2, Eq. (1.2.11), p. 64; Lee (2013), App. C, Thm C.14, p. 648

Let $I, J \subseteq \mathbb{R}$ be open intervals, let $f : I \times J \to \mathbb{R}$ be continuous with $\partial f/\partial x$ existing and continuous on $I \times J$ , and let $\alpha, \beta : I \to J$ be differentiable. Define

F(x) = \int_{\alpha(x)}^{\beta(x)} f(x,t)\, dt \qquad (x \in I).

Then $F$ is differentiable on $I$ and

F'(x) = f\bigl(x, \beta(x)\bigr)\, \beta'(x) \;-\; f\bigl(x, \alpha(x)\bigr)\, \alpha'(x) \;+\; \int_{\alpha(x)}^{\beta(x)} \frac{\partial f}{\partial x}(x,t)\, dt .

If in addition $\alpha$ and $\beta$ are of class $C^1$ , then $F$ is of class $C^1$ . The two boundary terms are the rate at which the interval of integration sweeps across the integrand at each end; the integral term is the previous statement with the limits frozen.

NoteOriginal

What it means. The question answered here is when a limit (the derivative in $x$ ) may be moved through an integral (over $t$ ). Continuity of $\partial f/\partial x$ on the closed strip makes the difference quotients in $x$ converge uniformly in $t \in [a,b]$ , which is one sufficient condition for a limit to pass under a sign of integration (uniform convergence); that is the whole content of the fixed-limit theorem. The variable-limit version is then the chain rule applied to $(x,u,v) \mapsto \int_u^v f(x,t)\,dt$ , whose three partial derivatives are, by the fundamental theorem of calculus and the fixed-limit theorem, $\int_u^v \partial f/\partial x\,dt$ , $-f(x,u)$ and $f(x,v)$ . On an unbounded interval of integration, or with an unbounded integrand, continuity is not enough and the hypothesis has to be a domination condition instead; see dominated convergence.

Idea of the proof(Fixed limits)Original

The derivative of $F$ at $x$ is a limit of averages of difference quotients of $f$ , and the integral of $\partial f/\partial x$ is the average of a limit. So the whole question is whether one limit, taken in $x$ , may be exchanged with an average taken in $t$ — and the one hypothesis that always licenses such an exchange over a bounded interval is that the convergence be uniform in $t$ .

That uniformity is bought by compactness. The mean value theorem rewrites each difference quotient in $x$ as a value of $\partial f/\partial x$ at an unknown nearby point; continuity of $\partial f/\partial x$ on a compact strip around the segment $\{x\}\times[a,b]$ makes that value close to $\partial f/\partial x(x,t)$ by an amount that does not depend on $t$ . Integrating an error smaller than $\varepsilon$ over an interval of length $b-a$ costs at most $(b-a)\varepsilon$ , which is the whole estimate. The same compactness argument, used once more, shows the resulting partial derivatives are continuous, which is what upgrades "differentiable" to $C^1$ .

Proof(Fixed limits)Original

Lee states the theorem on p. 648 without proving it: the paragraph that introduces it, at the foot of p. 647, says only that a proof may be found in [Apo74], [Rud76] or [Str00]. So the argument below is ours; it is the standard compactness route sketched in the remark.

Step 0 ( $F$ is defined). For each fixed $x \in U$ the map $t \mapsto f(x,t)$ is continuous on $[a,b]$ , hence Riemann integrable, so $F(x)$ is a well-defined real number.

Step 1 (localise). Fix $x_0 \in U$ and an index $i$ . Since $U$ is open there is $r > 0$ with the closed ball $K = \{x \in \mathbb{R}^n : \lvert x - x_0\rvert \le r\}$ contained in $U$ . The set $K \times [a,b] \subseteq \mathbb{R}^{n+1}$ is closed and bounded, hence compact by Heine-Borel.

Step 2 (uniform continuity of the derivative). By hypothesis $\partial f/\partial x_i$ is continuous on $U \times [a,b]$ , hence on the compact set $K \times [a,b]$ , hence uniformly continuous there (compactness, "what continuity gains on a compact space", part 2). So, given $\varepsilon > 0$ , there is $\delta$ with $0 < \delta \le r$ such that

\left\lvert \frac{\partial f}{\partial x_i}(y,s) - \frac{\partial f}{\partial x_i}(y',s) \right\rvert < \varepsilon \qquad\text{whenever } y, y' \in K,\ \lvert y - y'\rvert < \delta,\ s \in [a,b] .

This is the only place the continuity of $\partial f/\partial x_i$ is used; its mere existence at every point of $U \times [a,b]$ is used again in Step 3, to apply the mean value theorem. Existence alone would not give the display.

Step 3 (mean value theorem in the $i$ -th variable). Let $0 < \lvert h \rvert < \delta$ , so that $x_0 + h e_i \in K \subseteq U$ , and fix $t \in [a,b]$ . The one-variable function $s \mapsto f(x_0 + s e_i, t)$ is defined for $\lvert s \rvert \le r$ , is continuous there (as $f$ is continuous) and is differentiable with derivative $\partial f/\partial x_i (x_0 + s e_i, t)$ , because that partial derivative exists at every point of $U \times [a,b]$ . The mean value theorem therefore supplies a number $\theta \in (0,1)$ , depending on $t$ and $h$ , with

\frac{f(x_0 + h e_i, t) - f(x_0, t)}{h} = \frac{\partial f}{\partial x_i}\bigl(x_0 + \theta h e_i, t\bigr) .

Step 4 (the difference quotients converge uniformly in $t$ ). The point $x_0 + \theta h e_i$ lies in $K$ and is at distance $\lvert \theta h\rvert < \delta$ from $x_0$ , so Step 2 applied to $y = x_0 + \theta h e_i$ and $y' = x_0$ turns Step 3 into

\left\lvert \frac{f(x_0 + h e_i, t) - f(x_0, t)}{h} - \frac{\partial f}{\partial x_i}(x_0, t) \right\rvert < \varepsilon \qquad\text{for every } t \in [a,b] \text{ and every } 0 < \lvert h\rvert < \delta .

The bound is independent of $t$ : the difference quotients converge to $\partial f/\partial x_i(x_0, \cdot)$ uniformly on $[a,b]$ (uniform convergence).

Step 5 (integrate the estimate). By linearity of the integral and the elementary bound $\lvert\int_a^b g\rvert \le (b-a)\sup_{[a,b]}\lvert g\rvert$ , Step 4 gives, for every $0 < \lvert h \rvert < \delta$ ,

\left\lvert \frac{F(x_0 + h e_i) - F(x_0)}{h} - \int_a^b \frac{\partial f}{\partial x_i}(x_0,t)\,dt \right\rvert \le (b-a)\,\varepsilon .

As $\varepsilon > 0$ was arbitrary, the limit of the left-hand difference quotient as $h \to 0$ exists and equals the integral: $\partial F/\partial x_i(x_0)$ exists and is $\int_a^b \partial f/\partial x_i (x_0,t)\,dt$ . This is the displayed formula.

Nothing in Steps 1 to 5 singled out the point $x_0$ or the index $i$ : the hypotheses hold at every point of $U$ and for every $i$ , so the argument may be run again at any $x \in U$ . Hence every first partial derivative of $F$ exists at every point of $U$ and is given there by the displayed formula. That is what makes Step 6 legitimate as stated.

Step 6 (the partial derivatives are continuous). Fix again a point $x_0 \in U$ and the ball $K \subseteq U$ , $\varepsilon$ and $\delta$ of Steps 1 and 2. By the paragraph just above, $\partial F/\partial x_i$ exists at every point of $K$ and is the integral of $\partial f/\partial x_i$ there. So for $x, x' \in K$ with $\lvert x - x'\rvert < \delta$ , Step 2 bounds the integrands and the same elementary estimate gives

\left\lvert \frac{\partial F}{\partial x_i}(x) - \frac{\partial F}{\partial x_i}(x') \right\rvert \le \int_a^b \left\lvert \frac{\partial f}{\partial x_i}(x,t) - \frac{\partial f}{\partial x_i}(x',t) \right\rvert dt \le (b-a)\,\varepsilon .

Since $\varepsilon$ was arbitrary and $\delta$ depends only on $\varepsilon$ , this says $\partial F/\partial x_i$ is continuous at $x_0$ ; and $x_0 \in U$ was arbitrary. Every first partial derivative of $F$ therefore exists on $U$ and is continuous on $U$ , which is exactly the statement that $F$ is of class $C^1$ (partial derivatives and the chain rule).

Step 7 (iteration). Suppose $f$ and all its $x$ -partials up to order $k$ are continuous on $U \times [a,b]$ . Steps 1 to 6 apply to $f$ and give $\partial F/\partial x_i = \int_a^b \partial f/\partial x_i\,dt$ ; the integrand is a function to which the hypotheses again apply, with $k$ replaced by $k-1$ . Induction on $k$ therefore gives $F \in C^k$ and the interchange for every derivative of order at most $k$ . $\square$

Idea of the proof(Variable limits)Original

There is nothing new to prove. Freeze the two limits of integration as independent variables, so that the single function $\Phi(x,u,v) = \int_u^v f(x,t)\,dt$ contains all the information, and note that $F$ is the composite of $\Phi$ with the curve $x \mapsto (x, \alpha(x), \beta(x))$ .

The three partial derivatives of $\Phi$ are known separately: the derivative in $x$ is the fixed-limit theorem just proved, and the derivatives in $u$ and in $v$ are the fundamental theorem of calculus, with a minus sign at the lower end. The chain rule then assembles them, and the three terms of the Leibniz formula are just those three partial derivatives weighted by $1$ , $\alpha'(x)$ and $\beta'(x)$ .

Proof(Variable limits)Original

Write $\Phi(x,u,v) = \int_u^v f(x,t)\,dt$ for $(x,u,v) \in I \times J \times J$ , with the usual orientation convention $\int_u^v = -\int_v^u$ . Because $J$ is an interval, the closed segment between $u$ and $v$ lies in $J$ , so the integrand is defined and continuous there and $\Phi$ is well defined. By construction $F(x) = \Phi\bigl(x, \alpha(x), \beta(x)\bigr)$ , the composition being legitimate because $\alpha$ and $\beta$ take values in $J$ .

Step 1 (the two derivatives in the limits). Fix $x$ and $u$ . The map $v \mapsto \int_u^v f(x,t)\,dt$ is an indefinite integral of the continuous function $t \mapsto f(x,t)$ , so by the fundamental theorem of calculus $\partial \Phi/\partial v = f(x,v)$ . Exchanging the roles of the endpoints and using the orientation convention, $\partial \Phi/\partial u = -f(x,u)$ .

Step 2 (the derivative in the parameter). Fix $u, v \in J$ and choose a closed interval $[c,d] \subseteq J$ containing both $u$ and $v$ , and an open interval $I' \subseteq I$ containing the point $x$ of interest. Then $f$ and $\partial f/\partial x$ are continuous on $I' \times [c,d]$ , so the fixed-limit statement applies with $[c,d]$ in place of $[a,b]$ and gives

\frac{\partial \Phi}{\partial x}(x,u,v) = \int_u^v \frac{\partial f}{\partial x}(x,t)\,dt ,

the sign convention carrying over from $\int_u^v$ to $\int_u^v \partial f/\partial x$ unchanged.

Step 3 ( $\Phi$ is of class $C^1$ ). The two derivatives of Step 1 are continuous in $(x,u,v)$ because $f$ is continuous. For the third, let $(x,u,v)$ and $(x',u',v')$ both lie in $I' \times [c,d]^2$ with $I'$ and $[c,d]$ as in Step 2, and let $M = \max \lvert \partial f/\partial x\rvert$ over the compact set $\overline{I''} \times [c,d]$ for a compact subinterval $\overline{I''} \subseteq I'$ containing $x$ and $x'$ . Then

\left\lvert \frac{\partial \Phi}{\partial x}(x,u,v) - \frac{\partial \Phi}{\partial x}(x',u',v') \right\rvert \le \left\lvert \int_u^v \Bigl( \frac{\partial f}{\partial x}(x,t) - \frac{\partial f}{\partial x}(x',t) \Bigr) dt \right\rvert + M\bigl( \lvert u - u'\rvert + \lvert v - v'\rvert \bigr) ,

because changing the endpoints alters the integral only over two intervals of lengths $\lvert u - u'\rvert$ and $\lvert v-v'\rvert$ . The first term tends to $0$ as $x' \to x$ by Step 6 of the previous proof, and the second tends to $0$ as $(u',v') \to (u,v)$ . So all three partial derivatives of $\Phi$ exist and are continuous on $I \times J \times J$ , that is, $\Phi$ is of class $C^1$ , and in particular differentiable at every point (partial derivatives and the chain rule).

Step 4 (chain rule). The curve $\gamma(x) = \bigl(x, \alpha(x), \beta(x)\bigr)$ is differentiable on $I$ , since $\alpha$ and $\beta$ are differentiable by hypothesis, with $\gamma'(x) = \bigl(1, \alpha'(x), \beta'(x)\bigr)$ . As $\Phi$ is differentiable, the chain rule applies to $F = \Phi \circ \gamma$ and gives, at every $x \in I$ ,

F'(x) = \frac{\partial \Phi}{\partial x} \cdot 1 + \frac{\partial \Phi}{\partial u}\,\alpha'(x) + \frac{\partial \Phi}{\partial v}\,\beta'(x) ,

all three partials evaluated at $\gamma(x)$ . Substituting Steps 1 and 2 produces

F'(x) = \int_{\alpha(x)}^{\beta(x)} \frac{\partial f}{\partial x}(x,t)\,dt - f\bigl(x,\alpha(x)\bigr)\alpha'(x) + f\bigl(x,\beta(x)\bigr)\beta'(x) ,

which is the Leibniz formula.

Step 5 (the $C^1$ refinement). If in addition $\alpha, \beta \in C^1$ , then $\gamma$ is continuous with continuous derivative, and the right-hand side above is a continuous function of $x$ : the three partials of $\Phi$ are continuous by Step 3 and $\gamma$ , $\alpha'$ , $\beta'$ are continuous. A differentiable function with continuous derivative is of class $C^1$ , so $F \in C^1(I)$ . $\square$

Method(Evaluating an integral by introducing a parameter)Original

A definite integral with no visible parameter can often be computed by creating one and differentiating.

Embed the integral in a family. Write the wanted number as $F(x_0)$ for some $F(x) = \int_a^b f(x,t)\,dt$ , choosing the parameter $x$ so that $\partial f/\partial x$ is simpler than $f$ — typically by inserting a factor $e^{-xt}$ , $t^{x}$ or $\cos(xt)$ , or by promoting a constant already present to a variable.
Check the hypotheses on a region containing $x_0$ . Find an open $I \ni x_0$ on which $f$ and $\partial f/\partial x$ are continuous on $I \times [a,b]$ , endpoints included. If the $t$ -interval is infinite or $f$ is unbounded, this step must instead produce a function $g$ , integrable over the $t$ -range, with $\lvert \partial f/\partial x (x,t)\rvert \le g(t)$ for every $x \in I$ and every $t$ , so that dominated convergence applies.
Differentiate under the sign to get $F'(x) = \int_a^b \partial f/\partial x\,(x,t)\,dt$ , and evaluate that integral in closed form. If it is still hard, differentiate again, or combine with a reversal of the order of integration, the companion technique for a double integral (Arfken, Weber & Harris (2013), §1.10, Example 1.10.6, pp. 70–71).
Solve the resulting equation for $F$ . Usually $F' = \varphi(x)$ is integrated directly; sometimes step 3 yields a first-order differential equation for $F$ , for example $F' = c(x)F$ .
Fix the constant of integration from a value of $x$ at which $F$ is known, often an endpoint or a limit $x \to \infty$ ; frequently it is a value at which the integrand vanishes identically in $t$ , so that $F$ there is $0$ . The constant cannot be recovered from the differentiated equation alone, and this step is where the answer is pinned down.
Evaluate at $x_0$ , checking that $x_0$ lies in the interval $I$ on which steps 2–5 were justified.

How it fails. Omitting step 2 is the usual error: on a strip that is not closed, or with an infinite interval of integration, continuity of $f$ and $\partial f/\partial x$ does not justify step 3, and the identity can be false (see the remark). Omitting step 5 produces an answer right up to an additive constant. If $f$ or $\partial f/\partial x$ is discontinuous at some $t = t_0$ in $[a,b]$ , split the integral at $t_0$ first and use the Leibniz rule on each piece, as Stakgold and Holst do for a Green's function whose $t$ -derivative jumps at $t = x$ .

RemarkOriginal

Where the hypotheses are used, and what fails without them.

The continuity must hold on a closed strip $U \times [a,b]$ . Let $f(x,t) = \frac{2x^3}{(x^2+t^2)^2} ,$ which is continuous, together with $\partial f/\partial x$ , on $\mathbb{R}\times(0,\infty)$ — but on no strip $U \times [0,b]$ with $0 \in U$ , since $f(x,0) = 2/x$ has no limit as $x \to 0$ . With the full range $t \in (0,\infty)$ , $\int_0^\infty dt/(x^2+t^2)^2 = \pi/(4x^3)$ for $x > 0$ , so $F(x) = \int_0^\infty f(x,t)\, dt = \frac{\pi}{2}\ (x>0), \qquad F(0) = 0, \qquad F(x) = -\frac{\pi}{2}\ (x<0),$ and $F$ is not even continuous at $0$ , let alone differentiable; meanwhile $f(x,t)$ vanishes to third order in $x$ for each fixed $t > 0$ , so $\partial f/\partial x(0,t) = 0$ and the differentiated integral exists and equals $0$ . One side of the formula does not exist at all.
It is the open end at $t = 0$ , not the length of the $t$ -interval, that fails here. Shortening the range to $(0,1]$ changes nothing: $\int_0^1 f(x,t)\,dt = x/(x^2+1) + \arctan(1/x)$ , with one-sided limits $\pm\pi/2$ at $0$ while the value there is $0$ . Keeping the range infinite but closing it at the bad end repairs the identity: on $\mathbb{R}\times[1,\infty)$ one finds $\int_1^\infty f(x,t)\,dt = \arctan x - x/(1+x^2) = \tfrac23 x^3 + O(x^5)$ , whose derivative at $0$ is $0 = \int_1^\infty \partial f/\partial x(0,t)\,dt$ . What goes wrong in the first bullet is that the mass of $\partial f/\partial x$ escapes to $t = 0$ as $x \to 0$ : no function of $t$ integrable near $0$ dominates $\partial f/\partial x(x,t)$ for all small $x$ .
Existence of $\partial f/\partial x$ is not enough; its continuity is what is used. The fixed-limit statement is proved by controlling $(f(x+h,t)-f(x,t))/h$ uniformly in $t \in [a,b]$ , via the mean value theorem and the continuity of $\partial f/\partial x$ on the compact strip $K \times [a,b]$ for a compact $K \ni x$ . Drop that continuity and the difference quotients need not converge uniformly.
The boundary terms need only differentiability of the limits. In the Leibniz rule $\alpha$ and $\beta$ are required to be differentiable, not $C^1$ ; $F$ is then differentiable but perhaps not $C^1$ . They must also take values in $J$ , so that $f(x,\alpha(x))$ and $f(x,\beta(x))$ are defined, and $J$ must be an interval, so that the whole segment between $\alpha(x)$ and $\beta(x)$ lies in it.
A jump in the integrand across the diagonal. Stakgold and Holst apply the rule to a Green's-function integral $\int_0^1 g(x,\xi)\,q(\xi)\,d\xi$ — they call the source term $f$ , which we rename $q$ here to keep $f$ for the integrand of the statements above — where $\partial g/\partial x$ jumps at $\xi = x$ . The integral must first be split at $\xi = x$ , after which each piece has a continuous $x$ -derivative and a variable limit. The boundary terms from the two pieces cancel on the first differentiation precisely because $g$ is continuous across $\xi = x$ , and fail to cancel on the second, which is what produces the delta-like source term.
Conventions. Some texts call the fixed-limit identity "Leibniz's rule" and others reserve that name for the variable-limit formula with its two boundary terms, which is the usage here. The same name is also attached, in an unrelated context, to the formula for the $n$ th derivative of a product.
The general hypothesis. For an infinite interval or an unbounded integrand the right statement replaces continuity of $\partial f/\partial x$ by a bound $\lvert \partial f/\partial x(x,t)\rvert \le g(t)$ valid for all $x$ in a neighbourhood, with $g$ integrable over the $t$ -range; that belongs to dominated convergence. Domination is what replaces the compactness of $[a,b]$ used in the proof above.

Example(A parameter inserted where there was none)Original

Evaluate

N = \int_0^1 \frac{\ln(1+t)}{1+t^2}\, dt .

There is no parameter in sight, so we create one.

Embed in a family. Put $f(x,t) = \ln(1+xt)/(1+t^2)$ and $F(x) = \int_0^1 f(x,t)\,dt$ , so that $N = F(1)$ and $F(0) = 0$ . The point of this choice is that differentiating in $x$ removes the logarithm.
Check the hypotheses. Take $I = (-1,\infty)$ . For $x \in I$ and $t \in [0,1]$ we have $1 + xt \ge \min\{1, 1+x\} > 0$ , so $f$ is continuous on $I \times [0,1]$ , and so is $\frac{\partial f}{\partial x}(x,t) = \frac{t}{(1+xt)(1+t^2)} .$ The $t$ -interval is closed and bounded, so the fixed-limit theorem applies on the whole of $I$ and $F \in C^1(I)$ .
Differentiate under the sign and split the integrand. Partial fractions in $t$ give $\frac{t}{(1+xt)(1+t^2)} = \frac{1}{1+x^2}\left( \frac{-x}{1+xt} + \frac{t+x}{1+t^2} \right),$ which is checked by clearing denominators: $-x(1+t^2) + (t+x)(1+xt) = t + x^2 t = t(1+x^2)$ .
Evaluate the three pieces. For $x \neq 0$ , $\int_0^1 dt/(1+xt) = \ln(1+x)/x$ , while $\int_0^1 t\,dt/(1+t^2) = \tfrac12 \ln 2$ and $\int_0^1 dt/(1+t^2) = \pi/4$ . Hence $F'(x) = \frac{1}{1+x^2}\left( -\ln(1+x) + \tfrac12 \ln 2 + \frac{\pi x}{4} \right),$ and this also holds at $x = 0$ , where both sides equal $\tfrac12\ln 2$ .
Integrate from $0$ to $1$ . $F$ is $C^1$ on $I \supseteq [0,1]$ , so $F(1) = F(0) + \int_0^1 F'(x)\,dx$ with $F(0) = 0$ : $F(1) = -\int_0^1 \frac{\ln(1+x)}{1+x^2}\,dx + \frac{\ln 2}{2}\int_0^1 \frac{dx}{1+x^2} + \frac{\pi}{4}\int_0^1 \frac{x\,dx}{1+x^2} .$
Recognise the integral on the right. The first integral is $N = F(1)$ again — only the name of the variable has changed. The other two are $\pi/4$ and $\tfrac12\ln 2$ , so $F(1) = -F(1) + \frac{\ln 2}{2}\cdot\frac{\pi}{4} + \frac{\pi}{4}\cdot\frac{\ln 2}{2} = -F(1) + \frac{\pi \ln 2}{4} .$
Solve. $2F(1) = \pi \ln 2/4$ , so $N = \frac{\pi \ln 2}{8} = 0.2721982613\ldots$ Numerical quadrature of the original integral gives $0.2721982613$ , which is the check. Note that step 5, fixing the value of $F$ at a point where it is known, is what pins the answer down; step 4 alone determines $N$ only up to a constant.

Example(Variable limits: a Volterra integral and the equation it solves)Original

Let $q : \mathbb{R} \to \mathbb{R}$ be continuous and set

u(x) = \int_0^x \sin(x - s)\, q(s)\, ds .

We show that $u$ solves $u'' + u = q$ with $u(0) = u'(0) = 0$ , and then compute it for $q(s) = s$ .

Hypotheses. Here $f(x,s) = \sin(x-s)q(s)$ is continuous on $\mathbb{R}^2$ , and $\partial f/\partial x = \cos(x-s)q(s)$ is too; the limits $\alpha(x) = 0$ and $\beta(x) = x$ are of class $C^1$ with $\alpha' = 0$ , $\beta' = 1$ . The Leibniz rule applies, with $I = J = \mathbb{R}$ .
First differentiation. The boundary term at the lower limit carries the factor $\alpha' = 0$ , and the one at the upper limit is $f(x,x)\cdot 1 = \sin(0)q(x) = 0$ . So only the integral term survives: $u'(x) = \int_0^x \cos(x-s)\, q(s)\, ds .$ In particular $u(0) = u'(0) = 0$ , both integrals being over an interval of zero length.
Second differentiation. Apply the rule again, now to the integrand $\cos(x-s)q(s)$ , whose $x$ -derivative is $-\sin(x-s)q(s)$ . This time the upper boundary term is $\cos(0)q(x) = q(x)$ and does not vanish: $u''(x) = q(x) - \int_0^x \sin(x-s)\,q(s)\,ds = q(x) - u(x) .$ Hence $u'' + u = q$ . The whole content of the computation is that the boundary term is invisible at the first differentiation and decisive at the second.
A closed form. Take $q(s) = s$ . Substituting $\sigma = x - s$ , $u(x) = \int_0^x (x - \sigma)\sin\sigma\, d\sigma = x\bigl(1 - \cos x\bigr) - \bigl(\sin x - x\cos x\bigr) = x - \sin x .$
Check. $u'' = \sin x$ , so $u'' + u = \sin x + (x - \sin x) = x = q(x)$ , and $u(0) = u'(0) = 0$ . The value at $x = 1$ is $u(1) = 1 - \sin 1 = 0.1585290152\ldots,$ and numerical quadrature of $\int_0^1 \sin(1-s)\,s\,ds$ returns $0.1585290152$ .

ApplicationOriginal

Differentiating an integral in a parameter, and in its limits, is one of the standard moves of theoretical physics.

Impulse response and Duhamel's formula. For a mass on a spring driven from rest, $m\ddot{x} + kx = f(t)$ with $x(0) = \dot{x}(0) = 0$ and $\omega = \sqrt{k/m}$ , the solution is $x(t) = \frac{1}{m\omega}\int_0^t \sin\bigl(\omega(t-s)\bigr) f(s)\,ds ,$ and the verification is Example 2 word for word, with the factor $1/(m\omega)$ and $\omega$ inserted: the first differentiation loses its boundary term because the impulse response $\sin(\omega\tau)/(m\omega)$ vanishes at $\tau = 0$ (the system starts at rest), and the second keeps it because the derivative of that response is $1/m$ at $\tau = 0$ (an impulse delivers momentum). That surviving boundary term is the force at the present instant. For $m = 1$ , $\omega = 2$ and a step force $f \equiv 1$ the formula gives $x(t) = (1 - \cos 2t)/4$ , which oscillates about the static displacement $f/k = 1/4$ ; at $t = 1$ it equals $0.3540367091$ .
Green's functions with a kink. The same rule, applied twice, is how one checks that $u(x) = \int_0^1 g(x,\xi)\,q(\xi)\,d\xi$ solves $-u'' = q$ on $(0,1)$ with $u(0) = u(1) = 0$ . Because $\partial g/\partial x$ jumps across $\xi = x$ , the integral must first be split at $\xi = x$ , which turns a fixed-limit integral into two variable-limit ones. Stakgold and Holst carry this out (Stakgold & Holst (2011), §1.2, pp. 63–64): the boundary terms of the two pieces cancel on the first differentiation, because $g$ is continuous across the diagonal, and fail to cancel on the second, where their difference is exactly the unit jump of $\partial g/\partial x$ and delivers $-q(x)$ . Every delta-function source in a Green's-function calculation appears in this way.
Generating functions in statistical mechanics and quantum theory. Thermodynamic averages are parameter derivatives of integrals: from $Z(\beta) = \int e^{-\beta H}\,d\Gamma$ one gets $\langle H \rangle = -\partial(\ln Z)/\partial\beta$ , and differentiating again gives the energy variance, hence the heat capacity. Moments of any distribution come the same way from its moment generating function, and the Hellmann-Feynman theorem $\partial E/\partial\lambda = \langle \psi \vert \partial H/\partial\lambda \vert \psi\rangle$ is the same interchange applied to $E(\lambda) = \int \overline{\psi}\,H_\lambda \psi$ . In each case the hypothesis that has to be checked is the one of step 2 of the method: on an infinite range of integration continuity is not enough, and one needs a dominating integrable bound valid throughout a neighbourhood of the parameter value.

Method

Gaussian integrals

Adapted

Add to PDF

Collects the one-dimensional Gaussian integral and its moments, the effect of completing the square in the exponent, and the $n$ -dimensional formula $\int \exp(-\tfrac12 x^{\mathsf T}Ax + b^{\mathsf T}x)\,d^n x = (2\pi)^{n/2}(\det A)^{-1/2}\exp(\tfrac12 b^{\mathsf T}A^{-1}b)$ for symmetric positive definite $A$ , with the recipe for reducing an exponent to that form.

Prerequisites: Multiple integrals and change of variables, Differentiation under the integral sign, Quadratic forms and the Rayleigh quotient, Spectral theorem for Hermitian matrices, Fubini and Tonelli theorems, Cauchy's theorem

Definition(Notation)AdaptedAdapted from Arfken, Weber & Harris (2013), §13.1, Eq. (13.5), p. 600 and Eq. (13.8), p. 601

Throughout, $\alpha > 0$ and $n \in \mathbb{N}$ .

The gamma function is defined for $\operatorname{Re} z > 0$ by the Euler integral $\Gamma(z) = \int_0^\infty e^{-t} t^{z-1}\,dt$ , which converges there; it satisfies $\Gamma(z+1) = z\Gamma(z)$ , $\Gamma(k) = (k-1)!$ for $k \in \mathbb{N}$ , and $\Gamma(1/2) = \sqrt{\pi}$ . Only these facts are used below; for the full theory see the Gamma function.
The double factorial of an odd number is $(2m-1)!! = 1 \cdot 3 \cdot 5 \cdots (2m-1)$ for $m \ge 1$ , with the convention $(-1)!! = 1$ so that the formulas below hold at $m = 0$ . It is related to the gamma function by $\Gamma(m + \tfrac12) = (2m-1)!!\,\sqrt{\pi}/2^{m}$ .
A real symmetric $n \times n$ matrix $A$ is positive definite if $x^{\mathsf T} A x > 0$ for every $x \in \mathbb{R}^n \setminus \{0\}$ , equivalently if all its eigenvalues are positive (quadratic forms).
$d^n x$ denotes Lebesgue measure on $\mathbb{R}^n$ , and an integral written $\int_{\mathbb{R}^n}$ is over all of $\mathbb{R}^n$ . Only the symmetric part of a matrix in a quadratic exponent is ever visible, since $x^{\mathsf T} A x = x^{\mathsf T}\bigl(\tfrac12(A + A^{\mathsf T})\bigr)x$ for every $x$ , so $A$ may be assumed symmetric without loss.

Method(The one-dimensional Gaussian integral and its moments)AdaptedAdapted from Riley, Hobson & Bence (2006), §6.4.2, pp. 202–204; Arfken, Weber & Harris (2013), §13.1, Eq. (13.8), p. 601 and Exercise 13.1.10, p. 608

Let $\alpha > 0$ .

The integral. $\int_{-\infty}^{\infty} e^{-\alpha x^2}\, dx = \sqrt{\frac{\pi}{\alpha}} , \qquad\text{equivalently}\qquad \int_{-\infty}^{\infty} e^{-x^2/(2\sigma^2)}\, dx = \sigma\sqrt{2\pi} \quad (\sigma > 0).$ In gamma-function form the same statement is $\Gamma(1/2) = \sqrt{\pi}$ .
Even moments. For every integer $m \ge 0$ , $\int_{0}^{\infty} x^{2m} e^{-\alpha x^2}\, dx = \frac{\Gamma\!\left(m + \tfrac12\right)}{2\,\alpha^{\,m+1/2}} = \frac{(2m-1)!!}{2^{\,m+1}\alpha^{m}}\sqrt{\frac{\pi}{\alpha}},$ and, the integrand being even, twice that over the whole line: $\int_{-\infty}^{\infty} x^{2m} e^{-\alpha x^2}\, dx = \frac{\Gamma\!\left(m + \tfrac12\right)}{\alpha^{\,m+1/2}} = \frac{(2m-1)!!}{(2\alpha)^{m}}\sqrt{\frac{\pi}{\alpha}} .$
Odd moments. For every integer $m \ge 0$ , $\int_{0}^{\infty} x^{2m+1} e^{-\alpha x^2}\, dx = \frac{m!}{2\,\alpha^{\,m+1}}, \qquad \int_{-\infty}^{\infty} x^{2m+1} e^{-\alpha x^2}\, dx = 0 ,$ the second because the integrand is odd and absolutely integrable. The half-line formula here is part (a) of Arfken's Exercise 13.1.10; part (b) of that exercise is the half-line even moment of part 2 above. The full-line values in parts 2 and 3 are obtained from those two by the parity of the integrand.
Completing the square. For $\alpha > 0$ and any $\beta \in \mathbb{C}$ , $\int_{-\infty}^{\infty} e^{-\alpha x^2 + \beta x}\, dx = \sqrt{\frac{\pi}{\alpha}}\; e^{\beta^2/(4\alpha)} ,$ the integral converging absolutely for every complex $\beta$ . In particular, with $\beta = i b$ for real $b$ , $\int_{-\infty}^{\infty} e^{-\alpha x^2} e^{i b x}\, dx = \sqrt{\frac{\pi}{\alpha}}\; e^{-b^2/(4\alpha)} :$ the Fourier transform of a Gaussian is a Gaussian.

Method(The multidimensional Gaussian integral)AdaptedAdapted from Riley, Hobson & Bence (2006), §30.13, Eqs. (30.144)–(30.145), p. 1206 and §30.15.2, Eq. (30.149), pp. 1209–1210

Let $A$ be a real symmetric positive definite $n \times n$ matrix and let $b \in \mathbb{R}^n$ . Then the integral below converges absolutely and

\int_{\mathbb{R}^n} \exp\!\left( -\tfrac12\, x^{\mathsf T} A x + b^{\mathsf T} x \right) d^n x = \frac{(2\pi)^{n/2}}{\sqrt{\det A}}\; \exp\!\left( \tfrac12\, b^{\mathsf T} A^{-1} b \right).

Special cases and equivalent forms:

With $b = 0$ : $\displaystyle \int_{\mathbb{R}^n} e^{-\frac12 x^{\mathsf T} A x}\, d^n x = (2\pi)^{n/2} (\det A)^{-1/2}$ .
Writing $V = A^{-1}$ , which is also symmetric positive definite, and shifting by $\mu \in \mathbb{R}^n$ : $\int_{\mathbb{R}^n} \exp\!\left( -\tfrac12 (x-\mu)^{\mathsf T} V^{-1} (x-\mu) \right) d^n x = (2\pi)^{n/2} (\det V)^{1/2},$ the normalisation of the multivariate Gaussian density with covariance $V$ .
Without the factor $\tfrac12$ : $\displaystyle \int_{\mathbb{R}^n} e^{- x^{\mathsf T} A x}\, d^n x = \pi^{n/2}(\det A)^{-1/2}$ , from case 1 with $A$ replaced by $2A$ .

If $A$ is symmetric but not positive definite, the integral diverges for every $b$ : writing $A = Q\Lambda Q^{\mathsf T}$ with $Q$ orthogonal and substituting $x = Qy$ , with $c = Q^{\mathsf T}b$ , the integrand factorises and

\int_{\mathbb{R}^n} \exp\!\left( -\tfrac12 x^{\mathsf T}Ax + b^{\mathsf T}x \right) d^n x = \prod_{i=1}^{n} \int_{-\infty}^{\infty} e^{-\lambda_i y_i^2/2 + c_i y_i}\, dy_i ,

in which every factor is positive and any factor with $\lambda_i \le 0$ is $+\infty$ .

NoteOriginal

What it means. Everything here comes from two facts: the single number $\int e^{-x^2}dx = \sqrt{\pi}$ , which is not obtainable by antidifferentiation and is got instead by squaring the integral and passing to plane polar coordinates (hence this unit's place after change of variables); and the fact that a positive definite quadratic form becomes a sum of squares in suitable coordinates, which turns the $n$ -dimensional integral into a product of $n$ copies of the first fact. Every other formula is bookkeeping: completing the square removes a linear term at the cost of a constant, and differentiating the answer with respect to $\alpha$ or $b$ generates the moments.

Idea of the proof(One dimension)Original

The single value $\int_{\mathbb{R}} e^{-x^2}dx = \sqrt{\pi}$ cannot be reached by antidifferentiation, so it is reached sideways: the square of the integral is an integral over the plane, and in polar coordinates the radial factor $\rho$ supplied by the area element is exactly what makes $e^{-\rho^2}\rho$ elementary to integrate. Because the plane is unbounded, the squeeze between the disc inscribed in a large square and the disc circumscribing it does the work of a convergence argument.

Everything else in the statement is a change of variable away. The substitution $u = \alpha x^2$ turns a half-line moment into the Euler integral for $\Gamma$ , so the moments are gamma values; the odd full-line moments vanish by symmetry, which is legitimate only because the integrals converge absolutely. And a linear term in the exponent is removed by completing the square, which shifts the variable — by a real amount when $\beta$ is real, and by a complex one when it is not, in which case the shift is a deformation of the contour that costs nothing because $e^{-\alpha z^2}$ is entire.

Proof(One dimension)AdaptedAdapted from Riley, Hobson & Bence (2006), §6.4.2, pp. 202–204; Hassani (2013), Example 11.3.8, pp. 352–353

Part 1.

Convergence. The integrand is continuous and positive, and $e^{-\alpha x^2} \le e^{-\alpha \lvert x\rvert}$ for $\lvert x\rvert \ge 1$ , whose integral over $\lvert x \rvert \ge 1$ is finite. So the integral converges absolutely for every $\alpha > 0$ .

Reduction to $\alpha = 1$ . The substitution $x = t/\sqrt{\alpha}$ is a linear bijection of $\mathbb{R}$ with $dx = \alpha^{-1/2}dt$ , so $\int_{\mathbb{R}} e^{-\alpha x^2}dx = \alpha^{-1/2} I$ with $I = \int_{\mathbb{R}} e^{-t^2}dt$ . Only $I$ has to be found.

Truncate and square. Put $I(a) = \int_{-a}^{a} e^{-t^2}dt$ for $a > 0$ ; then $I(a) \to I$ as $a \to \infty$ . Since $e^{-x^2}e^{-y^2} = e^{-(x^2+y^2)}$ and the square $Q_a = [-a,a]^2$ is a product of intervals,

I(a)^2 = \int_{Q_a} e^{-(x^2+y^2)}\, dx\, dy ,

the double integral of a continuous function over a compact rectangle being the iterated one (iterated integrals).

Squeeze between two discs. Write $D_R$ for the closed disc of radius $R$ about the origin. The square $Q_a$ contains $D_a$ and is contained in $D_{a\sqrt2}$ , because its inradius is $a$ and its circumradius is $a\sqrt2$ . The integrand is positive, so

J(a) \le I(a)^2 \le J\bigl(a\sqrt2\bigr), \qquad J(R) = \int_{D_R} e^{-(x^2+y^2)}\,dx\,dy .

This is the step that needs care, and it is why Riley, Hobson and Bence treat the unbounded region explicitly (pp. 203–204).

Evaluate on a disc. In plane polar coordinates the area element is $\rho\,d\rho\,d\varphi$ (change of variables), so

J(R) = \int_0^{2\pi}\!\!\int_0^{R} e^{-\rho^2}\rho\, d\rho\, d\varphi = 2\pi \left[ -\tfrac12 e^{-\rho^2} \right]_0^{R} = \pi\bigl(1 - e^{-R^2}\bigr) .

Conclude. Hence $\pi(1 - e^{-a^2}) \le I(a)^2 \le \pi(1 - e^{-2a^2})$ . Both bounds tend to $\pi$ as $a \to \infty$ , so $I^2 = \pi$ , and $I > 0$ because the integrand is positive; therefore $I = \sqrt{\pi}$ and $\int_{\mathbb{R}}e^{-\alpha x^2}dx = \sqrt{\pi/\alpha}$ . Taking $\alpha = 1/(2\sigma^2)$ gives the second form, $\sqrt{2\pi\sigma^2} = \sigma\sqrt{2\pi}$ .

The gamma form. In $\Gamma(1/2) = \int_0^\infty e^{-t}t^{-1/2}dt$ substitute $t = x^2$ with $x > 0$ : then $dt = 2x\,dx$ and $t^{-1/2} = 1/x$ , so $\Gamma(1/2) = 2\int_0^\infty e^{-x^2}dx = I = \sqrt{\pi}$ , the last equality because the integrand is even.

Part 2 (even moments). On $(0,\infty)$ substitute $u = \alpha x^2$ , an increasing $C^1$ bijection of $(0,\infty)$ onto itself, with $x = (u/\alpha)^{1/2}$ and $dx = \tfrac12 \alpha^{-1/2}u^{-1/2}du$ . Then

\int_0^\infty x^{2m}e^{-\alpha x^2}dx = \int_0^\infty \Bigl(\frac{u}{\alpha}\Bigr)^{m} e^{-u}\,\tfrac12 \alpha^{-1/2}u^{-1/2}\,du = \frac{1}{2\,\alpha^{\,m+1/2}} \int_0^\infty u^{m-1/2}e^{-u}\,du = \frac{\Gamma(m+\tfrac12)}{2\,\alpha^{\,m+1/2}} ,

the last integral being the Euler integral for $\Gamma(m+\tfrac12)$ , convergent because $m - \tfrac12 > -1$ .

The double-factorial form follows from $\Gamma(z+1) = z\Gamma(z)$ applied $m$ times:

\Gamma\bigl(m+\tfrac12\bigr) = \bigl(m-\tfrac12\bigr)\bigl(m-\tfrac32\bigr)\cdots\tfrac12\,\Gamma\bigl(\tfrac12\bigr) = \frac{(2m-1)(2m-3)\cdots 1}{2^{m}}\sqrt{\pi} = \frac{(2m-1)!!}{2^{m}}\sqrt{\pi} ,

which at $m = 0$ reads $\Gamma(1/2) = \sqrt{\pi}$ because the product is empty and $(-1)!! = 1$ . Substituting and separating one factor $\alpha^{-1/2}$ ,

\frac{\Gamma(m+\tfrac12)}{2\alpha^{\,m+1/2}} = \frac{(2m-1)!!}{2^{\,m+1}\alpha^{m}}\sqrt{\frac{\pi}{\alpha}} .

The integrand $x^{2m}e^{-\alpha x^2}$ is even and, by the computation just made, absolutely integrable over $(0,\infty)$ ; hence the full-line integral is twice the half-line one, which gives $(2m-1)!!\,(2\alpha)^{-m}\sqrt{\pi/\alpha} = \Gamma(m+\tfrac12)\alpha^{-(m+1/2)}$ .

Part 3 (odd moments). On $(0,\infty)$ the same substitution $u = \alpha x^2$ , now in the form $du = 2\alpha x\,dx$ , absorbs one power of $x$ :

\int_0^\infty x^{2m+1}e^{-\alpha x^2}dx = \int_0^\infty \Bigl(\frac{u}{\alpha}\Bigr)^{m}e^{-u}\,\frac{du}{2\alpha} = \frac{\Gamma(m+1)}{2\,\alpha^{\,m+1}} = \frac{m!}{2\,\alpha^{\,m+1}} ,

using $\Gamma(k) = (k-1)!$ for positive integers. Over the whole line the integrand is odd, and it is absolutely integrable because $\int_{\mathbb{R}}\lvert x\rvert^{2m+1}e^{-\alpha x^2}dx = 2\,m!/(2\alpha^{m+1}) < \infty$ by what has just been computed; so the contributions of $(0,\infty)$ and $(-\infty,0)$ are negatives of one another and may be added, giving $0$ .

Part 4 (completing the square). Set $\gamma = \beta/(2\alpha) \in \mathbb{C}$ . Expanding, $-\alpha(x-\gamma)^2 + \alpha\gamma^2 = -\alpha x^2 + 2\alpha\gamma x = -\alpha x^2 + \beta x$ , and $\alpha\gamma^2 = \beta^2/(4\alpha)$ . The integral converges absolutely for every $\beta \in \mathbb{C}$ , since $\lvert e^{-\alpha x^2 + \beta x}\rvert = e^{-\alpha x^2 + (\operatorname{Re}\beta)x}$ and the exponent tends to $-\infty$ quadratically. So

\int_{-\infty}^{\infty} e^{-\alpha x^2 + \beta x}dx = e^{\beta^2/(4\alpha)} \int_{-\infty}^{\infty} e^{-\alpha (x-\gamma)^2}\,dx ,

and it remains to show the last integral is $\sqrt{\pi/\alpha}$ for every complex $\gamma$ .

Write $\gamma = p + iq$ with $p, q$ real. Translating the real variable by $p$ leaves the integral unchanged, so we may take $\gamma = iq$ , and the claim becomes that $e^{-\alpha z^2}$ has the same integral along the line $\operatorname{Im}z = -q$ as along the real axis. Integrate $e^{-\alpha z^2}$ , an entire function, around the boundary of the rectangle with vertices $\pm R$ and $\pm R - iq$ : the closed contour integral is zero by Cauchy's theorem. On the vertical side $z = R + is$ , with $s$ between $0$ and $-q$ , we have $\operatorname{Re}(z^2) = R^2 - s^2 \ge R^2 - q^2$ , so $\lvert e^{-\alpha z^2}\rvert \le e^{-\alpha(R^2 - q^2)}$ ; the side has length $\lvert q\rvert$ , so its contribution is at most $\lvert q\rvert e^{-\alpha(R^2-q^2)}$ , which tends to $0$ as $R \to \infty$ , and the same bound holds on the side at $-R$ . Letting $R \to \infty$ , the two horizontal sides therefore contribute equal and opposite amounts, which is the claim. This is the contour of Hassani's Example 11.3.8, used here for a general complex $\gamma$ rather than only a purely imaginary one. Finally, $\beta = ib$ with $b$ real gives $\beta^2/(4\alpha) = -b^2/(4\alpha)$ . $\square$

Idea of the proof(Several dimensions)Original

A positive definite quadratic form is a sum of squares in the right coordinates, and the right coordinates are reached by two volume-preserving moves: a translation that kills the linear term, and a rotation that diagonalises the matrix. After them the integrand is a product of $n$ one-variable Gaussians, so the integral is a product of $n$ copies of part 1, and the only trace of the matrix that survives is $\lambda_1\cdots\lambda_n = \det A$ under the square root. No new analysis is involved; this is why the $n$ -dimensional formula is as cheap as the one-dimensional one.

Proof(Several dimensions)AdaptedAdapted from Riley, Hobson & Bence (2006), §30.13, Eqs. (30.144)–(30.145), p. 1206 and §30.15.2, Eq. (30.149), pp. 1209–1210

Step 1 (diagonalise). $A$ is real symmetric, so $A = Q\Lambda Q^{\mathsf T}$ with $Q$ orthogonal and $\Lambda = \operatorname{diag}(\lambda_1,\dots,\lambda_n)$ real (spectral theorem). Positive definiteness of $A$ says exactly that every $\lambda_i > 0$ (eigenvalue test). Hence $\det A = \lambda_1\cdots\lambda_n > 0$ and $A$ is invertible.

Step 2 (complete the square). Put $c = A^{-1}b$ , so that $Ac = b$ . Expanding and using $A^{\mathsf T} = A$ ,

-\tfrac12 (x-c)^{\mathsf T}A(x-c) = -\tfrac12 x^{\mathsf T}Ax + \tfrac12 x^{\mathsf T}Ac + \tfrac12 c^{\mathsf T}Ax - \tfrac12 c^{\mathsf T}Ac .

Now $x^{\mathsf T}Ac = x^{\mathsf T}b = b^{\mathsf T}x$ and $c^{\mathsf T}Ax = (Ac)^{\mathsf T}x = b^{\mathsf T}x$ , the second using the symmetry of $A$ , so the two middle terms add to $b^{\mathsf T}x$ ; and $c^{\mathsf T}Ac = c^{\mathsf T}b = b^{\mathsf T}A^{-1}b$ . Rearranging,

-\tfrac12 x^{\mathsf T}Ax + b^{\mathsf T}x = -\tfrac12(x-c)^{\mathsf T}A(x-c) + \tfrac12\, b^{\mathsf T}A^{-1}b .

Step 3 (two volume-preserving substitutions). Substitute $y = x - c$ : a translation of $\mathbb{R}^n$ onto itself with Jacobian determinant $1$ . Then substitute $y = Qz$ : a linear bijection with $\lvert \det Q\rvert = 1$ , since $Q^{\mathsf T}Q = \mathrm{I}$ forces $(\det Q)^2 = 1$ . Both are covered by the change-of-variables theorem, and under the second

y^{\mathsf T}Ay = z^{\mathsf T}Q^{\mathsf T}AQz = z^{\mathsf T}\Lambda z = \sum_{i=1}^{n}\lambda_i z_i^2 .

Step 4 (factorise). The constant from Step 2 comes out of the integral, leaving

\int_{\mathbb{R}^n}\exp\Bigl(-\tfrac12 \sum_i \lambda_i z_i^2\Bigr) d^n z = \int_{\mathbb{R}^n} \prod_{i=1}^{n} e^{-\lambda_i z_i^2/2}\, d^n z = \prod_{i=1}^{n}\int_{-\infty}^{\infty} e^{-\lambda_i z_i^2/2}\,dz_i ,

the exchange being legitimate for a non-negative measurable integrand (Tonelli) and each one-variable factor being finite because $\lambda_i > 0$ . This also proves the absolute convergence asserted in the statement, the integrand being positive. By part 1 of the one-dimensional statement with $\alpha = \lambda_i/2$ ,

\prod_{i=1}^{n}\sqrt{\frac{2\pi}{\lambda_i}} = \frac{(2\pi)^{n/2}}{\sqrt{\lambda_1\cdots\lambda_n}} = \frac{(2\pi)^{n/2}}{\sqrt{\det A}} .

Step 5 (assemble). Multiplying by $e^{\frac12 b^{\mathsf T}A^{-1}b}$ from Step 2 gives the displayed formula.

Step 6 (the three special cases). Case 1 is $b = 0$ . For case 2, note that $V = A^{-1}$ is symmetric with eigenvalues $1/\lambda_i > 0$ , hence symmetric positive definite, and $\det V = 1/\det A$ ; applying case 1 to $V^{-1}$ after the translation $x \mapsto x - \mu$ gives $(2\pi)^{n/2}(\det V^{-1})^{-1/2} = (2\pi)^{n/2}(\det V)^{1/2}$ . For case 3, replace $A$ by $2A$ in case 1 and use $\det(2A) = 2^{n}\det A$ : the answer is $(2\pi)^{n/2}2^{-n/2}(\det A)^{-1/2} = \pi^{n/2}(\det A)^{-1/2}$ .

Step 7 (divergence without positive definiteness). Steps 1 and 3 used only that $A$ is real symmetric, so for any such $A$ the integral equals $\prod_i \int_{\mathbb{R}} e^{-\lambda_i y_i^2/2 + c_i y_i}dy_i$ with $c = Q^{\mathsf T}b$ , as recorded in the statement. Each factor is positive, and the $i$ -th is finite only when $\lambda_i > 0$ : it equals $\sqrt{2\pi/\lambda_i}\,e^{c_i^2/(2\lambda_i)}$ then, by part 4 of the one-dimensional statement; it is $\int_{\mathbb{R}} e^{c_i y}dy = +\infty$ when $\lambda_i = 0$ (whatever $c_i$ is), and $+\infty$ when $\lambda_i < 0$ , the exponent then tending to $+\infty$ at one end at least. Since all factors are positive, a single infinite factor makes the product infinite, so the integral diverges whenever some eigenvalue fails to be positive — that is, whenever $A$ is not positive definite. $\square$

Method(Reducing a Gaussian integral to standard form)AdaptedAdapted from Riley, Hobson & Bence (2006), §30.13, p. 1206 and §30.15.2, pp. 1209–1210; Hassani (2013), Example 11.3.8, pp. 352–353

To evaluate $\displaystyle \int_{\mathbb{R}^n} P(x)\, e^{Q(x)}\, d^n x$ , where $Q$ is a quadratic polynomial and $P$ a polynomial:

Separate the exponent by degree. Write $Q(x) = -\tfrac12 x^{\mathsf T} A x + b^{\mathsf T} x + c$ and replace $A$ by $\tfrac12(A + A^{\mathsf T})$ , which changes nothing. Pull $e^{c}$ outside.
Test positive definiteness of $A$ : all eigenvalues positive, or equivalently all $n$ leading principal minors of $A$ positive. If the test fails, stop; the integral diverges.
Complete the square. Using $A^{\mathsf T} = A$ , $-\tfrac12 x^{\mathsf T} A x + b^{\mathsf T} x = -\tfrac12 (x - A^{-1}b)^{\mathsf T} A\, (x - A^{-1}b) + \tfrac12\, b^{\mathsf T} A^{-1} b ,$ and substitute $y = x - A^{-1}b$ : a translation, with Jacobian determinant $1$ , so $d^n x = d^n y$ and the domain is still $\mathbb{R}^n$ .
Diagonalise. Write $A = Q \Lambda Q^{\mathsf T}$ with $Q$ orthogonal and $\Lambda = \operatorname{diag}(\lambda_1,\dots,\lambda_n)$ (spectral theorem) and substitute $y = Qz$ ; since $\lvert\det Q\rvert = 1$ , $d^n y = d^n z$ , and the integral factorises into $\prod_i \int e^{-\lambda_i z_i^2/2}\,dz_i = \prod_i \sqrt{2\pi/\lambda_i} = (2\pi)^{n/2}(\det A)^{-1/2}$ .
Handle a polynomial prefactor $P$ by differentiation. Each factor $x_j$ in $P$ is produced by $\partial/\partial b_j$ acting on $e^{b^{\mathsf T}x}$ , so $\int_{\mathbb{R}^n} x_{j_1}\cdots x_{j_k}\, e^{-\frac12 x^{\mathsf T}Ax + b^{\mathsf T}x}\, d^n x = \frac{\partial^k}{\partial b_{j_1} \cdots \partial b_{j_k}} \left[ \frac{(2\pi)^{n/2}}{\sqrt{\det A}} e^{\frac12 b^{\mathsf T}A^{-1}b} \right],$ evaluated at the wanted $b$ (usually $b = 0$ ). In one variable the same job is done by $\partial/\partial\alpha$ : differentiating $\int e^{-\alpha x^2}dx = \sqrt{\pi/\alpha}$ $m$ times in $\alpha$ produces the even moment of order $2m$ with a factor $(-1)^m$ . Differentiation under the integral sign is justified because every derivative is dominated, on a neighbourhood of the parameter value, by an integrable Gaussian times a polynomial (differentiation under the integral sign).
A linear term that is not real. The statement above assumes $b \in \mathbb{R}^n$ , so a complex $b$ is outside it; the identity nevertheless persists, and the following justification is ours. With $x$ still real, $\lvert \exp(-\tfrac12 x^{\mathsf T}Ax + b^{\mathsf T}x)\rvert = \exp(-\tfrac12 x^{\mathsf T}Ax + (\operatorname{Re}b)^{\mathsf T}x)$ , which is the real case already settled, so the integral converges absolutely. Diagonalise as in step 4 — $y = Qz$ , with $c = Q^{\mathsf T}b \in \mathbb{C}^n$ — and split the absolutely convergent integral into a product by Fubini's theorem (Fubini and Tonelli); each factor is $\int_{\mathbb{R}} e^{-\lambda_i z_i^2/2 + c_i z_i}dz_i = \sqrt{2\pi/\lambda_i}\,e^{c_i^2/(2\lambda_i)}$ by part 4 of the one-dimensional statement, which is proved there for every complex coefficient. Multiplying the factors and using $\sum_i c_i^2/\lambda_i = c^{\mathsf T}\Lambda^{-1}c = b^{\mathsf T}Q\Lambda^{-1}Q^{\mathsf T}b = b^{\mathsf T}A^{-1}b$ returns the same closed form. The one-dimensional ingredient is where the complex shift is paid for: there the real variable cannot simply be translated, and a rectangular contour is closed instead, using that $e^{-\alpha z^2}$ is entire — Hassani's route for $\int e^{iax-bx^2}dx$ .

How it fails. Step 2 is not a formality: symmetry of $A$ alone does not give convergence. Step 4 needs $A$ real symmetric, so that $Q$ can be taken orthogonal; a matrix with complex entries is outside the statement above, and extending the formula to one requires a choice of branch for $(\det A)^{-1/2}$ . Purely imaginary exponents, $\int e^{i x^{\mathsf T}Ax}d^nx$ (Fresnel-type), are not covered at all: those integrals do not converge absolutely and need an oscillatory-integral argument.

RemarkOriginal

Where the hypotheses are used, and conventions that differ.

Positive definiteness, not just symmetry. Riley, Hobson and Bence's worked evaluation of $\int \exp(-\tfrac12 (x-\mu)^{\mathsf T} V^{-1}(x-\mu))\,d^n x$ (p. 1210) assumes only that $V$ is symmetric, and its orthogonal diagonalisation then gives $(2\pi)^{n/2}(\lambda_1\cdots\lambda_n)^{1/2}$ . That step needs $\lambda_i > 0$ for every $i$ : each factor is $\int e^{-y_i^2/(2\lambda_i)}dy_i$ , which is $+\infty$ if $\lambda_i < 0$ , while if some $\lambda_i = 0$ then $V$ is singular and $V^{-1}$ does not exist, so the integral is not defined; and $(\det V)^{1/2}$ would not be real if an odd number of eigenvalues were negative. In the setting of that section $V$ is a covariance matrix, where positive definiteness is automatic for a non-degenerate distribution; the statement above makes the hypothesis explicit because the formula needs it.
Why the $n$ -dimensional case needs no new analysis. Convergence of the factorised integral is absolute, so Tonelli's theorem lets the $n$ -fold integral be written as a product (Fubini and Tonelli); the substitutions in steps 3 and 4 are linear changes of variables with $\lvert\det\rvert = 1$ , covered by the change-of-variables theorem.
The one-dimensional value is not an antiderivative. $e^{-x^2}$ has no elementary antiderivative, which is why the standard route squares the integral, converts the resulting double integral to plane polar coordinates, and controls the unbounded region by squeezing the square $[-a,a]^2$ between the discs of radius $a$ and $a\sqrt{2}$ — the care Riley, Hobson and Bence take on pp. 203–204 precisely because the region is unbounded.
Normalisation conventions. Three conventions for the exponent are in common use: $e^{-\alpha x^2}$ (giving $\sqrt{\pi/\alpha}$ ), $e^{-x^2/(2\sigma^2)}$ (giving $\sigma\sqrt{2\pi}$ ), and in $n$ dimensions $e^{-\frac12 x^{\mathsf T}Ax}$ versus $e^{-x^{\mathsf T}Ax}$ , which differ by $2^{n/2}$ as recorded in case 3 above. Statistics and quantum field theory keep the factor $\tfrac12$ ; many methods texts omit it. Check which is in force before quoting a result.
Double factorial, and ranges of integration. $(2m-1)!!$ denotes the product of the odd numbers up to $2m-1$ , with $(-1)!! = 1$ ; Arfken's Exercise 13.1.10 writes the even moment both with $\Gamma(m+\tfrac12)$ and with the double factorial, and the two agree because $\Gamma(m+\tfrac12) = (2m-1)!!\sqrt{\pi}/2^{m}$ . Both of his formulas are over $(0,\infty)$ , which is why part 2 states the half-line value first and the full-line value as twice it. Some authors write $(2m)!/(2^m m!)$ for the same number.
The odd moments on the half-line are not zero. Part 3 gives $m!/(2\alpha^{m+1})$ over $(0,\infty)$ and $0$ over $\mathbb{R}$ ; the second follows from the first only because the integral converges absolutely, so the two halves may be added. The symmetry argument is not available for a Gaussian multiplied by a non-integrable factor.

Example(One dimension: a moment, a linear term, and both at once)Original

Take $\alpha = 2$ throughout and compute three integrals.

A pure moment. With $m = 2$ and $\alpha = 2$ , part 2 of the one-dimensional statement gives $\int_{-\infty}^{\infty} x^{4}e^{-2x^2}dx = \frac{(2\cdot 2 - 1)!!}{(2\cdot 2)^{2}}\sqrt{\frac{\pi}{2}} = \frac{3}{16}\sqrt{\frac{\pi}{2}} = \frac{3\sqrt{2\pi}}{32} = 0.2349964007\ldots$ Independent check by the other route. Differentiating $\int_{\mathbb{R}}e^{-\alpha x^2}dx = \sqrt{\pi}\,\alpha^{-1/2}$ twice in $\alpha$ (step 5 of the method) gives $\int_{\mathbb{R}}x^4 e^{-\alpha x^2}dx = \tfrac34\sqrt{\pi}\,\alpha^{-5/2}$ , and at $\alpha = 2$ that is $\tfrac34\sqrt{\pi}\,2^{-5/2} = 0.2349964007$ . Numerical quadrature agrees to ten digits.
A linear term. With $\beta = 3$ , part 4 gives $\int_{-\infty}^{\infty} e^{-2x^2 + 3x}dx = \sqrt{\frac{\pi}{2}}\;e^{9/8} = 1.2533141373 \times 3.0802168489 = 3.8604793227\ldots$ The exponent is largest at $x = \beta/(2\alpha) = 3/4$ , where $-2x^2+3x = 9/8$ , which is where the factor $e^{9/8}$ comes from.
Both. For the second moment with the linear term present, differentiate the answer to 2 twice in $\beta$ rather than inserting $x^2$ by hand. With $G(\beta) = \sqrt{\pi/\alpha}\,e^{\beta^2/(4\alpha)}$ , $G''(\beta) = \sqrt{\frac{\pi}{\alpha}}\;e^{\beta^2/(4\alpha)}\left[ \left(\frac{\beta}{2\alpha}\right)^{2} + \frac{1}{2\alpha} \right],$ so at $\alpha = 2$ , $\beta = 3$ the bracket is $(3/4)^2 + 1/4 = 13/16 = 0.8125$ and $\int_{-\infty}^{\infty} x^{2}e^{-2x^2+3x}dx = 3.8604793227 \times 0.8125 = 3.1366394497\ldots$ Quadrature returns $3.1366394497$ . Sanity check on the bracket: it is $\langle x^2\rangle$ for the Gaussian with mean $3/4$ and variance $1/(2\alpha) = 1/4$ , namely $(3/4)^2 + 1/4$ , as it should be.

Example(Three dimensions with a linear term)Original

Evaluate

I = \int_{\mathbb{R}^3} \exp\!\left( -x_1^2 - x_2^2 - x_3^2 + x_1x_2 + x_2x_3 + x_1 + x_3 \right) d^3x .

Read off $A$ and $b$ . The quadratic part is $-\tfrac12 x^{\mathsf T}Ax$ and the linear part $b^{\mathsf T}x$ with $A = \begin{pmatrix} 2 & -1 & 0 \\ -1 & 2 & -1 \\ 0 & -1 & 2 \end{pmatrix}, \qquad b = \begin{pmatrix} 1 \\ 0 \\ 1 \end{pmatrix},$ since $-\tfrac12(2x_1^2+2x_2^2+2x_3^2 - 2x_1x_2 - 2x_2x_3)$ reproduces the quadratic terms. $A$ is already symmetric, so step 1 of the method changes nothing.
Test positive definiteness. The leading principal minors are $2$ , $\det\begin{pmatrix}2&-1\\-1&2\end{pmatrix} = 3$ and $\det A = 2(4-1) + 1(-2-0) = 4$ , all positive, so $A$ is positive definite by Sylvester's criterion (quadratic forms). As a cross-check the eigenvalues are $2-\sqrt2$ , $2$ and $2+\sqrt2$ , with product $(4-2)\cdot 2 = 4 = \det A$ .
Invert $A$ . Cofactor expansion gives $A^{-1} = \frac{1}{4}\begin{pmatrix} 3 & 2 & 1 \\ 2 & 4 & 2 \\ 1 & 2 & 3 \end{pmatrix},$ which is verified by $A A^{-1} = \mathrm{I}$ : the first column of the product is $\tfrac14(2\cdot 3 - 1\cdot 2,\; -3 + 4 - 1,\; -2 + 2)^{\mathsf T} = (1,0,0)^{\mathsf T}$ , and similarly for the others.
Evaluate the exponent. $A^{-1}b = \tfrac14(3+1,\,2+2,\,1+3)^{\mathsf T} = (1,1,1)^{\mathsf T}$ , so $b^{\mathsf T}A^{-1}b = 1 + 0 + 1 = 2$ and $\tfrac12 b^{\mathsf T}A^{-1}b = 1$ .
Assemble. $I = \frac{(2\pi)^{3/2}}{\sqrt{4}}\,e^{1} = \frac{(2\pi)^{3/2}e}{2} = \sqrt{2}\,e\,\pi^{3/2} = 21.4059392604\ldots$ Iterated numerical quadrature over $[-12,12]^3$ gives $21.40594$ , which is the check.
A moment of the same form. Setting $b = 0$ and using step 5 of the method, $\int x_i x_j\,e^{-\frac12 x^{\mathsf T}Ax}d^3x$ divided by $\int e^{-\frac12 x^{\mathsf T}Ax}d^3x$ equals $(A^{-1})_{ij}$ . So for the centred density proportional to $e^{-\frac12 x^{\mathsf T}Ax}$ the second moments are $\langle x_1^2\rangle = 3/4$ , $\langle x_2^2\rangle = 1$ and $\langle x_1 x_3\rangle = 1/4$ . The pair $x_1, x_3$ is the one pair that the exponent does not couple, yet its covariance is nonzero: correlations spread beyond the couplings, because it is $A^{-1}$ and not $A$ that gives the second moments. The end coordinates also fluctuate less than the middle one.

ApplicationOriginal

Gaussian integrals are the only multidimensional integrals that can be done in closed form in general, which is why so much of physics is arranged to look Gaussian.

Classical fluctuations of a harmonic system. Three equal masses on a line, joined to each other and to two fixed walls by identical springs of stiffness $\kappa$ , have potential energy $U = \tfrac{\kappa}{2}\bigl[x_1^2 + (x_2-x_1)^2 + (x_3-x_2)^2 + x_3^2\bigr] = \tfrac12\, x^{\mathsf T}Kx , \qquad K = \kappa A ,$ with $A$ the matrix of the three-dimensional example above. In the canonical ensemble at temperature $T$ the displacements have density proportional to $e^{-U/k_BT}$ , which is the Gaussian with $A_{\text{eff}} = K/k_BT$ , so the configurational partition function is $(2\pi k_BT/\kappa)^{3/2}(\det A)^{-1/2} = \tfrac12 (2\pi k_BT/\kappa)^{3/2}$ and the second moments are $\langle x_i x_j \rangle = k_BT\,(K^{-1})_{ij} = \frac{k_BT}{\kappa}\,(A^{-1})_{ij} ,$ so $\langle x_1^2\rangle = 0.75\,k_BT/\kappa$ , $\langle x_2^2\rangle = k_BT/\kappa$ and $\langle x_1x_3\rangle = 0.25\,k_BT/\kappa$ : the middle mass, held by softer effective restoring forces, wanders furthest. Summing the same formula recovers equipartition, $\langle U \rangle = \tfrac12\sum_{i,j}K_{ij}\langle x_ix_j\rangle = \tfrac12 k_BT\operatorname{tr}\bigl(KK^{-1}\bigr) = \tfrac32 k_BT ,$ that is, $\tfrac12 k_BT$ per quadratic degree of freedom, for any positive definite $K$ whatever. The free energy $-k_BT\ln Z$ contains the term $\tfrac12 k_BT \ln \det K$ , which is how the determinant of a stiffness matrix becomes a measurable quantity.
Wave packets and propagators. Part 4 with $\beta = ib$ says that the Fourier transform of a Gaussian is a Gaussian of reciprocal width: a wave packet of position spread $\sigma$ has momentum spread $\hbar/(2\sigma)$ , and the product is the minimum allowed by the uncertainty relation. The same formula, with $\alpha$ taken to the imaginary axis, is the free-particle propagator $\langle x\vert e^{-iHt/\hbar}\vert x_0\rangle$ ; and the Gaussian path integral of a quadratic action is the infinite-dimensional version of the formula of the second statement, with $(\det A)^{-1/2}$ becoming a functional determinant.
Error estimates and least squares. A measurement with Gaussian errors of covariance $V$ has likelihood proportional to $\exp(-\tfrac12 (x-\mu)^{\mathsf T}V^{-1}(x-\mu))$ , normalised by $(2\pi)^{n/2}(\det V)^{1/2}$ from case 2 of the second statement. Fitting parameters by maximising it is exactly minimising the quadratic form in the exponent, which is least squares; the confidence regions are its level sets, ellipsoids with axes along the eigenvectors of $V$ and semi-axes proportional to $\sqrt{\lambda_i}$ , so the diagonalisation of step 4 of the method is also the statement that correlated errors can always be reported as independent ones in rotated variables.