Stable Simulation: Avoiding Overshoot
Lecture, Dr. Lawlor
Almost all interesting simulations involve some sort of closed-loop
feedback, where a change in A causes a change in B which in turn
affects A again.
For example, consider the velocity of an object including drag
forces from air resistance.
F = -k * V
A = F/m (from F = mA)
dV / dt = A
Here F is the net force on the object, V is the object's velocity, k
the object's wind resistance (newtons of resistance per meter per
second of velocity), A is the object's acceleration, and m is the
If drag is the only force involved, and discretizing to first order,
dV = (-k*dt/m) * V
V = V + (-k*dt/m) * V = V + speed * V
where we've defined speed = (-k*dt/m).
That is, the change in V itself varies with V. If that speed
(-k*dt/m) is positive, then V will feed back on itself leading to
exponential growth. Luckily, speed is always negative because
real drag constants, timesteps, and masses are all positive.
the theoretical non-discretized equations, a negative speed term
V to drop exponentially, though possibly with a managably small
exponent. A small negative speed term causes the discrete V to
also drop exponentially.
But weird stuff starts happening if the speed term gets negative
Anytime speed<-2.0, V is going to spiral out of control fairly
quickly. This is similar to an inexperienced driver skidding
ice: they steer too hard, overcorrecting the car's orientation, and
up making the skid worse instead of better. The effect is that
the simulation "blows up", where the points have vanished off the
screen, or hit floating-point infinity or not-a-number. Not
- At speed=0, V doesn't change at all.
- At speed=-0.1, V shrinks by 10% per timestep.
- At speed=-0.5, V shrinks by 50% per timestep.
- At speed=-1.0, V = 0 after the first timestep.
- At speed=-1.5, V overshoots zero every timestep, ending each
step 50% shorter but facing the other way.
- At speed=-2.0, V overshoots zero completely, and ends up
facing exactly the opposite direction with the same magnitude.
- At speed=-2.1, V overshoots zero and gains 10% per timestep!
Note that speed = -k*dt/m, so speed is going to get out of control
a stability limitation on the timestep. If your steps get
than this, you can expect the simulation to explode.
- k, the drag coefficient, gets too big. It'd be weird
you could simulate air, but not the slower movements expected in
or tar, but it's quite possible.
- m, the object mass, gets too small. The drag coefficient
for small objects should probably get small as well.
- dt, the timestep, gets too big. This is common, but
luckily we can
adjust the simulation timestep easily, as long as our computer
Stability in Spring Systems
Try out this spring system stability
Generally, if the spring forces can push the mass to complete a
cycle in less than pi timesteps, the mass starts overshooting and
the simulation blows up. Here, a simple
harmonic oscillator has period 2*pi*sqrt(m/k). I've seen
programs that compute their own timestep by figuring out the minimum
stable period at startup.
- If you decrease the mass (m), the simulation blows up.
- If you increase the spring constant (k), the simulation blows
up. (Both of these make the mass-spring system move
- If you increase the timestep, the simulation blows up.
Note also that if you update P and V in the wrong order, the system
is unconditionally unstable. You can see this by updating cell
B5 to "=B4+A$2*C4" (using the *old* velocity, instead of the new
updated velocity), and choosing "fill down". Note it now
always blows up. Using the *new* velocity (C5) makes this "Leapfrog
integration"; using the old velocity (C4) makes it the less
accurate and unstable Euler.
The Speed of Sound
The other big limit on the timestep is the speed of sound.
If I take a steel bar, and break it into 100 small springs, the
only way a force exerted on one end will eventually move the other
end is via each of the 100 springs. This must take at least
100 timesteps, since each spring computes its forces based on the
*old* position of each neighbor.
The problem is the speed of sound in solid materials is very
fast; for example, the speed
of sound in steel is about 6 km/sec (just under orbital
velocity!). Thus, sound crosses a 1m steel bar in 1/6,000 of
a second. To simulate 100 smaller springs, we need 600,000
physics updates per second, which is usually way too many.
One saving grace for computer graphics is that it's not easy to
see a slower-than-correct speed of sound. And in 3D, using larger
(hence faster-reacting) tets inside an object can help reduce this
effect. But the net result is that accurate physics
typically takes at least several physics steps per animation