I couldn't find anything that really explained what accumulated impulses are, at first i tried to store the impulse of each contact inside the contact body its-self (without modifying the velocity of the object directly) and then applying this accumulated velocity in one big go. That made everything explode, spectacularly!
I think i'm onto something with my new iteration, i changed my physics loop from
- Find intersection
- Apply Impulses
- Integrate velocity & positoon
To
- Find intersections
- Integrate velocity
- Apply Impulses
- Integrate postiion
Which seems to have really helped. Is this change in update loop what's referred to as accumulated impulse, or am i doing something unrelated?
I also changed my ApplyImpulse function based on your feedback. I could not get the coefficient of restitution working with this new correction, so it's gone for now. If anyone is interested, the updated code looks like this:
void ApplyImpulse(RigidbodyVolume& A, RigidbodyVolume& B, const CollisionManifold& M, int c, float deltaTime) {
float invMass1 = A.InvMass();
float invMass2 = B.InvMass();
float invMassSum = invMass1 + invMass2;
if (invMassSum == 0.0f) {
return; // Both objects have infinate mass!
}
vec3 r1 = M.contacts[c] - A.position;
vec3 r2 = M.contacts[c] - B.position;
mat4 i1 = A.InvTensor();
mat4 i2 = B.InvTensor();
vec3 relativeVel = (B.velocity + Cross(B.angVel, r2)) - (A.velocity + Cross(A.angVel, r1));
vec3 relativeNorm = M.normal;
Normalize(relativeNorm);
// Moving away from each other? Do nothing!
/*if (Dot(relativeVel, relativeNorm) > 0.0f) {
return;
}*/
float bias = 0.2f;
float slop = 0.01f;
float correction = (bias / deltaTime) * fmaxf(M.depth - slop, 0.0f);
float e = fminf(A.cor, B.cor);
float numerator = (/*-(1.0f + e) **/ (-Dot(relativeVel, relativeNorm) + correction));
float d1 = invMassSum;
vec3 d2 = Cross(MultiplyVector(Cross(r1, relativeNorm), i1), r1);
vec3 d3 = Cross(MultiplyVector(Cross(r2, relativeNorm), i2), r2);
float denominator = d1 + Dot(relativeNorm, d2 + d3);
float j = (denominator == 0.0f) ? 0.0f : numerator / denominator;
/*if (M.contacts.size() > 0.0f && j != 0.0f) {
j /= (float)M.contacts.size();
}*/
//j = fmaxf(j, 0.0f);
vec3 impulse = relativeNorm * j;
A.velocity = A.velocity - impulse * invMass1;
B.velocity = B.velocity + impulse * invMass2;
A.angVel = A.angVel - MultiplyVector(Cross(r1, impulse), i1);
B.angVel = B.angVel + MultiplyVector(Cross(r2, impulse), i2);
// Friction
vec3 t = relativeVel - (relativeNorm * Dot(relativeVel, relativeNorm));
if (CMP(MagnitudeSq(t), 0.0f)) {
return;
}
Normalize(t);
numerator = -Dot(relativeVel, t);
d1 = invMassSum;
d2 = Cross(MultiplyVector(Cross(r1, t), i1), r1);
d3 = Cross(MultiplyVector(Cross(r2, t), i2), r2);
denominator = d1 + Dot(t, d2 + d3);
float jt = (denominator == 0.0f) ? 0.0f : numerator / denominator;
if (M.contacts.size() > 0.0f && jt != 0.0f) {
jt /= (float)M.contacts.size();
}
if (CMP(jt, 0.0f)) {
return;
}
float friction = sqrtf(A.friction * B.friction);
if (jt > j * friction) {
jt = j * friction;
}
else if (jt < -j * friction) {
jt = -j * friction;
}
vec3 tangentImpuse = t * jt;
A.velocity = A.velocity - tangentImpuse * invMass1;
B.velocity = B.velocity + tangentImpuse * invMass2;
A.angVel = A.angVel - MultiplyVector(Cross(r1, tangentImpuse), i1);
B.angVel = B.angVel + MultiplyVector(Cross(r2, tangentImpuse), i2);
}
At this point i've taken out all linear projection, and everything is running with the changes mentioned above. The simulation is pretty close. For reference,
This is a gif of a scene i build in unity: http://giphy.com/gifs/d1E2BSEJUEEz0dQ4
And this is my simulation: http://giphy.com/gifs/l0Ex6twOL4qhbX4Ig
I have a feeling the "crawling" i'm getting is due to either the removal of the coefficient restitution or maybe i'm doing friction wrong?
Thank you for taking the time to help!
