Open Educational Resources

Advanced Dynamics and Vibrations: Continuous Systems

Continuous Systems and Waves

The difference between discrete and continuous systems is that we go from a finite number to an infinite number of degrees of freedom. Mathematically, this means that the equations of motion go from ordinary to partial differential equations.

We will consider a number of continuous systems including strings, rods, and beams.

First, consider a flexible string of mass per unit length \rho(x) and stretched under a tension T.

Where y(x,t) is the shape. The string could also be loaded laterally with a f(x,t).

For small angles, \theta \approx \frac{\partial y}{\partial x}. Thus:

    \[ \begin{split} + \uparrow \sum F_y &= f(x,t) dx + \Big(T + \frac{\partial T}{\partial x} dx \Big) \Big( \frac{\partial y}{\partial x} + \frac{\partial^2 y(x)}{\partial x^2} dx \Big) - T \frac{\partial y}{\partial x} \\&= \rho dx\frac{\partial^2 y}{\partial t^2} \end{split} \]

    \[ \frac{\partial T}{\partial x} \frac{\partial y}{\partial x} dx + T \frac{\partial^2 y}{\partial x^2} dx + f(x,t)dx = \rho dx\frac{\partial^2 y}{\partial t^2} \]

Therefore:

    \[ \frac{\partial}{\partial x}\Big[ T(x) \frac{\partial y}{\partial x} \Big] + f(x,t) = \rho(x)\frac{\partial^2 y}{\partial t^2} \]

As we did for discrete systems, we look for the so-called synchronized motions and start with the free-vibration case. Therefore, assume:

    \[ y(x,t) = Y(x)G(t) \]

We are looking for special cases where all the particles are described by Y(x) and this varies in time.

    \[ \begin{split} \frac{1}{\rho (x) Y(x)} &= \frac{d}{dx} \Big[T(x) \frac{dY}{dx}\Big] \\&= \frac{1}{G(t)} \frac{d^2 G(t)}{dt} \end{split} \]

As the LHS is a function of only x and the RHS is a function of only time, each side must be equal to a constant. Therefore, without loss of generality, set:

    \[ -p^2 = \frac{1}{G(t)} \frac{d^2 G(t)}{d t^2} \]

Therefore, the solution for G(t) is:

    \[ G(t) = C \cos (pt - \phi) \]

Where C and \phi are constants. This leaves an O.D.E.

    \[ \frac{d}{dx}\Big[ T(x) \frac{d Y(x)}{dx} \Big] + p^2 \rho(x)Y(x) = 0 \quad (*) \]

instead of an algebraic problem that we had in the discrete case. Solving this, yields the natural frequencies and the mode shapes.

    \[ [p_1,Y_1],\ [p_2,Y_2]\ .... \]

However, there are an infinite number.

We can normalize the eigenfunctions, as they satisfy the orthogonality conditions.

    \[ \int_{0}^{L} \rho(x) Y_r(x) Y_s(x) dx = \delta rs \]

We will also see that

    \[ \int_{0}^{L} T(x) \frac{d Y_r}{d x} \frac{d Y_s}{d x} dx = p_r^2 \delta rs \]

And this leads to an expansion theorem analogous to what we saw in the discrete case. This is the normal mode solution and therefore not the only solution of the partial differential equation. We will consider the other types of solutions.

Consider an example of a uniform string, where \rho(x) = \rho (constant) and T(x) = T. Note that:

    \[ \beta^2 \equiv \frac{p^2 \rho}{T} \]

The equation (*) becomes:

    \[ \frac{d^2 Y}{d x^2} + \beta^2 Y(x) = 0 \]

Thus:

    \[ Y(x) = A \sin \beta x + B \cos \beta x \]

Where A and B are found from the boundary conditions.

    \[y(0,t) = 0,\ y(L,t) = 0 \]

Therefore:

    \[ Y(0) = 0,\ Y(L) = 0\]

    \[ B = 0 \]

Therefore:

    \[ \begin {split} Y(L) &= A \sin \beta L \\&= 0 \text{ (characteristic equation)} \end{split} \]

    \[ \beta L = n\pi \]

    \[ p_n \sqrt{\frac{\rho}{T}}L = n\pi \]

    \[p_n = \frac{n\pi}{L} \sqrt{\frac{T}{\rho}} \]

    \[ \sqrt{\frac{T}{\rho}} \equiv c \]

Thus:

    \[p_n = \frac{n\pi c}{L} \]

We can have other boundary conditions, for example

    \[ T\frac{\partial y}{\partial x}|_{x=0} = 0 \]

For pin supports, the mode shapes are simply:

    \[ A_n \sin \frac{n\pi x}{L} \]

1st Mode

    \[p_1 = \frac{\pi c}{L} \]

2nd Mode

    \[p_2 = \frac{2 \pi c}{L} \]

Check the orthogonality condition.

    \[ \int_{0}^{L} \rho(x) Y_r(x) Y_s(x) dx = \delta rs \]

If \rho(x) = \rho:

    \[ \begin{split} \rho \int_{0}^L A_r^2 \sin^2 \frac{r \pi x}{L} dx &= \frac{A^2_r \rho}{2} \int_{0}^L \Big(1 - \cos \frac{2r\pi x}{L} \Big) dx \\&= \frac{A_r^2 \rho}{2} \Big\{\Big[L - \frac{L}{2\pi r} \sin \frac{2 r \pi x}{L}\Big]_0^L \Big\} \\&= \frac{A_r^2 \rho L}{2} \\&= 1 \end{split} \]

Therefore:

    \[ A_r = \sqrt{\frac{2}{\rho L}},\ Y_r(x) = \sqrt{\frac{2}{\rho L}} \sin \frac{r \pi x}{L} \]

If r \neq s:

    \[ \int_0^L \rho Y_r Y_S dx = \frac{2}{L} \int_0^L \sin \frac{r \pi x}{L} \sin \frac{s \pi x}{L} dx \]

Use the identity \sin \alpha \sin \beta = \frac{1}{2} [ \cos ( \alpha - \beta) - \cos (\alpha + \beta) ]:

    \[ \begin{split} &= \frac{1}{L} \int_0^L \Big( \cos (r-s) \frac{\pi x}{L} - \cos (r+s) \frac{\pi x}{L} \Big) dx \\&= 0 \end{split} \]

If we consider the other boundary condition:

    \[ \frac{ dY } {dx}\Big|_{x=L} = 0,\ Y(x) = A \sin \beta x \]

    \[ \frac{ dY } {dx} = A\beta \cos \beta x \]

Therefore:

    \[ 0 = \cos \beta L \Rightarrow \beta L = \Big(2n + 1\Big) \frac{\pi}{2} \qquad n = 0, 1, 2, \ldots \]

    \[ Y(x) = A_n \sin (2n+1) \frac{\pi}{2} \frac{x}{L} \]

1st Mode

    \[p_1 = \frac{\pi c}{2L}, \ \sin \frac{\pi x}{2L} \]

2nd Mode

    \[p_2 = \frac{3 \pi c}{2L}, \ \sin \frac{3 \pi x}{2L} \]

Consider a more general solution of the equation when T is constant:

    \[ \frac{\partial^2 y}{\partial x^2} = \frac{\rho}{T} \frac{\partial^2 y}{\partial t^2} \]

And since we defined c = \sqrt{\frac{T}{\rho}}:

    \[ \frac{\partial^2 y}{\partial x^2} = \frac{1}{c^2} \frac{\partial^2 y}{\partial t^2} \]

If we take any functions of the form:

    \[ y = F(ct - x) + G(ct + x) \]

Where ct - x := \gamma and ct + x := \epsilon. Then F and G both satisfy the differential equation. Thus:

    \[ \begin{split} \frac{\partial y}{\partial x} &= \frac{\partial F}{\partial \gamma} \frac{\partial \gamma}{\partial x} + \frac{\partial G}{\partial \epsilon}\frac{\partial \epsilon}{\partial x}\\&= -\frac{\partial F}{\partial \gamma} + \frac{\partial G}{ \partial \epsilon} \end{split} \]

    \[ \frac{\partial^2 y}{\partial x^2} = \frac{\partial^2 F}{\partial \gamma^2} + \frac{\partial^2 G}{\partial \epsilon^2} \]

    \[ \frac{\partial y}{\partial t} = c\Big(\frac{\partial F}{\partial \gamma} + \frac{\partial G}{\partial \epsilon}\Big) \]

    \[ \frac{\partial^2 y}{\partial t^2} = c^2 \Big(\frac{\partial^2 F}{\partial \gamma^2} + \frac{\partial^2 G}{\partial \epsilon^2} \Big) \]

Therefore, both F and G satisfy the differential equation. Consider a physical interpretation of F and G. This leads to the idea of waves on the string. In this case, the particle motion is perpendicular to the wave motion (transverse wave) as opposed to longitudinal waves.

Consider any function F that might be:

If the ct_o - x_o = ct_1 - x_1, then the value of F must be the same. Therefore, x_1 - x_o = ct_1 - ct_o. Thus:

    \[ \begin{split} c &= \frac{x_1 - x_o}{t_1 - t_o} \\&= \frac{\Delta x}{\Delta t} \end{split} \]

This is a velocity and c is called the wave speed for a wave moving in the positive direction.

Similarly for G(ct + x), we can show that the wave is moving in the negative direction. This can be used in several applications as:

    \[ c = \sqrt{\frac{T}{\rho}} \]

Therefore, if we measure c, we can determine T!

What happens when this waves hits a boundary? As the disturbance reaches an end point it is reflected and the incident plus reflected wave must add to zero to satisfy the boundary condition.

As a result, there is a 180 degree phase shift and a wave of F(ct -x) becomes -F(ct+x).

So that the sum at x=\ell is zero. As the velocity is c, the time for the round trip is \frac{2\ell}{c}.

One of the functions that is a positive travelling wave F(ct-x) is a harmonic function – \sin(ct-x). Imagine a wave emanating from the origin in this form.

    \[\begin{split}\lambda &= c \cdot \text{period of oscillation} \\ &= c\tau = c\left(\frac{1}{f}\right) \\ \therefore c &= f\lambda \end{split} \]

We can then see the vibration occurs when the positive wave and negative wave at the same time.

This is then a standing wave.

This is what Naum Gabo called the “Kinetic Art” in 1919!

Example

Consider a guy wire 1 000 ft long made of 1 1/4″ diameter bridge strand cable that weighs 2.50 lbs/ft with breaking strength of 130 000 lbs, working at 50 000 lbs.

    \[ \begin{split} c &= \sqrt{\frac{(50000)(32.2)}{2.50}} \\&= 802 \text{ ft/sec} \end{split} \]

Therefore:

    \[ \begin{split} \frac{2 \ell}{c} &= \frac{2000}{802} \\&\approx 2.5 \text{ seconds} \end{split} \]

Note: Analogous to the string is the longitudinal and torsional vibration of a bar. Each differential equation looks the same however the speed c has a different interpretation.

Example: CBC Tower in Calgary

    \[ 1 \frac{1}{8}"\ \phi \text{ Bridge Strand} \]

    \[T @ 10^{\circ} C = 23400 \text{ lbs} \]

    \[\text{Cable Weighs } 2.61 \text{ lbs/ft} \]

    \[ c^2 = \frac{(23400)(32.2)}{2.61} \Rightarrow c = 537 \text{ ft/sec (366 MPH)} \]

Therefore, the time to travel 2 \cdot 1037 ft:

    \[ \begin{split} t &= \frac{1037(2)}{537} \\&= 3.86 \text{ seconds} \end{split} \]

On site measurements of wind-induced cable vibration showed a standing wave pattern

Which mode of vibration is this?

    \[ \begin{split} f_1 &= \frac{c}{2L} \\&= \frac{1}{3.86} \\&= 0.26 \text{ Hz} \end{split} \]

    \[ \begin{split} f_2 &= \frac{2c}{2L} \\&= 0.58 \text{ Hz} \end{split} \]

    \[ f_n = 8.7 \text{ Hz (measured)} \]

Where n = ?

    \[\begin{split} f_n &= \frac{nc}{2L} \\&= 8.7 \end{split} \]

    \[\begin{split} n &= 8.7(3.86) \\n &= 33 \sim 34 \text{ mode!} \end{split} \]

Longitudinal Vibration of a Bar

Assuming a 1-D stress field:

    \[ P = EA\frac{\partial u}{\partial x} \]

    \[ \Big(P + \frac{\partial P}{\partial x} dx \Big) + f(x,t) dx - P = \rho A dx \frac{\partial^2 u}{\partial t^2} \]

Therefore:

    \[ \frac{\partial}{\partial x} \Big(EA \frac{\partial u}{\partial x} \Big) + f(x,t) = \rho(x) A(x) \frac{\partial^2 u}{\partial t^2} \]

For a uniform bar:

    \[EA \frac{\partial^2 u}{\partial x^2} = \rho A \frac{\partial^2 u}{\partial t^2} \]

The equation is similar to the string:

    \[ \frac{ \partial^2 u}{\partial x^2} = \frac{1}{c_A^2} \frac{\partial^2 u}{\partial t^2} \]

For Steel:

    \[ \begin{split} c_A &= \sqrt{\frac{E}{\rho}} \\&= \sqrt{\frac{30 \cdot 10^6 \cdot 386}{286}} \\&= 16860 \text{ ft/s} \end{split} \]

Where the subscript A indicates axial wave velocity.

Again we can have boundary conditions which are similar to the string.

With:

C – damping coefficient

Therefore:

    \[ m_2\ddot{u} \Big|_{x=\ell} = -EA\frac{\partial u}{\partial x} (\ell, t) - ku(\ell, t) - C \frac{\partial u}{\partial t} (\ell t) \]

Is the general boundary condition. The solution follows the same procedure as the string:

    \[ u(x,t) = X(x)T(t) \]

    \[\begin{split} c_A^2 \frac{X''(x)}{X(x)} &= \frac{T''(t)}{T(t)} \\&= -p^2 \end{split} \]

Therefore:

    \[ X''(x) + \frac{p^2}{c_A^2}X(x) = 0 \]

    \[ X(x) = D \sin \frac{p}{c_A} x + B \cos \frac{p}{c_A} x \]

If at x = 0, u(0,t) = 0, therefore:

    \[ B = 0 \]

Thus:

    \[ u(x,t) = D \sin \frac{p}{c_A} x \cos( pt - \phi) \]

    \[ \ddot{u}(x,t) = -p^2 D \sin \frac{p}{c_A} x \cos(pt-\phi)\]

Consider the free boundary condition at x = L:

    \[ u(0,t) = 0,\ P(L,t) = 0 \Rightarrow \frac{\partial u}{\partial t}(L,t) = 0 \]

    \[ X(x) = D \sin \frac{p}{c_A} x + B \cos \frac{p}{c_A}x \]

    \[X(0) = 0 \Rightarrow B = 0 \]

    \[ \frac{\partial u}{\partial x}\Big|_{x=L} = 0 \Rightarrow \frac{dX}{dx} \Big|_{x=L} = 0 \]

    \[ D \frac{p}{c_A} \cos \frac{p}{c_A} L = 0 \]

Therefore:

    \[ p \frac{L}{c_A} = (2n + 1) \frac{\pi}{2} \quad n = 0,1,2,... \]

    \[ p_1 \frac{L}{c_A} = \frac{\pi}{2} \]

    \[ p_1 = \frac{c_A}{L} \frac{\pi}{2} \]

Therefore:

    \[ \begin{split} \text{for } p_1 &\Rightarrow D \sin\left( \frac{\pi c_A x}{2 L c_A}\right) = X_1 \quad \text{First Mode} \\X_1 &= D \sin\left( \frac{\pi x}{2L}\right) \end{split} \]

For the 2nd Mode:

    \[ \frac{p_2 L}{c_A} = \frac{3\pi}{2} \]

    \[ p_2 = \frac{3}{2} \pi \frac{c_A}{L} \]

    \[ X_2 = D \sin\left( \frac{3}{2} \pi \frac{c_A}{L} \frac{x}{c_A}\right) = D \sin\left( \frac{3 \pi x}{2L}\right) \]

We can now consider a mass m_2 at the end of the bar. The boundary condition is then:

    \[ m_2 \frac{d^2 u}{dx^2}\Big|_{x=L} = -EA \frac{\partial u}{\partial x} (L, t) \]

Therefore:

    \[ -EAD \frac{p}{c_A} \cos \frac{p}{c_A} x\Big|_{x=L} = -m_2p^2D \sin \frac{p}{c_A}x\Big|_{x=L} \]

    \[ \frac{EA}{m_2c_Ap} = \tan \frac{p}{c_A} L \]

Where c_A = \sqrt{\frac{E}{\rho}} and E = c_A^2\rho.

    \[ \frac{c_A^2 \rho A}{m_2c_Ap} = \tan \frac{p}{c_A} L \]

But \rho A L = m – mass of rod. Therefore:

    \[ \frac{pAL}{m_2} \frac{c_A}{pL} = \tan \frac{p}{c_A} L \]

    \[ \frac{m}{m_2} = \gamma \tan \gamma \]

Where \gamma = \frac{pL}{c_A}. We can consider the mass of a “spring” in the calculation of the system natural frequency.

    \[ \begin{split} \gamma &= \frac{pL}{c_A} \\&= pL\sqrt{\frac{\rho}{E}} \end{split} \]

To put this in a form we recognize, recall:

    \[ \frac{P}{A} = E\frac{\Delta L}{L} \]

Therefore:

    \[ \begin{split} P &= \frac{EA}{L} \Delta L \\&= "k"\Delta L \end{split} \]

Therefore:

    \[ \begin{split} \gamma &= pL\sqrt{\frac{A\rho}{kL}} \\&= p \sqrt{\frac{A \rho L}{k}} \\&= p\sqrt{\frac{m}{k}} \end{split} \]

Where m is the mass of the rod or spring. As a result, we must solve this transcendental equation to find p for a given \frac{m}{m_2} ratio.

Example

    \[ \frac{m}{m_2} = \frac{1}{2} \Rightarrow \gamma = 0.6534 \qquad \text{(for lowest }p, p_1\text{)}\]

    \[\begin{split} p_1 &= 0.6534\sqrt{\frac{k}{m}} \\&= 0.6534 \sqrt{\frac{2k}{m_2}} \\&= 0.924\sqrt{\frac{k}{m_2}} \end{split} \]

If we used the rule of thumb:

    \[ \begin{split} p_1 &= \sqrt{\frac{k}{m_2 + \frac{m}{3}}} \\&= \sqrt{\frac{k}{m_2 + \frac{m_2}{6}}} \\&= \sqrt{ \frac{6}{7}\frac{k}{m_2}} \\&= 0.926\sqrt{\frac{k}{m_2}} \end{split} \]

If \frac{m}{m_2} = 1, then:

    \[ \gamma = 0.8603 \]

    \[ p_1 = 0.860\sqrt{\frac{k}{m}} \qquad (m_2 = m)\]

Using R.O.T:

    \[\begin{split} p_1 &= \sqrt{\frac{3}{4} \frac{k}{m_2}} \\&= 0.866\sqrt{\frac{k}{m}} \end{split} \qquad (m_2 = m)\]

Torsional Vibration of Shafts & Rods

We assume that I(x) is the polar second moment of area of the x-section. J_o = I\rho is the moment of inertia (polar) per unit length of the rod, where \rho is the mass density.

    \[T = \mu I \frac{\partial \theta}{\partial x}\]

The equation of motion is

    \[\frac{\partial}{\partial x}[\mu I(x) \frac{\partial \theta}{\partial x}] + \tilde{f} (x,t) = J_o dx \frac{\partial ^{2} \theta}{\partial t^{2}}\]

if we consider a uniform rod for free vibration

    \[\mu I \frac{\partial ^{2} \theta}{\partial x^{2}} = \rho I \frac{\partial ^{2} \theta}{\partial t^{2}}\]

    \[\frac{\partial ^{2} \theta}{\partial x^{2}} = \frac{\rho}{\mu} \frac{\partial ^{2} \theta}{\partial t^{2}}\]

which is again the wave equation, therefore

    \[\theta (x,t) = \left(A \cos p \frac{x}{\tilde{c}} + B \sin p \frac{x}{\tilde{c}}\right)\cos\left(pt-\phi\right)\]

where

    \[\begin{split} \tilde{c} &= \sqrt{\frac{\mu}{\rho}} \\&= \sqrt{\frac{E}{2(1+\nu) \rho}} \\& \approx 0.62 C_A \end{split}\]

We can proceed in a completely analogous manner to the axial vibrations

The torsional stiffness \hat{k} = \frac{\mu I}{L}, J = J_o L = \rho LI, therefore

    \[\begin{split} \frac{\rho}{\mu} &= \frac{J/LI}{\hat{k} L/I} \\&= \frac{J}{\hat{k} L^{2}}\end{split}\]

Define \gamma = p\sqrt{\frac{J}{\hat{k}}} where J is the total moment of inertia of the rod and \hat{k} is the total torsional stiffness

    \[\theta (x,t) = \left(A \cos \frac{\gamma x}{L} + B \sin \frac{\gamma x}{L}\right) \cos (pt + \phi) \]

again we can solve the boundary value problem

The solution is again:

    \[\gamma \tan \gamma = \frac{J}{J_1}\]

if \gamma << 1

    \[\gamma ^{2} \approx \frac{J}{J_1}\]

    \[p^{2} \frac{J}{\hat{k}} = \frac{J}{J_1}\]

    \[p^{2} \approx \frac{\hat{k}}{J_1}\]

Bending Vibrations

The three problems discussed above are all of exactly the same form, the wave equation with two boundary conditions. For the bar in flexure, the bvp is a 4th order d.e.

we neglect the rotation of the element in the deformation and we neglect the shear deformation with respect to the bending deformation. (We will relax these restrictions later.)

    \[\begin{split} + \uparrow \sum F_y &= f(x,t)dx - V(x) + V(x,t) + \frac{\partial V}{\partial x}dx \\&= m(x) \frac{\partial ^{2} y}{\partial t^{2}} \end{split}\]

Since we neglect the rotary inertia

    \[\begin{split} \circlearrowleft + \sum M_x &= 0 \\&= M(x,t) + \frac{\partial M}{\partial x}dx - M(x,t) + V(x)dx + f(x,t) dx \frac{dx}{2} \end{split}\]

    \[\frac{\partial M}{\partial x} + V(x) = 0 \]

    \[-\frac{\partial ^{2} M}{\partial x^{2}} + f(x,t) = m(x) \frac{\partial ^{2} y}{\partial t^{2}}, 0<x<L\]

but the constitutive law gives

    \[M(x,t) = EI(x) \frac{\partial ^{2} y}{\partial x^{2}}\]

Therefore

    \[-\frac{\partial ^{2}}{\partial x^{2}}\bigg(EI(x) \frac{\partial ^{2} y}{\partial x^{2}}\bigg) + f(x,t) = m(x) \frac{\partial ^{2} y}{\partial t^{2}}\]

The boundary conditions are:

1.Clamped ends at x=0

    \[y(0,t) = 0 , \frac{\partial y}{\partial x}\bigg|_{x=0} = 0\]

2.Hinged end

    \[y(0,t) = 0 , EI(x)\frac{\partial ^{2} y}{\partial x^{2}}(x,t) \bigg|_{x=0} = 0 \]

3.Free end

    \[EI(x)\frac{\partial ^{2} y}{\partial x^{2}}(x,t) \bigg|_{x=0} = 0, \frac{\partial }{\partial x}\bigg(EI(x) \frac{\partial ^{2} y}{\partial x^{2}}\bigg)\bigg|_{x=0} = 0\]

We divide these boundary conditions into geometric and natural boundary conditions.

Geometric:

    \[y(0,t) = 0, \frac{\partial y}{\partial x}\bigg|_{x=0} = 0\]

Natural:

    \[EI(x)\frac{\partial ^{2} y}{\partial x^{2}}(x,t) \bigg|_{x=0} = 0, \frac{\partial }{\partial x}\bigg(EI(x) \frac{\partial ^{2} y}{\partial x^{2}}\bigg)\bigg|_{x=0} = 0\]

(This difference is important when dealing with approximate solutions.) Again we look for synchronous solutions.

    \[y(x,t) = Y(x)G(t)\]

This again leads to

    \[G(t) = C \cos(pt - \phi)\]

and an o.d.e

    \[\frac{d^{2}}{dx^{2}} \bigg(EI(x) \frac{d ^{2} Y(x)}{d x^{2}}\bigg) = p^{2}m(x)Y(x), 0<x<L \qquad (*)\]

This is an eigenvalue problem.

Additional b.c.’s

    \[\frac{\partial }{\partial x}\bigg(EI \frac{\partial ^{2} y}{\partial x^{2}}\bigg) = a[ky + c \frac{\partial y}{\partial t} + m \frac{\partial ^{2} w}{\partial t^{2}}]\]

where a = 1 for the right end, and a = -1 for the left end

    \[EI \frac{\partial ^{2} y}{\partial x^{2}} = a\bigg[\hat{k}\frac{\partial y}{\partial x} + \hat{c}\frac{\partial ^{2} y}{\partial x \partial t} + J_o \frac{\partial ^{3} w}{\partial x \partial t^{2}}\bigg]\]

This gives b.c.’s for a torsional spring damper & rotational inertia.

There is no general solution to (*) but there are solutions for certain special cases. For example, consider a uniform bar:

    \[\frac{d^{4}Y(x)}{dx^{4}} - \frac{p^{2}m}{EI} Y(x) = 0\]

let \beta ^{4} = \frac{p^{2}m}{EI}

The general solution of this equation is Y(x) = e^{sx}

    \[e^{sx}(s^{4} - \beta ^{4}) = 0\]

Therefore:

    \[(s^{2}+\beta ^{2})(s^{2}-\beta ^{2}) = 0\]

    \[s = \beta, -\beta, i \beta, -i \beta\]

Therefore:

    \[Y(x) = Ce^{\beta x} + De^{-\beta x} + Ee^{i \beta x} + Fe^{-i \beta x}\]

or equivalently

    \[Y(x) = C_1 \sin \beta x + C_2 \cos \beta x + C_3 \sinh \beta x + C_4 \cosh \beta x\]

Where C_i must be determined from the boundary conditions.

If we consider a simply supported beam then

    \[\textcircled{1} y(0,t) = 0\]

    \[\textcircled{2} \frac{\partial ^{2} y}{\partial x^{2}}(0,t) = 0\]

    \[\textcircled{3} y(L,t) = 0\]

    \[\textcircled{4} \frac{\partial ^{2} y}{\partial x^{2}}(L,t) = 0\]

\textcircled{1} gives:

    \[C_2 + C_4 = 0\]

\textcircled{2} gives:

    \[\beta ^{2} (-C_2 + C_4) = 0\]

Therefore:

    \[C_2 = C_4 = 0\]

\textcircled{3} gives:

    \[C_1 \sin \beta L + C_3 \sinh \beta L = 0\]

\textcircled{4} gives:

    \[ \beta ^{2}(-C_1 \sin \beta L + C_3 \sinh \beta L) = 0\]

Therefore:

    \[C_3 = 0\]

    \[C_1 sin \beta L = 0\]

    \[\beta L = r \pi \qquad r = 1, 2, ...\]

    \[\begin{split} \beta ^{2} &= \bigg(\frac{r \pi}{L}\bigg)^{2} \\&= p \sqrt{\frac{m}{EI}} \end{split}\]

    \[p = r^{2} \pi^{2} \sqrt{\frac{EI}{m L^{4}}}\]

    \[Y_r{x} = C_1 \sin \frac{r \pi x}{L}, r = 1,2, ...\]

If we use \int _0 ^{L} m Y_r ^{2}(x)dx = 1

then the normalized eigenvalue are Y_r(x) = \sqrt{\frac{2}{mL}} \sin \frac{r \pi x}{L}, r= 1,2,...

For the diagram below: p_1 = \pi^{2} \sqrt{\frac{EI}{mL^{4}}}

For the diagram below: p_2 = 4 \pi^{2} \sqrt{\frac{EI}{mL^{4}}}

However if we consider a cantilever beam

    \[\textcircled{1} y(0,t) = 0\]

    \[\textcircled{2} \frac{\partial y}{\partial x}(0,t) = 0\]

    \[\textcircled{3} \frac{\partial ^{2} y}{\partial x^{2}}(L,t) = 0\]

    \[\textcircled{4} \frac{\partial ^{3} y}{\partial x^{3}}(L,t) = 0\]

\textcircled{1} gives:

    \[C_2 + C_4 = 0\]

\textcircled{2} gives:

    \[C_1 + C_3 = 0\]

Therefore:

    \[Y(x) = C_1(\sin \beta x - \sinh \beta x) + C_2(\cos \beta x - \cosh \beta x)\]

Now, use \textcircled{3} and \textcircled{4} respectively

    \[C_1(\sin \beta L + \sinh \beta L) + C_2(\cos \beta L + \cosh \beta L) = 0\]

    \[C_1(\cos \beta L + \cosh \beta L) - C_2(\sin \beta L - \sinh \beta L) = 0\]

For a non-trivial solution the determinant of the coefficients is 0, therefore:

    \[\cos \beta _{r} L \cosh \beta _{r} L = -1\]

This is a transcendental equation for the \beta _{r} L which contains the natural frequencies. In general we can get:

    \[C_2 = \frac{-(\sin \beta _{r} L + \sinh \beta _{r} L)}{(\cos \beta _{r} L + \cosh \beta _{r} L)} C_1\]

or

    \[C_2 = \frac{(\cos \beta _{r} L + \cosh \beta _{r} L)}{(\sin \beta _{r} L - \sinh \beta _{r} L)} C_1\]

but these are equivalent, therefore the mode shapes are:

    \[Y_{r}(x) = C_1 \bigg[(\sin \beta _{r} x - \sinh \beta _{r} x) + \frac{(\cos \beta _{r} L + \cosh \beta _{r} L)}{(\sin \beta _{r} L - \sinh \beta _{r} L)}(\cos \beta _{r} x - \cosh \beta _{r} x)\bigg]\]

which is very complicated

    \[p_1 = (1.875)^{2} \sqrt{\frac{EI}{mL^{4}}}\]

    \[p_2 = (4.694)^{2} \sqrt{\frac{EI}{mL^{4}}}\]

    \[p_3 = (7.855)^{2} \sqrt{\frac{EI}{mL^{4}}}\]

This theory works well for the lower frequencies but not for the higher ones since shear and rotation now become important.

we have suggested that these eigenvalues are orthogonal in some sense. While this can be shown in general using operator theory we will show it in a limited way.

Consider two solutions of the general beam vibration problem.

    \[\textcircled{1} \frac{d^{2}}{dx^{2}}\bigg(EI(x)\frac{d^{2}Y_{r}(x)}{dx^{2}}\bigg) = p_{r}^{2} m(x) Y_{r}(x), 0<x<L\]

    \[\textcircled{2} \frac{d^{2}}{dx^{2}}\bigg(EI(x)\frac{d^{2}Y_{s}(x)}{dx^{2}}\bigg) = p_{s}^{2} m(x) Y_{s}(x), 0<x<L\]

multiply \textcircled{1} by Y_{s}(x) and integrate by parts to obtain:

    \[\int _{0}^{L} Y_s(x) \frac{d^{2}}{dx^{2}}\bigg[EI(x)\frac{d^{2}Y_{r}}{dx^{2}}\bigg] dx\]

    \[\begin{split} &= \bigg[Y_{s}(x) \frac{d}{dx}\bigg\{EI(x)\frac{d^{2}Y_{r}}{dx^{2}}\bigg\}\bigg]_{0}^{L} - \int_{0}^{L} \frac{dY_{s}(x)}{dx} \frac{d}{dx} \bigg\{EI(x)\frac{d^{2}Y_{r}}{dx^{2}}\bigg\} dx \\&= \bigg[Y_{s}(x) \frac{d}{dx}\bigg\{EI(x)\frac{d^{2}Y_{r}}{dx^{2}}\bigg\}\bigg]_{0}^{L} - \bigg[\frac{dY_{s}(x)}{dx} EI \frac{d^{2}Y_{r}}{dx^{2}}\bigg]_{0}^{L} + \int _{0}^{L} EI \frac{d^{2}Y_{s}(x)}{dx^{2}}\frac{d^{2}Y_{r}(x)}{dx^{2}} dx \\&= p_{r}^{2} \int _{0}^{L} m(x) Y_{r}(x) Y_{s}(x) dx \end{split}\]

Now do a similar thing for equation \textcircled{2} i.e multiply by Y_{r}(x) and integrate by parts twice. This gives:

    \[\int _{0}^{L} Y_r(x) \frac{d^{2}}{dx^{2}} \bigg\{EI\frac{d^{2}Y_{s}}{dx^{2}}\bigg\} dx\]

    \[\begin{split} &= \bigg[Y_r(x) \frac{d}{dx} \bigg\{ EI \frac{d^{2}Y_s}{dx^{2}} \bigg\} \bigg]_{0}^{L} - \bigg[\frac{dY_{r}}{dx} EI \frac{d^{2}Y_s}{dx^{2}} \bigg]_{0}^{L} + \int _{0}^{L} EI \frac{d^{2}Y_r}{dx^{2}}\frac{d^{2}Y_s}{dx^{2}} dx \\&= p_{s}^{2} \int _{0}^{L} m(x) Y_{r}(x) Y_{s}(x) dx \end{split}\]

Subtract these two to give:

    \[(p_r^{2} - p_s^{2}) \int _{0}^{L} m(x) Y_{r}(x) Y_{s}(x) dx = 0 ?\]

This is clearly true when we have any combination of clamped, simple support, or free ends. It can also be shown to be true when the ends are spring supported.

Therefore, since the eigenvalues are distinct

    \[\int _{0}^{L} m(x) Y_{r}(x) Y_{s}(x) dx = 0\]

and they are orthogonal with respect to the mass distribution.

While in the discrete case the eigenvectors were also orthogonal with respect to the stiffness matrix, the eigenfunctions are orthogonal with respect to the stiffness EI(x), only in a certain sense. To see this again consider the eigenvalue problem.

    \[\frac{d^2}{dx^2}\bigg(EI(x) \frac{d^{2}Y_{r}(x)}{dx^{2}}\bigg) = p_{r}^{2}m(x)Y_{r}(x)\]

multiply by Y_s(x) and integrate:

    \[\int _{0}^{L} Y_s(x) \frac{d^{2}}{dx^{2}} \bigg(EI(x)\frac{d^{2}Y_{r}(x)}{dx^{2}}\bigg) = p_{r}^{2} \int m(x)Y_{r}(x) Y_{s}(x) dx\]

But the RHS is zero. Therefore:

    \[\int_0^L Y_s(x)\frac{d^2}{dx^2}\bigg(EI\frac{d^2Y_r(x)}{dx^2}\bigg)dx = 0, \ r\neq s\]

This is somewhat complicated but again we integrate by parts twice to get:

    \[\int_0^L EI(x)\frac{d^2Y_r(x)}{dx^2}\frac{d^2Y_s(x)}{dx^2}dx=0\]

If we use the combination of natural combination of natural and/or geometric b.c.’s.

This says that the second derivative of the eigenfunctions are orthogonal with respect to the stiffness, not the eigenfunctions themselves

When r = s, we normalize the eigonfunctions as:

    \[\int_0^Lm(x)Y_r(x)Y_s(x)dx = \delta_{rs}\]

then:

    \[\int Y_s(x)\frac{d^2}{dx^2}\bigg(EI\frac{d^2Y_r(x)}{dx^2}\bigg)dx =p^2_r\delta_{rs}\]

If we have one end which has a concentrated mass M then these results must be modified.

For example the corresponding orthogonality becomes:

    \[\int_0^Lm(x)Y_r(x)Y_s(x)dx+MY_r(L)Y_s(L) = 0, \ r\neq s\]

NOTE:

1. The integrals are symmetric with respect to r & s. This is similar to the requirement that m_{rs} = m_{sr} and k_{ij} = k_{ji} for discrete systems.

2. Again there could be rigid body modes which means the system is positive semi-definite rather than positive definite. Again we would have to use conservation of linear momentum to exclude these modes.

All this leads to the expansion theorem for continuous systems.

Any function Y(x), satisfying the boundary conditions of the problem and such that \frac{d^2}{dx^2}\bigg[EI(x)\frac{d^2Y(x)}{dx^2}\bigg] is a continuous function can be represented by the absolutely and uniformly convergent series of the system eigenfunctions.

    \[Y(x) = \sum_{x=1}^\infty C_rY_r(x)\]

where:

    \[C_r = \int_0^Lm(x)Y(x)Y_r(x)dx, \ r = 1,2,...\]

This is essentially a “generalized” Fourier expansion. In fact when the eigenfunctions are harmonic the expansion theorem does reduce to a Fourier series representation.

The response of a system to initial excitation external excitation or both can be obtained conveniently by modal analysis. The method is entirely analogous to that for discrete systems which means we must first find the eigenvalues and eigenfunctions. This is of course not trivial.

As an illustration, consider a bar with both edges hinged and with initial conditions:

    \[y(x,0) = y_0(x)\]

    \[\frac{\partial y(x,t)}{\partial t}\bigg|_{t=0} v_0(x)\]

The original boundary value problem for a uniform bar is:

    \[-EI\frac{\partial^4y}{\partial x^4} +f(x,t) = m\frac{\partial^2y(x,t)}{\partial t^2}, \ 0<x<L\]

    \[y(0,t) = 0, \quad EI\frac{\partial^2y}{\partial x^2}\bigg|_{x = 0} = 0\]

    \[y(L,t) = 0, \quad EI\frac{\partial^2y}{\partial x^2}\bigg|_{x = L} = 0\]

We also know the solution to the eigenvalue problem for natural frequencies and mode shapes is:

    \[p_r = (r\pi)^2\sqrt{\frac{EI}{mL^4}}, \ r = 1,2,...\]

    \[Y_r(x) = \sqrt{\frac{2}{ml}}\sin\frac{r\pi x}{L}\]

We let the solution of the original equation be of the form:

    \[y(x,t) = \sum_{r = 1}^\infty Y_r(x)\eta_r(t) \quad (\star)\]

This is analogous to the situation in discrete systems in which we used the modal matrix [\mu] to transform the vector of coordinates q(t) into the normal coordinates [\eta(t)]. Here instead of the matrix [\mu] we use the eigenfunctions of the problem.

If we put (\star) into:

    \[-EI\frac{\partial^4y}{\partial x^2} + f(x,t) = m\frac{\partial^2y}{\partial t^2}\]

we get:

    \[\sum_{r=1}^\infty \ddot\eta_r(t)Y_r(x) + \sum_{r=1}^\infty\eta_r(t)EI\frac{d^4Y_r(x)}{dx^4} = f(x, t)\]

Multiply by Y_s(x) and integrate with respect to x to get:

    \[\ddot\eta_r(t) + p^2_r\eta_r(t) = Q_r(t)\]

Where Q_r(t) = \int_0^Lf(x,t)Y_r(x)dx are the generalized forces associated with the normal coordinates \eta_r(t). The resulting differential equation is again a SDOF system and its general solution is:

    \[\eta_r(t) = \frac{1}{p_r}\int Q_r(\tau)\sin p_r(t-\tau)d\tau + \eta_{r0}\cos p_rt + \frac{\dot\eta_{r0}}{p_r}\sin p_rt\]

where \eta_{r0} = \eta_r(0),\ \dot\eta_{r0} = \dot\eta_r(0), where again:

    \[\eta_{r0} = \int_0^Lmy_0(x)Y_r(x)dx,\ r = 1,2,...\]

    \[\dot\eta_{r0} = \int_0^Lm\dot y_0(x)Y_r(x)dx,\ r = 1,2,…\]

The general response is obtained by inserting \eta_r(t) back into (\star).

Example

    \[f(x,t) = F_0H(t)\]

    \[y(x,0) = 0, \ \dot y(x,0) = 0\]

    \[ \begin{split}Q_r(t) &= \int_0^LF_0Y_r(x)dx \\&= F_0H(t)\sqrt{\frac{2}{mL}}\int_0^L\sin\frac{r\pi x}{L} dx \\&= F_0H(t)\sqrt{\frac{2}{mL}}(1-\cos r\pi)\frac{L}{\pi r}, \ r = 1,2,3\\&= \frac{2LF_0H(t)}{\pi r}\sqrt{\frac{2}{mL}}, \ r = 1,3,5 \end{split}\]

Therefore only odd modes excited by uniform load.

    \[\begin{split}\eta_r &= \frac{1}{p_r}2\sqrt{\frac{2}{mL}}\left(\frac{LF_0}{\pi r}\right)\int_0^t\sin p_r(t-\tau)d\tau \\&= \frac{1}{p_r}2\sqrt{\frac{2}{mL}}\left(\frac{LF_0}{\pi r}\right)\frac{\cos p_r(t-\tau)}{p_r}\bigg]_0^t \\&= \frac{1}{p_r}2\sqrt{\frac{2}{mL}}\left(\frac{LF_0}{\pi r}\right)\bigg[1-\cos p_r t\bigg]\end{split}\]

