Integration error in high velocity

Question

I've implemented a simple simulation of two planets (simple 2D disks really) in which the only force is gravity and there is also collision detection/response (collisions are completely elastic). I can launch one planet into orbit of the other just fine.

The collision detection code though does not work so well. I noticed that when one planet hits the other in a free fall it speeds backward and goes much higher than its original position. Some poking around convinced me that the simplistic Euler integration is causing the error.

Consider this case. One object has a mass of 1kg and the other has a mass equal to earth. Say the object is 10 meters above ground. Assume that our dt (delta t) is 1 second. The object goes to the height of 9 meters at the end of the first iteration, 7 at the end of the second, 4 at the end of the third and 0 at the end of the fourth iteration.

At this points it hits the ground and bounces back with the speed of 10 meters per second. The problem is with dt=1, on the first iteration it bounces back to a height of 10. It takes several more steps to make the object change its course.

So my question is, what integration method can I use which fixes this problem. Should I split dt to smaller pieces when velocity is high? Or should I use another method altogether? What method do you suggest?

EDIT: You can see the source code here at github:https://github.com/elektito/diskworld/

Answer 1

Say the object is 10 meters above ground. Assume that our dt (delta t) is 1 second. The object goes to the height of 9 meters at the end of the first iteration

Here lies your problem. It is true that the velocity at the end of the first iteration is \$1 m.s^{-1}\$. However during that time the object has not travelled \$1 m\$.

In fact, since the acceleration is constant, the average object velocity is simply \$0.5 m.s^{-1}\$ and thus the object has only travelled \$0.5 m\$. The next frames will be:

$$ \begin{array}{c|c|c|l} \text{time} & \;\; v \;\; & v(avg) & \text{height} \\ \hline 0 & 0 & - & \;\;10.0 \\ 1 & 1 & 0.5 & \;\;9.5 \\ 2 & 2 & 1.5 & \;\;8.0 \\ 3 & 3 & 2.5 & \;\;5.5 \\ 4 & 4 & 3.5 & \;\;2.0 \\ 5 & 5 & 4.5 & -2.5 \text{(collision happened)} \end{array} $$

So, to know the new object position, use the average velocity instead of the new velocity. It will greatly improve your accuracy (in fact, in the case of constant acceleration, it will even give you exact results).

At this points it hits the ground and bounces back with the speed of 10 meters per second.

This is incorrect, too. Where does the value 10 come from? With your integration method, the velocity should be 4.

Answer 2

Just guessing from the description of your issue, since I don't have any actual code, I would also look in to the following if you haven't already:

What happens if you use a smaller time step, say 1/30 sec?

Check that you apply the time step variable (dt) in the appropriate parameters of your simulation calculation. Since you mentioned that the "combination causes the problem" and the physics are mangled with the collisions it might be necessary to apply the time step to both physics and collision detection/response code.

Make sure that the objects don't overlap during the collision. So the collision detection code will not run multiple times that may add the collision response force multiple times.

Make sure that physics calculations occur in fixed time step because if you use variable time step the integration produces small errors proportional to each time step that sum up during multiple frames.

Try to smooth the results of the physics calculations by averaging the results of the last two or three frames. That way you will also average the small errors and the end result will be perceived as more accurate.