Capsule Moment Of Inertia

Started by
8 comments, last by JohnnyCode 9 years, 6 months ago

Does anybody know a reference for the formula for the inertia tensor of a capped cylinder / capsule?

I always find the same equation but I actual want to understand the formula and not just use it, cause everybody does.

So the important part for me in the quote below:


  Ia = M1*(REAL(0.25)*radius*radius + (REAL(1.0)/REAL(12.0))*M1*length*length) +
    M2*(REAL(0.4)*radius*radius + REAL(0.375)*radius*length + REAL(0.25)*length*length);


I understood the following:

r = radius
l = length of cylinder

I = mass_of_cylinder * (1/4 * r * r + 1/12 * l* l)     //formula for cylinder inertia 
         +
    mass_of_sphere * (2/5 * r * r)  //formula for sphere inertia
         +
    mass_of_sphere * (3/8 * r * l + 1/4 * l * l) //what's this for?

I can only assume it maybe has to do something with the distance to the center , so the inertia will grow.

Thanks for your help smile.png


This is how it's done in Open Dynamics Engine (ODE), and it work well for me too =D


void dMassSetCappedCylinder (dMass *m, dReal density, int direction,
			     dReal radius, dReal length)
{
  dReal M1,M2,Ia,Ib;
  dAASSERT (m);
  dUASSERT (direction >= 1 && direction <= 3,"bad direction number");
  dMassSetZero (m);
  M1 = M_PI*radius*radius*length*density;			// cylinder mass
  M2 = (REAL(4.0)/REAL(3.0))*M_PI*radius*radius*radius*density;	// total cap mass
  m->mass = M1+M2;
  Ia = M1*(REAL(0.25)*radius*radius + (REAL(1.0)/REAL(12.0))*M1*length*length) +
    M2*(REAL(0.4)*radius*radius + REAL(0.375)*radius*length + REAL(0.25)*length*length);
  Ib = (M1*REAL(0.5) + M2*REAL(0.4))*radius*radius;
  m->_I(0,0) = Ia;
  m->_I(1,1) = Ia;
  m->_I(2,2) = Ia;
  m->_I(direction-1,direction-1) = Ib;

# ifndef dNODEBUG
  checkMass (m);
# endif
}
Advertisement

Physically, the moment of inertia is basically a rotational form of regular translational inertia. Newton's 2nd law (\(F = ma\)) describes how masses are accelerated due to forces, but mass can be seen as a resistance to acceleration. In that way, more mass is more resistant to being accelerated. For translations, \(F=ma\) describes this. For rotations, \(T = J \alpha\) is the equivalent.

The inertia tensor (matrix) describes in essence how "hard" any object is to be rotated about some point and some axis. To calculate it, you choose a reference point (can be on or off the body, but we like using either the centroid or center of mass), choose reference axes (can be arbitrary but should be a right-handed basis at least), break the objects up into small, infinitesimal pieces, and then calculate the moment of inertia of each piece, then sum them up (usually in a kind of weighted average fashion). This is why integrals are usually present in the matrix formulation. The inertia matrix is useful when you are rotating in 3D and not necessarily along an axis of the reference frame from which you calculated the matrix (i.e. not rotating about the X, Y, or Z axis).

However, a lot of times the objects of which you want the moment of inertia are easily broken down to more simple geometric forms that can be easily calculated. In your capsule example, they are breaking down the capsule into cylinders and spheres, something easily more calculable than using the integral method. We can use do a similar thing as above by calculating the moments of inertia of each geometric piece and summing them up with a kind of weighted average. If we are summing up the moments of inertia in this way, our "weighted average" is the parallel axis theorem. The parallel axis theorem allows us to calculate the moment of inertia of an object and compensate for the distance that object is away from the rotation axis. We can see this by spinning a pencil about its middle or about its ends. It's "harder" to spin about the ends than about the middle because the bulk of the pencil's mass (or we can also say the center of mass) is much farther away from the rotation axis. The parallel axis theorem gives us the correct adjustment factor to sum things up this way.

You labeled the cylinder and sphere terms correctly, but the "what's this for" is the term due to the parallel axis theorem by moving the half-spheres out away from the rotation axis, which is directly in the middle of the capsule. They seem to have combined the half-sphere moment of inertias and the parallel axis theorem transfer terms in the code.

Hope that helps!


//what's this for?

It accounts for the 2 half-spheres not being at the center of rotation.

EDIT: ph34r.png ninja'd

Please don't PM me with questions. Post them in the forums for everyone's benefit, and I can embarrass myself publicly.

You don't forget how to play when you grow old; you grow old when you forget how to play.

Look at the Parallel Axis Theorem. E.g. here: http://en.wikipedia.org/wiki/Parallel_axis_theorem

Look at the Parallel Axis Theorem. E.g. here: http://en.wikipedia.org/wiki/Parallel_axis_theorem

To help a little more, you add up the inertia tensors for two half-spheres and a cylinder. You use the parallel axis theorem on the half-spheres to move them up/down a long an axis. Look at the tensor generalization formula on the wiki page.

