Skip to content
Snippets Groups Projects
Commit 3257f78d authored by Jenni Rinker's avatar Jenni Rinker
Browse files

adding eigenanalysis

parent 07e29724
No related branches found
No related tags found
No related merge requests found
%% Cell type:markdown id:fa60cf2d tags:
# Eigenanalysis
When analysing a dynamical system, it is always beneficial to know the natural frequencies of the system. These are frequencies at which the system "likes" to vibrate, and an oscillatory force near this frequency will have result in a larger system response than at other frequencies.
For a dynamical system with a set of mass, stiffness, and damping matrices, the undamped natural frequencies are calculated from the following generalized eigenvalue problem:
\begin{equation}
\lambda[M]\boldsymbol{v} = [K]\boldsymbol{v}.
\end{equation}
If $[M]$ is invertible, this reduces to the standard eigenvalue problem,
\begin{equation}
\lambda\boldsymbol{v} = [A]\boldsymbol{v},
\end{equation}
where $[A]=[M]^{-1}[K]$. The eigenvalue problem can be solved using the `eig` (Matlab) and `np.linalg.eig` (Python) functions. The natural frequency in radians/sec is directly related to the eigenvalue $\lambda$:
\begin{equation}
\lambda = \omega^2.
\end{equation}
%% Cell type:markdown id:01046bff tags:
### Exercises for the reader
1. Create a function that returns Turbie's mass, stiffness and damping matrices. It is okay to hard-code Turbie's parameters (available in the Definition tab) in the function, or you can try to load them from `turbie_parameters.txt`.
2. Create a function that takes the mass and stiffness matrices as input and returns an array of the undamped natural frequencies **in Hz**. Be careful with units!
3. **Optional**. Reexamine the PSD plot you made considering the natural frequencies you calculated. Do you see anything of interest?
### Answers to the exercises
#### Exercise 1
%% Cell type:code id:ff463f58 tags:
``` python
import numpy as np
def get_turbie_system_matrices():
"""Assmble M, C and K for Turbie"""
turbie_path = 'turbie_parameters.txt'
mb, mn, mh, mt, c1, c2, k1, k2, fb, ft, drb, drt, Dr, rho = np.loadtxt(turbie_path, comments='%')
m1 = 3*mb # mass 1 is 3 blades
m2 = mh + mn + mt # mass 2 is hub, nacelle and tower
M = np.array([[m1, 0], [0, m2]]) # mass matrix
C = np.array([[c1, -c1], [-c1, c1+c2]]) # damping matrix
K = np.array([[k1, -k1], [-k1, k1+k2]]) # stiffness matrix
return M, C, K
M, C, K = get_turbie_system_matrices()
print('Mass matrix: \n', M)
print('Damping matrix: \n', C)
print('Stiffness matrix: \n', K)
```
%% Output
Mass matrix:
[[ 123000. 0.]
[ 0. 1179000.]]
Damping matrix:
[[ 4208. -4208.]
[-4208. 16938.]]
Stiffness matrix:
[[ 1711000. -1711000.]
[-1711000. 4989000.]]
%% Cell type:markdown id:2d4d9bc5 tags:
#### Exercise 2
%% Cell type:code id:7e11efc0 tags:
``` python
def get_undamped_natural_frequencies(M, K, in_hz=True):
"""Calculate the undamped natural frequencies of a system in Hz.
Arguments: mass, stiffness and damping matrices.
"""
lamda = np.linalg.eig(np.linalg.inv(M) @ K)[0] # solve eiegenvalue problem (only evals)
omega = np.sort(np.sqrt(lamda)) # nat freqs in rad/s (ascending order)
freqs = omega / 2 / np.pi # in Hz
if in_hz:
return freqs
else:
return omega
freqs = get_undamped_natural_frequencies(M, K)
print('Undamped natural frequencies (Hz): ', freqs)
```
%% Output
Undamped natural frequencies (Hz): [0.25000017 0.63011524]
%% Cell type:markdown id:6e7f244a tags:
#### Exercise 3
Below is a reprint of the PSD plot from the PSD exercise. We can see peaks at 0.25 Hz and at 0.63 Hz, which matches our natural frequencies! This is a nice demonstration of how the energy in the input signal, given by the blue line, is redistributed to new frequencies in the response.
<img src="figures/3-psd.png" alt="Power spectral density" style="width: 600px;"/>
%% Cell type:markdown id:d19df525 tags:
......@@ -18,4 +18,5 @@ Contents:
1-definition
2-loading-plotting
3-psd
\ No newline at end of file
3-psd
4-eigenanalysis
\ No newline at end of file
......@@ -5,6 +5,19 @@ import matplotlib.pyplot as plt
import numpy as np
def get_turbie_system_matrices(turbie_path=None):
"""Assmble M, C and K for Turbie"""
if turbie_path is None:
turbie_path = 'turbie_parameters.txt'
mb, mn, mh, mt, c1, c2, k1, k2, fb, ft, drb, drt, Dr, rho = np.loadtxt(turbie_path, comments='%')
m1 = 3*mb # mass 1 is 3 blades
m2 = mh + mn + mt # mass 2 is hub, nacelle and tower
M = np.array([[m1, 0], [0, m2]]) # mass matrix
C = np.array([[c1, -c1], [-c1, c1+c2]]) # damping matrix
K = np.array([[k1, -k1], [-k1, k1+k2]]) # stiffness matrix
return M, C, K
def load_results(path):
"""Load turbie results from text file.
Returns 4 [nt] arrays: t, u, x1 and x2.
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment