# Procedural approach to classical mechanics

## Table of Contents

**Note.** For more details see:

- Structure and Interpretation of Classical Mechanics (SICM) by G. J. Sussman and J. Wisdom (2001). The book revisits the fundamentals of classical mechanics, integrating theoretical concepts with formalized computational procedures.
- The Role of Programming in the Formulation of Ideas by G.J. Sussman and J. Wisdom 2002, AI Lab Technical Report, AIM-2002-018. The report focuses in details on the Euler-Lagrange equations example.

## 1 Traditional notation

Physics notation is often ambiguous. Unstated notational conventions
are effective tools if the context is clear. However, a formal
definition of procedures about *how* to automatically compute (as by a
computer) a given quantity requires unequivocal syntax.

As an example, let's consider Euler-Lagrange equations. In the traditional notation they read:

\begin{equation} \frac{d}{dt} \frac{\partial L}{\partial \dot{q}^i} - \frac{\partial L}{\partial q^i} = 0 \;. \end{equation}A dumb computer will have issues in accepting the above notation:

- We think of coordinates \(q^i\) and velocities \(\dot{q}^i\) as independent when taking partial derivatives. But we also now that \(\frac{d}{dt} q^i = \dot{q}^i\): here the interpretation changes and positions and velocities are not independent.
- It contains an obvious typo: the Lagrangian (and its partial derivatives) is a function of multiple variables, but \(d/dt\) only applies to functions of one variable.

The above notation works with the implicit assumption that we mean to substitute the path and its derivative (both are functions of time, so that the total time derivative makes sense) into the coordinate and velocity arguments of the partial derivative functions. If \(w^i\) is a path through the coordinate configuration space, then

\begin{equation} \frac{d}{dt} \left( \left. \frac{\partial L(t, q^i, \dot{q}^i)}{\partial \dot{q}^i} \right|_{q^i=w^i(t), \dot{q}^i=\frac{dw^i(t)}{dt}} \right) - \left. \frac{\partial L}{\partial q^i} \right|_{q^i=w^i(t), \dot{q}^i=\frac{dw^i(t)}{dt}} = 0 \;. \end{equation}The use of ambiguous notation is convenient in this simple case. Still, implicit and equivocal mathematical notation is an obstacle to interpret the equation above programmatically.

## 2 Functional notation

M. Spivak in *Calculus on Manifolds* argues in favor of a functional
notation. For instance, the chain rule would read \(D(f \circ g) =
((Df) \circ g) \cdot Dg\), so at \(a\):

Given the proper basic symbols definitions, this expression is unequivocal and it does not require the introduction of irrelevant letters (as is the case of the traditional notation).

The main advantage of functional notation over traditional notation is
that it is well adapted to be interpreted automatically, independently
from the context. Programming languages cannot deal in a robust,
general way with apparent contradicting statements that a human would
know how to interpret given the context. (There are exceptions,
though. See for instance the prototypical language discussed in
G. J. Sussman and A. Radul, *The Art of the Propagator*, 2009.)

### 2.1 Euler-Lagrange equations

Using a functional notation, Euler-Lagrange equations read as follows:

\begin{equation} D(\partial_2L \circ \Gamma[q]) - \partial_1L \circ \Gamma[q] = 0 \;. \end{equation}The Lagrangian is a function of the tuple composed by the time, coordinates and velocities \(L(t, x, v)\). The partial derivative suffix indicates the argument position in the \((t, x, v)\) tuple. The function \(\Gamma\) maps time to the local tuple \(\Gamma[q](t) = (t, q(t), Dq(t))\), given the coordinate path \(q\). This notation can be easily written in terms of a robust procedural syntax:

(define ((Lagrange-equations Lagrangian) q) (- (D (compose ((partial 2) Lagrangian) (Gamma q))) (compose ((partial 1) Lagrangian) (Gamma q))))

This defines a procedure called `Lagrange-equations`

that takes a
`Lagrangian`

and a coordinate path `q`

as inputs to compute the
Euler-Lagrange equations. Without going into the details of the syntax
definition, let's stress that the first argument of each parenthesis
names a procedure. For instance, `(+ 1 2)`

indicates the sum of the
integers 1 and 2, which evaluates to 3. Internal lists are evaluated
first: given the expression `(+ 1 (* 2 3))`

, first the product between
2 and 3 is evaluated to 6, which is then summed to 1.

Note the `compose`

procedure to compose the `Lagrangian`

procedure
with the one resulting from evaluating `(Gamma q)`

. As mentioned
above, given a coordinate path (here `up`

indicates an up-tuple,
think about it as a vector):

(define q (up (literal-function 'x) (literal-function 'y) (literal-function 'z)))

the procedure `Gamma`

returns the local tuple composed by the time,
the coordinate tuple and the velocity tuple:

(print-expression ((Gamma q) 't)) #| (up t (up (x t) (y t) (z t)) (up ((D x) t) ((D y) t) ((D z) t))) |#

(The result of evaluating the procedure is shown enclosed by
`#|...|#`

.)

The procedure defined for the Euler-Lagrange equation can be applied to any Lagrangian. Let's consider a couple of examples.

#### 2.1.1 The free particle

Let's apply the procedure defined above to the Lagrangian \(L=\frac{1}{2} m v^2\) for a free particle of mass \(m\). The corresponding procedure is:

```
(define ((L-free-particle mass) local)
(let ((v (velocity local)))
(* 1/2 mass (dot-product v v))))
```

Let's compose the Euler-Lagrange equations procedure with the free particle Lagrangian

(show-expression (((Lagrange-equations (L-free-particle 'm)) (literal-function 'x)) 't)) #| (* m (((expt D 2) x) t)) |#

This provides the expected equation of motion for a free particle \(mD^2x(t)=0\).

#### 2.1.2 The harmonic oscillator

The Lagrangian \(L(t, q, v)=\frac{1}{2}mv^2 - \frac{1}{2}kq^2\) of a particle of mass \(m\) constrained at the end of a spring (characterized by a spring constant \(k\)) is described by the following procedure :

```
(define ((L-harmonic m k) local)
(let ((q (coordinate local))
(v (velocity local)))
(- (* 1/2 m (square v)) (* 1/2 k (square q)))))
```

Let's check that the motion of a harmonic oscillator \(x(t)=a cos(\omega t)+\phi\):

(define (proposed-solution t) (* 'a (cos (+ (* 'omega t) 'phi))))

satisfies the Euler-Lagrange equations given the Lagrangian above:

(show-expression (((Lagrange-equations (L-harmonic 'm 'k)) proposed-solution) 't))

This gives

\begin{equation} a \cos\left( \omega t + \phi \right) (k - m {\omega}^{2}) \end{equation}As expected, the motion is harmonic and the solution allowed are those satisfying \(\omega^2=k/m\).

## 3 About the procedural notation

The procedures formalized in the examples above are robust, i.e., not bound to a specific context. The Euler-Lagrange equations procedure generalizes automatically to arbitrary Lagrangians or coordinate paths. It can be composed with other arbitrary functions, too (as long as their domains overlap).

The procedural syntax used in the examples is the Scheme dialect of
the Lisp family of programming languages. Lisp was designed by John
McCarthy in 1960 (which makes it one of the older programming
languages still in use, second only to Fortran). The name stands for
*List Processing*, because one of the key ideas is to use a simple
data structure called a list for both code and data. The Roots of Lisp
by P. Graham (2001) revisits McCarthy's original article illustrating
the code/data evaluation procedure.

Here we skipped several basic definitions and only aimed to give a flavor about the Scheme notation used in SICM. For more details, refer to:

- SICM, appendix
*Scheme*: an introduction to the Scheme programming language. - SICM, appendix
*Our notation*: details about the functional mathematical notation adopted by the book and how generic arithmetic extends basic operations, such as addition and multiplication, to a variety of mathematical types. The discussion about differentiation and of tuple (a generalization of tensors) arithmetic is also interesting.

## 4 Appendix: Procedures composition

**Note.** This section is rather technical and aimed to illustrate a
subtlety used above that shows the flexibility of the Lisp syntax
(despite its simplicity).

The procedure `L-free-particle`

was defined as a function of `mass`

and of the `local`

tuple as:

```
(define ((L-free-particle mass) local)
(let ((v (velocity local)))
(* 1/2 mass (dot-product v v))))
```

It may seem simpler to define it as ```
(define (L-free-particle mass
local) ...)
```

. However, we want to compose `(L-free-particle mass)`

with other functions by itself.

Let's consider a minimal example. The following definitions are
equivalent and bound the procedure `(lambda (x y) (+ x y))`

(where
`lambda`

denotes an anonymous function) to the symbol `foo`

:

(define foo (lambda (x y) (+ x y))) (define (foo x y) (+ x y)) ; Same as before, using syntactic sugar.

Let's consider also the following equivalent definitions:

(define bar (lambda (x) (lambda (y) (+ x y)))) (define (bar x) (lambda (y) (+ x y))) ; Syntactic sugar. (define ((bar x) y) (+ x y)) ; More syntactic sugar.

For the simple purpose of evaluating \(x+y\) (e.g., for \(x=2\) and \(y=3\)), both cases work:

(foo 2 3) ;Value: 5 ((bar 2) 3) ;Value: 5

However, the second case allows to compose `bar`

with other
functions of one argument. For instance, let's compute \(x^3+y\) (for
\(x=2\) and \(y=3\)) as:

(((compose bar cube) 2) 3) ;Value: 11

A similar composition of `foo`

with `cube`

is not possible because
`foo`

requires two arguments, while `cube`

one.

It is remarkable how a simple idea such as using lists for both
code and data naturally leads to (*a priori* non-trivial) procedure
composition.

## 5 Appendix: `scmutils`

The examples shown above are based on the powerful `scmutils`

library
for the Scheme programming language. The library is *libre* software,
developed for education and research at MIT.

Once installed, all the procedures illustrated above can be easily evaluated through the interactive Scheme REPL (read-eval-print loop), similar to (but way more powerful than) an interactive shell.

### 5.1 Installation

Official installation instructions are provided here. However, a local installation may be preferable:

- Download the system from the official website.
Extract the tar ball in a local folder (e.g.,

`~/scmutils`

):tar -xvzf [your-tarball.tar.gz]

This will extract two folders:

`bin`

and`scmutils`

.- Now you have two options: access the Scheme REPL with Edwin (an
Emacs-like editor written in Scheme) or Emacs.
Edwin. Set the

`SCMUTILS_ROOT`

environment variable to override the default location (`/usr/local`

) for the`scmutils`

root directory:SCMUTILS_ROOT=/path/to/scmutils export SCMUTILS_ROOT

Launch

`bin/mechanics`

script.Emacs. Copy the following

`mechanics`

Elisp command (e.g., in`~/.emacs`

or in the`*scratch*`

buffer) and evaluate it (change the root path below to the full one corresponding to the`scmutils`

folder extracted from the tar ball):(defun mechanics () (interactive) (let ((root "/path/to/scmutils/")) (run-scheme (concat root "/mit-scheme/bin/scheme --library " root "/mit-scheme/lib"))))

Launch the command as

`M-x mechanics`

from any Emacs buffer. This will open the Scheme REPL with the`scmutils`

library already available.

- If everything seems to work you can now remove the tarball.

### 5.2 Usage

After having launched the Scheme REPL in Ediwn/Emacs, open a new
buffer `C-x 2`

. Open a source code file `C-x C-f demo.scm`

. After
adding some code, evaluate and send sexps to the REPL by typing ```
C-x
C-e
```

at the end of the last parenthesis, or `C-M-x`

within a given
sexp. (This assumes that only default modes are installed in
Emacs. If, for instance, Geiser mode is also installed, see its
usage instructions or toggle `geiser-mode`

to disable it.)

To quit the REPL, type `(exit)`

.

Useful documents are included in the `scmutils/manual/`

folder. Scheme sources for the entire system are included in the
`scmutils/src/`

folder.

**Tip**: when editing Scheme source code in Emacs, install the
useful paredit extension.