Simple Pendulum

Parker-Sochacki: The Simple Pendulum



The next few blogs are applying the Parker-Sochacki method to a pendulum. I am starting with the Simple Pendulum and applied the simulation in Houdini. I will show that a higher order can handle larger time steps.

The simple pendulum is a nonlinear Ordinary Differential Equation but has no damping or driving forces acting on it. I will cover those next.

Here is the equation for the simple pendulum:simple_pendulum_2nd_order.PNG

This is a Second Order Ordinary Differential Equation where w_0 is a constant. For this example, w_0 is set to 1. To use Parker-Sochacki, this is first broken down into two first order ODEs.equations_simple_pendulum.PNG

Theta represents angular position in radians and omega represents angular velocity in radians per second. To use Parker-Sochacki, we make sin(theta(t)) an auxiliary variable. Doing this will require a derivative so we have another auxiliary variable for cos(theta(t)).

This requires the chain rule. Let a=sin(theta(t)). Taking the derivative of a leads to

a'(t) = cos(theta(t))*theta'(t)

We already established omega(t) as theta'(t). Here is our list of Ordinary Differential Equations:


Note that I switched a dot and b dot between equal signs. Having this setup makes it possible to turn theta and omega into a Maclaurin Series. Here is the recurrence relationship:


a and b require cauchy products. w does not because w_0^2 is a constant. You can set the initial conditions for theta and omega but note that it is in radians. a(0) = sin(theta(0)) and b(0) = cos(theta(0)).

The final form is a Maclaurin series. That Maclaurin series becomes the initial condition for the next iteration.

This render is my first use of HDRI images and I would like to thank Andrew Maynard for his help on that.

Everything was done in Houdini 15. This first simulation is at the default time step of 1/24 set to order 10 for Parker-Sochacki and took about an hour to render.

The second render is Parker-Sochacki at order 10 doubling the time step. The time step is now every other frame or 2/24:

The pendulum moves twice as fast.

Here is the point of this post: Higher order can handle larger time steps. I multiplied the time step by 10 so that now it is 10/24. Look at a simulation of the pendulum at order 2 and a time step of 10/24:

It quickly crashes. Look at a simulation of the pendulum at order 3 and a time step of 10/24:

It lasts longer than order 2 but it still crashes at taking such large time steps.

Now look at a  simulation of the pendulum at order 10 and a time step of 10/24:

It is stable for the entire 240 frames of the simulation. That is the power of using Parker-Sochacki at a higher order.