Collision Response

[goruka]

Hi! I purchased and rest most of the book "Physics for Game Developers" publisherd by O'Reilly.
By now I've read most of it and I think i get a good hang on implementing rigid bodies..

However, this complex algorithm which is used for collision response (bullet seems to use a smimilar formula) I can't get it to work nor figure it out:

Code: Select all

  pt1 = Collisions[i].vCollisionPoint - Bodies[b1].vPosition;

  pt2 = Collisions[i].vCollisionPoint - Bodies[b2].vPosition;

				

  // calculate impulse


  j = (-(1+fCr) * (Collisions[i].vRelativeVelocity*Collisions[i].vCollisionNormal)) /
       ( (1/Bodies[b1].fMass + 1/Bodies[b2].fMass) +

       (Collisions[i].vCollisionNormal * ( ( (pt1 ^ Collisions[i].vCollisionNormal)*Bodies[b1].mIeInverse )^pt1) ) +

       (Collisions[i].vCollisionNormal * ( ( (pt2 ^ Collisions[i].vCollisionNormal)*Bodies[b2].mIeInverse )^pt2) )

					);

  Vrt = Collisions[i].vRelativeVelocity * Collisions[i].vCollisionTangent;

	

  if(fabs(Vrt) > 0.0 && dofriction) {

	Bodies[b1].vVelocity += ( (j * Collisions[i].vCollisionNormal) + ((mu * j) * Collisions[i].vCollisionTangent) ) / Bodies[b1].fMass;										

	Bodies[b1].vAngularVelocityGlobal += (pt1 ^ ((j * Collisions[i].vCollisionNormal) + ((mu * j) * Collisions[i].vCollisionTangent)))*Bodies[b1].mIeInverse;


        //separate b2 in the opposite impulse
     } else {

             same thing but "mu" (friction coefficient) seems considered zero
     }

I'm sure this snippet looks familiar to those who implemented it. So Far I can guess:

pt1 and pt2 are the points in both rigid bodies' local coordinates (but not rotation).

vRelativeVelocity is both velocities at the collision point (with linear and angular processed), and then the difference, as in vel_obj1 - vel_obj2

but my doubt is.. what is exactly the collision normal vCollisionNormal ? is it a vector from the collision point, pointing orthogonally from the surfaces towards object 1 ?

Also, is there a place where this algorithm is properly explained? or should i just use it "as is" ?

[Erwin Coumans]

I'm sure this snippet looks familiar to those who implemented it.
but my doubt is.. what is exactly the collision normal vCollisionNormal ? is it a vector from the collision point, pointing orthogonally from the surfaces towards object 1 ?

It is the separating normal, pointing from one object to the other indeed. Most iterative constraint solvers use this approach, although it might look obfuscated. See Bullet collision response source code line 89 resolveSingleCollision, or slighly optimized in resolveSingleCollisionCombinedCacheFriendly line 166 of btSequentialImpulseConstraintSolver.cpp. In Open Dynamics Engine (ODE, quickstep) is looks a bit different,see SOR_LCP around line 500, but it is essentially the same as sequentially applying impulses when you analyse it further.

Also, is there a place where this algorithm is properly explained? or should i just use it "as is" ?

Check out equation 9 in Chris Heckers Collision Response article.

Hope this helps,
Erwin

[Dirk Gregorius]

I like this paper: http://www.essentialmath.com/AngularImpulse.pdf

Note that in the final equation 'n' ( the contact normal) is usually unit length so you can drop the n*n term in the denominator in this case

[goruka]

Thanks enormously for the replies! i've read the linked pdf and the source code of both bullet and ode and I think I have a much better understanding now. I think my solver is working much better now!

It is the separating normal, pointing from one object to the other indeed.

However, I still have a little doubt, what should I actually use as the impact normal? I am testing collision+response with spheres, however, the impact normal between them is always pointing to the center of mass of the spheres, thus colliding them does not produce impulses that change angular velocity. Is this to be expected? Does it mean I should be adding friction to the equation for angular velocity to happen?

Thanks!

[Dirk Gregorius]

For spheres you have only central collisions, basically the same like point masses. For convex polyhedra or general meshes you get also angular impulses.

[goruka]

For spheres you have only central collisions, basically the same like point masses. For convex polyhedra or general meshes you get also angular impulses.

I see, but Imagine I have a sphere rolling over a plane (billard ball) or two rubber spheres colliding.. I imagine they'll be getting angular impulses due to friction.. how do I add that to the model?

[Dirk Gregorius]

You can compute the friction impulse in the tangent direction as you did for the normal. Just switch the normal with the tangent direction and also use the relative tangent velocity. Then you test whether the imoulse lies in the friction cone and clamp if not.

To find the tangent direction at the contact point you can use:

v_rel = v2 + cross(omega2, r2) - v1 - cross(omega1, r1)
vn_rel = dot(v_rel, n) * n
vt_rel = v_rel - vn_rel
t = vt_rel / length(vt_rel)

