waggy wrote:
> I got this alternative version working, and as expected it is much slower and
> less accurate than yours for simple shapes (with a given number of x and y
> divisions, and a similar number of z divisions for my version).  I'll test it
> with more gnarly objects, see if I can find some cases where the alternative
> gives better performance, and post the results and the finished code here.
> 
> Thanks for providing your test objects.  Since I make so many mistakes, I always
> prefer to verify known working cases.  :)

I've now written my moment of inertia macro, and couldn't find an 
optimized way to do it in n^2, so here it is in n^3. My issue here is 
that I'm getting the wrong values for the moments for everything but a 
box, and I can't figure out why. Spheres and Toruses are giving almost 
exactly half the official value. This macro would have even less of an 
audience than my last, but again any math help would be appreciated.
I need to read up more on how to use moment of inertial tensors for 
calculating the moment about an arbitrary axis, and I left my dynamics 
textbooks at home this semester...
Also, this macro requires the center of mass to already be computed.

Chris

(I'm interchanging volume and mass here in the terminology, sorry. It 
will be better once density is involved.)

// Macro, description, and test values:

#declare X_Resolution=20;
#declare Y_Resolution=20;
#declare Z_Resolution=20;
// Finds the moment of inertia for an arbitrary object.
// Return value is <Ix,Iy,Iz> where each is the moment around that axis 
where it passes
// through the center of mass of the object, allowing for "easy" 
calculation of any moment
// with the parallel axis theorem
#macro FindMomentsOfInertia(Object_Identifier, Center_Of_Mass)
   // Set up reference vars
   #local maxcorner = max_extent(Object_Identifier);
   #local mincorner = min_extent(Object_Identifier);
   #local max_x = max(maxcorner.x,mincorner.x);
   #local min_x = min(maxcorner.x,mincorner.x);
   #local max_y = max(maxcorner.y,mincorner.y);
   #local min_y = min(maxcorner.y,mincorner.y);
   #local max_z = max(maxcorner.z,mincorner.z);
   #local min_z = min(maxcorner.z,mincorner.z);
   #local dx = abs(max_x-min_x)/X_Resolution;
   #local dy = abs(max_y-min_y)/Y_Resolution;
   #local dz = abs(max_z-min_z)/Z_Resolution;
   #local dv = dx*dy*dz;
   #local Moments = <0,0,0>;
   #local Volume = 0;
   // Go through entire 3D matrix
   #local curr_x = min_x+dx/2;
   #while (curr_x < max_x)
     #local curr_y = min_y+dy/2;
     #while (curr_y < max_y)
       #local curr_z = min_z+dz/2;
       #while (curr_z < max_z)
         #if(inside(Object_Identifier,<curr_x,curr_y,curr_z>))
           // I = sum(m*r^2) where m is each point mass, r is distance 
to axis
           #local Volume = Volume+dv;
           #local Moments = Moments+dv*<pow(abs(Center_Of_Mass.x-curr_x),2),
                                        pow(abs(Center_Of_Mass.y-curr_y),2),
 
pow(abs(Center_Of_Mass.z-curr_z),2)>;
         #end
         #local curr_z = curr_z+dz;
       #end
       #local curr_y = curr_y+dy;
     #end
     #local curr_x = curr_x+dx;
   #end
   //(Moments)
   DebugVector3("Calculated moments",Moments)
   #debug concat("Calculated volume: ",str(Volume,0,2),"\n")
#end

#declare BoxObj=box{-1,1}
#declare SphereObj=sphere{0,1.5}
#declare TorusObj=torus{1,0.5}

#debug "Box: \n"
#debug concat("Real moments: ",str(1/12*4*(pow(2,2)+pow(2,2)),0,2),"\n")
#debug concat("Real volume: ",str(pow(2,3),0,2),"\n")
FindMomentsOfInertia(BoxObj,<0,0,0>)

#debug "Sphere: \n"
#local sphereMass = 4/3*pi*pow(1.5,3);
#debug concat("Real moments: ",str(2/5*sphereMass*pow(1.5,2),0,2),"\n")
#debug concat("Real volume: ",str(sphereMass,0,2),"\n")
FindMomentsOfInertia(SphereObj,<0,0,0>)

#debug "Torus: \n"
#local torusMass = 2*pow(pi,2)*1*pow(0.5,2);
#debug concat("Real moments: 
",str(torusMass*(pow(1,2)+3/4*pow(0.5,2)),0,2),0,2),"\n")
#debug concat("Real volume: ",str(torusMass,0,2),"\n")
FindMomentsOfInertia(TorusObj,<0,0,0>)

Post a reply to this message