The
first time I implemented a spring model, it fascinated me for hours. The
effect itself is amazingly realistic, and its implementation is fairly simple.
Fortunately, I found a lot of articles, references, and source code to help
with my research. Nevertheless, as I went further down the road, I noticed
that in most cases these references limited to the standard applications
of spring models -- string, cloth and jelly.

This article
reviews the implementation of a spring model from its simplest form to
more sophisticated applications, taking the subject a step beyond the
material available in most references.

**Spring Basics**

Before modeling
springs with the computer, you need to understand the basic principles
of springs from your classic physics book.

As you compress
or extend a spring, it creates a force that is opposed to direction fo
the force that you are applying. This force is mathematically equated
by the formula:

**F
= -k*Dx**

Dx
= xc-xi

Where F is the resultant force, k is the spring coefficient, and Dx
the distance from its inertial position (xc = current distance, and xi
= distance at the inertial position).

Mathematically,
you can think of a spring as being two points separated from a distance
x. This is its inertial position. In other words, if you don't move these
points, no force will be generated because Dx
= 0. If you try to compress the distance between these points the Dx
will be negative, generating a positive force (expansion). If you try
to separate these points the Dx
will be positive, generating a negative force (compression). The picture
below illustrates these cases.

Pa
-/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\- Pb (Dx
= 0, F = 0)

Pa -/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\- Pb
(Dx
< 0, F > 0)

Pa -/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\/\-
Pb (Dx
> 0, F < 0)

**Springs
(using "ascii art")**

The k coefficient
is called the elasticity coefficient, and weights the spring's final force.
In other words, a spring with bigger k is stiffer because it creates a
larger force from its inertial state. Conversely, a spring with a smaller
k is more flexible because it creates a smaller force from its inertial
state.

**Spring Force Implementation**

Implement
the spring force is a direct application of the first two equations. In
other words, first compute the distance of the spring's extremities (Pa
and Pb). This is done through Pythagoras' theorem. Once you have the distance,
you subtract from the inertial distance getting the Dx.
Then, multiply this scalar value to the k coefficient and finally use
the inverse of this value to compute the force.

As simple
as it sounds, the implementation in C code has its details. First, the
distance Dx,
and the k coefficient are scalars, but the force itself is a vector. Therefore,
it requires a direction, which can be from Pa to Pb, or from Pb to Pb.
As of right now, the code below always use the direction Pb to Pa, or
mathematically Fdir = Pb-Pa.

As the force's
direction is normalized, you can use the distance you already computed
for the Dx
for the normalization. You must be careful is to avoid division by zero
when normalizing the vector. This is a nasty condition and it usually
does not happen, but you have to do something about it in case it does.
You have two options in this particular case: assign a random force to
the points, or null out the force. Theoretically, the force should be
big enough that as the points get closer they will never occupy the same
space. Therefore, it seems to be a good option to randomly assign a massive
force in this case. However, this can cause undesirable artifacts, because
some points inside the spring model might fly a part with extreme velocity
when compared with others. Although nulling out seems wrong, it doesn't
cause any major problems and preserves the spring's stability. The code
below shows the implementation; it's a little redundant for sake of clarity.
I null out the force in case of potential division by zero.

**Listing
1: Computing A Spring's Force**

#define SPRING_TOLERANCE
(0.0000000005f)

void SpringComputeSingleForce
(VECTOR *force, VECTOR *pa, VECTOR *pb, float k, float delta_distance)

{

VECTOR force_dir;

float intensity;

float distance;

float delta;

//get
the distance

float dx = pb->x-pa->x;

float dy = pb->y-pa->y;

float dz = pb->z-pa->z;

distance
= (float) sqrt ((dx*dx)+(dy*dy)+(dz*dz));

if
(distance < SPRING_TOLERANCE)

{

force->x =

force->y =

force->z = 0.0f;

return;

}

force_dir.x
= dx;

force_dir.y = dy;

force_dir.z = dz;

//normalize

force_dir.x /= distance;

force_dir.y /= distance;

force_dir.z /= distance;

//force

delta = distance-delta_distance;

intensity = k*delta;

//store

force->x = force_dir.x*intensity;

force->y = force_dir.y*intensity;

force->z = force_dir.z*intensity;

}