Skip to content

Dynamics

Acceleration of a rigid body

Linear Velocity

(i) Since from equation in general case

\[ \begin{align} \frac{d}{dt}({^A_BR}{^BQ})={^AV_Q}={^A_BR}{^BV_Q}+{^A\Omega_B}\times {^A_BR}{^BQ}+{^AV_B}.\label{linear-vel} \end{align} \]

take derivative on both sides of the second "\(=\)" we have

\[ ^A\dot{V}_Q=\frac{d}{dt}({^A_BR}{^BV_Q})+{^A\dot{\Omega}_B}\times ({^A_BR {^BQ}})+{^A\Omega_B}\times \frac{d}{dt}({^A_BR}{^BQ})+{^A\dot{V}_B}, \]

apply first "\(=\)" of equation \(\ref{linear-vel}\) to the brackets, we have

\[ \begin{align} \nonumber{^A\dot{V}_Q}=&{^A_BR}{^B\dot{V}_Q}+{^A\Omega_B}\times{^A_BR}{^BV_Q}+{^A\dot{\Omega}_B}\times {^A_BR {^BQ}}\\ \nonumber&+{^A\Omega_B}\times({^A_BR}{^BV_Q}+{^A\Omega_B}\times {^A_BR}{^BQ})+{^A\dot{V}_B}\\ \nonumber=&{^A_BR}{^B\dot{V}_Q}+2{^A\Omega_B}\times{^A_BR}{^BV_Q}+{^A\dot{\Omega}_B}\times {^A_BR {^BQ}}\\ &+{^A\Omega_B}\times{^A\Omega_B}\times {^A_BR}{^BQ}+{^A\dot{V}_B}.\label{linear-vel-general} \end{align} \]

For revolute joint, we have \(^BV_Q={^B\dot{V}_Q}=0\), so the above equation could be simplified as

\[ \begin{align} {^A\dot{V}_Q}={^A\dot{\Omega}_B}\times {^A_BR {^BQ}}+{^A\Omega_B}\times{^A\Omega_B}\times {^A_BR}{^BQ}+{^A\dot{V}_B}.\label{linear-vel-revolute} \end{align} \]

Angular Velocity

Consider the case in which \(\{B\}\) is rotating relative to \(\{A\}\) with \(^AQ_B\) and \(\{C\}\) is rotating relative to \(\{B\}\) with \(^B\Omega_C\). Then the angular velocity of \(\{C\}\) in frame \(\{A\}\) is

\[ ^A\Omega_C={^A\Omega_B}+{^AR_B}{^B\Omega_C}, \]

take derivative on both sides we have

\[ \begin{align} \nonumber {^A\dot{\Omega}_C}=&{^A\dot{\Omega}_B}+\frac{d}{dt}({^A_BR}{^B\Omega_C})\\ =&{^A\dot{\Omega}_B}+{^A\Omega_B}\times{^A_BR}{^B\Omega_C}+{^A_BR}{^B\dot{\Omega}_C}.\label{angular-vel} \end{align} \]

supplementary info: Fixed point motion of a rigid body

Reference

  • 理论力学, 马尔契夫

  • 理论物理学教程, 朗道

Here we give a deduction of motion of a rigid body in terms of physics.

To easily get the result, we view the rigid body as a set of discrete points. As for a continuous body, we change the mass \(m_i\) of a point \(i\) into \(dV\) with unit mass \(\rho dV\) and then integrate over it.

We have two basic frames, the fixed one and the moving one. The former one fr1: \(\{X,Y,Z\}\) is viewed as a global frame in our discussion and the latter one fr2: \(\{x_1,x_2,x_3\}\) is moving with the rigid, whose origin is usually set in the center of mass. The origin vector of fr2 is expressed in fr1 as \(\pmb{R}\) and the orientation of fr2 expressed in fr1 is expressed as a rotation matrix, which we would later discuss.

Then the motion can be decomposed into two parts, the first is the translation of center of mass, and the second is the rotation along the center of mass.

Assume any point \(P\) on the rigid body, we have its position \(\pmb{r}\) expressed in fr2 and \(\pmb{\tau}\) expressed in fr1. Assume is small time, it has a rotation \(d\pmb{\phi}\times \pmb{r}\), then

\[ d\pmb{\tau} = d\pmb{R} + d\pmb{\phi}\times \pmb{r} \]

divide both sides with \(dt\) we have

\[ \frac{d\pmb{\tau}}{dt}=\frac{d\pmb{R}}{dt}+ \frac{d\pmb{\phi}}{dt}\times \pmb{r} \]

introduce velocity

\[ \pmb{v}=\frac{d\pmb{\tau}}{dt}, \quad \pmb{V}=\frac{d\pmb{R}}{dt},\quad \pmb{\Omega}=\frac{d\pmb{\phi}}{dt} \]

and we have

\[ \begin{equation} \pmb{v}=\pmb{V}+\pmb{\Omega}\times \pmb{r}.\label{vel} \end{equation} \]

interpretation 1

Here \(\pmb{V}\) is the velocity of the origin of center of mass, \(\pmb{\Omega}\) is the rotation velocity of the rigid body.

fr2 do not coincide with center of mass

Assume center of mass if \(O\) and the origin of fr2 is \(O'\), with \(O'O=\pmb{a}\). Assume we have velocity of \(O'\) is \(\pmb{V}'\), its rotation velocity is \(\pmb{\Omega}'\).

Then the same point \(P\) has its position expressed in fr2 is \(\pmb{r}'\). So we have \(\pmb{r}=\pmb{r}'+\pmb{a}\). From equation \(\ref{vel}\) we have

\[ \pmb{v}=\pmb{V}+\pmb{\Omega}\times \pmb{r}'+\pmb{\Omega}\times \pmb{a}. \]

However we the same definition we still should have

\[ \pmb{v}=\pmb{V}'+\pmb{\Omega}'\times \pmb{r}' \]

combining the above two we have

\[ \begin{align} \pmb{V}'&=\pmb{V}+\pmb{\Omega}\times \pmb{a}\label{vel-translate-relation}\\ \pmb{\Omega}'&=\pmb{\Omega}.\label{vel-angular-relation} \end{align} \]

The above deduction is important, since equation \(\ref{vel-angular-relation}\) tells us the rotation velocity is independent of the frame 2. For all the frames, the rigid nody rotates around the axes which are parallel. Equation \(\ref{vel-translate-relation}\) tells us the following property.

Relationship of \(\pmb{V}\) and \(\pmb{\Omega}\)

If origin \(O\) is chosen, \(\pmb{V}\) and \(\pmb{\Omega}\) are perpendicular, then for all \(O'\), its \(\pmb{V}'\) and \(\pmb{\Omega}'\) are also perpendicular, since

\[ \pmb{V}'\cdot\pmb{\Omega}'=(\pmb{V}+\pmb{\Omega}\times \pmb{a})\cdot \pmb{\Omega}=0. \]

Easy to check that, in this situation, we the velocity of all points are located on the same plane. Then by equation \(\ref{vel-translate-relation}\), we could always find another origin \(O'\) such that \(\pmb{V}'=0\), which means the rigid body is only making rotation motion.

Mass Distribution

We introduce inertia tensor, which can be thought of as a generalization of the scalar moment of inertia of an object.

The kinetics of the rigid body expressed in discrete points, is

\[ T = \sum_{i}\frac{m_i \pmb{v}_i^2}{2}=\sum \frac{m\pmb{v}^2}{2} \]

to simplify we assume the summation is over \(i\).

Substitute equation \(\ref{vel}\) in and we have

\[ \begin{align*} T &= \sum\frac{m}{2}(\pmb{V}+\pmb{\Omega}\times \pmb{r})^2\\ &=\sum \frac{m}{2}\pmb{V}^2+\sum m\pmb{V}\cdot (\pmb{\Omega}\times \pmb{r})+\sum\frac{m}{2}(\pmb{\Omega}\times \pmb{r})^2 \end{align*} \]

the first is easy to be written as \(\frac{\mu\pmb{V}^2}{2}\), the translation of rigid boty. The second item could be derived by mixed product \(a\cdot(b\times c)=c\cdot(a\times b)\)

\[ \sum m\pmb{V}\cdot (\pmb{\Omega}\times \pmb{r})=\sum m\pmb{r}\cdot (\pmb{V}\times\pmb{\Omega}) \]

for all the point in the rigid body, \(\pmb{V}\) and \(\pmb{\Omega}\) are of the same, then by definition of center of mass

\[ \sum m\pmb{V}\cdot (\pmb{\Omega}\times \pmb{r})=(\pmb{V}\times\pmb{\Omega})\cdot (\sum m\pmb{r})=0. \]

As for the third part, we have

\[ \sum\frac{m}{2}(\pmb{\Omega}\times \pmb{r})^2=\sum\frac{m}{2} (\pmb{\Omega}^2\pmb{r}^2-(\pmb{\Omega}\cdot\pmb{r})^2), \]

which is called the rotational part of kinetics. To simplify the notation, we use \(x_i\), \(\Omega_k\) to be the coordinate of \(\pmb{r}\) and \(\pmb{\Omega}\), respectively. So

\[ \begin{align*} T_{bp}&=\frac{1}{2}\sum m [\Omega_i^2x_l^2-\Omega_ix_i\Omega_kx_k]\\ &=\frac{1}{2}\sum m [\Omega_i\Omega_k\delta_{ik}x_l^2-\Omega_ix_i\Omega_kx_k]\\ &=\frac{1}{2}\Omega_i\Omega_k\sum m (x_l^2\delta_{ik}-x_ix_k)\\ \end{align*} \]

If we introduce tensor

\[ I_{ik}=\sum m (x_l^2\delta_{ik}-x_ix_k) \]

then

\[ T = \frac{\mu \pmb{V}^2}{2} + \frac{1}{2}I_{ik}\Omega_i\Omega_k \]

Define a function called Lagrange function by subtracting \(T\) with potential energy \(U\)

\[ \begin{align} L=T-U=\frac{\mu \pmb{V}^2}{2} + \frac{1}{2}I_{ik}\Omega_i\Omega_k-U.\label{Lagrange-func} \end{align} \]

We define the inertia od moment of the system with respect to an axis \(u\),

\[ J_u=\sum_{v=1}^N m_v\rho_v^2, \]

where \(\rho_v\) equals the distance of point \(P_v\) of the system. Then we focus on passing on difference axes which pass the same point.

Assume \(u\) pass the origin of the coordinate \(Oxyz\), with its cosine along \(Ox\), \(Oy\), \(Oz\) is denoted by \(\alpha\), \(\beta\), \(\gamma\), respectively. Then denote \(e_u=(\alpha, \beta, \gamma)\), we have

\[ \begin{align} \nonumber J_u&=\sum_{v=1}^N m_v \rho_v^2=\sum_{v=1}^N m_v [r_v^2-(r_v\cdot e_u)^2]\\ \nonumber &=\sum_{v=1}^N m_v[x_v^2+y_v^2+z_v^2- (x_v\alpha+y_v\beta+z_v\gamma)^2]\\ \nonumber &=\sum_{v=1}^N m_v [(1-\alpha^2)x_v^2+(1-\beta^2)y_v^2+(1-\gamma^2)z_v^2-2\alpha \beta x_vy_v -2\alpha \gamma x_vz_v-2\gamma \beta z_vy_v]\\ \nonumber &=\sum_{v=1}^N m_v[(\beta^2+\gamma^2)x_v^2+(\alpha^2+\gamma^2)y_v^2+(\alpha^2+\beta^2)z_v^2-2\alpha \beta x_vy_v -2\alpha \gamma x_vz_v-2\gamma \beta z_vy_v]\\ \nonumber &=\sum_{v=1}^N m_v [(y_v^2+z_v^2)\alpha^2+(x_v^2+z_v^2)\beta^2+(z_v^2+y_v^2)\gamma^2 -2\alpha \beta x_vy_v -2\alpha \gamma x_vz_v-2\gamma \beta z_vy_v]\\ &=J_x\alpha^2+J_y\beta^2+J_z\gamma^2-2J_{xy}\alpha\beta-2J_{xz}\alpha\gamma-2J_{yz}\beta\gamma\label{expand-expresstion} \end{align} \]

where \(J_x, J_y, J_z\) are the same definition as the following definition.

Properties of Inertia Tensor

(i) \(I_{ik}=I_{ki}\).

(ii) we have the following tensor matrix

\[ \pmb{J}=I_{ik}=\left[\begin{array}{ccc} \sum m (y^2+z^2) & -\sum m xy & -\sum m xz\\ -\sum m yx & \sum m (x^2+z^2) & -\sum m yz\\ -\sum m zx & -\sum m zy & \sum m (x^2+y^2) \end{array}\right]=\left[\begin{array}{ccc} I_{xx} & -I_{xy} & -I_{xz}\\ -I_{yx} & I_{yy} & -I_{yz}\\ -I_{zz} & -I_{zy} & I_{zz} \end{array}\right]. \]

where \(I_{xx}\), \(I_{yy}\) and \(I_{zz}\) are called the mass moments of inertia, or axis moment of inertia. And the other elements, whose integrals over mixed indices are called the mass products of inertia, or centrifugal moment of inertia.

(iii) for each element, we could integrate

\[ I_{ik}=\int \rho (x_l^2\delta_{ik}-x_ix_k)dV \]

(iv) Equation \(\ref{expand-expresstion}\) has a geometric meaning. That is, if we cut down a segment on axis \(u\), denoted by \(ON\), satisfies \(|ON|=1/\sqrt{J_u}\), then \(N(x,y,z)\) satisfies

\[ \alpha=\sqrt{J_u}x,\quad \beta=\sqrt{J_u}y,\quad \gamma=\sqrt{J_u}z. \]

substitute in the equation \(\ref{expand-expresstion}\), we have

\[ J_x x^2+J_y y^2+J_z z^2-2J_{xy}xy-2J_{xz} xz-2J_{yz} yz=1. \]

which is an ellipsoid. If we choose an appropriate coordinate, we have its expression simplified as

\[ Ax_*^2+By_*^2+Cz_*^2=1. \]

that is, the mass product of moments equal zero. In this case, we call \(x_*\), \(y_*\), \(z_*\) the prinpiple axes, the corresponding coefficients \(A\), \(B\), \(C\) are called principle moments of inertia.

Actually, \(A,B,C\) are the eigenvalue of the matrix \(\pmb{J}\). If they are distinct, then the principle axis is uniquely determined.

Example. Some typical inertia tensor.

(i) For planar system, we have \(x_1,x_2\) in the plane, so \(x_3=0\) for all particles,

\[ I_1=\sum m x_2^2,\quad I_2=\sum m x_1^2,\quad I_3=I_1+I_2. \]

(ii) For linear system, we have \(x_3\) in the line, so \(x_1=x_2=0\), so

\[ I_1=I_2=\sum m x_3^2, \quad I_3=0. \]

Moment of Momentum

Multiple cross product

Show that

\[ \pmb{a}\times (\pmb{b}\times \pmb{c})=(\pmb{a}\cdot \pmb{c})\pmb{b}-(\pmb{a}\cdot\pmb{b})\pmb{c}. \]

Notice that \(e_i\times (e_j\times e_i)=e_j\), so

\[ \begin{align*} \pmb{a}\times (\pmb{b}\times \pmb{c})&=\sum a_i b_j c_k(e_i\times(e_j\times e_k))\\ &=\sum a_i b_j c_i e_j - \sum a_i b_i c_j e_j\\ &=(\sum a_i c_i)\sum b_je_j-(\sum a_i b_i)\sum c_je_j\\ &=(\pmb{a}\cdot \pmb{c})\pmb{b}-(\pmb{a}\cdot\pmb{b})\pmb{c} \end{align*} \]

Still, we have the origin of moving frame coincides with the center of mass. Then by definition of moment of momentum

\[ \begin{align*} \pmb{M}&=\sum m \pmb{r}\times \pmb{v}\\ &=\sum m \pmb{r}\times (\pmb{\Omega}\times \pmb{r})\\ &=\sum m (\pmb{r}\cdot\pmb{r})\pmb{\Omega}-(\pmb{r}\cdot\pmb{\Omega})\pmb{r}. \end{align*} \]
\[ \begin{align} \nonumber M_i&=\sum m (x_l^2\Omega_i-x_ix_k\Omega_k)\\ \nonumber&=\Omega_k \sum m (x_l^2\delta_{ik}-x_ix_k)\\ &=\Omega_kI_{ik}.\label{moment-tensor} \end{align} \]

here we have to sum over \(k\).

If \(\pmb{r}=(x_1,x_2,x_3)^T\), \(\pmb{\Omega}=(\Omega_1, \Omega_2, \Omega_3)^T\), then

\[ \begin{align*} M_i&=\sum_v m_v (\sum_j x_j^2)\Omega_i - \sum_v m_v (\sum_j x_j\Omega_j)x_i\\ &=\sum_v m_v (\sum_{j\neq i} x_j^2)\Omega_i- \sum m_v [\sum_{j\neq i} x_ix_j \Omega_j]\\ &=J_i \Omega_i - \sum_{j\neq i}J_{ij}\Omega_j. \end{align*} \]

or in language of matrix

\[ \pmb{M}=\pmb{J}\pmb{\Omega} \]

If \(Ox\), \(Oy\), \(Oz\) are the principle axes, then \(\pmb{J}\) is diagonal, so \(M_1=J_x\Omega_1\), \(M_2=J_y\Omega_2\), \(M_3=J_z\Omega_3\). Generally speaking, the moment of momentum \(\pmb{M}\) is not parallel to angular velocity \(\pmb{\Omega}\), unless \(\pmb{J}\) is disgomal, or the coordinate is the principle axes, meaning the three directions are decoupled.

Dynamic Equations

If we know the center of mass and the inertia of the link, then its mass distribution is completely characterized. To accelerate or decelerate it, we have to apply forces and torque.

To get the equations, we take derivative of \(\pmb{p}\) and \(\pmb{M}\) with respect to \(t\).

Translation ODE

Translation ODE

(i) For each particle, we have \(\dot{\pmb{p}}_i=\pmb{f}_i\), then sum over \(i\), and we have

\[ \begin{align} \pmb{F}=\mu \dot{V}.\label{eq1} \end{align} \]

the force is from the outside of system.

(ii) And the force \(\pmb{F}\) could be expressed by the potential energy of the rigid body in the outer field

\[ \pmb{F}=-\frac{\partial U}{\partial \pmb{R}}. \]

(iii) so we have equation \(\ref{eq1}\) using Lagrange function

\[ \frac{d}{dt}\frac{\partial L}{\partial \pmb{V}}=\frac{\partial L}{\partial \pmb{R}}. \]

Since \(\pmb{P}=\sum m\pmb{p}=\mu \pmb{V}\), we have

\[ \begin{align*} \pmb{F}&=\sum\pmb{f}=\sum \frac{\pmb{p}}{dt}=\mu \dot{\pmb{V}}. \end{align*} \]

\(\square\)

When the rigid body translates \(\delta\pmb{R}\), then the position vector \(\tau\) of each particle changes, the change of the corresponding potential

\[ \begin{align*} \delta U&=\sum \frac{\partial U}{\partial \pmb{\tau}}\cdot \pmb{\tau}\\ &=\delta \pmb{R}\cdot \sum \frac{\partial U}{\partial \pmb{\tau}}\quad (\delta\pmb{R}=d\pmb{\tau})\\ &=-\delta \pmb{R}\cdot\sum \pmb{f}\\ &=-\pmb{F}\cdot \delta \pmb{R}. \end{align*} \]

\(\square\)

By equation \(\ref{Lagrange-func}\), we have

\[ \frac{\partial L}{\partial \pmb{V}}=\mu \pmb{V}. \]

and using (ii) we have the result.

Rotation ODE

Rotation ODE

Assume the center of mass does not move. Then

(i)

\[ \frac{d\pmb{M}}{dt}=\pmb{K}=\sum \pmb{r}\times \pmb{f}. \]

where \(\pmb{K}\) is called the torque.

(ii) Using Lagrange function

\[ \frac{d}{dt}\frac{\partial L}{\partial \pmb{\Omega}}=\frac{\partial L}{\partial \pmb{\varphi}}. \]

(iii) when the origin of the frame translates, we have \(\pmb{r}=\pmb{r}'+\pmb{a}\), so

\[ \pmb{K}=\sum \pmb{r}\times \pmb{f}=\sum \pmb{r}'\times \pmb{f}+\sum \pmb{a}\times \pmb{f}=\pmb{K}'+\pmb{a}\times \pmb{F}. \]

Taking derivative of moment of momentum \(\pmb{M}\)

\[ \dot{\pmb{M}}=\frac{d}{dt} \sum \pmb{r}\times \pmb{p}=\sum \dot{\pmb{r}}\times \pmb{p} + \sum \pmb{r}\times \dot{\pmb{p}}. \]

Because \(\pmb{V}=0\), we have \(\dot{\pmb{r}}=\pmb{v}=\dot{\pmb{\tau}}\). Since \(\pmb{v}\) and \(\pmb{p}\) share the same direction, we have \(\dot{\pmb{r}}\times \pmb{p}=0.\) Replacing \(\dot{\pmb{p}}\) by \(\pmb{f}\), we have

\[ \frac{d\pmb{M}}{dt}=\pmb{K}=\sum \pmb{r}\times \pmb{f}. \]

\(\square\)

Take derivative of Lagrange function, or with respect to the coordinate of \(\pmb{\Omega}\)

\[ \frac{\partial L}{\partial \Omega_i}=I_{ik}\Omega_k=M_i. \]

the same as equation \(\ref{moment-tensor}\).

Then when the rigid body rotates about a small angle \(\delta \varphi\), the change of potential

\[ \begin{align*} \delta U&=-\sum \pmb{f}\cdot \delta\pmb{\tau}\\ &=-\sum \pmb{f}\cdot (\delta\pmb{\varphi}\times\pmb{r})\\ &=-\delta\pmb{\varphi}\sum \pmb{r}\times \pmb{f}\quad \delta\pmb{\varphi}\text{ are the same}\\ &=-\pmb{K}\cdot \delta\pmb{\varphi}. \end{align*} \]

and we are done.

\(\square\)

When \(\pmb{F}\) and \(\pmb{K}\) are perpendicular, then we could find a vector \(\pmb{a}\) such that

\[ \pmb{K}=\pmb{a}\times \pmb{F}, \quad \pmb{K}'=0. \]

The \(\pmb{a}+\pmb{F}\) is also an acceptable solution.

Example. Uniform force field. The force acting on particle is \(\pmb{f}=e\pmb{E}\), then

\[ \pmb{F}=\sum \pmb{f}=\pmb{E}\sum e,\quad \pmb{K}=\sum e\pmb{r}\times \pmb{E}. \]

assume \(\sum e\neq 0\), then let the center of 'mass' to be

\[ \pmb{r}_0=\frac{\sum e\pmb{r}}{\sum e} \]

so

\[ \pmb{K}=\pmb{r}_0\times \pmb{F}. \]

The above deduction is based on a fixed frame (either not translation or rotation).

Newton-Euler's Equation

Newton-Euler's Equation

For a fixed point motion of a rigid body, we have translation part ODE and rotation part ODE

\[ \mu\left[\frac{d\pmb{V}}{dt}+\pmb{\Omega}\times \pmb{V}\right]=\pmb{F}, \quad \pmb{J}\frac{d\pmb{\Omega}}{dt}+\pmb{\Omega}\times \pmb{J}\pmb{\Omega}=\pmb{K}. \]

Any vector \(\pmb{A}\) has a changing velocity with respect to a frame which is moving.

Assume \(d\pmb{A}/dt\) is the changing velocity of \(\pmb{A}\) with respect to a fixed frame. If \(\pmb{A}\) do not change in a rotating frame, then

\[ \frac{d\pmb{A}}{dt}=\pmb{\Omega}\times \pmb{A}. \]

If \(\pmb{A}\) moves in the rotating frame, then we add \(d'\pmb{A}/{dt}\) (we just use the results from fixed frame here)

\[ \frac{d\pmb{A}}{dt}=\frac{d'\pmb{A}}{dt}+\pmb{\Omega}\times \pmb{A}. \]

So apply this regulation to \(\pmb{P}\) and \(\pmb{M}\), we have

\[ \frac{d'\pmb{P}}{dt}+\pmb{\Omega}\times \pmb{P}=\pmb{F}, \quad \frac{d'\pmb{M}}{dt}+\pmb{\Omega}\times \pmb{M}=\pmb{K}. \]

Project to the coordinate of the moving frame, we have \(\left(\frac{d'\pmb{P}}{dt}\right)=\frac{dP_1}{dt},\cdots, \left(\frac{d'\pmb{M}}{dt}\right)=(\pmb{J}\frac{d\pmb{\Omega}}{dt})_1\), so the translation part ODE is

\[ \begin{align*} \mu\left(\frac{dV_1}{dt}+\Omega_2V_3-\Omega_3V_2\right)&=F_1,\\ \mu\left(\frac{dV_2}{dt}+\Omega_3V_1-\Omega_1V_3\right)&=F_2,\\ \mu\left(\frac{dV_3}{dt}+\Omega_1V_2-\Omega_2V_1\right)&=F_3.\\ \end{align*} \]

and the rotation part ODE is (a little complicated, only present the first line)

\[ I_x \frac{d\Omega_1}{dt}-I_{xy}\frac{d\Omega_2}{dt}-I_{xz}\frac{d\Omega_3}{dt}+ (I_z-I_y)\Omega_2\Omega_3+ I_{yz}(\Omega_3^2-\Omega_2^2)+\Omega_1(I_{xy}\Omega_3-I_{xz}\Omega_2)=K_1, \]

if the coordinate \(Ox,Oy,Oz\) of moving system are principle axes, then \(M_i=I_i\Omega_i (i=1,2,3)\), and \(I_{ik}=0 (i\neq k)\), so

\[ \begin{align*} I_1\frac{d\Omega_1}{dt}+(I_3-I_2)\Omega_2\Omega_3=K_1,\\ I_2\frac{d\Omega_2}{dt}+(I_1-I_3)\Omega_3\Omega_1=K_2,\\ I_3\frac{d\Omega_3}{dt}+(I_2-I_1)\Omega_1\Omega_2=K_3.\\ \end{align*} \]

Iterative Dynamic Formulation

Starting with link \(1\) and moving successively, link by link, outward to link \(n\).

Deduction: Outward

(i) Rotational velocity.

\[ {^{i+1}\omega_{i+1}}={^{i+1}_i R}{^i\omega_i}+\dot{\theta}_{i+1} {^{i+1}\hat{Z}_{i+1}}. \]

from equation \(\ref{angular-vel}\), we have

\[ \begin{align*} \dot{\omega}_{i+1}&=\dot{\omega}_{i}+{^0_iR} {^i \dot{\Omega}_{i+1}}+\omega_i\times {^0_iR}{^i\Omega_{i+1}}\\ &=\dot{\omega}_{i}+{^0_iR}(\ddot{\theta}_{i+1} {^i_{i+1}R}{^{i+1}\hat{Z}}_{i+1})+\omega_i\times {^0_iR}(\dot{\theta}_{i+1} {^i_{i+1}R}{^{i+1}\hat{Z}}_{i+1})\\ &=\dot{\omega}_{i}+\ddot{\theta}_{i+1}{^0_{i+1}R}{^{i+1}\hat{Z}_{i+1}}+\omega_i\times (\dot{\theta}_{i+1} {^0_{i+1}R}{^{i+1}\hat{Z}_{i+1}})\\ \Rightarrow \quad {^{i+1}\dot{\omega}_{i+1}}&={^{i+1}_i R}{^i\dot{\omega}_{i}}+\ddot{\theta}_{i+1}{^{i+1}\hat{Z}_{i+1}}+{^{i+1}_i R}{^i\omega_i}\times \dot{\theta}_{i+1}{^{i+1}\hat{Z}_{i+1}}. \end{align*} \]

For prismatic, we have \(\dot{\theta}_{i+1}=\ddot{\theta}_{i+1}=0\), so we only have the first item.

(ii) Linear velocity.

from equation \(\ref{linear-vel-general}\), we have

\[ \begin{align*} \dot{v}_{i+1}&=\dot{v}_i + {^0_iR} {^i\dot{V}_{i+1}}+2\omega_i\times {^0_iR}{^iV_{i+1}} + \dot{\omega}_i\times {^0_iR} {^iO_{i+1}}+\omega_i\times (\omega_i \times {^0_iR}{^iO_{i+1}})\\ &={^{i+1}_i R}[{^i\dot{v}_i}+{^i\dot{\omega}_i}\times {^iO_{i+1}}+{^i\omega_i}\times ({^i\omega_i} \times {^iO_{i+1}})]+{^{i+1}_i R}[{^i\dot{V}_{i+1}}+2{^i\omega_i}\times {^iV_{i+1}}]. \end{align*} \]

For prismatic, we have the second item with its corresonding change

\[ \begin{align*} {^{i+1}_i R}[{^i\dot{V}_{i+1}}+2{^i\omega_i}\times {^iV_{i+1}}]&= \ddot{d}_{i+1}{^{i+1}\hat{Z}_{i+1}}+{^{i+1}_i R}{^i\omega_i}\times \dot{d}_{i+1}{^{i+1}\hat{Z}_{i+1}}\\ &=\ddot{d}_{i+1}{^{i+1}\hat{Z}_{i+1}}+ {^{i+1}\omega_{i+1}}\times \dot{d}_{i+1}{^{i+1}\hat{Z}_{i+1}}. \end{align*} \]

For revolute, \(^iV_{i+1}={^i\dot{V}_{i+1}}=0\), so we only have the first item.

(iii) Acceleration of center of mass of link

Still from equation \(\ref{linear-vel-general}\), we have the second and third item equal zero

\[ \begin{align*} \dot{v}_{C_i}&=\dot{v}_i + {^0_iR} {^i\dot{V}_{C_i}}+2\omega_i\times {^0_iR}{^iV_{C_i}} + \dot{\omega}_i\times {^0_iR} {^iO_{C_i}}+\omega_i\times (\omega_i \times {^0_iR}{^iO_{C_i}})\\ &=\dot{v}_i +\dot{\omega}_i\times {^0_iR} {^iO_{C_i}}+\omega_i\times (\omega_i \times {^0_iR}{^iO_{C_i}})\\ \Rightarrow \quad {^i\dot{v}_{C_i}}&={^i\dot{v}_i} +{^i\dot{\omega}_i}\times {^iO_{C_i}}+{^i\omega_i}\times ({^i\omega_i} \times {^iO_{C_i}}). \end{align*} \]

(iv) force and torque

From Newton-Euler's Equation, we have

\[ \begin{align*} F_i&=m\dot{v}_{C_i},\\ {^{C_i}N_i}&={^{C_i}I}{^{C_i}\dot{\omega}_i}+ {^{C_i}{\omega}_i}\times {^{C_i}I}{^{C_i}{\omega}_i}\\ &={^{C_i}I}{^{i}\dot{\omega}_i}+ {^{i}{\omega}_i}\times {^{C_i}I}{^{i}{\omega}_i}. \end{align*} \]

Deduction: Inward

With forces and torques acting on links, we have to calculate joint torques that would result in these net forces and torques on the links. Here we denote \(f_i\) as force exerted on link \(i\) by link \(i-1\), \(n_i\) as torque exerted on link \(i\) by link \(i-1\).

Apparently net force on link \(i\) is

\[ \begin{align*} {^iF_i}&={^if_i}-{^i_{i+1}R}{^{i+1}f_{i+1}}.\\ \Rightarrow \quad {^if_i}&={^iF_i}+{^i_{i+1}R}{^{i+1}f_{i+1}}. \end{align*} \]

and for torque we substitute the above force relationship in

\[ \begin{align*} {^iN_i}&={^in_i}-{^in_{i+1}} +(-{^iO_{C_i}})\times {^if_i}+({^iO_{i+1}}-{^iO_{C_i}})\times (-{^if_{i+1}})\\ \Rightarrow {^in_i} &={^iN_i}+{^i_{i+1}R}{^{i+1}n_{i+1}}+{^iO_{C_i}}\times {^if_i}+({^iO_{i+1}}-{^iO_{C_i}})\times {^if_{i+1}}\\ &={^iN_i}+{^i_{i+1}R}{^{i+1}n_{i+1}}+{^iO_{C_i}}\times {^iF_i}+{^iO_{i+1}}\times {^i_{i+1}R}{^{i+1}f_{i+1}}. \end{align*} \]

For revolute, \(\tau_i={^in_i^T}{^i\hat{Z}_i}\), and for prismatic, \(\tau_i={^if_i^T}{^i\hat{Z}_i}\).

Structure of a manipulator's dynamic equations

State-space equation

The above Newton-Euler equations apply for any manipulators, and this yields a dynamic equation which could be written in the following form

\[ \begin{align} \tau=M(\Theta)\ddot{\Theta}+V(\Theta,\dot{\Theta})\dot{\Theta}+G(\Theta),\label{state-space} \end{align} \]

or for every row with \(\Theta=(\theta_1,\cdots,\theta_n)^T\),

\[ \tau_i=\sum_{j=1}^N m_{ij} \ddot{\theta_j}+\sum_{j,k=1}^Nc_{kji}\dot{\theta}_{k}\dot{\theta}_{j}-\sum_{j=1}^N m_j {^0\pmb{g}^T}\frac{\partial O_{C_j}}{\partial \theta_i} \]

where

\[ c_{kji}=\frac{1}{2}\left(\frac{\partial m_{ij}}{\partial \theta_k} + \frac{\partial m_{ik}}{\partial \theta_j}-\frac{\partial m_{kj}}{\partial \theta_i} \right)=c_{jki} \]

are called the first Christoffel sign, \(j\) and \(k\) are of equal state.

equation \(\ref{state-space}\) are called state-space equation, \(M(\Theta)_{N\times N}\) is called mass matrix of the manipulator, \(V(\Theta,\dot{\Theta})_{n\times 1}\) is called vector of centrifugal and Coriolis terms, which is where the main idea of the formula come from.

Configuration-space equation

Rewrite the velocity-dependent term \(V(\Theta, \dot{\Theta})\), we have

\[ \begin{align} \tau=M(\Theta)\ddot{\Theta}+B(\Theta)[\dot{\Theta}\dot{\Theta}]+C(\Theta)\dot{\Theta}^2+G(\Theta),\label{configuration-space} \end{align} \]

where \(B(\Theta)_{n\times n(n-1)/2}\) is a matrix of Coriolis coefficients, \([\dot{\Theta}\dot{\Theta}]_{n(n-1)/2\times 1}\) is a vector of joint velocity products given by

\[ [\dot{\Theta}\dot{\Theta}]=[\dot{\theta}_1\dot{\theta}_2, \dot{\theta}_1\dot{\theta}_3,\cdots,\dot{\theta}_{n-1}\dot{\theta}_n]^T. \]

\(C(\Theta)_{n\times n}\) is an matrix of centrifugal coefficients, and \(\dot{\Theta}^2_{n\times 1}\) is a vector given by

\[ \dot{\Theta}^2=[\dot{\theta}_1^2, \dot{\theta}_2^2,\cdots, \dot{\theta}_n^2]. \]

Equation \(\ref{configuration-space}\) is called configuration-space equation because of \(B\) and \(C\).