# Center of Mass and Moment of Inertia¶

Renato Naville Watanabe

In :
import numpy as np
import sympy as sym


### Center of mass¶

The center of mass of a set of particles is formally defined as the point where the sum of the vectors linking this point to each particle, weighted by its mass, is zero.

#### Example¶

Two masses with mass $m_1$ and $m_2$

\begin{equation} \vec{\bf{r_{cm}}} = \frac{m_1\vec{\bf{r_1}}+m_2\vec{\bf{r_2}}}{m_1+m_2} = \frac{m_1\vec{\bf{r_1}}+m_2\vec{\bf{r_1}}-m_2\vec{\bf{r_1}} +m_2\vec{\bf{r_2}}}{m_1+m_2}=\frac{(m_1+m_2)\vec{\bf{r_1}}}{m_1+m_2} + \frac{m_2(\vec{\bf{r_2}}-\vec{\bf{r_1}})}{m_1+m_2}=\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, = \vec{\bf{r_1}} + \frac{m_2}{m_1+m_2}(\vec{\bf{r_2}}-\vec{\bf{r_1}}) \end{equation}

For a rigid body, the center of mass is defined as the point where the integral of the vectors linking this point to each differential part of mass, weighted by this differential mass, is zero.

#### Example¶

bar

Now, we will compute the center of mass of the triangle below.

In :
x,y = sym.symbols('x,y')
rho, h, p, b = sym.symbols('rho, h, p, b')

In :
xcm = 2*(sym.integrate(sym.integrate(x, (y, 0, x*h/p)), (x, 0, p)) +
sym.integrate(sym.integrate(x, (y, 0, (b-x)*h/(b-p))), (x, p, b)))/(b*h)
print('x_cm = ',sym.simplify(xcm))

x_cm =  b/3 + p/3

In :
ycm = 2*(sym.integrate(sym.integrate(y, (y, 0, x*h/p)), (x, 0, p)) +
sym.integrate(sym.integrate(y, (y, 0, (b-x)*h/(b-p))), (x, p, b)))/(b*h)
print('y_cm = ',sym.simplify(ycm))

y_cm =  h/3


So, the center of mass of the triangle is $\vec{\bf{r_{cm}}} = \left(\frac{b+p}{3}\right)\hat{\bf{i}} + \frac{h}{3}\hat{\bf{j}}$

### Center of gravity¶

Center of gravity of a body is the point where the moment caused by the gravitational force in the whole body relative to this point is zero.

\begin{equation} \vec{\bf{M_0}} = \int\limits_{B} \vec{\bf{r_B}}\, \times\,\vec{\bf{g}}\,dm \end{equation}

If the acceleration of gravity being applied to the whole body is the same (for all practical purposes, in Biomechanics we can consider it the same in the whole body), the gravity vector can go out of the integral:

\begin{equation} \vec{\bf{M_0}}= \int\limits_{B}\vec{\bf{r_B}}\,dm\, \times\,\vec{\bf{g}} = \int\limits_{B}(\vec{\bf{r}} - \vec{\bf{r_{cm}}})\,dm\, \times\,\vec{\bf{g}} = \left(\int\limits_{B}\vec{\bf{r}}\,dm -\int\limits_{B}\vec{\bf{r_{G}}}\,dm\,\right) \times\,\vec{\bf{g}} \end{equation}

Now, we equal this moment to zero and isolate $\vec{\bf{r_{G}}}$:

\begin{equation} \left(\int\limits_{B}\vec{\bf{r}}\,dm -\int\limits_{B}\vec{\bf{r_{G}}}\,dm\right) \times\,\vec{\bf{g}} = 0 \longrightarrow \int\limits_{B}\vec{\bf{r}}\,dm -\int\limits_{B}\vec{\bf{r_{G}}}\,dm = 0 \longrightarrow \int\limits_{B}\vec{\bf{r}}\,dm -\vec{\bf{r_{G}}}\int\limits_{B}\,dm = 0\,\,\,\,\,\, \longrightarrow \vec{\bf{r_{G}}} = \frac{ \int\limits_{B}\vec{\bf{r}}\,dm}{\int\limits_{B}\,dm} = \frac{ \int\limits_{B}\vec{\bf{r}}\,dm}{m_B} \end{equation}

where $m_B$ is the mass of the body.

Note that in this case, that the gravity acceleration is constant, the center of gravity $\vec{\bf{r_{G}}}$ is equal to the center of mass $\vec{\bf{r_{cm}}}$.

### Geometric center¶

For a rigid body, the geometric center is defined as the point where the integral of the vectors linking this point to each differential part of mass is zero. It

This integral is:

