2.10 HW9

  2.10.1 problem description
  2.10.2 problem 1
  2.10.3 problem 2
  2.10.4 Problem 3
  2.10.5 Problem 4
  2.10.6 Key solution for HW 9

2.10.1 problem description

PDF

PDF (letter size)
PDF (legal size)

2.10.2 problem 1

   2.10.2.1 verification using Matlab
   2.10.2.2 Sketch of each mode

After solving problem \(4.3\) in text, Sketch the deformation of the system when it moves in each of the modes.

pict

There are 2 degrees of freedom, \(\theta _{1}\) and \(\theta _{2}\) as shown in this diagram, using anticlock wise rotation as positive

pict

We solved this problem in HW8, using classical Lagrangian method. This problem will now be solved using power balance method.  The static equilibrium position must be chosen so that all generalized coordinates have value zero. Hence, using the above diagram as the static equilibrium, we take \(\theta _{1}=\theta _{2}=0\) in this position.

Now, as in Lagrangian method, we always start by finding kinetic energy \(T\)\[ T=\frac {1}{2}I_{1}\dot {\theta }_{1}^{2}+\frac {1}{2}I_{2}\dot {\theta }_{2}^{2}\]

Where \(I_{1}=\frac {m_{1}L^{2}}{3}\) and \(I_{2}=\frac {m_{2}L^{2}}{12}+m_{2}\left (\frac {L}{4}\right ) ^{2}\) (using parallel axis theorem). Hence \(I_{2}=\frac {m_{2}L^{2}}{12}+m_{2}\frac {L^{2}}{16}= \)\(\frac {7}{48}L^{2}m_{2}\)

Now we compare the above expression to the quadratic form

\[ T=\frac {1}{2}\left (M_{11}\dot {\theta }_{1}^{2}+M_{22}\dot {\theta }_{2}^{2}+2M_{12}\dot {\theta }_{1}\dot {\theta }_{2}\right ) \]

Hence we see that \(M_{11}=I_{1},M_{22}=I_{2}\),\(M_{12}=M_{21}=0\), therefore the mass matrix is\[ \left [ M\right ] =\begin {pmatrix} I_{1} & 0\\ 0 & I_{2}\end {pmatrix} \]

We now find the potential energy due to springs. For this, we need to write down the relative displacement between end points of each spring. Let \(\Delta _{1}\) be the relative displacement in the first spring \(k_{1}\) and let \(\Delta _{2}\) be the relative displacement in the second spring \(k_{2}\). Hence (and assuming springs remain straight, since we are assuming very stiff springs and small angles) then\begin {align*} \Delta _{1} & =\frac {3L}{4}\theta _{1}+\frac {3L}{4}\theta _{2}\\ \Delta _{2} & =L\theta _{1}+\frac {L}{2}\theta _{2} \end {align*}

Then\begin {align*} V_{spring} & =\frac {1}{2}k_{1}\Delta _{1}^{2}+\frac {1}{2}k_{1}\Delta _{2}^{2}\\ & =\frac {1}{2}k_{1}\left (\frac {3L}{4}\theta _{1}+\frac {3L}{4}\theta _{2}\right ) ^{2}+\frac {1}{2}k_{2}\left (L\theta _{1}+\frac {L}{2}\theta _{2}\right ) ^{2}\\ & =\frac {1}{2}k_{1}\left (\frac {9}{16}L^{2}\theta _{1}^{2}+\frac {9}{8}L^{2}\theta _{1}\theta _{2}+\frac {9}{16}L^{2}\theta _{2}^{2}\right ) +\frac {1}{2}k_{2}\left (L^{2}\theta _{1}^{2}+L^{2}\theta _{1}\theta _{2}+\frac {1}{4}L^{2}\theta _{2}^{2}\right ) \\ & =\frac {9}{32}L^{2}\theta _{1}^{2}k_{1}+\frac {1}{2}L^{2}\theta _{1}^{2}k_{2}+\frac {9}{32}L^{2}\theta _{2}^{2}k_{1}+\frac {1}{8} L^{2}\theta _{2}^{2}k_{2}+\frac {9}{16}L^{2}\theta _{1}\theta _{2}k_{1}+\frac {1}{2}L^{2}\theta _{1}\theta _{2}k_{2}\\ & =\left (\frac {9}{32}L^{2}k_{1}+\frac {1}{2}L^{2}k_{2}\right ) \theta _{1}^{2}+\left (\frac {9}{32}L^{2}k_{1}+\frac {1}{8} L^{2}k_{2}\right ) \theta _{2}^{2}+\left (\frac {9}{16}L^{2}k_{1}+\frac {1}{2}L^{2}k_{2}\right ) \theta _{1}\theta _{2} \end {align*}

Now we compare the above to quadratic form\[ V_{spring}=\frac {1}{2}\left (K_{11}\theta _{1}^{2}+K_{22}\theta _{2}^{2}+2K_{12}\theta _{1}\theta _{2}\right ) \]

We see that\begin {align*} K_{11} & =\frac {9}{16}L^{2}k_{1}+L^{2}k_{2}\\ K_{22} & =\frac {9}{16}L^{2}k_{1}+\frac {1}{4} L^{2}k_{2}\\ K_{12} & =\frac {9}{16}L^{2}k_{1}+\frac {L^{2}}{2}k_{2} \end {align*}

Now we need to find \(V_{gravity}.\)Taking the static equilibrium position as the datum, then upward displacement of center of gravity will be positive and downward displacement is negative. This means the left bar will add positive potential energy due to gravity and the right bar will add negative potential energy, hence\[ V_{g}=m_{1}g\frac {L}{2}\sin \theta _{1}-m_{2}g\frac {L}{4}\sin \theta _{2}\]

Now we need to find the components of the gravity potential energy stiffness matrix. Notice that each term is evaluated at static equilibrium\begin {align*} V_{g_{11}} & =\left (\frac {\partial V_{g}^{2}}{\partial \theta _{1}^{2}}\right ) _{\substack {\theta _{1}=0\\\theta _{2}=0}}=\left (-m_{1}g\frac {L}{2}\sin \theta _{1}\right ) _{\theta _{1}=0}=0\\ V_{g_{22}} & =\left (\frac {\partial V_{g}^{2}}{\partial \theta _{2}^{2}}\right ) _{\substack {\theta _{1}=0\\\theta _{2}=0}}=\left (-m_{1}g\frac {L}{2}\sin \theta _{2}\right ) _{\theta _{2}=0}=0\\ V_{g_{12}} & =\left (\frac {\partial V_{g}^{2}}{\partial \theta _{1}\partial \theta _{2}}\right ) _{\substack {\theta _{1}=0\\\theta _{2}=0}}=0 \end {align*}

Hence, no contribution from gravity is added to the stiffness matrix. All contribution comes from the springs potential energy. Therefore, the stiffness matrix is\[ \left [ K\right ] =\begin {pmatrix} \frac {9}{16}L^{2}k_{1}+L^{2}k_{2} & \frac {9}{16}L^{2}k_{1}+\frac {L^{2}}{2}k_{2}\\ \frac {9}{16}L^{2}k_{1}+\frac {L^{2}}{2}k_{2} & \frac {9}{16}L^{2}k_{1}+\frac {1}{4} L^{2}k_{2}\end {pmatrix} \]

Now since there is no damping, then \(P_{disp}=0\). To find \(P_{in}\) we need to find\[ P_{in}=Q_{1}\theta _{1}+Q_{2}\theta _{2}\]

The only external force is \(F\) which generates a torque \(F\frac {L}{2}\theta _{2}\), hence by comparing to the above\begin {align*} P_{in} & =F\frac {L}{2}\theta _{2}\\ Q_{2} & =F\frac {L}{2} \end {align*}

Now we can make the matrix of EOM\begin {align*} MX^{\prime \prime }+kX & =Q\\\begin {pmatrix} I_{1} & 0\\ 0 & I_{2}\end {pmatrix}\begin {pmatrix} \ddot {\theta }_{1}\\ \ddot {\theta }_{2}\end {pmatrix} +\begin {pmatrix} \frac {9}{16}L^{2}k_{1}+L^{2}k_{2} & \frac {9}{16}L^{2}k_{1}+\frac {L^{2}}{2}k_{2}\\ \frac {9}{16}L^{2}k_{1}+\frac {L^{2}}{2}k_{2} & \frac {9}{16}L^{2}k_{1}+\frac {1}{4} L^{2}k_{2}\end {pmatrix}\begin {pmatrix} \theta _{1}\\ \theta _{2}\end {pmatrix} & =\begin {pmatrix} 0\\ F\frac {L}{2}\end {pmatrix} \\ L^{2}\begin {pmatrix} \frac {m_{1}}{3} & 0\\ 0 & \frac {7}{48}m_{2}\end {pmatrix}\begin {pmatrix} \ddot {\theta }_{1}\\ \ddot {\theta }_{2}\end {pmatrix} +L^{2}\begin {pmatrix} \frac {9}{16}k_{1}+k_{2} & \frac {9}{16}k_{1}+\frac {1}{2}k_{2}\\ \frac {9}{16}k_{1}+\frac {1}{2}k_{2} & \frac {9}{16}k_{1}+\frac {1}{4} k_{2}\end {pmatrix}\begin {pmatrix} \theta _{1}\\ \theta _{2}\end {pmatrix} & =\begin {pmatrix} 0\\ F\frac {L}{2}\end {pmatrix} \end {align*}

Now we can solve the problem given.

When \(m_{1}=m,m_{2}=2m,k_{1}=k,k_{2}=\frac {k}{2}\) we obtain\[ mL^{2}\begin {pmatrix} \frac {1}{3} & 0\\ 0 & \frac {7}{24}\end {pmatrix}\begin {pmatrix} \ddot {\theta }_{1}\\ \ddot {\theta }_{2}\end {pmatrix} +kL^{2}\begin {pmatrix} \frac {17}{16} & \frac {13}{16}\\ \frac {13}{16} & \frac {11}{16}\end {pmatrix}\begin {pmatrix} \theta _{1}\\ \theta _{2}\end {pmatrix} =\begin {pmatrix} 0\\ F\frac {L}{2}\end {pmatrix} \]

To find modes of free vibration, let the RHS \(\left \{ 0\right \} \) then we write\[ \left [ \left [ K\right ] -\omega ^{2}\frac {m}{k}\left [ M\right ] \right ] \left \{ \Phi \right \} =\left \{ 0\right \} \]

Let \(\lambda =\omega ^{2}\frac {m}{k}\), hence\[ \left [ \left [ K\right ] -\lambda \left [ M\right ] \right ] \left \{ \Phi \right \} =\left \{ 0\right \} \]

Solving for eigenvalues\begin {align*} \det \left [ \begin {pmatrix} \frac {17}{16} & \frac {13}{16}\\ \frac {13}{16} & \frac {11}{16}\end {pmatrix} -\lambda \begin {pmatrix} \frac {1}{3} & 0\\ 0 & \frac {7}{24}\end {pmatrix} \right ] & =0\\ \det \begin {bmatrix} \frac {17}{16}-\frac {1}{3}\lambda & \frac {13}{16}\\ \frac {13}{16} & \frac {11}{16}-\frac {7}{24}\lambda \end {bmatrix} & =0\\ 9.7222\times 10^{-2}\lambda ^{2}-0.53906 \lambda +7.0313\times 10^{-2} & =0 \end {align*}

Hence, taking the positive square root only we find\begin {align*} \lambda _{1} & =0.13366\\ \lambda _{2} & =5.4110 \end {align*}

When \(\lambda _{1}=0.13366\) \[ \left [ \left [ k\right ] -\lambda _{1}\left [ M\right ] \right ] \left \{ \Phi \right \} _{1}=\left \{ 0\right \} \]\[\begin {bmatrix} \frac {17}{16}-\frac {1}{3}\lambda _{1} & \frac {13}{16}\\ \frac {13}{16} & \frac {11}{16}-\frac {7}{24}\lambda _{1}\end {bmatrix}\begin {Bmatrix} \Phi _{11}\\ \Phi _{21}\end {Bmatrix} =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \]

Let \(\Phi _{11}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} \frac {17}{16}-\frac {1}{3}\lambda _{1} & \frac {13}{16}\\ \frac {13}{16} & \frac {11}{16}-\frac {7}{24}\lambda _{1}\end {bmatrix}\begin {Bmatrix} 1\\ \Phi _{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ \frac {17}{16}-\frac {1}{3}\lambda _{1}+\frac {13}{16}\Phi _{21} & =0\\ \Phi _{21} & =\frac {16}{13}\left (\frac {1}{3}\lambda _{1}-\frac {17}{16}\right ) =\frac {16}{13}\left (\frac {1}{3}\left (0.137\right ) -\frac {17}{16}\right ) =-1.253 \end {align*}

Hence the first mode associated with \(\lambda _{1}=0.13366\) is \[\begin {Bmatrix} 1\\ -1.253 \end {Bmatrix} \]

When \(\lambda _{2}=5.4110\) \[ \left [ \left [ k\right ] -\lambda _{1}\left [ M\right ] \right ] \left \{ \Phi \right \} _{1}=\left \{ 0\right \} \]\[\begin {bmatrix} \frac {17}{16}-\frac {1}{3}\lambda _{2} & \frac {13}{16}\\ \frac {13}{16} & \frac {11}{16}-\frac {7}{24}\lambda _{2}\end {bmatrix}\begin {Bmatrix} \Phi _{12}\\ \Phi _{22}\end {Bmatrix} =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \]

Let \(\Phi _{12}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} \frac {17}{16}-\frac {1}{3}\lambda _{2} & \frac {13}{16}\\ \frac {13}{16} & \frac {11}{16}-\frac {7}{24}\lambda _{2}\end {bmatrix}\begin {Bmatrix} 1\\ \Phi _{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ \frac {17}{16}-\frac {1}{3}\lambda _{2}+\frac {13}{16}\Phi _{22} & =0\\ \Phi _{22} & =\frac {16}{13}\left (\frac {1}{3}\lambda _{2}-\frac {17}{16}\right ) =\frac {16}{13}\left (\frac {1}{3}\left (5.411\right ) -\frac {17}{16}\right ) =0.912 \end {align*}

Hence the second mode associated with \(\lambda _{2}=5.411\) is \[\begin {Bmatrix} 1\\ 0.912 \end {Bmatrix} \]

Summary\[\begin {tabular} [c]{|l|l|}\hline $\omega $ (rad/sec) & mode shape\\\hline $\lambda =\omega ^{2}\frac {m}{k}\Rightarrow \omega =\sqrt {\frac {k}{m}}\sqrt {0.137}=0.366\sqrt {\frac {k}{m}}$ & $\begin {Bmatrix} 1\\ -1.253 \end {Bmatrix}$ \\\hline $\lambda =\omega ^{2}\frac {m}{k}\Rightarrow \omega =\sqrt {\frac {k}{m}}\sqrt {5.411}=2.326\sqrt {\frac {k}{m}}$ & $\begin {Bmatrix} 1\\ 0.912 \end {Bmatrix} $\\\hline \end {tabular} \ \ \ \ \]

2.10.2.1 verification using Matlab
EDU>> M=[1/3 0;0 7/24]; K=[17/16 13/16;13/16 11/16];
EDU>> [phi,omega]=eig(K,M);
EDU>> sqrt(omega)

    0.3656         0
         0    2.3262

EDU>> phi(:,1)/abs(phi(1,1))

    1.0000
   -1.2529

EDU>> phi(:,2)/abs(phi(1,2))

   -1.0000
   -0.9122

2.10.2.2 Sketch of each mode

\(\begin {Bmatrix} 1\\ 0.912 \end {Bmatrix} \) means that \(\theta _{1}\) and \(\theta _{2}\) are in phase, and for each \(1\) unit rotation of \(\theta _{1}\) there will be \(0.912\) units of rotation of \(\theta _{2}\), while \(\begin {Bmatrix} 1\\ -1.253 \end {Bmatrix} \)means that \(\theta _{1}\) and \(\theta _{2}\) are out of phase, and for each \(1\) unit rotation of \(\theta _{1}\) there will be \(1.253\) units of rotation of \(\theta _{2}\) but in the opposite direction. This is a sketch of both modes

pict

2.10.3 problem 2

pict

pict

Using power balance method, we start by finding the kinetic energy.

pict

Since the top bar does not have a point that is fixed in inertial space as the lower bar does, then we take its moment of inertia around its center of mass, and add a translational kinetic energy due to the motion of its center of mass in space. For the lower bar, since it has a point that is fixed in space, then we take the moment of inertia around that point, and we do not need to account for translational kinetic energy for the lower bar. To find the speed of the center of mass of the top bar, we can either use its coordinates system \(x,y\) differentiate these w.r.t time, or we can use the angular motion of the base of the second bar and add it to the speed of the center of mass of the second bar relative to the base. This is what will be done next:

pict

Therefore, the speed components of the center of mass of the top bar is\begin {align*} v_{x} & =\frac {L}{2}\dot {\theta }_{2}\cos \theta _{2}+L\dot {\theta }_{1}\cos \theta _{1}\\ v_{y} & =-\frac {L}{2}\dot {\theta }_{2}\sin \theta _{2}-L\dot {\theta }_{1}\sin \theta _{1} \end {align*}

So the velocity of the center of mass is\[ v_{c.g.}=\sqrt {v_{x}^{2}+v_{y}^{2}}\]

Now that we have the translation velocity of the top bar, and we know its moment of inertia around its c.g. then we have all the terms needed to obtain the kinetic energy.\[ T=\frac {1}{2}I_{1}\dot {\theta }_{1}^{2}+\frac {1}{2}I_{2}\dot {\theta }_{2}^{2}+\frac {1}{2}mv_{c.g.}^{2}\]

Again, the important thing is to note that \(I_{1}\) is taken around the base of lower rod while \(I_{2}\) is taken around the center of mass of the top rod. Hence \begin {align*} T & =\frac {1}{2}\frac {mL^{2}}{3}\dot {\theta }_{1}^{2}+\frac {1}{2}\frac {mL^{2}}{12}\dot {\theta }_{2}^{2}+\frac {1}{2}m\left (v_{x}^{2}+v_{y}^{2}\right ) \\ & =\frac {1}{2}\frac {mL^{2}}{3}\dot {\theta }_{1}^{2}+\frac {1}{2}\frac {mL^{2}}{12}\dot {\theta }_{2}^{2}+\frac {1}{2}m\left (\left (\frac {L}{2}\dot {\theta }_{2}\cos \theta _{2}+L\dot {\theta }_{1}\cos \theta _{1}\right ) ^{2}+\left ( -\frac {L}{2}\dot {\theta }_{2}\sin \theta _{2}-L\dot {\theta }_{1}\sin \theta _{1}\right ) ^{2}\right ) \\ & =\frac {1}{2}\frac {mL^{2}}{3}\dot {\theta }_{1}^{2}+\frac {1}{2}\frac {mL^{2}}{12}\dot {\theta }_{2}^{2}+\frac {1}{2}m\\ & \left (\left (\frac {L^{2}}{4}\dot {\theta }_{2}^{2}\cos ^{2}\theta _{2}+L^{2}\dot {\theta }_{1}^{2}\cos ^{2}\theta _{1}+L^{2}\dot {\theta }_{2}\dot {\theta }_{1}\cos \theta _{2}\cos \theta _{1}\right ) +\left (\frac {L^{2}}{4}\dot {\theta }_{2}^{2}\sin ^{2}\theta _{2}+L^{2}\dot {\theta }_{1}^{2}\sin ^{2}\theta _{1}+L^{2}\dot {\theta }_{2}\dot {\theta }_{1}\sin \theta _{2}\sin \theta _{1}\right ) \right ) \end {align*}

Simplifying the last term, and using \(\cos ^{2}\theta _{1}+\sin ^{2}\sin ^{2}\theta _{1}=1\) we obtain\[ T=\frac {1}{2}\frac {mL^{2}}{3}\dot {\theta }_{1}^{2}+\frac {1}{2}\frac {mL^{2}}{12}\dot {\theta }_{2}^{2}+\frac {1}{2}m\left (\frac {L^{2}}{4}\dot {\theta }_{2}^{2}+L^{2}\dot {\theta }_{1}^{2}+L^{2}\dot {\theta }_{2}\dot {\theta }_{1}\left (\cos \theta _{2}\cos \theta _{1}+\sin \theta _{2}\sin \theta _{1}\right ) \right ) \]

To compare with the quadratic form, we collect all terms as follows\[ T=\dot {\theta }_{1}^{2}\left (\frac {1}{2}\frac {mL^{2}}{3}+\frac {1}{2}L^{2}m\right ) +\dot {\theta }_{2}^{2}\left (\frac {1}{2}\frac {mL^{2}}{12}+\frac {1}{8}mL^{2}\right ) +\dot {\theta }_{2}\dot {\theta }_{1}\left (\frac {1}{2}mL^{2}\left (\cos \theta _{2}\cos \theta _{1}+\sin \theta _{2}\sin \theta _{1}\right ) \right ) \]

Using \(\cos \theta _{2}\cos \theta _{1}+\sin \theta _{2}\sin \theta _{1}=\cos \left ( \theta _{2}-\theta _{1}\right ) \) the above becomes\[ T=\dot {\theta }_{1}^{2}\left (\frac {4}{6}mL^{2}\right ) +\dot {\theta }_{2}^{2}\left (\frac {1}{6}mL^{2}\right ) +\dot {\theta }_{2}\dot {\theta }_{1}\left ( \frac {1}{2}mL^{2}\cos \left (\theta _{2}-\theta _{1}\right ) \right ) \]

We now compare the above to \[ T=\frac {1}{2}M_{11}\dot {\theta }_{1}^{2}+\frac {1}{2}M_{22}\dot {\theta }_{2}^{2}+M_{12}\dot {\theta }_{2}\dot {\theta }_{1}\]

Therefore\begin {align*} M_{11} & =\frac {4}{3}mL^{2}\\ M_{22} & =\frac {1}{3}mL^{2}\\ M_{12} & =M_{21}=\frac {1}{2}mL^{2}\cos \left (\theta _{2}-\theta _{1}\right ) \end {align*}

I am not sure how to get the same answer given for the mass matrix. Even if I assume that \(\theta _{2}-\theta _{1}\) is very small, hence \(M_{12}=\frac {1}{2}mL^{2}\) then the mass matrix is\[ \left [ M\right ] =mL^{2}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix} \]

Now we find the \(V_{spring}\) the potential energy due to springs.\begin {align*} V_{spring} & =\frac {1}{2}\beta mgL\theta _{1}^{2}+\frac {1}{2}\beta mgL\left ( \theta _{2}-\theta _{1}\right ) ^{2}\\ & =\frac {1}{2}\beta mgL\theta _{1}^{2}+\frac {1}{2}\beta mgL\left (\theta _{2}^{2}+\theta _{1}^{2}-2\theta _{2}\theta _{1}\right ) \\ & =\theta _{1}^{2}\left (\beta mgL\right ) +\theta _{2}^{2}\left (\frac {1}{2}\beta mgL\right ) +\theta _{2}\theta _{1}\left (-\beta mgL\right ) \end {align*}

Comparing to quadratic form \[ V_{spring}=\frac {1}{2}K_{11}\theta _{1}^{2}+\frac {1}{2}K_{11}\theta _{1}^{2}+K_{12}\theta _{1}\theta _{2}\]

Then\begin {align*} K_{11} & =2\beta mgL\\ K_{22} & =\beta mgL\\ K_{12} & =K_{21}=-\beta mgL \end {align*}

Hence the stiffness matrix due to springs only is\[ \left [ K\right ] =mgL\begin {pmatrix} 2\beta & -\beta \\ -\beta & \beta \end {pmatrix} \]

We know need to find the gravity contribution to stiffness. We start by finding the \(V_{gravity}\). We take the datum as the horizontal line at the bottom the lower bar.\[ V_{gravity}=mg\frac {L}{2}\cos \theta _{1}+mg\left (L\cos \theta _{1}+\frac {L}{2}\cos \theta _{2}\right ) \]

Hence\[ V_{11}=\frac {\partial ^{2}V_{g}}{\partial ^{2}\theta _{1}}=-mg\frac {L}{2}\cos \theta _{1}-mg\left (L\cos \theta _{1}\right ) \]

evaluate at \(\theta _{1}=0\) gives\begin {align*} V_{11} & =-mg\frac {L}{2}-mgL\\ & =-\frac {3}{2}mgL \end {align*}

and\[ V_{22}=\frac {\partial ^{2}V_{g}}{\partial ^{2}\theta _{2}}=-mg\left (\frac {L}{2}\cos \theta _{2}\right ) \]

evaluate at \(\theta _{2}=0\) gives\[ V_{22}=-mg\frac {L}{2}\]

and \[ V_{12}=\frac {\partial ^{2}V_{g}}{\partial \theta _{1}\partial \theta _{2}}=0 \]

Hence the stiffness matrix due to gravity is\[ \left [ K\right ] =mgL\begin {pmatrix} -\frac {3}{2} & 0\\ 0 & -\frac {1}{2}\end {pmatrix} \]

Combine the above with the stiffness matrix due to springs we obtain\begin {align*} \left [ K\right ] & =mgL\begin {pmatrix} -\frac {3}{2} & 0\\ 0 & -\frac {1}{2}\end {pmatrix} +mgL\begin {pmatrix} 2\beta & -\beta \\ -\beta & \beta \end {pmatrix} \\ & =mgL\begin {pmatrix} 2\beta -\frac {3}{2} & -\beta \\ -\beta & \beta -\frac {1}{2}\end {pmatrix} \end {align*}

There is no \(P_{disp}\) and no \(P_{in}\) hence the equations of motion are\[ mL^{2}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix}\begin {pmatrix} \ddot {\theta }_{1}\\ \ddot {\theta }_{2}\end {pmatrix} +mgL\begin {pmatrix} 2\beta -\frac {3}{2} & -\beta \\ -\beta & \beta -\frac {1}{2}\end {pmatrix}\begin {pmatrix} \theta _{1}\\ \theta _{2}\end {pmatrix} =\begin {pmatrix} 0\\ 0 \end {pmatrix} \]

For \(\beta =4\)\[ mL^{2}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix}\begin {pmatrix} \ddot {\theta }_{1}\\ \ddot {\theta }_{2}\end {pmatrix} +mgL\begin {pmatrix} \frac {13}{2} & -4\\ -4 & \frac {7}{2}\end {pmatrix}\begin {pmatrix} \theta _{1}\\ \theta _{2}\end {pmatrix} =\begin {pmatrix} 0\\ 0 \end {pmatrix} \]

To find modes of free vibration, we write\[ \left [ \left [ k\right ] -\omega ^{2}\left [ M\right ] \right ] \left \{ \Phi \right \} =\left \{ 0\right \} \]

Solving for eigenvalues\begin {align*} \det \left [ mgL\begin {pmatrix} \frac {13}{2} & -4\\ -4 & \frac {7}{2}\end {pmatrix} -\omega ^{2}mL^{2}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix} \right ] & =0\\ \det \left [ \begin {pmatrix} \frac {13}{2} & -4\\ -4 & \frac {7}{2}\end {pmatrix} -\omega ^{2}\frac {L}{g}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix} \right ] & =0 \end {align*}

Let \(\omega ^{2}\frac {L}{g}=\eta \), hence\begin {align*} \det \left [ \begin {pmatrix} \frac {13}{2} & -4\\ -4 & \frac {7}{2}\end {pmatrix} -\eta \begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix} \right ] & =0\\ \det \begin {bmatrix} \frac {13}{2}-\frac {4}{3}\eta & -4-\frac {1}{2}\eta \\ -4-\frac {1}{2}\eta & \frac {7}{2}-\frac {1}{3}\eta \end {bmatrix} & =0\\ \frac {7}{36}\eta ^{2}-\frac {65}{6}\eta +\frac {27}{4} & =0 \end {align*}

Hence \(\eta =55.084,\eta =0.63023\)

When \(\eta =55.084\)\[ \left [ \left [ k\right ] -\eta \left [ M\right ] \right ] \left \{ \Phi \right \} _{1}=\left \{ 0\right \} \]\begin {align*} \begin {bmatrix} \frac {13}{2}-\frac {4}{3}\eta & -4-\frac {1}{2}\eta \\ -4-\frac {1}{2}\eta & \frac {7}{2}-\frac {1}{3}\eta \end {bmatrix}\begin {Bmatrix} \Phi _{11}\\ \Phi _{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} -66.945 & -31.542\\ -31.542 & -14.861 \end {bmatrix} \begin {Bmatrix} \Phi _{11}\\ \Phi _{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

Let \(\Phi _{11}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} -66.945 & -31.542\\ -31.542 & -14.861 \end {bmatrix}\begin {Bmatrix} 1\\ \Phi _{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ -66.945-31.542\Phi _{21} & =0\\ \Phi _{21} & =-\frac {66.945}{31.542}=-2.1224 \end {align*}

Hence the first mode associated with \(\eta =55.084\) is \[\begin {Bmatrix} 1\\ -2.1224 \end {Bmatrix} \]

When \(\eta =0.63023\)\[ \left [ \left [ k\right ] -\omega _{2}^{2}\left [ M\right ] \right ] \left \{ \Phi \right \} _{2}=\left \{ 0\right \} \]\begin {align*} \begin {bmatrix} \frac {13}{2}-\frac {4}{3}\eta & -4-\frac {1}{2}\eta \\ -4-\frac {1}{2}\eta & \frac {7}{2}-\frac {1}{3}\eta \end {bmatrix}\begin {Bmatrix} \Phi _{12}\\ \Phi _{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} 5.6597 & -4.3151\\ -4.3151 & 3.2899 \end {bmatrix} \begin {Bmatrix} \Phi _{12}\\ \Phi _{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

Let \(\Phi _{12}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} 5.6597 & -4.3151\\ -4.3151 & 3.2899 \end {bmatrix}\begin {Bmatrix} 1\\ \Phi _{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ 5.6597-4.3151\Phi _{22} & =0\\ \Phi _{22} & =\frac {-5.6597}{-4.3151}=1.3116 \end {align*}

Hence the second mode associated with \(\eta =0.630\) is \[\begin {Bmatrix} 1\\ 1.3116 \end {Bmatrix} \]

Summary, \(\omega ^{2}\frac {L}{g}=\eta \ \) hence \(\omega =\sqrt {\eta }\sqrt {\frac {g}{L}}\) \[\begin {tabular} [c]{|l|l|}\hline $\omega _{n}$ (rad/sec) & mode shape\\\hline $7.422\sqrt {\frac {g}{L}}$ & $\begin {Bmatrix} 1\\ -2.1224 \end {Bmatrix} $\\\hline $0.794\sqrt {\frac {g}{L}}$ & $\begin {Bmatrix} 1\\ 1.3116 \end {Bmatrix} $\\\hline \end {tabular} \ \ \ \ \]

For \(\beta =2\)\[ mL^{2}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix}\begin {pmatrix} \ddot {\theta }_{1}\\ \ddot {\theta }_{2}\end {pmatrix} +mgL\begin {pmatrix} 4-\frac {3}{2} & -2\\ -2 & 2-\frac {1}{2}\end {pmatrix}\begin {pmatrix} \theta _{1}\\ \theta _{2}\end {pmatrix} =\begin {pmatrix} 0\\ 0 \end {pmatrix} \]

To find modes of free vibration, we write\[ \left [ \left [ k\right ] -\omega ^{2}\left [ M\right ] \right ] \left \{ \Phi \right \} =\left \{ 0\right \} \]

Solving for eigenvalues\begin {align*} \det \left [ mgL\begin {pmatrix} \frac {5}{2} & -2\\ -2 & \frac {3}{2}\end {pmatrix} -\omega ^{2}mL^{2}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix} \right ] & =0\\ \det \left [ \begin {pmatrix} \frac {5}{2} & -2\\ -2 & \frac {3}{2}\end {pmatrix} -\omega ^{2}\frac {L}{g}\begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix} \right ] & =0 \end {align*}

Let \(\omega ^{2}\frac {L}{g}=\eta \), hence\begin {align*} \det \left [ \begin {pmatrix} \frac {5}{2} & -2\\ -2 & \frac {3}{2}\end {pmatrix} -\eta \begin {pmatrix} \frac {4}{3} & \frac {1}{2}\\ \frac {1}{2} & \frac {1}{3}\end {pmatrix} \right ] & =0\\ \det \begin {bmatrix} \frac {5}{2}-\frac {4}{3}\eta & -2-\frac {1}{2}\eta \\ -2-\frac {1}{2}\eta & \frac {3}{2}-\frac {1}{3}\eta \end {bmatrix} & =0\\ \frac {7}{36}\eta ^{2}-\frac {29}{6}\eta -\frac {1}{4} & =0 \end {align*}

Hence \(\eta =24.909,\eta =-5.162\times 10^{-2}\)

Since \(\omega ^{2}\frac {L}{g}=\eta \ \) hence when \(\eta =-5.162\times 10^{-2}\ \) then \(\omega =\sqrt {\eta }\sqrt {\frac {g}{L}}\) which means there will a complex number for \(\omega \) which is not possible as the frequency must be positive. This means such a system is not stable.It is not possible to obtain the shape functions when \(\omega \) is complex.

2.10.4 Problem 3

problem 4.11 in text: the mass and stiffness matrices of a system are \(\left [ M\right ] =\begin {bmatrix} 4 & 1\\ 1 & 3 \end {bmatrix} \)kg, \(\left [ K\right ] =\begin {bmatrix} 300 & 1\\ 1 & 200 \end {bmatrix} \)kN/m, determine the system natural frequencies and normal vibration modes. (hint, normal modes means mass normalized modes).

Answer:

To find modes of free vibration, we write\[ \left [ \left [ k\right ] -\omega ^{2}\left [ M\right ] \right ] \left \{ v\right \} =\left \{ 0\right \} \]

Solving for eigenvalues\begin {align*} \det \left [ \begin {bmatrix} 300\times 10^{3} & 1000\\ 1000 & 200\times 10^{3}\end {bmatrix} -\omega ^{2}\begin {bmatrix} 4 & 1\\ 1 & 3 \end {bmatrix} \right ] & =0\\ \det \begin {bmatrix} 300\times 10^{3}-4\omega ^{2} & 1000-\omega ^{2}\\ 1000-\omega ^{2} & 200\times 10^{3}-3\omega ^{2}\end {bmatrix} & =0\\ 11\omega ^{4}-1698000\omega ^{2}+59999000000 & =0 \end {align*}

Hence the positive roots are \(\omega =234.02,\omega =315.59\)

When \(\omega _{1}=234.02\) rad/sec then \[ \left [ \left [ k\right ] -\omega _{1}^{2}\left [ M\right ] \right ] \left \{ v\right \} _{1}=\left \{ 0\right \} \]\begin {align*} \begin {bmatrix} 300\times 10^{3}-4\omega _{1}^{2} & 1000-\omega _{1}^{2}\\ 1000-\omega _{1}^{2} & 200\times 10^{3}-3\omega _{1}^{2}\end {bmatrix}\begin {Bmatrix} v_{11}\\ v_{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} 300\times 10^{3}-4\left (234.02\right ) ^{2} & 1000-\left (234.02\right ) ^{2}\\ 1000-\left (234.02\right ) ^{2} & 200\times 10^{3}-3\left (234.02\right ) ^{2}\end {bmatrix}\begin {Bmatrix} v_{11}\\ v_{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} 80939. & -53765.\\ -53765. & 35704. \end {bmatrix} \begin {Bmatrix} v_{11}\\ v_{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

Let \(v_{11}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} 80939. & -53765.\\ -53765. & 35704. \end {bmatrix}\begin {Bmatrix} 1\\ v_{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ 80939.-53765v_{21} & =0\\ v_{21} & =-\frac {80939}{53765}=1.5054 \end {align*}

Hence \[ \left \{ v\right \} _{1}=\begin {Bmatrix} 1\\ 1.5054 \end {Bmatrix} \]

When \(\omega _{2}=315.59\) rad/sec then \[ \left [ \left [ k\right ] -\omega _{2}^{2}\left [ M\right ] \right ] \left \{ v\right \} _{2}=\left \{ 0\right \} \]\begin {align*} \begin {bmatrix} 300\times 10^{3}-4\omega _{2}^{2} & 1000-\omega _{2}^{2}\\ 1000-\omega _{2}^{2} & 200\times 10^{3}-3\omega _{2}^{2}\end {bmatrix}\begin {Bmatrix} v_{12}\\ v_{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} 300\times 10^{3}-4\left (315.59\right ) ^{2} & 1000-\left ( 315.59\right ) ^{2}\\ 1000-\left (315.59\right ) ^{2} & 200\times 10^{3}-3\left ( 315.59\right ) ^{2}\end {bmatrix}\begin {Bmatrix} v_{12}\\ v_{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} -98388. & -98597.\\ -98597. & -98791. \end {bmatrix} \begin {Bmatrix} v_{12}\\ v_{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

Let \(v_{12}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} -98388. & -98597.\\ -98597. & -98791. \end {bmatrix}\begin {Bmatrix} 1\\ v_{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ -98388-98597v_{22} & =0\\ v_{22} & =-\frac {98388}{98597}=-0.998 \end {align*}

Hence\[ \left \{ v\right \} _{2}=\begin {Bmatrix} 1\\ -0.998 \end {Bmatrix} \]

To obtain the mass normalized shape functions:\begin {align*} \mu _{1} & =\left \{ v\right \} _{1}^{T}\left [ M\right ] \left \{ v\right \} _{1}\\ & =\begin {Bmatrix} 1\\ 1.5054 \end {Bmatrix} ^{T}\begin {bmatrix} 4 & 1\\ 1 & 3 \end {bmatrix}\begin {Bmatrix} 1\\ 1.5054 \end {Bmatrix} \\ & =\begin {bmatrix} 5.5054 & 5.5162 \end {bmatrix}\begin {Bmatrix} 1\\ 1.5054 \end {Bmatrix} \\ & =13.809 \end {align*}

And\begin {align*} \mu _{2} & =\left \{ v\right \} _{2}^{T}\left [ M\right ] \left \{ v\right \} _{2}\\ & =\begin {Bmatrix} 1\\ -0.999 \end {Bmatrix} ^{T}\begin {bmatrix} 4 & 1\\ 1 & 3 \end {bmatrix}\begin {Bmatrix} 1\\ -0.998 \end {Bmatrix} \\ & =\begin {bmatrix} 3.002 & -1.994 \end {bmatrix}\begin {Bmatrix} 1\\ -0.998 \end {Bmatrix} \\ & =4.992 \end {align*}

Hence\[ \left \{ \Phi \right \} _{1}=\frac {\left \{ v\right \} _{1}}{\sqrt {\mu _{1}}}=\frac {\begin {Bmatrix} 1\\ 1.505 \end {Bmatrix} }{\sqrt {13.809}}=\begin {Bmatrix} 0.269\\ 0.405 \end {Bmatrix} \]

and\[ \left \{ \Phi \right \} _{2}=\frac {\left \{ v\right \} _{2}}{\sqrt {\mu _{2}}}=\frac {\begin {Bmatrix} 1\\ -0.999 \end {Bmatrix} }{\sqrt {4.992}}=\begin {Bmatrix} 0.446\\ -0.447 \end {Bmatrix} \]

Summary\[\begin {tabular} [c]{|l|l|l|}\hline $\omega _{n}$ (rad/sec) & original mode shape & normal mode shapes\\\hline $234.02$ & $\begin {Bmatrix} 1\\ 1.5054 \end {Bmatrix} $ & $\begin {Bmatrix} 0.269\\ 0.405 \end {Bmatrix} $\\\hline $315.59$ & $\begin {Bmatrix} 1\\ 0.999 \end {Bmatrix} $ & $\begin {Bmatrix} 0.44759\\ -0.44665 \end {Bmatrix} $\\\hline \end {tabular} \ \]

Hence\[ \left [ \Phi \right ] =\begin {bmatrix} 0.2691 & 0.448\\ 0.40511 & -0.447 \end {bmatrix} \]

To verify\begin {align*} \left [ \Phi \right ] ^{T}\left [ M\right ] \left [ \Phi \right ] & =\begin {bmatrix} 0.269 & 0.448\\ 0.405 & -0.447 \end {bmatrix} ^{T}\begin {bmatrix} 4 & 1\\ 1 & 3 \end {bmatrix}\begin {bmatrix} 0.269 & 0.446\\ 0.405 & -0.447 \end {bmatrix} \\ & =\begin {bmatrix} 0.269 & 0.405\\ 0.448 & -0.447 \end {bmatrix}\begin {bmatrix} 4 & 1\\ 1 & 3 \end {bmatrix}\begin {bmatrix} 0.269 & 0.448\\ 0.405 & -0.447 \end {bmatrix} \\ & =\begin {bmatrix} 1.0 & 8.840\times 10^{-5}\\ 8.840\times 10^{-5} & 1.0 \end {bmatrix} \end {align*}

Which is approximately \(\begin {bmatrix} 1.0 & 0\\ 0 & 1.0 \end {bmatrix} \) as expected. calculations were not done with high enough accuracy, so that is why the off-diagonal numerical values were not an exact zeros.

To verify with the \(\left [ K\right ] \) matrix\begin {align*} \left [ \Phi \right ] ^{T}\left [ K\right ] \left [ \Phi \right ] & =\begin {bmatrix} 0.269 & 0.448\\ 0.405 & -0.447 \end {bmatrix} ^{T}\begin {bmatrix} 300\times 10^{3} & 1000\\ 1000 & 200\times 10^{3}\end {bmatrix}\begin {bmatrix} 0.269 & 0.448\\ 0.405 & -0.447 \end {bmatrix} \\ & =\begin {bmatrix} 54765. & 6.594\\ 6.594 & 99600.0 \end {bmatrix} \end {align*}

Note \(\omega _{1}^{2}=234.02^{2}=54765.\) and \(\omega _{2}^{2}=315.59^{2}=99597\) and these are the values on the diagonal as expected. The values off the diagonal should be an exact zero, since the \(\left [ K\right ] \) matrix should be decoupled. Due to low precision in the above calculations, the values did not come out to be zero.

Verify using Matlab. Note that Matlab eig() returns the shape function that are mass normalized

EDU>> M=[4 1;1 3];
EDU>> K=[300*10^3 1000;1000 200*10^3];
EDU>> [eig,lam]=eig(K,M)

eig =

   -0.2691   -0.4476
   -0.4051    0.4467

lam =

   1.0e+04 *

    5.4764         0
         0    9.9600

EDU>> eig'*M*eig

    1.0000         0
         0    1.0000

EDU>> eig'*K*eig

   1.0e+04 *

    5.4764         0
         0    9.9600

2.10.5 Problem 4

problem 4.29 in text.

pict

The mapping between the generalized coordinates \(y_{1},y_{2}\) and \(y_{g},\theta \) is given by

pict

\begin {align*} y_{g}\left (A+B\right ) & =y_{1}A+y_{2}B\\ \theta \left (A+B\right ) & =y_{2}-y_{1} \end {align*}

In our care, \(A=B=\frac {L}{2}\), hence the above becomes\begin {align*} y_{g} & =\frac {y_{1}\frac {L}{2}+y_{2}\frac {L}{2}}{L}=\frac {y_{1}+y_{2}}{2}\\ \theta & =\frac {y_{2}-y_{1}}{L} \end {align*}

Hence taking derivative\begin {align*} \dot {y}_{g} & =\frac {\dot {y}_{1}+\dot {y}_{2}}{2}\\ \dot {\theta } & =\frac {\dot {y}_{2}-\dot {y}_{1}}{L} \end {align*}

Using the power balance method, we start by finding the kinetic energy \(T\)\begin {align*} T & =\frac {1}{2}m\dot {y}_{g}^{2}+\frac {1}{2}I_{cg}\dot {\theta }^{2}\\ & =\frac {1}{2}m\left (\frac {\dot {y}_{1}+\dot {y}_{2}}{2}\right ) ^{2}+\frac {1}{2}\left (mr_{G}^{2}\right ) \left (\frac {\dot {y}_{2}-\dot {y}_{1}}{L}\right ) ^{2} \end {align*}

where \(r_{G}\) is the radius of gyration \(0.4L\), hence\begin {align*} T & =\frac {1}{8}m\left (\dot {y}_{1}^{2}+\dot {y}_{2}^{2}+2\dot {y}_{1}\dot {y}_{2}\right ) +\frac {1}{2}\left (m\left (\frac {4}{10}L\right ) ^{2}\right ) \frac {1}{L^{2}}\left (\dot {y}_{2}^{2}+\dot {y}_{1}^{2}-2\dot {y}_{1}\dot {y}_{2}\right ) \\ & =\frac {1}{8}m\left (\dot {y}_{1}^{2}+\dot {y}_{2}^{2}+2\dot {y}_{1}\dot {y}_{2}\right ) +\frac {8}{100}m\left (\dot {y}_{2}^{2}+\dot {y}_{1}^{2}-2\dot {y}_{1}\dot {y}_{2}\right ) \\ & =\dot {y}_{1}^{2}\left (\frac {1}{8}m+\frac {8}{100}m\right ) +\dot {y}_{2}^{2}\left (\frac {1}{8}m+\frac {8}{100}m\right ) +\dot {y}_{1}\dot {y}_{2}\left (\frac {2}{8}m-\frac {16}{100}m\right ) \\ & =\frac {41}{200}m\dot {y}_{1}^{2}+\frac {41}{200}m\dot {y}_{2}^{2}+\frac {9}{100}m\dot {y}_{1}\dot {y}_{2} \end {align*}

Comparing the above to quadratic form \(T=\frac {1}{2}M_{11}\dot {y}_{1}^{2}\) \(+\frac {1}{2}M_{22}\dot {y}_{2}^{2}+M_{12}\dot {y}_{1}\dot {y}_{2}\) then\begin {align*} M_{11} & =\frac {41}{100}m=0.41m\\ M_{22} & =\frac {41}{100}m=0.41m\\ M_{12} & =0.09m \end {align*}

Hence the mass matrix is\[ \left [ M\right ] =m\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix} \]

\begin {align*} V_{spring} & =\frac {1}{2}k_{left}y_{1}^{2}+\frac {1}{2}k_{right}y_{1}^{2}\\ & =\frac {1}{2}\left (\frac {3}{2}k\right ) y_{1}^{2}+\frac {1}{2}ky_{1}^{2} \end {align*}

Comparing to quadratic form \(\frac {1}{2}K_{11}y_{1}^{2}\) \(+\frac {1}{2}K_{22}y_{2}^{2}+K_{12}y_{1}y_{2}\) then\[ \left [ k_{spring}\right ] =k\begin {bmatrix} 1.5 & 0\\ 0 & 1 \end {bmatrix} \]

\[ V_{gravity}=mgy_{g}=mg\frac {y_{1}+y_{2}}{2}\]

Since this will be evaluated at \(y_{1}=y_{2}=0\) then we see right away that there is no contribution to potential energy to the stiffness matrix. Hence the EOM are\[ m\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {bmatrix} \ddot {y}_{1}\relax (t) \\ \ddot {y}_{2}\relax (t) \end {bmatrix} +k\begin {bmatrix} 1.5 & 0\\ 0 & 1 \end {bmatrix}\begin {bmatrix} y_{1}\relax (t) \\ y_{2}\relax (t) \end {bmatrix} =\begin {bmatrix} 0\\ 0 \end {bmatrix} \]

To convert to \(t^{\prime }\) space, given by \(t^{\prime }=\sqrt {\frac {k}{m}}t\) as required, then we see that \(\frac {dy}{dt}=\) \(\frac {dy}{dt^{\prime }}\frac {dt^{\prime }}{dt}=\frac {dy}{dt^{\prime }}\sqrt {\frac {k}{m}}\), and \(\frac {d^{2}y}{dt^{2}}=\) \(\frac {d^{2}y}{dt^{\prime 2}}\frac {k}{m}\)

Hence the ODE becomes\begin {align*} m\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {bmatrix} \ddot {y}_{1}\left (t^{\prime }\right ) \\ \ddot {y}_{2}\left (t^{\prime }\right ) \end {bmatrix} \frac {k}{m}+k\begin {bmatrix} 1.5 & 0\\ 0 & 1 \end {bmatrix}\begin {bmatrix} y_{1}\left (t^{\prime }\right ) \\ y_{2}\left (t^{\prime }\right ) \end {bmatrix} & =\begin {bmatrix} 0\\ 0 \end {bmatrix} \\ k\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {bmatrix} \ddot {y}_{1}\left (t^{\prime }\right ) \\ \ddot {y}_{2}\left (t^{\prime }\right ) \end {bmatrix} +k\begin {bmatrix} 1.5 & 0\\ 0 & 1 \end {bmatrix}\begin {bmatrix} y_{1}\left (t^{\prime }\right ) \\ y_{2}\left (t^{\prime }\right ) \end {bmatrix} & =\begin {bmatrix} 0\\ 0 \end {bmatrix} \end {align*}

Since \(k\neq 0\) we can divide by it, hence\[\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {bmatrix} \ddot {y}_{1}\left (t^{\prime }\right ) \\ \ddot {y}_{2}\left (t^{\prime }\right ) \end {bmatrix} +\begin {bmatrix} 1.5 & 0\\ 0 & 1 \end {bmatrix}\begin {bmatrix} y_{1}\left (t^{\prime }\right ) \\ y_{2}\left (t^{\prime }\right ) \end {bmatrix} =\begin {bmatrix} 0\\ 0 \end {bmatrix} \]

To find modes of free vibration, we write\[ \left [ \left [ K\right ] -\omega ^{2}\left [ M\right ] \right ] \left \{ v\right \} =\left \{ 0\right \} \]

Solving for eigenvalues\begin {align*} \det \left [ \begin {bmatrix} 1.5 & 0\\ 0 & 1 \end {bmatrix} -\omega ^{2}\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix} \right ] & =0\\ \det \begin {bmatrix} 1.5-0.41\omega ^{2} & -\omega ^{2}0.09\\ -\omega ^{2}0.09 & 1-\omega ^{2}0.41 \end {bmatrix} & =0\\ 0.16\omega ^{4}-1.025\omega ^{2}+1.5 & =0 \end {align*}

Hence the positive roots are \(\omega _{1}=2.0357,\omega _{2}=1.5041\)

When \(\omega _{1}=2.0357\) then \[ \left [ \left [ k\right ] -\omega _{1}^{2}\left [ M\right ] \right ] \left \{ v\right \} _{1}=\left \{ 0\right \} \]\begin {align*} \begin {bmatrix} 1.5-0.41\left (2.0357\right ) ^{2} & -\left (2.0357\right ) ^{2}0.09\\ -\left (2.0357\right ) ^{2}0.09 & 1-\left (2.0357\right ) ^{2}0.41 \end {bmatrix}\begin {Bmatrix} v_{11}\\ v_{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} -0.199 & -0.373\\ -0.373 & -0.699 \end {bmatrix} \begin {Bmatrix} v_{11}\\ v_{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

Let \(v_{11}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} -0.199 & -0.373\\ -0.373 & -0.699 \end {bmatrix}\begin {Bmatrix} 1\\ v_{21}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ -0.199-0.373v_{21} & =0\\ v_{21} & =-\frac {0.199}{0.373}=-0.534 \end {align*}

Hence \[ \left \{ v\right \} _{1}=\begin {Bmatrix} 1\\ -0.534 \end {Bmatrix} \]

When \(\omega _{2}=1.5041\) then\[ \left [ \left [ k\right ] -\omega _{2}^{2}\left [ M\right ] \right ] \left \{ v\right \} _{2}=\left \{ 0\right \} \]\begin {align*} \begin {bmatrix} 1.5-0.41\left (1.504\right ) ^{2} & -\left (1. 504\right ) ^{2}0.09\\ -\left (1.504\right ) ^{2}0.09 & 1-\left (1. 504\right ) ^{2}0.41 \end {bmatrix}\begin {Bmatrix} v_{12}\\ v_{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\\begin {bmatrix} 0.572 & -0.204\\ -0.204 & 0.072 \end {bmatrix} \begin {Bmatrix} v_{12}\\ v_{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

Let \(v_{12}\) be the arbitrary value \(1\) hence\begin {align*} \begin {bmatrix} 0.572 & -0.204\\ -0.204 & 0.072 \end {bmatrix}\begin {Bmatrix} 1\\ v_{22}\end {Bmatrix} & =\begin {Bmatrix} 0\\ \times \end {Bmatrix} \\ 0.572-0.204v_{22} & =0\\ v_{22} & =\frac {0.572}{0.204}=2.812 \end {align*}

Hence\[ \left \{ v\right \} _{2}=\begin {Bmatrix} 1\\ 2.812 \end {Bmatrix} \]

To obtain the mass normalized shape functions:\begin {align*} \mu _{1} & =\left \{ v\right \} _{1}^{T}\left [ M\right ] \left \{ v\right \} _{1}\\ & =\begin {Bmatrix} 1\\ -0.534 \end {Bmatrix} ^{T}\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {Bmatrix} 1\\ -0.534 \end {Bmatrix} \\ & =0.43073 \end {align*}

And\begin {align*} \mu _{2} & =\left \{ v\right \} _{2}^{T}\left [ M\right ] \left \{ v\right \} _{2}\\ & =\begin {Bmatrix} 1\\ 2.812 \end {Bmatrix} ^{T}\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {Bmatrix} 1\\ 2.812 \end {Bmatrix} \\ & =4.1569 \end {align*}

Hence\[ \left \{ \Phi \right \} _{1}=\frac {\left \{ v\right \} _{1}}{\sqrt {0.43073}}=\frac {\begin {Bmatrix} 1\\ -0.53374 \end {Bmatrix} }{\sqrt {0.43073}}=\begin {Bmatrix} 1.5237\\ -0.81326 \end {Bmatrix} \]

and\[ \left \{ \Phi \right \} _{2}=\frac {\left \{ v\right \} _{2}}{\sqrt {\mu _{2}}}=\frac {\begin {Bmatrix} 1\\ 2.812 \end {Bmatrix} }{\sqrt {4.1569}}=\begin {Bmatrix} 0.491\\ 1.379 \end {Bmatrix} \]

Summary\[\begin {tabular} [c]{|l|l|l|}\hline $\omega $ (rad/sec) & original mode shape & normal mode shapes\\\hline $2.0357$ & $\begin {Bmatrix} 1\\ -0.534 \end {Bmatrix} $ & $\begin {Bmatrix} 1.524\\ -0.813 \end {Bmatrix} $\\\hline $1.5041$ & $\begin {Bmatrix} 1\\ 2.812 \end {Bmatrix} $ & $\begin {Bmatrix} 0.491\\ 1.379 \end {Bmatrix} $\\\hline \end {tabular} \ \ \ \]

Hence\[ \left [ \Phi \right ] =\begin {bmatrix} 1.524 & 0.491\\ -0.813 & 1.379 \end {bmatrix} \]

To verify\begin {align*} \left [ \Phi \right ] ^{T}\left [ M\right ] \left [ \Phi \right ] & =\begin {bmatrix} 1.524 & 0.491\\ -0.813 & 1.379 \end {bmatrix} ^{T}\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {bmatrix} 1.524 & 0.491\\ -0.813 & 1.379 \end {bmatrix} \\ & =\begin {bmatrix} 1.0 & -1.9688\times 10^{-4}\\ -1.9688\times 10^{-4} & 1.0 \end {bmatrix} \end {align*}

Which is approximately \(\begin {bmatrix} 1.0 & 0\\ 0 & 1.0 \end {bmatrix} \) as expected. calculations were not done with high enough accuracy, so that is why the off-diagonal numerical values were not an exact zeros.

To verify with the \(\left [ K\right ] \) matrix\begin {align*} \left [ \Phi \right ] ^{T}\left [ K\right ] \left [ \Phi \right ] & =\begin {bmatrix} 1.5237 & 0.49047\\ -0.81326 & 1.3790 \end {bmatrix} ^{T}\begin {bmatrix} 1.5 & 0\\ 0 & 1 \end {bmatrix}\begin {bmatrix} 1.5237 & 0.49047\\ -0.81326 & 1.3790 \end {bmatrix} \\ & =\begin {bmatrix} 4.1439 & -4.9183\times 10^{-4}\\ -4.9183\times 10^{-4} & 2.2625 \end {bmatrix} \end {align*}

Verify using Matlab

EDU>> K=[1.5 0;0 1]; M=[0.41  0.09;0.09 0.41];
EDU>> [eig,lam]=eig(K,M)

eig =
   -0.4905   -1.5238
   -1.3789    0.8130

lam =
    2.2624         0
         0    4.1439

Now we can solve the problem. Using \(\left \{ x\right \} =\left [ \Phi \right ] \left \{ \eta \right \} \), where \begin {align*} \left \{ \eta \right \} & =\left [ \Phi \right ] ^{-1}\left \{ x\right \} \\ & =\left [ \Phi \right ] ^{T}\left [ M\right ] \left \{ x\right \} \end {align*}

Hence, initial conditions in the \(\left \{ \eta \right \} \) space is \begin {align*} \left \{ \eta \right \} _{0} & =\left [ \Phi \right ] ^{T}\left [ M\right ] \left \{ x\right \} _{0}\\ & =\left [ \Phi \right ] ^{T}\left [ M\right ] \begin {Bmatrix} y_{1}\relax (0) \\ y_{2}\relax (0) \end {Bmatrix} \\ & =\begin {bmatrix} 1.524 & 0.491\\ -0.813 & 1.379 \end {bmatrix} ^{T}\begin {bmatrix} 0.41 & 0.09\\ 0.09 & 0.41 \end {bmatrix}\begin {Bmatrix} \frac {mg}{k}\\ 0 \end {Bmatrix} \\ & =\begin {Bmatrix} 0.552\frac {g}{k}m\\ 0.325\frac {g}{k}m \end {Bmatrix} \end {align*}

and \begin {align*} \left \{ \eta ^{\prime }\right \} _{0} & =\left [ \Phi \right ] ^{T}\left [ M\right ] \left \{ x^{\prime }\right \} _{0}\\ & =\left [ \Phi \right ] ^{T}\left [ M\right ] \begin {Bmatrix} y_{1}^{\prime }\relax (0) \\ y_{2}^{\prime }\relax (0) \end {Bmatrix} \\ & =\left [ \Phi \right ] ^{T}\left [ M\right ] \begin {Bmatrix} 0\\ 0 \end {Bmatrix} \\ & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

So, we need to solve \[\begin {bmatrix} 1 & 0\\ 0 & 1 \end {bmatrix}\begin {Bmatrix} \ddot {\eta }_{1}\\ \ddot {\eta }_{2}\end {Bmatrix} +\begin {bmatrix} 4.1439 & 0\\ 0 & 2.2625 \end {bmatrix}\begin {Bmatrix} \eta _{1}\\ \eta _{2}\end {Bmatrix} =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \]

with the initial conditions\begin {align*} \begin {Bmatrix} \eta _{1}\relax (0) \\ \eta _{2}\relax (0) \end {Bmatrix} & =\begin {Bmatrix} 0.552\frac {g}{k}m\\ 0.325\frac {g}{k}m \end {Bmatrix} \\\begin {Bmatrix} \eta _{1}^{\prime }\relax (0) \\ \eta _{2}^{\prime }\relax (0) \end {Bmatrix} & =\begin {Bmatrix} 0\\ 0 \end {Bmatrix} \end {align*}

The solution is given by\[ \eta _{1}\relax (t) =A\cos \sqrt {4.144}t+B\sin \sqrt {4.144}t \]

When \(t=0\), \(\eta _{1}\relax (0) =0.55152\frac {g}{k}m=A\). Taking derivative gives \(\eta _{1}^{\prime }\relax (t) =-A\sin \sqrt {4.144}t+B\cos \sqrt {4.144}t\), hence when \(t=0\) we have \(0=B,\) therefore\[ \eta _{1}\relax (t) =0.552\frac {g}{k}m\cos \sqrt {4.144}t \]

Now we solve for \(\eta _{2}\relax (t) ,\)The solution is given by\[ \eta _{2}\relax (t) =A\cos \sqrt {2.263}t+B\sin \sqrt {2.263}t \]

When \(t=0\), \(\eta _{2}\relax (0) =0.3252\frac {g}{k}m=A\). and \(0=B,\) therefore\[ \eta _{2}\relax (t) =0.325\frac {g}{k}m\cos \sqrt {2.263}t \]

Now we obtain the solution in the \(y\) space\begin {align*} \left \{ y\right \} & =\left [ \Phi \right ] \left \{ \eta \right \} \\\begin {Bmatrix} y_{1}\left (t^{\prime }\right ) \\ y_{2}\left (t^{\prime }\right ) \end {Bmatrix} & =\begin {bmatrix} 1.524 & 0.491\\ -0.813 & 1.379 \end {bmatrix}\begin {Bmatrix} 0.552\frac {g}{k}m\cos \sqrt {4.1439}t\\ 0.325\frac {g}{k}m\cos \sqrt {2.2625}t \end {Bmatrix} \\ & =\begin {Bmatrix} 0.840\frac {g}{k}m\cos \left (2.036t^{\prime }\right ) +0.1595\frac {g}{k}m\cos \left (1.504t^{\prime }\right ) \\ 0.448\frac {g}{k}m\cos \left (1.504t^{\prime }\right ) -0.449\frac {g}{k}m\cos \left (2.036t^{\prime }\right ) \end {Bmatrix} \end {align*}

We are supposed to obtain the answer \[\begin {Bmatrix} y_{1}\relax (t) \\ y_{2}\relax (t) \end {Bmatrix} =\begin {Bmatrix} 0.16\cos \left (1.5t^{\prime }\right ) +0.84\cos \left (2t^{\prime }\right ) \\ 0.45\cos \left (1.5t^{\prime }\right ) -0.45\cos \left (2t^{\prime }\right ) \end {Bmatrix} \]

The answers agree. The scalar \(\frac {g}{k}m\) for some reason is not shown in the key solution.

2.10.6 Key solution for HW 9

PDF