To test for the friction cone:

if ( jt > mu * jn )
jt = mu * jn

[goruka]

You can compute the friction impulse in the tangent direction as you did for the normal. Just switch the normal with the tangent direction and also use the relative tangent velocity. Then you test whether the imoulse lies in the friction cone and clamp if not.

This worked great so far, thanks!!! Now the spheres in my test roll against each other and on the floor perfectly!

The only issue left i seem have so far is that sometimes (in some rare cases) they keep rolling very slow against an axis (not | or -- , more like \ -respective to the floor-) even though they seem to be standing still in position. I could understand -- but \ seems rare.. is this to be expected?
To clarify a little more, it's like it has no linear velocity but a small bit of angular velocity..

Again thanks a lot!

[Dirk Gregorius]

I don' understand? Are they slidinging to the side?

[goruka]

I don' understand? Are they slidinging to the side?

I managed to make an animated gif showing what happens.

This is a first test, without using Z (2D-like), which seems to work ok:
(EDIT: I forgot to mention, there's invisible floor and invisible walls on left and right)



And this second test, i added a little Z displacement to each ball. Notice the ball rolling strangely towards the camera, also other balls keep rolling without moving, but not in a parallel axis to the floor:



I'm not sure if this is to be expected, I think i implemented friction ok, as overall they seem to work, but I feel i'm missing something..

Actual code is:

Code: Select all

		Vector3 velt_rel = vel_rel - n * vel_rel.dot(n);
		float velt_rel_len = velt_rel.length();
		
		if (j>0 && velt_rel_len > CMP_EPSILON) {
		
			Vector3 t = velt_rel/velt_rel_len;
			
			float numt = velt_rel_len;
			float denomt = t.dot( t * ( inv_mass + c.against->get_inv_mass() ) );
			
			inertia_local = Vector3::cross( contact_point, t );
			inv_inertia_tensor.transform( inertia_local );
				
			inertia_foreign = Vector3::cross( contact_point_foreign, t );
			c.against->get_inv_inertia_tensor().transform( inertia_foreign );
			
			denomt += (Vector3::cross( inertia_local, contact_point ) + Vector3::cross( inertia_foreign, contact_point_foreign )).dot( t );
			
			float jt=numt/denomt;
			
			float mu =1.0;
			
			if ( jt > mu * j )
				jt = mu * j;
	
			if ( jt < -mu * j )
				jt = -mu * j;
	
			impulse += friction_impulse;
		}

          linear_velocity+= impulse * inv_mass;
          Vector3 angvel_rel=Vector3::cross( contact_point, impulse );
          inv_inertia_tensor.transform(angvel_rel);
          angular_velocity+=angvel_rel

[Dirk Gregorius]

Do you only apply positive incremental impulses? You can try the following and accumulate your impulse and clamp against the accumulated impulse. So instead of some test whether the relative velocity in the direction of the contact normal is negative (approaching) you always compute it (actually allowing attracting impulse) and then clamp like this:

float jn = "Normal impulse strength"

float jn_old = mAccumN;
mAccumN = Max( jn_old + jn, 0.0f )
jn = mAccumN - jn_old;

Now apply jn...


Then before you start the solve your system you choose two arbitrary orthorgonal axes in the contact plane for each contact point and accumulate the friction in these directions. Basically the same as you do for the normal. You need the "basis" for the accumulation. Then you clamp the accumulated friction against the accumulated normal impulse.


float jt1 = "First friction impulse"

float jt1_old = mAccumF1;
mAccumF1 = Clamp( jt1_old + jt1, -mu * mAccumN, mu * mAccumN );
jt1 = mAccumF1 - jt1_old;

Now apply jt1...



Repeat the same for the second friction direction.


You can improve the friction if you choose the first friction direction t1 in the direction of the relative tangent velocity. For now reset the accumulated impulses each frame.

For a good explanation of this see the presentations of Erin Catto here: http://www.gphysics.com


Cheers,
-Dirk

[goruka]

Hey, so far this worked great, thanks enormously for all the help!
It's kind of frustrating overall to work on physics simulation as i find there is a lot of information around, but i can't seem to find either good tutorials or somehow something that tells me how to integrate everything together, how each technique works, what performance/precision ratio to expect, etc. And I don't feel like spending a lot of time trying each only to find out it's limitations later.. so again, thanks for your help on this!

[Dirk Gregorius]

Personally I would recommend learning from Box2D. Start with the Lite version and port it to 3D - just sphere, box and ball-socket joint with a velocity only (Baumgarte) and a two phase velocity/position solver (see Erin's latest post). Then look into the newer release and extend your version to general convex polyhedra, all joint types and different friction models (per point vs central). Then you look into Bullet and implement GJK/EPA and an incremental manifold and compare it against the SAT version. Finally go continuous, add bouancy and try to couple particle with rigid body physics.

Should be some fun...