\begin{equation} \int\limits_{B} {\bf{\vec{r}_{/gc}}} dV = \int\limits_{B} {\bf{\vec{r}}} - {\bf{\vec{r}_{gc}}} dV = \int\limits_{B} {\bf{\vec{r}}}\,dV - \int\limits_{B}{\bf{\vec{r}_{gc}}} dV = \int\limits_{B} {\bf{\vec{r}}}\,dV - {\bf{\vec{r}_{gc}}}\int\limits_{B}\, dV \end{equation}

Now we equal this integral to zero and isolate ${\bf\vec{r_{gc}}}$:

\begin{equation} \int\limits_{B} {\bf{\vec{r}}}\,dV - {\bf{\vec{r}_{gc}}}\int\limits_{B}\, dV = 0 \longrightarrow {\bf{\vec{r}_{gc}}} = \frac{ \int\limits_{B} {\bf{\vec{r}}}\,dV}{\int\limits_{B}\, dV} = \frac{ \int\limits_{B} {\bf{\vec{r}}}\,dV}{V} \label{eq:gc} \end{equation}

where $V$ is the volume of the body. Note that when the body has a constant density $\rho$, the center of mass is equal to the geometric center:

\begin{equation} {\bf{\vec{r}_{cm}}} = \frac{ \int\limits_{B} {\bf{\vec{r}}}\,dm}{m_B} = \frac{ \int\limits_{B} {\bf{\vec{r}}}\rho\,dV}{\rho V} = \frac{ \rho\int\limits_{B} {\bf{\vec{r}}}\,dV}{\rho V} = \frac{ \int\limits_{B} {\bf{\vec{r}}}\,dV}{V} = {\bf{\vec{r}_{gc}}} \end{equation}

### Center of mass of a complex system¶

Now, we will consider a set of $n$ bodies. The center of mass of this set of bodies can be computed by integrating the Eq.~\eqref{eq:cm} over all bodies:

\begin{equation} \vec{\bf{r_{G}}} = \frac{ \int\limits_{B1,B2,...,Bn}\vec{\bf{r}}\,dm}{m_{B1}+m_{B2}+... + m_{Bn}} =\frac{\int\limits_{B1}\vec{\bf{r}}\,dm+\int\limits_{B2}\vec{\bf{r}}\,dm+...+\int\limits_{Bn}\vec{\bf{r}}\,dm}{m_{B1}+m_{B2}+... + m_{Bn}} = \frac{\frac{\int\limits_{B1}\vec{\bf{r}}\,dm}{m_{B1}}m_{B1} + \frac{\int\limits_{B2}\vec{\bf{r}}\,dm}{m_{B2}}m_{B2} +...+ \frac{\int\limits_{Bn}\vec{\bf{r}}\,dm}{m_{Bn}}m_{Bn}}{m_{B1}+m_{B2}+... + m_{Bn}} = \frac{\vec{\bf{r_{cm_1}}} m_{B1} + \vec{\bf{r_{cm_2}}}m_{B2} +...+ \vec{\bf{r_{cm_n}}}m_{Bn}}{m_{B1}+m_{B2}+... + m_{Bn}}=\frac{\vec{\bf{r_{cm_1}}} m_{B1} + \vec{\bf{r_{cm_2}}}m_{B2} +...+ \vec{\bf{r_{cm_n}}}m_{Bn}}{m_T} \end{equation}

where $\vec{\bf{r_{cm_i}}}$ is the center of mass of the body $i$, $m_{Bi}$ is the mass of the body $i$ and $m_T$ is the total mass of the bodies. The expression above shows that we can interpret each body as a particle with its mass and position.

Now we will compute the center of mass of the system shown below.

### Moment of inertia¶

The moment of inertia is a measure of how the mass of the body distributes relative to a given axis passing to a given point O. In other words, the moment of inertia relative to a point O for an axis a is defined how far the mass of the body is from this axis. If we define the versor $\hat{\bf{n_a}}$ as a versor in the direction of the axis, this distance is the croo product between the vector linking the point O to each part of the body and this versor.

Now we will compute the moment of inertia relative to the center of mass of the same bar we computed the center of mass.

The moment of inertia relative to the center of mass is:

\begin{equation} I_{zz}^{cm} = \int\limits_B x_{/cm}^2+y_{/cm}^2\,dm = \int\limits_{-l/2}^{l/2} x^2\rho\,dx = \rho\left.\left(\frac{x^3}{3}\right)\right|_{-l/2}^{l/2} = \rho\left(\frac{l^3}{24}+\frac{l^3}{24}\right)=\rho\frac{l^3}{12} = \frac{ml^2}{12} \end{equation}

Radius of gyration is defined as the distance that a particle would be from the a point O to have the same moment of inertia that the body have. So, the radius of gyration is defined as:

\begin{equation} r_{gyr} = \sqrt{\frac{I_{zz}^{cm}}{m}} \end{equation}

For example, for the homogenous bar with length $l$, the radius of gyration is:

\begin{equation} r_{gyr} = \sqrt{\frac{\frac{ml^2}{12}}{m}}=\sqrt{\frac{l^2}{12}} = \frac{l\sqrt{3}}{6} \end{equation}

### Parallel axis theorem¶

If we have computed the moment of inertia relative to the center of mass for an axis, for example the z-axis, and now want to compute the moment of inertia relative to another point O for an axis parallel to the first, there is an expression to aid the computation of this moment of inertia.

In the figure below, the axis is perpendicular to the figure plane.

Now we will compute the moment of inertia the moment of inertia of a bar relative to one of its extremities, by using the parallel axis theorem.

\begin{equation} I_{zz}^O = I_{zz}^{cm} + m\left(\frac{l}{2}\right)^2 = \frac{ml^2}{12} + \frac{ml^2}{4} = \frac{ml^2}{3} \end{equation}

### Moment of inertia of a complex system¶

To compute the moment of inertia of a set o $n$ bodies relative to a point O for a given axis (for example, the z axis), we must apply the Eq.~\eqref{eq:inmomzzo} for all bodies:

\begin{equation} I_{zz}^{O} = \int\limits_{B1,B2,...,Bn} x_{/O}^2+y_{/O}^2\,dm = \int\limits_{B1} x_{/O}^2+y_{/O}^2\,dm +\int\limits_{B2} x_{/O}^2+y_{/O}^2\,dm +...+\int\limits_{Bn} x_{/O}^2+y_{/O}^2\,dm \end{equation}

Now, using the parallel axis theorem, we can write the equation above in terms of the moment of inertia relative to the center of mass of each body:

\begin{equation} I_{zz}^{O} = I_{zz_{B1}}^{cm_1} + m_{B1}.||\vec{\bf{r_{cm_1/O}}}||^2 +I_{zz_{B2}}^{cm_2} + m_{B2}.||\vec{\bf{r_{cm_2/O}}}||^2 +...+I_{zz_{Bn}}^{cm_n} + m_{Bn}.||\vec{\bf{r_{cm_n/O}}}||^2 \label{eq:paral} \end{equation}

where $I_{zz_{Bi}}^{cm_i}$ is the moment of inertia of the body $i$ relative to its center of mass.

Now we will compute the moment of inertia of the set of eight equal and homogenous bars depicted in the figure below relative to its center of mass.

### Matrix of Inertia¶

For 3D movements we must consider the moments of inertia inertia in the three directions besides three other quantities called products of inertia. A product of inertia is defined as in Eq.~\eqref{eq:inmomgen}, but with each versor orthogonal to each other. So, the product of inertia relative to point O for axes a and b is:

\begin{equation} I_{ab}^O= \int\limits_B (\vec{\bf{r_{/O}}} \times \hat{\bf{n_a}}).(\vec{\bf{r_{/O}}} \times \hat{\bf{n_b}})\,dm \label{eq:prodinmomgen} \end{equation}

Note that $I_{ab}^O$ is equal to $I_{ba}^O$.

The matrix of inertia is defined as:

\begin{equation} I^O = \left[\begin{array}{ccc}I_{xx}^O&I_{xy}^O&I_{xz}^O\\ I_{yx}^O&I_{yy}^O&I_{yz}^O\\ I_{zx}^O&I_{zy}^O&I_{zz}^O\end{array}\right] \end{equation}

where each term is defined as follows:

