# Finding a possible Lagrangian corresponding to a differential equation

A typical problem in the calculus of variations is to extremize a functional

$$J[y]=\int f(x,y,y') \, \mathrm{d} x.$$

This usually involves solving the Euler-Lagrange equation

$$\frac{\mathrm{d}}{\mathrm{d}x} \frac{\partial f}{\partial y'}-\frac{\partial f}{\partial y}=0.$$

The inverse of this problem is given a differential equation, determine $$f$$ such that a function $$y$$ is a solution to the original differential equation if and only if $$y$$ is a solution to the Euler-Lagrange equation. For some cases, like second-order linear differential equations, there are known solutions to this problem. In fact, Wolfram|Alpha can solve this exact problem. For example, running

WolframAlpha["y''+2y'-y=0", {{"PossibleLagrangian", 1}, "ComputableData"}]


gives one such solution:

Hold[ℒ[y', y, x] == 1/2 (E^(2 x) y^2 + E^(2 x) (y')^2)]


Is there a way to do this directly in Mathematica, without using Wolfram|Alpha or resorting to known formulas?

I've written a short function based on a known formula (The Calculus of Variations by Brunt, section 3.4), which gives the same answer as WA:

(* y'' + P y' + Q y - G == 0 *)

PossibleLagrangian[P_, Q_, G_, x_] :=
Block[{p = Exp[Integrate[P[t], {t, 0, x}]], q = Q[x] p, g = G[x] p},
1/2 (p (y')^2 - q y^2 + 2 g y)]

(* y'' + 2y' - y == 0 *)
PossibleLagrangian[2 &, -1 &, 0 &, x]
1/2 (E^(2 x) y^2 + E^(2 x) (y')^2)


 Names["**Lagrang*"]
(*{"NDSolveFEMHangingNodeLagrangeMultipliers",
"NDSolveFEMLagrangeMultipliers", "RegionLagrangeNearest"}*)