To design efficient feedback control systems, we use mathematical models to understand the system, and then we design controllers that achieve the desired specifications. Developing suitable mathematical models is problem-specific and can be very complicated, requiring specialized knowledge. However, there are some well-known techniques for certain classes of systems. This book covers one general technique that is useful for a large class of rigid-body mechanical systems. The technique is known as the Euler-Lagrange method and requires a mathematical description of the kinetic and potential energy of the system, given its generalized position and velocity. The theory behind the Euler-Lagrange method is deep and well beyond the scope of this book. However, we will give a brief introduction that will enable the reader to apply Euler-Lagrange methods to a surprisingly large class of problems.

As motivation for the Euler-Lagrange equations, consider a particle of mass \(m\) in a gravity field, as shown in Fig. 1.5, where \(g\) is the force of gravity at sea level, \(f_a\) is an externally applied force, and the friction due to air resistance is modeled by a gain \(b\) times the velocity \(\dot{y}\).

Since all of the forces are in the vertical direction, the motion of the vehicle can be completely described by the vertical position of particle \(y(t)\). From Newton’s law, the equations of motion are \[\tag{1.1} m\ddot{y} = f_a-mg-b\dot{y},\] where dots over the variable will be used throughout the book to denote differentiation with respect to time, i.e., \(\dot{y}\stackrel{\triangle}{=}\frac{dy}{dt}\) and \(\ddot{y}\stackrel{\triangle}{=}\frac{d^2y}{dt^2}\). Rearranging 1.1 gives \[\tag{1.2} m\ddot{y} + mg = f_a -b\dot{y}.\] The kinetic energy for the particle is given by \(K(\dot{y})=\frac{1}{2}mv^2=\frac{1}{2}m\dot{y}^2\). The potential energy, which is induced by gravitational pull proportional to the distance from the center of the earth, is given by \(P(y)=P_0+mgy\), where \(P_0\) is the potential energy of the particle when \(y=0\). Note that if the velocity \(\dot{y}\) is thought of as independent of \(y\), then \[\frac{\partial K(\dot{y})}{\partial \dot{y}} = m\dot{y}\] and \[\frac{d}{dt}\left(\frac{\partial K(\dot{y})}{\partial \dot{y}}\right) = m\ddot{y}.\] Also note that \[\frac{\partial P(y)}{\partial y} = mg.\] Therefore, Equation 1.2 can be written as \[\tag{1.3} \frac{d}{dt}\left(\frac{\partial K(\dot{y})}{\partial \dot{y}}\right) + \frac{\partial P(y)}{\partial y} = f_a-b\dot{y}.\] Equation 1.3 is a particular form of the so-called Euler-Lagrange equation. The generalization to a more complicated system takes the form \[\tag{1.4} \frac{d}{dt}\left(\frac{\partial K(y,\dot{y})}{\partial \dot{y}}-\frac{\partial P(y)}{\partial \dot{y}}\right) - \left(\frac{\partial K(y,\dot{y})}{\partial y} - \frac{\partial P(y)}{\partial y}\right) = f_a-b\dot{y},\] which reduces to Equation 1.3 for our particular case since \(K\) is only a function of \(\dot{y}\). To simplify the notation, we define the so-called Lagrangian as \[L(y,\dot{y})\stackrel{\triangle}{=}K(y,\dot{y})-P(y),\] and write Equation 1.4 as \[\tag{1.5} \frac{d}{dt}\left(\frac{\partial L(y,\dot{y})}{\partial \dot{y}}\right) - \left(\frac{\partial L(y,\dot{y})}{\partial y}\right) = f_a-b\dot{y}.\] The position of the system can be generalized to both positions and angles, and the velocity can be generalized to both linear and angular velocities. Throughout the book we will use the notation \(\mathbf{q} = (q_1, q_2, \dots, q_n)^\top\) to denote the general configuration of the mechanical system that we are modeling. The vector \(\mathbf{q}\) is called the generalized-coordinate vector. The kinetic energy of the system is generally a function of both the generalized coordinates \(\mathbf{q}\) and the generalized velocity \(\dot{\mathbf{q}}\). Therefore, we denote the kinetic energy as \(K(\mathbf{q},\dot{\mathbf{q}})\). On the other hand, in this book we will only consider potential fields that are conservative and depend only on the configuration variable \(\mathbf{q}\). Under these assumptions, the general form of the Euler-Lagrange equations is given by \begin{align} \frac{d}{dt}\left(\frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial \dot{q}_1}\right) - \frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial q_1} &= \tau_1-b_{11} \dot{q}_1 - \cdots - b_{1n}\dot{q}_n \notag \\ \frac{d}{dt}\left(\frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial \dot{q}_2}\right) - \frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial q_2} &= \tau_2-b_{21} \dot{q}_1 - \cdots - b_{2n}\dot{q}_n \tag{1.6} \\ \vdots \\ \frac{d}{dt}\left(\frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial \dot{q}_n}\right) - \frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial q_n} &= \tau_n-b_{n1} \dot{q}_1 - \cdots - b_{nn}\dot{q}_n, \notag \end{align} where \(\tau_i\) are the generalized forces and torques that are applied to the system. Or in vector form we can write \[\tag{1.8} \frac{d}{dt}\left(\frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial \dot{\mathbf{q}}}\right) - \frac{\partial L(\mathbf{q},\dot{\mathbf{q}})}{\partial \mathbf{q}} = \boldsymbol{\tau} -B \dot{\mathbf{q}},\] where \(B\) is a matrix with positive elements and \(\boldsymbol{\tau}=(\tau_1,\dots,\tau_n)^\top\) are the generalized forces. The Euler-Lagrange equation will be discussed in more detail in Chapter 6. However, before showing how to use Equation 1.8 to derive the equations of motion for practical systems, we need to show how to compute the kinetic and potential energy for general (holonomic) rigid body systems.

2 Kinetic Energy of Mechanical Systems

2.1 Theory

Consider a point mass with mass \(m\), traveling at velocity \(\mathbf{v}\), as shown in Fig. 2.1. If the velocity is resolved with respect to a coordinate axis, then \(\mathbf{v}=(v_x, v_y, v_z)^\top\). The kinetic energy of the point mass is given by \[\begin{aligned} K &= \frac{1}{2} m \left|\mathbf{v}\right|^2 \\ &= \frac{1}{2} m \mathbf{v}^\top \mathbf{v} \\ & = \frac{1}{2} m (v_x^2 + v_y^2 + v_z^2). \end{aligned}\]

Now consider two point masses as shown in Fig. 2.2.

The total kinetic energy of the system is given by summing the kinetic energy of each of the point masses as \[\tag{2.1} K = \frac{1}{2} m_1 \mathbf{v}_1^\top \mathbf{v}_1 + \frac{1}{2} m_2 \mathbf{v}_2^\top \mathbf{v}_2.\] If the two masses are connected by a massless rod and the particles are spinning about a point on the rod, as shown in Fig. 2.2, where the angular velocity vector \(\boldsymbol{\omega}\) points out of the page, and the magnitude is given by \(\dot{\theta}\), i.e., \(\boldsymbol{\omega}=(0, 0, \dot{\theta})^\top\), then from elementary physics, the velocity of point \(i\) is given by \[\mathbf{v}_i = \boldsymbol{\omega}\times\mathbf{p}_i,\] where in our case \(\mathbf{p}_1 = (L_1\cos\theta, L_1\sin\theta, 0)^\top\), and \(\mathbf{p}_2 = (-L_2\cos\theta, -L_2\sin\theta, 0)^\top\). Therefore \[\begin{aligned} \mathbf{v}_1 &= \begin{pmatrix} 0 \\ 0 \\ \dot{\theta}\end{pmatrix}\times\begin{pmatrix}L_1\cos\theta \\ L_1\sin\theta \\ 0\end{pmatrix} = \begin{pmatrix} -L_1\dot{\theta}\sin\theta \\ L_1\dot{\theta}\cos\theta \\ 0 \end{pmatrix} \\ \mathbf{v}_2 &= \begin{pmatrix} 0 \\ 0 \\ \dot{\theta} \end{pmatrix}\times\begin{pmatrix}-L_2\cos\theta \\ -L_2\sin\theta \\ 0\end{pmatrix} = \begin{pmatrix} L_2\dot{\theta}\sin\theta \\ -L_2\dot{\theta}\cos\theta \\ 0 \end{pmatrix}. \end{aligned}\] An alternative method for computing the velocities \(\mathbf{v}_1\) and \(\mathbf{v}_2\) is to differentiate \(\mathbf{p}_1\) and \(\mathbf{p}_2\) with respect to time to obtain \[\begin{aligned} \mathbf{v}_1 &\stackrel{\triangle}{=}\dot{\mathbf{p}}_1 = \frac{d\mathbf{p}_1}{dt} = \frac{d}{dt}(L_1\cos\theta, L_1\sin\theta, 0)^\top \\ &= \begin{pmatrix}-L_1\left(\frac{d\theta}{dt}\right)\sin\theta, & L_1\left(\frac{d\theta}{dt}\right)\cos\theta, & 0 \end{pmatrix}^\top \\ &= \begin{pmatrix}-L_1\dot{\theta}\sin\theta, & L_1\dot{\theta}\cos\theta, & 0 \end{pmatrix}^\top \\ \mathbf{v}_2 &\stackrel{\triangle}{=}\dot{\mathbf{p}}_2 = \frac{d\mathbf{p}_2}{dt} = \frac{d}{dt}(-L_2\cos\theta, -L_2\sin\theta, 0)^\top \\ &= \begin{pmatrix}L_2\left(\frac{d\theta}{dt}\right)\sin\theta, & -L_2\left(\frac{d\theta}{dt}\right)\cos\theta, & 0 \end{pmatrix}^\top \\ &= \begin{pmatrix}L_2\dot{\theta}\sin\theta, & -L_2\dot{\theta}\cos\theta, & 0 \end{pmatrix}^\top. \end{aligned}\]

Therefore, using Equation 2.1 the kinetic energy is \[\begin{aligned} K &= \frac{1}{2} m_1 \mathbf{v}_1^\top \mathbf{v}_1 + \frac{1}{2} m_2 \mathbf{v}_2^\top \mathbf{v}_2 \\ &= \frac{1}{2} m_1 \left[ (-L_1\dot{\theta}\sin\theta)^2 + (L_1\dot{\theta}\cos\theta)^2 \right] \\ &\qquad + \frac{1}{2} m_2 \left[ (L_2\dot{\theta}\sin\theta)^2 + (-L_2\dot{\theta}\cos\theta)^2 \right] \\ &= \frac{1}{2} m_1 \left[ L_1^2\dot{\theta}^2\sin^2\theta + L_1^2\dot{\theta}^2\cos^2\theta \right] % + \frac{1}{2} m_2 \left[ L_2^2\dot{\theta}^2\sin^2\theta + L_2^2\dot{\theta}^2\cos^2\theta \right] \\ &= \frac{1}{2} (m_1 L_1^2 + m_2 L_2^2)\dot{\theta}^2. \end{aligned}\]

Now consider the case where there are \(N\) point masses spinning about a general angular velocity vector \(\boldsymbol{\omega}\), where the position of the \(i^{th}\) point mass, \(i=1,\dots,N\), in a specified coordinate system is given by \(\mathbf{p}_i\), as shown in Fig. 2.3.

We will see in the discussion below that it is most convenient to pick the coordinate system to be at the center of mass. The kinetic energy of the \(i^{th}\) particle is \(K_i = \frac{1}{2}m_i\mathbf{v}_i^\top\mathbf{v}_i\). Therefore, the kinetic energy of the total system is given by \[K = \sum_{i=1}^N \frac{1}{2}m_i\mathbf{v}_i^\top\mathbf{v}_i.\] Using the fact that \(\mathbf{v}_i = \boldsymbol{\omega}\times\mathbf{p}_i\) gives \[K = \sum_{i=1}^N \frac{1}{2}m_i\left(\boldsymbol{\omega}\times\mathbf{p}_i\right)^\top\left(\boldsymbol{\omega}\times\mathbf{p}_i\right).\] Recalling the cross product rule \[(\mathbf{a}\times\mathbf{b})^\top(\mathbf{c}\times\mathbf{d}) = (\mathbf{a}^\top\mathbf{c})(\mathbf{b}^\top\mathbf{d})-(\mathbf{a}^\top\mathbf{d})(\mathbf{b}^\top\mathbf{c})\] gives \[\begin{aligned} K &= \frac{1}{2} \sum_{i=1}^N m_i \left[ (\boldsymbol{\omega}^\top\boldsymbol{\omega})(\mathbf{p}_i^\top\mathbf{p}_i)-(\boldsymbol{\omega}^\top\mathbf{p}_i)(\mathbf{p}_i^\top\boldsymbol{\omega})\right] \\ &= \frac{1}{2} \sum_{i=1}^N m_i \boldsymbol{\omega}^\top \left[ (\mathbf{p}_i^\top\mathbf{p}_i) I_3 - \mathbf{p}_i \mathbf{p}_i^\top \right] \boldsymbol{\omega} \\ &= \frac{1}{2} \boldsymbol{\omega}^\top \left[ \sum_{i=1}^N m_i \left((\mathbf{p}_i^\top\mathbf{p}_i) I_3 - \mathbf{p}_i \mathbf{p}_i^\top \right) \right] \boldsymbol{\omega} \\ &\stackrel{\triangle}{=}\frac{1}{2} \boldsymbol{\omega}^\top J \boldsymbol{\omega}. \end{aligned}\] where \(I_3\) is the \(3 \times 3\) identity matrix. The matrix \[J = \sum_{i=1}^N m_i \left((\mathbf{p}_i^\top\mathbf{p}_i) I_3 - \mathbf{p}_i \mathbf{p}_i^\top \right)\] is called the inertia matrix of the mass system. If the mass system is solid then there is an infinite number of particles, and the sum becomes an integral resulting in \[J = \int_{\mathbf{p}\in\text{body}} \left((\mathbf{p}^\top\mathbf{p}) I_3 - \mathbf{p} \mathbf{p}^\top \right)dm(\mathbf{p}),\] where \(dm(\mathbf{p})\) is the elemental mass located at position \(\mathbf{p}\) and the integral is over all elemental masses in the system.

Suppose that \(\mathbf{p}=(x, y, z)^\top\) is expressed with respect to a coordinate system fixed in a rigid body. Then \[\begin{aligned} J &= \int_{(x,y,z)\in\text{body}} \Bigg[ \begin{pmatrix} x & y & z \end{pmatrix}\begin{pmatrix} x \\ y \\ z \end{pmatrix} \begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix} % - \begin{pmatrix} x \\ y \\ z \end{pmatrix} \begin{pmatrix} x & y & z \end{pmatrix} \Bigg] dm(x,y,z) \\ &= \int_{(x,y,z)\in\text{body}} \Bigg[ (x^2 + y^2 + z^2)\begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix} % - \begin{pmatrix} x^2 & xy & x z \\ yx & y^2 & yz \\ zx & zy & z^2 \end{pmatrix} \Bigg] dm(x,y,z) \\ &= \int_{(x,y,z)\in\text{body}} \begin{pmatrix} y^2 + z^2 & -xy & -x z \\ -yx & x^2+z^2 & -yz \\ -zx & -zy & x^2+y^2 \end{pmatrix} dm(x,y,z). \end{aligned}\]

As an example, to compute the inertia matrix of the thin rod of mass \(m\) about its center of mass, as shown in Fig. 2.4, the first step is to define a coordinate system at the center of mass, as shown in Fig. 2.4, where \(\hat{i}\) and \(\hat{j}\) are orthogonal unit vectors pointing parallel and perpendicular to the rod, respectively. The unit vector \(\hat{k}\) forms a right-handed coordinate system and points out of the page.

Assuming that the rod is infinitely thin and that the mass is evenly distributed throughout the rod, then the elemental mass is given by \[dm(x) = \frac{m}{\ell}dx.\] Since the rod is infinitely thin in the \(\hat{j}\) and \(\hat{k}\) directions, the inertia matrix is given by \begin{align} J &= \int_{x=-\ell/2}^{\ell/2}\int_{y=0}^0\int_{z=0}^0 \begin{pmatrix} y^2 + z^2 & -xy & -x z \\ -yx & x^2+z^2 & -yz \\ -zx & -zy & x^2+y^2 \end{pmatrix} dz \, dy \, \frac{m}{\ell}dx \notag\\ &= \frac{m}{\ell} \int_{x=-\ell/2}^{\ell/2} \begin{pmatrix} 0 & 0 & 0 \\ 0 & x^2 & 0 \\ 0 & 0 & x^2 \end{pmatrix} dx \notag\\ &= \frac{m}{\ell} \begin{pmatrix} 0 & 0 & 0 \\ 0 & \frac{x^3}{3}\big|_{x=-\frac{\ell}{2}}^{\frac{\ell}{2}} & 0 \\ 0 & 0 & \frac{x^3}{3}\big|_{x=-\frac{\ell}{2}}^{\frac{\ell}{2}} \end{pmatrix} = \begin{pmatrix} 0 & 0 & 0 \\ 0 & \frac{m\ell^2}{12} & 0 \\ 0 & 0 & \frac{m\ell^2}{12} \end{pmatrix}. \tag{2.2} \end{align} Note that the infinitely thin rod assumption implies that there is only inertia, or rotational mass, about the \(\hat{j}\) and \(\hat{k}\) axes.

As another example, to compute the inertia matrix of the thin plate of mass \(m\) about its center of mass, as shown in Fig. 2.5, the first step is to define a coordinate system at the center of mass, as shown in Fig. 2.5, where \(\hat{i}\) and \(\hat{j}\) are orthogonal unit vectors pointing along the length and width of the plate, respectively, and \(\hat{k}\) forms a right-handed coordinate system and points perpendicular to the plate.

Assuming that the plate is infinitely thin, and that the mass is evenly distributed throughout the plate, then the elemental mass is given by \[dm(x,y) = \frac{m}{\ell w}dx\,dy.\] Since the plate is infinitely thin in the \(\hat{k}\) direction, the inertia matrix is given by \[\begin{aligned} J &= \int_{x=-\ell/2}^{\ell/2}\int_{y=-w/2}^{w/2}\int_{z=0}^0 \begin{pmatrix} y^2 + z^2 & -xy & -x z \\ -yx & x^2+z^2 & -yz \\ -zx & -zy & x^2+y^2 \end{pmatrix} dz \, dy \, \frac{m}{\ell w}\,dx \\ &= \frac{m}{\ell w} \int_{x=-\ell/2}^{\ell/2} \int_{y=-w/2}^{w/2} \begin{pmatrix} y^2 & -xy & 0 \\ -yx & x^2 & 0 \\ 0 & 0 & x^2 +y^2\end{pmatrix} dy\,dx \\ % &= \frac{m}{\ell w} \int_{x=-\ell/2}^{\ell/2} \begin{pmatrix} \frac{w^3}{12} & 0 & 0 \\ 0 & x^2 w & 0 \\ 0 & 0 & x^2 w +\frac{w^3}{12}\end{pmatrix} dx \\ % &= \frac{m}{\ell w} \begin{pmatrix} \frac{w^3}{12}\ell & 0 & 0 \\ 0 & \frac{\ell^3}{12} w & 0 \\ 0 & 0 & \frac{\ell^3}{12}w +\frac{w^3}{12}\ell\end{pmatrix} = \begin{pmatrix} \frac{mw^2}{12} & 0 & 0 \\ 0 & \frac{m\ell^2}{12} & 0 \\ 0 & 0 & \frac{m\ell^2}{12} +\frac{mw^2}{12}\end{pmatrix}. \end{aligned}\]

The diagonal elements of \(J\) are called the moments of inertia and represent the rotational mass about the \(\hat{i}\), \(\hat{j}\), and \(\hat{k}\) axes respectively.

