This web page contains solved examples from statics, dynamics, linear and nonlinear vibrations, and acoustics. All examples are solved using the Python programming language.
Solved Examples - Free Vibration of Single Degree of Freedom Systems
Example 1 The system shown in Figure 1, consist of a body with mass of 1 [kg], and spring with stfness of 2500 [N/m]. Calcuate the magnitude, and the phase of the response of the system when the initial displacement is -5 mm and initial velocity of 200 [mm/s]. Using Python plot the system response in time interval of 20 [s].
Figure 1 Spring-mass system in vertical position.
Solution The amplitude and the phase of the response can be calculated using the expressions:
\begin{eqnarray}
A &=& \sqrt{x_0^2 + \left(\frac{\dot{x}_0}{\omega_n}\right)^2},\nonumber\\
\phi &=& \tan^{-1}\left(\frac{\dot{x}_0}{x_0\omega_n}\right) \nonumber
\end{eqnarray}
To calculate the amplitude and the phase of the respnse it is necessary to calculate the natural frequency. The natural frequency is equal to:
\begin{eqnarray}
\omega_n &=& \sqrt{\frac{k}{m}} = \sqrt{\frac{2500}{1}} = 50 \left[\mathrm{rad/s}\right]\nonumber\\
\end{eqnarray}
The amplitude and the pase are equal to:
\begin{eqnarray}
A &=& \sqrt{x_0^2 + \left(\frac{\dot{x}_0}{\omega_n}\right)^2} = \sqrt{(-5)^2 + \left(\frac{200}{50}\right)^2}\nonumber\\\
A &=& \sqrt{25+16} = 6.40312 \left[\mathrm{mm}\right]\nonumber\\ \
\phi &=& \tan^{-1}\left(\frac{\dot{x}_0}{x_0\omega_n}\right) = \tan^{-1}\left(\frac{200}{-5\cdot 50}\right) = -38.6598°\nonumber\\
\phi &=& 141.34°\nonumber\\
\phi &=& 2.46684 \left[\mathrm{rad}\right]\nonumber
\end{eqnarray}
The next step is to investigate if \(\frac{\dot{x}_0}{\omega_n}\) is positive or negative and if \(x_0\) is also postive or negative. This must be done to see in which quadrant \(\phi\) belogs to i.e. determine the value of \(\phi\). The value of ratio \(\frac{\dot{x}_0}{\omega_n}\) is equal to:
\begin{eqnarray}
\frac{\dot{x}_0}{\omega_n} &=& \frac{200}{50} = 4.\nonumber
\end{eqnarray}
The value of ratio \(\frac{\dot{x}_0}{\omega_n}\) is greater than 0 i.e. positive, and the \(x_0\) is negative which means that \(\phi\) belongs to second quadrant. So the \(\phi\) is equal to 141.34° or \(2.46684\) [rad]. The response of the system can be written in general form as:
\begin{eqnarray}
x(t) &=& A \cos (\omega_n t -\phi).\nonumber
\end{eqnarray}
By inserting the amplitude, natural freuqnecy and phase angle value the response of the system can be written as:
\begin{eqnarray}
x(t) &=& 6.40312\cos(50t - 2.46684) \left[\mathrm{mm}\right].\nonumber
\end{eqnarray}
To plot the system response using Python programming langulage two libraries are required numpy and matplotlib. To import those two libraries simply type in following code.
import numpy as np import matplotlib.pyplot as plt
Since system response will be shown for first 20 [s] using numpy the list of values will be generated that represent afformentioned time interval. However, to generate smooth curve of system response the step between list elements in range from 0 to 20.01 will be 0.01. The code for generating time interval si given below.
t = np.arange(0,20.01, 0.01)
The next step is to calculate the value of system resposne for every value in the time list. This will be done using for loop and the value of system response will be saved in new lest that is named SR.
SR = [] for i in range(len(t)): SR.append(6.40312*np.cos(50*t[i] - 2.46684))
Now that the system response is calculated for all the values of \(t\) the next step is to plot system response using matplotlib library. However, before plotting the values the font stlye will be set to "Times New Roman", and the font size for axes, ticks, legend and title will be set. The code for fonty style and font size is given below.
plt.rcParams["font.family"] = "Times New Roman" SMALL_SIZE = 18 MEDIUM_SIZE=20 BIGGER_SIZE=21 plt.rc('font', size=SMALL_SIZE) # DEFAULT TEXT SIZE plt.rc('axes', titlesize=SMALL_SIZE) # AXES TITLES FONTSIZE plt.rc('axes', labelsize=MEDIUM_SIZE) # x AND y LABELS FONTSIZES plt.rc('xtick', labelsize=SMALL_SIZE) # x-TICK FONTSIZE plt.rc('ytick', labelsize=SMALL_SIZE) # y-TICK FONTSIZE plt.rc('legend', fontsize=SMALL_SIZE) # LEGEND FONTSIZE plt.rc('figure', titlesize=BIGGER_SIZE) # FIGURE TITLE FONTSIZE
The size of the figure is set to \(12\times 8\left[\mathrm{in}\right]\). So now everything is ready for plotting the system response.
In previous lines of code the plot function is used to plot the system response, the grid function is used to show grid in graph, the xlabel and ylabel are used to create the labels on \(x\) and \(y\) axes, and finally to show the graph the show function is used. The graph is shown in Figure 2.
Figure 2 The system response in interval of 20 [s].
The entire code necessary to show the system response in time interval of 20 [s] is given below.
import numpy as np import matplotlib.pyplot as plt t = np.arange(0,20.01,0.01) SR = [] for i in range(len(t)): SR.append(6.40312*np.cos(50*t[i] - 2.46684)) plt.rcParams["font.family"] = "Times New Roman" SMALL_SIZE = 18 MEDIUM_SIZE = 20 BIGGER_SIZE = 21 plt.rc('font', size=SMALL_SIZE) # DEFAULT TEXT SIZE plt.rc('axes', titlesize=SMALL_SIZE) # AXES TITLES FONTSIZE plt.rc('axes', labelsize=MEDIUM_SIZE) # x AND y LABELS FONTSIZES plt.rc('xtick', labelsize=SMALL_SIZE) # x-TICK FONTSIZE plt.rc('ytick', labelsize=SMALL_SIZE) # y-TICK FONTSIZE plt.rc('legend', fontsize=SMALL_SIZE) # LEGEND FONTSIZE plt.rc('figure', titlesize=BIGGER_SIZE) # FIGURE TITLE FONTSIZE plt.figure(figsize=(12,8)) plt.plot(t,SR) plt.grid(True) plt.xlabel("t [s]") plt.ylabel("x(t) [mm]") plt.show()
Example 2 A spring-mass system has a natural period of 0.21 [s]. Calculate the new period of the system if the spring constant is: I) increased by 50%, and II) decreased by 50%. Solution: The formula for the natural period can be written as:
\begin{eqnarray}
\tau &=& 2\pi \sqrt{\frac{m}{k}}.
\end{eqnarray}
From the original natural period the variable \(\sqrt{m}\) will be expressed:
\begin{eqnarray}
\tau_o &=& 2\pi \sqrt{\frac{m}{k}}\nonumber \\
0.21 &=& 2\pi\sqrt{\frac{m}{k}}\Big/\cdot \frac{\sqrt{k}}{2\pi}\nonumber \\
\sqrt{m} &=& \frac{0.21\sqrt{k}}{2\pi}\nonumber
\end{eqnarray}
In I) case the spring constant is increased by 50% i.e.
\begin{eqnarray}
k_I &=& 1.5 k\nonumber
\end{eqnarray}
\begin{eqnarray}
\tau_I &=& \frac{2\pi\sqrt{m}}{\sqrt{k_I}}\nonumber\\
\tau_I &=& \frac{2\pi\sqrt{m}}{1.5k}\nonumber\\
\tau_I &=& \frac{2\pi\frac{0.21\sqrt{k}}{2\pi}}{\sqrt{1.5k}}\nonumber\\
\tau_I &=& 0.1715 \left[\mathrm{s}\right]\nonumber
\end{eqnarray}
In II) case the spring constant is decreased by 50% i.e.
\begin{eqnarray}
k_{II} &=& 0.5 k\nonumber
\end{eqnarray}
\begin{eqnarray}
\tau_{II} &=& \frac{2\pi\sqrt{m}}{\sqrt{k_{II}}}\nonumber\\
\tau_{II} &=& \frac{2\pi\sqrt{m}}{0.5k}\nonumber\\
\tau_{II} &=& \frac{2\pi\frac{0.21\sqrt{k}}{2\pi}}{\sqrt{0.5k}}\nonumber\\
\tau_{II} &=& 0.2970 \left[\mathrm{s}\right]\nonumber
\end{eqnarray}
Example 3 A simply supported beam with square cross section with sides \(5\left[\mathrm{mm}\right]\times 5\left[\mathrm{mm}\right] \times 1\left[\mathrm{m}\right]\). The body of the mass of 2.3 [kg] is positioned at the middle on the beam. The system shown in Figure 3 has the natural freuqency of 30 [rad/s]. Determine the Young's modulus of elasticty of the beam.
Figure 3 - Simply supported beam with body at the cenetr of the beam.
Solution - the weight of the beam is neglected. The natural freuqnecy of transverse vibraiton of the beam can be written as:
\begin{eqnarray}
\omega_n &=& \sqrt{\frac{k}{m}}\nonumber\\
k &=& \frac{192EI}{l^3}\nonumber\\
\end{eqnarray}
In previous equations the \(l\) is the length, \(I\) is the area moment of inertia and is calculated as:
\begin{eqnarray}
I &=& \frac{1}{12}(5\times 10^{-3})^3(5\times 10^{-3})^3 = 0.5208 \times 10^{-10} \left[\mathrm{m}^4\right]\nonumber
\end{eqnarray}
The equation for natural frequency can be used to express the stifness k:
\begin{eqnarray}
\omega_n &=& \sqrt{\frac{k}{m}}\Big/^2\nonumber\\
\omega_n^2 &=& \frac{k}{m}\Big/\cdot m\nonumber\\
k &=& \omega_n^2\cdot m\nonumber
\end{eqnarray}
By equating the previous equation with the equation \(k=\frac{192EI}{l^3}\) the Youngs modulus can be calculated.
\begin{eqnarray}
\omega_n^2\cdot m &=& \frac{192EI}{l^3}\Big/\cdot \frac{l^3}{192I}\nonumber\\
E &=& \frac{\omega_n^2\cdot m\cdot l^3}{192I}\nonumber\\
E &=& \frac{30^2\cdot 2.3 \cdot 1^3}{192 \cdot 0.5208 \times 10^{-10}}\nonumber\\
E &=& 207.0132 \times 10^9 \left[\frac{\mathrm{N}}{\mathrm{m}^2}\right]
\end{eqnarray}
Example 4 Determine the natural frequency of the system shown in Figure 4 assuming that the pulleys to be frictonless and have negligible mass.
Figure 4 - Pulley system
Solution To calculate the natural freuqency first the equivalent stiffness of the system must be obtained in order to sove the example as single-degree-of-freedom problem. In the description of the example it was stated that pulleys are frictionless and massless so the rope's tension is constant and equal to weight \(W\) of the mass \(m\). Using static equilibrium the total movement of the mass m can be determined. The upward force acting on a pulley 1 is 2W and the downward force acting on the pully 2 is also 2W. The pulley 1 moves up by a distance \(\frac{2W}{k_1}\), and the pulley moves down by \(\frac{2W}{k_2}\). The total movement of the body with mass \(m\) is then equal to
\begin{eqnarray}
&2&\left(\frac{2W}{k_1} + \frac{2W}{k_2}\right)\nonumber
\end{eqnarray}
The coefficient \(k_{eq}\) will denote the equivalent spring constant then weight of the mass/ equivalent spring constant is equal to net displacemetn of the mass.
\begin{eqnarray}
\frac{W}{k_{eq}} &=& 4W\left(\frac{1}{k_1} + \frac{1}{k_2}\right) = \frac{4W(k_1 + k_2)}{k_1k_2}\nonumber\\
\frac{1}{k_{eq}} &=& \frac{4(k_1+k_2)}{k_1k_2}\nonumber\\
k_{eq} &=& \frac{k_1k_2}{4(k_1+k_2)}\nonumber
\end{eqnarray}
The differential equation that describes the motion of the mass \(m\) can be written as:
\begin{eqnarray}
m\ddot{x} + k_{eq}x &=& 0\Big/:m\nonumber\\
\ddot{x} + \omega_nx &=& 0.\nonumber
\end{eqnarray}
Then the natural frequency is equal to:
\begin{eqnarray}
\omega_n &=& \sqrt{\frac{k_{eq}}{m}} = \sqrt{\frac{k_1k_2}{4m(k_1+k_2)}}.\nonumber
\end{eqnarray}
Usin the relation \(f_n = \frac{\omega_n}{2\pi}\) the natural freuqency can be written as:
\begin{eqnarray}
f_n &=& \frac{\omega_n}{2\pi} = \frac{1}{4\pi}\sqrt{\frac{k_1k_2}{m(k_1+k_2)}}.
\end{eqnarray}
Example 5 Determine the effect of the mass of the spring on the nautral frequency of the spring-mass system shown in Figure 5.
Figure 5 Spring mass system
Solution: The effect of the mass of the spring on the natural freuqnecy of the spring-mass system will be investigated using energy method to determine the natural frequency. The lenght of the spring is denoted with variable \(l\). The \(x\) denotes the displacemnet of the body with mass \(m\), the displacement at distance \(y\) from the support is equal \(y(x/l)\). The similar approch is addopted for springs velocity i.e. \(y(\dot{x}/l)\).
The kinetic energy of the spring element of length \(dy\) is equal to:
\begin{eqnarray}
dT_s &=& \frac{1}{2} \left(\frac{m_s}{l}dy\right)\left(\frac{y\dot{x}}{l}\right)^2\nonumber
\end{eqnarray}
where \(m_s\) is the mass of the spring. The total kinetic energy of the system consist of kinetic energy of the spring and kinetic energy of the body with mass \(m\).
\begin{eqnarray}
T &=& \frac{1}{2}m\dot{x}^2 + \int_{y=0}^l \frac{1}{2}\left(\frac{m_s}{l}dy\right)\left(\frac{y\dot{x}}{l}\right)^2\nonumber\\
T &=& \frac{1}{2}m\dot{x}^2 + \frac{1}{2}\frac{m_s}{3}\dot{x}^2\nonumber
\end{eqnarray}
The total potential energy of the system is equal to:
\begin{eqnarray}
U &=& \frac{1}{2}kx^2.\nonumber
\end{eqnarray}
The assumption is that the system moves in form of harmonic function which can be written as:
\begin{eqnarray}
x(t) &=& X\cos\omega_n t,\nonumber
\end{eqnarray}
where \(X\) is the maximum displacement of the mass and \(\omega_n\) is the natural freuqency. The maximum kinetic and potential energies can be written as:
\begin{eqnarray}
T_{max}&=& \frac{1}{2}\left(m+\frac{m_s}{3}\right)X^2\omega_n^2\nonumber\\
U_{max} &=& \frac{1}{2}kX^2.\nonumber
\end{eqnarray}
By equating \(T_{max}\) and \(U_{max}\) the expression for natural frequency is obtained.
\begin{eqnarray}
\frac{1}{2}\left(m+\frac{m_s}{3}\right)X^2 \omega_n^2 &=& \frac{1}{2}kX^2\nonumber\\
\left(m+\frac{m_s}{3}\right)\omega_n^2 &=& k\nonumber\\
\omega_n&=& \sqrt{\left(\frac{k}{m+\frac{m_s}{3}}\right)}.
\end{eqnarray}
Example 6 The parameters of a SDOF (single degree of freedom) system shown in Figure 6 are \(m=1 \left[\mathrm{kg}\right]\), \(c = 5 \left[\frac{\mathrm{N-s}}{\mathrm{m}}\right]\), and \(k = 16 \left[\frac{\mathrm{N}}{\mathrm{s}}\right]\). Find the response of the system for the initial conditions:
\(x(0) = 0.1 \left[\mathrm{m}\right]\), and \(\dot{x}(0) = 2 \left[\frac{\mathrm{m}}{\mathrm{s}}\right]\)
\(x(0) = -0.1 \left[\mathrm{m}\right]\), and \(\dot{x}(0) = 2 \left[\frac{\mathrm{m}}{\mathrm{s}}\right]\)
Figure 6 - Single degree of freedom system with cisous damper.
Solution:The response of the system can be generally written as:
\begin{eqnarray}
x(t) &=& X_0e^{\zeta \omega_n t} \sin (\omega_d t + \phi_0)\nonumber
\end{eqnarray}
where \(X_0\) and \(\phi_0\) are equal to:
\begin{eqnarray}
X_0 &=& \frac{\sqrt{x_0^2 \omega_n^2 + \dot{x}_0^2 + 2\zeta \omega_n x_0 \dot{x}_0}}{\omega_d}\nonumber\\
\phi_0 &=& \tan^{-1}\left(\frac{x_0 \omega_d}{\dot{x}_0 + \zeta\omega_n x_0}\right)\nonumber
\end{eqnarray}
Before determining the system response for different boundary conditions the parameters such as natural frequency, damping coefficinet etc have to be determined.
\begin{eqnarray}
\omega_n &=& \sqrt{\frac{k}{m}} = \sqrt{\frac{16}{1}} = 4 \left[\frac{\mathrm{rad}}{\mathrm{s}}\right]\nonumber\\
c_c &=& 2\sqrt{km} = 2 \sqrt{16} = 8 \left[\frac{\mathrm{Ns}}{\mathrm{m}}\right]\nonumber\\
\zeta &=& \frac{c}{c_c} = \frac{5}{8} = 0.625\nonumber\\
\omega_d &=& \sqrt{1-\zeta^2}\omega_n =4\cdot\sqrt{1-0.625^2} = 3.1225 \left[\frac{\mathrm{rad}}{\mathrm{s}}\right]\nonumber
\end{eqnarray}
For the first case of boundary conditions the \(X_0\) and \(\phi_0\) are equal to:
\begin{eqnarray}
X_0 &=& \frac{\sqrt{x_0^2 \omega_n^2 + \dot{x}_0^2 + 2\zeta \omega_n x_0 \dot{x}_0}}{\omega_d}\nonumber\\
X_0 &=& \frac{\sqrt{0.1^2\cdot 4^2 + 2^2 + 2\cdot 0.625 \cdot 4 \cdot 0.1 \cdot 2}}{3.1225} \nonumber\\
X_0 &=& \frac{2.2716}{3.1225} = 0.7275 \left[\mathrm{m}\right]\nonumber\\
\phi_0 &=& \tan^{-1}\left(\frac{x_0 \omega_d}{\dot{x}_0 + \zeta\omega_n x_0}\right)\nonumber\\
\phi_0 &=& \tan^{-1}\left(\frac{0.1\cdot 3.1225}{2+0.625\cdot 4 \cdot 0.1}\right)\nonumber\\
\phi_0 &=& \tan^{-1}\left(\frac{0.3122}{2.25}\right)\nonumber\\
\phi_0 &=& 7.8997° (0.1379 \mathrm{rad})\nonumber
\end{eqnarray}
Now that \(X_0\) and \(\phi_0\) are determined the response of the system can be written as:
\begin{eqnarray}
x(t) &=& 0.7275e^{-0.25t}\sin(3.1225t + 0.1379)\left[\mathrm{m}\right]\nonumber
\end{eqnarray}
For the second case of boundary conditions \(X_0\) and \(\phi_0\) are equal to:
\begin{eqnarray}
X_0 &=& \frac{\sqrt{(-0.1)^2\cdot 4^2 + 2^2 + 2\cdot 0.625\cdot 4\cdot (-0.4)\cdot 2}}{3.1225}\nonumber\\
X_0 &=& \frac{1.7776}{3.1225} = 0.5693\left[\mathrm{m}\right]\nonumber\\
\phi_0 &=& \tan^{-1}\left(\frac{-0.1\cdot (3.1225)}{2+0.625\cdot 4 \cdot (-0.1)}\right) \nonumber\\
\phi_0 &=& \tan^{-1}(-0.1784)\nonumber\\
\phi_0 &=& -10.1167°(-0.1766\mathrm{rad})\nonumber\\
\phi_0 &=& 349.8833°(6.1066 \mathrm{rad})\nonumber
\end{eqnarray}
Since \(\sin \phi_0 = \frac{x_0}{X_0} = \frac{-0.1}{0.5693} = -0.1756 < 0\) and \(\cos\phi_0 = \frac{\dot{x}_0 + \zeta\omega_n x_0}{X_0\omega_d} = 0.9867 > 0\) \(\phi_0\) falls in the fourth quadrant with \(\phi_0 = 349.8833° (6.1066 \mathrm{rad})\). The response of the system can be expressed as
\begin{eqnarray}
x(t) &=& 0.5693 e^{-0.25t}\sin(3.1225t + 6.1066)\left[\mathrm{m}\right]\nonumber
\end{eqnarray}
Now both of system responses will be graphiclly shown using Python. First import numpy and matplotlib library. Then define the list which represent time interval of 20 s, the two empty lists that will be filled with system response values for different values of time.
import numpy as np import matplotlib.pyplot as plt t = np.arange(0,20.001, 0.001) SR1 = []; SR2 = []; for i in range(len(t)): SR1.append(0.7275*np.exp(-0.25*t[i])*np.sin(3.1225*t[i]+0.1379)) SR2.append(0.5693*np.exp(-0.25*t[i])*np.sin(3.1225*t[i]+6.1066))
Before typing code for plotting the system responses the font style, and font size have to be defined which is given below.
plt.rcParams["font.family"] = "Times New Roman" SMALL_SIZE = 18 MEDIUM_SIZE=20 BIGGER_SIZE=21 plt.rc('font', size=SMALL_SIZE) # DEFAULT TEXT SIZE plt.rc('axes', titlesize=SMALL_SIZE) # AXES TITLES FONTSIZE plt.rc('axes', labelsize=MEDIUM_SIZE) # x AND y LABELS FONTSIZES plt.rc('xtick', labelsize=SMALL_SIZE) # x-TICK FONTSIZE plt.rc('ytick', labelsize=SMALL_SIZE) # y-TICK FONTSIZE plt.rc('legend', fontsize=SMALL_SIZE) # LEGEND FONTSIZE plt.rc('figure', titlesize=BIGGER_SIZE) # FIGURE TITLE FONTSIZE
The size of the diagram will be set \(12\cdot 8\left[\mathrm{in}\right]\). Both of the system responses will be plot using plt.plot function with different labels that are denoted as \(x_1(t)\), and \(x_2(t)\). The grid will be set with command plt.grid(True), the xlabel is defined "t [s]" while ylabel is named "x [m]". The legend is activated with plt.legned() and finally to show the graph is plt.show(). The code is given below.
Figure 7 - The responses of the system with different initial conditions.
Example 7 The metal block is placed on rough surface and is attached to a spring with initial displacement of 10 [cm] form its equilibrium position. After 10 cycles of oscillation in 3 [s], the final position of the metal block is 2 [cm] from its equilibrium position. Determine the fiction coefficient between the surface and the metal block.
Solution The system makse 10 cycles in 3 [s] so the period is equal to:
\begin{eqnarray}
\tau_n &=& \frac{3}{10} = 0.3\left[\mathrm{s}\right],\nonumber\\
\end{eqnarray}
the frequency of oscillation is equal to:
\begin{eqnarray}
\omega_n &=& \frac{2\pi}{\tau_n} = \frac{2\pi}{0.3}\nonumber\\
\omega_n &=& 20.944 \left[\frac{\mathrm{rad}}{\mathrm{s}}\right]\nonumber
\end{eqnarray}
The reduction of amplitude in 10 cycles is equal to:
\begin{eqnarray}
10\left(\frac{4\mu m g}{k}\right) &=& 0.1 - 0.02 = 0.08 \left[\mathrm{m}\right]\nonumber\\
\frac{4\mu m g}{k} &=& \frac{0.08}{10}\Big/:4\nonumber\\
\frac{\mu m g}{\omega_n^2 m} &=& \frac{0.08}{40}\Big/\cdot \frac{\omega_n^2}{g}\nonumber\\
\mu &=& \frac{0.2 \cdot 20.944^2}{40\cdot9.81}\nonumber\\
\mu &=& 0.0894294\nonumber
\end{eqnarray}