CIE4140_Lecture_13_Part_2_Python

Let us calculate the natural and resonanse frequencies of a string with arbitrary (though linear) boundary conditions

The universal recipe is to substitute $w(x,t)=W(x) e^{i \omega t})

Strictly speaking, the natural frequencies are defined for undamped systems only. Therefore, we assume that there is no damping along the string nor in the boundaries.

import sympy as sp

The equation of motion in this case reads

w = sp.symbols('w',cls=sp.Function,real=True)
x,t = sp.symbols('x,t',real=True)
c = sp.symbols('c',real=True)
EQM = sp.diff(w(x,t),t,2)-c**2*sp.diff(w(x,t),x,2)
display(EQM)

\[\displaystyle - c^{2} \frac{\partial^{2}}{\partial x^{2}} w{\left(x,t \right)} + \frac{\partial^{2}}{\partial t^{2}} w{\left(x,t \right)}\]

Substituting \(w(x,t)=W(x) e^{i \omega t}\), we obtain the following ordinary differential equation, in which \(\omega\) is one of the natural frequencies

W = sp.symbols('W',cls=sp.Function,real=True)
omega = sp.symbols('omega',real=True)
EQM_Fr = sp.simplify(EQM.subs(w(x,t),W(x)*sp.exp(sp.I*omega*t)))
display(EQM_Fr)

\[\displaystyle - c^{2} e^{i \omega t} \frac{d^{2}}{d x^{2}} W{\left(x \right)} - \omega^{2} W{\left(x \right)} e^{i \omega t}\]

The general solution to this equation reads

C1, C2 = sp.symbols('C1, C2',real=True)
W_sol = sp.dsolve(EQM_Fr,W(x)).rhs
display(W_sol)
W_sol = C1*sp.sin(omega*x/c)+C2*sp.cos(omega*x/c)
display(W_sol)

\[\displaystyle C_{1} \sin{\left(x \left|{\frac{\omega}{c}}\right| \right)} + C_{2} \cos{\left(\frac{\omega x}{c} \right)}\]

\[\displaystyle C_{1} \sin{\left(\frac{\omega x}{c} \right)} + C_{2} \cos{\left(\frac{\omega x}{c} \right)}\]

Now we will consider different boundary conditions

1: Fixed-fixed string: \(W(0)=0\), \(W(L)=L\)

L = sp.symbols('L',real=True)
BC_left =  W_sol.subs(x,0)
BC_right = W_sol.subs(x,L)
display(BC_left)
display(BC_right)

\[\displaystyle C_{2}\]

\[\displaystyle C_{1} \sin{\left(\frac{L \omega}{c} \right)} + C_{2} \cos{\left(\frac{L \omega}{c} \right)}\]

When \(\omega\) equals one of the natural frequencies, the determinant of the coefficient matrix composed of the two boundary conditions (i.e., the frequency equation) should be zero.

Let us compose the coefficient matrix

Coeff_Matrix= sp.Matrix([[BC_left.diff(C1),BC_left.diff(C2)],
                         [BC_right.diff(C1),BC_right.diff(C2)]])
display(Coeff_Matrix)

\[\begin{split}\displaystyle \left[\begin{matrix}0 & 1\\\sin{\left(\frac{L \omega}{c} \right)} & \cos{\left(\frac{L \omega}{c} \right)}\end{matrix}\right]\end{split}\]

Frequency_Equation = sp.det(Coeff_Matrix)
display(Frequency_Equation)

\[\displaystyle - \sin{\left(\frac{L \omega}{c} \right)}\]

Now we search for the roots of the frequency equation, which, by definition, are the natural frequencies (in the undamped case these frequencies coincide with the resonance frequencies)

Often, especially in the case of more complex boundary conditions, we introduce \(\Omega_{dl}= {\Omega \over L} C\). For this dimensionless parameter, the roots of the frequency equation can be always found numerically

Omega_dl = sp.symbols('Omega_dl',real=True)
Frequency_Equation_DL = Frequency_Equation.subs(omega,Omega_dl/L*c)
display(Frequency_Equation_DL)

\[\displaystyle - \sin{\left(\Omega_{dl} \right)}\]

Obviously, independently of the wave speed and the string length \(\Omega_{dl}=\pi n\), where \(n\) is an integer

sol = sp.solveset(sp.Eq(Frequency_Equation_DL,0),Omega_dl)
display(sol)

\[\displaystyle \left\{2 n \pi\; \middle|\; n \in \mathbb{Z}\right\} \cup \left\{2 n \pi + \pi\; \middle|\; n \in \mathbb{Z}\right\}\]

Leta us find values of \(\beta\) also numerically

Dimensionless_Natural_Frequency = []
step = 0.01
for i in range(5000):
    xi = i * step
    left_value = Frequency_Equation_DL.subs(Omega_dl,xi)
    right_value = Frequency_Equation_DL.subs(Omega_dl,xi+step)
    if left_value*right_value < 0:
        Dimensionless_Natural_Frequency.append(xi+step/2)
        
display(Dimensionless_Natural_Frequency)

[3.145,
 6.285,
 9.425,
 12.565000000000001,
 15.705000000000002,
 18.845,
 21.995,
 25.134999999999998,
 28.275,
 31.415,
 34.55500000000001,
 37.695,
 40.845000000000006,
 43.98500000000001,
 47.125]

2: Fixed-free string: \(W(0)=0, {\left. {{{dW\left( x \right)} \over {dx}}} \right|_{x = L}} = 0\)

BC_left =  W_sol.subs(x,0)
BC_right = W_sol.diff(x).subs(x,L)
display(BC_left)
display(BC_right)

\[\displaystyle C_{2}\]

\[\displaystyle \frac{C_{1} \omega \cos{\left(\frac{L \omega}{c} \right)}}{c} - \frac{C_{2} \omega \sin{\left(\frac{L \omega}{c} \right)}}{c}\]

Coeff_Matrix= sp.Matrix([[BC_left.diff(C1),BC_left.diff(C2)],
                         [BC_right.diff(C1),BC_right.diff(C2)]])
display(Coeff_Matrix)

\[\begin{split}\displaystyle \left[\begin{matrix}0 & 1\\\frac{\omega \cos{\left(\frac{L \omega}{c} \right)}}{c} & - \frac{\omega \sin{\left(\frac{L \omega}{c} \right)}}{c}\end{matrix}\right]\end{split}\]

Frequency_Equation = sp.det(Coeff_Matrix)
display(Frequency_Equation)

\[\displaystyle - \frac{\omega \cos{\left(\frac{L \omega}{c} \right)}}{c}\]

Frequency_Equation_DL = Frequency_Equation.subs(omega,Omega_dl/L*c)
display(Frequency_Equation_DL)

\[\displaystyle - \frac{\Omega_{dl} \cos{\left(\Omega_{dl} \right)}}{L}\]

sol = sp.solveset(sp.Eq(Frequency_Equation_DL,0),Omega_dl)
display(sol)

\[\displaystyle \left\{0\right\} \cup \left\{2 n \pi + \frac{\pi}{2}\; \middle|\; n \in \mathbb{Z}\right\} \cup \left\{2 n \pi + \frac{3 \pi}{2}\; \middle|\; n \in \mathbb{Z}\right\}\]

3: Fixed-{Elastically Supported Inertial End}:

• \({w\left( {0,t} \right) = 0}\)

• \({\left. {T{{\delta w\left( {x,t} \right)} \over {\delta x}}} \right|_{x = L}} = - m{\left. {{{{\delta ^2}w\left( {x,t} \right)} \over {\delta {x^2}}}} \right|_{x = L}} + kw\left( {x,t} \right)\)

T, m, k = sp.symbols('T, m, k',real=True)
BC_left =  W_sol.subs(x,0)
BC_right = T*sp.diff(w(x,t),x)+m*sp.diff(w(x,t),t,2)+k*w(x,t)
display(BC_left)
display(BC_right)

\[\displaystyle C_{2}\]

\[\displaystyle T \frac{\partial}{\partial x} w{\left(x,t \right)} + k w{\left(x,t \right)} + m \frac{\partial^{2}}{\partial t^{2}} w{\left(x,t \right)}\]

BC_right = sp.simplify(BC_right.subs(w(x,t),W_sol*sp.exp(sp.I*omega*t))).subs(x,L)
display(BC_right)

\[\displaystyle T \left(\frac{C_{1} \omega \cos{\left(\frac{L \omega}{c} \right)}}{c} - \frac{C_{2} \omega \sin{\left(\frac{L \omega}{c} \right)}}{c}\right) e^{i \omega t} + k \left(C_{1} \sin{\left(\frac{L \omega}{c} \right)} + C_{2} \cos{\left(\frac{L \omega}{c} \right)}\right) e^{i \omega t} - m \omega^{2} \left(C_{1} \sin{\left(\frac{L \omega}{c} \right)} + C_{2} \cos{\left(\frac{L \omega}{c} \right)}\right) e^{i \omega t}\]

\(e^{i \Omega t}\) is not dropped, so can dropped manually, not necessary:

BC_right = BC_right.subs(sp.exp(sp.I*omega*t),1)
display(BC_right)

\[\displaystyle T \left(\frac{C_{1} \omega \cos{\left(\frac{L \omega}{c} \right)}}{c} - \frac{C_{2} \omega \sin{\left(\frac{L \omega}{c} \right)}}{c}\right) + k \left(C_{1} \sin{\left(\frac{L \omega}{c} \right)} + C_{2} \cos{\left(\frac{L \omega}{c} \right)}\right) - m \omega^{2} \left(C_{1} \sin{\left(\frac{L \omega}{c} \right)} + C_{2} \cos{\left(\frac{L \omega}{c} \right)}\right)\]

Coeff_Matrix= sp.Matrix([[BC_left.diff(C1),BC_left.diff(C2)],
                         [BC_right.diff(C1),BC_right.diff(C2)]])
display(Coeff_Matrix)

\[\begin{split}\displaystyle \left[\begin{matrix}0 & 1\\\frac{T \omega \cos{\left(\frac{L \omega}{c} \right)}}{c} + k \sin{\left(\frac{L \omega}{c} \right)} - m \omega^{2} \sin{\left(\frac{L \omega}{c} \right)} & - \frac{T \omega \sin{\left(\frac{L \omega}{c} \right)}}{c} + k \cos{\left(\frac{L \omega}{c} \right)} - m \omega^{2} \cos{\left(\frac{L \omega}{c} \right)}\end{matrix}\right]\end{split}\]

Frequency_Equation = sp.det(Coeff_Matrix)
display(Frequency_Equation)

\[\displaystyle \frac{- T \omega \cos{\left(\frac{L \omega}{c} \right)} - c k \sin{\left(\frac{L \omega}{c} \right)} + c m \omega^{2} \sin{\left(\frac{L \omega}{c} \right)}}{c}\]

This equation should be solved numerically. The natural frequencies depend both on the string parameters and the mass and stiffness of the boundary

Frequency_Equation_subs = Frequency_Equation.subs([(L,1),(c,2),(T,10),(k,20),(m,7)])

sol = sp.solveset(sp.Eq(Frequency_Equation_subs,0),omega)
display(sol)

\[\displaystyle \left\{\omega\; \middle|\; \omega \in \mathbb{C} \wedge 7 \omega^{2} \sin{\left(\frac{\omega}{2} \right)} - 5 \omega \cos{\left(\frac{\omega}{2} \right)} - 20 \sin{\left(\frac{\omega}{2} \right)} = 0 \right\}\]

Dimensionless_Natural_Frequency = []
step = 0.01
for i in range(5000):
    xi = i * step
    left_value = Frequency_Equation_subs.subs(omega,xi)
    right_value = Frequency_Equation_subs.subs(omega,xi+step)
    if left_value*right_value < 0:
        Dimensionless_Natural_Frequency.append(xi+step/2)
        
display(Dimensionless_Natural_Frequency)

[1.9449999999999998,
 6.515,
 12.685,
 18.925,
 25.185,
 31.465,
 37.73500000000001,
 44.015]