Exercise 1.12: Lagrange's equations (code)

(require '[sicmutils.env :refer :all])

This exercise asks us to write Clojure implementations for each of the three systems described in Exercise 1.9.

Before we begin, here is a function that will display an up-tuple of:

  • inline_formula not implemented, the generalized force

  • inline_formula not implemented, the generalized momenta

  • inline_formula not implemented, the derivative of our momenta

  • The Lagrange equations for the system.

(defn lagrange-equation-steps [L q]
  (let [p1 (compose ((partial 1) L) (Gamma q))
        p2 (compose ((partial 2) L) (Gamma q))
        dp2 (D p2)]
    ((up p1 p2 dp2 (- dp2 p1))
     't)))

These are the steps required on the road to a derivation of Lagrange's equations.

Preliminary Notes

I found this exercise to be challenging because I was searching for a particular elegant form of the Lagrange equations for each system. Because the Lagrange equations are residuals, any linear combination of the equations also has to equal 0. In a few of the exercises below, I reached a solution that was technically correct, but cluttered.

If I were using Lagrangian mechanics to develop a game, or to simulate motion in some virtual world, I would just move on. But it seems that one of the tasks for the learner in modern physics is to take transferable lessons from the equations, and this means that it's important to try and unmask terms that might appear in different systems under superficially different forms.

Exercise 1.14 has an example of this problem that took me a long time to notice. The systems we analyze here happen to have yielded nice, familiar solutions. But note now that this is not a gimme.

Part A: Ideal Planar Pendulum

From the book:

An ideal planar pendulum consists of a bob of mass inline_formula not implemented connected to a pivot by a massless rod of length inline_formula not implemented subject to uniform gravitational acceleration inline_formula not implemented. A Lagrangian is inline_formula not implemented. The formal parameters of inline_formula not implemented are inline_formula not implemented, inline_formula not implemented, and inline_formula not implemented; inline_formula not implemented measures the angle of the pendulum rod to a plumb line and inline_formula not implemented is the angular velocity of the rod.

Here is the Lagrangian described by the exercise:

(defn L-pendulum [m g l]
  (fn [[_ theta thetadot]]
    (+ (* (/ 1 2) m (square l) (square thetadot))
       (* m g l (cos theta)))))

And the steps that lead us to Lagrange's equations:

(lagrange-equation-steps
 (L-pendulum 'm 'g 'l)
 (literal-function 'theta))

The final entry is the Lagrange equation, equal to inline_formula not implemented. Divide out the shared factors of inline_formula not implemented and inline_formula not implemented:

(let [L (L-pendulum 'm 'g 'l)
      theta (literal-function 'theta)
      eqs ((Lagrange-equations L) theta)]
  ((/ eqs (* 'm 'l))
   't))

This is the familiar equation of motion for a planar pendum.

Part B: 2D Potential

The next problem is in rectangular coordinates. This means that we'll end up with two Lagrange equations that have to be satisfied.

From the book:

A particle of mass inline_formula not implemented moves in a two-dimensional potential inline_formula not implemented, where inline_formula not implemented and inline_formula not implemented are rectangular coordinates of the particle. A Lagrangian is inline_formula not implemented.

I have no intuition for what this potential is, by the way. One term, inline_formula not implemented, looks like half the square of the distance of the particle away from 0, or inline_formula not implemented. What are the other terms? I've been so well trained that I simply start calculating.

Define the Lagrangian to be the difference of the kinetic energy and some potential inline_formula not implemented that has access to the coordinates:

(defn L-2d-potential [m]
  (fn [V]
    (fn [local]
      (- (* (/ 1 2) m (square (velocity local)))
         (V (coordinate local))))))

Thanks to the tuple algebra of scmutils, This form of the Lagrangian is general enough that it would work for any number of dimensions in rectangular space, given some potential inline_formula not implemented. square takes a dot product, so we end up with a kinetic energy term for every spatial dimension.

Note this for later, as this idea will become useful when the book reaches the discussion of coordinate transformations.

Next define the potential from the problem description:

(defn V [[x y]]
  (- (+ (/ (+ (square x)
              (square y))
           2)
        (* (square x) y))
     (/ (cube y) 3)))

Our helpful function generates the Lagrange equations, along with each intermediate step:

(lagrange-equation-steps
 ((L-2d-potential 'm) V)
 (up (literal-function 'x)
     (literal-function 'y)))

The final down-tuple gives us the Lagrange equations that inline_formula not implemented and inline_formula not implemented (respectively) must satisfy.

Part C: Particle on a Sphere

This problem is slightly more clear. From the book:

A Lagrangian for a particle of mass inline_formula not implemented constrained to move on a sphere of radius inline_formula not implemented is inline_formula not implemented. The angle inline_formula not implemented is the colatitude of the particle and inline_formula not implemented is the longitude; the rate of change of the colatitude is inline_formula not implemented and the rate of change of the longitude is inline_formula not implemented.

So the particle has some generalized kinetic energy with terms for:

  • its speed north and south, and

  • its speed east and west, scaled to be strongest at 0 longitude along the inline_formula not implemented axis and fall off to nothing at the inline_formula not implemented axis.

Here is the Lagrangian:

(defn L-sphere [m R]
  (fn [[_ [theta] [alpha beta]]]
    (* (/ 1 2) m (square R)
       (+ (square alpha)
          (square (* beta (sin theta)))))))

Here is the full derivation:

(lagrange-equation-steps
 (L-sphere 'm 'R)
 (up (literal-function 'theta)
     (literal-function 'phi)))

The final Lagrange residuals have a few terms that we can divide out. Scheme doesn't know that these are meant to be residuals, so it won't cancel out factors that we can see by eye are missing.

Isolate the Lagrange equations from the derivation and manually simplify each equation by dividing out, respectively, inline_formula not implemented and inline_formula not implemented:

(let [L     (L-sphere 'm 'R)
      theta (literal-function 'theta)
      q     (up theta (literal-function 'phi))
      le    ((Lagrange-equations L) q)
      eq1   (ref le 0)
      eq2   (ref le 1)]
  ((up (/ eq1 (* 'm (square 'R)))
       (/ eq2 (* (sin theta) 'm (square 'R))))
   't))

These are the Lagrange equations for inline_formula not implemented and inline_formula not implemented, respectively.

Runtimes (1)