Hi All,

Been implementing a 2d sequential impulse engine (mostly for the academics of it). The main pipeline of my current engine is:

1) Detect Collisions

2) Resolve % of each penetration

3) Repeat 1 & 2 'n' times

3) Resolve velocities

4) Integrate

5) Repeat

This seems to give fairly descent results. However, no matter what I try when objects stack in vertex-edge configurations, or edge-edge in pyramid-like structure, then the angled objects will push others out of the way regardless of how many objects are in the pile being pushed. Here are some videos showing the results:

Initially I thought there was something wrong with my friction implementation but I've tried increasing friction to no avail. Next, for edge-edge cases I only resolve a single contact, so I've experimented with different contact points also to no avail.

I can only think of two things which may be causing this, either:

1) This is un-avoidable for a sequential impulse solver -- which doesn't feel entirely correct. sliding from sequential solvers should be gradual drifts; not this dramatic but perhaps I am wrong.

2) Steps 1 & 2 in my loop place a large emphasis on resolving penetrations (via position corrections) rather than giving more weight to #3. As I've been reading around online, it seems many don't even do #2 but rather only resolve the velocity by applying the appropriate impulse plus some small portion of the penetration (baumgaurte stabilization?).

As far as my penetration resolution implementation: it resolves the penetrating objects in the direction of the contact normal according to the objects' mass. Each pair is resolved in the order of greatest separating velocity.

Simple answer seems to be put everything to sleep.... though I feel like that's just masking a stability problem.

I'd appreciate any feedback on the ordering of my main loop, and any suggestions on what may be causing my stacks to push apart / how I might be able to mitigate that? Would a TOI solver fix this, or is simultaneous resolution the only answer? Anyways, I feel like I have lots of ideas, but no clear direction; so some guidance would be very helpful.

Let me know if you need more info from me.

Thanks

Are you warmstarting the friction? You can get pretty bad friction otherwise.

I'm not. Any recommended ways for doing this? Or would simply feeding in the previous frame's friction impulse work?

Here's my current friction implementation:

Geometry::Vector2D tangent;					
tangent.x = relVel->x - normal->x  * relativeVelOnNormal;
tangent.y = relVel->y - normal->y  * relativeVelOnNormal;
Geometry::Normalize(tangent, tangent);
//tangent.x = -normal->y;
//tangent.y =  normal->x;

flttype planarVel = (relVel->x * tangent.x + relVel->y * tangent.y); 
flttype changeInPlanarVel = planarVel;
flttype relAPT = (posRelA->x * tangent.y - posRelA->y * tangent.x);
flttype relBPT = (posRelB->x * tangent.y - posRelB->y * tangent.x);

flttype relACrossTangent = entityAInertiaTensorInv * relAPT;
flttype relBCrossTangent = entityBInertiaTensorInv * relBPT;
					
flttype angularMassAT = relACrossTangent * relAPT;
flttype angularMassBT = relBCrossTangent * relBPT;
flttype frictionImpulse = changeInPlanarVel / (curContact->totalLinearMass + angularMassAT + angularMassBT);


// clip friction mag to min (impulseMag * frictionCoefficient, frictionMag)
flttype frictionMag;
if (frictionImpulse < curContact->friction * impulseMag)
{
    // static friction
    frictionMag = -frictionImpulse;
}
else
{
    // dynamic friction
    frictionMag = -curContact->friction * frictionImpulse;
}	


/////////////////////////////
// 
//   Update velocities
//
//

// linear vel
entityA->velocity.x -= (impulseX + frictionMag * tangent.x) * entityAInvMass;
entityA->velocity.y -= (impulseY + frictionMag * tangent.y) * entityAInvMass;

entityB->velocity.x += (impulseX + frictionMag * tangent.x) * entityBInvMass;
entityB->velocity.y += (impulseY + frictionMag * tangent.y) * entityBInvMass;

// angular vel
entityA->angularVelocity -= (relACrossNormal * impulseMag + relACrossTangent * frictionMag);
entityB->angularVelocity += (relBCrossNormal * impulseMag + relBCrossTangent * frictionMag);

First I would make sure that your tangent doesn't go to zero when the relative velocity vanishes. What is the impulseMag? Is this the incremental or the accumulated impulse? You need to clamp against the accumulated impulse.

I would clamp the friction impulse like this (assuming impulseMag > 0) :

frictionImpulse = clamp( frictionImpulse, -curContact->friction * impulseMag, curContact->friction * impulseMag );

Clamping the friction as you mentioned seems to have helped a bit, though at first there is definitely this odd jarring to right of the pyramid before it falls. Once it falls though, the objects do settle down pretty quickly without pushing apart as much as before, though that could just be due to the configuration of how the blocks lie once they fall:

impulseMag is supposed to be the total impulse required to cause the appropriate change in velocity. I assume it is the accumulated value because I solve for this impulse in one go rather than in several small iterations (per time step). Or does accumulated mean something else? Here is how I compute impulseMag:

flttype relativeVelOnNormal = (relVel->x * normal->x + relVel->y * normal->y);

if (relativeVelOnNormal > 0.0)
      continue;

// get applied velocity from previous step along normal (i.e. gravity)		
flttype velFromAcc = (entityB->appliedVel.x - entityA->appliedVel.x) * normal->x +
	             (entityB->appliedVel.y - entityA->appliedVel.y) * normal->y;

// determine desired change in velocity (w/out applied vel from prev step)
flttype changeInVel = -velFromAcc + (-1.0 - curContact->coefficientOfRestitution) * (relativeVelOnNormal - velFromAcc);

// determine impulse required to cause desired change in vel
flttype relAPN = (posRelA->x * normal->y - posRelA->y * normal->x);
flttype relBPN = (posRelB->x * normal->y - posRelB->y * normal->x);

flttype relACrossNormal = entityAInertiaTensorInv * relAPN;
flttype relBCrossNormal = entityBInertiaTensorInv * relBPN;					
flttype angularMassA = relACrossNormal * relAPN;					
flttype angularMassB = relBCrossNormal * relBPN;

flttype totalLinearMass = entityAInvMass + entityBInvMass;
flttype impulseMag = changeInVel / (totalLinearMass + angularMassA + angularMassB);

// keep impulse positive
impulseMag = impulseMag < 0.0 ? 0.0 : impulseMag;

1) Detect Collisions

2) Resolve % of each penetration

3) Repeat 1 & 2 'n' times

You don't re-detect collisions. You iterate only on the current contact points. I recommend looking at Box2D Lite to see how this works. Erin also published a bunch of papers that should help you.

Isn't it possible though that by resolving one contact pair, new contacts or collisions are created?

No, you only change velocities. New contacts would be only created if you also update the position/orientation after each impulse. That's not how it works

Ah, okay that makes sense. I was thinking in this manner: fix overlap (via position corrections) then resolve velocity (via vel corrections or impulse). But from what you're saying it sounds like I should not fix overlap directly via position corrections, but only through impulses. Is that correct?

Yes! You first solve the velocities. The projection steo give you the 'best' velocities to advance to the positions. Then you correct positions. I would not re-run the collision detection though, but keep contact points local and measure the penetration along the normal.

Again, I recommend having look at Box2D and Box2D Lite.