Suspicious results from MassProperties algebra

[avl]

I think there may be a bug somewhere in how MassProperties 'add' operations are handled in rapier3d 0.8.0.

See the following program:

use rapier3d::na::{Point3, Isometry3, Vector3};
use rapier3d::dynamics::MassProperties;
use rapier3d::math::AngVector;

#[test]
pub fn test_rapier_inertial_tensor_calculator() {

    for p in vec![0.01, 100.0] {
        let mut mass_points = Vec::new();
        for z in vec![-50.0,50.0f32] {
            mass_points.push(Point3::new(p*0.5, p*0.5, z));
            mass_points.push(Point3::new(p*0.5, p*-0.5, z));
            mass_points.push(Point3::new(p*-0.5, p*0.5, z));
            mass_points.push(Point3::new(p*-0.5, p*-0.5, z));
        }

        let mut mass_prop = MassProperties::new(Point3::new(0.0,0.0,0.0), 1.0, AngVector::new(1e-3,1e-3,1e-3));

        for p in mass_points.iter().copied() {
            let mass_per_point = 1.0;
            let mut point_mass_prop = MassProperties::new(Point3::new(0.0,0.0,0.0), mass_per_point, AngVector::new(1e-3, 1e-3, 1e-3));
            point_mass_prop = point_mass_prop.transform_by(&Isometry3::new(p.coords, Vector3::new(0.0,0.0,0.0)));
            mass_prop += point_mass_prop;
        }
        println!("p={}: MassProp: {:?}", p, mass_prop.inv_principal_inertia_sqrt); //The small difference for p = 0.01 and 100.0 is implausible
    }
}

It outputs:

p=0.01: MassProp: Matrix { data: [[0.0056241076, 0.0069408445, 0.0060162893]] }
p=100: MassProp: Matrix { data: [[0.0034270522, 0.0037224356, 0.0034776148]] }

As the first object is like a very long thin rod, and the other is a fat box, I think there should be a larger difference between the two.

This is a simplified example. The background to this bug is that objects in my application are behaving in 'obviously' wrong ways, as if they have too high moments of inertia. I'm not an expert in mechanics, so it's entirely possible I've made a mistake somewhere.

[wyattjsmith1]

I know this is quite old, but did you ever resolve this issue? I am noticing something similar where the force required to move an object laterally is substantially higher than the force required to rotate the object.

[avl]

What I finally did is I did the calculations myself, and constructed the MassProperties object like this:

let mass_prop = MassProperties::new(mass_center, total_mass, AngVector::new(Ix,Iy,Iz));

But I did some testing just now. Take a look at this:

let mut point_mass_prop = MassProperties::new(nalgebra::Point3::new(0.0,0.0,0.0), 1.0, AngVector::new(1e-3, 1e-3, 1e-3));
let a = point_mass_prop.transform_by(&Isometry3::new(Vector3::new(0.001,0.0,0.0), Vector3::new(0.0,0.0,0.0)));
let b = point_mass_prop.transform_by(&Isometry3::new(Vector3::new(10000.0,0.0,0.0), Vector3::new(0.0,0.0,0.0)));
println!("A: {:?}", a.reconstruct_inertia_matrix());
println!("B: {:?}", b.reconstruct_inertia_matrix());
assert_eq!(a.reconstruct_inertia_matrix(),b.reconstruct_inertia_matrix()); //This holds!

The two inertia matrices become identical! Despite a huge difference in leverage. The function transform_by has this comment:

// NOTE: we don't apply the parallel axis theorem here
// because the center of mass is also transformed.

I'm starting to think the "bug" is simply that: a) transform_by does not do what I thought it would do b) the struct MassProperties is not simply the InertiaTensor, for the given center of mass, it is a center of mass and an inertia tensor? And using transform_by on it just moves the object including center of mass, and thus does not change the value of reconstruct_inertia_matrix/inv_principal_inertia_sqrt.

[avl]

Actually I think we can close this bug. Possibly one could clarify the docs regarding the field inv_principal_inertia_sqrt in MassProperties, adding something like this:

This is the angular moment of inertia around the object's local center of mass, it's not the moment of inertia with respect to the origin.

[wyattjsmith1]

I also realized my issue seems to be a bug with pixels_per_meter as setting that to 1 resolved my issue.

[sebcrozet]

The MassProperties will indeed give the angular inertia around the object’s center of mass.

It is also worth noting that the difference between the two bodies in the original comment are larger than they can seem. inv_principal_inertia_sqrt is the square root of the inverse angular inertia. So the actual angular inertia value, e.g., along the second principal axis, is 1 / 0.0069408445² = 20757 for p = 0.01 and is 1 / 0.0037224356² = 72168 for p = 100.

I will add the suggested comment (which will happen in parry).