# Equation of motion through the Lagrangian with Lagrange multipliers

We have such a Lagrangian:

$$L=\frac{1}{2}m(\dot{x}^2+\dot{y}^2)-\lambda(x+xy+y-1)$$

Here $$\lambda(x+xy+y-1)$$ is the constraint on the phase variables.

I need to derive the equation of motion given the constraints and solve them numerically with the help of NDSolve.

We do this in accordance with the classic formula:

$$\frac{d}{dt}(\frac{dL}{d\dot{q}})-\frac{dL}{dq}=0$$

Where $$q=[x,y]$$ are generalized coordinates. I'm not sure about the Lagrange multiplier as a generalized coordinate.

Clear["Derivative"]

ClearAll["Global*"]

T = 1/2 m (x'[t]^2 + y'[t]^2);(*Kinetic Energy*)

f = \[Lambda] (x[t] + x[t] y[t] +y[t] - 1);(*Constraint*)

L = T - f;(*Lagrangian*)

D[D[L, x'[t]], t] - D[L, x[t]];

D[D[L, y'[t]], t] - D[L, y[t]];

D[D[L, \[Lambda]'[t]], t] - D[L, \[Lambda][t]];
`

Question: how are the Lagrange multipliers included in this system when compiling the ODE system and numerically solving it?

Maybe this help?

https://farside.ph.utexas.edu/teaching/336k/lectures/node90.html

https://www.sciencedirect.com/science/article/abs/pii/0045782588900850

• You're trying to move along the curve $x+xy+y=1$, then? If so then your function doesn't specify how the velocity should evolve. Alternately, I could understand this as there being a potential energy penalizing your distance from there? If it is the latter then you don't really have a Lagrange multiplier, you just have an intensity coefficient for the potential energy.
– Ian
Oct 16, 2021 at 3:29
• @lan Yes, and I need to understand the principle of obtaining a system of differential equations (i.e., obtain a system of equations of motion) for different types of constraints and solve it numerically.
– dtn
Oct 16, 2021 at 3:31
• Well I don't think that constraint will tell you how the evolution will proceed, if you remain on the curve all the time. You need some other source of force otherwise the system will just sit at one point on the curve forever.
– Ian
Oct 16, 2021 at 3:32
• I don't think you can just put the system in motion and subject it to the constraint and then just have it "decide" how to stay on the curve without a force field, but I could be wrong, this isn't my area of expertise.
– Ian
Oct 16, 2021 at 3:38
• Well then you're back to ordinary Lagrangian dynamics essentially, there's no "Lagrange multiplier", there's just a potential energy pushing you towards the desired curve.
– Ian
Oct 16, 2021 at 3:40

If the underlying Lagrangian is just kinetic energy (i.e. there is no potential energy in the system) and you need to stay on the curve $$F=0$$, then the Newton equations for the constrained Lagrangian dynamics become just $$m\ddot{x}-\lambda(t) F_x=0,m\ddot{y}-\lambda(t) F_y=0$$, so the force points parallel to $$\nabla F$$ but with a coefficient that depends on time, and the requirement that the constraint is satisfied all the time determines the evolution of $$\lambda$$. In order for these equations to be consistent, the initial velocity will need to be tangent to the curve, but in this case you can solve the two differential equations together with the algebraic equation to describe the evolution.

Note that this constrained Lagrangian dynamics physically corresponds to "making up the right force that keeps you on the specified curve" (so it is a bit less magical than it looks at first glance). For a classic example you can think of moving on the circle $$x^2+y^2=r^2$$ starting at $$(r,0)$$ with initial velocity $$(0,v)$$. From this you can recover the familiar $$\mathbf{a}=-\frac{v^2}{r} \mathbf{e}_r$$.

A perhaps more physical view of a situation like this in the absence of any meaningful potential energy is to just take some strictly convex function $$G$$ minimized at zero (e.g. $$G(x)=x^2$$). and run the evolution throughout with $$m(\ddot{x},\ddot{y})^T=-\nabla G(F(x,y))$$.

• Thank you for your answer. I see that the question is not as simple as it seems and there are several nuances. 1. What will change if potential energy is added to the system? 2. And how to make the Lagrangian evolve correctly, i.e. to choose a change in the $\lambda$ coefficient?
– dtn
Oct 16, 2021 at 4:15
• @dtn 1. You can follow your first link to see how another term is introduced if there is another force like gravity on the system. This term takes the usual $\frac{\partial L}{\partial q}$ form. 2. From the math side this is called a differential algebraic equation (DAE)...how they're handled numerically depends a lot on the actual problem you're faced with. I'm not an expert on this either.
– Ian
Oct 16, 2021 at 4:46
• Please see mathematica.stackexchange.com/questions/256885/… maybe it will be useful in the future
– dtn
Oct 16, 2021 at 11:18

This is a DAE system of index 2 or 3. You obtain an ODE system by adding derivatives of the equation, mainly the constraint, to the system. $$x+xy+y=1,\\ \dot x(1+y)+\dot y(1+x)=0,\\ \ddot x(1+y)+\ddot y(1+x)+2\dot x\dot y=0$$ From the last equation eliminate the second derivatives to get an equation for the Lagrange multiplier. With another derivative (thus index 3) you get an ODE for the multiplier. $$m\ddot x=-λ(1+y)\\ m\ddot y=-λ(1+y)\\ \implies -λ[(1+x)^2+(1+y)^2]+2m\dot x\dot y=0\\ λ=\frac{2m\dot x\dot y}{(1+x)^2+(1+y)^2}$$

• Thank you for your answer! Please see mathematica.stackexchange.com/questions/256885/… maybe it will be useful in the future
– dtn
Oct 16, 2021 at 11:18
• I added the result one should get for $λ$. This is sufficient to get the equations of motion in $x,y,\dot x,\dot y$ only. For the derivative of $λ$ apply the quotient rule to get a more complicated rational expression. Oct 17, 2021 at 6:44
• This approach didn't work. I tried both that and the one suggested here. The numerical solution of the obtained systems of equations shows that the dynamics of the Lagrangian multiplier does not fulfill the indicated restrictions. Please see mathematica.stackexchange.com/questions/256885/… It seems that I programmed everything correctly. I don't see the error yet.
– dtn
Oct 18, 2021 at 8:42
• Your initial conditions there do no satisfy the constraint equation. $\frac12+\frac14·\frac12+\frac14=\frac78$, and the numerical integration should approximately stay at this constant for $x+xy+y$. Oct 18, 2021 at 9:12
• Yes, I have already fixed this error. And the principle is now more or less completely clear to me, it seems.
– dtn
Oct 18, 2021 at 9:13