Here is a common, yet designer-unfriendly, implementation of numeric springing. Let x_i and v_i denote the numeric value and velocity, respectively, at frame i; and let k, d, h, x_t, denote the spring stiffness, damping factor, time step between frames, and target value, respectively.

v_{i+1} = d * v_i + h * k * (x_t - x_i)

x_{i+1} = x_i + h * v_{i+1}

Let's dissect this implementation.

The difference between the target value and the current value (x_t - x_i) contributes to the velocity change after multiplied by the spring stiffness k. The larger the difference and spring stiffness, the larger the change in velocity.

The new velocity v_{i+1} inherits the value from the old velocity v_i scaled by the damping factor d, which should be a value between zero and one. If the damping factor is equal to one, then no scaling happens, and the new velocity completely inherits the value from the old velocity. If the damping factor is a fraction of one, then the magnitude of oscillation will decrease gradually.

That seems all reasonable and good, but just what exactly are the values for k and d do we need to achieve a specific springing effect, say an oscillation at 5Hz frequency and the oscillation magnitude decreases 90% every one second? Many people just decide to bite the bullet and spend a large amount of time hand-tweaking the values for k and d, until the result is acceptably close to the desired effect.

Now let's take a look at the differential equation for a damped spring system centered around x = x_t.

x'' + 2 * zeta * omega * x' + omega^2 * (x - x_t) = 0,

I won't go into details as to how this equation is derived. All you really need to know is that omega and governs how the system behaves and what they mean.

Omega is the angular frequency of the oscillation. An angular frequency of 2 pi (radians per second) means the oscillation completes one full period over one second, i.e. 1Hz.

Zeta is the damping ratio. A damping ratio of zero means there is no damping at all, and the oscillation just continues indefinitely. A damping ratio between zero and one means the spring system is underdamped; oscillation happens, and the magnitude of oscillation decreases exponentially over time. A damping ratio of 1 signifies a critically damped system, where this is the point the system stops showing oscillation, but only converging to the target value exponentially. Any damping ratio above 1 means the system is overdamped, and the effect of springing becomes more draggy as the damping ratio increases.

Here's a figure I borrowed form Erin Catto's presentation on soft constraints, to show you the comparison of undamped, underdamped, critically damped, and overdamped systems.

Below are the equations for simulating a damped spring system using the implicit Euler method.

x_{i+1} = x_i + h * v_{i+1}

v_{i+1} = v_i - h * (2 * zeta * omega v_{i+1}) - h omega^2 * (x_{i + 1} - x_t)

Solving x_{i+1} and v_{i+1} in terms of x_i and v_i using Cramer's rule, we get:

x_{i+1} = Delta_x / Delta

v_{i+1} = Delta_v / Delta

where:

Delta = (1 + 2 * h * zeta * omega) + h^2 * omega^2

Delta_x = (1 + 2 h * zeta * omega) x_i + h * v_i + h^2 * omega^2 * x_t

Delta_v = v_i + h * omega^2 * (x_t - x_i)

Below is a sample implementation in C++. The variables x and v are initialized once and then passed into the function by reference every frame, where the function keeps updating their values every time it's called.

Now we have our formula and implementation all worked out, what parameters should we expose to the designers from our spring system so that they can easily tweak the system?

Inspired by the aforementioned presentation by Erin Catto, I propose exposing the oscillation frequency f in Hz, and the fraction of oscillation magnitude reduced p_d over a specific duration t_d due to damping.

Mapping f to omega is easy, you just multiply it by 2 pi:

omega = 2 * pi * f

Mapping p_d and t_d to zeta is a little bit more tricky.

First you need to understand that the oscillation magnitude decreases exponentially with this curve:

y(t) = e^{-zeta * omega * t}

After plugging in p_d and t_d, we get:

y(t_d) = e^{-zeta * omega * t_d} = p_d

So we can solve for zeta:

zeta = ln(p_d)} / (-omega * t_d)

For example, if we want the oscillation to reduce 90% (p_d = 0.1) every 0.5 second (t_d = 0.5) with oscillation frequency 2Hz (omega = 4 pi), then we get:

zeta = ln(0.1)} / (-4 * pi * 0.5) = 0.367

Numeric springing is a very powerful tool for procedural animation, and now you know how to precisely control it. Given the desired oscillation frequency and percentage of oscillation magnitude reduction over a specific duration, you can now compute the exact angular frequency omega and damping ratio zeta necessary to create the desired springing effect.

By the way, let's look at the animated example shown above again:

This animation was created with omega = 8 pi (4 full oscillation periods per second) and zeta = 0.23 (oscillation magnitude reduces 99.7% every second).