Non-trivial functional derivatives

Question

Looking for ways to calculate functional derivatives in Mathematica, I found this question. The answer by Jens points to a documented application of DiracDelta exactly for this purpose. The documentation suggests

\[ScriptCapitalD][functional_, f_[y_]] := 
 Assuming[Element[y, Reals], 
  Limit[((functional /. 
        f :> Function[x, f[x] + ε DiracDelta[x - y]]) - 
      functional)/ε, ε -> 0]]

However, this implementation as well the other answers to that question are very limited in scope. As pointed out in the comments, they fail for an example as simple as $\delta G[f(x)]/\delta f(y)$ (or as code \[ScriptCapitalD][G[f[x]], f[y]]), which should output

$$G^\prime[f(x)] \, \delta(x-y).$$

In my case I'd like to compute functional derivatives of the form

$$\frac{\delta^n \, \Gamma[\rho(\varphi)]}{\delta \varphi_{i_1}(p_1) \dots \delta \varphi_{i_n}(p_n)},$$

where $\rho(\varphi) = \frac{1}{2} \sum_{i = 1}^n \varphi_i^2$ is a regular function of the $\varphi_i$, $i \in \{1,\dots,n\}$.

Is there a way to equip Mathematica with the ability to calculate such derivatives with justifiable effort?

Answer 1

Update 2

My first update added support for derivatives of the functionals. This second update adds support for multiple functional derivatives, as long as they are all the same functional. Adding support for multiple, different functionals is harder, but could be done.

FunctionalD[expr_, v:(f_[_]|{f_[_],_Integer}) .., OptionsPattern[]] := Internal`InheritedBlock[{f},
    f /: D[f[x_], f[y_], NonConstants->{f}] := DiracDelta[x-y];
    f /: D[f, f[y_], NonConstants->{f}] := DiracDelta[#-y]&;

    D[expr, v, NonConstants->{f}]
]

Here is FunctionalD on the OP's simple example:

FunctionalD[G[f[x]], f[y]]
% //TeXForm

DiracDelta[x - y] Derivative[1][G][f[x]]

$\delta (x-y) G'(f(x))$

After the update, the following example of a Lagrangian and action for a simple 1D particle in a potential now works:

L = 1/2 m x'[t]^2 - V[x[t]];
S = Integrate[L, {t, t0, t1}];
S //TeXForm

$\int_{\operatorname{t0}}^{\operatorname{t1}} \left(\frac{1}{2} m x'(t)^2-V(x(t))\right) \, dt$

The functional derivative of this action is:

Assuming[t0<τ<t1, FunctionalD[S, x[τ]]] //TeXForm

$-m x''(\tau )-V'(x(\tau ))$

and setting the above to zero yields the usual equations of motion for a 1D particle.

And finally, the example from the comments with multiple functional derivatives:

FunctionalD[Exp[-(1/2) f[x]^2], {f[y],2}, f[z]] //Factor //TeXForm

$-e^{-\frac{1}{2} f(x)^2} f(x) \left(f(x)^2-3\right) \delta (x-y)^2 \delta (x-z)$