Therefore:

    \[ y(x,t) = \frac{4F_0L^4}{\pi^5EI}\sum_{r = 1}^\infty \frac{1}{(2r-1)^5}\sin\frac{(2r-1)\pi x}{L}\bigg[1-\cos(2r-1)^2\pi^2\sqrt{\frac{EI}{mL^4}}t\bigg] \]

(NOTE: r = 1, 2, \dots \infty now).

It is easy to see that the first mode is predominant.

1st Term:

    \[y_1 = \frac{4F_0L^4}{\pi^5EI}\frac{1}{1^5}\sin\frac{\pi x}{L}\bigg[ 1 - \cos\pi^2\sqrt{\frac{EI}{mL^4}}t\bigg]\]

    \[y_{1\text{(MAX)}} = \frac{4F_0L^4}{\pi^5EI}(2) \]

    \[\begin{split}y_{2\text{(MAX)}} &= \frac{4F_0L^4}{\pi^5EI}\bigg(\frac{1}{3^5}\bigg)(2) \\&= 0.004y_{1\text{(MAX)}}\end{split}\]

    \[y_{\text{STATIC(MAX)}} = \frac{5}{384}\frac{F_0L^4}{EI}\]

    \[\frac{y_{\text{DYNAMIC}}}{y_{\text{STATIC}}}\bigg|_{\text{MAX}} = \frac{2\cdot 4}{\pi^5}\left(\frac{384}{5}\right) =\underline{\underline{2.01}}\]

