# Legendre Transform of a function of a 3-vector

I am trying to implement the Legendre transform of a function in Mathematica with the purpose of calculating the Hamiltonian of a system starting from the Lagrangian.

I have found this page which explained how to perform the transform on a function of a single variable. But I need to evaluate the transform of a function of a vector in 3D.

Here it is my clumsy attempt at it:

Off[Solve::ifun];
Off[InverseFunction::ifun];
Clear[legendreTransform];
legendreTransform[f_, v_, p_] :=
Module[{d1 = D[f, v[[1]]], d2 = D[f, v[[2]]], d3 = D[f, v[[3]]],transform},
transform /; (transform = (-f + v[[1]]*d1 + v[[2]]*d2 + v[[3]]*d3)
/.Solve[{p[[1]] == d1 && p[[2]] == d2 && p[[3]] == d3}, v[[1]], v[[2]], v[[3]]}] // Last;
FreeQ[transform, Solve] && FreeQ[transform, InverseFunction])] /; ! (v === p);
vel = {vel1, vel2, vel3};
mom = {mom1, mom2, mom3};
legendreTransform[Exp[Total[vel]], vel, mom]


This piece of code is not working. I was not able to troubleshoot it also because I don't fully understand the code in the page I linked in the first place. I looked up every word in the documentation, but I still cannot understand code flow.

If you could give my any hint or suggestion about where to look I would be grateful.

EDIT (SOLUTION):

The VariationalMethods path seemed me more appealing at first and in fact I used it to do my calculations. But it produced slightly wrong results in my case. It took me some time to realize it but there was a problem with a relative sign.

Eventually I wrote this function that takes the best of the two approaches (the direct and the alternate answer) and it is working as intended in my case:

legendreTransform[l_, q_, p_] :=
Module[{mom = Flatten[{p}], pos = Flatten[{q}], h},
First[h /.Quiet[Solve[h == D[pos,t].VariationalD[l[pos], D[pos, t], t] - l[pos] &&
mom == VariationalD[l[pos], D[pos, t], t], Append[D[pos, t], h]],
{Solve::incnst, Solve::ifun}]]]


More insight on what I am doing: I am doing calculations in primordial cosmology and, as you may know, in General Relativity the Hamiltonian is null in value (but not in form). This is sometimes called "frozen formalism" because with a null Hamiltonian is trickier to define the evolution of the system in time. Maybe that is the reason the standard methods in VariationalMethods are failing me. More precisely, there is the wrong sign between the kinetical term and the potential term (the spatial curvature term) in the Hamiltonian produced by Legendre transform of the Lagrangian (I am working in the ADM formalism). This at first was not a big issue since all I did was introducing the correct sign by hand.

• Your code has syntax errors. It is not bracketed correctly. Either you did not copy the code correctly from its source or you made clumsy modifications or you did copy and paste in correctly into your question. No matter which, we can't help you until you correct the syntax errors. Also if have taken code from from someone else, you should give a reference to the source in your question. – m_goldberg Jul 15 '17 at 1:56
• Apparently the linked code is from @Bob Hanlon, and he may be able to provide the most direct answer. I just posted an alternative route. – Jens Jul 15 '17 at 6:23
• @m_goldberg I agree, it is very probable that the code in the question is not correctly copied and pasted. Anyway it is not working even on my notebook in Mathematica. This is why I am asking here. I already included the link to the reference, maybe you didn't notice. – LastStarDust Jul 16 '17 at 9:24
• @Jens Thank you very much. Your direct answer is very readable. I think I got it. You are right, my test function is not a scalar. That is a blunder on my side. Anyway I had already tested my code with a scalar function and it wasn't working as well. – LastStarDust Jul 16 '17 at 9:27

The function you were trying to transform in the question is not a scalar and can therefore not be used as an example.

Here is a version of legendreTransform that works for one or multiple variables, with only the bare minimum of code so as to make it more readable (hopefully):

legendreTransform[f_, v_, p_] :=
Module[{vel = Flatten[{v}], mom = Flatten[{p}], h},
First[h /.
Append[vel, h]], {Solve::incnst, Solve::ifun}]]]


It leaves the hard work to Solve, with the relevant warnings turned off only inside the function, not globally. The arguments for the variables can be symbols (in one dimension) or lists of symbols (in two or more dimensions). Both cases are reduced to the same code by wrapping even a 1D variable in a List (using Flatten) so that it can be used in Grad. Of course the names for the components of velocity and momentum must be different, but are otherwise arbitrary. For the sign of the result, I used the convention of Wikipedia, which agrees with the one in the question.

The first example reproduces the one-dimensional test case:

legendreTransform[Exp[v], v, p]


$p (\log (p)-1)$

The second example is a toy Lagrangian in three dimensions:

legendreTransform[
1/2 (v1^2 + v2^2 + v3^2) - (v1 - v2)/2, {v1, v2, v3}, {p1, p2, p3}]


$$\frac{1}{4} \left(2 \text{p1}^2+2 \text{p1}+2 \text{p2}^2-2 \text{p2}+2\text{p3}^2+1\right)$$

Given that you're interested in going from a Lagrangian to a Hamiltonian, you could use a generalization of this answer:

Needs["VariationalMethods"]

Clear[pos, mom, r, p, t, h, V]
pos[t_] := Table[r[i][t], {i, 3}]
mom[t_] := Table[p[i][t], {i, 3}]

L = 1/2 m  pos'[t].pos'[t] - V[pos[t]];

Expand[h /. First@Solve[
mom[t] == VariationalD[L, pos'[t], t] && h == FirstIntegral[t] /.
FirstIntegrals[L, pos[t], t], Append[pos'[t], h]]]


This uses the VariationalMethods package and therefore I had to define positions pos and momenta mom as functions of the time t.

Then I introduced a standard Lagrangian L from which you get to the Hamiltonian by Legendre transformation. But instead of doing the transformation by hand, I let FirstIntegrals find the Hamiltonian in terms of the positions and velocities. Then I just have to substitute the momenta for the velocities, using the fact that momenta are the variational derivative of L with respect to velocities pos'[t].

The substitution is done in Solve. I ask it to solve for the unknown Hamiltonian h, and also for the velocities. Although I then discard the velocities, this strategy makes sure that Solve finds a unique answer in terms of the variables I want to keep, i.e., the positions and momenta.

• Thank you so much. I adapted your direct answer to the particular calculation I am doing in cosmology and it worked as a charm. In the next days I am going to study your alternative answer as well. Thank you again. Actually I wasn't expecting such a quick and exhaustive answer. You swept me off my feet, so to say. – LastStarDust Jul 16 '17 at 10:53
• Glad to hear that. The VariationalMethods approach is especially useful if you want to automatically get equations of motion from a Lagrangian, see for example the double pendulum where a manual approach is very ugly. The analogous thing for Hamilton's equations of motion is addressed in another post. – Jens Jul 16 '17 at 17:17
• You are right. I ended up using the VariationalMethods` approach too. – LastStarDust Jul 20 '17 at 6:47