After two posts full of boring text and graphs I feel like I have to give away a crowd pleaser. Again, not a realtime simulation. About 120.000 particles, running at almost one FPS.
Wednesday, August 18, 2010
For a long time I've wondered how impulses develop inside a rigid body solver. Mainly for trying to come up with a prediction of where the final magnitude would be as early as possible. I dumped the data for all impulses out to a CSV file and graphed it for a couple of simple test scenarios. Let's start with a small test case, a box pyramid of two layers (three boxes). All dirty tricks are switched off, so this is a plain vanilla solver. The Y axis is impulse magnitude for the normal impulses and the X axis is iterations.
At a first glance, there seems to be no suprises here. Impulses form asymptotic curves towards their solution and there is not a lot of action after about twelve iterations. But look again, especially at the bright yellow curve, that seems to steadily decrease and form a non-horizontal asymptote, and what about the deviation around iteration eight and nine? Let's try a more complex exmple, with a box pyramid of five rows (15 boxes). There is around 90 contacts, and it would be too messy to show them all in the same graph, so I chose ten of them from the middle of the stack:
This is quite interesting - most of the curves are asymtotic to a non-horizontal line. Even at 64 iterations they have not found a stable configuration. Most interesting is maybe contact 45, which first goes up and then steadily decreases, while most othe contacts increase. Somehow, the solver is shifting around the weight. Maybe just because it can. Another graph we could study is the relative velocity along the normal direction (v from my previous post) and how it changes while iterating:
This graph is from the same scenario, but it's not from the same run, so contact numbering might not match exactly. Interesting here is that all curves steadily decrease after just a few iterations and seem asymptotic to a horisontal line, so our odd graph fellow above doesn't actually mean that the velocity is shifting similarly. Most likely it is an effect of an overdetermined system (with four co-planar contact points, we can shift the impulses around a bit and still get the same result. Not that the solver really has a good reason to do this, but it seems to be doing it anyway).
I also tried a variant of the solver that doesn't clamp accumulated normal impulses. Hence, the delta impulse is just clamped to be positive at each iteration. This prevents the accumulated impulse to ever have a negative inclination:
This is more what I expected. Each impulse is asymptotic to a horizontal line. Now if there only was a way to predict that asymptote and use it as initial guess, we would end up with a really stiff solver without having to keep track of contact points from the previous frame. I'm personally quite sceptical that it can be done from just looking at the curves, but it might be worth a shot. What i find most interesting is how good visual results you can get from just terminating the whole proces after three or four iterations, considering how far off the curves are.
Sunday, August 1, 2010
Following my last post about scientific papers not being written for engineers I will do an attempt at explaining a rigid body solver without equations:
Even though a rigid body scene may consists of hundreds of objects and thousands of contact points, a popular way to solve the problem is to solve each contact point in sequential order, one at a time. It sounds kind of lame, and compared to other methods it is, but if iterated a couple of times it gives really good results, and this is what most games are actually using, so let's focus on solving one contact without friction first:
So you have two objects and one contact point with a contact normal. Start by computing the velocity at the contact point for both objects and compute the difference between those vectors. Project that difference onto the contact normal (dot product). This is the contact's relative velocity along the contact normal and it indicates how much the objects are moving towards each other or away from each other at the contact point. Let's call this velocity v. If v is positive, the objects are moving away from each other and we're done. If v is negative we need to proceed onto computing and applying an impulse.
This is the key computation in the solver. The ultimate question we want to answer is - how big of an impulse do we need to apply to make the objects stop moving towards each other? The direction of the impulse is going to be the contact normal (since there is no friction yet), so we're only looking for the magnitude - a scalar quantity. The velocity v that we just computed is also a scalar quantity. The impulse magnitude will be proportinonal to v - the more the objects are approaching each other, the bigger impulse we need to apply. Easy! So what we're really looking for is the correlation between those two (another scalar quantity).
Let's assume we're applying an impulse of magnitude 1.0 (unit impulse) and see how big of a velocity change that would cause. The beauty of linear systems is that everything scales.. well.. linearly, so if we know how much of a velocity change a unit impulse causes we can just compare it to the desired velocity change, v, and adjust it accordingly. Say a unit impulse would cause a velocity change of 0.25 and v is -0.5, then we just apply two unit impulses (magnitude 2.0) and we'll reach our target relative normal velocity (zero). So make a copy of the linear and angular velocities for the two objects, apply an equal but opposite unit impulse and measure the relative contact velocity again following the exact same procedure as described above. Now that you know how much velocity change a unit impulse causes, the final computation boils down to a simple division. That's it folks. Apply the impulse you just computed on both objects in opposite directions and they will not move towards each other any more at the contact point.
Friction can be handled in the exact same way as described above, but in the direction of tangential relative velocity at the contact point. Friction impulses need to be capped to the magnitude of the normal impulse scaled by the friction coefficient, so that if the normal impulse is 2.0 and the friction coefficient is 0.9 your maximum friction impulse is 1.8. After you compute the friction impulse magnitude using above method, if it's more than 1.8 just apply 1.8.
Now to make the solver stable you need to go over all contacts several times, and there is also a method called "accumulated impulses" that improves accuracy a lot which is not covered here, but most importantly, you need to compensate for penetration. This is usually done in a very pragmatic way - if the objects are in penetration, adjust the target relative contact velocity so that it is not zero, but slightly positive (the more penetration the more positive). This means that after we're done solving, the objects will move slightly away from each other along the contact normal instead of not moving at all.