Effect of Axial Force

    \[\theta = \frac{\partial y}{\partial x}\]

    \[\bigg(V+\frac{\partial V}{\partial x}dx\bigg) - V + f(x,t)dx + p+\frac{\partial p}{\partial x}dx\bigg(\frac{\partial y}{\partial x} +\frac{\partial^2 y}{\partial x^2}dx\bigg) - p\frac{\partial y}{\partial x} = mdx\frac{\partial^2 y}{\partial x^2}\]

From the moment equation:

    \[\frac{\partial M}{\partial x} = -V\]

    \[\theta = \frac{\partial y}{\partial x} \]

    \[ \theta + N\theta = \frac{\partial y}{\partial x} + \frac{\partial^2y}{\partial x^2}dx\]

    \[M = EI\frac{\partial^2y}{\partial x^2}\]

    \[\frac{\partial^2}{\partial x^2}\bigg( EI(x)\frac{\partial^2 y}{\partial x^2}\bigg) + m\frac{\partial^2 y}{\partial t^2} - p\frac{\partial^2y}{\partial x^2} = f(x,t)\]

Therefore:

Again we look for the eigenvalues and eigenfunctions by looking at a variables separable solution.

    \[ y(x,t) = Y(x)G(t)\]

Therefore:

    \[y(x,t) = Y(x)\cos(pt - \phi)\]