\begin{equation} I_{xx}^O = \int\limits_B (\vec{\bf{r_{/O}}} \times \hat{\bf{i}}).(\vec{\bf{r_{/O}}} \times \hat{\bf{i}})\,dm = \int\limits_B y_{/O}^2+z_{/O}^2\,dm \end{equation}\begin{equation} I_{yy}^O = \int\limits_B (\vec{\bf{r_{/O}}} \times \hat{\bf{j}}).(\vec{\bf{r_{/O}}} \times \hat{\bf{j}})\,dm = \int\limits_B x_{/O}^2+z_{/O}^2\,dm \end{equation}\begin{equation} I_{zz}^O = \int\limits_B (\vec{\bf{r_{/O}}} \times \hat{\bf{k}}).(\vec{\bf{r_{/O}}} \times \hat{\bf{k}})\,dm = \int\limits_B x_{/O}^2+y_{/O}^2\,dm \end{equation}\begin{equation} I_{xy}^O = I_{yx}^O =\int\limits_B (\vec{\bf{r_{/O}}} \times \hat{\bf{i}}).(\vec{\bf{r_{/O}}} \times \hat{\bf{j}})\,dm = \int\limits_B \left[(x_{/O}\hat{\bf{i}}+y_{/O}\hat{\bf{j}}+z_{/O}\hat{\bf{k}}) \times \hat{\bf{i}}\right].\left[(x_{/O}\hat{\bf{i}}+y_{/O}\hat{\bf{j}}+z_{/O}\hat{\bf{k}}) \times \hat{\bf{j}}\right]\,dm = \int\limits_B (-y_{/O}\hat{\bf{k}}+z_{/O}\hat{\bf{j}}).(x_{/O}\hat{\bf{k}}-z_{/O}\hat{\bf{i}})\,dm \longrightarrow I_{xy}^O = \int\limits_B - y_{/O}x_{/O}\,dm \end{equation}\begin{equation} I_{xz}^O = I_{zx}^O =\int\limits_B (\vec{\bf{r_{/O}}} \times \hat{\bf{i}}).(\vec{\bf{r_{/O}}} \times \hat{\bf{k}})\,dm = \int\limits_B \left[(x_{/O}\hat{\bf{i}}+y_{/O}\hat{\bf{j}}+z_{/O}\hat{\bf{k}}) \times \hat{\bf{i}}\right].\left[(x_{/O}\hat{\bf{i}}+y_{/O}\hat{\bf{j}}+z_{/O}\hat{\bf{k}}) \times \hat{\bf{k}}\right]\,dm = \int\limits_B (-y_{/O}\hat{\bf{k}}+z_{/O}\hat{\bf{j}}).(-x_{/O}\hat{\bf{j}}+y_{/O}\hat{\bf{i}})\,dm \longrightarrow I_{xz}^O = \int\limits_B - z_{/O}x_{/O}\,dm \end{equation}\begin{equation} I_{yz}^O = I_{zy}^O =\int\limits_B (\vec{\bf{r_{/O}}} \times \hat{\bf{j}}).(\vec{\bf{r_{/O}}} \times \hat{\bf{k}})\,dm = \int\limits_B \left[(x_{/O}\hat{\bf{i}}+y_{/O}\hat{\bf{j}}+z_{/O}\hat{\bf{k}}) \times \hat{\bf{j}}\right].\left[(x_{/O}\hat{\bf{i}}+y_{/O}\hat{\bf{j}}+z_{/O}\hat{\bf{k}}) \times \hat{\bf{k}}\right]\,dm = \int\limits_B (x_{/O}\hat{\bf{k}}-z_{/O}\hat{\bf{i}}).(-x_{/O}\hat{\bf{j}}+y_{/O}\hat{\bf{i}})\,dm \longrightarrow I_{yz}^O = \int\limits_B - z_{/O}y_{/O}\,dm \end{equation}

Now, we will consider a cylinder and compute the moment of inertia for the three axes, relative to the center of mass. The mass of this cylinder is $m = \rho \pi R^2$ and, by symmetry, it is in the middle of the cylinder.

In :
h, R, z, theta, r, rho, m1 = sym.symbols('h, R, z, theta, r, rho, m')
m = rho*h*sym.pi*R**2
Ixx = rho*sym.integrate(sym.integrate(sym.integrate(r*z**2+(r**3)*sym.sin(theta)**2, (theta, 0, 2*sym.pi)), (r, 0, R)), (z,-h/2,h/2))/m
sym.simplify(m1*Ixx)

Out:
$$\frac{m}{12} \left(3 R^{2} + h^{2}\right)$$

So, the moment of inertia of the cylinder for the axis $x$ is:

$I_{xx}^{cm} = \frac{m}{12}(3R^2+h^2)$

By symmetry, the moment of inertia for the axis $y$ is equal to the moment of inertia for the axis $x$.

$I_{yy}^{cm} = \frac{m}{12}(3R^2+h^2)$

#### Principal axes¶

It is always possible to choose the direction of the axes in a way that the products of inertia are zero. In this situation, the inertia matrix is:

\begin{equation} I^{cm} = \left[\begin{array}{ccc}I_{1}^O&0&0\\ 0&I_{2}^O&0\\ 0&0&I_{3}^O\end{array}\right] \end{equation}

where $I_1$, $I_2$ and $I_3$ are the moments of inertia for each of the principal axes.

### Problems¶

Solve the problems 3.2.5, 3.2.7, 15.3.11, 15.3.15 from Ruina and Pratap book.

In [ ]: