6DoF Rigid Body Dynamics

If you throw an arbitrarily-shaped rigid object into the air with some random rotational motion, the motion can proceed semi-chaotically, unless it happens to be spinning purely around one of its “principle axes”. (You could try your mobile phone, for example – give most spin around its tall axis for interesting results.) General rigid body rotational motion is not simple, but it is a universal scenario, with applications ranging from submarine simulators to satellite design.

The Analogy

There is an analogy between translational dynamics and rotational dynamics. Force is analogous to torque, velocity to angular velocity, and so on. This analogy can easily mislead one into believing that rotational motion can be described as simply as translational motion. So why, in spite of this simple analogy, is rotational motion so much more complex than translational motion?

The signal diagram below summarizes the entire process: Euler’s laws of motion (derived from Newton’s laws of motion), as well as the relevant kinematic principles.

All the shaded blocks in this diagram are 3 dimensional vectors, except Orientation, which is some other 3-DOF representation (See Expressing Rotation).

• $\int dt$ denotes vector integration over time.
• “rotational integral” denotes some form of (non-linear) integration that depends on the representation used for orientation.
• $m^{-1}$ denotes scalar division by the body’s mass.
• $\bf{I}^{-1}$ denotes matrix left-multiplication by the matrix inverse of the body’s inertia tensor.
• The arrow from Orientation to $\bf{I}^{-1}$ denotes the fact that $\bf{I}$‘s orientation is attached to the body’s frame of reference.

Loosely thinking, we can say (keeping in mind that we are dealing with vectors): Force tells momentum what to do. Momentum implies velocity. Velocity tells position what to do.

This holds in some sense for both translational motion and rotational motion. The difference between them arises in the italicized relationship above. Translational momentum and velocity are simply related by a scale factor (the mass), whereas their angular analogs are related by a linear transformation (the inertia tensor).

Inertia Tensor

The inertia tensor is a symmetric matrix, and therefore has orthogonal eigenspaces. Its eigenvalues are the moments of inertia about its eigenvectors (called the principle axes). So the angular velocity vector is a distorted version of the angular momentum vector, biased away from dominating principle axes. For example, a disc-shaped body will distort angular velocity towards the disc’s plane, while a rod-shaped body will distort angular velocity towards the rod’s line.

So, returning to our object in the air, we know that for all its convoluted gymnastics, its angular momentum remains constant. Therefore its angular velocity is determined by its orientation. That is to say, the body is turning about an axis that is influenced by its instantaneous orientation. It is easy to see how this can lead to some degree of chaos in the general case.

Dynamic Model

From the diagram above, we can straightforwardly derive a dynamic model. The translational part of the model is defined by these equations:

$\large{\mathbf{\dot{x}} = \mathbf{v}}$

$\large{\mathbf{v} = \mathbf{p}/m}$

$\large{\mathbf{\dot{p}} = \mathbf{F}}$

We will choose the matrix representation of rotation (see Expressing Rotation), so that $\mathbf{R}\mathbf{I}\mathbf{R^{-1}}$ is the oriented inertia tensor. The rotational part of the model is then defined by these equations:

$\large{\mathbf{\dot{R}} = \boldsymbol{[\omega]_{\times}\mathbf{R}}}$

$\large{\boldsymbol{\omega} = \mathbf{R}\mathbf{I^{-1}}\mathbf{R^{-1}}\mathbf{L}}$

$\large{\mathbf{\dot{L}} = \boldsymbol{\tau}}$

Where

• $\mathbf{x}$ is the 3D position vector
• $\mathbf{v}$ is the 3D velocity vector
• $\mathbf{p}$ is the 3D momentum vector
• $m$ is the body’s mass
• $\mathbf{F}$ is the 3D force vector (sum of forces acting on the body)

and

• $\mathbf{R}$ is the 3×3 rotation matrix
• $\boldsymbol{\omega}$ is the 3D angular velocity vector
• $\mathbf{L}$ is the 3D angular momentum vector
• $\boldsymbol{\tau}$ is the 3D torque vector (sum of torques acting on the body)
• $\boldsymbol{[\quad]_{\times}}$ is what we might call a 3×3 “cross-product matrix”: the skew-symmetric matrix such that $\mathbf{[a]_{\boldsymbol{\times}}}\mathbf{b}=\mathbf{a}\boldsymbol{\times}\mathbf{b}$ (for any 3D vectors $\mathbf{a}$ and $\mathbf{b}$).

$\boldsymbol{[a]_{\times}}$ is given by $\begin{bmatrix}\,\,0&\!-a_3&\,\,\,a_2\\\,\,\,a_3&0&\!-a_1\\\!-a_2&\,\,a_1&\,\,0\end{bmatrix}$

The diagram below is the Simulink implementation of the model above, with a trivial difference: It is convenient to integrate acceleration rather than force, so that gravitational effects can be introduced as an acceleration without an additional dependence on the body’s mass.

The rotational kinetic energy (given by $\mathbf{L}\cdot\boldsymbol{\omega}$) is conserved under zero input-torque conditions. This serves as a useful sanity check for the model’s correctness.

The rotation matrix is supposed to remain an orthogonal matrix. However, since numerical round-off errors would cause it to drift over time, it needs to be re-normalized and re-orthogonalized occasionally. For most applications this issue can be ignored. I originally developed this model for a helicopter flight simulator, and after hours of test-flying the rotation matrix still had not drifted detectably.

6 Responses to “6DoF Rigid Body Dynamics”

• Erik says:

Thanks for this great example. It’s exactly what I’m looking for. I am trying to build this within Simulink, but I am struggling with the Inertia Tensor & Cross Product Matrix blocks. Would it be possible to spell out exactly what is in them and/or send me (or set up a link to) a runnable copy of the model please? Thanks in advance.

I’m glad you found it useful Erik. What are you modelling?

The Inertia Tensor (principle) block is a constant 3×3 matrix. The matrix itself depends on the shape of the rigid body, but for any rigid body it can be expressed as a diagonal matrix by choosing a convenient reference orientation. So in the end it is determined by only 3 numbers, corresponding to the body’s moment of inertia about each of 3 principle axes. Symmetry in the body can simplify things further. For example, if the body is a sphere or a cube, the three diagonal elements are equal.

For a Simulink implementation of the Cross Product Matrix block, take a look at section 9.2.6 of this document. There may also be some discussion of inertia tensors in there, if my memory serves.

• Erik says:

Thanks for the quick reply. I’ve just skim read your document; it looks very interesting and I will read it in more detail later. I did an implementation of Cross Product Matrix using Embedded Matlab (which may not have existed when you did this report) as follows:
 function a_cross = fcn(a) %#codegen

 

a_cross = [ 0, -a(3), a(2); a(3), 0, -a(1); -a(2), a(1), 0 ]; 

I’m still working on the Inertia Tensor constant. I am trying to model a car.

Yes that’s equivalent, and clearer.

For a back-of-the-envelope approximation I would go with a cuboid plus a point-mass in front representing the engine. Note that the point-mass would shift the cuboid’s center of mass, necessitating a few calculations.

For best accuracy one would derive the tensor directly from a 3D description of the car’s density over space, and one could also include in the model gyroscopic torques from the spinning wheels and flywheel, like I had to include the helicopter’s rotors in that report. But perhaps those effects are negligible for a car.

• Erik says:

Cuboid is fine for now. Thanks for the link.

I think the engine can have an effect, depending on type, but a bit beyond the stage I’m at.

Thanks again for all your help.

• omer says:

thank you for great explanation, very important topic for engineering students