And for uniform beams:

    \[IE\frac{d^4Y}{dx^4} - mp^2Y(x) - p\frac{d^2Y}{dx^2} = 0\]

    \[Y(x) = e^{sx}\]

    \[s^4 - \frac{p}{EI}s^2 - \frac{mp^2}{EI} = 0\]

    \[s^2 = \frac{p}{2EI}\pm\sqrt{\frac{p^2}{4E^2I^2}+ \frac{mp^2}{EI}}\]

Therefore:

    \[S^2_1 = \frac{p}{2EI} + \sqrt{\frac{p^2}{4R^2I^2} + \frac{mp^2}{EI}} = \delta^2 \ \text{(positive)}\]

    \[S^2_2 = \frac{p}{2EI} + \sqrt{\frac{p^2}{4R^2I^2} - \frac{mp^2}{EI}} = -\gamma^2 \ \text{(negative)}\]

Therefore:

    \[ S_1 = -\delta^2, \quad S_2 = \pm i\gamma\]

Therefore:

    \[Y(x) = C_1\cosh \delta x + C_2\sinh \delta x + C_3\cos\gamma x + C_4\sin \gamma x \]

For the S.S case:

    \[Y(0) = Y(L) = 0, \ Y''(0) = Y''(L) = 0\]

    \[C_1 +C_3 = 0, \ C_1\delta^2 - C_3\gamma^ 2 = 0, \ C_3 = -C_1\]

Therefore:

    \[C_1(\delta^2 +\gamma^2) = 0\]

Therefore:

    \[C_1\bigg[ \frac{p}{2EI} + \sqrt{\frac{p^4}{4E^2I^2} + \frac{mp^2}{EI}} - \frac{p}{2EI} + \sqrt{\frac{p^4}{4E^2I^2} + \frac{mp^2}{EI}}\bigg] = 2C_1\sqrt{\frac{p^4}{4E^2I^2} + \frac{mp^2}{EI}} = 0\]

Therefore:

    \[C_1 = 0, \ C_3 = 0\]

    \[Y(L) = C_2\sinh\delta L + C_4\sin \gamma L = 0\]

    \[Y''(L) = C_2\delta^2\sinh\delta L - C_4\gamma^2\sin\gamma L = 0\]

Multiply by \gamma^2 and add to get:

    \[C_2(\gamma^2 + \delta^2)\sinh\delta L = 0\]

    \[ \implies C_2 = 0\]

Therefore:

    \[\sin\gamma L = 0\]

    \[ \gamma L = r\pi, \quad r = 1,2,3...\]

    \[ \gamma^2 = \bigg(\frac{r\pi}{L}\bigg)^2 = -\frac{p}{2EI} + \sqrt{\frac{p^2}{4E^2I^2} + \frac{mp^2}{EI}}\]

    \[\frac{p}{2EI} + \bigg(\frac{r\pi}{L}\bigg)^2 = \sqrt{\frac{p^2}{4E^2I^2} + \frac{mp^2}{EI}}\]

    \[\frac{p^2}{4E^2I^2} + \frac{2p}{2EI}\bigg(\frac{r\pi}{L}\bigg)^2 + \frac{r^4\pi^4}{L^4} = \frac{p^2}{4E^2I^2} + \frac{mp^2}{EI}\]

Therefore :

    \[\begin{split}p^2 &= \frac{EI}{m}\bigg[\frac{r^4\pi^2}{L^4} +\frac{r^2\pi^2}{L^4}\frac{p}{EI}\bigg] \\&= \frac{EI\pi^4}{mL^4}\bigg[ r^4 + \frac{r^2pl^2}{\pi^2EI}\bigg]\end{split}\]

    \[p = \pi^2\sqrt{\frac{EI}{ml^4}}\bigg[r^4 + \frac{r^2pL^2}{\pi^2EI}\bigg] ^{\frac{1}{2}}\]

Special Cases:

  1. if p = 0:

        \[ p = r^2\pi^2\sqrt{\frac{EI}{ml^4}} \ \text{(S.S beam as before)}\]

  2. if EI \rightarrow 0:

        \[ p = \pi^2\sqrt{\frac{EI}{ml^4}}\frac{r}{\pi}\sqrt{\frac{pL^2}{EI}}\]


        \[p = \frac{r\pi}{L}\sqrt{\frac{p}{m}} = \frac{r\pi}{L}\sqrt{\frac{p}{\rho}} \ \text{(string result)}\]

  3. p>0, then the natural frequencies are higher than without axial force.
  4. p<0, then the natural frequencies decrease from that without axial force.
    as p \rightarrow -\frac{r^2\pi^2EI}{L^2}, p \rightarrow 0.

This is the state buckling load of the beam. The effect of axial load can be used as a technique to estimate the load in a particular member.

Approximate Methods for Beams

It is obvious that only a few problems can be solved analytically for beam vibrations. As a result, approximate methods, including FEA, are necessary. There are a few approximate methods that can be used to give a crude approximation and require relatively little work.

We well consider the Rayleigh’s Quotient and from it Dunkerley’s Method. These allow bounding of the lowest natural frequency and can be used as a starting point for more sophisticated techniques. They are based on an energy approach similar to that used for finite dimensional problems.

Rayleigh’s Quotient

Again we find the kinetic and potential energy assuming SSHM for the particles of a beam.

The kinetic energy of dx is

    \[ \frac{1}{2}mdx\Big(\frac{\partial y}{\partial t}\Big)^2 \]

Therefore the total kinetic energy of the beam is

    \[ T = \int_L \frac{1}{2}m\Big(\frac{\partial y}{\partial t}\Big)^2dx \]

where we integrate over the entire beam.

To find the potential energy (stored in the beam) consider the bending of an element from its static shape.

The work done is \frac{1}2Md\theta

    \[ \theta \approx \frac{\partial y}{\partial x} \]

Therefore:

    \[ d\theta = \frac{\partial^2y}{\partial x^2}dx \]

    \[ M = EI \frac{\partial^2y}{\partial x^2} \qquad \text{(as before)} \]

Therefore work done on the element

    \[ = \frac{1}{2}EI\Big(\frac{\partial^2 y}{\partial x^2}\Big)^2dx \]

and the energy stored is the total work done on the entire beam

    \[ U = \int_L \frac{1}{2}EI\Big(\frac{\partial^2 y}{\partial x^2}\Big)^2dx \]

Rayleigh again assume SSHM and T_{\text{MAX}} = U_{\text{MAX}}

Therefore assume y(x,t) = Y(x)\sin(pt+\phi)

    \[ \begin{split} \frac{\partial^2y}{\partial x^2} &= \frac{d^2Y}{dx^2}\sin(pt+\phi) \\ \frac{\partial y}{\partial t} &= -Y(x)\cos(pt+\phi) \end{split} \]

    \[ \begin{split} U_{\text{MAX}} &= V_{\text{MAX}} \\ U_{\text{MAX}} &= \frac{1}{2}\int_0^L EI \Big(\frac{d^2Y}{dx^2}\Big)^2dx \\ V_{\text{MAX}} &= \frac{1}{2}p^2\int_0^L m(Y)^2dx \\ p^2 &= \frac{\int_0^L EI \Big(\frac{d^2Y}{dx^2}\Big)^2dx}{\int_0^L m(Y)^2dx} \end{split} \]

but we must have the exact mode shape to get the right answer. It can be shown that we can get an upper bound for p_1^2 if we select a shape function (mode shape) which satisfies the geometric b.c.

    \[ p^2 \leq \frac{\int_L EI \big(\frac{d^2Y}{dx^2}\big)^2dx}{\int_L m \{Y\}^2dx} \]

e.g. Cantilever Beam

    \[ \begin{split} Y(x) &= Ax^2 \\ Y'(x) &= 2Ax \\ Y''(x) &= 2A \end{split} \]

    \[ \begin{split} p_1^2 \leq \frac{\int_0^L EI(4A)^2dx}{\int_0^L mA^2x^4dx} &= \frac{4EIL}{mx^5/5\big]_0^L} \\ &= \frac{20EI}{mL^4} \end{split} \]

Therefore:

    \[ \begin{split} p \leq \sqrt{20} \sqrt{\frac{EI}{mL^4}} &= 4.47 \sqrt{\frac{EI}{mL^4}} \\ p_{\text{EXACT}} &= 3.52 \sqrt{\frac{EI}{mL^4}} \end{split} \]

    \[ \begin{split} Y &= Ax^3 \\ Y' &= 3Ax^2 \\ Y'' &= 6Ax \end{split} \]

    \[ \begin{split} p_1^2 \leq \frac{EI\int_0^L 36x^2 dx}{m\int_0^L x^6dx} &= \frac{EI 12L^3}{m L^7/7} \\ p_1^2 \leq \frac{84EI}{mL^4} \\qquad \text{Worse} \end{split} \]

We can extend the idea of Rayleigh to the case in which there are point masses and springs attached at a point on a beam.

Consider first a point mass, M_1 at position x_1

It has kinetic energy T where

    \[ \begin{split} T &= \frac{1}{2}m\Big(\frac{\partial y}{\partial t}\Big|_{x=\hat{x}}\Big)^2 \\ T_{\text{MAX}} &= \frac{1}{2}mp^2\{Y(\hat{x}\}^2 \end{split} \]

when we assume

    \[ y(x,t,) = Y(x)\sin(pt+\phi) \]

For a spring at position x^*

The potential energy = \frac{1}{2}k\{y(x^*)\}^2

Again if we assume y(x) = Y(x)\sin(pt+\phi)

    \[U_{\text{MAX}} = \frac{1}{2}k\{Y(x^*)\}^2 \]

We simply add these energies to the T_{\text{MAX}} and U_{\text{MAX}} of the beam to get

    \[ p_1^2 \leq \frac{\int_L EI\big(\frac{d^2Y}{dx^2}\big)^2dx + \sum_{i=1}^n k_i{Y(x_i)}^2}{\int_L m{Y(x)}^2dx + \sum^m_{j=1} M_j {Y(x_j)}^2} \]

for n springs and m masses

Therefore the integrations are over the entire length of the beam.

We can consider beams width varying X section and varying mass distribution.

Dunkerley’s Method for Beams is an extension of Rayleigh’s Quotient where there are no additional springs.

Consider a beam carrying several masses

The Rayleigh Quotient is

    \[ p_1^2 = \frac{\int_L EI\big(\frac{d^2Y}{dx^2}\big)^2dx}{\int_L m\{Y\}dx + \sum_{i=1}^r M_i \{Y^*(x_i)\}^2} \]

Assume that Y^*(x) is the actual shape function for the combined system so that p_1^2 is the actual natural frequency.

Therefore:

    \[ \frac{1}{p_1^2} = \frac{\int_L m(Y^*)^2 dx}{\int_L EI\big(\frac{d^2Y^*}{dx^2}\big)^2dx} + \frac{\sum^r_{i=1}M_i\{Y^*(x_i)\}^2}{\int_L EI\big(\frac{d^2Y^*}{dx^2}\big)^2dx} \]

The first term is an approximation to the fundamental frequency of the beam alone

Therefore

    \[ p^2_{1(B)} \leq \frac{\int_L EI\Big(\frac{d^2Y^*}{dx^2}\Big)^2dx}{\int_L m\{Y^*\}^2dx} \]

    \[ \frac{1}{p^2_{1(B)}} \geq \frac{\int_0^L m{Y^*}^2dx}{\int_0^L EI\Big(\frac{d^2Y^*}{dx^2}\Big)^2dx} \]

since Y^* is the exact answer for the combined system.

Similarly to each of the added masses.

    \[ \frac{1}{(p_{1i})^2} \geq \frac{M_i\{Y^*(x)\}^2}{\int_0^L EI\big(\frac{d^2Y^*}{dx^2}\big)^2dx} \]

Therefore:

    \[ \frac{1}{p_1^2} \leq \frac{1}{p_{1(B)}^2} + \sum^r_{i=1}\frac{1}{(p_1i)^2} \]

This allows the calculation of the lower bound for the natural frequency p_1

Example

S.S beam of the mass/length m carrying M at midspan

    \[ p^2_{1(B)} = \frac{\pi^4EI}{mL^4} \qquad \qquad \Delta_{L/2} = \frac{M_gL^3}{48EI} \quad (\text{static}) \]

Therefore

    \[ p^2_1 = \frac{g}{\Delta} \]

Therefore

    \[ \begin{split} \frac{1}{p_1^2} & \leq \frac{mL^4}{\pi^4EI} + \frac{ML^3}{48EI} \\ &= \frac{48mL^4 + \pi^4ML^3}{48\pi^4EI} \end{split} \]

    \[ p_1^2 > \frac{48EI}{L^3\big[M + \frac{48}{\pi^4}mL\big]} = \frac{48EI}{L^3\big[M + 0.493mL\big]} \]

Use Rayleigh to get upper bound. Try:

    \[ \begin{split} Y(x) &= Ax(L-x) \\ &= AxL - Ax^2 \\ Y'(x) &= AL - 2Ax \\ Y''(x) &= -2A \end{split} \]

    \[ \begin{split} Y\Big(\frac{L}{2}\Big) &= A\frac{L}{2}\Big(L- \frac{L}{2}\Big) \\ &= \frac{AL^2}{4} \end{split} \]

    \[ \begin{split} p^2_{1(B)} & \leq \frac{\int^L_0 EI(4A^2)dx}{\int_0^L m(A^2x^2L^2 + A^2x^4 - 2A^2x^3L)dx + M\frac{A^2L^4}{16}} \\ &= \frac{EI4L}{m\big[ \frac{x^3L^2}{3} + \frac{x^5}{5} - \frac{2x^4L}{4} \big]^L_0 + \frac{ML^4}{16}} \\ &= \frac{4EI}{mL^4\big[\frac{1}{3} + \frac{1}{5} - \frac{1}{2}\big] + \frac{ML^3}{16}} \\ &= \frac{4EI}{mL^4 \big[ \frac{10+6-15}{30} \big] + \frac{ML^3}{16}} \\ &= \frac{4EI}{\frac{mL^4}{30} + \frac{ML^3}{16}} \\ &= \frac{48EI}{L^3\big[\frac{3}{4}M + \frac{2}{5}mL\big]} \end{split} \]

Let M=mL

    \[ \begin{split} p_1^2 &> \frac{48EI}{1.493mL^4} = 32.15 \frac{EI}{mL^4} \\ p_1 &> 5.67\sqrt{\frac{EI}{mL^4}} \end{split} \]

    \[ \begin{split} p_1^2 &< \frac{48EI}{\frac{23}{20}mL^4} = 41.74 \frac{EI}{mL^4} \\ p_1 &< 6.46 \sqrt{\frac{EI}{mL^4}} \end{split} \]

Vibration of a Free-Free Beam

The bending vibration of a beam are described by the following equation:

    \[ EI \frac{\partial ^4 y}{\partial x^4} + \rho A \frac{\partial ^2 y}{\partial t^2} = 0 \quad (1) \]

E,I,\rho, A are respectively the Young Modulus, second moment of area of the cross section, density and cross section area of the beam. L is the length of the beam. The solution of Eq.(1) can be written as a standing wave^1 y(x,t) = w(x)u(t), separating the spatial and temporal component. This leads to the following characteristic equation that relates the circular frequency \omega to the wavenumber k:

    \[\omega ^2 = \frac{EI}{\rho A}k^4 \quad (2)\]

The spatial part can be written as:

    \[w(x) = C_1 \sin(kx) + C_2 \cos(kx) + C_3\sinh(kx) + C_4\cosh(kx) \quad (3)\]

For a Free-Free Beam the boundary conditions are (vanishing of force and moment):

      \begin{equation*} \begin{cases} w''(0)=0\\ w''(0)=0\\ w''(L)=0\\ w''(L)=0 \end{cases} \text{to get} \begin{cases} -C_2 + C_4 = 0\\ -C_1 + C_3 = 0\\ -C_1 \sin(kL) - C_2 \cos(kL) + C_3 \sinh(kL) + C_4 \cosh(kL)=0\\ -C_1 \cos(kL) + C_2 \sin(kL) + C_3 \cosh(kL) + C_4 \sinh(kL)=0\\ \end{cases} \end{equation*}

——————————————————————

^1 A standing wave or stationary wave is a wave ‘frozen’ in the space and vibrating in time. It result by the sum of two identical waves travelling in opposite directions: 

    \[\begin{split} y(x,t) &= y_+ (x,t) + y_-(x,t) \\&= Y_o \sin (\omega t -kx) + Y_o \sin (\omega t + kx) \\&= 2Y_o \sin (kx) \cos(\omega t) \\&= w(x)u(t) \end{split}\]

Using the first two equations the 3^{rd} and 4^{th} can be arranged in matrix form:

      \begin{equation*} \begin{bmatrix} \sinh(kL) - \sin(kL) & \cosh(kL) - \cos(kL) \\ \cosh(kL) - \cos(kL) & \sin(kL) + \sinh(kL) \\ \end{bmatrix} \begin{bmatrix} C_1 \\ C_2 \\ \end{bmatrix} = \begin{bmatrix} 0 \\ 0 \\ \end{bmatrix} \quad (4) \end{equation*}

For a non trivial solution the determinant of the matrix has to vanish to get:

    \[\cosh(kL) \cos(kL) = 1 \quad (5)\]

The transcendental Eq.(5) has infinite solutions, it can be solved numerically, the first five values are reported here:

Mode order nknL
00
14.7300
27.8532
310.9956
414.1371
517.2787

Putting these values back in Eq. (5) gives the modeshapes corresponding to the natural frequencies \omega _n, that can be calculated from the characteristic Eq. (2). The mode shapes are given by the following: 

    \[\omega _n (x) = [\sinh(k_nx) +\sin(k_nx)] + \frac{\sin(k_nL) - \sinh(k_nL)}{\cosh(k_nL) - \cos(k_nL)}[\cosh(k_nx + \cos(k_nx)] \quad (6)\]

Figure 1: First 5 mode shapes for a free-free beam

While many musical instruments use vibrations as the source of the sound they produce, one of simplest that is very commonly used outdoors is the wind chimes. As illustrated the classical version uses a number (usually 5 or 6) tubes of various lengths and held at various positions that are struck with a central clapper near the center of the tubes. To accomplish this various length tubes are positioned at different lengths relative to the supporting structure. The illustrated chimes are composed of 5 tubes corresponding to a piano’s black keys, making up the so-called minor pentatonic scale. These were first studied by the ancient Greeks (including Pythagoras) who realized the link between the length of a vibrating body and the notes of a musical scale. Actually, the tubes can be treated as free-free beams which are supported at certain positions. A do-it-yourself magazine suggested building your own using 5 copper tubes (3/4 inch type M with an OD = 0.875 in., ID = 0.811 in., weight of 0.328 lbs./ft and Modulus of Elasticity 23 X 10 psi) with: 

Length(in.)Hang Point (from top) (in.)Musical pitch
11 1/22 9/16C-sharp
10 7/82 7/16D-sharp
102 1/4F-sharp
9 7/162 1/8G-sharp
8 7/82A-sharp

a) Calculate the first natural frequency of the five tubes of various lengths

b) Calculate the ratio of the position of the hang point to the total length of each tube. Why is this ratio essentially the same?

c) An 88 key piano keyboard (numbered from 1-88) has standard fundamental frequencies (available on the internet). Find the key numbers that correspond to the 5 tube free-free beam frequencies and show their position on the standard piano keyboard. 

How to Make Wind Chimes in 5 Easy Steps 

1.) Gather Materials

