How does one correctly calculate the second-order variation of an action? I have started an attempt at the calculation (restricting the scalar fields for simplicity), but I'm unsure how to proceed.
Starting with the action for a free scalar field $$S[\phi]=\int\;d^{4}x\mathcal{L}=\frac{1}{2}\int\;d^{4}x\left(\partial_{\mu}\phi(x)\partial^{\mu}\phi(x)-m^{2}\phi^{2}(x)\right)$$
with Minkowski sign convention $(+,-,-,-)$. Naively, if I expand this to second-order, I get $$S[\phi+\delta\phi]=S[\phi]+\int\;d^{4}x\frac{\delta S[\phi(x)]}{\delta\phi(x)}\delta\phi(x)+\int\;d^{4}x d^{4}y\frac{\delta^{2} S[\phi(x)]}{\delta\phi(x)\delta\phi(y)}\delta\phi(x)\delta\phi(y)$$ Now, assuming that $\phi(x)$ satisfies the equations of motion (EOM), then the first-order term vanishes, however, I'm unsure how to calculate the second-order variation. So far, my attempt is $$\delta^{2}S=\int\;d^{4}x d^{4}y\frac{\delta^{2} S[\phi(x)]}{\delta\phi(x)\delta\phi(y)}\delta\phi(x)\delta\phi(y)=\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\;\;\;\\=\int\;d^{4}x d^{4}y\left(\frac{\partial^{2}\mathcal{L}}{\partial\phi(x)\partial\phi(y)}\delta\phi(x)\delta\phi(y)+2\frac{\partial^{2}\mathcal{L}}{\partial\phi(x)\partial(\partial_{\mu}\phi(y)}\delta\phi(x)\delta(\partial_{\mu}\phi(y))\\+\frac{\partial^{2}\mathcal{L}}{\partial(\partial_{\mu}\phi(x))\partial(\partial_{\mu}\phi(y))}\delta(\partial_{\mu}\phi(x))\delta(\partial_{\mu}\phi(y))\right)$$ However, I am unsure how to progress (integration by parts doesn't seem to work as nicely in this case), as naively it seems as though the only term that would survive is $\frac{\partial^{2}\mathcal{L}}{\partial\phi(x)\partial\phi(y)}$, but I've seen references stating that $\frac{\delta^{2} S[\phi(x)]}{\delta\phi(x)\delta\phi(y)}$ is of the form $\frac{\delta^{2} S[\phi(x)]}{\delta\phi(x)\delta\phi(y)}\sim\Box +m^{2}$.
Any help would be much appreciated.
Update
I have though about it a bit more and have come up with a general formula (for a Lagrangian with up to first-order derivatives in the fields) that I hope is correct: $$\frac{\delta^{2} S[\phi(x)]}{\delta\phi(x)\delta\phi(y)}=\frac{\delta^{2}}{\delta\phi(x)\delta\phi(y)}\int\;d^{4}z\,\mathcal{L}\left(\phi(z),\partial_{\mu}\phi(z)\right)\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\\=\frac{\delta}{\delta\phi(x)}\int\;d^{4}z\left[\frac{\partial\mathcal{L}}{\partial\phi(z)}\delta^{4}(z-y)+\frac{\partial\mathcal{L}}{\partial(\partial_{\mu}\phi(z))}\partial_{\mu}\left(\delta^{4}(z-y)\right)\right]\qquad\qquad\qquad\qquad\\=\int\;d^{4}z\left[\frac{\partial^{2}\mathcal{L}}{\partial\phi(z)^{2}}\delta^{4}(z-x)\delta^{4}(z-y)+2\frac{\partial^{2}\mathcal{L}}{\partial\phi(z)\partial(\partial_{\mu}\phi(z))}\partial_{\mu}\left(\delta^{4}(z-x)\right)\delta^{4}(z-y)\\+\frac{\partial^{2}\mathcal{L}}{\partial(\partial_{\mu}\phi(z))\partial(\partial_{\nu}\phi(z))}\partial_{\mu}\left(\delta^{4}(z-x)\right)\partial_{\nu}\left(\delta^{4}(z-x)\right)\right]\\ =\frac{\partial^{2}\mathcal{L}}{\partial\phi(x)^{2}}\delta^{4}(x-y)+2\frac{\partial^{2}\mathcal{L}}{\partial\phi(x)\partial(\partial_{\mu}\phi(x))}\partial_{\mu}\left(\delta^{4}(x-y)\right)\qquad\qquad\qquad\qquad\qquad\qquad\qquad-\frac{\partial^{2}\mathcal{L}}{\partial(\partial_{\mu}\phi(x))\partial(\partial_{\nu}\phi(x))}\partial_{\mu}\partial_{\nu}\left(\delta^{4}(x-y)\right)$$ which, upon further integrations by parts (neglecting boundary terms), gives $$\frac{\delta^{2} S[\phi(x)]}{\delta\phi(x)\delta\phi(y)}=\delta^{4}(x-y)\left[\frac{\partial^{2}\mathcal{L}}{\partial\phi(x)^{2}}-2\partial_{\mu}\left(\frac{\partial^{2}\mathcal{L}}{\partial\phi(x)\partial(\partial_{\mu}\phi(x))}\partial_{\mu}\right)\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad-\partial_{\mu}\partial_{\nu}\left(\frac{\partial^{2}\mathcal{L}}{\partial(\partial_{\mu}\phi(x))\partial(\partial_{\nu}\phi(x))}\right)\right]$$ This seems to give the expected answer in the case of a free scalar field. Indeed, $$\frac{\partial^{2}\mathcal{L}}{\partial\phi^{2}}=-m^{2}\, ,\qquad \frac{\partial^{2}\mathcal{L}}{\partial(\partial_{\mu}\phi)\partial(\partial_{\nu}\phi)}=\eta^{\mu\nu}\, , \qquad\frac{\partial^{2}\mathcal{L}}{\partial\phi\partial(\partial_{\mu}\phi)}=0$$ and hence, $$\frac{\delta^{2} S[\phi(x)]}{\delta\phi(x)\delta\phi(y)}=-\delta^{4}(x-y)\left[\Box+m^{2}\right]$$ Any feedback on whether this is correct or not would be much appreciated.
Best Answer
A proper treatment (and how you should usually go about these things if you forget) is to remember the definition of the functional derivative. It is linear, defined to obey a chain rule, a product rule, and has the fundamental feature
$$\frac{\delta\phi(y)}{\delta\phi(x)}=\delta(x-y)$$
Thus, in painstaking detail, we have
$$\frac{\delta S[\phi]}{\delta\phi(x)}=\frac{1}{2}\int\mathrm{d}^dy\left[\frac{\delta}{\delta\phi(x)}\left(\partial\phi(y)\cdot\partial\phi(y)\right)-m^2\frac{\delta}{\delta\phi(x)}\phi(y)^2\right]\\ =\int\mathrm{d}^dy\left[\partial_{\mu}\delta(x-y)\partial^{\mu}\phi(y)-m^2\delta(x-y)\phi(y)\right]\\ =-(\square+m^2)\phi(x)$$
Thus, we can simply differentiate again to obtain
$$\frac{\delta^2S[\phi]}{\delta\phi(x)\delta\phi(y)}=-\frac{\delta}{\delta\phi(y)}\left[(\square_x+m^2)\phi(x)\right]=-(\square_x+m^2)\delta(x-y)$$
Which is the desired result (note that $\square_x$ simply means that the derivative is only with respect to $x$ -- sometimes this matters)! Note that the delta function comes after the Klein-Gordon operator.
And that's it! No need to expand to second order or pull your hair out deciding whether you have to integrate by parts and when you can.
I hope this helps!
B-B-B-BONUS ROUND
This type of manipulation is actually extremely useful! For instance, in the path integral formulation, we have
$$\langle\mathcal{F}[\phi](x)\rangle=\int\mathcal{D}\phi\,\mathcal{F}[\phi](x)\,e^{iS[\phi]}$$
With this, we can use the above manipulations to find correlation functions! The key is to note that the path integral of a total functional derivative is zero. Thus, we have
$$\int\mathcal{D}\phi\,\frac{\delta^2}{\delta\phi(x)\delta\phi(y)}e^{iS[\phi]}=i\int\mathcal{D}\phi\left[\frac{\delta^2S}{\delta\phi(x)\delta\phi(y)}+i\frac{\delta S}{\delta\phi(x)}\frac{\delta S}{\delta\phi(y)}\right]e^{iS[\phi]}\\ =i\bigg\langle\frac{\delta^2S}{\delta\phi(x)\delta\phi(y)}+i\frac{\delta S}{\delta\phi(x)}\frac{\delta S}{\delta\phi(y)}\bigg\rangle=0$$
This holds for any action $S[\phi]$. In particular, in your free theory, this gives us
$$\left(\square_y+m^2\right)\left(\square_x+m^2\right)\langle\phi(x)\phi(y)\rangle=-i\left(\square_y+m^2\right)\delta(x-y)$$
Eliminating $\square_y+m^2$ from each side tells you that the two point function for a free theory is the Green's function of the Klein-Gordon operator. No need for generating functionals or all that messy second quantization.