In this article, I hopefully clear the most important parts in rotational dynamics related to 3D rigid body simulation and introduce them in an intuitive way, this is mainly for non-physicists (I'm not even a physicist myself, so that's that...), so if you're a physicist or already familiar with moments of inertia well, then I don't think you're going to get anything new out of this, otherwise then this is probably for you!
It's mainly focused on engineers and programmers, so there's not a whole lot of walls of math text and giggles, mostly focused on the very fundamentals and applications of them
The article is also mainly to make sure I have gotten certain concepts right for the most part, so there may be some mistakes here and there, be warned
Now the way to think of a regular rigid body is not as a totally solid shape per se, but as tons of point masses tied together to form a solid shape (box, sphere, cylinder etc.), and they are all constrained to each other almost perfectly that they literally become a solid shape altogether, this mindset is helpful since we want to approximate the orientation of the body more precisely that way, and apply rotations from specific points on a rigid body, we will derive things out of this idea alone so keep that in mind
You may also realize that this means if you want to take a bunch of point masses and constraint their distance to be hard individually towards each other, then you've still gotten a rigid body that way and you can still approximate the orientation of the body in some way, so really there are other ways you could represent a rigid body, however taking a known static shape and deriving things analytically is much faster which is also the common (preferred) way on most physics engines
This is also exactly the reason why most physics engines has the concept of contact points when detecting collisions between bodies, they are helpful in terms of determining which exact point mass on the rigid body got overlapped with the other body so we apply a rotational force from that specific point mass relative to the nearby distribution of the body mass
In physics, cross products don't work the way you would normally think, so things start to get interesting from here
Say you have two vectors A and B, crossing them together (A x B) would give you a vector that's perpendicular to the oriented plane with it's length equal to the parallelogram area between the two vectors, it's orientation would depend on the product order as going from one vector to another, it's not associative, switching the product order to B x A would give an opposite orientation and the normal vector would be on the opposite side of the plane for instance
This is AB oriented plane and it's orientation order is counter-clockwise like in the image (the orientation goes from A towards B), with it's correspondence purple vector that's normal to the plane
You could also think of this oriented plane as something similar to Euler's formula, where it's an exponential complex number e^ix, the plane AB is where the rotation happen and it's area x being the orientation magnitude on that plane
To make it more meaninigful, we could visualize Euler's formula in a complex plane, keep increasing the orientation magnitude what you would notice is that you get a circular motion on that plane
We could use it as well like the following to describe orientations on a complex plane relative to an identity complex number which would be 1
As you increase pi's magnitude, you would get to i, -1, -i, 1, i, and so on, so it's constrained to be on the circle regardless of orientation magnitude
And we could also simulate Euler's formula in a matrix and use it as a linear map operation on regular vectors
So what we're really doing here is mapping regular vector's A space to some other space B with a matrix used as a linear map operator, so again, it can also be used to describe some oriented space
And you could generalize Euler's formula to 3D as well, and what you would come up with is a quaternion, instead of rotating in a single plane, you would rotate in three planes (in which case the complex parts are three too, i being XY, j being XZ, k being YZ, and the scalar part, it's total area being the orientation magnitude, so it's essentially being a rotation on all three planes simultaneously in 3D rather than on one plane in 2D)
For more info about Euler's formula, you could look up this link, it's excellent: http://www.songho.ca/math/euler/euler.html
Also the reason behind this cross product "magically" working with rotation forces would be explained later, it much has to do with center of masses and forces and their areas which is something that would be called "torque"
So when it comes to physics, cross products are just the resultant of an oriented plane spanned by the two vectors and it's area parallelogram being the magnitude orientation on that plane
Angular Velocity is the rate of rotation expressed in radians per second (or rad/s), it's similar to linear velocity (or commonly referred to as just "velocity") which is only translational rather being rotational, it's defined as omega_w = (r x v) / r^2, r is the distance from the "center of mass" of the body to the point mass, v is the point mass velocity, x here is the cross product, I will describe center of mass later
(Quick note, the omega vector symbol is going to be represented as w for the rest of the article)
If you put a point on a circle and start moving from one place to another at a speed of 1m/s, it would be equavielent to rotating at a speed of 2pi/s which is in radians 1 full revolution instead of one meter
Now in case where the body is only rotating and not moving, you could derive one of the point masses pure/linear velocity through the angular velocity with v = w x r
And you can see when you take the cross product between the point mass velocity and the r vector, you get the angular veloity (in blue) with the same direction as the original one (in orange), implying that it's how it indeed rotates!
However that's not always the case, sometimes the angular velocity goes in a different direction and as to why is going to be explained later
Angular Acceleration accelerates angular velocity, and it's unit is rad/s^2, and is defined as alpha_a = r x a, a is the point mass acceleration
Angular Momentum is somewhat similar from the concept of linear momentum, it's defined as L = r x (m * v) or L = r x P, m is mass, P is momentum, this one is going to be interesting to talk about later
Now in case you have a collision from a specific point on a rigid body, there's going to be a torque to accelerate body's angular velocity
From there we need to apply torque at that specific point, in which case the torque would be something like this
Torque is essentially a rotational force, unlike the "normal" force which has an origin to rotate around but at infinity which is undefined, it's so far that it essentially becomes a translational force (like say gravity that acts on a rigid body), torque has an origin to rotate around, so it makes intuitive sense. in which case you could also say both are equivalent forces, but it's better to distinguish them since it's not really helpful to put them both at the same time
Torque (in this case) is defined as T = r x F, where r is a vector pointing from the center of mass of the rigid body to a specific point on the body, now center of mass is quite literally the center of the mass distribution of the rigid body, the box for instance it's center of mass distribution would be right at the origin of the body's axis, and in case of non-uniform shapes however that's not always true like so
Also another important thing to note is that when you pull the body at a certain point mass of the rigid body and the closer you get to the center of rotation, the harder it is to rotate the object, and the opposite is true, and you could imagine the following case
The further point has greater area which in this case has more magnitude of torque than the one closer to the center of rotation, this is the primary reason behind the use of cross products, it's useful as a way to tell how much we want to orient from the center of mass as we get further away or how parallel is the force's direction with the r vector, so in case of the vector r it's increasing & decreasing the magnitude of the rotational force, so the following image should hopefully make everything clear
And In the context of masses, you would probably notice that force is F = m * a, now m (mass) here is where makes everything more complicated in the context of rotations
The mass here is the moment of inertia tensor, it's basically a rotational mass, in case of 2D it's just a scalar, in 3D it's a 3x3 matrix, which seems to explode by just moving up to one higher dimension, why is that?
In 2D, the oriented plane is just XY (relative to the two dimensional orthogonal coordinates which are just x and y), so there's nothing outside of the oriented plane to take into account, which is really easy to deal with
In 3D however, there's now three oriented planes and they are XY, XZ, YZ, and you can imagine the point masses can be outside one of the three oriented planes, which complicates things even further because rotations only happen inside an oriented plane, not outside of it!
Now in this case we can try to derive moment of inertia for a point mass in a rigid body using the angular momentum equation for our scenario, which is L = m * r x (w x r) (w here is angular velocity)
Now as explained earlier, (w x r) is really just point's velocity direction relative to the center of mass which is constrained to the surface of a circle/sphere in our case like so
So as the time goes, you would notice the point mass is still moving and constrained to the circle/sphere, so essentially it's still rotating
The goal is to try and isolate the angular velocity from the products and derive moment of inertia, you could skip to the result if you don't care much about the derivation
Since r x (w x r) is vector triple product, we can expand it to the following
And turning them into matrices would look like the following
We then take those separated matrices and decrement them
And the result is the moment of inertia for a single point mass out of many, this is also one of many ways of deriving it, I found this to be the easiest I've seen
Now I don't think this would make intuitive sense to you whatsoever at first, but I'll give my best effort to make this a little more intuitive (don't worry, it's not fully intuitive to me either or to most people in fact)
So for the easier part, in case all the points lie on the oriented plane, then the resultant of the vector triple product r x (w x r) would be in the same direction as the angular velocity, in which case the orientation on that oriented plane would be as you would expect
However, in case points don't lie on the oriented plane, this is where things get interesting, the resultant vector out of the vector triple product would not give a vector that's parallel to the angular velocity anymore like so
When you start rotating on that oriented plane, it would also start rotating on another plane simultaneously, as you have bigger mass somewhere that's not uniform across the body, then it would start rotating on another plane that's not the one you're primarily rotating on like so
The sphere here that's attached to the box is mostly the cause for the orientation on the blue part, and orange is mostly caused by combination of the box and sphere
And you could hopefully see that on the final matrix too, the purple parts are the ones right on the planes, and the orange parts are the ones that are outside of the planes that didn't get canceled out by the other sides of the planes, as you notice the orange patterns, there are equivalent parts that appear on two rows at once, as you rotate on one it would resist and start adding some acceleration on the other row
Hopefully this all made sense... and if not, it's still fine I guess, it's important to know what it's purpose for at least
Also, mainly for games that don't rely on realism, it's common to store the moments of inertia as a 3D vector, not a 3x3 matrix
struct Inertia
{
float yz; // resistance on the yz oriented plane
float xz; // resistance on the xz oriented plane
float xy; // resistance on the xy oriented plane
};It's the common way to avoid weird behaviors caused by moments of inertia on non-uniform mass distributed shapes that are not boxes or spheres, it's less physically accurate but it works, this also saves some memory
And also if you have a sphere, then you can just store the inertia as a single float
Now when you rotate the body, the moment of inertia would be invalid immediately since it's axis don't match the body's world axis anymore
To fix this you have to do a special kind of transformation that's specialized to tensors which is R x I x inverse(R), R is the body's orientation matrix, I is the moment of inertia, this would update the inertia's axis to the new axis that matches the body, or you could just put the body in it's local coordinates and go on without doing nasty transformations with the tensor if you wish, both ways work, also tensors have their own way of transformation rules that's not similar to what you're already familiar with, my explanation of this is quite vague, I recommend this video series on tensors in case you're curious
Now the derivation for moments of inertia analytically is a little complex for known shapes, it does require calculus for integration over the area/volume for certain shapes, but there's already a list on wikipedia for various known shapes, you could also use the above resultant moment of inertia tensor and use that for every point that's inside the shape iteratively and sum them up altogether, you could imagine this could get expensive depending on how small the error is and there are methods to speed this process up
The concept of parallel axis theorem would also be handy, here's a link for the formula, it's essentially offsetting the mass distribution of the inertia tensor over it's principle axis (think of it as "translation" like in homogenous algebra, but not exactly), this helps in case you have a shape that's made up of many shapes tied together for a single body and you want to move those around and sum these tensors together in case you did derive the moments at the origin of each of the shapes axis
Now back to torque, to apply torque on the body, you would need something like this (described in C/C++)
// NOTE: this is an example, it's not the most efficient way of doing this!
struct Rigid_Body
{
Vector3 position;
Quaternion orientation;
Vector3 center_of_mass;
Vector3 linear_velocity;
float mass;
Matrix3x3 inertia; // moment of inertia
Vector3 angular_velocity;
};
void apply_torque(
Rigid_Body* body, Vector3 force, Vector3 point,
float delta_time)
{
Vector3 world_com = body->position + body->center_of_mass; // world space center of mass
Vector3 r = point - world_com; // from center of mass to the point
Vector3 torque = cross(r, force);
Matrix3x3 orientation_matrix = get_quaternion_matrix(body->orientation);
Matrix3x3 world_inertia = orientation_matrix * body->inertia * matrix_inverse(orientation_matrix);
Vector3 J = (matrix_inverse(world_inertia) * torque) * delta_time;
body->angular_velocity += J;
}And to update body's orientation every tick, you do the following
void update_body_orientation(Rigid_Body* body, float delta_time)
{
// Internal torque caused by moments of inertia
Matrix3x3 orientation_matrix = get_quaternion_matrix(body->orientation);
Matrix3x3 world_inertia = orientation_matrix * body->inertia * matrix_inverse(orientation_matrix);
body->angular_velocity += (matrix_inverse(world_inertia) * cross(body->angular_velocity, world_inertia * body->angular_velocity)) * delta_time;
Vector3 world_COM = body->position + quaternion_rotate_vector(body->orientation, body->center_of_mass);
Vector3 COM_to_pos = body->position - world_COM; // from world center of mass to the body's origin
Vector3 oriented_plane = body->angular_velocity;
float orientation_magnitude = length(body->angular_velocity) * delta_time;
Quaternion theta = quaternion_angle_axis(orientation_magnitude, oriented_plane);
body->orientation = theta * body->orientation; // orientation composition
// rotate the body origin with respect to the updated center of mass
body->position = world_COM + quaternion_rotate_vector(theta, COM_to_pos);
}And that's all there is to it.
Since we're pretty much done here, there are various other resources I would recommend, since this article is supposed to give you a small push into rotational dynamics in general with the fundamentals, it doesn't contain enough applications and examples, but here's one that does, it's also part of the entire Pixar's physically based modeling course, and here's another one that goes much deeper into the rotational dynamics specifically in 3D, which is also part of Chris Hecker's magazine series
That's it, thank you for reading the article.