Tom Lee

Software development geekery.

Jakobsen’s Advanced Character Physics (Annotated)

This is a reproduction of Thomas Jakobsen’s Advanced Character Physics, which I’m slowly annotating with my own notes.

I think it was after reading this post (by the guy who wrote the completely awesome Soldat) maybe 5-10 years back that I first took a look at Jakobsen’s paper. I struggled to absorb much first time through. It’s relatively simple as far as physics papers go, I guess, but my math & physics skills are more than a little wanting. :)

So now, years later, I’m taking another, more leisurely pass at trying to absorb this stuff.

Please note this is a work in progress & it’s a long way from being complete.

You can find my notes indented, in italics.

Advanced Character Physics

by Thomas Jakobsen
(Annotated by Tom Lee)


1. Introduction

The use of physically-based modeling to produce nice-looking animation has been considered for some time and many of the existing techniques are fairly sophisticated. Different approaches have been proposed in the literature [Baraff, Mirtich, Witkin, and others] and much effort has been put into the construction of algorithms that are accurate and reliable. Actually, precise simulation methods for physics and dynamics have been known for quite some time from engineering. However, for games and interactive use, accuracy is really not the primary concern (although it’s certainly nice to have) – rather, here the important goals are believability (the programmer can cheat as much as he wants if the player still feels immersed) and speed of execution (only a certain time per frame will be allocated to the physics engine). In the case of physics simulation, the word believability also covers stability; a method is no good if objects seem to drift through obstacles or vibrate when they should be lying still, or if cloth particles tend to “blow up”.

When building games, we care more about immersion and performance than accuracy: we can blatantly fake the “physics” so long as it *feels* believable to the user. In a game physics simulation, “believable” implies “stable”: a physics simulation is not believable if stupid shit happens (i.e. if it is “unstable”).

The methods demonstrated in this paper were created in an attempt to reach these goals. The algorithms were developed and implemented by the author for use in IO Interactive’s computer game Hitman: Codename 47, and have all been integrated in IO’s in-house game engine Glacier. The methods proved to be quite simple to implement (compared to other schemes at least) and have high performance.

The techniques in Jakobsen’s paper can be used to build performant, stable physics systems. They have been used in published AAA game titles & used in IO Interactive’s in-house game engine.

The algorithm is iterative such that, from a certain point, it can be stopped at any time. This gives us a very useful time/accuracy trade-off: If a small source of inaccuracy is accepted, the code can be allowed to run faster; this error margin can even be adjusted adaptively at run-time. In some cases, the method is as much as an order of magnitude faster than other existing methods. It also handles both collision and resting contact in the same framework and nicely copes with stacked boxes and other situations that stress a physics engine.

This system has a bunch of cool properties I only partially understand. It gracefully deals with several situations that tax other physics engines.

In overview, the success of the method comes from the right combination of several techniques that all benefit from each other:

  • A so-called Verlet integration scheme.
  • Handling collisions and penetrations by projection.
  • A simple constraint solver using relaxation.
  • A nice square root approximation that gives a solid speed-up.
  • Modeling rigid bodies as particles with constraints.
  • An optimized collision engine with the ability to calculate penetration depths.

The algorithm described in this paper is a combination of several simpler, mutually beneficial techniques.

Each of the above subjects will be explained shortly. In writing this document, the author has tried to make it accessible to the widest possible audience without losing vital information necessary for implementation. This means that technical mathematical explanations and notions are kept to a minimum if not crucial to understanding the subject. The goal is demonstrating the possibility of implementing quite advanced and stable physics simulations without dealing with loads of mathematical intricacies.

Jakobsen did his best to keep this paper accessible.

I started feeling dumb about two paragraphs back, making this a Document of -1 Self Worth.

The content is organized as follows. First, in Section 2, a “velocity-less” representation of a particle system will be described. It has several advantages, stability most notably and the fact that constraints are simple to implement.

Okay, this is important: given that he’s talking about a *physics* system, when Jakobsen refers to a “particle system” you should understand that he’s not talking about the kind of particle system often associated with computer graphics. Rather, he’s talking about a system of particles in the “physics” sense. The phrase system of particles at a minimum seems to yield more relevant stuff from a quick google.

It seems you can substitute “system of particles” for “particle system” for the rest of the paper & it will still make perfect sense. :)

Section 3 describes how collision handling takes place. Then, in Section 4, the particle system is extended with constraints allowing us to model cloth.

Collision handling: what should our physics system do when an object runs into something?

When he talks about “model[ing] cloth”, he’s literally talking about simulating a piece (or pieces!) of cloth. Still not sure what I mean? Check out this great little cloth simulation in JavaScript. If you lack a physics brain, a quick skim of the text on that page is very helpful. The explanation of constraints with the visual aid is really intuitive.

Section 5 explains how to set up a suitably constrained particle system in order to emulate a rigid body.

A rigid body is just a solid object.

Next, in Section 6, it is demonstrated how to further extend the system to allow articulated bodies (that is, systems of interconnected rigid bodies with angular and other constraints). Section 7 contains various notes and shares some experience on implementing friction etc. Finally, in Section 8 a brief conclusion.

This should otherwise all be pretty self explanatory with the exception of maybe “articulated bodies”, but he gives a pretty good definition I think.

In the following, bold typeface indicates vectors. Vector components are indexed by using subscript, i.e.,     ![$\begin{displaymath}\mathbf{x} = ( x{1}, x{2}, x_{3} )\end{displaymath}$][9]


2. Verlet integration

The heart of the simulation is a particle system. Typically, in implementations of particle systems, each particle has two main variables: Its position x and its velocity v. Then in the time-stepping loop, the new position ![$\mathbf{x{\prime}}$][10] and velocity ![$\mathbf{v{\prime}}$][11] are often computed by applying the rules

    ![ $$ \begin{displaymath} \begin{array}{lcl} \mathbf{x{\prime}} = \mathbf{x} + \mathbf{v} \cdot \Delta t \ \mathbf{v{\prime}} = \mathbf{v} + \mathbf{a} \cdot \Delta t \end{array} \end{displaymath} $$ ][12]

where ![$\Delta t$][13] is the time step, and ![$a$][14] is the acceleration computed using Newton’s law ![$f = ma$][15] (where ![$f$][16] is the accumulated force acting on the particle). This is simple Euler integration.

Here Jakobsen is stating the typical approach people tend to take when simulating movement in a game engine: the new position of a particle (![$\mathbf{x{\prime}}$][10]) is calculated by its current position (![$\mathbf{x}$][17]) plus its current velocity (![$\mathbf{v}$][18]) multiplied by the amount of time elapsed since the last physics “tick” (![$\Delta t$][13]).

Similarly, the new velocity of the same particle (![$\mathbf{v{\prime}}$][11]) is calculated by its current velocity (![$\mathbf{v}$][18]) plus its acceleration (![$\mathbf{a}$][19]) multiplied by the amount of time elapsed (![$\Delta t$][13]).

If you’re more familiar with the C programming language, this might look a little more familiar:

GameObject::update(float elapsed)
  position  = velocity * elapsed;
  velocity  = acceleration * elapsed;

Pretty standard stuff so far.

Here, however, we choose a velocity-less representation and another integration scheme: Instead of storing each particle’s position and velocity, we store its current position ![$\mathbf{x}$][17] and its previous position ![$\mathbf{x{*}}$][20]. Keeping the time step fixed, the update rule (or integration step) is then

    ![ $$ \begin{displaymath} \begin{array}{lcl} \mathbf{x{\prime}} = 2\mathbf{x} - \mathbf{x{*}} + \mathbf{a} \cdot \Delta t{2} \ \mathbf{x{*}} = \mathbf{x} \end{array} \end{displaymath} $$ ][21]

This is called Verlet integration (see [Verlet]) and is used intensely when simulating molecular dynamics. It is quite stable since the velocity is implicitly given and consequently it is harder for velocity and position to come out of sync. (As a side note, the well-known demo effect for creating ripples in water uses a similar approach.)

And bam, just like that we have Verlet integration: a way to express the movement of a particle without explicit velocity. Let’s look at this a little closer:

![$\mathbf{x{\prime}} = 2\mathbf{x} - \mathbf{x{*}} + \mathbf{a} \cdot \Delta t{2}$][22]

You’ll recognize most of the variables here from the Euler integration stuff we saw just before, with the exception of ![$\mathbf{x{*}}$][20] which is the particle’s previous position. Jakobsen says this explicitly, but you can see it expressed on the second line too:

![$\mathbf{x{*}} = \mathbf{x}$][23]

With this, ![$\mathbf{x{*}}$][20] is ready for the next physics tick.

Here’s what it might look like in C source code:

GameObject::update(float elapsed)
  const Vector3D current_position = position;
  position  = position - last_position   acceleration * (elapsed * elapsed);
  last_position = current_position;

Note he also implies that a fixed time-step should be used. I don’t fully understand why, but I recall hearing elsewhere that this is important.

TODO more reading RE: fixed time step verlet.

Lastly, note that Jakobsen claims that this approach is “quite stable” because “it is harder for velocity and position to come out of sync”. Given his earlier rant about stability, I’m going to go ahead and assume that he means “more stable than Euler integration” since he doesn’t go into more detail. Perhaps it’s obvious to mathlier geeks than I. :)