Round up at least 5 feet of Type M 3/4-inch copper tubing, seven eye screws, five No. 6 1-inch machine screws and nuts, nylon twine, and 1 x 6 lumber.

2.) Map Pipe Mounts

Center a 4 1/2-inch-diameter circle in a 5 1/2-inch-square cut of lumber. Mark the circle at five equidistant points. Insert eye screws at the circle’s center and at all five points.

3.) Cut Pythagoras 

Cut five pieces of tubing to the lengths in the table below and deburr. The chime’s five notes, which correspond to a piano’s black keys, make up the minor pentatonic scale. The notes are pleasing in any order. The ancient Greeks such as Pythagoras were the first to study the link between the length of a vibrating body and the notes of a musical scale. 

Musical PitchLength 3/4-inch copperHang Point from top
C-sharp11 1/2 inches2 9/16 inches
D-sharp10 7/8 inches2 7/16 inches
F-sharp10 inches2 1/4 inches
G-sharp9 7/16 inches2 1/8 inches
A-sharp8 7/8 inches2 inches
4.) Tie With Twine 

Drill a 5/32-inch hole through each pipe as listed in the table. These hanging points produce the best chime resonance. Insert a machine screw through the hole and fit a nut onto the screw shank. Tie a 7-inch length of twine from the circle of eye screws to the screw shank in each pipe. 

5.) Hang the Chime 

Use a 2 1/2-inch holesaw to cut a clapper from a 1x scrap. Use more 1x waste to make a V-shaped, 3-inch-long wind scoop. Hang each from the center eye screw. Cut and glue two smaller 1x squares to the top of the first square. Top-center-mount an eye screw in the smallest square. Hang the chime in the breeze and enjoy. 

Table 1: Dimensions and Physical Characteristics of Copper Tube: TYPE M
Table 2: Dimensions and Physical Characteristics of Copper Tube: DWV (Drain, Waste and Vent)

Piano key frequencies 

This is a list of the fundamental frequencies in hertz (cycles per second) of the keys of a modern 88-key standard or 108-key extended piano in twelve-tone equal temperament, with the 49th key, the fifth A (called A4), tuned to 440 Hz (referred to as A440). 12 Since every octave is made of twelve steps and equals two times the frequency (for example, the fifth A is 440 Hz and the higher octave A is 880 Hz), each successive pitch is derived by multiplying (ascending) or dividing (descending) the previous by the twelfth root of two (approximately 1.059463). For example, to get the frequency a semitone up from A4 (A#4), multiply 440 by the twelfth root of two. To go from A4 to B4 (up a whole tone, or two semitones), multiply 440 twice by the twelfth root of two (or just by the sixth root of two, approximately 1.122462). To go from A4 to C5 (which is a minor third), multiply 440 three times by the twelfth root of two, (or just by the fourth root of two, approximately 1.189207). For other tuning schemes refer to musical tuning. 

This list of frequencies is for a theoretically ideal piano. On an actual piano the ratio between semitones is slightly larger, especially at the high and low ends, where string stiffness causes inharmonicity, i.e., the tendency for the harmonic makeup of each note to run sharp. To compensate for this, octaves are tuned slightly wide, stretched according to the inharmonic characteristics of each instrument. This deviation from equal temperament is called the Railsback curve. The following equation gives the frequency f of the nth key, as shown in the table:

The following equation gives the frequency f of the nth key, as shown in the table:

    \[f(n) = (\sqrt[12]{2})^{n-49} \times 440 \text{Hz}\]

(a’ = A4 = A440 is the 49th key on the idealized standard piano)

Alternatively, this can be written as:

    \[f(n) = 2^{\frac{n-49}{12}} \times 440 \text{Hz}\]

Conversely, starting from a frequency on the idealized standard piano tuned to A440, one obtains the key number by: 

    \[\begin{split} n &= 12 \log_2(\frac{f}{440\text{Hz}}) + 49 \\&= 39.86 \log _{10} (\frac{f}{440\text{Hz}}) + 49 \end{split}\]

An 88-key piano, with the octaves numbered and Middle C (cyan) and A440 (yellow) highlighted 

Values in bold are exact on an ideal piano. Keys shaded gray are rare and only appear on extended pianos. The normally included 88 keys have been numbered 1-88, with the extra low keys numbered 89-97 and the extra high keys numbered 98-108. (A 108-key piano that extends from Co to B8 was first built in 2018 by Stuart & Sons.

interactive piano frequency table (http://shakahara.com/pianopitch2.php) – A PHP script allowing the reference pitch of A4 to be altered from 440 Hz.

PySynth (https://web.archive.org/web/20190510200224/https://mdoege.github.io/PySynth/) – A simple Python-based software synthesizer that prints the key frequencies table and then creates a few demo songs based on that table.

“Keyboard and frequencies (http://www.sengpielaudio.com/calculator-notenames.htm)”, SengpielAudio.com. 

Notefreqs (http://www.deimos.ca/notefreqs) – A complete table of note frequencies and ratios for midi, piano, guitar, bass, and violin. Includes fret measurements (in cm and inches) for building instruments. 

Solution:

For a free-free beam the natural frequencies are found from:

    \[\cosh(kL)\cos(kL) = 1\]

and therefore kL = 0 is a solution corresponding to a rigid body mode with no oscillation. The first mode is with (kL)_1 = 4.7300 and the frequency in Hz is:

    \[f_1 = \frac{(4.73)^2}{2\pi L^2} \sqrt{\frac{EI}{m}}\]

where E = modulus of elasticity, I section moment of area of the X section (bending), m = mass/unit length & L is the length of the beam for annulus \times sections I = \frac{\pi}{64}(d_o^4 - d_i^4)

a) For the copper tubing

b)

    \[m = \frac{0.328}{12 \times 386} \frac{\text{lb} -\text{sec}^2}{\text{in}^2}\]

    \[E = 23 \times 10^6\]

c)

    \[I = \frac{\pi}{64}((0.875)^4 - (0.811)^4) = 5.297 \times 10^{-3} \text{m}^4\]

For the 10 inch tube:

    \[\begin{split} f_1 &= \frac{22.373}{2\pi (10)^2} \sqrt{\frac{23 \times 10^6 \times 5.297 \times 10^{-3} \times 12 \times 386}{0.328}} \\&= 1477 \text{Hz} \end{split}\]

as the other tubes are the dimension except for length

Length (m)Hang Pt / Length RatioMusical Pitch(Hz)Musical Pitch (Desig)Key
8 7/80.2251875A#74
9 7/160.2251658G#72
100.2251477F#70
10 7/80.2241248D#67
11 1/20.2231117C#65

The hang point corresponds to the nodal point of the fundamental frequency mode shape as shown.