Suppose we want to find the electrostatic potential $\phi$(r) inside of some volume $V$ with known boundary conditions. The physical field configuration should minimize the electrostatic potential energy function, which is defined as
$U[\phi] = \frac{\epsilon_0}{2} \int_v ( \nabla \phi)^2$d r
Show that requiring $\phi$ to be a stationary point of $U$ is equivalent to Laplace's equation, $\nabla^2\phi = 0$.
A hint I was given: take a small deformation of $\phi$, that is $\phi + \delta\phi$, such that $\delta\phi$ vanishes on the boundary of V. We must show that $\delta U = U[\phi + \delta\phi] - U[\phi] = 0$, to linear order in $\delta\phi$.
I have tried limits from first principles and first order taylor expansions. I always seem to just move in a loop however, constantly reaching no conclusion. If anyone is able to show me how this can be done, that would be fantastic.
$\newcommand{\angles}[1]{\left\langle\, #1 \,\right\rangle} \newcommand{\braces}[1]{\left\lbrace\, #1 \,\right\rbrace} \newcommand{\bracks}[1]{\left\lbrack\, #1 \,\right\rbrack} \newcommand{\ceil}[1]{\,\left\lceil\, #1 \,\right\rceil\,} \newcommand{\dd}{{\rm d}} \newcommand{\ds}[1]{\displaystyle{#1}} \newcommand{\expo}[1]{\,{\rm e}^{#1}\,} \newcommand{\fermi}{\,{\rm f}} \newcommand{\floor}[1]{\,\left\lfloor #1 \right\rfloor\,} \newcommand{\half}{{1 \over 2}} \newcommand{\ic}{{\rm i}} \newcommand{\iff}{\Longleftrightarrow} \newcommand{\imp}{\Longrightarrow} \newcommand{\pars}[1]{\left(\, #1 \,\right)} \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}} \newcommand{\pp}{{\cal P}} \newcommand{\root}[2][]{\,\sqrt[#1]{\vphantom{\large A}\,#2\,}\,} \newcommand{\sech}{\,{\rm sech}} \newcommand{\sgn}{\,{\rm sgn}} \newcommand{\totald}[3][]{\frac{{\rm d}^{#1} #2}{{\rm d} #3^{#1}}} \newcommand{\ul}[1]{\underline{#1}} \newcommand{\verts}[1]{\left\vert\, #1 \,\right\vert}$ \begin{align} \verts{\nabla\pars{\phi + \delta\phi}}^{2} - \verts{\nabla\phi}^{2} &=\pars{\nabla\phi + \nabla\delta\phi}\cdot\pars{\nabla\phi + \nabla\delta\phi} -\pars{\nabla\phi}^{2} \\[5mm]&=2\nabla\phi\cdot\nabla\delta\phi + \pars{\nabla\delta\phi}^{2} \end{align}
Then, $$ 0=\delta{\rm U}=\int\bracks{2\nabla\cdot\pars{\delta\phi\,\nabla\phi} - 2\delta\phi\,\nabla^{2}\phi}\,\dd{\bf r} =2\int_{S}\delta\phi\,\nabla\phi\cdot\dd\vec{\rm S} -2\int\delta\phi\,\nabla^{2}\phi\ \dd{\bf r} $$