It works due to the fact that ![$2 \boldsymbol{\mathbf{x}} - \boldsymbol{\mathbf{x{*}}} = \boldsymbol{\mathbf{x}} + (\boldsymbol{\mathbf{x}} - \boldsymbol{\mathbf{x{*}}})$][24] and ![$\boldsymbol{\mathbf{x}}-\boldsymbol{\mathbf{x{*}}}$][25] is an approximation of the current velocity (actually, it’s the distance traveled last time step).

Verlet integration works because we can calculate an approximation of our current velocity from our previous and current position.

For example, if our previous position was ![$\boldsymbol{\mathbf{x{*}}} = (0, 1, 0)$][26] and our current position is ![$\boldsymbol{\mathbf{x}} = (1, 1, 0)$][27], we can estimate our instantaneous velocity by ![$\boldsymbol{\mathbf{x}} - \boldsymbol{\mathbf{x{*}}} = (1, 1, 0) - (0, 1, 0) = (1, 0, 0)$][28]. Thus, we can use ![$(1, 0, 0)$][29] as an approximation of the velocity for our current tick.

What Jakobsen is trying to say here is if you look carefully at the Verlet equations again, you’ll see that the implicit velocity is there in plain sight:

    ![ $$ \begin{displaymath} \begin{array}{lcl} \boldsymbol{\mathbf{x{\prime}}} & = & 2\boldsymbol{\mathbf{x}} - \boldsymbol{\mathbf{x{*}}} + \boldsymbol{\mathbf{a}} \cdot \Delta t{2} \ & = & \boldsymbol{\mathbf{x}} + (\boldsymbol{\mathbf{x}} - \boldsymbol{\mathbf{x{*}}}) + \boldsymbol{\mathbf{a}} \cdot \Delta t{2} \end{array} \end{displaymath} $$ ][30]

So using Verlet integration, our new position ![$\boldsymbol{\mathbf{x{\prime}}}$][31] is our current position ![$\boldsymbol{\mathbf{x}}$][32] plus an approximation of our velocity ![$\boldsymbol{\mathbf{x}} - \boldsymbol{\mathbf{x{*}}}$][33] plus some acceleration ![$\boldsymbol{\mathbf{a}}$][34] multiplied by the squared time step ![$\Delta t{2}$][35].

It is not always very accurate (energy might leave the system, i.e., dissipate) but it’s fast and stable.

Reading between the lines: weird shit should happen less.

By changing the update rule to something like ![$\boldsymbol{\mathbf{x{\prime}}} = 1.99\boldsymbol{\mathbf{x}}-0.99\boldsymbol{\mathbf{x{*}}+\mathbf{a{*}}} \cdot \Delta t{2}$][36] a small amount of drag can also be introduced to the system.

XXX ??? “Viscous Drag”

Note that ![$\boldsymbol{\mathbf{a{*}}}$][37] is introduced here.

At the end of each step, for each particle the current position ![$\boldsymbol{\mathbf{x}}$][32] gets stored in the corresponding variable ![$\boldsymbol{\mathbf{x{*}}}$][38]. Note that when manipulating many particles, a useful optimization is possible by simply swapping array pointers.

TODO Sounds like this might save you a vector copy or two per particle, but need to hack out an implementation to understand fully what he’s on about here.

The resulting code would look something like this (the Vector3 class should contain the
appropriate member functions and overloaded operators for manipulation of vectors):

// Sample code for physics simulation
class ParticleSystem {
  Vector3 m_x[NUM_PARTICLES];    // Current positions
  Vector3 m_oldx[NUM_PARTICLES]; // Previous positions
  Vector3 m_a[NUM_PARTICLES];    // Force accumulators
  Vector3 m_vGravity;            // Gravity
  float  m_fTimeStep;

  void TimeStep();

  void Verlet();
  void SatisfyConstraints();
  void AccumulateForces();

// (constructors, initialization etc. omitted)

// Verlet integration step
ParticleSystem::Verlet() {
  for(int i=0; i Several existing systems have a hard time dealing with collisions. When collisions do occur, nasty hacks are often involved. Further, handling them can be slow.

Here, we use yet another strategy. Offending points are simply projected out of the
obstacle. By projection, loosely speaking, we mean moving the point as little as possible
until it is free of the obstacle. Normally, this means moving the point perpendicularly out
towards the collision surface.

Let’s examine an example. Assume that our world is the inside of the cube (0,0,0)-
(1000,1000,1000) and assume also that the particles’ restitution coefficient is zero (that is,
particles do not bounce off surfaces when colliding). To keep all positions inside the valid
interval, the corresponding projection code would be:

// Implements particles in a box
void ParticleSystem::SatisfyConstraints() {
  for(int i=0; i