Copyright © Michael Richmond.
This work is licensed under a Creative Commons License.
Physicists face many complicated situations which can only be solved by numerical methods; one example is the motion of more than two bodies due to their mutual gravitational attraction. As you have seen, numerical methods are not foolproof -- in fact, they can easily lead to incorrect, even absurd, results. If you know the correct answer ahead of time (as in the case of a skydiver falling without air resistance), then you can recognize errors as "wrong." But what if you don't know the right answer; how can you be sure that your calculations aren't being led astray by roundoff or truncation errors?
Fortunately, in a number of important cases (such as bodies in space moving due to gravitational forces), there is a way to check that one's calculations haven't gone awry: use integrals of motion, also known as conserved quantities. A truly conserved quantity should always have the same value, at the beginning, middle, and end of a long calculation. By monitoring one of the integrals of motion, one can detect that errors are creeping into the simulation, or gain some confidence that they are not.
Just what are these "integrals of motion?" Here are a few examples:
Figuring out which quantities ought to be conserved in a particular situation might be considered the Physics part of the problem; using those quantities as checks on the numerical calculations could be considered the Math/computer part of the problem. One must be both a physicist and a computer scientist to do the best job.
Warning: preserving the integrals of motion is a necessary condition for numerical simulations -- but is not sufficient to prove that they are correct. Imagine a simulation of the solar system which calculates the initial value of energy, say, and then leaves the planets fixed in place as time moves forward. Clearly, that isn't an accurate depiction of the true behavior ... but it will yield exactly the same energy at all times.
Suppose that your job is to calculate the motion of planets in the solar system. You are given the position and velocity of each planet at a starting time, and decide to use Euler's method to calculate positions and velocities at intervals one day. Okay, sounds like a good plan.
But first, before you start the integration, you calculate the kinetic energy of each body (the sum of which will be positive) and the gravitational potential energy between each pair of bodies (the sum of which will be negative). You add these together to find the initial total mechanical energy:
Initial Energy = Sum of KE + Sum of GPE
Next, you use Euler's method to calculate the motions of the planets for ten years, with a timestep of one day. At the end of the entire simulation, you once again add up the kinetic and gravitational potential energies:
Final Energy = Sum of KE + Sum of GPE
The Final Energy should be exactly the same as the Initial Energy. One measure of the quality of your simulation is the fractional change in Energy:
Final Energy - Initial Energy
fractional change in Energy = -----------------------------
Initial Energy
If this fractional change is small enough (say, 0.01 percent), then your results are probably accurate. On the other hand, if the fractional change is significant (say, 3 percent), then your calculations are introducing significant errors, and the results are very likely to be wrong.
What should you do if the fractional change in this conserved quantity is too large? You have two choices:
The first choice is easy: cut the time step in half, run from the beginning, and once again calculate the change in Energy. If the change has shrunk to an acceptably small level, you're done. If not, you can either decrease the time step further, or switch to a different method of integration.
Copyright © Michael Richmond.
This work is licensed under a Creative Commons License.