Unit 13 - Notes: Examples of Lagrangian Mechanics

We’ve seen the Lagrangian formulation of mechanics, and we have seen how to use it to derive the equations of motion for a system. We will cover a few common examples of Lagrangian mechanics and point out some of the subtleties that arise in each case.

Example: The Atwood Machine¶

The Atwood machine consists of two masses, $M$ and $m$ , connected by a massless string that passes over a massless pulley. The system is subject to gravity alone. The figure below shows the system along with a choice of coordinates $y_1$ and $y_2$ for the two masses.

These coordinates are measured from the center of the pulley and positive $y_1$ and $y_2$ are taken to be upward. Let’s try to use the Lagrangian formalism to find the equations of motion for this system.

V = + Mgy_1 + mgy_2

(1)

T = \dfrac{1}{2}M\dot{y}_1^2 + \dfrac{1}{2}m\dot{y}_2^2

(2)

Equation of Constraint¶

But notice that $y_1$ and $y_2$ are no independent coordinates. If we unravel the string that is over the pulley, we find that (assume length of string is $l$ ):

y_1 + \pi R + y_2 = l

(3)

where $R$ is the radius of the pulley. That is shown in the figure below.

The equation above is called an equation of constraint. It relates the coordinates $y_1$ and $y_2$ to each other. We can use this equation to eliminate one of the coordinates. Let’s eliminate $y_2$ :

l = y_1 + \pi R + y_2 \rightarrow y_1 = (l - \pi R) - y_2

(4)

This constraint has implications for velocities,

\dfrac{dy_1}{dt} = \dfrac{d}{dt}\left((l-\pi R) - y_2\right) = -\dfrac{dy_2}{dt}

(5)

\dot{y}_1 = -\dot{y}_2.

(6)

This is likely what we could have expected, that the two masses move in opposite directions at the same speed.

Constructing the Lagrangian¶

Let’s use this constraint to reduce the number of coordinates in the Lagrangian.

\mathcal{L}(y_1, y_2, \dot{y}_1, \dot{y}_2) \rightarrow \mathcal{L}(y_1, \dot{y}_1)

(7)

We can do this by substituting $y_2$ in terms of $y_1$ into the energy equations:

T(y_1, \dot{y}_1) = \dfrac{1}{2}M\dot{y}_1^2 + \dfrac{1}{2}m\dot{y}_2^2 = \dfrac{1}{2}(M+m)\dot{y}_1^2 = T(\dot{y}_1)

(8)

V(y_1, y_2) = +Mgy_1 + mgy_2

(9)

V(y_1, y_2) = Mgy_1 + mg((l-\pi R) - y_1)

(10)

V(y_1, y_2) = (M-m)g y_1 + mg(l-\pi R)

(11)

V(y_1, y_2) = (M-m)g y_1 + U_0

(12)

where $U_0 = mg(l-\pi R)$ is a constant and will not affect the equations of motion.

\mathcal{L}(y_1, \dot{y}_1) = T(\dot{y}_1) - V(y_1) = \dfrac{1}{2}(M+m)\dot{y}_1^2 - (M-m)g y_1

(13)

\dfrac{\partial \mathcal{L}}{\partial y_1} = -(M-m)g

(14)

\dfrac{\partial \mathcal{L}}{\partial \dot{y}_1} = (M+m)\dot{y}_1

(15)

These derivatives give the following equation of motion:

-(M-m) g - \dfrac{d}{dt}\left((M+m)\dot{y}_1\right) = 0

(16)

Generalized Force¶

Notice that the first term in the above equation is the force on the mass $M$ in the Newtonian picture: the weight of $M$ minus the weight of $m$ . That makes sense because the Lagrangian formalism is supposed to reproduce Newton’s laws, and the spatial derivative of the Lagrangian produces a generalized force.

\dfrac{\partial \mathcal{L}}{\partial q_i} = -\dfrac{\partial V}{\partial q_i} = F_i

(17)

The kinetic term has no spatial dependence, so it does not contribute to the generalized force.

Generalized Momentum¶

