The Problem

The standard projectile motion problem becomes compounded when one considers the rotational dynamics of the projectile. In order to simplify the problem one considers the translation of a point on the surface of the projectile to be a function of two different translations. The first translation results from the translation of the relative axis. The second translation results from the rotation of the relative axis. This problem concentrates only on the second translation.

We fix the center of mass at the orgin, and give the object an intial rotation vector [wx(0) wy(0) wz(0)] where wx is the angular velocity about the x-axis, wy about the y-axis, and wz about the z-axis. Our problem is to solve [wx(t) wy(t) wz(t)] over a time interval, at discrete time-steps.

We choose the our coordinate system such that Ixy=Iyz=Izx=0 (tensors are zero). We can then find the moments of inertia about each axis by using equations of the form Ixx= integral (y^2+z^2) w.r.t. mass. There are formulas for these quantities in any standard dynamics textbook

Once we have found our Ixx, Iyy, and Izz. We can start to solve our problem. We want to know the solution to the vector w(t) at discrete time-steps.

Steps to the Solution

The solution to this problem can be found by solving the following non-linear system of equations, commonly known as the Euler Equations:

Ixx*wdotx = (Iyy - Izz) * wy * wz

Iyy*wdoty = (Izz - Ixx) * wz * wx

Izz*wdotz = (Ixx - Iyy) * wx * wy

Putting this into matrix form, one gets the following skew-symmetric system:

Sdot(w) = J(w) * S(w)

Where Sdot(w) is the lefthand side of the system of equations,

S(w) = [Ixx*wx] a vector,     and     J(w) = [ 0    wz  -wy]  a 3x3 Skew Matrix
       [Iyy*wy]                              [-wz   0    wx]
       [Izz*wz] (Angular Momentum)           [ wy  -wx   0 ] 

Using Trapezoidal Rule, each iteration is the linear solution of:

S(k+1,tn+1)=S(k,tn) + h * J( (w(k,tn)+w(k,tn+1))/2 ) * ( S(k+1,tn)+S(k+1,tn+1) ) /2

for S(k+1,tn+1) where k is the iteration number, tn is the value of time at the nth grid point, and tn+1 is the value of time at the n+1st grid point. h=stepsize (h=tn+1-tn)

This iteration is rather slow, but it does converge. For a uniform-mass block, the solution is rather periodic. The exact solution to this particular case, is shown below.