We have shown that if the motion of the system of \(N\) masses is caused purely by spinning about the angular velocity vector \(\boldsymbol{\omega}\), then the kinetic energy is given by \(K=\frac{1}{2}\boldsymbol{\omega}^\top J \boldsymbol{\omega}\). If the system of masses is also translating with velocity \(\mathbf{v}\), then the velocity of the \(i^{th}\) particle is given by \(\mathbf{v}_i = \mathbf{v} + \boldsymbol{\omega}\times\mathbf{p}_i.\) In that case, the kinetic energy is given by \[\begin{aligned} K &= \frac{1}{2}\sum_{i=1}^N m_i \mathbf{v}_i^\top\mathbf{v}_i = \frac{1}{2}\sum_{i=1}^N m_i (\mathbf{v} + \boldsymbol{\omega}\times\mathbf{p}_i)^\top(\mathbf{v} + \boldsymbol{\omega}\times\mathbf{p}_i) \\ &= \frac{1}{2}\sum_{i=1}^N m_i \Big(\mathbf{v}^\top\mathbf{v} + \mathbf{v}^\top(\boldsymbol{\omega}\times\mathbf{p}_i)+(\boldsymbol{\omega}\times\mathbf{p}_i)^\top\mathbf{v} % + (\boldsymbol{\omega}\times\mathbf{p}_i)^\top(\boldsymbol{\omega}\times\mathbf{p}_i) \Big) \\ &= \frac{1}{2}\left(\sum_{i=1}^N m_i\right)\mathbf{v}^\top\mathbf{v} + \frac{1}{2} \sum_{i=1}^N m_i (\boldsymbol{\omega}\times\mathbf{p}_i)^\top(\boldsymbol{\omega}\times\mathbf{p}_i) \\ &\qquad\qquad + \frac{1}{2}\sum_{i=1}^N m_i\left[\mathbf{v}^\top(\boldsymbol{\omega}\times\mathbf{p}_i) + (\boldsymbol{\omega}\times\mathbf{p}_i)^\top\mathbf{v}\right]. \end{aligned}\] Using the cross product property \[\mathbf{a}^\top(\mathbf{b}\times\mathbf{c}) = \mathbf{b}^\top(\mathbf{c}\times\mathbf{a}) = \mathbf{c}^\top(\mathbf{a}\times\mathbf{b})\] we get that \[\begin{gathered} \sum_{i=1}^N m_i\left[\mathbf{v}^\top(\boldsymbol{\omega}\times\mathbf{p}_i) + (\boldsymbol{\omega}\times\mathbf{p}_i)^\top\mathbf{v}\right] = \\ (\sum_{i=1}^N m_i\mathbf{p}_i)^\top(\mathbf{v}\times\boldsymbol{\omega}) + (\sum_{i=1}^N m_i\mathbf{p}_i)^\top(\mathbf{v}\times\boldsymbol{\omega}). \end{gathered}\] If we choose the center of the coordinate system to be the center of mass, then \(\sum_{i=1}^N m_i\mathbf{p}_i=0\). Defining the total mass as \(m=\sum_{i=1}^N m_i\), the kinetic energy is given by \[\begin{aligned} K&= \frac{1}{2}\left(\sum_{i=1}^N m_i\right)\mathbf{v}_{cm}^\top\mathbf{v}_{cm} + \frac{1}{2} \sum_{i=1}^N m_i (\boldsymbol{\omega}\times\mathbf{p}_i)^\top(\boldsymbol{\omega}\times\mathbf{p}_i) \notag \\ &= \frac{1}{2}m\mathbf{v}_{cm}^\top\mathbf{v}_{cm} + \frac{1}{2} \boldsymbol{\omega}^\top J_{cm} \boldsymbol{\omega}. \end{aligned}\]

Therefore, the kinetic energy of a rigid body moving with velocity \(\mathbf{v}_{cm}\) and spinning at angular velocity \(\boldsymbol{\omega}\) is \[\tag{2.3} K=\frac{1}{2}m\mathbf{v}_{cm}^\top\mathbf{v}_{cm} + \frac{1}{2} \boldsymbol{\omega}^\top J_{cm} \boldsymbol{\omega}.\] If the system consists of \(n\) rigid bodies, then the total kinetic energy is computed by summing the kinetic energy for each body as \[K = \sum_{j=1}^n \left[\frac{1}{2}m_j\mathbf{v}_{cm,j}^\top\mathbf{v}_{cm,j} + \frac{1}{2} \boldsymbol{\omega}_j^\top J_{cm,j} \boldsymbol{\omega}_j\right].\]

2.1.1 Example: Mass Spring System

As an example, consider the mass spring system shown in Fig. 2.6 where two masses are moving along a surface and are connected by springs with constants \(k_1\) and \(k_2\) and dampers with constants \(b_1\) and \(b_2\). The kinetic energy of the system is determined by the velocity of the two masses and is therefore given by \[K = \frac{1}{2} m_1 \dot{z}_1^2 + \frac{1}{2} m_2 \dot{z}_2^2.\]

2.1.2 Example: Spinning Dumbbell

As another simple example, consider the spinning dumbbell system shown in Fig. 2.7, where two masses are spinning about a point on the adjoining line between the masses. The angular velocity vector is pointing out of the page with magnitude \(\dot{\theta}\) and the dumbbell system is moving to the right at velocity \(\dot{z}\).

The first step in computing the kinetic energy is to define an inertial coordinate system as shown in Fig. 2.7 where \(\hat{i}\) is a unit vector along the direction of travel, \(\hat{j}\) is perpendicular to the direction of travel, and where \(\hat{k}\) forms a right-handed coordinate system and points out of the page. In this coordinate system, the position of the point of rotation is given by \(\mathbf{p}(t)=(z(t), 0, 0)^\top\). Accordingly, the positions of \(m_1\) and \(m_2\) are given by \[\begin{aligned} \mathbf{p}_1(t) &= (z(t)+\ell_1\cos\theta(t), \ell_1\sin\theta(t), 0)^\top, \\ \mathbf{p}_2(t) &= (z(t)-\ell_2\cos\theta(t), -\ell_2\sin\theta(t), 0)^\top. \end{aligned}\] Since \(z\) and \(\theta\) are functions of time, using the chain rule \[\frac{d}{dt}f(g(t)) = \frac{\partial f}{\partial g}\frac{d g}{dt}\] and differentiating with respect to time, we get the velocities of the masses as \[\begin{aligned} \mathbf{v}_1 &= (\dot{z} - \ell_1\dot{\theta}\sin\theta, \ell_1\dot{\theta}\cos\theta, 0)^\top,\\ \mathbf{v}_2 &= (\dot{z} + \ell_2\dot{\theta}\sin\theta, -\ell_2\dot{\theta}\cos\theta, 0)^\top. \end{aligned}\] The kinetic energy is therefore given by \[\begin{aligned} K &= \frac{1}{2}m_1\mathbf{v}_1^\top \mathbf{v}_1 + \frac{1}{2}m_2\mathbf{v}_2^\top \mathbf{v}_2 \\ &= \frac{1}{2}m_1 \left[(\dot{z} - \ell_1\dot{\theta}\sin\theta)^2 + (\ell_1\dot{\theta}\cos\theta)^2 \right] \\&\qquad + \frac{1}{2}m_2\left[(\dot{z} + \ell_2\dot{\theta}\sin\theta)^2 + (-\ell_2\dot{\theta}\cos\theta)^2 \right] \\ &= \frac{1}{2}m_1 \left[\dot{z}^2 - 2\ell_1\dot{z}\dot{\theta}\sin\theta + \ell_1^2\dot{\theta}^2\sin^2\theta+\ell_1^2\dot{\theta}^2\cos^2\theta \right] \\ &\quad + \frac{1}{2}m_2\left[\dot{z}^2 + 2\ell_2\dot{z}\dot{\theta}\sin\theta + \ell_2^2\dot{\theta}^2\sin^2\theta+\ell_2^2\dot{\theta}^2\cos^2\theta \right] \\ &= \frac{1}{2}(m_1+m_2)\dot{z}^2 + \frac{1}{2} (m_1\ell_1^2+m_2\ell_2^2)\dot{\theta}^2 + (m_2\ell_2-m_1\ell_1)\dot{z}\dot{\theta}\sin\theta. \end{aligned}\] Note that if the pivot point is located at the center of mass, then \(m_1\ell_1=m_2\ell_2\), and the kinetic energy is simply given by \[K = \frac{1}{2}(m_1+m_2)\dot{z}^2 + \frac{1}{2} (m_1\ell_1^2+m_2\ell_2^2)\dot{\theta}^2,\] which is in the form \[K = \frac{1}{2}m\mathbf{v}_{cm}^\top\mathbf{v}_{cm} + \frac{1}{2}\boldsymbol{\omega}^\top J \boldsymbol{\omega},\] where the total mass is \(m=m_1+m_2\), the velocity of the center of mass is \(\mathbf{v}_{cm}=(\dot{z}, 0, 0)^\top\), the angular velocity vector is \(\boldsymbol{\omega} = (0, 0, \dot{\theta})^\top\) and the \((3,3)\) element of \(J\) is given by \(J_z = m_1\ell_1^2 + m_2\ell_2^2\).

2.2 Design Study A: Single Link Robot Arm

A definition of the single-link robot arm is given in Design Study A.

2.2.1 Example Problem 3.1.2

(a) Using the configuration variable \(\theta\), write an expression for the kinetic energy of the system.

(b) Referring to Appendices P.1, P.2, and P.3, write a Python or MATLAB class, or a MATLAB function that creates an animation of the single-link robot arm. Simulate the animation to display a sinusoidal variation on the configuration variable \(q=\theta\).

Design Study A guided exercise

Derive the arm kinetic energy

Work from the center-of-mass motion and the rotational inertia of the slender link.

Use this geometry to identify the center-of-mass position and the angular velocity before deriving kinetic energy.

Step 1. Center-of-mass position

Enter the horizontal and vertical components of \(p_{cm}\).

\(p_x\) \(p_y\)

Check Step 1

Step 2. Link inertia

Enter the moment of inertia about the center of mass for the slender link.

\(J_z\)

Check Step 2

Step 3. Total kinetic energy

Enter the multiplier on \(\dot{\theta}^2\) in the final kinetic energy expression.

Multiplier

Check Step 3

Pseudocode

The complete solution includes an animation class and a simulation script. Before opening it, outline each component separately.

These prompts are ungraded. They are intended to connect the derivation to the implementation code in the complete solution.

Show complete solution

2.3 Design Study B. Pendulum on Cart

A definition of the pendulum on a cart (or inverted pendulum) is given in Design Study 3.2.

2.3.1 Example Problem 3.2.2

(a) Using the configuration variables \(z\) and \(\theta\), write an expression for the kinetic energy of the system.

(b) Referring to Appendices P.1, P.2, and P.3, write a Python or MATLAB class, or a MATLAB function that creates an animation of the inverted pendulum. Simulate the animation to display sinusoidal variations on the configuration variables \(q=(z, \theta)^\top\).

Design Study B guided exercise

Build the pendulum kinetic energy

Use the geometry to identify positions, velocities, and energy terms for the cart and rod.

Use this geometry to identify the cart position, rod center-of-mass position, and pendulum angular velocity.

Step 1. Position variables

Rod center of mass includes Cart position

Check Step 1

Step 2. Velocity variables

Rod velocity depends on Cart velocity

Check Step 2

Step 3. Energy terms

Energy types included

Check Step 3

Pseudocode

The complete solution includes an animation class and a simulation script. Before opening it, outline each component separately.

These prompts are ungraded. They are intended to connect the derivation to the implementation code in the complete solution.

Show complete solution

2.4 Design Study C. Satellite Attitude Control

A definition of the satellite attitude control problem is given in Design Study C.

2.4.1 Example Problem 3.3.2

(a) Using the configuration variables \(\theta\) and \(\phi\), write an expression for the kinetic energy of the system.

(b) Referring to Appendices P.1, P.2, and P.3, write a Python or MATLAB class, or a MATLAB function that creates an animation of the satellite. Simulate the animation to display sinusoidal variations on the configuration variables \(q=(\theta, \phi)^\top\).

Design Study C guided exercise

Build the satellite kinetic energy

The satellite model has two rotational inertias: the body inertia \(J_s\) associated with \(\theta\), and the panel inertia \(J_p\) associated with \(\phi\).

Use inertial angular rates for kinetic energy. The relative deflection \(\phi-\theta\) appears in spring and damping terms later.

Step 1. Match inertias to coordinates

Body inertia coordinate Panel inertia coordinate

Check Step 1

Step 2. Identify the kinetic-energy terms

Number of rotational terms Relative rate used directly?

Check Step 2

Step 3. Enter the total kinetic-energy structure

Body term coefficient Panel term coefficient

Check Step 3

Pseudocode

The complete solution includes satellite animation code. Before opening it, outline how an animation should update the body and panel orientations from \(\theta(t)\) and \(\phi(t)\).

This prompt is ungraded. It connects the energy model to the animation code in the complete solution.

Show complete solution

2.5 Interactive Inertia Check

2.6 Knowledge Check

2.7 Notes and References

The inertia matrix in Equation 2.3 is computed about the center of mass of the object. When an object is composed of several connected rigid bodies and the inertia matrix is known for each body, the parallel axis theorem provides a convenient way to find the total inertia about the center of mass. The theorem states that if \(J_{cm}\) is the inertia matrix about the center of mass, then the inertia matrix about a point located at vector \(\mathbf{r}\) from the center of mass is \[J(\mathbf{r}) = J_{cm} + m\left(\mathbf{r}^\top\mathbf{r} I_3 - \mathbf{r}\mathbf{r}^\top\right).\] For the moment of inertia about a single axis, where the axis is shifted from the center of mass by distance \(d\), the parallel axis theorem implies that \[I(d) = I_{cm} + md^2.\]