The second term in the above equation is the time derivative of the momentum of the system using $y_1$ as the coordinate:

(M+m)\dot{y}_1 = M\dot{y}_1 - m\dot{y}_2 = p_{y_1}

(18)

Again, that is a sensible result because the Lagrangian formalism is supposed to reproduce Newton’s laws, and the generalized force is related to the time derivative of the generalized momentum.

\dfrac{d}{dt}\left(\dfrac{\partial \mathcal{L}}{\partial \dot{q}_i}\right) = \dfrac{d}{dt}\left(\dfrac{\partial T}{\partial \dot{q}_i}\right) = \dfrac{dp_{q_i}}{dt}

(19)

Equation of Motion¶

The equation of motion can be written as

(M+m)\ddot{y}_1 = -(M-m)g

(20)

\ddot{y}_1 = -\dfrac{(M-m)}{(M+m)}g.

(21)

With $M > m$ , this acceleration is downward as the larger mass $M$ accelerates down. With $\dot{y}_1 = -\dot{y}_2$ , we know

\ddot{y}_2 = -\ddot{y}_1

(22)

, so that:

\ddot{y}_2 = \dfrac{(M-m)}{(M+m)}g

(23)

Again, with $M > m$ , this acceleration is upward as the smaller mass $m$ accelerates up.

Example: Atwood Machine with Rotating Pulley¶

In the previous example, we didn’t take into account the energy needed to rotate the pulley. Let’s do that now. Because the rope cannot slip, any small rotation $R d\phi$ of the pulley gives a change $dy_1$ in the position of mass $M$ . This is the no slip constraint.

If the pulley has a mass $M_p$ and radius $R$ , then we must introduce it’s kinetic energy:

T_{pulley} = \dfrac{1}{2} I \omega^2

(24)

where $I$ is the moment of inertia of the pulley and $\omega$ is the angular velocity of the pulley, $\omega = \dot{\phi}$ . This angular velocity is related to the linear velocities of the masses.

I = \dfrac{1}{2}M_p R^2

(25)

T_{pulley} = \dfrac{1}{2} \left(\dfrac{1}{2}M_p R^2\right) \dot{\phi}^2

(26)

T_{pulley} = \dfrac{1}{4}M_p R^2 \dot{\phi}^2

(27)

We now map this additional kinetic energy into the problem.

T(\dot{y}_1, \dot{\phi}) = \dfrac{1}{2}(M+m)\dot{y}_1^2 + \dfrac{1}{4}M_p R^2 \dot{\phi}^2

(28)

But the constraint is such that,

dy_1 = Rd\phi \rightarrow y_1 = R\phi +\underbrace{R\phi_0}_{\textrm{const.}}

(29)

\dot{y}_1 = R\dot{\phi}

(30)

We work these back into $T$ , $V$ , and $\mathcal{L}$ .

T(\dot{\phi}) = \dfrac{1}{2}(M+m)\dot{y}_1^2 + \dfrac{1}{4}M_p R^2 \dot{\phi}^2

(31)

T(\dot{\phi}) = \dfrac{1}{2}(M+m+\frac{1}{2}M_p)R^2 \dot{\phi}^2

(32)

V(\phi) = (M-m)gy_1 + U_0

(33)

V(\phi) = (M-m)g(R\phi + R\phi_0) + U_0

(34)

V(\phi) = (M-m)gR\phi + \tilde{U}_0

(35)

where $\tilde{U}_0 = (M-m)gR\phi_0 + U_0$ is another constant that will not affect the equations of motion.

\mathcal{L}(\phi, \dot{\phi}) = T(\dot{\phi}) - V(\phi)

(36)

\mathcal{L}(\phi, \dot{\phi}) = \dfrac{1}{2}(M+m+\frac{1}{2}M_p)R^2 \dot{\phi}^2 - (M-m)gR\phi

(37)

Torque and Angular Momentum¶

Our generalized coordinate is $\phi$ and our generalized velocity is $\dot{\phi}$ . We apply the Euler-Lagrange equation. We obtain the generalized force:

\dfrac{\partial \mathcal{L}}{\partial \phi} = -(M-m)gR = F_{\phi}

(38)

Notice that in this case, the generalized force is not the same as the force on mass $M$ in the Newtonian picture. It’s a torque around the pulley.

We can find the generalized momentum in a similar way:

\dfrac{\partial \mathcal{L}}{\partial \dot{\phi}} = (M+m+\frac{1}{2}M_p)R^2\dot{\phi} = p_{\phi}

(39)

This equals the generalized momentum associated with the pulley angle. The hanging masses move vertically, not tangentially around the pulley, so their contributions should be computed from their actual vertical position vectors and velocities rather than from $R\hat r$ and $R\dot\phi\hat\phi$ .

\vec{L}_{total} = \vec{L}_{disk} + \vec{L}_{M} + \vec{L}_{m}

(40)

\vec{L}_{disk} = I_{disk}\omega = \dfrac{1}{2}M_p R^2 \dot{\phi}\,\text{(out of the page)}

(41)

\vec{L}_{M} = M\vec{r}_{M}\times \vec{v}_{M}\quad\text{(computed using the actual vertical position and velocity)}

(42)

\vec{L}_{m} = m\vec{r}_{m}\times \vec{v}_{m}\quad\text{(computed using the actual vertical position and velocity)}

(43)

Add them up:

\vec{L}_{total} = \left(\dfrac{1}{2}M_p + M + m\right)R^2 \dot{\phi}\,\text{(out of the page)}

(44)

Or the magnitude:

L_{total} = \left(\dfrac{1}{2}M_p + M + m\right)R^2 \dot{\phi}

(45)

p_{\phi} = L_{total}

(46)

Equation of Motion¶

We return to the differential equation of motion:

$-(M-m)gR - \dfrac{d}{dt}\left((M+m+\frac{1}{2}M_p)R^2\dot{\phi}\right)=0$

-(M-m)gR - (M+m+\frac{1}{2}M_p)R^2\ddot{\phi} = 0

(47)

which produces the following equation of motion:

\ddot{\phi} = -\dfrac{g}{R}\dfrac{(M-m)}{(M+m+\frac{1}{2}M_p)}

(48)

which is a constant acceleration.

Comparing this to the massless-pulley result $\ddot{y}_1 = -\frac{(M-m)}{(M+m)}g$ from before, we see the pulley’s moment of inertia acts like extra inertia in the denominator, slowing the whole system down. The plot below shows the distance the heavier mass $M$ has descended after time $t$ for a few choices of pulley mass $M_p$ .

import numpy as np
import matplotlib.pyplot as plt
plt.style.use('seaborn-v0_8-colorblind')

g = 9.8
M, m = 2.0, 1.0
t = np.linspace(0, 2, 200)

fig, ax = plt.subplots(figsize=(7, 5))
for Mp, color in zip([0, 1, 3, 6], ['C0', 'C1', 'C2', 'C3']):
    a = (M - m) / (M + m + 0.5 * Mp) * g
    y1 = 0.5 * a * t**2
    label = 'Massless pulley' if Mp == 0 else fr'$M_p={Mp}$ kg'
    ax.plot(t, y1, color=color, label=label)

ax.set_xlabel('Time (s)')
ax.set_ylabel(r'Distance descended by $M$ (m)')
ax.set_title('Effect of Pulley Inertia on the Atwood Machine')
ax.legend()
ax.grid(True)
plt.tight_layout()
plt.show()

Example: Bead in a Parabolic Bowl¶

A bead of mass $m$ is constrained to move along a parabolic bowl. There is a gravitational force acting on the bead. The bowl is symmetric about the $z$ -axis and the bead is constrained to move along the surface without friction or rolling. The bowl is described by the equation:

z = \dfrac{1}{2}c(x^2 + y^2)

(49)

where $c$ is a constant that describes the curvature of the bowl. The figure below shows the system.

In this case the system is better solved in cylindrical coordinates. The coordinates are $(r, \phi, z)$ , where $r$ is the distance from the $z$ -axis, $\phi$ is the angle around the $z$ -axis, and $z$ is the height above the $xy$ -plane as shown above.

With $\langle \rho, \phi, z\rangle$ as the coordinates, equation for the constraint is:

z = \dfrac{1}{2}c(x^2 + y^2) = \dfrac{1}{2}c\rho^2.

(50)

Note that $c$ has units.

[z] = m \quad [\rho^2] = m^2 \quad [c] = \dfrac{1}{m}

(51)

The speed in cylindrical coordinates can be derived from the expressions in Cartesian coordinates, but we quote the result here:

v^2 = \dot{\rho}^2 + \rho^2\dot{\phi}^2 + \dot{z}^2.

(52)

Constructing the Lagrangian¶

We write the kinetic and potential energy of the bead in terms of the coordinates $(r, \phi, z)$ .

T = \dfrac{1}{2}m\left(\dot{\rho}^2 + \rho^2\dot{\phi}^2 + \dot{z}^2\right)

(53)

V = mgz

(54)

In principle, the Lagrangian can depend on all three coordinates and all three velocities.

\mathcal{L}(\rho, \dot{\rho}, \phi, \dot{\phi}, z, \dot{z}, t) = T - V

(55)

There is no explicit time dependence, so we can ignore $t$ . When the Lagrangian has no explicit time dependence, we should expect the energy to be conserved. This form of Lagrangian analysis does not account for dissipation.

\mathcal{L}(\rho, \dot{\rho}, \phi, \dot{\phi}, z, \dot{z})

(56)

Moreover, with the symmetry of the problem, we can expect no $\phi$ dependence. That is the gravitational potential energy only depends on $z$ and not on $\phi$ . This is a consequence of the symmetry of the problem and indicates that $\phi$ is a cyclic coordinate. We thus expect angular momentum to be conserved about the $z$ -axis.

\mathcal{L}(\rho, \dot{\rho}, \dot{\phi}, z, \dot{z})

(57)

Lastly, the constraint equation gives us $z$ in terms of $\rho$ :

z = \dfrac{1}{2}c\rho^2

(58)

And thus we can find the time derivative of $z$ :

\dot{z} = c\rho\dot{\rho}.

(59)

So the Lagrangian can be simplified three variables:

\mathcal{L}(\rho, \dot{\rho}, \dot{\phi}) = T - V

(60)

\mathcal{L}(\rho, \dot{\rho}, \dot{\phi}) = \dfrac{1}{2}m\left((1+c^2\rho^2)\dot{\rho}^2 + \rho^2\dot{\phi}^2\right) - mg\dfrac{1}{2}c\rho^2

(61)

Equations of Motion¶

We can now apply the Euler-Lagrange equations to find the equations of motion. We will do this for each coordinate. Let’s start with $\phi$ because there is only one term in the Lagrangian that depends on $\dot{\phi}$ .

\underbrace{\dfrac{\partial \mathcal{L}}{\partial \phi}}_0 - \dfrac{d}{dt}\left(\dfrac{\partial \mathcal{L}}{\partial \dot{\phi}}\right) = 0

(62)

\dfrac{d}{dt}\left(\dfrac{\partial \mathcal{L}}{\partial \dot{\phi}}\right) = \dfrac{d}{dt}\underbrace{\left(m\rho^2\dot{\phi}\right)}_{L_z} = 0

(63)

This equation of motion indicates that the angular momentum about the $z$ -axis is conserved, as we expected from the symmetry of the problem.

For the coordinate $\rho$ , we have:

\dfrac{\partial \mathcal{L}}{\partial \rho} - \dfrac{d}{dt}\left(\dfrac{\partial \mathcal{L}}{\partial \dot{\rho}}\right) = 0

(64)

\dfrac{\partial \mathcal{L}}{\partial \rho} = m\rho\dot{\phi}^2 + mc^2\rho\dot{\rho}^2 - mgc\rho

(65)

\dfrac{\partial \mathcal{L}}{\partial \dot{\rho}} = m(1+c^2\rho^2)\dot{\rho}

(66)

\dfrac{d}{dt}\left(\dfrac{\partial \mathcal{L}}{\partial \dot{\rho}}\right) = m(1+c^2\rho^2)\ddot{\rho} + 2mc^2\rho\dot{\rho}^2

(67)

We can now write the equation of motion:

m\rho\dot{\phi}^2 + mc^2\rho\dot{\rho}^2 - mgc\rho - m(1+c^2\rho^2)\ddot{\rho} - 2mc^2\rho\dot{\rho}^2 = 0

(68)

We can clean this up a little bit:

\ddot{\rho}(1+c^2\rho^2) + c^2\rho\dot{\rho}^2 = \rho\dot{\phi}^2 - gc\rho

(69)

\ddot{\rho}(1+c^2\rho^2) + c^2\rho\dot{\rho}^2 - \rho\dot{\phi}^2 + gc\rho = 0

(70)

We try to write both equations in terms of their accelerations. We have:

\ddot{\rho} = -\dfrac{c^2\rho\dot{\rho}^2 - \rho\dot{\phi}^2 + gc\rho}{1+c^2\rho^2}

(71)

\ddot{\phi} = - \dfrac{2\rho\dot{\rho}\dot{\phi}}{\rho^2}

(72)

For which we can develop a solution anywhere away from the origin ( $\rho \neq 0$ ).

Preparing for Numerical Solution¶

We need to write these equations in a form that is ready for numerical solution. We can do this by writing the equations in terms of the first derivatives. Let $\omega = \dot{\phi}$ and $v = \dot{\rho}$ . We get 4 1st order equations:

\dot{\rho} = v

(73)

\dot{\phi} = \omega

(74)

\dot{v} = -\dfrac{c^2\rho v^2 - \rho\omega^2 + gc\rho}{1+c^2\rho^2}

(75)

\dot{\omega} = - \dfrac{2 v\omega}{\rho}

(76)

Let’s integrate these equations numerically, launching the bead at $\rho_0=1$ with $\dot\rho_0=0$ and a nonzero angular velocity $\dot\phi_0=1.2$ . We plot $\rho(t)$ to see the radial oscillation, and reconstruct the 3D path $(x,y,z)$ to see the bead trace out its orbit on the bowl.

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
plt.style.use('seaborn-v0_8-colorblind')

g = 9.8
c = 1.0

def bead(t, state):
    rho, phi, v, omega = state
    drho = v
    dphi = omega
    dv = -(c**2 * rho * v**2 - rho * omega**2 + g * c * rho) / (1 + c**2 * rho**2)
    domega = -2 * v * omega / rho
    return [drho, dphi, dv, domega]

rho0, omega0 = 1.0, 1.2
state0 = [rho0, 0.0, 0.0, omega0]
t_span = (0, 15)
t_eval = np.linspace(*t_span, 3000)
sol = solve_ivp(bead, t_span, state0, t_eval=t_eval)

rho, phi, v, omega = sol.y
z = 0.5 * c * rho**2
x = rho * np.cos(phi)
y = rho * np.sin(phi)

fig = plt.figure(figsize=(11, 5))
ax1 = fig.add_subplot(1, 2, 1)
ax1.plot(sol.t, rho, 'C0')
ax1.set_xlabel('Time (s)')
ax1.set_ylabel(r'$\rho$ (m)')
ax1.set_title(r'Radial Motion $\rho(t)$')
ax1.grid(True)

ax2 = fig.add_subplot(1, 2, 2, projection='3d')
ax2.plot(x, y, z, 'C1')
ax2.set_xlabel('x')
ax2.set_ylabel('y')
ax2.set_zlabel('z')
ax2.set_title('Bead Trajectory on the Paraboloid')

plt.tight_layout()
plt.show()

The bead’s radius oscillates steadily between a minimum and maximum value rather than settling down, because there is nothing to dissipate energy in this frictionless system, while conservation of $L_z = m\rho^2\dot\phi$ forces $\dot\phi$ to speed up whenever $\rho$ shrinks, tracing out the looping path on the bowl’s surface shown at right.