Welcome to my tour of Structure and Interpretation of Classical Mechanics. I’m working on this book to develop my sense of the best way to do research in public; this book is heavy on math, programming and visualization, and should stress the normal tools.
I’m attempting to take notes in on org-mode file, and generate all my code from there.
I don’t think I have the heart, or the time, to really do high-class notes of every single section; but I am going to do each of the exercises, and explore some of the code-based concepts in each section.
The preface is already intriguing. A tour through the new notation, plus some discussion of why a programming language is the best route in to this stuff. Both of these are extremely powerful ideas, and why I was pulled to this book in the first place.
The functional notation is:
\begin{equation} D (∂_2 L ˆ Γ[q]) - (∂_1 L ˆ Γ[q]) = 0 \end{equation}
Compare that to the traditional notation:
\begin{equation} \frac{d}{dt} \frac{∂ L}{∂ \dot q^i} -\frac{∂ L}{∂ q^i}= 0 \end{equation}
They have a nice riff on how this is totally ambiguous.
This exercise is designed to get you thinking about the idea of configuration space. One goal of classical mechanics is to predict the future of a system, based on a description of the current state of the system. So how to describe the system?
One way would be to track, for every instant of time, a 3-dimensional position for every particle in a system, along with velocities in each direction, for a total of 6 numbers - or dimensions - per particle.
That’s fine for particles moving in straight lines through space, not affecting each other. But remember, this is all an accounting trick that we’re using to represent reality, not reality itself. If there is some more convenient way to track the state of the system, we’re free to use that as well, provided we can recover the original positions and velocities after evolving the system.
Maybe two particles are attached to each other, so their positions in space are the same, or offset by a tiny amount. Then it’s enough to track:
- the position of one of the particles (3 numbers)
- the angle and distance of the offset (2 numbers)
What seemed like a system that required 6 numbers actually required 5.
Theo Jansen’s incredible Strandbeesten are built out of copies of Jansen’s Linkage. Each of these legs has 11 pipes that flex and bend, but only one degree of freedom.
The first exercise gives us some practice thinking about the redundancy in different physical systems.
For each of the mechanical systems described below, give the number of degrees of freedom of the configuration space. (SICM, ex1)
Exercise 1.2 asks about the “generalized coordinates” of each system. What are
the actual numbers that we want to track for each system, if not the
For each of the systems in exercise 1.1, specify a system of generalized coordinates that can be used to describe the behavior of the system.
- Three juggling pins.
The system has **18 degrees of freedom**. Each pin requires 3 coordinates to specify its center of mass, and 3 angles for each pin. If you assume that each pin is symmetric about its central axis, then it doesn’t matter how far around the pin has rotated and you can make do with **15 degrees of freedom**. 3 positions and 2 angles for each.
- A spherical pendulum consisting of a point mass (the pendulum bob) hanging
from a rigid massless rod attached to a fixed support point. The pendulum bob
may move in any direction subject to the constraint imposed by the rigid rod.
The point mass is subject to the uniform force of gravity.
This system has only **2 degrees of freedom**. One for the latitude of the pendulum, and one for the longitude.
- A spherical double pendulum, consisting of one point mass hanging from a
rigid massless rod attached to a second point mass hanging from a second
massless rod attached to a fixed support point. The point masses are subject
to the uniform force of gravity.
**4 degrees of freedom**; two angles from previous, plus two more angles for the second pendulum off of the first.
- A point mass sliding without friction on a rigid curved wire.
**1 degree of freedom**; the distance along the wire.
- A top consisting of a rigid axisymmetric body with one point on the symmetry
axis of the body attached to a fixed support, subject to a uniform
gravitational force.
This system seems to have **2 degrees of freedom**, for the angles off of vertical. It’s like a spherical pendulum, but upside down. What I find strange about this answer is that the top does have a rotation speed, which is a measure of how far the top has rotated in time. How can we track this velocity if we don’t track the top’s spin angle? This may be wrong.
- The same as 5, but not axisymmetric.
A non-symmetric top has **3 degrees of freedom**. 2 from before, and an additional angle to measure how far around the top has rotated.
This problem has us exploring some consequences for optics of the principle of least time. Exercise 1.3 states:
Fermat observed that the laws of reflection and refraction could be accounted for by the following facts: Light travels in a straight line in any particular medium with a velocity that depends upon the medium. The path taken by a ray from a source to a destination through any sequence of media is a path of least total time, compared to neighboring paths. Show that these facts imply the laws of reflection and refraction.
The law of reflection is described in the footnote:
For reflection the angle of incidence is equal to the angle of reflection.
Here’s the setup. The horizontal line is a mirror. The law states that
We have to show that if we consider all possible paths from a given starting point to a given endpoint, the path of minimum time will give us the law of reflection.
The actual path of minimum time is the straight line that avoids the mirror,
of course. If we force the light to bounce off of the mirror, then we have to
figure out where it will hit, where
There are two ways to solve this problem. We can use geometry and visual intuition, or we can use calculus.
First, recall this fact from the problem text:
Light travels in a straight line in any particular medium with a velocity that depends upon the medium.
There’s no medium change, so if there were no mirror in its path, the light beam would continue in a straight line. Instead of figuring out what the beam will do when it hits the mirror, reflect the endpoint across the mirror and draw a straight line between the start and “end” points:
The angle that the beam makes with the plane of the mirror is the same on both sides of the mirror.
Now reflect the the “end” point and the segment of the beam that’s crossed the
mirror back up. By symmetry,
We can also solve this with calculus. Because the beam doesn’t change media, its
speed
Set
\begin{equation} d(x_p) = \sqrt{y_1^2 + x_p^2} + \sqrt{(x_2 - x_p)^2 + y_2^2} \end{equation}
For practice, we can also define this function in Scheme.
(define ((total-distance x1 y1 x2 y2) xp)
(+ (sqrt (+ (square (+ x1 xp))
(square y1)))
(sqrt (+ (square (- x2 (+ x1 xp)))
(square y2)))))Here’s the function again, generated from code, with general
(->tex-equation
((total-distance 'x_1 'y_1 'x_2 'y_2) 'x_p))#+RESULTS[084acf42d4fe771c97db9cf39e92c75383662d30]: \begin{equation} \sqrt{{{x}1}2 + 2 {x}1 {x}p + {{x}p}2 + {{y}1}2} + \sqrt{{{x}1}2 - 2 {x}1 {x}2 + 2 {x}1 {x}p + {{x}2}2 - 2 {x}2 {x}p + {{x}p}2 + {{y}2}2} \end{equation}
To find the
- take the derivative with respect to
$x_p$ , - set it equal to 0 and
- solve for
$x_p$ .
The derivative will look cleaner in code if we keep the components of the sum separate and prevent Scheme from “simplifying”. Redefine the function to return a tuple:
(define ((total-distance* x1 y1 x2 y2) xp)
(up (sqrt (+ (square (+ x1 xp))
(square y1)))
(sqrt (+ (square (- x2 (+ x1 xp)))
(square y2)))))Here are the sum components:
(->tex-equation
((total-distance* 0 'y_1 'x_2 'y_2) 'x_p))#+RESULTS[8080e49ee342b7a2a69c9c84337c37bc473a3c58]: \begin{equation} \begin{pmatrix} \displaystyle{ \sqrt{{{x}p}2 + {{y}1}2}} \cr \cr \displaystyle{ \sqrt{{{x}2}2 - 2 {x}2 {x}p + {{x}p}2 + {{y}2}2}}\end{pmatrix} \end{equation}
Taking a derivative is easy with scmutils. Just wrap the function in D:
(let* ((distance-fn (total-distance* 0 'y_1 'x_2 'y_2))
(derivative (D distance-fn)))
(->tex-equation
(derivative 'x_p)))#+RESULTS[5bbf36ca4a362ee2f2d2423071a6f818c8c93cab]: \begin{equation} \begin{pmatrix} \displaystyle{ {{{x}p}\over {\sqrt{{{x}p}2 + {{y}1}2}}}} \cr \cr \displaystyle{ {{ - {x}2 + {x}p}\over {\sqrt{{{x}2}2 - 2 {x}2 {x}p + {{x}p}2 + {{y}2}2}}}}\end{pmatrix} \end{equation}
The first component is the base of base
The bottom component is
\begin{equation} \label{eq:reflect-laws} cos θ_1 = cos θ_2 \implies θ_1 = θ_2 \end{equation}
This description in terms of the two incident angles isn’t so obvious from the Scheme code. Still, you can use Scheme to check this result.
If the two angles are equal, then the left and right triangles are similar, and the ratio of each base to height is equal:
\begin{equation} \label{eq:reflect-ratio} {x_p \over y_1} = {{x_2 - x_p} \over y_2} \end{equation}
Solve for
\begin{equation} \label{eq:reflect-ratio2} x_p = {{y_1 x_2} \over {y_1 + y_2}} \end{equation}
Plug this in to the derivative of the original total-distance function, and we
find that the derivative equals 0, as expected:
(let* ((distance-fn (total-distance 0 'y_1 'x_2 'y_2))
(derivative (D distance-fn)))
(->tex-equation
(derivative (/ (* 'y_1 'x_2) (+ 'y_1 'y_2)))))#+RESULTS[535d1b50ac55ba86347a21920c8bbf87153148eb]: \begin{equation} 0 \end{equation}
If a beam of light travels in a way that minimizes total distance (and therefore time in a constant medium), then it will reflect off of a mirror with the same angle at which it arrived. The law of reflection holds.
The law of refraction is also called Snell’s law. Here’s the description from the footnote:
Refraction is described by Snell’s law: when light passes from one medium to another, the ratio of the sines of the angles made to the normal to the interface is the inverse of the ratio of the refractive indices of the media. The refractive index is the ratio of the speed of light in the vacuum to the speed of light in the medium.
First we’ll tackle this with calculus.
The setup here is slightly different. We have a light beam traveling from one
medium to another and changing speeds at a boundary located
The refractive index
Time is distance over speed, so the total time that the beam spends between the
start and end points as a function of
\begin{equation} \begin{aligned} t(y_p) & = {c \sqrt{a^2 + y_p^2}\over v_1} + {c \sqrt{(x_2 - x_p)^2 + y_2^2} \over v_2} \cr & = {n_1 \over c} \sqrt{a^2 + y_p^2} + {n_2 \over c} \sqrt{(x_2 - x_p)^2 + y_2^2} \end{aligned} \end{equation}
Take the derivative:
\begin{equation} Dt(y_p) = {1 \over c} \left({n_1 y_p \over \sqrt{a^2 + y_p^2}} - {n_2 (x_2 - x_p) \over \sqrt{(x_2 - x_p)^2 + y_2^2}}\right) \end{equation}
Set the derivative equal to 0 and split terms:
\begin{equation} \label{eq:almost-snell} {n_1 y_p \over \sqrt{a^2 + y_p^2}} = {n_2 (x_2 - x_p) \over \sqrt{(x_2 - x_p)^2 + y_2^2}} \end{equation}
Similar to the law of reflection’s result, each term (up to its
Equation \eqref{eq:almost-snell} simplifies to:
\begin{equation} n_1 sin θ_1 = n_2 sin θ_2 \end{equation}
Rearranging yields Snell’s law:
\begin{equation} {n_1 \over n_2} = {sin θ_2 \over sin θ_1} \end{equation}
I won’t recreate this here, but the Feynman Lectures on Physics, in Lecture 26, has a fantastic discussion about, and derivation of, the law of refraction using no calculus, just geometry. I highly recommend you check out that lecture. Feynman lays out a number of examples of how the principle of least time is not just a restatement of the optical rules we already knew.
You can use the idea to guess what shape of mirror you’d want to build to focus many light rays on a single point (a parabola), or how you might force all light rays coming out of a single point to meet up again at another point (build a converging lens).
This whole area of optics and least time has obsessed scientists for hundreds of years. Spend a few minutes poking around and see what you find.
:header-args+: :tangle ch1/sec1-4.scm :comments orgI don’t plan on doing this for every section in the book, but section 1.4 is the first place where we’re introduced to Scheme, so I followed along and made a few notes.
This is the first demo of how any of this stuff works, starting on page 15. Here’s our first Lagrangian, super simple.
(define ((L-free-particle mass) local)
(let ((v (velocity local)))
(* 1/2 mass (dot-product v v))))L-free-particle is a function that takes some mass and returns a new
function. The new function takes an instance of a “local tuple” and returns the
value of the “Lagrangian”. This is the function that you query at every point
along some evolving path in configuration space. For realizable physical paths,
the integral of this function should by minimized, or stationary.
Why? That’s what we’re trying to develop here.
Suppose we let
(define q
(up (literal-function 'x)
(literal-function 'y)
(literal-function 'z)))The value
(->tex-equation
((Gamma q) 't))#+RESULTS[1f4aaac455bf48bd20965b4268009969bd7fd58e]: \begin{equation} \begin{pmatrix} \displaystyle{ t} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ x\left( t \right)} \cr \cr \displaystyle{ y\left( t \right)} \cr \cr \displaystyle{ z\left( t \right)}\end{pmatrix}} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ Dx\left( t \right)} \cr \cr \displaystyle{ Dy\left( t \right)} \cr \cr \displaystyle{ Dz\left( t \right)}\end{pmatrix}}\end{pmatrix} \end{equation}
This is just
Composing the Langrangian with
(->tex-equation
((compose (L-free-particle 'm) (Gamma q)) 't))#+RESULTS[49b01dd30b3d679e70016f72a5e51c78ecbf6c38]: \begin{equation} {{1}\over {2}} m {\left( Dz\left( t \right) \right)}2 + {{1}\over {2}} m {\left( Dy\left( t \right) \right)}2 + {{1}\over {2}} m {\left( Dx\left( t \right) \right)}2 \end{equation}
This particular formula is written in terms of
This function calculates the action
(define (Lagrangian-action L q t1 t2)
(definite-integral (compose L (Gamma q)) t1 t2))Here’s an example path that a particle might take, moving along a straight line
as a function of
(define (test-path t)
(up (+ (* 4 t) 7)
(+ (* 3 t) 5)
(+ (* 2 t) 1)))Calculate the action for a particle of mass 3, between
(Lagrangian-action (L-free-particle 3) test-path 0.0 10.0)#+RESULTS[78be3ee817bba9527071b0d765c261f632735187]:
#| 435. |#
This happens to be the minimal action, since the path we provided was a uniform path and the Lagrangian was for a free particle. If we’d provided a different path, we would still get an action. Just not a stationary action. Infinitesimal wiggles would change the action.
This exercise has us calculating the actual value of the action along some realizable path taken by a free particle.
For a free particle an appropriate Lagrangian is
\begin{equation} \label{eq:14lagrangian} L(t, x, v) = {1 \over 2}mv^2 \end{equation}
Suppose that x is the constant-velocity straight-line path of a free particle, such that
$x_a = x(t_a)$ and$x_b = x(t_b)$ . Show that the action on the solution path is
\begin{equation} \label{eq:14result} {m \over 2}{{(x_b - x_a)^2} \over {t_b - t_a}} \end{equation}
I’m not sure I see the point of this exercise, for developing intuition about
Langrangian mechanics. I think it may be here to make sure we understand that
we’re not minimizing the function
The velocity is constant between the two points, so it must be equal to the difference in position over the difference in time:
\begin{equation} \label{eq:constant-v} v = {{x(t_b) - x(t_a)} \over {t_b - t_a}} = {{x_b - x_a} \over {t_b - t_a}} \end{equation}
The action is equal to the integral of
\begin{equation} \label{eq:2} \begin{aligned} S[q](t_a, t_b) & = ∫t_at_b L(t, x, v) dx \cr & = ∫t_at_b {1 \over 2}mv(t)^2 dx \cr & = {m \over 2}{v(t)^2 t} \Bigr|t_at_b \cr & = {m \over 2}{v(t_b)^2 t_b - v(t_a)^2 t_a} \end{aligned} \end{equation}
The velocity is constant, so substitute in equation \eqref{eq:constant-v} and simplify:
\begin{equation} \label{eq:4} \begin{aligned} S[q](t_a, t_b) & = {m \over 2}{({{x_b - x_a} \over {t_b - t_a}})^2 (t_b - t_a)} \cr & = {m \over 2}{(x_b - x_a)^2 \over {t_b - t_a}} \end{aligned} \end{equation}
Boom, solution achieved.
:header-args+: :tangle ch1/min-action-paths.scm :comments orgThis section takes us through an example action calculation on a path with an adjustable “variation”, or wiggle. We should see that, if we consider a “realizable path”, then any wiggle we add will increase the calculated action.
The only restriction on a variation is that it can’t affect the endpoints of the realizable path. the times and positions of the start and end of the path are pinned.
make-eta takes some path
(define ((make-eta nu t1 t2) t)
(* (- t t1) (- t t2) (nu t)))Next, define a function that calculates the Lagrangian for a free particle, like
before, but adds in the path variation
(define ((varied-free-particle-action mass q nu t1 t2) epsilon)
(let ((eta (make-eta nu t1 t2)))
(Lagrangian-action (L-free-particle mass)
(+ q (* epsilon eta))
t1
t2)))Consider some variation like make-eta to pin its endpoints)
is larger (top entry) than the action of the non-varied path (bottom entry), as
expected:
<<test-path>>
(let ((action-fn (varied-free-particle-action 3.0 test-path
(up sin cos square)
0.0 10.0)))
(->tex-equation
(up (action-fn 0.001)
(action-fn 0))))#+RESULTS[10c4c8e6a40701b7176b4b1b2e9dd30b80185e6f]: #| test-path |#
\begin{equation} \begin{pmatrix} \displaystyle{ 436.2912142857153} \cr \cr \displaystyle{ 435.}\end{pmatrix} \end{equation}
What value of
We can search over values of minimize function:
(let ((action-fn (varied-free-particle-action
3.0 test-path
(up sin cos square)
0.0 10.0)))
(->tex-equation
(car (minimize action-fn -2.0 1.0))))#+RESULTS[87d808676b8bc1d49dd05d4ce8126fb3f223c1ce]: \begin{equation} 5.134781488891349e-15 \end{equation}
The result shows that the minimum action occurs at
Is it possible to use this principle to actually find a path, instead of simply checking it?
First we need a function that builds a path. This version generates a path of
individual points, bracketed by the supplied start and end points
(define (make-path t0 q0 t1 q1 qs)
(let ((n (length qs)))
(let ((ts (linear-interpolants t0 t1 n)))
(Lagrange-interpolation-function
(append (list q0) qs (list q1))
(append (list t0) ts (list t1))))))The next function sort-of-composes make-path and Lagrangian-action into a
function that takes
(define ((parametric-path-action L t0 q0 t1 q1) qs)
(let ((path (make-path t0 q0 t1 q1 qs)))
(Lagrangian-action L path t0 t1)))Finally, find-path takes the previous function’s arguments, plus a parameter
(define (find-path L t0 q0 t1 q1 n)
(let ((initial-qs (linear-interpolants q0 q1 n)))
(let ((minimizing-qs
(multidimensional-minimize
(parametric-path-action L t0 q0 t1 q1)
initial-qs)))
(make-path t0 q0 t1 q1 minimizing-qs))))Let’s test it out with a Lagrangian for a one dimensional harmonic oscillator
with spring constant
(define ((L-harmonic m k) local)
(let ((q (coordinate local))
(v (velocity local)))
(- (* 1/2 m (square v))
(* 1/2 k (square q)))))Now we invoke the procedures we’ve built, and plot the final, path-minimizing trajectory.
(define harmonic-path
(find-path (L-harmonic 1.0 1.0) 0.0 1.0 :pi/2 0.0 3))
(define win2 (frame 0.0 :pi/2 0 1))
(plot-function win2 harmonic-path 0 :pi (/ :pi 100))The path looks like a harmonic oscillator that starts high and bounces down,
after
The goal of this The goal of this exercise is to watch the minimization process that we just discussed proceed, from the initial guess of a straight-line path to the final, natural looking harmonic oscillation.
The exercise states:
We can watch the progress of the minimization by modifying the procedure parametric-path-action to plot the path each time the action is computed.
The functions the authors provide in the exercise define a window, and then a
version of parametric-path-action that updates the graph as it minimizes:
(define win2 (frame 0.0 :pi/2 0.0 1.2))
<<L-harmonic>>
(define ((parametric-path-action Lagrangian t0 q0 t1 q1)
intermediate-qs)
(let ((path (make-path t0 q0 t1 q1 intermediate-qs)))
;; display path
(graphics-clear win2)
(plot-function win2 path t0 t1 (/ (- t1 t0) 100))
;; compute action
(Lagrangian-action Lagrangian path t0 t1)))Run the minimization with the same parameters as in the previous section:
(find-path (L-harmonic 1.0 1.0) 0.0 1.0 :pi/2 0.0 2)and watch the plot update:
The authors have lightly demonstrated that plausible-looking paths have stationary action between fixed endpoints. What happens if we overconstrain the problem?
The exercise asks:
Suppose we try to obtain a path by minimizing an action for an impossible problem. For example, suppose we have a free particle and we impose endpoint conditions on the velocities as well as the positions that are inconsistent with the particle being free. Does the formalism protect itself from such an unpleasant attack? You may find it illuminating to program it and see what happens.
I spent a good amount of time thinking about this exercise. When I attempted to read this book in 2015, I found it very confusing.
Let’s say you take, as the authors suggest, some path, and impose velocity constraints on the endpoints in addition to the required position constraints.
Usually, you constrain the coordinates at each endpoint and force a path that minimizes the action between two times. What does it mean to impose velocity conditions?
The key is to realize that on the computer, you’re forcing a path to be composed of a bunch of discrete points. If you can force a point into the path that is NOT controlled by the optimizer, then you can force a velocity at some point in the path that makes no sense for minimal action.
Let’s define a new version of parametric-path-action that also takes an offset
for the initial and final points. We’ll force the first and last intermediate
point to be equal to the start and end points, plus some offset we can supply to
the function.
Then, we can try to find an action-minimizing path, but force the optimizer to deal with not just our endpoint conditions, but these two extra points as well. Forcing two points on each end will force an initial velocity condition. An offset of 0 would be equivalent to imposing a velocity of 0 at the start.
I simply proceeded with the implementation, but I’d recommend you take a minute to consider what you think will happen here. A hint is that the code is attempting to minimize action, given the constraint of the actual Lagrangian. What will it do when it’s forced to battle with a new exterior constraint, not captured in the Lagrangian?
Here’s the implementation of the modification described earlier:
(define (((parametric-path-action* win)
Lagrangian t0 q0 offset0 t1 q1 offset1)
intermediate-qs)
;; See the two new points?
(let ((intermediate-qs* (append (list (+ q0 offset0))
intermediate-qs
(list (+ q1 offset1)))))
(let ((path (make-path t0 q0 t1 q1 intermediate-qs*)))
;; display path
(graphics-clear win)
(plot-function win path t0 t1 (/ (- t1 t0) 100))
;; compute action
(Lagrangian-action Lagrangian path t0 t1))))You might try a similar trick by modifying the first and last entries of
intermediate-qs instead of appending a point, but I suspect that the optimizer
would be able to figure out how to undo your offset. (Try this as an exercise.)
Next, a new version of find-path that passes the offsets through to the new
parametric-path-action*:
(define ((find-path* win) L t0 q0 offset0 t1 q1 offset1 n)
(let ((initial-qs (linear-interpolants q0 q1 n)))
(let* ((action (parametric-path-action* win))
(minimizing-qs
(multidimensional-minimize
(action L t0 q0 offset0 t1 q1 offset1)
initial-qs)))
(make-path t0 q0 t1 q1 minimizing-qs))))And finally, a function that can execute runs of our formalism-killing experiment.
(define (one-six offset0 offset1 n)
(let* ((tmax 10)
(win (frame -1 (+ tmax 1) -0.2 (+ 1.2 offset0 offset1)))
(find (find-path* win))
(L (L-free-particle 3.0))
(path (find L
0. 1. offset0
tmax 0. offset1
n)))
(Lagrangian-action L path 0 tmax)))one-six takes two offsets and runs the minimization routine against
L-free-particle, moving from position 1 to 0 over 10 seconds. n controls the
number of interpolation points that the system will use.
Internally, remember, parametric-path-action* will append two extra fixed
offset points to the n intermediate points that the optimizer gets to control.
Let’s run the code with 0 offsets and 3 interpolation points. Note that this should still distort the path, since we now have two fixed points at the start and end. This is effectively imposing a 0 velocity constraint at the beginning and end.
Here’s the code, and its output:
(one-six 0 0 3)
The path ends up looking almost sinusoidal, and takes a while to converge. This is the best polynomial that the system can come up with that matches the 7 points (3 interpolated, 2 offsets, 1 start and 1 end).
The actual realizable path should be a straight line between the two points. The initial velocity of
Here’s a small positive velocity imposed at the beginning, and 0 at the end:
(one-six 0.2 0 3)
The system takes longer to converge. Here’s a larger impulse of 0.5 at the beginning:
(one-six 0.5 0 3)
And a moderate negative velocity, just for fun:
(one-six -0.5 0 3)
The process __does__ converge, but this is only because we only used 3 intermediate points. If you bump up to 10 points, with this code:
(one-six -0.5 0 10)The optimization freezes.
What is going on here? Why does the minimizer converge?
With velocity constraints imposed, we’re no longer minimizing the action with respect to some Lagrangian. We’re minimizing the action given two constraints. You have the Lagrangian, and then the warring goal of the polynomial interpolation forcing a certain shape on the path. At some point, the minimizer breaks; internally it ends up pinned between two tugging constraints.
If you make the impulse too big or force too many intermediate points, then the war is too hardcore and the process never converges. But it’s important to note here that these are details of the computational process. This detail doesn’t break reality. (It would break your model of reality, as it did here, if you have constraints or forces that you don’t capture in the Lagrangian.)
If you do need to impose velocity conditions, it turns out you can use a Lagrangian that takes acceleration into account. This is discussed in Exercise 1.10.
This exercise asks us to prove various products of the variation operator
The product rule for variations states that:
\begin{equation} \label{eq:var-prod} δ_η (f g)[q] = δ_η f[q] g[q] + f[q] δ_η g[q] \end{equation}
Write out the left side explicitly, using the definition of
\begin{equation} \label{eq:var-prod-proof} δ_η (f g)[q] = limε → 0 \left( {f[q + ε\eta]g[q + ε\eta] - f[q]g[q]} \over ε \right) \end{equation}
Make the inspired move to add and subtract
\begin{equation} \label{eq:var-prod-proof2} δ_η (f g)[q] = limε → 0 \left( {g[q + ε\eta](f[q + ε\eta] - f[q])} \over ε \right) + f[q] limε → 0 \left( {(g[q + ε η] - g[q])} \over ε \right) \end{equation}
You might recognize that we’ve now isolated terms that look like
\begin{equation} \label{eq:var-prod2} δ_η (f g)[q] = δ_η f[q]\,g[q] + f[q]\,δ_η g[q] \end{equation}
The sum rule is easier. Our goal is to show that:
\begin{equation} \label{eq:var-sum} δ_η (f + g)[q] = δ_η f[q] + δ_η g[q] \end{equation}
Expand out the definition of the variation operator, regroup terms, allow
\begin{equation} \label{eq:var-sum-proof} \begin{aligned} δ_η (f + g)[q] & = limε → 0 \left( {(f[q + ε\eta] + g[q + ε\eta]) - (f[q] + g[q])} \over ε \right) \cr & = limε → 0 \left( {f[q + ε\eta] - f[q]} \over ε \right) + limε → 0 \left( {g[q + ε\eta] - g[q]} \over ε \right) \cr & = δ_η f[q] + δ_η g[q] \end{aligned} \end{equation}
Done!
We want to show that
\begin{equation} \label{eq:var-scalar} δ_η (c g)[q] = c δ_η g[q] \end{equation}
Expand out the definition of the variation operator:
\begin{equation} \label{eq:var-scalar-proof} \begin{aligned} δ_η (c g)[q] & = limε → 0 \left( {c f[q + ε\eta] - c f[q]} \over ε \right) \cr & = c limε → 0 \left( {f[q + ε\eta] - f[q]} \over ε \right) \cr & = c δ_η f[q] \end{aligned} \end{equation}
Done, since the limit operator preserves scalar multiplication.
The chain rule for variations states that:
\begin{equation} \label{eq:var-chain} δ_η h[q] = (DF ˆ g[q])\, δ_η g[q] \textrm{ with } h[q] = F ˆ g[q] \end{equation}
Expand this out using the definition of
\begin{equation} \label{eq:var-chain-proof} δ_η (F ˆ g[q]) = limε → 0 \left( {(F ˆ g[q + ε\eta]) - (F ˆ g[q])} \over ε \right) \end{equation}
Now multiply the term inside the limit by
\begin{equation} \label{eq:var-chain-proof2} \begin{aligned} δ_η (F ˆ g[q]) & = limε → 0 \left( {((F ˆ g[q + ε\eta]) - (F ˆ g[q]))({g[q + ε\eta] - g[q]})} \over {({g[q + ε\eta] - g[q]}) ε} \right) \cr & = limε → 0 \left( {(F ˆ g[q + ε\eta]) - (F ˆ g[q])} \over {g[q + ε\eta] - g[q]} \right) δ_η g[q] \end{aligned} \end{equation}
The remaining term inside the limit has the form of a derivative of some
function
\begin{equation} \label{eq:var-chain-proof3} Df(a) = limb → a \left( {f(b) - f(a)} \over {b - a} \right) \end{equation}
Where
Remember that that this is all function algebra, so composition here is
analogous to function application; so
\begin{equation} \label{eq:var-chain-proof4} δ_η (F ˆ g[q]) = (DF ˆ g[q])\, δ_η g[q] \end{equation}
We need to show the derivative can commute with a normal derivative of the
function that
\begin{equation} \label{eq:var-commute} D δ_η f[q] = δ_η g[q] \textrm{ with } g[q] = D(f[q]) \end{equation}
Expand the left side by the definition of
\begin{equation} \label{eq:var-commute-proof} D (δ_η f[q]) = D limε → 0 \left( {(f[q + ε\eta]) - (f[q])} \over ε \right) \end{equation}
The derivative
\begin{equation} \label{eq:var-commute-proof2} \begin{aligned} D (δ_η f[q]) & = limε → 0 \left( {D(f[q + ε\eta]) - D(f[q])} \over ε \right) \cr & = δ_η(D(f[q])) \end{aligned} \end{equation}
Our goal is achieved.
:header-args+: :tangle ch1/ex1-8.scm :comments orgThe goal here is to implement
Suppose we have a procedure
fthat implements a path-dependent function: for pathqand timetit has the value((f q) t). The procedure delta computes the variation$δ_η f[q](t)$ as the value of the expression((((delta eta) f) q) t). Complete the definition ofdelta:
After laboriously proving all of the properties above, the actual implementation feels so simple.
The key is equation 1.22 in the book:
\begin{equation} \label{eq:1-22} δ_η f[q] = limε → 0 \left( {g(ε) - g(0)} \over ε \right) = Dg(0) \end{equation}
Given
(define (((delta eta) f) q)
(let ((g (lambda (eps)
(f (+ q (* eps eta))))))
((D g) 0)))It’s almost spooky, that
Part B’s problem description gave us a path-dependent function similar to this one:
(define ((fn sym) q)
(let* ((Local (UP Real (UP* Real) (UP* Real)))
(F (literal-function sym (-> Local Real))))
(compose F (Gamma q))))I’ve modified it slightly to take in a symbol, since we’ll need to generate multiple functions for a few of the rules.
UP*?) – and returns a generic expression for a
path dependent function
The textbook also gives us this function from
(define q (literal-function 'q (-> Real (UP Real Real))))
(define eta (literal-function 'eta (-> Real (UP Real Real))))These weren’t easy to write down, but they really do constitute proofs. If you trust the system managing the algebra, then the equalities here are general. This is an area of computing I haven’t worked with much, but I’m left with the eery feeling that these are more powerful than any tests I might have decided to write, if I weren’t guided by this exercise.
Equation \eqref{eq:var-prod} states the product rule for variations. Here it is in code. I’ve implemented the right and left sides and subtracted them. As expected, the result is 0:
(let* ((f (fn 'f))
(g (fn 'g))
(de (delta eta)))
(let ((left ((de (* f g)) q))
(right (+ (* (g q) ((de f) q))
(* (f q) ((de g) q)))))
(->tex-equation
((- left right) 't))))#+RESULTS[a1c8ed29c32a4b30412af5c2bbd632edd7a63a42]: \begin{equation} 0 \end{equation}
The sum rule is similar. Here’s the Scheme implementation of equation \eqref{eq:var-sum}:
(let* ((f (fn 'f))
(g (fn 'g))
(de (delta eta)))
(let ((left ((de (+ f g)) q))
(right (+ ((de f) q)
((de g) q))))
(->tex-equation
((- left right) 't))))#+RESULTS[2109df50bd5a522b7140b5656b20b90ee52e5e63]: \begin{equation} 0 \end{equation}
Here’s equation \eqref{eq:var-scalar} in code. The sides are equal, so their difference is 0:
(let* ((g (fn 'g))
(de (delta eta)))
(let ((left ((de (* 'c g)) q))
(right (* 'c ((de g) q))))
(->tex-equation
((- left right) 't))))#+RESULTS[9feb0bdde1e870101cd1e1e41df0774ae2be5d6b]: \begin{equation} 0 \end{equation}
To compute the chain rule we’ll need a version of fn that takes the derivative
of the inner function:
(define ((Dfn sym) q)
(let* ((Local (UP Real (UP* Real) (UP* Real)))
(F (literal-function sym (-> Local Real))))
(compose (D F) (Gamma q))))For the Scheme implementation, remember that both fn and Dfn have
Here’s a check that the two sides of equation \eqref{eq:var-chain} are equal:
(let* ((h (fn 'F))
(dh (Dfn 'F))
(de (delta eta)))
(let ((left (de h))
(right (* dh (de Gamma))))
(->tex-equation
(((- left right) q) 't))))#+RESULTS[539d27f8dd21a8838b24cf4b3c24a157597e0a42]: \begin{equation} 0 \end{equation}
Our final test. Here’s equation \eqref{eq:var-commute} in code, showing that the derivative commutes with the variation operator:
(let* ((f (fn 'f))
(g (compose D f))
(de (delta eta)))
(let ((left (D ((de f) q)))
(right ((de g) q)))
(->tex-equation
((- left right) 't))))#+RESULTS[5faf1227403cf244932447535a4366a066995299]: \begin{equation} 0 \end{equation}
This exercise asks us to derive Lagrange’s equations in steps for three different systems. Is this conscionable, given that we have the computer available to do the algebra for us? I think that this exercise is good practice for understanding the syntax, and maybe for refreshing your memory of how to symbolically manipulate derivatives.
But I’m feeling more and more that we’re in the middle of a thicket of exercises that are smugly attempting to show us how bad life is with pencil and paper, and how nice a computer algebra system can be. These paper-and-pencil problems are going to get harder and harder, while they stay trivial in Scheme.
Decide for yourself. Exercise 1.12 solves implements each Lagrangian in Scheme and demonstrates the steps required to get to Lagrange’s equations; I do buy that it would be difficult to do this without a good handle on the syntax.
Give it a try, then go examine exercise 1.12.
:header-args+: :tangle ch1/ex1-10.scm :comments orgThe only reason that the Lagrangians we’ve been considering don’t take any local tuple components beyond velocity is that the physics of our universe seems concerned with updating velocities and nothing beyond. Newton’s second law gives us an update rule for the velocities in a system, and we picked the Lagrangian so that Lagrange’s equations would match Newton’s second law.
But the formula for action works just as well if the Lagrangian takes many, or infinite, derivatives of the original coordinates. This exercise asks us to:
Derive Lagrange’s equations for Lagrangians that depend on accelerations. In particular, show that the Lagrange equations for Lagrangians of the form
$L(t, q, \dot{q}, \ddot{q})$ with$\ddot{q}$ terms are\begin{equation} D^2(∂_3L ˆ Γ[q]) - D(∂_2L ˆ Γ[q]) + (∂_1L ˆ Γ[q]) = 0 \end{equation}
In other words, find the constraint that has to be true of the Lagrangian for the action to be stationary along the supplied path.
This derivation follows the derivation of the Lagrange equations from the book,
starting on page 28. Begin with the equation for action with an acceleration
argument to
\begin{equation} \begin{aligned} S[q](t_a, t_b) & = ∫t_at_b L(t, x(t), v(t), a(t)) dx \cr & = ∫t_at_b (L ˆ Γ[q]) \end{aligned} \end{equation}
apply the variation operator,
\begin{equation} δ_η S[q](t_a, t_b) = ∫t_at_b δ_η (L ˆ Γ[q]) \end{equation}
Expand the right side out using the chain rule for variations, equation \eqref{eq:var-chain}:
\begin{equation} ∫t_at_b δ_η (L ˆ Γ[q]) =∫t_at_b ((DL) ˆ Γ[q]) δ_η\Gamma[q] \end{equation}
From equations 1.20 and 1.21 in the book, we know that
\begin{equation} δ_η\Gamma[q] = (0, η, Dη, D^2η, \ldots) \end{equation}
Expand the chain rule up to the $n$th derivative of the coordinate:
\begin{equation} \begin{aligned} ∫t_at_b & ((DL) ˆ Γ[q]) δ_η\Gamma[q] = \cr & ∫t_at_b 0 + (∂_1L ˆ Γ[q])η + (∂_2L ˆ Γ[q])Dη + \ldots + (∂n + 1L ˆ Γ[q])D^nη \end{aligned} \end{equation}
Our goal now is to find some quantity inside the integral that doesn’t depend on
\begin{equation} ∫t_at_b (∂_2L ˆ Γ[q])Dη \end{equation}
Integrate by parts:
\begin{equation} ∫t_at_b (∂_2L ˆ Γ[q])Dη = (∂_2L ˆ Γ[q])η \Bigr|t_at_b - ∫t_at_b D(∂_2L ˆ Γ[q])η \end{equation}
The first of the two terms disappears, since, by definition,
\begin{equation} ∫t_at_b (∂_2L ˆ Γ[q])Dη = ∫t_at_b D(∂_2L ˆ Γ[q])η \end{equation}
And reducing the full equation to:
\begin{equation} \begin{aligned} ∫t_at_b & ((DL) ˆ Γ[q]) δ_η\Gamma[q] = \cr & ∫t_at_b ((∂_1L ˆ Γ[q]) - D(∂_2L ˆ Γ[q]))η + (∂_3L ˆ Γ[q])D^2η \cr & + \ldots + (∂n + 1L ˆ Γ[q])D^nη \end{aligned} \end{equation}
The original Lagrange equations are peeking at us, multiplied by
Can we get another term into that form and move it in with the original Lagrange equation terms? Take the next term and integrate by parts twice:
\begin{equation} \begin{aligned} ∫t_at_b (∂_3L ˆ Γ[q])D^2η & = (∂_3L ˆ Γ[q])Dη \Bigr|t_at_b - ∫t_at_b D(∂_2L ˆ Γ[q])Dη \cr & = (∂_3L ˆ Γ[q])Dη \Bigr|t_at_b - D(∂_2L ˆ Γ[q])η \Bigr|t_at_b + ∫t_at_b D^2(∂_2L ˆ Γ[q])η \end{aligned} \end{equation}
The second of the two definite evaluations disappears, since, as before,
If we require
The term collapses to:
\begin{equation} \begin{aligned} ∫t_at_b (∂_3L ˆ Γ[q])D^2η & = (∂_3L ˆ Γ[q])Dη \Bigr|t_at_b - D(∂_2L ˆ Γ[q])η \Bigr|t_at_b + ∫t_at_b D^2(∂_2L ˆ Γ[q])η \cr & = ∫t_at_b D^2(∂_2L ˆ Γ[q])η \end{aligned} \end{equation}
Fold this back in to the full equation:
\begin{equation} \begin{aligned} ∫t_at_b & ((DL) ˆ Γ[q]) δ_η\Gamma[q] = \cr & ∫t_at_b ((∂_1L ˆ Γ[q]) - D(∂_2L ˆ Γ[q]) + D^2(∂_3L ˆ Γ[q]))η \cr & + \ldots + (∂n + 1L ˆ Γ[q])D^nη \end{aligned} \end{equation}
For a Lagrangian of the form
\begin{equation} \begin{aligned} δ_η S[q](t_a, t_b) & = ∫t_at_b ((DL) ˆ Γ[q]) δ_η\Gamma[q] \cr & = ∫t_at_b ((∂_1L ˆ Γ[q]) - D(∂_2L ˆ Γ[q]) + D^2(∂_3L ˆ Γ[q]))η \end{aligned} \end{equation}
The goal was to find conditions under which the action is stationary, ie,
\begin{equation} (∂_1L ˆ Γ[q]) - D(∂_2L ˆ Γ[q]) + D^2(∂_3L ˆ Γ[q]) = 0 \end{equation}
This is the result we were seeking.
To keep going, we have to integrate by parts once more for each new term of the local tuple that the Lagrangian depends on. For each new term we gain a new constraint:
\begin{equation} Dn-1η(t_a) = Dn-1η(t_b) = 0 \end{equation}
And a new term in the ever-higher-dimensional Lagrange’s equations:
\begin{equation} (-1)n Dn(∂n+1L ˆ Γ[q]) \end{equation}
The fully general Lagrange’s equations are, for a Lagrangian that depends on
the local tuple up to the $n$th derivative of
\begin{equation} 0 = ∑i = 0^n(-1)^i Di(∂i+1L ˆ Γ[q]) \end{equation}
Constrained by, for all
\begin{equation} D^iη(t_a) = D^iη(t_b) = 0 \end{equation}
Equivalently, the constraint is that all derivatives of
Exercise 1.13 implements a procedure that generates the residual required by these higher-dimensional Lagrange equations in Scheme.
:header-args+: :tangle ch1/ex1-11.scm :comments orgThis exercise asks us to derive Kepler’s third law by considering a Langrangian that describes two particles rotating in a circular orbit around their center of mass at some rate.
Here’s the Lagrangian for “central force”, in polar coordinates. This is
rotational kinetic energy, minus some arbitrary potential
(define ((L-central-polar m V) local)
(let ((q (coordinate local))
(qdot (velocity local)))
(let ((r (ref q 0)) (phi (ref q 1))
(rdot (ref qdot 0)) (phidot (ref qdot 1)))
(- (* 1/2 m
(+ (square rdot) (square (* r phidot))))
(V r)))))This function defines gravitational potential energy:
(define ((gravitational-energy G m1 m2) r)
(- (/ (* G m1 m2) r)))What is the mass
(define (reduced-mass m1 m2)
(/ (* m1 m2)
(+ m1 m2)))If you want to see why the reduced mass has the form it does, check out this derivation.
The Lagrangian is written in terms of some angle
(define ((q r omega) t)
(let ((phi (* omega t)))
(up r phi)))Write the Lagrange equations, given
(let ((eqfn (Lagrange-equations
(L-central-polar (reduced-mass 'm1 'm2)
(gravitational-energy 'G 'm1 'm2)))))
(->tex-equation
((eqfn (q 'a 'n)) 't)))#+RESULTS[e8565d29067487b8460bcb96e1a427d9eef22a0c]: \begin{equation} \begin{bmatrix} \displaystyle{ {{ - {a}3 ⋅ m1 ⋅ m2 ⋅ {n}2 + G {m1}2 ⋅ m2 + G ⋅ m1 ⋅ {m2}2}\over {{a}2 ⋅ m1 + {a}2 ⋅ m2}}} \cr \cr \displaystyle{ 0}\end{bmatrix} \end{equation}
These two entries are residuals, equal to zero. Stare at the top residual and you might notice that you can can factor out:
- the reduced mass, and
- a factor of
$1 \over a^2$
Manually factor these out:
(let ((eqfn (Lagrange-equations
(L-central-polar (reduced-mass 'm1 'm2)
(gravitational-energy 'G 'm1 'm2)))))
(->tex-equation
(* ((eqfn (q 'a 'n)) 't)
(/ (square 'a)
(reduced-mass 'm1 'm2)))))#+RESULTS[847a32694f3f07ce8f6d0e332bcded15b1e9abc5]: \begin{equation} \begin{bmatrix} \displaystyle{ - {a}3 {n}2 + G ⋅ m1 + G ⋅ m2} \cr \cr \displaystyle{ 0}\end{bmatrix} \end{equation}
And, boom, with some cleanup, we see Kepler’s third law:
\begin{equation} \label{eq:kepler3} n^2a^3 = G(m_1 + m_2) \end{equation}
:header-args+: :tangle ch1/ex1-12.scm :comments orgThis exercise asks us to write Scheme 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:
-
$∂_1 L ˆ Γ[q]$ , the generalized force -
$∂_2 L ˆ Γ[q]$ , the generalized momenta -
$D(∂_2 L ˆ Γ[q])$ , the derivative of our momenta - The Lagrange equations for the system.
(define (lagrange-equation-steps L q)
(let* ((p1 (compose ((partial 1) L) (Gamma q)))
(p2 (compose ((partial 2) L) (Gamma q)))
(dp2 (D p2)))
(->tex-equation
((up p1 p2 dp2 (- dp2 p1))
't))))These are the steps required on the road to a derivation of Lagrange’s equations.
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.
From the book:
An ideal planar pendulum consists of a bob of mass
$m$ connected to a pivot by a massless rod of length$l$ subject to uniform gravitational acceleration$g$ . A Lagrangian is$L(t, θ, \dot{θ}) = {1 \over 2} ml^2\dot{θ}^2 + mglcos θ$ . The formal parameters of$L$ are$t$ ,$θ$ , and$\dot{θ}$ ;$θ$ measures the angle of the pendulum rod to a plumb line and$\dot{θ}$ is the angular velocity of the rod.
Here is the Lagrangian described by the exercise:
(define ((L-pendulum m g l) local)
(let ((theta (coordinate local))
(theta_dot (velocity local)))
(+ (* 1/2 m (square l) (square theta_dot))
(* 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))#+RESULTS[aaf5812bee20b3464cf996ad648f7000e663545b]: \begin{equation} \begin{pmatrix} \displaystyle{ \left( - 1 \right) g l m sin\left( θ\left( t \right) \right)} \cr \cr \displaystyle{ {l}2 m Dθ\left( t \right)} \cr \cr \displaystyle{ {l}2 m {D}2θ\left( t \right)} \cr \cr \displaystyle{ g l m sin\left( θ\left( t \right) \right) + {l}2 m {D}2θ\left( t \right)}\end{pmatrix} \end{equation}
The final entry is the Lagrange equation, equal to
(let* ((L (L-pendulum 'm 'g 'l))
(theta (literal-function 'theta))
(eqs ((Lagrange-equations L) theta)))
(->tex-equation
((/ eqs (* 'm 'l))
't)))#+RESULTS[4342304f5c2ad636bb5fde728172b8977ea32776]: \begin{equation} g sin\left( θ\left( t \right) \right) + l {D}2θ\left( t \right) \end{equation}
This is the familiar equation of motion for a planar pendum.
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
$m$ moves in a two-dimensional potential$V(x, y) = {(x^2 + y^2) \over 2} + x^2 y - {y^3 \over 3}$ , where$x$ and$y$ are rectangular coordinates of the particle. A Lagrangian is$L(t;x, y; v_x, v_y) = {1 \over 2} m (v_x^2 + v_y^2) - V(x, y)$ .
I have no intuition for what this potential is, by the way. One term,
Define the Lagrangian to be the difference of the kinetic energy and some
potential
(define (((L-2d-potential m) V) 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 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:
(define (V q)
(let ((x (ref q 0))
(y (ref q 1)))
(- (+ (/ (+ (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)))#+RESULTS[934ef9c2a8a96a3dd4ad506e022deb8a6886de5c]: \begin{equation} \begin{pmatrix} \displaystyle{ \begin{bmatrix} \displaystyle{ - 2 y\left( t \right) x\left( t \right) - x\left( t \right)} \cr \cr \displaystyle{ {\left( y\left( t \right) \right)}2 - {\left( x\left( t \right) \right)}2 - y\left( t \right)}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ m Dx\left( t \right)} \cr \cr \displaystyle{ m Dy\left( t \right)}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ m {D}2x\left( t \right)} \cr \cr \displaystyle{ m {D}2y\left( t \right)}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ m {D}2x\left( t \right) + 2 y\left( t \right) x\left( t \right) + x\left( t \right)} \cr \cr \displaystyle{ m {D}2y\left( t \right) - {\left( y\left( t \right) \right)}2 + {\left( x\left( t \right) \right)}2 + y\left( t \right)}\end{bmatrix}}\end{pmatrix} \end{equation}
The final down-tuple gives us the Lagrange equations that
This problem is slightly more clear. From the book:
A Lagrangian for a particle of mass
$m$ constrained to move on a sphere of radius$R$ is$L(t; θ, φ; α, β) = {1 \over 2} m R^2(α^2+(β sin\theta)^2)$ . The angle$θ$ is the colatitude of the particle and$φ$ is the longitude; the rate of change of the colatitude is$α$ and the rate of change of the longitude is$β$ .
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
$x$ axis and fall off to nothing at the$y$ axis.
Here is the Lagrangian:
(define ((L-sphere m R) local)
(let* ((q (coordinate local))
(qdot (velocity local))
(theta (ref q 0))
(alpha (ref qdot 0))
(beta (ref qdot 1)))
(* 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)))#+RESULTS[5d05668e1d73bcd0f884eb8ba013c4eb56e72fda]: \begin{equation} \begin{pmatrix} \displaystyle{ \begin{bmatrix} \displaystyle{ {R}2 m sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) {\left( Dφ\left( t \right) \right)}2} \cr \cr \displaystyle{ 0}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ {R}2 m Dθ\left( t \right)} \cr \cr \displaystyle{ {R}2 m {\left( sin\left( θ\left( t \right) \right) \right)}2 Dφ\left( t \right)}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ {R}2 m {D}2θ\left( t \right)} \cr \cr \displaystyle{ 2 {R}2 m sin\left( θ\left( t \right) \right) Dθ\left( t \right) cos\left( θ\left( t \right) \right) Dφ\left( t \right) + {R}2 m {\left( sin\left( θ\left( t \right) \right) \right)}2 {D}2φ\left( t \right)}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ - {R}2 m sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) {\left( Dφ\left( t \right) \right)}2 + {R}2 m {D}2θ\left( t \right)} \cr \cr \displaystyle{ 2 {R}2 m sin\left( θ\left( t \right) \right) Dθ\left( t \right) cos\left( θ\left( t \right) \right) Dφ\left( t \right) + {R}2 m {\left( sin\left( θ\left( t \right) \right) \right)}2 {D}2φ\left( t \right)}\end{bmatrix}}\end{pmatrix} \end{equation}
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,
(let* ((L (L-sphere 'm 'R))
(theta (literal-function 'theta))
(q (up theta (literal-function 'phi)))
(le ((Lagrange-equations L) q)))
(let ((eq1 (ref le 0))
(eq2 (ref le 1)))
(->tex-equation
((up (/ eq1 (* 'm (square 'R)))
(/ eq2 (* (sin theta) 'm (square 'R))))
't))))#+RESULTS[db89e66528e605f07cf1f57c845568eb0c67777c]: \begin{equation} \begin{pmatrix} \displaystyle{ - sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) {\left( Dφ\left( t \right) \right)}2 + {D}2θ\left( t \right)} \cr \cr \displaystyle{ 2 Dθ\left( t \right) cos\left( θ\left( t \right) \right) Dφ\left( t \right) + sin\left( θ\left( t \right) \right) {D}2φ\left( t \right)}\end{pmatrix} \end{equation}
These are the Lagrange equations for
This exercise completes exercise 1.10 by asking for implementations of the higher-order Lagrange equations that we derived. This was a nice Scheme exercise; I would argue that this implementation should exist in the standard library. But that would ruin the fun of the exercise…
From the book:
Write a procedure to compute the Lagrange equations for Lagrangians that depend upon acceleration, as in exercise 1.10. Note that Gamma can take an optional argument giving the length of the initial segment of the local tuple needed. The default length is 3, giving components of the local tuple up to and including the velocities.
Now that we know the math, the implementation is a straightforward extension of
the Lagrange-equations procedure presented in the book:
(define ((Lagrange-equations3 L) q)
(let ((state-path (Gamma q 4)))
(+ ((square D) (compose ((partial 3) L) state-path))
(- (D (compose ((partial 2) L) state-path)))
(compose ((partial 1) L) state-path))))Now it’s time to use the new function. From the book:
Use your procedure to compute the Lagrange equations for the Lagrangian
\begin{equation} L(t, x, v, a) = - {1 \over 2}mxa - {1 \over 2}kx^2 \end{equation}
Do you recognize the resulting equation of motion?
Here is the Lagrangian described in the problem:
(define ((L-1-13 m k) local)
(let ((x (coordinate local))
(a (acceleration local)))
(- (* -1/2 m x a)
(* 1/2 k (square x)))))Use the new function to generate the Lagrange equations. This call includes a
factor of
(->tex-equation
(- (((Lagrange-equations3 (L-1-13 'm 'k))
(literal-function 'x)) 't)))#+RESULTS[d9cb178c8c9c14831b8950938516b881e5599d78]: \begin{equation} k x\left( t \right) + m {D}2x\left( t \right) \end{equation}
This looks like the equation of motion for a classical harmonic oscillator. Again, I leave this problem with no new physical intuition for what is going on here, or what type of system would need an acceleration dependent Lagrangian. I suspect that we could build a harmonic oscillator for Lagrange equations of any order by properly tuning the Lagrangian. But I don’t know why this would be helpful.
Now, some more fun with Scheme. It just feels nice to implement Scheme procedures. From the book:
For more fun, write the general Lagrange equation procedure that takes a Lagrangian that depends on any number of derivatives, and the number of derivatives, to produce the required equations of motion.
As a reminder, this is the equation that we need to implement for each coordinate:
\begin{equation} 0 = ∑i = 0^n(-1)^i Di(∂i+1L ˆ Γ[q]) \end{equation}
There are two ideas playing together here. Each term takes an element from:
- an alternating sequence of
$1$ and$-1$ - a sequence of increasing-index
$D^i(∂_i L ˆ Γ[q])$ terms
The alternating
If you need to take
alternate generates a list of elems:
(define (cycle n elems)
(apply append (make-list n elems)))
(define (alternating n elems)
(let* ((l (length elems))
(times (quotient (+ n (-1+ l)) l)))
(list-head (cycle times elems) n)))Now, the general Lagrange-equations* implementation.
This function defines an internal function term that generates the $i$th term
of the derivative combination described above. This sequence is zipped with the
sequence of fold-left generates the sum.
(define ((Lagrange-equations* L n) q)
(let ((state-path (Gamma q (1+ n))))
(define (term i)
((expt D i)
(compose ((partial (1+ i)) L) state-path)))
(let ((terms (map term (iota n))))
(fold-left (lambda (acc pair)
(+ acc (apply * pair)))
0
(zip (alternating n '(1 -1))
(reverse terms))))))Generate the Lagrange equations from part b once more to check that we get the same result:
(->tex-equation
(- (((Lagrange-equations* (L-1-13 'm 'k) 3)
(literal-function 'x)) 't)))#+RESULTS[1393bf1b0b9cbfbd5c386cb8fcdd91c86cf38e32]: \begin{equation} k x\left( t \right) + m {D}2x\left( t \right) \end{equation}
There it is again, the harmonic oscillator. I don’t have any intuition for higher order Lagrangians, so I can’t cook up any further examples to test the implementation.
:header-args+: :tangle ch1/ex1-14.scm :comments orgLook carefully at what this exercise is asking us to do:
Check that the Lagrange equations for central force motion in polar coordinates and in rectangular coordinates are equivalent. Determine the relationship among the second derivatives by substituting paths into the transformation equations and computing derivatives, then substitute these relations into the equations of motion.
The punchline that we’ll encounter soon is that a coordinate transformation of
applied to some path function
\begin{equation} C ˆ Γ[q] = Γ[q’] \end{equation}
Because function composition is associative, instead of ever transforming the
path, you can push the coordinate transformation into the Lagrangian to generate
a new Lagrangian:
The section of textbook just before the exercise has given us two Lagrangians in
different coordinates – L-central-polar and L-rectangular – and generated
Lagrange equations from each.
Our task is to directly transform the Lagrange equations by substituting the first and second derivatives of the coordinate transformation expression into one of the sets of equations, and looking to see that it’s equivalent to the other.
Fair warning: this is much more painful than transforming the Lagrangian before generating the Lagrange equations. This exercise continues the theme of devastating you with algebra as a way to show you the horror that the later techniques were developed to avoid. Let us proceed.
Here are the two Lagrangians from the book:
(define ((L-central-rectangular m U) local)
(let ((q (coordinate local))
(v (velocity local)))
(- (* 1/2 m (square v))
(U (sqrt (square q))))))
(define ((L-central-polar m U) local)
(let ((q (coordinate local))
(qdot (velocity local)))
(let ((r (ref q 0)) (phi (ref q 1))
(rdot (ref qdot 0)) (phidot (ref qdot 1)))
(- (* 1/2 m
(+ (square rdot)
(square (* r phidot))) )
(U r)))))Here are the rectangular equations of motion:
(->tex-equation
(((Lagrange-equations
(L-central-rectangular 'm (literal-function 'U)))
(up (literal-function 'x)
(literal-function 'y)))
't)
"eq:rect-equations")#+RESULTS[3d6b30672d7d2fe77035ee8a5ae5b0f3046f2f90]: \begin{equation} \begin{bmatrix} \displaystyle{ {{m {D}2x\left( t \right) \sqrt{{\left( x\left( t \right) \right)}2 + {\left( y\left( t \right) \right)}2} + x\left( t \right) DU\left( \sqrt{{\left( x\left( t \right) \right)}2 + {\left( y\left( t \right) \right)}2} \right)}\over {\sqrt{{\left( x\left( t \right) \right)}2 + {\left( y\left( t \right) \right)}2}}}} \cr \cr \displaystyle{ {{m \sqrt{{\left( x\left( t \right) \right)}2 + {\left( y\left( t \right) \right)}2} {D}2y\left( t \right) + y\left( t \right) DU\left( \sqrt{{\left( x\left( t \right) \right)}2 + {\left( y\left( t \right) \right)}2} \right)}\over {\sqrt{{\left( x\left( t \right) \right)}2 + {\left( y\left( t \right) \right)}2}}}}\end{bmatrix} \label{eq:rect-equations} \end{equation}
And the polar Lagrange equations:
(->tex-equation
(((Lagrange-equations
(L-central-polar 'm (literal-function 'U)))
(up (literal-function 'r)
(literal-function 'phi)))
't))#+RESULTS[557893e62f632f0900e9ece883c176eaf5bcfd05]: \begin{equation} \begin{bmatrix} \displaystyle{ - m r\left( t \right) {\left( Dφ\left( t \right) \right)}2 + m {D}2r\left( t \right) + DU\left( r\left( t \right) \right)} \cr \cr \displaystyle{ m {D}2φ\left( t \right) {\left( r\left( t \right) \right)}2 + 2 m r\left( t \right) Dφ\left( t \right) Dr\left( t \right)}\end{bmatrix} \end{equation}
Once again, our goal is to show that, if you can write down coordinate transformations for the coordinates, velocities and accelerations and substitute them in to one set of Lagrange equations, the other will appear.
To do this by hand, take the coordinate transformation described in 1.64 in the book:
\begin{equation} \begin{split} x &= r cos φ \cr y &= r sin φ \end{split} \end{equation}
Note that
\begin{equation} \begin{split} Dx(t) &= Dr(t) cos φ(t) - r(t) Dφ(t) sin φ(t) \cr Dy(t) &= Dr(t) sin φ(t) + r(t) Dφ(t) cos φ(t) \end{split} \end{equation}
The rectangular equations of motion have second derivatives, so we need to keep going. This is too devastating to imagine doing by hand. Let’s move to Scheme.
Write the coordinate transformation for polar coordinates to rectangular in Scheme:
(define (p->r local)
(let* ((polar-tuple (coordinate local))
(r (ref polar-tuple 0))
(phi (ref polar-tuple 1))
(x (* r (cos phi)))
(y (* r (sin phi))))
(up x y)))Now use F->C, first described on page 46. This is a function that takes a
coordinate transformation like p->r and returns a new function that can
convert an entire local tuple from one coordinate system to another; the
The version that the book presents on page 46 can only generate a velocity
transformation given a coordinate transformation, but scmutils contains a more
general version that will convert as many path elements as you pass to it.
Here are the rectangular positions, velocities and accelerations, written in polar coordinates:
(let ((convert-path (F->C p->r))
(polar-path (up 't
(up 'r 'phi)
(up 'rdot 'phidot)
(up 'rdotdot 'phidotdot))))
(->tex-equation
(convert-path polar-path)))#+RESULTS[1bac37835829a8c17e06699ef20f2e42676725b9]: \begin{equation} \begin{pmatrix} \displaystyle{ t} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ r cos\left( φ \right)} \cr \cr \displaystyle{ r sin\left( φ \right)}\end{pmatrix}} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ - \dot{φ} r sin\left( φ \right) + \dot{r} cos\left( φ \right)} \cr \cr \displaystyle{ \dot{φ} r cos\left( φ \right) + \dot{r} sin\left( φ \right)}\end{pmatrix}} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ - {\dot{φ}}2 r cos\left( φ \right) - 2 \dot{φ} \dot{r} sin\left( φ \right) - \ddot{φ} r sin\left( φ \right) + \ddot{r} cos\left( φ \right)} \cr \cr \displaystyle{ - {\dot{φ}}2 r sin\left( φ \right) + 2 \dot{φ} \dot{r} cos\left( φ \right) + \ddot{φ} r cos\left( φ \right) + \ddot{r} sin\left( φ \right)}\end{pmatrix}}\end{pmatrix} \end{equation}
Ordinarily, it would be too heartbreaking to substitute these in to the rectangular equations of motion. The fact that we have Scheme on our side gives me the strength to proceed.
Write the rectangular Lagrange equations as a function of the local tuple, so we can call it directly:
(define (rect-equations local)
(let* ((q (coordinate local))
(x (ref q 0))
(y (ref q 1))
(v (velocity local))
(xdot (ref v 0))
(ydot (ref v 1))
(a (acceleration local))
(xdotdot (ref a 0))
(ydotdot (ref a 1))
(U (literal-function 'U)))
(up (/ (+ (* 'm xdotdot (sqrt (+ (square x) (square y))))
(* x ((D U) (sqrt (+ (square x) (square y))))))
(sqrt (+ (square x) (square y))))
(/ (+ (* 'm ydotdot (sqrt (+ (square x) (square y))))
(* y ((D U) (sqrt (+ (square x) (square y))))))
(sqrt (+ (square x) (square y)))))))Verify that these are, in fact, the rectangular equations of motion by passing in a symbolic rectangular local tuple:
(let ((rect-path (up 't
(up 'x 'y)
(up 'xdot 'ydot)
(up 'xdotdot 'ydotdot))))
(->tex-equation
(rect-equations rect-path)))#+RESULTS[5e06ee310a8c8e3036c17c5863647c80960dc568]: \begin{equation} \begin{pmatrix} \displaystyle{ {{m \ddot{x} \sqrt{{x}2 + {y}2} + x DU\left( \sqrt{{x}2 + {y}2} \right)}\over {\sqrt{{x}2 + {y}2}}}} \cr \cr \displaystyle{ {{m \ddot{y} \sqrt{{x}2 + {y}2} + y DU\left( \sqrt{{x}2 + {y}2} \right)}\over {\sqrt{{x}2 + {y}2}}}}\end{pmatrix} \end{equation}
Now use the p->r conversion to substitute each of the rectangular values above
with their associated polar values:
(let* ((convert-path (F->C p->r))
(polar-path (up 't
(up 'r 'phi)
(up 'rdot 'phidot)
(up 'rdotdot 'phidotdot)))
(local (convert-path polar-path)))
(->tex-equation
(rect-equations local)))#+RESULTS[577a75ce68ba364856b29d9b49e8f95c130ec027]: \begin{equation} \begin{pmatrix} \displaystyle{ - m {\dot{φ}}2 r cos\left( φ \right) - 2 m \dot{φ} \dot{r} sin\left( φ \right) - m \ddot{φ} r sin\left( φ \right) + m \ddot{r} cos\left( φ \right) + DU\left( r \right) cos\left( φ \right)} \cr \cr \displaystyle{ - m {\dot{φ}}2 r sin\left( φ \right) + 2 m \dot{φ} \dot{r} cos\left( φ \right) + m \ddot{φ} r cos\left( φ \right) + m \ddot{r} sin\left( φ \right) + DU\left( r \right) sin\left( φ \right)}\end{pmatrix} \end{equation}
Oh no. This looks quite different from the polar Lagrange equations above. What is the problem?
I had to stare at this for a long time before I saw what to do. Notice that the
terms we want from the polar Lagrange equations all seem to appear in the first
equation with a
\begin{equation} (cos φ)^2 + (sin φ)^2 = 1 \end{equation}
I realized that I could recover the first equation through a linear combination
of both terms. Multiply the first by
A similar trick recovers the second equation,given an extra factor of
(let* ((convert-path (F->C p->r))
(polar-path (up 't
(up 'r 'phi)
(up 'rdot 'phidot)
(up 'rdotdot 'phidotdot)))
(local (convert-path polar-path))
(eq (rect-equations local)))
(->tex-equation
(up (+ (* (cos 'phi) (ref eq 0))
(* (sin 'phi) (ref eq 1)))
(- (* 'r (cos 'phi) (ref eq 1))
(* 'r (sin 'phi) (ref eq 0))))))#+RESULTS[4a5dd5a22afc855770bba85efe4a6bffecdc6228]: \begin{equation} \begin{pmatrix} \displaystyle{ - m {\dot{φ}}2 r + m \ddot{r} + DU\left( r \right)} \cr \cr \displaystyle{ 2 m \dot{φ} r \dot{r} + m \ddot{φ} {r}2}\end{pmatrix} \end{equation}
This was a powerful lesson. We’re allowed to take a linear combination here
because each equation is a residual, equal to zero.
There is something important going on here, with the way we were able to remove
This is one of the more important exercises in the chapter. The problem asks for a proof that it’s possible to absorb a coordinate transformation directly into the Lagrangian. If you can do this, you can express your paths and your forces in whatever coordinates you like, so long as you can transition between them.
I also found that this exposed, and repaired, my weakness with the functional notation that Sussman and Wisdom have used in the book.
The problem states:
Show by direct calculation that the Lagrange equations for
$L’$ are satisfied if the Lagrange equations for$L$ are satisfied.
Imagine a coordinate transformation
\begin{equation} x = F(t, x’) \end{equation}
You might imagine that
Now imagine some Lagrangian
\begin{equation} L’ ˆ Γ[q’] = L ˆ Γ[q] \end{equation}
The function
\begin{equation} \label{eq:1-15-relations} L ˆ Γ[q] = L ˆ C ˆ Γ[q’] = L’ ˆ Γ[q’] \end{equation}
Function composition is associative, so two facts stare at us from \eqref{eq:1-15-relations}:
\begin{equation} \begin{aligned} L’ & = L ˆ C \cr Γ[q] & = C ˆ Γ[q’] \end{aligned} \end{equation}
This implies that if we want describe our path in some primed coordinate system
(imagine polar coordinates), but we only have a Lagrangian defined in an
unprimed coordinate system (like rectangular coordinates), if we can just write
down
The goal is to show that if the Lagrange equations hold in the unprimed coordinate system:
\begin{equation} \label{eq:1-15-lagrange} D(∂_2L ˆ Γ[q]) - (∂_1L ˆ Γ[q]) = 0 \end{equation}
Then they hold in the primed coordinate system:
\begin{equation} \label{eq:lagrange-prime} D(∂_2L’ ˆ Γ[q’]) - (∂_1L’ ˆ Γ[q’]) = 0 \end{equation}
The approach we’ll take will be to:
- Write down the form of
$C$ , given some arbitrary$F$ - start calculating the terms of the Lagrange equations in the primed coordinate
system using
$L’ = L ˆ C$ - Keep an eye out for some spot where we can use our assumption that the Lagrange equations hold in the unprimed coordinates.
I’ll walk through a solution using “pen and paper”, then show how we can run this derivation using Scheme to help us along.
First, some Scheme tools that will help us in both cases.
Equation (1.77) in the book describes how to implement F->C on page 46.
The following function is a slight redefinition that allows us to use an
(define ((F->C* F) local)
(let ((t (time local))
(x (coordinate local))
(v (velocity local)))
(up t
(F t x)
(+ (((partial 0) F) t x)
(* (((partial 1) F) t x)
v)))))#+RESULTS[8ea33eda9107b7b0d6ce890f316eb453b2d96fca]:
#| F->C* |#
Next we define qprime, a
function that can represent our unprimed coordinate path function.
The types here all imply that the path has one real coordinate. I did this to make the types easier to understand; the derivation applies equally well to paths with many dimensions.
(define F
(literal-function 'F (-> (X Real Real) Real)))
(define C (F->C* F))
(define L
(literal-function 'L (-> (UP Real Real Real) Real)))
(define qprime
(literal-function 'qprime))When we apply
(->tex-equation
((compose C (Gamma qprime)) 't))#+RESULTS[ecef7fe872de385af827c689dbd0db76aa5319f0]: \begin{equation} \begin{pmatrix} \displaystyle{ t} \cr \cr \displaystyle{ F\left( t, {q}^′\left( t \right) \right)} \cr \cr \displaystyle{ D{q}^′\left( t \right) {∂}1F\left( t, {q}^′\left( t \right) \right) + {∂}0F\left( t, {q}^′\left( t \right) \right)}\end{pmatrix} \end{equation}
This looks correct. We can also transform the path before passing it to
(define ((to-q F qp) t)
(F t (qp t)))Subtract the two forms to see that they’re equivalent:
(->tex-equations
((- (compose C (Gamma qprime))
(Gamma (to-q F qprime)))
't))#+RESULTS[10aea4a5b833aa00e0a1c63049d844dc37528289]: \begin{equation} \begin{pmatrix} \displaystyle{ 0} \cr \cr \displaystyle{ 0} \cr \cr \displaystyle{ 0}\end{pmatrix} \end{equation}
Now that we know Lprime:
(define q (to-q F qprime))
(define Lprime (compose L C))Begin by calculating the components of the Lagrange equations in equation
\eqref{eq:lagrange-prime}. Examine the
As we discussed above, function composition is associative, so:
\begin{equation} \label{eq:c-l} (L ˆ C) ˆ Γ[q’] = L’ ˆ Γ[q’] \implies L’ = L ˆ C \end{equation}
Substituting
\begin{equation} ∂_2L’ = ∂_2(L ˆ C) = ((DL) ˆ C) ∂_2 C \end{equation}
I found the next step troubling until I became more comfortable with the functional notation.
The tuple algebra described in Chapter 9 defines multiplication between an up and down tuple as a dot product, or a “contraction” in the book’s language. This means that we can expand out the product above:
\begin{equation} (DL ˆ C)∂_2 C = (∂_0L ˆ C)(I_0 ˆ ∂_2 C) + (∂_1L ˆ C)(I_1 ˆ ∂_2 C) + (∂_2L ˆ C)(I_2 ˆ ∂_2 C) \end{equation}
Example the value of
(->tex-equation
(((partial 2) C) (up 't 'xprime 'vprime))
"eq:p2c")#+RESULTS[2ad70f16eca982284368a2f70f31de3b2de41025]: \begin{equation} \begin{pmatrix} \displaystyle{ 0} \cr \cr \displaystyle{ 0} \cr \cr \displaystyle{ {∂}1F\left( t, {x}^′ \right)}\end{pmatrix} \label{eq:p2c} \end{equation}
The first two components are 0, leaving us with:
\begin{equation} ∂_2 L’ = (DL ˆ C)∂_2 C = (∂_2L ˆ C)(I_2 ˆ ∂_2 C) \end{equation}
Compose this quantity with
\begin{equation} \begin{aligned} ∂_2L’ ˆ Γ[q’] & = (∂_2L ˆ C)(I_2 ˆ ∂_2 C) ˆ Γ[q’] \cr & = (∂_2L ˆ C ˆ Γ[q’])(I_2 ˆ ∂_2 C ˆ Γ[q’]) \cr & = (∂_2L ˆ Γ[q])(I_2 ˆ ∂_2 C ˆ Γ[q’]) \end{aligned} \end{equation}
Take the derivative (with respect to time, remember, from the types):
\begin{equation} D(∂_2L’ ˆ Γ[q’]) = D\left[(∂_2L ˆ Γ[q])(I_2 ˆ ∂_2 C ˆ Γ[q’])\right] \end{equation}
Substitute the second term using \eqref{eq:p2c}:
\begin{equation} D(∂_2L’ ˆ Γ[q’]) = D\left[(∂_2L ˆ Γ[q])∂_1F(t, q’(t))\right] \end{equation}
Expand using the product rule:
\begin{equation} \label{eq:ex1-15-dp2l} D(∂_2L’ ˆ Γ[q’]) = \left[ D(∂_2L ˆ Γ[q]) \right]∂_1F(t, q’(t)) + (∂_2L ˆ Γ[q])D\left[ ∂_1F(t, q’(t)) \right] \end{equation}
A term from the unprimed Lagrange’s equation is peeking. Notice this, but don’t make the substitution just yet.
Next, expand the
\begin{equation} ∂_1L’ = ∂_1(L ˆ C) = ((DL) ˆ C) ∂_1 C \end{equation}
Calculate
(->tex-equation
(((partial 1) C) (up 't 'xprime 'vprime)))#+RESULTS[40670631406e70c2545a7426039a01214769d4fd]: \begin{equation} \begin{pmatrix} \displaystyle{ 0} \cr \cr \displaystyle{ {∂}1F\left( t, {x}^′ \right)} \cr \cr \displaystyle{ {v}^′ {{∂}1}2\left( F \right)\left( t, {x}^′ \right) + \left( {∂}0 {∂}1 \right)\left( F \right)\left( t, {x}^′ \right)}\end{pmatrix} \end{equation}
Expand the chain rule out and remove 0 terms, as before:
\begin{equation} \begin{aligned} ∂_1L’ & = ((DL) ˆ C) ∂_1 C \cr & = (∂_0L ˆ C)(I_0 ˆ ∂_1 C) + (∂_1L ˆ C)(I_1 ˆ ∂_1 C) + (∂_2L ˆ C)(I_2 ˆ ∂_1 C) \cr & = (∂_1L ˆ C)(I_1 ˆ ∂_1 C) + (∂_2L ˆ C)(I_2 ˆ ∂_1 C) \end{aligned} \end{equation}
Compose
\begin{equation} \begin{aligned} ∂_1L’ ˆ Γ[q’] & = (∂_1L ˆ Γ[q])(I_1 ˆ ∂_1 C ˆ Γ[q’]) + (∂_2L ˆ Γ[q])(I_2 ˆ ∂_1 C ˆ Γ[q’]) \cr & = (∂_1L ˆ Γ[q])(∂_1F(t, q’(t))) + (∂_2L ˆ Γ[q])D(∂_1F(t, q’(t))) \end{aligned} \end{equation}
We now have both components of the primed Lagrange equations from \eqref{eq:lagrange-prime}.
Subtract the two terms, extract common factors and use our assumption \eqref{eq:1-15-lagrange} that the original Lagrange equations hold:
\begin{equation} \begin{aligned} D(∂_2L’ ˆ Γ[q’]) - (∂_1L’ ˆ Γ[q’]) & = \left[ D(∂_2L ˆ Γ[q]) - (∂_1L ˆ Γ[q]) \right] ∂_1F(t, q’(t)) \cr & + \left[ (∂_2L ˆ Γ[q]) - (∂_2L ˆ Γ[q]) \right] D(∂_1F(t, q’(t))) \cr & = \left[ D(∂_2L ˆ Γ[q]) - (∂_1L ˆ Γ[q]) \right] ∂_1F(t, q’(t)) \cr & = 0 \end{aligned} \end{equation}
And boom, we’re done! The primed Lagranged equations \eqref{eq:lagrange-prime} hold if the unprimed equations \eqref{eq:1-15-lagrange} hold.
I’m not sure what to make of the new constant terms. The new Lagrange equations
are scaled by
Can we use Scheme to pursue the same derivation? If we can write the relationships of the derivation in code, then we’ll have a sort of computerized proof that the primed Lagrange equations are valid.
First, consider
(->tex-equation
((compose ((partial 1) Lprime) (Gamma qprime))
't))#+RESULTS[5eac70f3edfcf5feeb3a2eb9b98d89646f21ade3]: \begin{equation} D{q}^′\left( t \right) {∂}2L\left( \begin{pmatrix} \displaystyle{ t} \cr \cr \displaystyle{ F\left( t, {q}^′\left( t \right) \right)} \cr \cr \displaystyle{ D{q}^′\left( t \right) {∂}1F\left( t, {q}^′\left( t \right) \right) + {∂}0F\left( t, {q}^′\left( t \right) \right)}\end{pmatrix} \right) {{∂}1}2\left( F \right)\left( t, {q}^′\left( t \right) \right) + {∂}1F\left( t, {q}^′\left( t \right) \right) {∂}1L\left( \begin{pmatrix} \displaystyle{ t} \cr \cr \displaystyle{ F\left( t, {q}^′\left( t \right) \right)} \cr \cr \displaystyle{ D{q}^′\left( t \right) {∂}1F\left( t, {q}^′\left( t \right) \right) + {∂}0F\left( t, {q}^′\left( t \right) \right)}\end{pmatrix} \right) + {∂}2L\left( \begin{pmatrix} \displaystyle{ t} \cr \cr \displaystyle{ F\left( t, {q}^′\left( t \right) \right)} \cr \cr \displaystyle{ D{q}^′\left( t \right) {∂}1F\left( t, {q}^′\left( t \right) \right) + {∂}0F\left( t, {q}^′\left( t \right) \right)}\end{pmatrix} \right) \left( {∂}0 {∂}1 \right)\left( F \right)\left( t, {q}^′\left( t \right) \right) \end{equation}
This is completely insane, and already unhelpful. The argument to
(define (->eq expr)
(write-string
(replace-all (->tex-equation* expr)
(->tex* ((Gamma q) 't))
"\\circ \\Gamma[q]")))Try again:
(->eq
((compose ((partial 1) Lprime) (Gamma qprime))
't))#+RESULTS[ee6c612cadd91730dd0142a4dd7476fda80d6e07]: \begin{equation} D{q}^′\left( t \right) {∂}2L\left( ˆ Γ[q] \right) {{∂}1}2\left( F \right)\left( t, {q}^′\left( t \right) \right) + {∂}1F\left( t, {q}^′\left( t \right) \right) {∂}1L\left( ˆ Γ[q] \right) + {∂}2L\left( ˆ Γ[q] \right) \left( {∂}0 {∂}1 \right)\left( F \right)\left( t, {q}^′\left( t \right) \right) \end{equation}
Ignore the parentheses around
The
(let* ((factor (((partial 1) F) 't (qprime 't))))
(->eq
((* factor (compose ((partial 1) L) (Gamma q)))
't)))#+RESULTS[31b49f8349da58b715e99a7cbcdf983dc4e12d1e]: \begin{equation} {∂}1F\left( t, {q}^′\left( t \right) \right) {∂}1L\left( ˆ Γ[q] \right) \end{equation}
Next, consider the
(->eq
((D (compose ((partial 2) Lprime) (Gamma qprime)))
't))#+RESULTS[9813f05b0fc4d5e7ad9ce035bca6e494cb087d84]: \begin{equation} {\left( D{q}^′\left( t \right) \right)}2 {∂}1F\left( t, {q}^′\left( t \right) \right) {{∂}1}2\left( F \right)\left( t, {q}^′\left( t \right) \right) {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) + D{q}^′\left( t \right) {\left( {∂}1F\left( t, {q}^′\left( t \right) \right) \right)}2 \left( {∂}1 {∂}2 \right)\left( L \right)\left( ˆ Γ[q] \right) + 2 D{q}^′\left( t \right) {∂}1F\left( t, {q}^′\left( t \right) \right) \left( {∂}0 {∂}1 \right)\left( F \right)\left( t, {q}^′\left( t \right) \right) {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) + {\left( {∂}1F\left( t, {q}^′\left( t \right) \right) \right)}2 {D}2{q}^′\left( t \right) {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) + {∂}0F\left( t, {q}^′\left( t \right) \right) {∂}1F\left( t, {q}^′\left( t \right) \right) \left( {∂}1 {∂}2 \right)\left( L \right)\left( ˆ Γ[q] \right) + D{q}^′\left( t \right) {∂}2L\left( ˆ Γ[q] \right) {{∂}1}2\left( F \right)\left( t, {q}^′\left( t \right) \right) + {∂}1F\left( t, {q}^′\left( t \right) \right) {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) {{∂}0}2\left( F \right)\left( t, {q}^′\left( t \right) \right) + {∂}1F\left( t, {q}^′\left( t \right) \right) \left( {∂}0 {∂}2 \right)\left( L \right)\left( ˆ Γ[q] \right) + {∂}2L\left( ˆ Γ[q] \right) \left( {∂}0 {∂}1 \right)\left( F \right)\left( t, {q}^′\left( t \right) \right) \end{equation}
This, again, is total madness. We really want some way to control how Scheme expands terms.
But we know what we’re looking for. Expand out the matching term of the unprimed Lagrange equations:
(->eq
((D (compose ((partial 2) L) (Gamma q)))
't))#+RESULTS[9937f7c80dcdc8744c8f2a1a57505172c55bee17]: \begin{equation} {\left( D{q}^′\left( t \right) \right)}2 {{∂}1}2\left( F \right)\left( t, {q}^′\left( t \right) \right) {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) + D{q}^′\left( t \right) {∂}1F\left( t, {q}^′\left( t \right) \right) \left( {∂}1 {∂}2 \right)\left( L \right)\left( ˆ Γ[q] \right) + 2 D{q}^′\left( t \right) \left( {∂}0 {∂}1 \right)\left( F \right)\left( t, {q}^′\left( t \right) \right) {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) + {∂}1F\left( t, {q}^′\left( t \right) \right) {D}2{q}^′\left( t \right) {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) + {∂}0F\left( t, {q}^′\left( t \right) \right) \left( {∂}1 {∂}2 \right)\left( L \right)\left( ˆ Γ[q] \right) + {{∂}2}2\left( L \right)\left( ˆ Γ[q] \right) {{∂}0}2\left( F \right)\left( t, {q}^′\left( t \right) \right) + \left( {∂}0 {∂}2 \right)\left( L \right)\left( ˆ Γ[q] \right) \end{equation}
Staring at these two equations, it becomes clear that the first contains the
second, multiplied by
Try writing out the primed Lagrange equations, and subtracting the unprimed Lagrange equations, scaled by this factor:
(let* ((primed-lagrange
(- (D (compose ((partial 2) Lprime) (Gamma qprime)))
(compose ((partial 1) Lprime) (Gamma qprime))))
(lagrange
(- (D (compose ((partial 2) L) (Gamma q)))
(compose ((partial 1) L) (Gamma q))))
(factor
(compose coordinate ((partial 1) C) (Gamma qprime))))
(->tex-equation
((- primed-lagrange (* factor lagrange))
't)))#+RESULTS[ae8a3f6f564ba9a64eff102d29316faebfb7e1f6]: \begin{equation} 0 \end{equation}
Done! It seems that the extra terms on each side exactly cancel. As with the pen and paper derivation, we’ve shown that the primed Lagranged equations \eqref{eq:lagrange-prime} hold if the unprimed equations \label{eq:1-15-lagrange} hold.
I’m troubled by my lack of intuition around two ideas:
- what is the meaning of the
$∂_1F(t, q’(t))$ scaling factor? - Both sides acquire a constant
$(∂_2L ˆ Γ[q]) ⋅ D(∂_1F(t, q’(t)))$ .
The right factor is a total time derivative. Is this meaningful? We know from later discussion that if we add a total time derivative to the Lagrangian we don’t affect the shape of the realizable path.
I learned quite a bit about functional notation from this exercise, and I think that test that the result represents is critical for leaning hard on the coordinate transformations that we’ll continue to explore. But I do feel like I’m leaving intuitive cake on the table.
:header-args+: :tangle ch1/ex1-16.scm :comments org(load "ch1/utils.scm")This exercise gives you practice in writing Lagrangians as compositions.
Find Lagrangians for central force motion in three dimensions in rectangular coordinates and in spherical coordinates. First, find the Lagrangians analytically, then check the results with the computer by generalizing the programs that we have presented.
Analytically in 3d coordinates, we have a straightforward extension:
\begin{equation} L = T - V = {1 \over 2} m (v_x^2 + v_y^2 + v_z^2) - U(\sqrt{x^2 + y^2 + z^2}) \end{equation}
How would you find the spherical Lagrangian analytically?
The potential part is easy. We know that
To transform the kinetic energies:
- write down the coordinate transformation from spherical to rectangular:
- get the velocity components by taking derivatives, chain rule
- substitute these in to the rectangular Lagrangian
This is equivalent to generating
How do you get the velocities?
Remember, from equation 1.77 in the book:
\begin{equation} C(t, x’, v’) = (t,\, F(t, x’),\, ∂_0F(t, x’) + ∂_1F(t, x’)v’) \end{equation}
\begin{equation} \begin{aligned} x & = r sin θ cos φ \cr y & = r sin θ sin φ \cr z & = r cos θ \end{aligned} \end{equation}
It’s important to remember that the partials in equation 1.77 are taken as if the arguments were static. The chain rule has already been applied. It’s more san
\begin{equation} ∂_0F(t, x’) = (0, 0, 0) \end{equation}
This has three components; one for each of the components in the return value.
The inner tuples collapse upon multiplication by
\begin{equation} \begin{aligned} v_x & = \dot{r} sin θ cos φ + r \dot{θ} cos θ cos φ + r \dot{φ} sin θ sin φ \cr v_y & = \dot{r} sin θ sin φ + r \dot{θ} cos θ sin φ + r \dot{φ} sin θ cos φ \cr v_z & = \dot{r} cos θ + r \dot{θ} sin θ \end{aligned} \end{equation}
Substituting these in results, eventually, in a waterfall of cancellations and leaves:
\begin{equation} L’ = {1 \over 2} m (r^2 \dot{φ}^2 sin^2φ + r^2 \dot{θ}^2 + \dot{r}^2) - U(r) \end{equation}
To show the rectangular Lagrangian, get the procedure from page 41:
(define ((L-central-rectangular m U) local)
(let ((q (coordinate local))
(v (velocity local)))
(- (* 1/2 m (square v))
(U (sqrt (square q))))))This is already written in a form that can handle an arbitrary number of coordinates. Confirm the rectangular Lagrangian by passing in a local tuple with 3 dimensional coordinates and velocities:
(->tex-equation
((L-central-rectangular 'm (literal-function 'U))
(up 't
(up 'x 'y 'z)
(up 'v_x 'v_y 'v_z))))#+RESULTS[9bf8352b0d3e667e7149e5519cc568efbbf2f331]: \begin{equation} {{1}\over {2}} m {{v}x}2 + {{1}\over {2}} m {{v}y}2 + {{1}\over {2}} m {{v}z}2 - U\left( \sqrt{{x}2 + {y}2 + {z}2} \right) \end{equation}
Next, the spherical. Write down the coordinate transformation from spherical to rectangular coordinates as a Scheme procedure:
(define (spherical->rect local)
(let* ((q (coordinate local))
(r (ref q 0))
(theta (ref q 1))
(phi (ref q 2)))
(up (* r (sin theta) (cos phi))
(* r (sin theta) (sin phi))
(* r (cos theta)))))Here are the velocities calculated above by hand:
(->tex-equation
(velocity
((F->C spherical->rect)
(up 't
(up 'r 'theta 'phi)
(up 'rdot 'thetadot 'phidot)))))#+RESULTS[b41cbe257b2a859243c9b36052b54e533d08484e]: \begin{equation} \begin{pmatrix} \displaystyle{ - \dot{φ} r sin\left( φ \right) sin\left( θ \right) + r \dot{θ} cos\left( φ \right) cos\left( θ \right) + \dot{r} cos\left( φ \right) sin\left( θ \right)} \cr \cr \displaystyle{ \dot{φ} r cos\left( φ \right) sin\left( θ \right) + r \dot{θ} sin\left( φ \right) cos\left( θ \right) + \dot{r} sin\left( φ \right) sin\left( θ \right)} \cr \cr \displaystyle{ - r \dot{θ} sin\left( θ \right) + \dot{r} cos\left( θ \right)}\end{pmatrix} \end{equation}
Now that we have
(define (L-central-spherical m U)
(compose (L-central-rectangular m U)
(F->C spherical->rect)))Confirm that this is equivalent to the analytic solution:
(->tex-equation
((L-central-spherical 'm (literal-function 'U))
(up 't
(up 'r 'theta 'phi)
(up 'rdot 'thetadot 'phidot))))#+RESULTS[e51d842361ff78bbab3e1c959cf534c69411b061]: \begin{equation} {{1}\over {2}} m {\dot{φ}}2 {r}2 {\left( sin\left( θ \right) \right)}2 + {{1}\over {2}} m {r}2 {\dot{θ}}2 + {{1}\over {2}} m {\dot{r}}2 - U\left( r \right) \end{equation}
Langrangian coordinate transformation from spherical -> rectangular on paper, which of course is a total nightmare, writing vx^2 + vy^2 + vz^2 and simplifying. BUT then, of course, you write down the spherical => rectangular position change…
the explicit link to function composition, and how the new lagrangian is (Lagrangian A + A<-B + B<-C)… really drives home how coordinate transforms can stack associatively through function composition. the lesson is, prove that the code works, then trust the program to go to crazy coordinate systems.
Later, the authors add in a very simple-to-write coordinate transform that has one of the angles depend on t. and then compose that in, and boom, basically for free you’re in rotating spherical coords.
:header-args+: :tangle ch1/ex1-17.scm :comments org(load "ch1/utils.scm")This, and the next three exercises, are here to give you practice in the real art, of difficulty, of any dynamics problem. It’s easy to change coordinates. So what coordinates do you use?
A bead of mass
$m$ is constrained to move on a frictionless helical wire. The helix is oriented so that its axis is horizontal. The diameter of the helix is$d$ and its pitch (turns per unit length) is$h$ . The system is in a uniform gravitational field with vertical acceleration$g$ . Formulate a Lagrangian that describes the system and find the Lagrange equations of motion.
I’ll replace this with a better picture later, but this is the setup:
(define ((turns->rect d h) local)
(let* ((turns (coordinate local))
(theta (* turns 2 'pi)))
(up (/ turns h)
(* (/ d 2) (cos theta))
(* (/ d 2) (sin theta)))))Or you could do this. Remember, these transformations need to be functions of a
local tuple, so if you’re going to compose them, remember to put coordinate at
the beginning of the composition.
(define ((turns->x-theta h) q)
(up (/ q h)
(* q 2 'pi)))
(define ((x-theta->rect d) q)
(let* ((x (ref q 0))
(theta (ref q 1)))
(up x
(* (/ d 2) (cos theta))
(* (/ d 2) (sin theta)))))
(define (turns->rect* d h)
(compose (x-theta->rect d)
(turns->x-theta h)
coordinate))The transformations are identical:
(->tex-equation
((- (turns->rect 'd 'h)
(turns->rect* 'd 'h))
(up 't 'n 'ndot)))#+RESULTS[3fa8a307ddde828e0f3b9b6c573e5935c577f76a]: \begin{equation} \begin{pmatrix} \displaystyle{ 0} \cr \cr \displaystyle{ 0} \cr \cr \displaystyle{ 0}\end{pmatrix} \end{equation}
Define the Lagrangian:
(define ((L-rectangular m U) local)
(let ((q (coordinate local))
(v (velocity local)))
(- (* 1/2 m (square v))
(U q))))
(define (L-turns m d h U)
(compose (L-rectangular m U)
(F->C (turns->rect d h))))The potential is a uniform gravitational acceleration:
(define ((U-grav m g) q)
(* m g (ref q 2)))Final Lagrangian:
(->tex-equation
((L-turns 'm 'd 'h (U-grav 'm 'g))
(up 't 'n 'ndot)))#+RESULTS[a486b9871810f8a74d1883be47666fa323a6de81]: \begin{equation} {{{{1}\over {2}} {d}2 {h}2 m {\dot{n}}2 {π}2 - {{1}\over {2}} d g {h}2 m sin\left( 2 n π \right) + {{1}\over {2}} m {\dot{n}}2}\over {{h}2}} \end{equation}
Lagrange equations of motion:
(let* ((L (L-turns 'm 'd 'h (U-grav 'm 'g)))
(n (literal-function 'n)))
(->tex-equation
(((Lagrange-equations L) n) 't)))#+RESULTS[92605835c702fcc963cc25ca8ccd8832270e11f2]: \begin{equation} {{{d}2 {h}2 m {π}2 {D}2n\left( t \right) + d g {h}2 m π cos\left( 2 π n\left( t \right) \right) + m {D}2n\left( t \right)}\over {{h}2}} \end{equation}
:header-args+: :tangle ch1/ex1-16.scm :comments org(load "ch1/utils.scm")A bead of mass
$m$ moves without friction on a triaxial ellipsoidal surface. In rectangular coordinates the surface satisfies\begin{equation} {x^2 \over a^2} + {y^2 \over b^2} + {z^2 \over c^2} = 1 \end{equation}
for some constants
$a$ ,$b$ , and$c$ . Identify suitable generalized coordinates, formulate a Lagrangian, and find Lagrange’s equations.
The transformation to elliptical coordinates is very similar to the spherical
coordinate transformation, but with a fixed
(define ((elliptical->rect a b c) local)
(let* ((q (coordinate local))
(theta (ref q 0))
(phi (ref q 1)))
(up (* a (sin theta) (cos phi))
(* b (sin theta) (sin phi))
(* c (cos theta)))))Next, the Lagrangian:
(define ((L-free-particle m) local)
(* 1/2 m (square
(velocity local))))
(define (L-central-triaxial m a b c)
(compose (L-free-particle m)
(F->C (elliptical->rect a b c))))Final Lagrangian:
(let ((local (up 't
(up 'theta 'phi)
(up 'thetadot 'phidot))))
(->tex-equation
((L-central-triaxial 'm 'a 'b 'c) local)))#+RESULTS[96029124f8a649771c860516d1c6e668422de93e]: \begin{equation} {{1}\over {2}} {a}2 m {\dot{φ}}2 {\left( sin\left( φ \right) \right)}2 {\left( sin\left( θ \right) \right)}2 - {a}2 m \dot{φ} \dot{θ} sin\left( φ \right) sin\left( θ \right) cos\left( φ \right) cos\left( θ \right) + {{1}\over {2}} {a}2 m {\dot{θ}}2 {\left( cos\left( φ \right) \right)}2 {\left( cos\left( θ \right) \right)}2 + {{1}\over {2}} {b}2 m {\dot{φ}}2 {\left( sin\left( θ \right) \right)}2 {\left( cos\left( φ \right) \right)}2 + {b}2 m \dot{φ} \dot{θ} sin\left( φ \right) sin\left( θ \right) cos\left( φ \right) cos\left( θ \right) + {{1}\over {2}} {b}2 m {\dot{θ}}2 {\left( sin\left( φ \right) \right)}2 {\left( cos\left( θ \right) \right)}2 + {{1}\over {2}} {c}2 m {\dot{θ}}2 {\left( sin\left( θ \right) \right)}2 \end{equation}
I’m sure there’s some simplification in there for us. But why?
Lagrange equations of motion:
(let* ((L (L-central-triaxial 'm 'a 'b 'c))
(theta (literal-function 'theta))
(phi (literal-function 'phi)))
(->tex-equation
(((Lagrange-equations L) (up theta phi))
't)))#+RESULTS[a4de48608ac69178397f5dbb9854dd6e4815faf5]: \begin{equation} \begin{bmatrix} \displaystyle{ - 2 {a}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) Dθ\left( t \right) Dφ\left( t \right) {\left( cos\left( θ\left( t \right) \right) \right)}2 - {a}2 m {\left( cos\left( φ\left( t \right) \right) \right)}2 {\left( Dθ\left( t \right) \right)}2 sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) - {a}2 m {\left( cos\left( φ\left( t \right) \right) \right)}2 {\left( Dφ\left( t \right) \right)}2 sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) - {b}2 m {\left( sin\left( φ\left( t \right) \right) \right)}2 {\left( Dθ\left( t \right) \right)}2 sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) - {b}2 m {\left( sin\left( φ\left( t \right) \right) \right)}2 {\left( Dφ\left( t \right) \right)}2 sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) + 2 {b}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) Dθ\left( t \right) Dφ\left( t \right) {\left( cos\left( θ\left( t \right) \right) \right)}2 - {a}2 m {D}2φ\left( t \right) sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) + {a}2 m {\left( cos\left( φ\left( t \right) \right) \right)}2 {\left( cos\left( θ\left( t \right) \right) \right)}2 {D}2θ\left( t \right) + {b}2 m {D}2φ\left( t \right) sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) + {b}2 m {\left( sin\left( φ\left( t \right) \right) \right)}2 {\left( cos\left( θ\left( t \right) \right) \right)}2 {D}2θ\left( t \right) + {c}2 m {\left( Dθ\left( t \right) \right)}2 sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) + {c}2 m {\left( sin\left( θ\left( t \right) \right) \right)}2 {D}2θ\left( t \right)} \cr \cr \displaystyle{ 2 {a}2 m {\left( sin\left( φ\left( t \right) \right) \right)}2 Dθ\left( t \right) Dφ\left( t \right) sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) + {a}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) {\left( Dθ\left( t \right) \right)}2 {\left( sin\left( θ\left( t \right) \right) \right)}2 + {a}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) {\left( Dφ\left( t \right) \right)}2 {\left( sin\left( θ\left( t \right) \right) \right)}2 - {b}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) {\left( Dθ\left( t \right) \right)}2 {\left( sin\left( θ\left( t \right) \right) \right)}2 - {b}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) {\left( Dφ\left( t \right) \right)}2 {\left( sin\left( θ\left( t \right) \right) \right)}2 + 2 {b}2 m {\left( cos\left( φ\left( t \right) \right) \right)}2 Dθ\left( t \right) Dφ\left( t \right) sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) + {a}2 m {D}2φ\left( t \right) {\left( sin\left( φ\left( t \right) \right) \right)}2 {\left( sin\left( θ\left( t \right) \right) \right)}2 - {a}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) {D}2θ\left( t \right) + {b}2 m {D}2φ\left( t \right) {\left( cos\left( φ\left( t \right) \right) \right)}2 {\left( sin\left( θ\left( t \right) \right) \right)}2 + {b}2 m sin\left( φ\left( t \right) \right) cos\left( φ\left( t \right) \right) sin\left( θ\left( t \right) \right) cos\left( θ\left( t \right) \right) {D}2θ\left( t \right)}\end{bmatrix} \end{equation}
This is fairly horrifying. This really demands animation, as I bet it looks cool, but it’s not comprehensible in this form.
:header-args+: :tangle ch1/ex1-16.scm :comments orgDouble pendulum, sort of, except the whole thing can fly around the plane.
The system description is:
(load "ch1/utils.scm")The two-bar linkage shown in figure 1.3 is constrained to move in the plane. It is composed of three small massive bodies interconnected by two massless rigid rods in a uniform gravitational field with vertical acceleration g. The rods are pinned to the central body by a hinge that allows the linkage to fold. The system is arranged so that the hinge is completely free: the members can go through all configurations without collision. Formulate a Lagrangian that describes the system and find the Lagrange equations of motion. Use the computer to do this, because the equations are rather big.
This is new. Now we have multiple bodies:
We can handle this by treating our coordinate space as having new dimensions
for, say,
Without any constraints, we have six degrees of freedom.
-
$x, y$ for the central body -
$θ$ and$φ$ for the angles off center.
(Sketch these out on the picture for the final version.)
\begin{equation} \begin{aligned} x_2(t) & = x_2(t) \cr y_2(t) & = y_2(t) \cr x_1(t) & = x_2(t) + l_1 sin θ \cr y_1(t) & = y_2(t) - l_1 cos θ \cr x_3(t) & = x_2(t) + l_2 sin φ \cr y_3(t) & = y_2(t) - l_2 cos φ \end{aligned} \end{equation}
Sketch out why this makes sense. Each angle is positive CCW for consistency, since they can swing all the way around.
Write the coordinate transformation in scheme.
(define ((double-linkage->rect l1 l2) local)
(let* ((q (coordinate local))
(theta (ref q 0))
(phi (ref q 1))
(x2 (ref q 2))
(y2 (ref q 3)))
(up (+ x2 (* l1 (sin theta)))
(- y2 (* l1 (cos theta)))
x2
y2
(+ x2 (* l2 (sin phi)))
(- y2 (* l2 (cos phi))))))Next, the Lagrangian given rectangular coordinates, assuming no constraints. Remember, we have a uniform gravitational field pointing down; this means that each of the components has a potential dragging on it.
(define ((L-double-linkage-rect m1 m2 m3 U) local)
(let* ((v (velocity local))
(vx1 (ref v 0))
(vy1 (ref v 1))
(vx2 (ref v 2))
(vy2 (ref v 3))
(vx3 (ref v 4))
(vy3 (ref v 5)))
(- (+ (* m1 (+ (square vx1)
(square vy1)))
(* m2 (+ (square vx2)
(square vy2)))
(* m3 (+ (square vx3)
(square vy3))))
(U (coordinate local)))))And the composition:
(define (L-double-linkage l1 l2 m1 m2 m3 U)
(compose (L-double-linkage-rect m1 m2 m3 U)
(F->C (double-linkage->rect l1 l2))))Gravitational potential:
(define ((U-gravity g m1 m2 m3) q)
(let* ((y1 (ref q 1))
(y2 (ref q 3))
(y3 (ref q 5)))
(* g (+ (* m1 y1)
(* m2 y2)
(* m3 y3)))))(let ((local (up 't
(up 'theta 'phi 'x_2 'y_2)
(up 'thetadot 'phidot 'xdot_2 'ydot_2)))
(U (U-gravity 'g 'm_1 'm_2 'm_3)))
(->tex-equation
((L-double-linkage 'l_1 'l_2 'm_1 'm_2 'm_3 U) local)))#+RESULTS[871defa58c4637289b4aa4236e470eceed35fc24]: \begin{equation} {{l}1}2 {m}1 {\dot{θ}}2 + 2 {l}1 {m}1 \dot{θ} {\dot{x}}2 cos\left( θ \right) + 2 {l}1 {m}1 \dot{θ} {\dot{y}}2 sin\left( θ \right) + {{l}2}2 {m}3 {\dot{φ}}2 + 2 {l}2 {m}3 \dot{φ} {\dot{x}}2 cos\left( φ \right) + 2 {l}2 {m}3 \dot{φ} {\dot{y}}2 sin\left( φ \right) + g {l}1 {m}1 cos\left( θ \right) + g {l}2 {m}3 cos\left( φ \right) - g {m}1 {y}2 - g {m}2 {y}2 - g {m}3 {y}2 + {m}1 {{\dot{x}}2}2 + {m}1 {{\dot{y}}2}2 + {m}2 {{\dot{x}}2}2 + {m}2 {{\dot{y}}2}2 + {m}3 {{\dot{x}}2}2 + {m}3 {{\dot{y}}2}2 \end{equation}
Lagrange equations of motion:
(let* ((U (U-gravity 'g 'm_1 'm_2 'm_3))
(L (L-double-linkage 'l_1 'l_2 'm_1 'm_2 'm_3 U))
(theta (literal-function 'theta))
(phi (literal-function 'phi))
(x2 (literal-function 'x_2))
(y2 (literal-function 'y_2)))
(->tex-equation
(((Lagrange-equations L) (up theta phi x2 y2))
't)))#+RESULTS[891b221c072e1a0a84d8a553f44953d7e1708fec]: \begin{equation} \begin{bmatrix} \displaystyle{ g {l}1 {m}1 sin\left( θ\left( t \right) \right) + 2 {{l}1}2 {m}1 {D}2θ\left( t \right) + 2 {l}1 {m}1 sin\left( θ\left( t \right) \right) {D}2{y}2\left( t \right) + 2 {l}1 {m}1 cos\left( θ\left( t \right) \right) {D}2{x}2\left( t \right)} \cr \cr \displaystyle{ g {l}2 {m}3 sin\left( φ\left( t \right) \right) + 2 {{l}2}2 {m}3 {D}2φ\left( t \right) + 2 {l}2 {m}3 sin\left( φ\left( t \right) \right) {D}2{y}2\left( t \right) + 2 {l}2 {m}3 cos\left( φ\left( t \right) \right) {D}2{x}2\left( t \right)} \cr \cr \displaystyle{ - 2 {l}1 {m}1 sin\left( θ\left( t \right) \right) {\left( Dθ\left( t \right) \right)}2 - 2 {l}2 {m}3 sin\left( φ\left( t \right) \right) {\left( Dφ\left( t \right) \right)}2 + 2 {l}1 {m}1 {D}2θ\left( t \right) cos\left( θ\left( t \right) \right) + 2 {l}2 {m}3 {D}2φ\left( t \right) cos\left( φ\left( t \right) \right) + 2 {m}1 {D}2{x}2\left( t \right) + 2 {m}2 {D}2{x}2\left( t \right) + 2 {m}3 {D}2{x}2\left( t \right)} \cr \cr \displaystyle{ 2 {l}1 {m}1 {\left( Dθ\left( t \right) \right)}2 cos\left( θ\left( t \right) \right) + 2 {l}2 {m}3 {\left( Dφ\left( t \right) \right)}2 cos\left( φ\left( t \right) \right) + 2 {l}1 {m}1 {D}2θ\left( t \right) sin\left( θ\left( t \right) \right) + 2 {l}2 {m}3 {D}2φ\left( t \right) sin\left( φ\left( t \right) \right) + g {m}1 + g {m}2 + g {m}3 + 2 {m}1 {D}2{y}2\left( t \right) + 2 {m}2 {D}2{y}2\left( t \right) + 2 {m}3 {D}2{y}2\left( t \right)}\end{bmatrix} \end{equation}
Kill some clear factors:
(let* ((U (U-gravity 'g 'm_1 'm_2 'm_3))
(L (L-double-linkage 'l_1 'l_2 'm_1 'm_2 'm_3 U))
(theta (literal-function 'theta))
(phi (literal-function 'phi))
(x2 (literal-function 'x_2))
(y2 (literal-function 'y_2))
(eqs (((Lagrange-equations L) (up theta phi x2 y2))
't)))
(->tex-equation
(up (/ (ref eqs 0) 'l_1 'm_1)
(/ (ref eqs 1) 'l_2 'm_3)
(/ (ref eqs 2) 2)
(ref eqs 3))))#+RESULTS[9aef049c71af2604c4b6e5d179f815d38d262792]: \begin{equation} \begin{pmatrix} \displaystyle{ g sin\left( θ\left( t \right) \right) + 2 {l}1 {D}2θ\left( t \right) + 2 sin\left( θ\left( t \right) \right) {D}2{y}2\left( t \right) + 2 cos\left( θ\left( t \right) \right) {D}2{x}2\left( t \right)} \cr \cr \displaystyle{ g sin\left( φ\left( t \right) \right) + 2 {l}2 {D}2φ\left( t \right) + 2 sin\left( φ\left( t \right) \right) {D}2{y}2\left( t \right) + 2 cos\left( φ\left( t \right) \right) {D}2{x}2\left( t \right)} \cr \cr \displaystyle{ - {l}1 {m}1 sin\left( θ\left( t \right) \right) {\left( Dθ\left( t \right) \right)}2 - {l}2 {m}3 sin\left( φ\left( t \right) \right) {\left( Dφ\left( t \right) \right)}2 + {l}1 {m}1 {D}2θ\left( t \right) cos\left( θ\left( t \right) \right) + {l}2 {m}3 {D}2φ\left( t \right) cos\left( φ\left( t \right) \right) + {m}1 {D}2{x}2\left( t \right) + {m}2 {D}2{x}2\left( t \right) + {m}3 {D}2{x}2\left( t \right)} \cr \cr \displaystyle{ 2 {l}1 {m}1 {\left( Dθ\left( t \right) \right)}2 cos\left( θ\left( t \right) \right) + 2 {l}2 {m}3 {\left( Dφ\left( t \right) \right)}2 cos\left( φ\left( t \right) \right) + 2 {l}1 {m}1 {D}2θ\left( t \right) sin\left( θ\left( t \right) \right) + 2 {l}2 {m}3 {D}2φ\left( t \right) sin\left( φ\left( t \right) \right) + g {m}1 + g {m}2 + g {m}3 + 2 {m}1 {D}2{y}2\left( t \right) + 2 {m}2 {D}2{y}2\left( t \right) + 2 {m}3 {D}2{y}2\left( t \right)}\end{pmatrix} \end{equation}
This was not as gnarly as the previous problem. Perhaps I did something wrong there. We’ll see when we get animation.
:header-args+: :tangle ch1/ex1-16.scm :comments org(load "ch1/utils.scm")Consider a pendulum of length
$l$ attached to a support that is free to move horizontally, as shown in figure 1.4. Let the mass of the support be$m1$ and the mass of the pendulum bob be$m2$ . Formulate a Lagrangian and derive Lagrange’s equations for this system.
This is interesting, and totally not-obvious how to represent with Newtonian mechanics. Here it is pretty simple. The setup:
We can use 2 coordinates:
- the horizontal position of the cart
- the angle
$θ$ of the bob.
Here’s the conversion to rectangular:
\begin{equation} \begin{aligned} x_1(t) & = x_1(t) \cr y_1(t) & = l \cr x_2(t) & = x_1(t) + l sin θ \cr y_2(t) & = l(1 - cos θ) \end{aligned} \end{equation}
Draw these on the picture to make it clearer.
Write the coordinate transformation in scheme.
(define ((sliding-pend->rect l) local)
(let* ((q (coordinate local))
(x1 (ref q 0))
(theta (ref q 1)))
(up x1
l
(+ x1 (* l (sin theta)))
(* l (- 1 (cos theta))))))Next, the Lagrangian given rectangular coordinates, assuming no constraints:
(define ((L-sliding-pend-rect m1 m2 U) local)
(let* ((v (velocity local))
(vx1 (ref v 0))
(vy1 (ref v 1))
(vx2 (ref v 2))
(vy2 (ref v 3)))
(- (+ (* m1 (+ (square vx1)
(square vy1)))
(* m2 (+ (square vx2)
(square vy2))))
(U (coordinate local)))))And the composition:
(define (L-sliding-pend l m1 m2 U)
(compose (L-sliding-pend-rect m1 m2 U)
(F->C (sliding-pend->rect l))))Gravitational potential. I could include the cart here, but since we know it’s fixed gravitationally it wouldn’t change the equations of motion.
(define ((U-gravity g m2) q)
(let* ((y2 (ref q 3)))
(* m2 g y2)))(let ((local (up 't
(up 'x_1 'theta)
(up 'xdot_1 'thetadot)))
(U (U-gravity 'g 'm_2)))
(->tex-equation
((L-sliding-pend 'l 'm_1 'm_2 U) local)))#+RESULTS[a0d77887cb7ea10a807f5718c1383453f676c148]: \begin{equation} {l}2 {m}2 {\dot{θ}}2 + 2 l {m}2 \dot{θ} {\dot{x}}1 cos\left( θ \right) + g l {m}2 cos\left( θ \right) - g l {m}2 + {m}1 {{\dot{x}}1}2 + {m}2 {{\dot{x}}1}2 \end{equation}
Lagrange equations of motion:
(let* ((U (U-gravity 'g 'm_2))
(L (L-sliding-pend 'l 'm_1 'm_2 U))
(x1 (literal-function 'x_1))
(theta (literal-function 'theta)))
(->tex-equation
(((Lagrange-equations L) (up x1 theta))
't)))#+RESULTS[a556eaea3f66fe4d94bd82a34c9f2e786a4187a1]: \begin{equation} \begin{bmatrix} \displaystyle{ - 2 l {m}2 sin\left( θ\left( t \right) \right) {\left( Dθ\left( t \right) \right)}2 + 2 l {m}2 {D}2θ\left( t \right) cos\left( θ\left( t \right) \right) + 2 {m}1 {D}2{x}1\left( t \right) + 2 {m}2 {D}2{x}1\left( t \right)} \cr \cr \displaystyle{ g l {m}2 sin\left( θ\left( t \right) \right) + 2 {l}2 {m}2 {D}2θ\left( t \right) + 2 l {m}2 {D}2{x}1\left( t \right) cos\left( θ\left( t \right) \right)}\end{bmatrix} \end{equation}
Cleaner:
(let* ((U (U-gravity 'g 'm_2))
(L (L-sliding-pend 'l 'm_1 'm_2 U))
(x1 (literal-function 'x_1))
(theta (literal-function 'theta))
(eqs (((Lagrange-equations L) (up x1 theta))
't)))
(->tex-equation
(up (ref eqs 0)
(/ (ref eqs 1) 'l 'm_2))))#+RESULTS[f3d6bb9e68a92d9b0f1c99d8aef6c9b4028d4578]: \begin{equation} \begin{pmatrix} \displaystyle{ - 2 l {m}2 sin\left( θ\left( t \right) \right) {\left( Dθ\left( t \right) \right)}2 + 2 l {m}2 {D}2θ\left( t \right) cos\left( θ\left( t \right) \right) + 2 {m}1 {D}2{x}1\left( t \right) + 2 {m}2 {D}2{x}1\left( t \right)} \cr \cr \displaystyle{ g sin\left( θ\left( t \right) \right) + 2 l {D}2θ\left( t \right) + 2 {D}2{x}1\left( t \right) cos\left( θ\left( t \right) \right)}\end{pmatrix} \end{equation}
:header-args+: :tangle ch1/ex1-21.scm :comments orgThe uneven dumbbell. We’ve just made it through four exercises which embrace the idea that you can bake constraints into the coordinate transformation. But why should we believe that this is allowed?
This exercise comes after a section called “Why it Works”.
The next exercise tries to do a coordinate change that is really careful about not changing the dimension of the configuration space, so that we can show that this move is allowed. Here’s the setup:
(load "ch1/utils.scm")The idea here is to take the distance between the particles
Goal is to assume that Newtonian mechanics’ approach to constraints, shown here,
I think, is correct, and then show that the resulting equations of motion let us
treat the distance coordinate
Takes in any number of up tuples and zips them into a new list of up-tuples by taking each element.
Many exercises have been dealing with multiple particles so far. Let’s introduce some functions that let us pluck the appropriate coordinates out of the local tuple.
If we have the velocity and mass of a particle, its kinetic energy is easy to define:
(define (KE-particle m v)
(* 1/2 m (square v)))This next function, extract-particle, takes a number of components – 2 for a
particle with 2 components, 3 for a particle in space, etc – and returns a
function of local and i, a particle index. This function can be used to
extract a sub-local-tuple for that particle from a flattened list.
(define ((extract-particle pieces) local i)
(let* ((indices (apply up (iota pieces (* i pieces))))
(extract (lambda (tuple)
(vector-map (lambda (i)
(ref tuple i))
indices))))
(up (time local)
(extract (coordinate local))
(extract (velocity local)))))Write Newton’s equations for the balance of forces for the four rectangular coordinates of the two particles, given that the scalar tension in the rod is
$F$ .
TODO: Write these down from the notebook.
Write the formal Lagrangian
\begin{equation} L(t; x_0, y_0, x_1, y_1, F; \dot{x}_0, \dot{y}_0, \dot{x}_1, \dot{y}_1, \dot{F}) \end{equation}
such that Lagrange’s equations will yield the Newton’s equations you derived in part a.
Here is how we model constraint forces. Each pair of particles has some constraint potential acting between them:
(define (U-constraint q0 q1 F l)
(* (/ F (* 2 l))
(- (square (- q1 q0))
(square l))))And here’s a Lagrangian for two free particles, subject to a constraint
potential
(define ((L-free-constrained m0 m1 l) local)
(let* ((extract (extract-particle 2))
(p0 (extract local 0))
(q0 (coordinate p0))
(qdot0 (velocity p0))
(p1 (extract local 1))
(q1 (coordinate p1))
(qdot1 (velocity p1))
(F (ref (coordinate local) 4)))
(- (+ (KE-particle m0 qdot0)
(KE-particle m1 qdot1))
(U-constraint q0 q1 F l))))Finally, the path. This is rectangular coordinates for each piece, plus
(define q-rect
(up (literal-function 'x_0)
(literal-function 'y_0)
(literal-function 'x_1)
(literal-function 'y_1)
(literal-function 'F)))This shows the lagrangian itself, which answers part b:
(let* ((L (L-free-constrained 'm_0 'm_1 'l))
(f (compose L (Gamma q-rect))))
(->tex-equation
(f 't)))#+RESULTS[9e535179734061e62fa19be6d57ad8f8846008d9]: \begin{equation} {{{{1}\over {2}} l {m}0 {\left( D{x}0\left( t \right) \right)}2 + {{1}\over {2}} l {m}0 {\left( D{y}0\left( t \right) \right)}2 + {{1}\over {2}} l {m}1 {\left( D{x}1\left( t \right) \right)}2 + {{1}\over {2}} l {m}1 {\left( D{y}1\left( t \right) \right)}2 + {{1}\over {2}} {l}2 F\left( t \right) - {{1}\over {2}} F\left( t \right) {\left( {x}1\left( t \right) \right)}2 + F\left( t \right) {x}1\left( t \right) {x}0\left( t \right) - {{1}\over {2}} F\left( t \right) {\left( {x}0\left( t \right) \right)}2 - {{1}\over {2}} F\left( t \right) {\left( {y}1\left( t \right) \right)}2 + F\left( t \right) {y}1\left( t \right) {y}0\left( t \right) - {{1}\over {2}} F\left( t \right) {\left( {y}0\left( t \right) \right)}2}\over {l}} \end{equation}
Here are the Lagrange equations, which, if you squint, are like Newton’s equations from part a.
(let* ((L (L-free-constrained 'm_0 'm_1 'l))
(f ((Lagrange-equations L) q-rect)))
(->tex-equation
(f 't)))#+RESULTS[2d4f6bae712013ae3a23025d941fa644e6cbf412]: \begin{equation} \begin{bmatrix} \displaystyle{ {{l {m}0 {D}2{x}0\left( t \right) - F\left( t \right) {x}1\left( t \right) + F\left( t \right) {x}0\left( t \right)}\over {l}}} \cr \cr \displaystyle{ {{l {m}0 {D}2{y}0\left( t \right) - F\left( t \right) {y}1\left( t \right) + F\left( t \right) {y}0\left( t \right)}\over {l}}} \cr \cr \displaystyle{ {{l {m}1 {D}2{x}1\left( t \right) + F\left( t \right) {x}1\left( t \right) - F\left( t \right) {x}0\left( t \right)}\over {l}}} \cr \cr \displaystyle{ {{l {m}1 {D}2{y}1\left( t \right) + F\left( t \right) {y}1\left( t \right) - F\left( t \right) {y}0\left( t \right)}\over {l}}} \cr \cr \displaystyle{ {{ - {{1}\over {2}} {l}2 + {{1}\over {2}} {\left( {x}1\left( t \right) \right)}2 - {x}1\left( t \right) {x}0\left( t \right) + {{1}\over {2}} {\left( {x}0\left( t \right) \right)}2 + {{1}\over {2}} {\left( {y}1\left( t \right) \right)}2 - {y}1\left( t \right) {y}0\left( t \right) + {{1}\over {2}} {\left( {y}0\left( t \right) \right)}2}\over {l}}}\end{bmatrix} \end{equation}
Make a change of coordinates to a coordinate system with center of mass coordinates $xCM$, $yCM
$, angle $ θ$, distance between the particles$c$ , and tension force$F$ . Write the Lagrangian in these coordinates, and write the Lagrange equations.
This is a coordinate change that is very careful not to reduce the degrees of freedom.
First, the coordinate change:
(define ((cm-theta->rect m0 m1) local)
(let* ((q (coordinate local))
(x_cm (ref q 0))
(y_cm (ref q 1))
(theta (ref q 2))
(c (ref q 3))
(F (ref q 4))
(total-mass (+ m0 m1))
(m0-distance (* c (/ m1 total-mass)))
(m1-distance (* c (/ m0 total-mass))))
(up (- x_cm (* m0-distance (cos theta)))
(- y_cm (* m0-distance (sin theta)))
(+ x_cm (* m1-distance (cos theta)))
(+ y_cm (* m1-distance (sin theta)))
F)))Then the coordinate change applied to the local tuple:
(let ((local (up 't
(up 'x_cm 'y_cm 'theta 'c 'F)
(up 'xdot_cm 'ydot_cm 'thetadot 'cdot 'Fdot)))
(C (F->C (cm-theta->rect 'm_0 'm_1))))
(->tex-equation
(C local)))#+RESULTS[2c67e287b90be4d75ab8c396222d968e56d71383]: \begin{equation} \begin{pmatrix} \displaystyle{ t} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ {{ - c {m}1 cos\left( θ \right) + {m}0 {x}cm + {m}1 {x}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{ - c {m}1 sin\left( θ \right) + {m}0 {y}cm + {m}1 {y}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{c {m}0 cos\left( θ \right) + {m}0 {x}cm + {m}1 {x}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{c {m}0 sin\left( θ \right) + {m}0 {y}cm + {m}1 {y}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ F}\end{pmatrix}} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ {{c {m}1 \dot{θ} sin\left( θ \right) - \dot{c} {m}1 cos\left( θ \right) + {m}0 {\dot{x}}cm + {m}1 {\dot{x}}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{ - c {m}1 \dot{θ} cos\left( θ \right) - \dot{c} {m}1 sin\left( θ \right) + {m}0 {\dot{y}}cm + {m}1 {\dot{y}}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{ - c {m}0 \dot{θ} sin\left( θ \right) + \dot{c} {m}0 cos\left( θ \right) + {m}0 {\dot{x}}cm + {m}1 {\dot{x}}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{c {m}0 \dot{θ} cos\left( θ \right) + \dot{c} {m}0 sin\left( θ \right) + {m}0 {\dot{y}}cm + {m}1 {\dot{y}}cm}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ \dot{F}}\end{pmatrix}}\end{pmatrix} \end{equation}
Then the Lagrangian in the new coordinates;
(define (L-free-constrained* m0 m1 l)
(compose (L-free-constrained m0 m1 l)
(F->C (cm-theta->rect m0 m1))))This shows the lagrangian itself, after the coordinate transformation:
(let* ((q (up (literal-function 'x_cm)
(literal-function 'y_cm)
(literal-function 'theta)
(literal-function 'c)
(literal-function 'F)))
(L (L-free-constrained* 'm_0 'm_1 'l))
(f (compose L (Gamma q))))
(->tex-equation
(f 't)))#+RESULTS[d36e7f3c5c5e5a8663d599a055face7a7f6ddb61]: \begin{equation} {{l {m}0 {m}1 {\left( c\left( t \right) \right)}2 {\left( Dθ\left( t \right) \right)}2 + l {{m}0}2 {\left( D{x}cm\left( t \right) \right)}2 + l {{m}0}2 {\left( D{y}cm\left( t \right) \right)}2 + 2 l {m}0 {m}1 {\left( D{x}cm\left( t \right) \right)}2 + 2 l {m}0 {m}1 {\left( D{y}cm\left( t \right) \right)}2 + l {m}0 {m}1 {\left( Dc\left( t \right) \right)}2 + l {{m}1}2 {\left( D{x}cm\left( t \right) \right)}2 + l {{m}1}2 {\left( D{y}cm\left( t \right) \right)}2 + {l}2 {m}0 F\left( t \right) + {l}2 {m}1 F\left( t \right) - {m}0 F\left( t \right) {\left( c\left( t \right) \right)}2 - {m}1 F\left( t \right) {\left( c\left( t \right) \right)}2}\over {2 l {m}0 + 2 l {m}1}} \end{equation}
Here are the Lagrange equations:
(let* ((q (up (literal-function 'x_cm)
(literal-function 'y_cm)
(literal-function 'theta)
(literal-function 'c)
(literal-function 'F)))
(L (L-free-constrained* 'm_0 'm_1 'l))
(f ((Lagrange-equations L) q)))
(->tex-equation
(f 't)))#+RESULTS[c6a3572c8b6115d086c9f8ee85939b79942967ed]: \begin{equation} \begin{bmatrix} \displaystyle{ {m}0 {D}2{x}cm\left( t \right) + {m}1 {D}2{x}cm\left( t \right)} \cr \cr \displaystyle{ {m}0 {D}2{y}cm\left( t \right) + {m}1 {D}2{y}cm\left( t \right)} \cr \cr \displaystyle{ {{2 {m}0 {m}1 Dc\left( t \right) c\left( t \right) Dθ\left( t \right) + {m}0 {m}1 {\left( c\left( t \right) \right)}2 {D}2θ\left( t \right)}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{ - l {m}0 {m}1 c\left( t \right) {\left( Dθ\left( t \right) \right)}2 + l {m}0 {m}1 {D}2c\left( t \right) + {m}0 F\left( t \right) c\left( t \right) + {m}1 F\left( t \right) c\left( t \right)}\over {l {m}0 + l {m}1}}} \cr \cr \displaystyle{ {{ - {{1}\over {2}} {l}2 + {{1}\over {2}} {\left( c\left( t \right) \right)}2}\over {l}}}\end{bmatrix} \end{equation}
That final equation states that
You may deduce from one of these equations that
$c(t) = l$ . From this fact we get that$Dc = 0$ and$D^2c = 0$ . Substitute these into the Lagrange equations you just computed to get the equation of motion for $xCM$, $yCM$, $ θ$.
We can substitute the constant value of
(let* ((q (up (literal-function 'x_cm)
(literal-function 'y_cm)
(literal-function 'theta)
(lambda (t) 'l)
(literal-function 'F)))
(L (L-free-constrained* 'm_0 'm_1 'l))
(f ((Lagrange-equations L) q)))
(->tex-equation
(f 't)))#+RESULTS[a51bdf5bbe673771262278b5dca048646c257752]: \begin{equation} \begin{bmatrix} \displaystyle{ {m}0 {D}2{x}cm\left( t \right) + {m}1 {D}2{x}cm\left( t \right)} \cr \cr \displaystyle{ {m}0 {D}2{y}cm\left( t \right) + {m}1 {D}2{y}cm\left( t \right)} \cr \cr \displaystyle{ {{{l}2 {m}0 {m}1 {D}2θ\left( t \right)}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ {{ - l {m}0 {m}1 {\left( Dθ\left( t \right) \right)}2 + {m}0 F\left( t \right) + {m}1 F\left( t \right)}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ 0}\end{bmatrix} \end{equation}
This is saying that the acceleration on the center of mass is 0.
The fourth equation, the equation of motion for the
(define (reduced-mass m1 m2)
(/ (* m1 m2)
(+ m1 m2)))If we let
\begin{equation} \label{eq:constraint-force} F(t) = m l \dot{θ}^2 \end{equation}
We can verify this with Scheme by subtracting the two equations:
(let* ((F (literal-function 'F))
(theta (literal-function 'theta))
(q (up (literal-function 'x_cm)
(literal-function 'y_cm)
theta
(lambda (t) 'l)
F))
(L (L-free-constrained* 'm_0 'm_1 'l))
(f ((Lagrange-equations L) q))
(m (reduced-mass 'm_0 'm_1)))
(->tex-equation
(- (ref (f 't) 3)
(- (F 't)
(* m 'l (square ((D theta) 't)))))))#+RESULTS[98183c269755d12e4a77e511208ba05dd8872966]: \begin{equation} 0 \end{equation}
Make a Lagrangian (
$= T − V$ ) for the system described with the irredundant generalized coordinates $xCM$, $yCM$, $ θ$ and compute the Lagrange equations from this Lagrangian. They should be the same equations as you derived for the same coordinates in part d.
For part e, I wrote this in the notebook - it is effectively identical to the substitution that is happening on the computer, so I’m going to ignore this. You just get more cancellations.
But let’s go at it, for fun.
Here’s the Lagrangian of 2 free particles:
(define ((L-free2 m0 m1) local)
(let* ((extract (extract-particle 2))
(p0 (extract local 0))
(q0 (coordinate p0))
(qdot0 (velocity p0))
(p1 (extract local 1))
(q1 (coordinate p1))
(qdot1 (velocity p1)))
(+ (KE-particle m0 qdot0)
(KE-particle m1 qdot1))))Then a version of cm-theta->rect where we ignore
(define ((cm-theta->rect* m0 m1 l) local)
(let* ((q (coordinate local))
(x_cm (ref q 0))
(y_cm (ref q 1))
(theta (ref q 2))
(total-mass (+ m0 m1))
(m0-distance (* l (/ m1 total-mass)))
(m1-distance (* l (/ m0 total-mass))))
(up (- x_cm (* m0-distance (cos theta)))
(- y_cm (* m0-distance (sin theta)))
(+ x_cm (* m1-distance (cos theta)))
(+ y_cm (* m1-distance (sin theta))))))#+RESULTS[aee3cbebd833662ba2610e81cef9e42936dcb703]: #| cm-theta->rect* |#
The Lagrangian:
(define (L-free-constrained2 m0 m1 l)
(compose (L-free2 m0 m1)
(F->C (cm-theta->rect* m0 m1 l))))Equations:
(let* ((q (up (literal-function 'x_cm)
(literal-function 'y_cm)
(literal-function 'theta)
(lambda (t) 'l)
(literal-function 'F)))
(L1 (L-free-constrained* 'm_0 'm_1 'l))
(L2 (L-free-constrained2 'm_0 'm_1 'l)))
(->tex-equation
((- ((Lagrange-equations L1) q)
((Lagrange-equations L2) q))
't)))#+RESULTS[132abd940897cd8f9d1e34d11fb3c415fce45a61]: \begin{equation} \begin{bmatrix} \displaystyle{ 0} \cr \cr \displaystyle{ 0} \cr \cr \displaystyle{ 0} \cr \cr \displaystyle{ {{ - l {m}0 {m}1 {\left( Dθ\left( t \right) \right)}2 + {m}0 F\left( t \right) + {m}1 F\left( t \right)}\over {{m}0 + {m}1}}} \cr \cr \displaystyle{ 0}\end{bmatrix} \end{equation}
The only remaining equation is \eqref{eq:constraint-force} from above. This
remains because the simplified Lagrangian ignores the
(load "ch1/utils.scm")Show that the Lagrangian (1.93) can be used to describe the driven pendulum (section 1.6.2), where the position of the pivot is a specified function of time: Derive the equations of motion using the Newtonian constraint force prescription, and show that they are the same as the Lagrange equations. Be sure to examine the equations for the constraint forces as well as the position of the pendulum bob.
I might be slightly misunderstanding this question. The problem states that we should use equations 1.91 and 1.92 in the book to express the equations of motion of the driven pendulum, and then to rederive the equations of motion using the Lagrangian.
The final note about examining the constraint forces means that we’ll need to
follow the approach of equation 1.21 and include a coordinate transformation
from some
Step one is to use 1.91 and 1.92 in the book to express
\begin{equation} V(t, θ, \dot{θ}) = m g (y_s(t) - l cos θ) \end{equation}
So equation 1.91 becomes, for the single pendulum bob:
\begin{equation} \begin{aligned} m D^2x(t) & = F(t) {{-x(t)} \over l} \cr m D^2y(t) & = -mgl + F(t) {{y_s(t) - y(t)} \over l} \end{aligned} \end{equation}
The assumption here is that the pendulum support sits at
Now write the Lagrangian for the driven pendulum in rectangular coordinates. The constraint force takes the same shape as in exercise 1.21:
<<U-constraint>>The Lagrangian is similar, but only involves a single particle – the pendulum bob. We can generate the constraint force by directly building the support’s coordinates, rather than extracting them from the local tuple.
(define ((L-driven-free m l y U) local)
(let* ((extract (extract-particle 2))
(bob (extract local 0))
(q (coordinate bob))
(qdot (velocity bob))
(F (ref (coordinate local) 2)))
(- (KE-particle m qdot)
(U q)
(U-constraint (up 0 (y (time local)))
q
F
l))))Here is the now-familiar equation for a uniform gravitational potential, acting
on the
(define ((U-gravity g m) q)
(let* ((y (ref q 1)))
(* m g y)))Now use the new Lagrangian to generate equations of motion for the three
coordinates
(let* ((q (up (literal-function 'x)
(literal-function 'y)
(literal-function 'F)))
(U (U-gravity 'g 'm))
(y (literal-function 'y_s))
(L (L-driven-free 'm 'l y U))
(f ((Lagrange-equations L) q)))
(->tex-equation
(f 't)))#+RESULTS[a7c6fc40bbcba84c2c0ef9058455e719c4ac8045]: \begin{equation} \begin{bmatrix} \displaystyle{ {{l m {D}2x\left( t \right) + F\left( t \right) x\left( t \right)}\over {l}}} \cr \cr \displaystyle{ {{g l m + l m {D}2y\left( t \right) - F\left( t \right) {y}s\left( t \right) + F\left( t \right) y\left( t \right)}\over {l}}} \cr \cr \displaystyle{ {{ - {{1}\over {2}} {l}2 + {{1}\over {2}} {\left( {y}s\left( t \right) \right)}2 - {y}s\left( t \right) y\left( t \right) + {{1}\over {2}} {\left( y\left( t \right) \right)}2 + {{1}\over {2}} {\left( x\left( t \right) \right)}2}\over {l}}}\end{bmatrix} \end{equation}
The first two equations of motion match the equations we derived in part A, using Newton’s equations. The third states that
\begin{equation} l^2 = x(t)^2 + (y_s(t) - y(t))^&2 \end{equation}
Verified, with some extra terms to force the simplification:
(let* ((q (up (literal-function 'x)
(literal-function 'y)
(literal-function 'F)))
(U (U-gravity 'g 'm))
(y (literal-function 'y_s))
(L (L-driven-free 'm 'l y U))
(f ((Lagrange-equations L) q))
(eq (ref (f 't) 2)))
(->tex-equation
(- eq
(/ (* 1/2 (- (+ (square ((literal-function 'x) 't))
(square ((- y (literal-function 'y)) 't)))
(square 'l)))
'l))))#+RESULTS[2bf6d50380d98e5aff9861c349b4e289b31fa168]: \begin{equation} 0 \end{equation}
Now we want to verify that we get the same Lagrangian and equations of motion as in 1.88 from the book. We also want to analyze the constraint forces. To do this we need to introduce a coordinate change.
To analyze the constraint forces, we have to do the same trick as in exercise
1.21 and use a coordinate
(define ((driven-polar->rect y) local)
(let* ((q (coordinate local))
(theta (ref q 0))
(c (ref q 1))
(F (ref q 2)))
(up (* c (sin theta))
(- (y (time local)) (* c (cos theta)))
F)))Compose the coordinate change with the rectangular Lagrangian:
(define (L-driven-pend m l y U)
(compose (L-driven-free m l y U)
(F->C (driven-polar->rect y))))Examine the Lagrangian itself, after the coordinate transformation. (Notice that
we’re using a constant function for
(let* ((q (up (literal-function 'theta)
(lambda (t) 'l)
(literal-function 'F)))
(y (literal-function 'y_s))
(L (L-driven-pend 'm 'l y (U-gravity 'g 'm)))
(f (compose L (Gamma q))))
(->tex-equation
(f 't)))#+RESULTS[306c957e09bfd22ebe85e5c6fb9380ff977df1e9]: \begin{equation} {{1}\over {2}} {l}2 m {\left( Dθ\left( t \right) \right)}2 + l m D{y}s\left( t \right) sin\left( θ\left( t \right) \right) Dθ\left( t \right) + g l m cos\left( θ\left( t \right) \right) - g m {y}s\left( t \right) + {{1}\over {2}} m {\left( D{y}s\left( t \right) \right)}2 \end{equation}
Looks just like equation 1.88.
Next, examine the Lagrange equations, using the same substitution of
(let* ((q (up (literal-function 'theta)
(lambda (t) 'l)
(literal-function 'F)))
(y (literal-function 'y_s))
(L (L-driven-pend 'm 'l y (U-gravity 'g 'm)))
(f ((Lagrange-equations L) q)))
(->tex-equation
(f 't)))#+RESULTS[f9c380046af367934c5d46dd0d3793a1082a018e]: \begin{equation} \begin{bmatrix} \displaystyle{ g l m sin\left( θ\left( t \right) \right) + {l}2 m {D}2θ\left( t \right) + l m {D}2{y}s\left( t \right) sin\left( θ\left( t \right) \right)} \cr \cr \displaystyle{ - l m {\left( Dθ\left( t \right) \right)}2 - g m cos\left( θ\left( t \right) \right) - m {D}2{y}s\left( t \right) cos\left( θ\left( t \right) \right) + F\left( t \right)} \cr \cr \displaystyle{ 0}\end{bmatrix} \end{equation}
The third equation is 0 because of the substitution of constant
The second equation describes the constraint force on the driven pendulum as a function of the other coordinates and the support position.
:header-args+: :tangle ch1/ex1-23.scm :comments org(load "ch1/utils.scm")TODO: Expand out the explicit Lagrangian, using a coordinate transformation, and do the manual substitution…
:header-args+: :tangle ch1/ex1-24.scm :comments orgThis is a special case of a solution we found in exercise 1.22. In that
exercise, we found the constraint forces on a driven pendulum. By setting
(load "ch1/utils.scm")Take some definitions that we need:
<<L-driven-free>>
<<U-gravity>>
<<driven-polar->rect>>
<<L-driven-pend>>The second equation of motion, for the (lambda (t) 'l):
(let* ((q (up (literal-function 'theta)
(lambda (t) 'l)
(literal-function 'F)))
(y (lambda (t) 'l))
(L (L-driven-pend 'm 'l y (U-gravity 'g 'm)))
(f ((Lagrange-equations L) q)))
(->tex-equation
(ref (f 't) 1)))#+RESULTS[ef32abeb390f21f1abca982d55c01e5ea78b4d80]: \begin{equation}
- l m {\left( Dθ\left( t \right) \right)}2 - g m cos\left( θ\left( t \right) \right) + F\left( t \right)
\end{equation}
Solve for
\begin{equation} F(t) = m (g cos θ + l \dot{θ}^2) \end{equation}
:header-args+: :tangle ch1/ex1-25.scm :comments org(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")nice, easy to guess.
(define ((FA m) local)
(let ((x (coordinate local)))
(* m x)))Show the function of t, and confirm that both methods are equivalent.
(check-f (FA 'm)
(literal-function 'x))NOT a total time derivative.
Define G directly:
(define ((GB m) local)
(let* ((t (time local))
(v_x (velocity local))
(GB0 0)
(GB1 (* m (cos t))))
(+ GB0 (* GB1 v_x))))And show the full G, for fun:
(let ((f (compose (GB 'm) (Gamma (literal-function 'x)))))
(se (f 't)))It’s easier to confirm that this is not a total time derivative by checking the partials.
(define (GB-properties m)
(let ((GB0 (lambda (local) 0))
(GB1 (lambda (local)
(* m (cos (time local))))))
(G-properties GB0 GB1 (literal-function 'x))))It’s clear here that the second and third tuple entries aren’t equal, so we don’t have a total time derivative.
(se (GB-properties 'm))no problem, we’ve got a total time derivative on our hands.
(define (FC local)
(let ((t (time local))
(x (coordinate local)))
(* x (cos t))))
(check-f FC (literal-function 'x))
(define GC-properties
(let ((GC0 (lambda (local)
(* -1
(coordinate local)
(sin (time local)))))
(GC1 (lambda (local)
(cos (time local)))))
(G-properties GC0 GC1 (literal-function 'x))))Boom, the second and third entries are equal, as we’d expect.
(se GC-properties)This is NOT a total time derivative; you can tell by taking the partials of each side, G0 and G1, as we’ll see here.
(define GD-properties
(let ((GD0 (lambda (local)
(* (coordinate local)
(sin (time local)))))
(GD1 (lambda (local)
(cos (time local)))))
(G-properties GD0 GD1 (literal-function 'x))))The partials for each side don’t match.
(se GD-properties)This is strange to me, because I thought that this thing had to produce a tuple.
OH, but the secret is that Qdot is also a tuple, so you contract them together.
Here’s the function F that we can use to derive it:
(define (FE local)
(let* ((t (time local))
(q (coordinate local))
(x (ref q 0))
(y (ref q 1)))
(* (+ (square x) (square y))
(cos t))))Boom, total time derivative!
(check-f FE (up (literal-function 'x)
(literal-function 'y)))And let’s show that we pass the tests by decomposing this into G0 and G1:
(define GE-properties
(let (
;; any piece of the function without a velocity multiplied.
(GE0 (lambda (local)
(let* ((t (time local))
(q (coordinate local))
(x (ref q 0))
(y (ref q 1)))
(* -1
(+ (square x) (square y))
(sin t)))))
;; The pieces multiplied by velocities, split into a down tuple of
;; components, one for each of the coordinate components.
(GE1 (lambda (local)
(let* ((t (time local))
(q (coordinate local))
(x (ref q 0))
(y (ref q 1)))
(down
(* 2 x (cos t))
(* 2 y (cos t)))))))
(G-properties GE0 GE1 (up (literal-function 'x)
(literal-function 'y)))))BOOM!
We’ve recovered F; the partials are equal, and the final matrix is symmetric.
(se GE-properties)This one is interesting, since the second partial is a tuple. This is not so obvious to me, so first let’s check the properties:
(define GF-properties
(let (
;; any piece of the function without a velocity multiplied.
(GF0 (lambda (local)
(let* ((t (time local))
(q (coordinate local))
(x (ref q 0))
(y (ref q 1)))
(* -1
(+ (square x) (square y))
(sin t)))))
;; The pieces multiplied by velocities, split into a down tuple of
;; components, one for each of the coordinate components.
(GF1 (lambda (local)
(let* ((t (time local))
(q (coordinate local))
(x (ref q 0))
(y (ref q 1)))
(down
(+ (cube y) (* 2 x (cos t)))
(+ x (* 2 y (cos t))))))))
(G-properties GF0 GF1 (up (literal-function 'x)
(literal-function 'y)))))AND it looks like we DO have a total time derivative, maybe. We certainly pass the first test here, since the second and third tuple entries are equal.
BUT we fail the second test; the hessian that we get from ((partial 1) G1) is not symmetric.
(se GF-properties)I’ll do this for a single particle, since it’s annoying to get the sum going for many; and the lagrangian is additive, so no problem.
(define (uniform-translate-shift->rect local)
(let* ((t (time local))
(q (coordinate local))
(xprime (ref q 0))
(delta_x (ref q 1))
(delta_v (ref q 2)))
(+ xprime delta_x (* t delta_v))))
(define (L-translate-shift m)
(compose (L-free-particle m)
(F->C uniform-translate-shift->rect)))First, confirm that if we have a constant, we get what we expected from paper.
(let* ((q (up (literal-function 'xprime)
(lambda (t) 'Delta_x)
(lambda (t) 'Delta_v)))
(f (compose (L-translate-shift 'm) (Gamma q))))
(->tex-equation (f 't)))#+RESULTS[5d2b4de08cfab4779bf7cdab31d518191b40a4d2]: \[\begin{equation} {{1}\over {2}} {{Δ}v}2 m + {Δ}v m D{x}^′\left( t \right) + {{1}\over {2}} m {\left( D{x}^′\left( t \right) \right)}2 \end{equation}\]
We can change this a little to see the extra terms; substract off the free particle lagrangian, to see the extra stuff.
(let* ((q (up (literal-function 'xprime)
(lambda (t) 'Delta_x)
(lambda (t) 'Delta_v)))
(L (- (L-translate-shift 'm)
(L-free-particle 'm)))
(f (compose L (Gamma q))))
(->tex-equation (f 't)))#+RESULTS[c17004e61fec7edb3835203cdc99c562940bee7c]: \[\begin{equation} {{1}\over {2}} {{Δ}v}2 m + {Δ}v m D{x}^′\left( t \right) \end{equation}\]
Here’s the gnarly version with both entries as actual functions. Can this be a
total time derivative? It CANNOT be, because we have a
(let* ((q (up (literal-function 'xprime)
(literal-function 'Delta_x)
(literal-function 'Delta_v)))
(L (- (L-translate-shift 'm)
(L-free-particle 'm)))
(f (compose L (Gamma q))))
(->tex-equation (f 't)))#+RESULTS[ded4f6dec25954c9b7536153e1db8db0315cb399]: \[ \begin{equation} {{1}\over {2}} m {t}2 {\left( D{Δ}v\left( t \right) \right)}2 + m t D{x}^′\left( t \right) D{Δ}v\left( t \right) + m t D{Δ}v\left( t \right) {Δ}v\left( t \right) + m t D{Δ}v\left( t \right) D{Δ}x\left( t \right) + m D{x}^′\left( t \right) {Δ}v\left( t \right) + m D{x}^′\left( t \right) D{Δ}x\left( t \right) - {{1}\over {2}} m {\left( D{Δ}v\left( t \right) \right)}2 + {{1}\over {2}} m {\left( {Δ}v\left( t \right) \right)}2 + m {Δ}v\left( t \right) D{Δ}x\left( t \right) \end{equation} \]
Let’s simplify by making the
We know that we have a total derivative when
(let* ((q (lambda (dx)
(up (literal-function 'xprime)
dx
(lambda (t) 'Delta_v))))
(L (- (L-translate-shift 'm)
(L-free-particle 'm)))
(f (lambda (dx)
(compose L (Gamma (q dx))))))
(->tex-equation
((- (f (literal-function 'Delta_x))
(f (lambda (t) 'Delta_x)))
't)))#+RESULTS[1a9463beb2f26c1661f1978633ca830ba12f73ec]: \[\begin{equation} {Δ}v m D{Δ}x\left( t \right) + m D{x}^′\left( t \right) D{Δ}x\left( t \right) \end{equation}\]
Take a look. there is a quadratic velocity term in here! We have
SO, only if the shift and uniform translation are constant do we not affect the Lagrangian value.
:header-args+: :tangle ch1/ex1-30.scm :comments org(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")(load "ch1/utils.scm")Notation Appendix. This is all about getting cozy with scheme, and with the various idiosyncracies of the tuple and functional notation.
:header-args+: :tangle ch9/ex9-1.scm :comments orgYou’re supposed to do these by hand, so I’ll do that in the textbook. But here, let’s redo them on the machine.
First, let’s define the functions we need.
(define (F x y)
(* (square x)
(cube y)))
(define (G x y)
(up (F x y) y))
(define (H x y)
(F (F x y) y))You can do this with explicit partials:
(let ((f (down ((partial 0) F) ((partial 1) F))))
(->tex-equation
(f 'x 'y)))#+RESULTS[b8eaf52d98e5903b52306509dcdc8f8eeb97144c]: \begin{equation} \begin{bmatrix} \displaystyle{ 2 x {y}3} \cr \cr \displaystyle{ 3 {x}2 {y}2}\end{bmatrix} \end{equation}
Or with the
(->tex-equation
((D F) 'x 'y))#+RESULTS[f3fba605ac97a3ebd30b4a96aca31eed921e2e93]: \begin{equation} \begin{bmatrix} \displaystyle{ 2 x {y}3} \cr \cr \displaystyle{ 3 {x}2 {y}2}\end{bmatrix} \end{equation}
Or, we could show that they’re equivalent this way:
(let ((f (down ((partial 0) F) ((partial 1) F))))
(->tex-equation
(- ((D F) 'x 'y)
(f 'x 'y))))#+RESULTS[bbfc31a98ddca1b434403a34cefb730e354f1be8]: \begin{equation} \begin{bmatrix} \displaystyle{ 0} \cr \cr \displaystyle{ 0}\end{bmatrix} \end{equation}
(->tex-equation
((D H) 'x 'y))#+RESULTS[22a0dfcbcf713d36b0f899b6baac6dbf1ec4b56d]: \begin{equation} \begin{bmatrix} \displaystyle{ 4 {x}3 {y}9} \cr \cr \displaystyle{ 9 {x}4 {y}8}\end{bmatrix} \end{equation}
(->tex-equation
((D G) 'x 'y))#+RESULTS[548f447f81ffe817f686965fb5fdc1d0cbecc5f9]: \begin{equation} \begin{bmatrix} \displaystyle{ \begin{pmatrix} \displaystyle{ 2 x {y}3} \cr \cr \displaystyle{ 0}\end{pmatrix}} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ 3 {x}2 {y}2} \cr \cr \displaystyle{ 1}\end{pmatrix}}\end{bmatrix} \end{equation}
(->tex-equation
(up ((D F) 'a 'b)
((D G) 3 5)
((D H) (* 3 (square 'a)) (* 5 (cube 'b)))))#+RESULTS[e0ef4bfc15551f9d05baeb3970cd8dafaf02db65]: \begin{equation} \begin{pmatrix} \displaystyle{ \begin{bmatrix} \displaystyle{ 2 a {b}3} \cr \cr \displaystyle{ 3 {a}2 {b}2}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ \begin{pmatrix} \displaystyle{ 750} \cr \cr \displaystyle{ 0}\end{pmatrix}} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ 675} \cr \cr \displaystyle{ 1}\end{pmatrix}}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ 210937500 {a}6 {b}27} \cr \cr \displaystyle{ 284765625 {a}8 {b}24}\end{bmatrix}}\end{pmatrix} \end{equation}
:header-args+: :tangle ch9/ex9-2.scm :comments orgA further exercise is to try defining the functions so that they use explicit tuples, so you can compose them:
(define (F* v)
(let ((x (ref v 0))
(y (ref v 1)))
(* (square x) (cube y))))
(define (G* v)
(let ((x (ref v 0))
(y (ref v 1)))
(up (F* v) y)))
(define H* (compose F* G*))to be really pro, I’d make a function that takes these as arguments and prints a nice formatted exercise output. Let’s do the final exercise, for fun:
(->tex-equation
(up ((D F*) (up 'a 'b))
((D G*) (up 3 5))
((D H*) (up (* 3 (square 'a)) (* 5 (cube 'b))))))#+RESULTS[1e43354828c8ce0ba497bcc6bd9e64c4f4e20419]: \begin{equation} \begin{pmatrix} \displaystyle{ \begin{bmatrix} \displaystyle{ 2 a {b}3} \cr \cr \displaystyle{ 3 {a}2 {b}2}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ \begin{pmatrix} \displaystyle{ 750} \cr \cr \displaystyle{ 0}\end{pmatrix}} \cr \cr \displaystyle{ \begin{pmatrix} \displaystyle{ 675} \cr \cr \displaystyle{ 1}\end{pmatrix}}\end{bmatrix}} \cr \cr \displaystyle{ \begin{bmatrix} \displaystyle{ 210937500 {a}6 {b}27} \cr \cr \displaystyle{ 284765625 {a}8 {b}24}\end{bmatrix}}\end{pmatrix} \end{equation}
This is an example of how we might structure an org-mode file that can export out to Github flavored Markdown, or to a PDF.
First, let’s get some code loaded up and written. Here’s a block that converts polar coordinates to rectangular coordinates.
(define (p->r local)
(let* ((polar-tuple (coordinate local))
(r (ref polar-tuple 0))
(phi (ref polar-tuple 1))
(x (* r (cos phi)))
(y (* r (sin phi))))
(up x y)))This is some good stuff.
(load "ch1/utils.scm")
<<p->r>>
<<spherical->rect>>And another, that gets us from spherical to rectangular.
(define (spherical->rect local)
(let* ((spherical-tuple (coordinate local))
(r (ref spherical-tuple 0))
(theta (ref spherical-tuple 1))
(phi (ref spherical-tuple 2)))
(up (* r (sin theta) (cos phi))
(* r (sin theta) (sin phi))
(* r (cos theta)))))#+RESULTS[f4f039075baf66ba4fe071844815bfcffe281eaa]:
;Loading "ch1/utils.scm"... done #| "" |#
This block will generate a LaTeX version of the code I’ve supplied:
(->tex-equation
((+ (literal-function 'c)
(D (literal-function 'z)))
't)
"eq:masterpiece")#+RESULTS[b383d2f5d6c252ac04a5f44aaeaec678132b8449]: \begin{equation} c\left( t \right) + Dz\left( t \right) \label{eq:masterpiece} \end{equation}
You can even reference these with equation numbers, like Equation \eqref{eq:masterpiece} above.
(up 1 2 't)#| (up 1 2 t) |#
Here’s (a test) of
And again this is a thing:
\[ eiπ = -1 \]
\[ ∫_0^∞ e-x^2 dx = \frac{\sqrt{π}}{2} \]