Update modal extraction

This commit is contained in:
Thomas Dehaeze 2019-07-16 11:09:14 +02:00
parent c0641f7865
commit e4d8712c94
8 changed files with 119 additions and 63 deletions

Binary file not shown.

Before

Width:  |  Height:  |  Size: 136 KiB

After

Width:  |  Height:  |  Size: 114 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 52 KiB

Binary file not shown.

View File

@ -43,6 +43,29 @@ Theses matrices will be used to tune the Simscape (multi-body) model.
The modes we want to identify are those in the frequency range between 0Hz and 150Hz. The modes we want to identify are those in the frequency range between 0Hz and 150Hz.
#+name: tab:terminology_numbers
#+caption: Terminology for further analysis
| <c> | <l> | <c> |
| Symbol | Meaning | Value |
|--------+---------------------------------+-------|
| $p$ | Number of solid body considered | 6 |
| $m$ | Number of accelerometers | 23 |
| $n$ | Number of identified modes | 21 |
| $q$ | Number of frequency points | 801 |
| $s$ | Number of excitation | 3 |
#+name: tab:terminology_elements
#+caption: Terminology for further analysis
| <c> | <l> |
| Symbol | Meaning |
|-----------------+--------------------------------------|
| $[\Lambda]$ | Complex eigen value matrix |
| $[\Psi]$ | Complex eigen vector matrix |
| $\omega_r$ | Eigen frequency of mode $r$ [rad/s] |
| $\xi_r$ | Modal damping for mode $r$ |
| $\{\psi\}_r$ | Complex mode shape of mode $r$ |
| $[M], [C], [K]$ | Mass, damping and stiffness matrices |
| $a_r$ | "Modal A" for mode $r$ |
The modal analysis of the ID31 Micro-station thus consists of several parts: The modal analysis of the ID31 Micro-station thus consists of several parts:
- [[file:measurement.org][Frequency Response Measurements]] - [[file:measurement.org][Frequency Response Measurements]]

Binary file not shown.

View File

@ -41,6 +41,9 @@
#+PROPERTY: header-args:latex+ :output-dir figs #+PROPERTY: header-args:latex+ :output-dir figs
:END: :END:
* Introduction :ignore:
The goal is to measure the dynamic of the Micro-Station and to extract Frequency Response Functions.
* ZIP file containing the data and matlab files :ignore: * ZIP file containing the data and matlab files :ignore:
#+begin_src bash :exports none :results none #+begin_src bash :exports none :results none
if [ matlab/modal_frf_coh.m -nt data/modal_frf_coh.zip ]; then if [ matlab/modal_frf_coh.m -nt data/modal_frf_coh.zip ]; then
@ -67,9 +70,6 @@
<<matlab-init>> <<matlab-init>>
#+end_src #+end_src
* Goal
The goal is to measure the dynamic of the Micro-Station and to extract Frequency Response Functions.
* Instrumentation Used * Instrumentation Used
In order to perform to *Modal Analysis* and to obtain first a *Response Model*, the following devices are used: In order to perform to *Modal Analysis* and to obtain first a *Response Model*, the following devices are used:
- An *acquisition system* (OROS) with 24bits ADCs (figure [[fig:oros]]) - An *acquisition system* (OROS) with 24bits ADCs (figure [[fig:oros]])
@ -956,3 +956,7 @@ Below, it could be due to:
As it is usually very important to measure point response $\frac{X_j}{F_j}$, that is to say exciting and measuring the *same* DOF, should the measurements be redone? As it is usually very important to measure point response $\frac{X_j}{F_j}$, that is to say exciting and measuring the *same* DOF, should the measurements be redone?
On the modal software that is used to extract modal parameters, it is suppose that we are exciting the *exact* same DOF as the one measured by the accelerometer 11. As is it not the case in practice, could this induce large errors in the modal model? On the modal software that is used to extract modal parameters, it is suppose that we are exciting the *exact* same DOF as the one measured by the accelerometer 11. As is it not the case in practice, could this induce large errors in the modal model?
#+end_warning #+end_warning
* Bibliography :ignore:
bibliographystyle:unsrt
bibliography:ref.bib

