Solving the Pendulum ODE with Python

Introduction

Pendulums are a classic example of physics in motion, embodying the elegance of mechanics and mathematics. From the simple swing of a playground swing to the precise ticking of a clock, pendulums have fascinated scientists and engineers for centuries. Their motion is governed by a second-order ordinary differential equation (ODE), which encapsulates the forces acting on the system.

In this post, we’ll explore the mathematical framework of the pendulum’s motion and show how to solve the governing equation using Python. By leveraging powerful libraries like SciPy for computation and Matplotlib for visualization, we’ll demonstrate how numerical methods can transform complex equations into meaningful results. Whether you’re an engineering student, a Python enthusiast, or someone curious about mechanics, this tutorial will guide you step by step.

What Makes the Pendulum Special?

The pendulum is more than a theoretical construct—it’s a cornerstone of mechanics. Its study has led to groundbreaking discoveries in fields like timekeeping, seismology, and even space exploration. For small displacements, the motion of a pendulum approximates simple harmonic motion, but for larger angles, the dynamics become nonlinear, revealing a richer and more complex behavior.

Understanding the pendulum’s motion is not just an academic exercise; it provides insight into the fundamental principles of oscillatory systems, which appear in everything from bridges and skyscrapers to electronic circuits and molecular vibrations.

The Pendulum Equation

The motion of a simple pendulum is described by the equation:

\(
\frac{d^2\theta}{dt^2} + \frac{g}{l}\sin\theta = 0
\)

Where:

  • \(\theta\): Angle of the pendulum from the vertical (radians).
  • \(g\): Acceleration due to gravity (m/s²).
  • \(l\): Length of the pendulum (m).

Small Angle Approximation

For small oscillations (\(\sin\theta \approx \theta\)), the equation simplifies to:

\(
\frac{d^2\theta}{dt^2} + \frac{g}{l}\theta = 0
\)

This is a second-order linear ODE, representing simple harmonic motion.

The Full Equation

For larger displacements, the approximation does not hold, and the full nonlinear equation must be solved. This requires numerical methods, as an exact solution is not possible for arbitrary initial conditions.

Want to dive deeper into ODEs?
👉 Introduction to Differential Equations in Mechanics

Step 1: Define the ODE in Python

Numerical solvers like SciPy require ODEs to be expressed as a system of first-order equations. The second-order pendulum equation:

\(
\frac{d^2\theta}{dt^2} + \frac{g}{l}\sin\theta = 0
\)

can be rewritten as two first-order equations:

\(
\frac{d\theta}{dt} = \omega
\) \( \frac{d\omega}{dt} = -\frac{g}{l} \sin\theta \)

Here, \(\omega\) is the angular velocity. This transformation allows us to solve the system numerically.

Defining the System in Python

import numpy as np

# Parameters of the pendulum
g = 9.81  # Gravity (m/s^2)
l = 1.0   # Length of the pendulum (m)

# Function defining the ODEs
def pendulum(t, y):
    theta, omega = y
    dtheta_dt = omega
    domega_dt = -(g / l) * np.sin(theta)
    return [dtheta_dt, domega_dt]

Step 2: Solve the ODE Numerically

We use scipy.integrate.solve_ivp to compute the solution. This function allows us to specify:

  • The time span of the simulation.
  • Initial conditions (angle and angular velocity).
  • The numerical method for solving the ODE.

Code for Solving the ODE

from scipy.integrate import solve_ivp

# Initial conditions: [initial angle, initial angular velocity]
theta0 = np.pi / 6  # 30 degrees in radians
omega0 = 0.0        # No initial velocity
y0 = [theta0, omega0]

# Time span for the simulation
t_span = (0, 10)  # 10 seconds
t_eval = np.linspace(t_span[0], t_span[1], 1000)

# Solve the ODE
solution = solve_ivp(pendulum, t_span, y0, t_eval=t_eval)

If you’re new to Python, start with this guide:
👉 Getting Started with Python for Engineers

Step 3: Visualize the Results

Visualization is a powerful way to understand the behavior of a system. Using Matplotlib, we can plot how the angle of the pendulum changes over time.

Code for Visualization

pythonCopyEditimport matplotlib.pyplot as plt

# Extract results
time = solution.t
theta = solution.y[0]

# Plot the results
plt.figure(figsize=(8, 4))
plt.plot(time, theta)
plt.title("Pendulum Motion")
plt.xlabel("Time (s)")
plt.ylabel("Angle (rad)")
plt.grid()
plt.show()

The resulting graph shows the oscillatory motion of the pendulum. Notice how the amplitude remains constant, reflecting energy conservation in the absence of damping.

Why Use Numerical Methods?

While analytical solutions provide elegance and insight, they are not always possible, especially for nonlinear systems. Numerical methods:

  • Handle complex equations like the full pendulum equation.
  • Allow for arbitrary initial conditions and parameters.
  • Provide a practical approach for real-world problems where exact solutions are unavailable.

These methods are foundational for simulations in engineering, physics, and applied mathematics.

Applications of Pendulum Dynamics

  1. Timekeeping: Pendulums revolutionized time measurement, leading to the development of accurate clocks.
  2. Seismology: Studying pendulum motion helps in understanding and designing instruments to measure ground movement.
  3. Education: The pendulum is a perfect example for teaching the principles of oscillatory motion and energy conservation.
  4. Engineering Systems: Pendulum-like dynamics appear in cranes, robotic arms, and even suspension bridges.

Key Takeaways

  1. Understanding the motion of a pendulum provides a foundation for exploring more complex systems.
  2. Numerical methods like solve_ivp allow us to solve equations that cannot be handled analytically.
  3. Python, with its rich ecosystem of libraries, is a powerful tool for engineers and scientists to analyze and visualize data.

Ready to Explore More?

This is just the beginning! Dive deeper into the fascinating world of numerical mechanics by exploring related topics:

Together, we’ll unlock the power of computational mechanics, one equation at a time.