Skip to main content

The Ordinary Differential Equations Project

Thomas W. Judson

Section 5.2 Hamiltonian Systems

An undamped harmonic oscillator, \(m y'' + k y = 0\text{,}\) can be written as the system
\begin{align*} y' \amp = v\\ v' \amp = -\frac{k}{m} y. \end{align*}
Now suppose that \((y(t), v(t))\) is a solution curve in the \(yv\)-plane. We will calculate the slope of the solution curve, \(dv/dy\text{.}\) Using the fact from calculus that the derivative of an inverse function is
\begin{equation*} \frac{d}{dx} f^{-1}(x) = \frac{1}{f'(x)}, \end{equation*}
we have
\begin{equation*} \frac{dv}{dy} = \frac{dv}{dt} \cdot \frac{dt}{dy} = \frac{dv/dt}{dy/dt} = - \frac{ky}{mv}. \end{equation*}
In general for the system
\begin{align*} \frac{dx}{dt} \amp = f(x, y)\\ \frac{dy}{dt} \amp = g(x, y), \end{align*}
we have
\begin{equation*} \frac{dy}{dx} = \frac{g(x, y)}{f(x, y)}. \end{equation*}
Using separation of variables to solve the equation
\begin{equation*} \frac{dv}{dy} = - \frac{ky}{mv}, \end{equation*}
we obtain the solution
\begin{equation} \frac{1}{2} m v^2 + \frac{1}{2} k y^2 = C,\tag{5.2.1} \end{equation}
where \(C\) is a constant. The first term of (5.2.1) is the kinetic energy function of the harmonic oscillator
\begin{equation*} K = \frac{1}{2} m v^2, \end{equation*}
while the second term is the potential energy function
\begin{equation*} U = \int_0^y k s \, ds = \frac{1}{2} k y^2. \end{equation*}
For this reason, we call the function
\begin{equation*} E = K + U = \frac{1}{2} m v^2 + \frac{1}{2} k y^2 \end{equation*}
the total energy function of the harmonic oscillator. Equation (5.2.1) tells us that energy is conserved. That is, the sum of the potential energy and the kinetic energy is constant.

Subsection 5.2.1 The Nonlinear Pendulum

While pendulums have long been used in clocks to keep time, they have also been used to measure gravity as well as used in early seismometers to measure the effect of earthquakes. One of the more interesting uses of a pendulum has been to measure the rotation of the earth. In 1851, the French physicist, Léon Foucault, used a pendulum (Figure 5.2.1) to demonstrate that the earth actually rotated on its axis. The fact that the earth rotates had been known for a long time, but Foucault’s experiment gave the first simple proof of this phenomena.
a photo of the Foucault Pendulum Clock at the California Academy of Sciences
Figure 5.2.1. Foucault Pendulum Clock—California Academy of Sciences
Let us consider a pendulum made of a light rod of length \(L\) called the arm of the pendulum with a mass \(m\) on the end of the rod called the bob (Figure 5.2.2). We will ignore the mass of the arm in our system. The position of the bob is given by \(\theta(t)\text{,}\) which we will measure in a counterclockwise direction. We will assume that \(\theta(0)\) is when the bob is in the vertical position.
a diagram of a pendulum with an arm of length L and displacement of theta and a bob of mass m
Figure 5.2.2. The nonlinear pendulum
There are two forces acting on the pendulum—gravity and friction. The position of the bob at time \(t\) is
\begin{equation*} (L \sin \theta(t), - L \cos \theta(t)), \end{equation*}
and the velocity of the bob is \(L (d \theta / dt)\text{,}\) the length of the velocity vector (Figure 5.2.3). The component of the acceleration that points along the direction of motion of the bob is
\begin{equation*} L \frac{d^2 \theta}{dt^2}. \end{equation*}
We can take the force due to friction to be proportional to the velocity,
\begin{equation*} - b L \frac{d \theta}{dt}, \end{equation*}
where \(b \geq 0\text{.}\) Thus, Newton’s second law tells us that
\begin{equation*} m L \frac{d^2 \theta}{dt^2} = -b L \frac{d \theta}{dt} - mg \sin \theta \end{equation*}
or
\begin{equation*} \frac{d^2 \theta}{dt^2} + \frac{b}{m} \frac{d \theta}{dt} + \frac{g}{L} \sin \theta = 0. \end{equation*}
As a system, we can model the pendulum as
\begin{align*} \frac{d \theta}{dt} \amp = v\\ \frac{dv}{dt} \amp = - \frac{b}{m} v - \frac{g}{L} \sin \theta. \end{align*}
a diagram of a pendulum with an arm of length L and displacement of theta and a bob of mass m and the forces acting on the mass
Figure 5.2.3. Forces on the nonlinear pendulum
If there is no damping for our pendulum, we say that we have an ideal pendulum. In this case,
\begin{align*} \frac{d \theta}{dt} \amp = v\\ \frac{dv}{dt} \amp = - \frac{g}{L} \sin \theta. \end{align*}
Since
\begin{equation*} \frac{dv}{d \theta} = - \frac{g \sin \theta}{Lv}, \end{equation*}
we know that the energy function for the pendulum is
\begin{equation*} E(\theta, v) = \frac{1}{2} v^2 - \frac{g}{L} \cos \theta. \end{equation*}
Of course, \(E(\theta, v)\) is constant on the ideal pendulum.

Subsection 5.2.2 Hamiltonian Systems

The ideal pendulum and the undamped harmonic oscillator are examples of Hamiltonian systems. More specifically, a Hamiltonian system is a system of the form
\begin{align*} \frac{dx}{dt} \amp = \frac{\partial H}{\partial y} (x, y)\\ \frac{dy}{dt} \amp = -\frac{\partial H}{\partial x}(x, y), \end{align*}
where \(H : {\mathbb R}^2 \rightarrow {\mathbb R}\) is a smooth function, sometimes referred to as a Hamiltonian function.

Example 5.2.4.

The system
\begin{align*} \frac{dx}{dt} \amp = -2x - 3y^2\\ \frac{dy}{dt} \amp = -3x^2 + 2y \end{align*}
is Hamiltonian since for \(H(x, y) = x^3 - 2xy - y^3\)
\begin{align*} \frac{\partial H}{\partial y} \amp = -2x - 3y^2\\ \frac{\partial H}{\partial x} \amp = -(-3x^2 + 2y). \end{align*}

Example 5.2.5.

In the case of the harmonic oscillator and ideal pendulum, \(H\) is just the energy function. For the harmonic oscillator, the Hamiltonian function is
\begin{equation*} E = K + U = \frac{1}{2} m v^2 + \frac{1}{2} k y^2. \end{equation*}
In this case,
\begin{align*} \frac{\partial E}{\partial v} \amp = y' = mv\\ -\frac{\partial E}{\partial y} \amp = v' = -k y. \end{align*}
A change of variables (\(y = my\)) tells us that this system is equivalent to the system
\begin{align*} y' \amp = v\\ v' \amp = -\frac{k}{m} y. \end{align*}
For the ideal pendulum
\begin{align*} \frac{\partial E}{\partial v} \amp =\frac{d \theta}{dt} = v\\ -\frac{\partial E}{\partial \theta} \amp =\frac{dv}{dt} = - \frac{g}{L} \sin \theta, \end{align*}
where
\begin{equation*} E(\theta, v) = \frac{1}{2} v^2 - \frac{g}{L} \cos \theta. \end{equation*}
The following theorem tells the importance of Hamiltonian systems. That is, the solution curves of the system are simply the level sets of the Hamiltonian function.

Proof.

Let \((x(t), y(t))\) be a solution curve for the system. Using the chain rule
\begin{align*} \frac{dH}{dt} \amp = \frac{\partial H}{\partial x} \frac{dx}{dt} + \frac{\partial H}{\partial y} \frac{dy}{dt}\\ \amp = \frac{\partial H}{\partial x} \left(\frac{\partial H}{\partial y} \right)+ \frac{\partial H}{\partial y} \left( -\frac{\partial H}{\partial x} \right)\\ \amp = 0. \end{align*}
Theorem 5.2.6 tells us how to draw the solution curves in the phase plane without solving the system. Assuming that the Hamiltonian function \(H\) is not constant on any open set in \({\mathbb R}^2\text{,}\) we simply need to plot the level curves, \(H(x, y) = C\text{.}\) The solutions of the system live on these level sets, and all we need to do is find the direction of the solution curve. However, this is quite easy since we know the vector field of the system. Furthermore, the equilibrium points of the Hamiltonian system occur at the critical points of \(H\) (where the partials of \(H\) vanish). For example, we can see the solution curves of
\begin{align*} \frac{\partial H}{\partial y} \amp = -2x - 3y^2\\ \frac{\partial H}{\partial x} \amp = -(-3x^2 + 2y). \end{align*}
in Figure 5.2.7.
solution curves on a direction field of slope arrows
Figure 5.2.7. Solution curves of \(x' = x - 3 y^2\) and \(y' = -y\)

Subsection 5.2.3 The Ideal Pendulum Revisited

The ideal pendulum
\begin{align*} \frac{d \theta}{dt} & = v\\ \frac{dv}{dt} & = - \frac{g}{L} \sin \theta. \end{align*}
is Hamiltonian, since the total energy function
\begin{equation} H(\theta, v) = E(\theta, v) = \frac{1}{2} v^2 - \frac{g}{L} \cos \theta.\tag{5.2.2} \end{equation}
is the Hamiltonian function for the system. If we set \(g/L = 1\text{,}\) then the solution curves of the system are just the level curves of (5.2.2). In Figure 5.2.8, the closed ellipses correspond to the normal motion of a pendulum, while cosine curves correspond to a pendulum that always rotates in the same direction. The curves that join the equilibrium points correspond to the pendulum that rotates exactly to the top of the arc and then rotates back in the other direction.
solution curves on a direction field of slope arrows for the ideal pendulum
Figure 5.2.8. The ideal pendulum
Hamiltonian systems are rather rare. Given a system, we need a test to see if it is indeed Hamiltonian. For the system
\begin{align*} \frac{dx}{dt} & = f(x, y)\\ \frac{dy}{dt} & = g(x, y), \end{align*}
we wish to find a function \(H(x, y)\) such that
\begin{align*} \frac{\partial H}{\partial y} & = f(x, y)\\ \frac{\partial H}{\partial x} & = -g(x, y). \end{align*}
If such a function exists, then
\begin{equation*} \frac{\partial f}{\partial x} = \frac{\partial^2 H}{\partial x \partial y} = \frac{\partial^2 H}{\partial y \partial x} = - \frac{\partial g}{\partial y}. \end{equation*}
Thus,
\begin{equation*} \frac{\partial f}{\partial x} = - \frac{\partial g}{\partial y}, \end{equation*}
if our system is Hamiltonian.

Example 5.2.9.

We already know that the system
\begin{align*} \frac{dx}{dt} \amp = -2x - 3y^2\\ \frac{dy}{dt} \amp = -3x^2 + 2y \end{align*}
is Hamiltonian. Indeed,
\begin{equation*} \frac{\partial f}{\partial x} = - \frac{\partial g}{\partial y} = -2. \end{equation*}
The condition that
\begin{equation*} \frac{\partial f}{\partial x} = - \frac{\partial g}{\partial y}, \end{equation*}
is also sufficient for
\begin{align*} \frac{dx}{dt} & = f(x, y)\\ \frac{dy}{dt} & = g(x, y), \end{align*}
to be a Hamiltonian system. That is, we can construct a Hamiltonian function if \(\partial f/ \partial x = - \partial g/ \partial y\text{.}\) If
\begin{equation*} f(x, y) = \frac{\partial H}{\partial y}(x, y), \end{equation*}
then
\begin{equation*} H(x, y) = \int f(x, y) \, dy + \phi(x). \end{equation*}
Differentiating \(H\) with respect to \(x\text{,}\) tells us that
\begin{equation*} \frac{\partial H}{\partial x}(x, y) = \frac{\partial}{\partial x} \int f(x, y) \, dy + \phi'(x) = - g(x, y). \end{equation*}
Now we can determine \(\phi(x)\) by solving the first-order differential equation
\begin{equation*} \phi'(x) = - g(x, y) - \frac{\partial}{\partial x} \int f(x, y) \, dy. \end{equation*}
We summarize our results in the following theorem.

Example 5.2.11.

As an example, we will show how we constructed the Hamiltonian function given in Example 5.2.4 for the system
\begin{align*} \frac{dx}{dt} \amp = -2x - 3y^2\\ \frac{dy}{dt} \amp = -3x^2 + 2y \end{align*}
First,
\begin{equation*} H(x, y) = \int f(x, y) \, dy + \phi(x) = \int (-2x - 3y^2) \, dy + \phi(x) = -2xy - y^3 + \phi(x). \end{equation*}
Consequently,
\begin{align*} \phi'(x) \amp = - g(x, y) - \frac{\partial}{\partial x} \int f(x, y) \, dy\\ \amp = -(-3x^2 + 2y) - \frac{\partial}{\partial x}(-2xy - y^3)\\ \amp = 3x^2 - 2y + 2y\\ \amp = 3x^2, \end{align*}
and \(\phi(x) = x^3 + C\text{,}\) where \(C\) is any constant. If we choose \(C = 0\text{,}\) our Hamiltonian function is
\begin{equation*} H(x, y) = -2xy - y^3 + \phi(x) = -2xy - y^3 + x^3. \end{equation*}

Example 5.2.12.

As a second example, suppose
\begin{align*} \frac{dx}{dt} & = -x^2 \cos y \sin x\\ \frac{dy}{dt} & = x^2 \cos x \sin y + 2x(\sin x \sin y - 1). \end{align*}
Since
\begin{equation*} \frac{\partial f}{\partial x} = x^2 \cos x \cos y - 2x \cos y \sin x = - \frac{\partial g}{\partial y} \end{equation*}
our system is Hamiltonian. To find a Hamiltonian, we compute
\begin{align*} H(x, y) \amp = \int f(x, y) \, dy + \phi(x)\\ \amp = \int (-x^2 \cos y \sin x) \, dy + \phi(x)\\ \amp = -x^2 \sin x \sin y + \phi(x). \end{align*}
Thus,
\begin{align*} \phi'(x) & = - g(x, y) - \frac{\partial}{\partial x} \int f(x, y) \, dy\\ & = -(x^2 \cos x \sin y + 2x(\sin x \sin y - 1))\\ \amp - ( 2x \sin x \sin y - x^2 \cos x \sin y)\\ & = 2x, \end{align*}
and \(\phi(x) = x^2 + C\text{.}\) Letting \(C = 0\text{,}\) our Hamiltonian function is
\begin{equation*} H(x, y) = x^2 - x^2 \sin x \sin y = x^2(1 - \sin x \sin y). \end{equation*}
Solution curves for this system are given in Figure 5.2.13.
solution curves as level sets on a direction field of slope arrows
Figure 5.2.13. Solution curves for Example 5.2.12

Activity 5.2.1. Hamiltonian Systems.

For each of the systems,
  1. Show that the system is Hamiltonian.
  2. Find a Hamiltonian function \(H(x, y)\) for the system.
  3. Plot the level sets of \(H(x, y)\text{.}\)
  4. Plot the phase portrait of the system, describing any equilibrium solutions.
     1 
    Use technology.
(a)
\begin{align*} \frac{dx}{dt} & = -e^x \sin y\\ \frac{dy}{dt} & = -e^x \cos y \end{align*}
(b)
\begin{align*} \frac{dx}{dt} & = x^2 \cos x \sin y\\ \frac{dy}{dt} & = -x^2 \cos y \sin x + 2x(\cos x \cos y - 1) \end{align*}

Subsection 5.2.4 Equilibrium Solutions of Hamiltonian Systems

The system that we studied in Example 5.2.4,
\begin{align*} \frac{dx}{dt} \amp = -2x - 3y^2\\ \frac{dy}{dt} \amp = -3x^2 + 2y \end{align*}
has equilibrium solutions at \((0,0)\) and \((-2/3, 2/3)\text{.}\) Referring to the phase portrait for the system (Figure 5.2.7), we might guess that \((0,0)\) is a saddle, and \((-2/3, 2/3)\) is either a center or a spiral sink. Since solution curves must follow contours for the Hamiltonian function \(H(x, y) = x^3 - 2xy - y^3\text{,}\) spiral sinks do not make sense.
Let us examine the possible types of equilibrium solutions for a Hamiltonian system. Suppose that \((x_0, y_0)\) is an equilibrium solution for the system
\begin{align*} \frac{dx}{dt} & = f(x, y) = \frac{\partial H}{\partial y} (x, y)\\ \frac{dy}{dt} & = g(x, y) = -\frac{\partial H}{\partial x}(x, y). \end{align*}
In order to determine the nature of the equilibrium solution, we will compute the Jacobian matrix of the system at \((x_0, y_0)\text{,}\)
\begin{equation*} J = \begin{pmatrix} f_x(x_0, y_0) & f_y(x_0, y_0) \\ g_x(x_0, y_0) & g_y(x_0, y_0) \end{pmatrix}= \begin{pmatrix} H_{yx}(x_0, y_0) & H_{yy}(x_0, y_0) \\ -H_{xx}(x_0, y_0) & -H_{xy}(x_0, y_0) \end{pmatrix}. \end{equation*}
If we let
\begin{align*} \alpha & = H_{xy}(x_0, y_0) = H_{yx}(x_0, y_0)\\ \beta & = H_{yy}(x_0, y_0)\\ \gamma & = -H_{xx}(x_0, y_0), \end{align*}
then \(J\) becomes
\begin{equation*} \begin{pmatrix} \alpha & \beta \\ \gamma & - \alpha \end{pmatrix}. \end{equation*}
The characteristic polynomial of this matrix is
\begin{equation*} \lambda^2 - \alpha^2 - \beta \gamma. \end{equation*}
Therefore, our matrix has eigenvalues
\begin{equation*} \lambda = \pm \sqrt{\alpha^2 + \beta \gamma}, \end{equation*}
and we have the following possibilities.
  • If \(\alpha^2 + \beta \gamma \gt 0\text{,}\) we have two real eigenvalues of opposite signs. Therefore, our equilibrium solution is a saddle.
  • If \(\alpha^2 + \beta \gamma = 0\text{,}\) the only eigenvalue is zero.
  • If \(\alpha^2 + \beta \gamma \lt 0\text{,}\) we have two purely imaginary eigenvalues.
In particular, a Hamiltonian system has no spiral sinks or sources.
Returning to Example 5.2.4, the Jacobian matrix is
\begin{equation*} J = \begin{pmatrix} H_{yx}(x_0, y_0) \amp H_{yy}(x_0, y_0) \\ -H_{xx}(x_0, y_0) \amp -H_{xy}(x_0, y_0) \end{pmatrix} = \begin{pmatrix} -2 \amp -6y_0 \\ -6x_0 \amp 2 \end{pmatrix}. \end{equation*}
For the equilibrium solution \((0,0)\text{,}\) we have eigenvalues \(\lambda = \pm 2\text{,}\) which tells us that the origin is a nodal saddle. The Jacobian matrix corresponding to \((x_0, y_0) = (-2/3, 2/3)\) is
\begin{equation*} J = \begin{pmatrix} -2 \amp -4 \\ -4 \amp 2 \end{pmatrix}. \end{equation*}
Since \(J\) has eigenvalues \(\lambda = \pm 2 \sqrt{5} \, i\text{.}\) Therefore, solution curves near \((-2/3, 2/3)\) should look like closed orbits.

Subsection 5.2.5 Important Lessons

  • A pendulum can be modeled by the equation
    \begin{equation*} \frac{d^2 \theta}{dt^2} + \frac{b}{m} \frac{d \theta}{dt} - \frac{g}{L} \sin \theta = 0. \end{equation*}
    or the system
    \begin{align*} \frac{d \theta}{dt} & = v\\ \frac{dv}{dt} & = - \frac{b}{m} v - \frac{g}{L} \sin \theta. \end{align*}
  • A Hamiltonian system is a system of the form
    \begin{align*} \frac{dx}{dt} & = \frac{\partial H}{\partial y} (x, y)\\ \frac{dy}{dt} & = -\frac{\partial H}{\partial x}(x, y), \end{align*}
    where \(H : {\mathbb R}^2 \rightarrow {\mathbb R}^2\) is a smooth function.
  • If
    \begin{align*} \frac{dx}{dt} & = \frac{\partial H}{\partial y} (x, y)\\ \frac{dy}{dt} & = -\frac{\partial H}{\partial x}(x, y), \end{align*}
    is a Hamiltonian system, then \(H\) is constant along any solution curve. In particular, the solution curves of a Hamiltonian system are the level sets of \(H\text{.}\)
  • If the system
    \begin{align*} \frac{dx}{dt} & = f(x, y)\\ \frac{dy}{dt} & = g(x, y), \end{align*}
    is Hamiltonian, then
    \begin{equation*} \frac{\partial f}{\partial x} = - \frac{\partial g}{\partial y}. \end{equation*}
  • If \(f_x = - g_y\) is true on a rectangle \(R\text{,}\) then the system
    \begin{align*} \frac{dx}{dt} & = f(x, y)\\ \frac{dy}{dt} & = g(x, y), \end{align*}
    is Hamiltonian on \(R\) and
    \begin{equation*} H(x, y) = \int f(x, y) \, dy + \phi(x), \end{equation*}
    where
    \begin{equation*} \phi'(x) = -\frac{\partial}{\partial x}\left( \int f(x, y) \, dy \right) - g(x, y). \end{equation*}
  • A Hamiltonian system has no spiral sinks or spiral sources.

Reading Questions 5.2.6 Reading Questions

1.

Explain in your own words what a Hamiltonian system is.

2.

What types of equilibrium solutions are possible in a Hamiltonian system?

Exercises 5.2.7 Exercises

Hamiltonian Systems.

For each of the systems in Exercise Group 5.2.7.1–4,
  1. Show that the system is Hamiltonian.
  2. Find a Hamiltonian function \(H(x, y)\) for the system.
  3. Plot the level sets of \(H(x, y)\text{.}\)
  4. Plot the phase portrait of the system, describing any equilibrium solutions.
1.
\begin{align*} \frac{dx}{dt} & = 3y^2\\ \frac{dy}{dt} & = -3x^2 \end{align*}
2.
\begin{align*} \frac{dx}{dt} & = 3y^2 - 4x\\ \frac{dy}{dt} & = 3x^2 + 4y \end{align*}
3.
\begin{align*} \frac{dx}{dt} & = x \cos y + 2y\\ \frac{dy}{dt} & = -\sin y \end{align*}
4.
\begin{align*} \frac{dx}{dt} & = -x \sin y - \cos y\\ \frac{dy}{dt} & = - \cos y \end{align*}

5.

Consider an ideal pendulum with bob mass \(m\) and arm length \(l\) given by
\begin{align*} \frac{d\theta}{dt} & = v\\ \frac{dv}{dt} & = - \frac{g}{l} \sin \theta. \end{align*}
  1. What is the linearization of the ideal pendulum system at the equilibrium point \((0, 0)\text{?}\)
  2. Using \(g = 9.8\) m/s\(^2\text{,}\) how should \(l\) and \(m\) be chosen so that small swings of the pendulum have period 1 second?
  3. For the linearization of the ideal pendulum at \((0, 0)\text{,}\) the period of oscillation is independent of the amplitude. Does the same statement hold for the ideal pendulum itself? Is the period of oscillation the same no matter how high the pendulum swings? If not, will the period be shorter or longer for high swings?
  4. An ideal pendulum clock— a clock containing an ideal pendulum that “ticks” once for each swing of the pendulum arm— keeps perfect time when the pendulum makes very high swings. Will the clock rum fast or slow if the amplitude of the swings is very small?
  5. If the arm length of the ideal pendulum is doubled from \(l\) to \(2l\text{,}\) what is the effect on the period of small amplitude swinging solutions?
  6. What is the rate of change of the period of small amplitude swings as \(l\) varies?
  7. Will an ideal pendulum clock that keeps perfect time on earth run fast or slow on the moon?

6.

Describe the phase portrait for
\begin{align*} x' & = x^2 - 1\\ y' & = -xy + a(x^2 -1), \end{align*}
when \(a \lt 0\text{.}\)

7.

Consider the autonomous system
\begin{align*} x' & = y^2 + \cos x\\ y' & =1 + 2x + y \sin x. \end{align*}
  1. Show that this system is a Hamiltonian system.
  2. Find a Hamiltonian function for this system.
PrevTopNext