[proofplan]
If the two coordinate directions are the same, the assertion is only an equality of a quantity with itself. Otherwise, we study the rectangular second difference of $f$ in the $x_i$- and $x_j$-directions near $a$. Applying the one-variable [mean value theorem](/theorems/186) first in one coordinate and then in the other writes the same rectangular quotient as a nearby value of $\partial_{x_j}\partial_{x_i}f$, while reversing the order writes it as a nearby value of $\partial_{x_i}\partial_{x_j}f$. Letting the rectangle shrink to $a$ and using the assumed continuity of both mixed partials gives the equality.
[/proofplan]
[step:Dispose of the identical-coordinate case]
If $i=j$, then both sides of the desired identity are the same iterated [partial derivative](/page/Partial%20Derivative):
\begin{align*}
\partial_{x_i}\partial_{x_i}f(a)=\partial_{x_i}\partial_{x_i}f(a).
\end{align*}
Hence the conclusion holds in this case. For the rest of the proof, assume $i\neq j$.
[/step]
[step:Choose a coordinate rectangle contained in the differentiability neighbourhood]
Let $e_1,\ldots,e_m$ denote the standard basis of $\mathbb{R}^m$. For each $\rho>0$, let $B(a,\rho):=\{x\in\mathbb{R}^m:|x-a|<\rho\}$ denote the open Euclidean ball of radius $\rho$ centred at $a$. Since $V$ is an open neighbourhood of $a$, there exists $\rho>0$ such that
\begin{align*}
B(a,\rho)\subset V.
\end{align*}
Choose $\delta>0$ such that $\delta\sqrt{2}<\rho$. Then, whenever $s,t\in\mathbb{R}$ satisfy $0<|s|<\delta$ and $0<|t|<\delta$, every point of the coordinate rectangle
\begin{align*}
a+r e_i+q e_j
\end{align*}
with $r$ between $0$ and $s$ and $q$ between $0$ and $t$ lies in $V$, because
\begin{align*}
|r e_i+q e_j|\leq \sqrt{r^2+q^2}<\delta\sqrt{2}<\rho.
\end{align*}
Thus all one-variable functions obtained by restricting $f$, $\partial_{x_i}f$, and $\partial_{x_j}f$ to these coordinate line segments are defined on the relevant intervals.
[/step]
[step:Express the rectangular quotient through $\partial_{x_j}\partial_{x_i}f$]
Define the rectangular second-difference map $D:(( -\delta,\delta)\setminus\{0\})\times(( -\delta,\delta)\setminus\{0\})\to\mathbb{R}$ by
\begin{align*}
D(s,t)=f(a+s e_i+t e_j)-f(a+s e_i)-f(a+t e_j)+f(a).
\end{align*}
Define the one-variable function $F_{s,t}:[0,s]\to\mathbb{R}$ if $s>0$, and $F_{s,t}:[s,0]\to\mathbb{R}$ if $s<0$, by
\begin{align*}
F_{s,t}(r)=f(a+r e_i+t e_j)-f(a+r e_i).
\end{align*}
The function $F_{s,t}$ is continuous on its closed interval and differentiable on its interior, because the partial derivative $\partial_{x_i}f$ exists at every point of the rectangle. Its derivative is
\begin{align*}
F_{s,t}'(r)=\partial_{x_i}f(a+r e_i+t e_j)-\partial_{x_i}f(a+r e_i).
\end{align*}
By the one-variable mean value theorem applied to $F_{s,t}$, there exists a real number $\theta_{s,t}$ strictly between $0$ and $1$ such that
\begin{align*}
D(s,t)=s\bigl(\partial_{x_i}f(a+\theta_{s,t}s e_i+t e_j)-\partial_{x_i}f(a+\theta_{s,t}s e_i)\bigr).
\end{align*}
Now define the one-variable function $G_{s,t}:[0,t]\to\mathbb{R}$ if $t>0$, and $G_{s,t}:[t,0]\to\mathbb{R}$ if $t<0$, by
\begin{align*}
G_{s,t}(q)=\partial_{x_i}f(a+\theta_{s,t}s e_i+q e_j).
\end{align*}
The function $G_{s,t}$ is continuous on its closed interval and differentiable on its interior, because the iterated partial derivative $\partial_{x_j}(\partial_{x_i}f)$ exists at every point of the rectangle. Its derivative is
\begin{align*}
G_{s,t}'(q)=\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+q e_j).
\end{align*}
Applying the one-variable mean value theorem to $G_{s,t}$ gives a real number $\eta_{s,t}$ strictly between $0$ and $1$ such that
\begin{align*}
\partial_{x_i}f(a+\theta_{s,t}s e_i+t e_j)-\partial_{x_i}f(a+\theta_{s,t}s e_i)=t\,\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+\eta_{s,t}t e_j).
\end{align*}
Therefore
\begin{align*}
\frac{D(s,t)}{st}=\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+\eta_{s,t}t e_j).
\end{align*}
[guided]
The purpose of this step is to turn a two-coordinate difference quotient into an ordinary one-variable problem twice. Fix nonzero $s,t$ with $|s|<\delta$ and $|t|<\delta$, and define
\begin{align*}
D(s,t)=f(a+s e_i+t e_j)-f(a+s e_i)-f(a+t e_j)+f(a).
\end{align*}
This is the change in the $x_j$-increment of $f$ as the $x_i$-coordinate moves from $0$ to $s$.
To make the first application of the mean value theorem legitimate, define $F_{s,t}$ on the interval with endpoints $0$ and $s$ by
\begin{align*}
F_{s,t}(r)=f(a+r e_i+t e_j)-f(a+r e_i).
\end{align*}
All points appearing here lie in $V$ by the rectangle choice. Since $\partial_{x_i}f$ exists on $V$, the restrictions of $f$ to the two $x_i$-coordinate line segments are differentiable in the one-variable sense. Hence $F_{s,t}$ is continuous on the closed interval and differentiable on the open interval between its endpoints. Its derivative is
\begin{align*}
F_{s,t}'(r)=\partial_{x_i}f(a+r e_i+t e_j)-\partial_{x_i}f(a+r e_i).
\end{align*}
The one-variable mean value theorem therefore gives a number $\theta_{s,t}$ strictly between $0$ and $1$ such that
\begin{align*}
D(s,t)=F_{s,t}(s)-F_{s,t}(0)=s\bigl(\partial_{x_i}f(a+\theta_{s,t}s e_i+t e_j)-\partial_{x_i}f(a+\theta_{s,t}s e_i)\bigr).
\end{align*}
We now apply the same one-variable idea to the remaining difference, this time in the $x_j$-direction and to the function $\partial_{x_i}f$. Define $G_{s,t}$ on the interval with endpoints $0$ and $t$ by
\begin{align*}
G_{s,t}(q)=\partial_{x_i}f(a+\theta_{s,t}s e_i+q e_j).
\end{align*}
The existence of $\partial_{x_j}(\partial_{x_i}f)$ on $V$ says exactly that this one-variable restriction is differentiable in the $q$ variable at interior points. Thus
\begin{align*}
G_{s,t}'(q)=\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+q e_j).
\end{align*}
A second application of the one-variable mean value theorem gives a number $\eta_{s,t}$ strictly between $0$ and $1$ such that
\begin{align*}
\partial_{x_i}f(a+\theta_{s,t}s e_i+t e_j)-\partial_{x_i}f(a+\theta_{s,t}s e_i)=t\,\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+\eta_{s,t}t e_j).
\end{align*}
Substituting this into the previous formula for $D(s,t)$ and dividing by the nonzero product $st$ yields
\begin{align*}
\frac{D(s,t)}{st}=\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+\eta_{s,t}t e_j).
\end{align*}
[/guided]
[/step]
[step:Express the same quotient through $\partial_{x_i}\partial_{x_j}f$]
Define the one-variable function $H_{s,t}:[0,t]\to\mathbb{R}$ if $t>0$, and $H_{s,t}:[t,0]\to\mathbb{R}$ if $t<0$, by
\begin{align*}
H_{s,t}(q)=f(a+s e_i+q e_j)-f(a+q e_j).
\end{align*}
The function $H_{s,t}$ is continuous on its closed interval and differentiable on its interior, with
\begin{align*}
H_{s,t}'(q)=\partial_{x_j}f(a+s e_i+q e_j)-\partial_{x_j}f(a+q e_j).
\end{align*}
By the one-variable mean value theorem, there exists a real number $\alpha_{s,t}$ strictly between $0$ and $1$ such that
\begin{align*}
D(s,t)=t\bigl(\partial_{x_j}f(a+s e_i+\alpha_{s,t}t e_j)-\partial_{x_j}f(a+\alpha_{s,t}t e_j)\bigr).
\end{align*}
Define the one-variable function $K_{s,t}:[0,s]\to\mathbb{R}$ if $s>0$, and $K_{s,t}:[s,0]\to\mathbb{R}$ if $s<0$, by
\begin{align*}
K_{s,t}(r)=\partial_{x_j}f(a+r e_i+\alpha_{s,t}t e_j).
\end{align*}
The function $K_{s,t}$ is continuous on its closed interval and differentiable on its interior, with
\begin{align*}
K_{s,t}'(r)=\partial_{x_i}\partial_{x_j}f(a+r e_i+\alpha_{s,t}t e_j).
\end{align*}
Applying the one-variable mean value theorem to $K_{s,t}$ gives a real number $\beta_{s,t}$ strictly between $0$ and $1$ such that
\begin{align*}
\partial_{x_j}f(a+s e_i+\alpha_{s,t}t e_j)-\partial_{x_j}f(a+\alpha_{s,t}t e_j)=s\,\partial_{x_i}\partial_{x_j}f(a+\beta_{s,t}s e_i+\alpha_{s,t}t e_j).
\end{align*}
Therefore
\begin{align*}
\frac{D(s,t)}{st}=\partial_{x_i}\partial_{x_j}f(a+\beta_{s,t}s e_i+\alpha_{s,t}t e_j).
\end{align*}
[/step]
[step:Shrink the rectangle and use continuity at $a$]
For every $0<|s|<\delta$ and $0<|t|<\delta$, the two formulas for the same quotient give
\begin{align*}
\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+\eta_{s,t}t e_j)=\partial_{x_i}\partial_{x_j}f(a+\beta_{s,t}s e_i+\alpha_{s,t}t e_j).
\end{align*}
Since $\theta_{s,t},\eta_{s,t},\alpha_{s,t},\beta_{s,t}\in(0,1)$, we have
\begin{align*}
|\theta_{s,t}s e_i+\eta_{s,t}t e_j|\leq \sqrt{s^2+t^2}
\end{align*}
and
\begin{align*}
|\beta_{s,t}s e_i+\alpha_{s,t}t e_j|\leq \sqrt{s^2+t^2}.
\end{align*}
Thus both evaluation points tend to $a$ as $(s,t)\to(0,0)$ with $s\neq 0$ and $t\neq 0$.
The function $\partial_{x_j}\partial_{x_i}f:V\to\mathbb{R}$ is continuous at $a$, so
\begin{align*}
\partial_{x_j}\partial_{x_i}f(a+\theta_{s,t}s e_i+\eta_{s,t}t e_j)\to \partial_{x_j}\partial_{x_i}f(a)
\end{align*}
as $(s,t)\to(0,0)$ through nonzero $s$ and $t$. Similarly, the continuity at $a$ of $\partial_{x_i}\partial_{x_j}f:V\to\mathbb{R}$ gives
\begin{align*}
\partial_{x_i}\partial_{x_j}f(a+\beta_{s,t}s e_i+\alpha_{s,t}t e_j)\to \partial_{x_i}\partial_{x_j}f(a).
\end{align*}
Passing to the limit in the equality of the two expressions therefore yields
\begin{align*}
\partial_{x_j}\partial_{x_i}f(a)=\partial_{x_i}\partial_{x_j}f(a).
\end{align*}
This is the desired equality of mixed partial derivatives.
[/step]