Real-Time Simulation of Rigid Body Dynamics
Final Project: CSE 557
Charles Gordon
Overview
A great deal of work has been done on creating realistic simulations of rigid body dynamics at frame rates. This
project attempts to combine the results of a number of different, state of the art algorithms for modelling realistic
physics and collisions. The goal was to create realistic approximations to real-world physics while maintaining at
least the 30 frames/second needed for animation purposes.
References
Chris Hecker has done a series of articles for Game Developer magazine and put up a wonderful
web page detailing his results and the library-full
of references he used.
One of those references is to David Baraff's set of papers on
physicaly based modelling,
which were used as a tutorial for SIGGRAPH '95. These include a simple algorithm for collision detection
and response which turned out to be monstrously difficult to implement (see my notes below).
The collision detection algorithm used in my program is the same one used for
Q-COLLIDE. This algorithm combines some
interesting results from computation geometry with some (fairly difficult to understand) algorithms for
finding closest points on arbitrary convex polymeshes. Kelvin Leung's master's thesis on QCOLLIDE can be found
here.
Details
The most common cause of aliasing problems in physics modelling is the numerical integrator. The easiest
way to approximate the closed form solutions of classical physics is to use an Euler integration step, which is
really just a first order Taylor polynomial. In order to increase the accuracy of the integrator and to also
increase the size of the time step allowable, I used an adaptive midpoint approximation, which is a second order
Taylor approximation that detects when the error is beyond some tolerance and subdivides the time step (this is a
surprisingly simple operation). Where the Euler integrator required time steps on the order of hundredths of a
second to be anywhere near accurate, the midpoint evaluator uses time steps right around the 1/30th of a second
necessary for real-time animation.
The program simulates both linear and angular dynamics. Quaternions are used to specify the orientation of the
bodies in space (and also because they are much easier to handle with the integrator). Initial conditions can
be set and the bodies follow all of the normal constraints of gravity. Special objects I called "attractors" can
be placed in scenes to act like small planetary bodies, giving out an attractive force that is inversely proportional
to the distance squared and directly proportional to the mass of the object. These can be set to any multiple of
the gravitational constant so that they either act as strong gravity wells or, if the scalar is negative, as
strongly repulsive forces.
Collision detection turned out to be an immensely more difficult problem than I had anticipated. This is
especially true of trying to make it accurate in real-time, as the calculations involved can be very
time consuming. In order to reduce the number of calculations there must be a simple, high granularity step
using bounding boxes or spheres. When the bounding boxes of two objects collide there is another step to
determine if the objects themselves are actually in contact. This is followed by a determination of exactly
which points on each object are in contact and then the calculation of an impulse force to be applied there.
The algorithms I used for each of these steps is summarized below.
Bounding Boxes: The first step in the CD algorithm is to put axis-aligned bounding boxes around all of
the objects in the scene. To avoid the overhead of determining the maximum extent of each object in three-space,
each object also had an object-aligned bounding box that underwent the same rotations and translations as the object.
The axis aligned bounding box was determined from the eight points of this object aligned box. To determine
collisions, the axis-aligned boxes were intially sorted along each axis. At every time step the coordinate lists
were traversed and resorted using a simple insertion sort (so that only objects that had actually changed position
were considered). Any objects that had overlapping boxes were added to a list of objects to be checked in the
next phase of collision detection. The bounding box algorithm required a small overhead at every time step to
keep the object-aligned boxes in the correct orientation, followed by an expected linear time to re-sort the axis
aligned boxes.
Separating Vectors: The separating plane algorithm used by Baraff can take up to O(mn) time to execute,
where m and n are the number of vertices of the colliding polymeshes. A better algorithm using a separating
vector is used by the QCOLLIDE library and takes time O(m+n). By caching the vector (or plane) between time
steps this can be reduced to nearly constant time. In this program a separating vector is found between polymeshes
everytime their bounding boxes overlap. This vector can be found very quickly by using incremental updates to
its position based on the results of a set of dot products at the last position. It is detailed in the thesis
mentioned above. Unfortunately, due to time constraints, I was unable to implement the caching algorithm, so
this step always takes O(m+n) time every time bounding boxes are found to intersect.
Closest Features: This is by far the most difficult part of collision detection and it took the bulk of my
time (this is an immense understatement, I have over 30 hours into this). The first approach I tried was an
O(mn) algorithm that simply compared the distance between every pair of features on both polymeshes until it found
the closest pair. This approach is neither effective nor fast. It could only find that, say, a point intersected
a face, but not where on the face they intersected. This led to very unrealistic looking collisions between objects
with large faces (such as squares). Although it worked well for spheres with a large number of small faces. It also
ran very slowly, bringing the simulation to its knees. The next approach I tried used the results from the previous
step along with an internal graph of neighboring vertices for each vertex in the polymesh. This allowed a much
faster lookup of the colliding features, but did not fix the problem of not knowing exactly which points were
colliding. It was at this point that I came across a description of the Gilber-Johnson-Keerthi algorithm for
collision detection. It is explained in the master's thesis mentioned above as well as in:
Bergen, Gino Van Den. A Fast and Robust GJK Implementation for Collision Detection of Convex Objects. Department of Mathematics and Computer Science, Eindhoven University of Technology.
This algorithm does an amazing job of finding the closest points between any two arbitray polymeshes in expected
linear time (using temporal coherence and caching). The basic idea is to do a pairwise subtraction of all the points
in both polymeshes and then find the points that are closest to the origin. Rather than actually performing the
pairwise subtraction, they devised a method for incrementally descending to the correct pair of points starting from
the convex hull of any four points on the surface. The Johnson algorithm, which is used as a subalgorithm, finds
the closest point from these four points to the origin. This point is used to find a new point on the original
surface and the process is repeated. Empirical results in the paper (and repeated in my program) show that the
algorithm generally terminates after 4 to 6 of these iterations. With some selected use of caching this can be
reduced even further.
The results of using this algorithm in my project leave a lot to be desired. Besides the mathematical difficulty
of the method, very little was said in either of the papers I found about parameter tuning, which is the most
difficult part to making it work correctly. The GJK algorithm includes a number of small error thresholds that
must be tuned very precisely to insure that the algorithm terminates correctly. After 15 hours of watching it
go into infinite loops I can assure you that this is no triviality. I am fairly confident that my current
implementation is correct, but it does not interface with the collision response algorithm very well, which means
that objects will occasionally get "glued together" because the physics engine can't generate the correct force
to get them unstuck. This is a fairly difficult problem, but not unsolvable. In the meantime I have fudged a
small repulsive force to guarantee that they do not stick. This makes the result less than realistic, but it
was all I had time for.