Particle dynamics simulation and conservation of energy

Question

Suppose you want to simulate the dynamics of a particle subject to Newton's laws of motion. We could just fix a discretization of the time dt and solve approximate differential equations like

    y+=vy*dt;
    vy+=ay*dt;

and so on for the other coordinates. If we are unsatisfied with the accuracy we could make dt smaller or use higher order approximation methods, at the cost of slowing the simulation and this would reduce the error, but the error can still sum up and become big in after a long time.

Now my point is: suppose that you want to run the simulation for long time so you don't mind short run accuracy but you want the system to preserve long run properties such as conservation laws.

Concrete example: suppose the particle is subject to gravity with the acceleration ay constantly equal to g. We have the conservation law 1/2 m v^2 + m g h =constant.

Now after the approximation step I could compute the energy and discover it is slightly changed, and this lead to a degeneration after long times (for example it can lead to overflow errors or particle slowing down). At each step I could then adjust v and h in order to make the energy conserved, but I have a lot of choice here... how do I decide for the better in this constrained choice?

In general what is the best way to deal with situation where you can sacrifice short time accuracy to save speed and preserve the long time behavior?

Answer 1

You're mentioning a very good point! The simple answer is: In some cases simple Euler integration is not good enough, and no increase in accuracy will help.

Clearly, you have to distinguish the goal of your particular simulation. If you're simulating a whole bunch particles bouncing, interacting with each other and/or fields and gravity en masse, then the Euler method should be good enough. Your system might diverge from the "real" solution, but for all practical purposes no one will notice.

That being said, there are cut cases where one intuitively knows the solution or some property of it. If you have a plane orbiting a star, it should do so on an elliptical orbit. If you have a weight attached to a spring, the elongation of the spring should be sinusoidal (without dampening). Then simple Euler integration will not just produce some inaccuracy, it will produce a solution that completely diverges from the closed path. The spring will get ridiculously elongated overtime, a planet will be escape the solar system...

When decreasing the timestep of an Euler integration, it will take longer time for the divergence to become noticeable, but if the system is prone to it, it will always occur.

There are other integrations scheme that inherently perform much better in conservative systems. Most notable are the Verlet integration and the Runge-Kutta methods. Verlet is much easier to understand and already quite good, especially in systems with conservation of energy. The Runge-Kutta method is usually the best in accuracy, but it is quite complicated and it has some inherit dampening which most of the time helps it to converge, but you might find a case where it performs worse.

Here is a nice comparison I just found: http://codeflow.org/entries/2010/aug/28/integration-by-example-euler-vs-verlet-vs-runge-kutta/

Sorry, for not providing implementations to compare. You should be able to find plenty when searching for these keywords.

Answer 2

Pretty sure you can just use symplectic euler integration and call it good. Over time numerical accuracies should be insignificant enough it won't matter.

Often times in a physics simulation you'll have a simple way of measuring error, and add energy back into the system. However with just simple particle dynamics conserving energy like this usually isn't necessary. Usually particle dynamics is only done to animate tons and tons of particles in an efficient manner, whereas accuracy isn't really the goal.

Answer 3

This

where you can sacrifice short time accuracy to save speed and preserve the long time behavior?

from the request is nonsensical as it stands. As you extend the simulation, and particularly if you do not carefully guard and preserve the accuracy of every result, you will end up with the rounding error of your calculations accumulating to overwhelm everything else.

The only way to preserve long-term accuracy is to ensure even greater accuracy at every shorter interval.