Here's my code for a 'capsule' with different radii. You can simplify it for your case - the idea is as mentioned above - compute the inertia tensors for the half spheres and cylinder, then apply the parallel axis theorem to transform and add the inertia tensors for all parts in the same coordinate frame.


        const Real r1r1 = radius1*radius1;
	const Real r1r2 = radius1*radius2;
	const Real r2r2 = radius2*radius2;
	
	// Compute some initial values needed for the inertia tensor calculation.
	const Real cap1Volume = (Real(2)/Real(3))*math::pi<Real>()*r1r1*radius1;
	const Real cap2Volume = (Real(2)/Real(3))*math::pi<Real>()*r2r2*radius2;
	const Real shaftVolume = (Real(1)/Real(3))*math::pi<Real>()*height*( r1r1 + r1r2 + r2r2 );
	
	volume = cap1Volume + cap2Volume + shaftVolume;
	
	Vector3 cap1CenterOfMass = endpoint1 + axis*(-(Real(3)/Real(8))*radius1);
	Vector3 cap2CenterOfMass = endpoint2 + axis*((Real(3)/Real(8))*radius2);
	Vector3 shaftCenterOfMass = endpoint1 + axis*( (height*(r1r1 + Real(2)*r1r2 + Real(3)*r2r2)) /
												(Real(4)*(r1r1 + r1r2 + r2r2)) );
	
	centerOfMass = (cap1Volume*cap1CenterOfMass + 
					cap2Volume*cap2CenterOfMass + 
					shaftVolume*shaftCenterOfMass) / volume;
	
	const Real cap1InertiaXXYY = Real(0.259)*cap1Volume*r1r1;
	const Real cap1InertiaZZ = Real(0.4)*cap1Volume*r1r1;
	const Real cap2InertiaXXYY = Real(0.259)*cap2Volume*r2r2;
	const Real cap2InertiaZZ = Real(0.4)*cap2Volume*r2r2;
	
	
	Matrix3 Icap1 = Matrix3( cap1InertiaXXYY, Real(0), Real(0),
							Real(0), cap1InertiaXXYY, Real(0),
							Real(0), Real(0), cap1InertiaZZ );
	
	Matrix3 Icap2 = Matrix3( cap2InertiaXXYY, Real(0), Real(0),
							Real(0), cap2InertiaXXYY, Real(0),
							Real(0), Real(0), cap2InertiaZZ );
	
	// Compute the shaft's moment of inertia in the X and Y directions.
	Real shaftInertiaXXYY = ((math::pi<Real>()*height)/Real(30))*
							(height*height*(r1r1 + Real(3)*r1r2 + Real(6)*r2r2) +
							Real(1.5)*(r1r1*r1r1 + r1r1*r1r2 + r1r1*r2r2 + r1r2*r2r2 + r2r2*r2r2));
	
	// Compute the shaft's moment of inertia in the Z direction.
	Real shaftInertiaZZ = ((math::pi<Real>()*height)/Real(10))*
							(r1r1*r1r1 + r1r1*r1r2 + r1r1*r2r2 + r1r2*r2r2 + r2r2*r2r2);
	
	// Compute an inertia tensor for the capsule shaft about its center of mass.
	Matrix3 Ishaft = Matrix3( shaftInertiaXXYY, Real(0), Real(0),
							Real(0), shaftInertiaXXYY, Real(0),
							Real(0), Real(0), shaftInertiaZZ );
	
	// Compute an orientation matrix for this capsule's axis.
	Matrix3 shaftOrientation = Matrix3::planeBasis(axis);
	
	// Compute the total capsule-space inertia tensor.
	Matrix3 Icapsule = util::transformInertiaTensor( Ishaft, shaftVolume, shaftOrientation, shaftCenterOfMass - centerOfMass ) + 
						util::transformInertiaTensor( Icap1, cap1Volume, shaftOrientation, cap1CenterOfMass - centerOfMass ) + 
						util::transformInertiaTensor( Icap2, cap2Volume, shaftOrientation, cap2CenterOfMass - centerOfMass );
	
	// Transform the inertia tensor into shape space.
	volumeInertiaTensor = util::transformInertiaTensor( Icapsule, shaftOrientation );
	

Nice :)
Thanks.

I will definitely have a closer look at the axis theorem.

The lack of consistent and correct derivations, both on the internet and in game programming books, for moments of inertia of shapes like capsules motivated me to write "Moments of Inertia for Common Shapes" in Game Engine Gems, Volume 1. It has full derivations for a box, cylinder, pyramid, cone, ellipsoid, dome, capsule, truncated pyramid, and truncated cone. (It generalizes anything round to the elliptical case.) This chapter also discusses the parallel axis theorem.

I agree with Eric. You find a bunch of scary solutions out there. So I quickly wrote my derivation up for you (if you don't have Eric's book handy) and also summarized some common gotchas to watch out for when using the parallel axis theorem. It is very short and on the point, but hopefully this helps:

https://dl.dropboxusercontent.com/u/55377219/D.%20Gregorius%20-%20Capsule%20Inertia.docx

moment of innertia exist only against its very collision distribution. You can very well reduce this to a bounding volume of a sphere or of a box to aproximate it. It would be grewsomne to distribute homogenous mass against IT tensor even

This topic is closed to new replies.

Advertisement