Mathcad Library

Differential Algebraic Equations: Pendulum

by Inutech, Inc.

This file is part of an electronic book that demonstrates advanced Differential Equation solving techniques in Mathcad. If this topic is of use to you, download the whole collection .

A Differential Algebraic Equation (DAE) is an ordinary differential equation whose solution is subjected to additional equality constraints, not just initial or boundary conditions. A simple example is that of the motion of a pendulum. Typically, this problem is recast in polar coordinates, because of the difficulty of including the algebraic constraint on the values of motion in x and y in rectangular coordinates.

Simple Pendulum in Rectangular Coordinates

Mass of the pendulum's mass point:

Gravitational acceleration:

Length of the pendulum:

Final time for integration:

In addition to equations for Newton's Law (force = mass x acceleration), there is a constraint that these positions must sum to the length of the pendulum:

Gives the x component of the force.

Gives the y-component of the force, and the additional downward force due to gravity acting on the pendulum.

The pendulum begins at the bottom of its trajectory, with no downward velocity and a horizontal velocity of 1. Note that these initial conditions must satisfy the algebraic constraint equation.

The period of this pendulum, in the small angle approximation, is

Plot the solution

The pendulum swings back and forth in space between the extremes of its x displacement, which is related to the fixed period and the length of the bob. Note the slight variation in the maximum values of y caused by the numerical simulation.

Error for the solution to the DAE is given by how closely the algebraic condition, in this case, the equation for x and y position in terms of the pendulum length, is met.

Index of a DAE

Mathcad's DAE algorithm solves DAE's up to index 3. Roughly speaking, the index of a DAE is the number of differentiations needed until we get a system consisting only of ordinary differential equations with no algebraic constraints. In the case above, we'll need to get an equation for the change in force over time, without the length constraint. If we consider the plain pendulum equations above, differentiating the algebraic constraint with respect to time yields

Differentiation 1

Differentiation 2

Substituting in the original expressions for acceleration in x and y, reduces to

equation that results after differentiating the algebraic constraint twice.

Observe how the calculated value for y grows further and further from the value required to meet the algebraic constraint. This is a result of the index reduction (index 1 here).

One more integration will produce a system of ordinary differential equations — three differentiations to reach this point means that the original DAE is of index 3. This procedure can be used to reduce the index of a DAE so that it may be solved by Mathcad, recognizing that the original algebraic constraints are still valid and must be fulfilled by any solution, as well as by any initial value. In our example of the plain pendulum, cast with an Index 2 formulation, the constraint would now be

Even though the constraint

is no longer visible, it must be satisfied by initial values, which must also fulfill the new constraint equation.

Index Reduction and Drift Effect

Since it is possible to reduce a system DAEs to a system of ODEs by index reduction, why consider the algebraic equations explicitly? One reason is that the numerical solution drifts away from the algebraic constraint proportional to the amount of index reduction (as demonstrated with the index 1 solution shown). If we reduce the index by one, the error increases linearly; if we reduce the index by two, it increases quadratically. The complete transformation to a system of ODEs could destroy the algebraic constraint completely. The goal is to only reduce a DAE system to the point where it is tractable: in Mathcad, systems of index 3 or smaller will be solvable.

Reference

Hairer E., Wanner G. (1991): Solving Ordinary Differential Equations II. Stiff and Differential-Algebraic Problems, Springer Series Computational Mathematics, Vol. 14, Springer