Binary file not shown.

View File

@ -273,41 +273,41 @@ They are all exported in a text file named =modes.asc=.
Its first 20 lines as shown below. Its first 20 lines as shown below.
#+begin_src bash :results output :exports results :eval no-export #+begin_src bash :results output :exports results :eval no-export
sed 20q mat/modes.asc | sed $'s/\r//' sed 20q mat/modes_2019_07_16.asc | sed $'s/\r//'
#+end_src #+end_src
#+RESULTS: #+RESULTS:
#+begin_example #+begin_example
Created by N-Modal Created by N-Modal
Estimator: bbfd Estimator: bbfd
01-Jul-19 16:44:11 16-Jul-19 09:35:42
Mode 1 Mode 1
freq = 11.41275Hz freq = 11.86509Hz
damp = 8.72664% damp = 12.20318%
modal A = -4.50556e+003-9.41744e+003i modal A = 4.13559e+003+6.22828e+003i
modal B = -7.00928e+005+2.62922e+005i modal B = 4.98475e+005-2.49344e+005i
Mode matrix of local coordinate [DOF: Re IM] Mode matrix of local coordinate [DOF: Re IM]
1X-: -1.04114e-001 3.50664e-002 1X+: -1.08546e-001 3.17769e-002
1Y-: 2.34008e-001 5.04273e-004 1Y+: 2.32355e-001 1.71539e-003
1Z+: -1.93303e-002 5.08614e-003 1Z+: -1.98455e-002 5.79676e-003
2X-: -8.38439e-002 3.45978e-002 2X+: -8.59861e-002 2.84442e-002
2Y-: 2.42440e-001 0.00000e+000 2Y+: 2.40594e-001 0.00000e+000
2Z+: -7.40734e-003 5.17734e-003 2Z+: -7.65144e-003 4.22952e-003
3Y-: 2.17655e-001 6.10802e-003 3Y+: 2.14274e-001 7.53040e-003
3X+: 1.18685e-001 -3.54602e-002 3X-: 1.22352e-001 -2.90328e-002
3Z+: -2.37725e-002 -1.61649e-003 3Z+: -2.12103e-002 -3.81368e-003
#+end_example #+end_example
We split this big =modes.asc= file into sub text files using =bash=. The obtained files are described one table [[tab:modes_files]]. We split this big =modes.asc= file into sub text files using =bash=. The obtained files are described one table [[tab:modes_files]].
#+begin_src bash :results none #+begin_src bash :results none
sed '/^\s*[0-9]*[XYZ][+-]:/!d' mat/modes.asc > mat/mode_shapes.txt sed '/^\s*[0-9]*[XYZ][+-]:/!d' mat/modes_2019_07_16.asc > mat/mode_shapes.txt
sed '/freq/!d' mat/modes.asc | sed 's/.* = \(.*\)Hz/\1/' > mat/mode_freqs.txt sed '/freq/!d' mat/modes_2019_07_16.asc | sed 's/.* = \(.*\)Hz/\1/' > mat/mode_freqs.txt
sed '/damp/!d' mat/modes.asc | sed 's/.* = \(.*\)\%/\1/' > mat/mode_damps.txt sed '/damp/!d' mat/modes_2019_07_16.asc | sed 's/.* = \(.*\)\%/\1/' > mat/mode_damps.txt
sed '/modal A/!d' mat/modes.asc | sed 's/.* =\s\+\([-0-9.e]\++[0-9]\+\)\([-+0-9.e]\+\)i/\1 \2/' > mat/mode_modal_a.txt sed '/modal A/!d' mat/modes_2019_07_16.asc | sed 's/.* =\s\+\([-0-9.e]\++[0-9]\+\)\([-+0-9.e]\+\)i/\1 \2/' > mat/mode_modal_a.txt
sed '/modal B/!d' mat/modes.asc | sed 's/.* =\s\+\([-0-9.e]\++[0-9]\+\)\([-+0-9.e]\+\)i/\1 \2/' > mat/mode_modal_b.txt sed '/modal B/!d' mat/modes_2019_07_16.asc | sed 's/.* =\s\+\([-0-9.e]\++[0-9]\+\)\([-+0-9.e]\+\)i/\1 \2/' > mat/mode_modal_b.txt
#+end_src #+end_src
#+name: tab:modes_files #+name: tab:modes_files
@ -366,32 +366,29 @@ The obtained mode frequencies and damping are shown below.
data2orgtable([(1:length(freqs_m))', freqs_m, damps_m], {}, {'Mode number', 'Frequency [Hz]', 'Damping [%]'}, ' %.1f '); data2orgtable([(1:length(freqs_m))', freqs_m, damps_m], {}, {'Mode number', 'Frequency [Hz]', 'Damping [%]'}, ' %.1f ');
#+end_src #+end_src
#+name: tab:obtained_modes_freqs_damps
#+caption: Obtained eigen frequencies and modal damping
#+RESULTS: #+RESULTS:
| Mode number | Frequency [Hz] | Damping [%] | | Mode number | Frequency [Hz] | Damping [%] |
|-------------+----------------+-------------| |-------------+----------------+-------------|
| 1.0 | 11.4 | 8.7 | | 1.0 | 11.9 | 12.2 |
| 2.0 | 18.5 | 11.8 | | 2.0 | 18.6 | 11.7 |
| 3.0 | 37.6 | 6.4 | | 3.0 | 37.8 | 6.2 |
| 4.0 | 39.4 | 3.6 | | 4.0 | 39.1 | 2.8 |
| 5.0 | 54.0 | 0.2 | | 5.0 | 56.3 | 2.8 |
| 6.0 | 56.1 | 2.8 | | 6.0 | 69.8 | 4.3 |
| 7.0 | 69.7 | 4.6 | | 7.0 | 72.5 | 1.3 |
| 8.0 | 71.6 | 0.6 | | 8.0 | 84.8 | 3.7 |
| 9.0 | 72.4 | 1.6 | | 9.0 | 91.3 | 2.9 |
| 10.0 | 84.9 | 3.6 | | 10.0 | 105.5 | 3.2 |
| 11.0 | 90.6 | 0.3 | | 11.0 | 106.6 | 1.6 |
| 12.0 | 91.0 | 2.9 | | 12.0 | 112.7 | 3.1 |
| 13.0 | 95.8 | 3.3 | | 13.0 | 124.2 | 2.8 |
| 14.0 | 105.4 | 3.3 | | 14.0 | 145.3 | 1.3 |
| 15.0 | 106.8 | 1.9 | | 15.0 | 150.5 | 2.4 |
| 16.0 | 112.6 | 3.0 | | 16.0 | 165.4 | 1.4 |
| 17.0 | 116.8 | 2.7 |
| 18.0 | 124.1 | 0.6 | #+name: tab:obtained_modes_freqs_damps
| 19.0 | 145.4 | 1.6 | #+caption: Obtained eigen frequencies and modal damping
| 20.0 | 150.1 | 2.2 | #+RESULTS:
| 21.0 | 164.7 | 1.4 |
** Theory ** Theory
It seems that the modal analysis software makes the *assumption* of viscous damping for the model with which it tries to fit the FRF measurements. It seems that the modal analysis software makes the *assumption* of viscous damping for the model with which it tries to fit the FRF measurements.
@ -472,6 +469,9 @@ With $a_i$ is the "Modal A" parameter linked to mode i.
modal_b_ext_M = blkdiag(modal_b_M, conj(modal_b_M)); modal_b_ext_M = blkdiag(modal_b_M, conj(modal_b_M));
#+end_src #+end_src
"Modal A" and "modal B" are linked through the following formula:
\[ B = - A \Lambda \]
* Obtained Mode Shapes animations * Obtained Mode Shapes animations
<<sec:mode_shape_display>> <<sec:mode_shape_display>>
From the modal parameters, it is possible to show the modal shapes with an animation. From the modal parameters, it is possible to show the modal shapes with an animation.
@ -511,8 +511,8 @@ There are two main ways to verify the validity of the modal model
From the modal model, we want to synthesize the Frequency Response Functions that has been used to build the modal model. From the modal model, we want to synthesize the Frequency Response Functions that has been used to build the modal model.
Let's recall that: Let's recall that:
- $M$ is the number of measured DOFs ($23 \times 3 = 69$) - $M$ is the number of measured DOFs ($3 \times n_\text{acc}$)
- $N$ is the number of modes identified ($21$) - $N$ is the number of modes identified
We then have that the FRF matrix $[H_{\text{syn}}]$ can be synthesize using the following formula: We then have that the FRF matrix $[H_{\text{syn}}]$ can be synthesize using the following formula:
#+begin_important #+begin_important
@ -529,12 +529,34 @@ with:
- $\psi_{pr}$: scaled modal coefficient for output DOF $p$, mode $r$ - $\psi_{pr}$: scaled modal coefficient for output DOF $p$, mode $r$
- $\lambda_r$: complex modal frequency - $\lambda_r$: complex modal frequency
From the modal software documentation:
#+begin_quote
*Modal A*
Scaling constant for a complex mode. It has the same properties as modal mass for normal modes (undamped or proportionally damped cases). Assuming
- $\psi_{pr}$ = Modal coefficient for measured degree of freedom p and mode r
- $\psi_{qr}$ = Modal coefficient for measured degree of freedom q and mode r
- $A_{pqr}$ = Residue for measured degree of freedom p, measured degree of q and mode r
- $M_{Ar}$ = Modal A of mode r
Then
\[ A_{pqr} = \frac{\psi_{pr}\psi_{qr}}{M_{Ar}} \]
*Modal B*
Scaling constant for a complex mode. It has the same properties as modal stiffness for normal modes (undamped or proportionally damped cases). Assuming
- $M_{Ar}$ = Modal A of mode r
- $\lambda_r$ = System pole of mode r
Then
\[ M_{Br} = - \lambda_r M_{Ar} \]
#+end_quote
** Matlab Implementation ** Matlab Implementation
#+begin_src matlab #+begin_src matlab
Hsyn = zeros(69, 69, 801); Hsyn = zeros(acc_n*dir_n, acc_n*dir_n, length(freqs));
for i = 1:length(freqs) for i = 1:length(freqs)
Hsyn(:, :, i) = eigen_vec_ext_M*(inv(modal_a_ext_M)/(diag(j*2*pi*freqs(i) - diag(eigen_val_ext_M))))*eigen_vec_ext_M.'; Hsyn(:, :, i) = eigen_vec_ext_M*((j*2*pi*freqs(i)).^2*inv(modal_a_ext_M)/(diag(j*2*pi*freqs(i) - diag(eigen_val_ext_M))))*eigen_vec_ext_M.';
end end
#+end_src #+end_src
@ -542,13 +564,13 @@ Because the synthesize frequency response functions are representing the displac
#+begin_src matlab #+begin_src matlab
for i = 1:size(Hsyn, 1) for i = 1:size(Hsyn, 1)
Hsyn(i, :, :) = squeeze(-Hsyn(i, :, :)).*(j*2*pi*freqs).^2; Hsyn(i, :, :) = squeeze(Hsyn(i, :, :)).*(j*2*pi*freqs).^2;
end end
#+end_src #+end_src
** Original and Synthesize FRF matrix comparison ** Original and Synthesize FRF matrix comparison
#+begin_src matlab :exports none #+begin_src matlab :exports none
acc_o = 15; dir_o = 1; dir_i = 1; acc_o = 1; dir_o = 1; dir_i = 1;
figure; figure;
ax1 = subplot(2, 1, 1); ax1 = subplot(2, 1, 1);
@ -556,9 +578,10 @@ Because the synthesize frequency response functions are representing the displac
plot(freqs, abs(squeeze(FRFs(3*(acc_o-1)+dir_o, dir_i, :))), 'DisplayName', 'Original'); plot(freqs, abs(squeeze(FRFs(3*(acc_o-1)+dir_o, dir_i, :))), 'DisplayName', 'Original');
plot(freqs, abs(squeeze(Hsyn(3*(acc_o-1)+dir_o, 3*(11-1)+dir_i, :))), 'DisplayName', 'Synthesize'); plot(freqs, abs(squeeze(Hsyn(3*(acc_o-1)+dir_o, 3*(11-1)+dir_i, :))), 'DisplayName', 'Synthesize');
hold off; hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log'); set(gca, 'yscale', 'log');
set(gca, 'XTickLabel',[]); set(gca, 'XTickLabel',[]);
ylabel('Magnitude [$\frac{m/s^2}{N}$]'); ylabel('Magnitude [$\frac{m/s^2}{N}$]');
title(sprintf('From acc %i %s to acc %i %s', 11, dirs(dir_i), acc_o, dirs(dir_o)))
legend('location', 'northwest'); legend('location', 'northwest');
ax2 = subplot(2, 1, 2); ax2 = subplot(2, 1, 2);
@ -567,7 +590,6 @@ Because the synthesize frequency response functions are representing the displac
plot(freqs, mod(180/pi*phase(squeeze(Hsyn(3*(acc_o-1)+dir_o, 3*(11-1)+dir_i, :)))+180, 360)-180); plot(freqs, mod(180/pi*phase(squeeze(Hsyn(3*(acc_o-1)+dir_o, 3*(11-1)+dir_i, :)))+180, 360)-180);
hold off; hold off;
yticks(-360:90:360); ylim([-180, 180]); yticks(-360:90:360); ylim([-180, 180]);
set(gca, 'xscale', 'log');
xlabel('Frequency [Hz]'); ylabel('Phase [deg]'); xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
linkaxes([ax1,ax2],'x'); linkaxes([ax1,ax2],'x');
@ -583,10 +605,17 @@ Because the synthesize frequency response functions are representing the displac
#+CAPTION: Comparison of the Original and Synthesize FRF matrix element #+CAPTION: Comparison of the Original and Synthesize FRF matrix element
[[file:figs/compare_synthesize_original_frf.png]] [[file:figs/compare_synthesize_original_frf.png]]
#+name: fig:from11xto1x
#+caption: Original FRF and synthesize FRF using the modal software. From force applied in the X direction to the acceleration of accelerometer 11 in the X direction
[[file:img/modal_software/from11xto1x.jpg]]
#+begin_warning
The synthesize FRFs from the modal software do not match the synthesize FRFs computed here.
It is possible that it uses residues that are not exported? Nothing is said about that in the documentation.
#+end_warning
** Synthesize FRF that has not yet been measured ** Synthesize FRF that has not yet been measured
#+begin_src matlab :exports none #+begin_src matlab :exports none
dir_names = {'X', 'Y', 'Z'};
accs = [1]; dirs = [1:3]; accs = [1]; dirs = [1:3];
figure; figure;
@ -594,7 +623,7 @@ Because the synthesize frequency response functions are representing the displac
hold on; hold on;
for acc_i = accs for acc_i = accs
for dir_i = dirs for dir_i = dirs
plot(freqs, abs((1./(j*2*pi*freqs').^2).*squeeze(Hsyn(3*(acc_i-1)+dir_i, 3*(acc_i-1)+dir_i, :))), 'DisplayName', sprintf('Acc %i - %s', acc_i, dir_names{dir_i})); plot(freqs, abs((1./(j*2*pi*freqs').^2).*squeeze(Hsyn(3*(acc_i-1)+dir_i, 3*(acc_i-1)+dir_i, :))), 'DisplayName', sprintf('Acc %i - %s', acc_i, dirs(dir_i)));
end end
end end
hold off; hold off;
@ -846,9 +875,9 @@ This inversion is equivalent to a mean square problem.
** Matlab Implementation ** Matlab Implementation
The obtained mode shapes matrix that gives the mode shapes of each solid bodies with respect to the fixed frame $\{O\}$, =mode_shapes_O=, is an $n \times p \times q$ with: The obtained mode shapes matrix that gives the mode shapes of each solid bodies with respect to the fixed frame $\{O\}$, =mode_shapes_O=, is an $n \times p \times q$ with:
- $n$ is the number of modes: 21 - $n$ is the number of identified modes
- $p$ is the number of DOFs for each solid body: 6 - $p$ is the number of DOFs for each solid body (6)
- $q$ is the number of solid bodies: 6 - $q$ is the number of solid bodies
#+begin_src matlab #+begin_src matlab
mode_shapes_O = zeros(mod_n, 6, length(solid_names)); mode_shapes_O = zeros(mod_n, 6, length(solid_names));
@ -885,7 +914,7 @@ with
\psi_{1, x} & \psi_{1, y} & \psi_{1, z} & \psi_{1, \theta_x} & \psi_{1, \theta_y} & \psi_{1, \theta_z} & \psi_{2, x} & \dots & \psi_{6, \theta_z} \psi_{1, x} & \psi_{1, y} & \psi_{1, z} & \psi_{1, \theta_x} & \psi_{1, \theta_y} & \psi_{1, \theta_z} & \psi_{2, x} & \dots & \psi_{6, \theta_z}
\end{Bmatrix}^T \] \end{Bmatrix}^T \]
With $M = 6 \times 6$ is the new number of DOFs and $N=21$ is the number of modes. With $M = 6 \times n_\text{solid}$ is the new number of DOFs and $N$ is the number of modes.
#+begin_src matlab #+begin_src matlab
eigen_vec_O = reshape(mode_shapes_O, [mod_n, 6*length(solid_names)]).'; eigen_vec_O = reshape(mode_shapes_O, [mod_n, 6*length(solid_names)]).';
@ -1013,7 +1042,7 @@ Because the synthesize frequency response functions are representing the displac
We want to have the two eigen matrices. We want to have the two eigen matrices.
They should have the same size $n \times n$ where $n$ is the number of mode_shapes as well as the number of degrees of freedom. They should have the same size $n \times n$ where $n$ is the number of mode_shapes as well as the number of degrees of freedom.
Thus, if we consider 21 modes, we should restrict our system to have only 21 DOFs. Thus, if we consider $N$ modes, we should restrict our system to have only $N$ DOFs.
Actually, we are measured 6 DOFs of 6 solids, thus we have 36 DOFs. Actually, we are measured 6 DOFs of 6 solids, thus we have 36 DOFs.
@ -1140,7 +1169,7 @@ with
\psi_{1, x} & \psi_{1, y} & \psi_{1, z} & \psi_{1, \theta_x} & \psi_{1, \theta_y} & \psi_{1, \theta_z} & \psi_{2, x} & \dots & \psi_{6, \theta_z} \psi_{1, x} & \psi_{1, y} & \psi_{1, z} & \psi_{1, \theta_x} & \psi_{1, \theta_y} & \psi_{1, \theta_z} & \psi_{2, x} & \dots & \psi_{6, \theta_z}
\end{Bmatrix}^T \] \end{Bmatrix}^T \]
With $M = 6 \times 6$ is the new number of DOFs and $N=21$ is the number of modes. With $M = 6 \times n_\text{solids}$ is the new number of DOFs and $N$ is the number of modes.
Each eigen vector is normalized. Each eigen vector is normalized.