In the previous post, we have covered a technique to resolve constraint in our physics simulation. Because we are working in the real-time domain, we need to sacrifice accuracy for speed by using an iterative method to resolve our constraints. Sequential Impulse doesn’t give the exact solution for the problems, but the solution it offers is good enough for games where accuracy is not the priority, as long as it still feels believable. Today, we will learn and implement two techniques that are widely used to help Sequential Impulse and to improve the stability of your physics.

Warm Starting

There is one property of games that we can exploit to help sequential impulse compute a better result. In games, the property of objects like velocity, and position doesn’t change much from one frame to the next frame. This property is called “frame coherence”. If the state of objects doesn’t differ much, then the solution for the constraints should not differ much either. So, instead of recomputing the impulse every frame, we can use the solution from the previous frame as the starting point and refine it instead. This technique is called Warm Starting, Below is a source code I copy from my physics engine. In the previous post, we keep the total impulse from one iteration to another. For warm starting to work we need this value, so we need to keep this value across the frame (In the sample code, nIAcc for contact impulse and tIAcc for friction impulse). We apply this total impulse to each body before the resolution logic takes place.

void ftConstraintSolver::warmStart()
{
    for (int32 i = 0; i < m_constraintGroup.nConstraint; ++i)
    {
        ftContactConstraint *constraint = &(m_constraintGroup.constraints[i]);

        int32 bodyIDA = constraint->bodyIDA;
        int32 bodyIDB = constraint->bodyIDB;
        ftVector2 normal = constraint->normal;
        ftVector2 tangent = normal.tangent();

        ftManifold *manifold = &(constraint->contact->manifold);

        for (int32 j = 0; j < constraint->numContactPoint; ++j)
        {

            ftContactPointConstraint *pointConstraint = &(constraint->pointConstraint[j]);

            pointConstraint->nIAcc = manifold->contactPoints[j].nIAcc;
            pointConstraint->tIAcc = manifold->contactPoints[j].tIAcc;

            ftVector2 impulse = pointConstraint->nIAcc * normal;
            impulse += pointConstraint->tIAcc * tangent;

            m_constraintGroup.velocities[bodyIDA] -= constraint->invMassA * impulse;
            m_constraintGroup.velocities[bodyIDB] += constraint->invMassB * impulse;

            ftVector2 r1 = pointConstraint->r1;
            ftVector2 r2 = pointConstraint->r2;
            m_constraintGroup.angularVelocities[bodyIDA] -= constraint->invMomentA * r1.cross(impulse);
            m_constraintGroup.angularVelocities[bodyIDB] += constraint->invMomentB * r2.cross(impulse);
        }
    }
}

As previously mentioned, Sequential Impulse is an iterative method, the more iteration it does, the more accurate the solution will be. Warm Starting works as if we do more iteration by continuing the previous iteration.

Slop

Baumgarte Stabilization is a hacky technique that we apply to fix positional drift. With this technique even if there is only a little penetration, our physics will add energy to the system and push the objects apart. Imagine if we have a stack of boxes, The lowermost boxes will penetrate the surface because of gravity and the boxes above pushing it down. The solver will resolve this and push it up, creating penetration between the boxes. This penetration between boxes in return will push the box back again to the surface creating a cycle of small upward and downward movement. This jittery effect could be avoided by using slop. Slop allows a small penetration or positional drift to exist. So, instead of resolving every penetration. we only resolve penetration that exceeds the slop value. The slop value should be small enough that it is visually insignificant but large enough that it is numerically significant. Box2d use 0.005f (or 5cm).

$b=–\frac{\beta }{△t} \left(max\left(d – slop, 0\right)\right)$

real bias = ftMax(manifold->penetrationDepth[i] - m_option.linearSlop, 0);
pointConstraint->positionBias = m_option.baumgarteCoef * bias;

Beside using slop at contact, we can use this technique on another resolution like joint and restitution.

Final Result: