Thomas Dehaeze
2caf80a292
- accelerometers - shaker Also added top accelereomters and device under test
729 lines
28 KiB
Org Mode
729 lines
28 KiB
Org Mode
#+TITLE: Vibration Table
|
|
:DRAWER:
|
|
#+LANGUAGE: en
|
|
#+EMAIL: dehaeze.thomas@gmail.com
|
|
#+AUTHOR: Dehaeze Thomas
|
|
|
|
#+HTML_LINK_HOME: ../index.html
|
|
#+HTML_LINK_UP: ../index.html
|
|
|
|
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="https://research.tdehaeze.xyz/css/style.css"/>
|
|
#+HTML_HEAD: <script type="text/javascript" src="https://research.tdehaeze.xyz/js/script.js"></script>
|
|
|
|
#+BIND: org-latex-image-default-option "scale=1"
|
|
#+BIND: org-latex-image-default-width ""
|
|
#+BIND: org-latex-bib-compiler "biber"
|
|
|
|
#+LaTeX_CLASS: scrreprt
|
|
#+LaTeX_CLASS_OPTIONS: [a4paper, 10pt, DIV=12, parskip=full]
|
|
#+LaTeX_HEADER_EXTRA: \input{preamble.tex}
|
|
#+LATEX_HEADER_EXTRA: \addbibresource{ref.bib}
|
|
|
|
#+PROPERTY: header-args:matlab :session *MATLAB*
|
|
#+PROPERTY: header-args:matlab+ :comments org
|
|
#+PROPERTY: header-args:matlab+ :exports both
|
|
#+PROPERTY: header-args:matlab+ :results none
|
|
#+PROPERTY: header-args:matlab+ :eval no-export
|
|
#+PROPERTY: header-args:matlab+ :noweb yes
|
|
#+PROPERTY: header-args:matlab+ :tangle no
|
|
#+PROPERTY: header-args:matlab+ :mkdirp yes
|
|
#+PROPERTY: header-args:matlab+ :output-dir figs
|
|
|
|
#+PROPERTY: header-args:latex :headers '("\\usepackage{tikz}" "\\usepackage{import}" "\\import{$HOME/Cloud/tikz/org/}{config.tex}")
|
|
#+PROPERTY: header-args:latex+ :imagemagick t :fit yes
|
|
#+PROPERTY: header-args:latex+ :iminoptions -scale 100% -density 150
|
|
#+PROPERTY: header-args:latex+ :imoutoptions -quality 100
|
|
#+PROPERTY: header-args:latex+ :results file raw replace
|
|
#+PROPERTY: header-args:latex+ :buffer no
|
|
#+PROPERTY: header-args:latex+ :tangle no
|
|
#+PROPERTY: header-args:latex+ :eval no-export
|
|
#+PROPERTY: header-args:latex+ :exports results
|
|
#+PROPERTY: header-args:latex+ :mkdirp yes
|
|
#+PROPERTY: header-args:latex+ :output-dir figs
|
|
#+PROPERTY: header-args:latex+ :post pdf2svg(file=*this*, ext="png")
|
|
:END:
|
|
|
|
#+begin_export html
|
|
<hr>
|
|
<p>This report is also available as a <a href="./vibration-table.pdf">pdf</a>.</p>
|
|
<hr>
|
|
#+end_export
|
|
|
|
* Introduction
|
|
|
|
This document is divided as follows:
|
|
- Section [[sec:experimental_setup]]: the experimental setup and all the instrumentation are described
|
|
- Section [[sec:meas_transformation]]: the mathematics used to compute the 6DoF motion of a solid body from several inertial sensor is derived
|
|
- Section [[sec:simscape_model]]: a Simscape model of the vibration table is developed
|
|
- Section [[sec:table_dynamics]]: the table dynamics is identified and compared with the Simscape model
|
|
|
|
* Experimental Setup
|
|
<<sec:experimental_setup>>
|
|
** Introduction :ignore:
|
|
** CAD Model
|
|
|
|
#+name: fig:vibration-table-cad-view
|
|
#+caption: CAD View of the vibration table
|
|
#+attr_latex: :width 0.8\linewidth
|
|
[[file:figs/vibration-table-cad-view.png]]
|
|
|
|
** Instrumentation
|
|
|
|
#+begin_note
|
|
Here are the documentation of the equipment used for this vibration table:
|
|
- Modal Shaker: Watson and Gearing
|
|
- Inertial Shaker: [[file:doc/inertial_shakers.pdf][IS20]]
|
|
- Viscoelastic supports: [[file:doc/810002_doc.pdf][810002]]
|
|
- Spring supports: [[file:doc/9129fcb6ec46bb52925bb16155a850f3be01c479.pdf][MV803-12CC]]
|
|
- Optical Table: [[https://www.thorlabs.com/thorproduct.cfm?partnumber=B4545A][B4545A]]
|
|
- Triaxial Accelerometer: [[https://www.pcb.com/products?model=356b18][356B18]]
|
|
- OROS
|
|
#+end_note
|
|
|
|
** Suspended table
|
|
|
|
- Dimensions :: 450 mm x 450 mm x 60 mm
|
|
- Mass :: 21.30 kg
|
|
|
|
#+name: fig:compliance_optical_table
|
|
#+caption: Compliance of the B4545A optical table
|
|
#+attr_latex: :width 0.8\linewidth
|
|
[[file:figs/B4545A_Compliance_inLb-780.png]]
|
|
|
|
** Inertial Sensors
|
|
|
|
| Equipment |
|
|
|----------------------------------|
|
|
| (2x) 1D accelerometer [[https://www.pcbpiezotronics.fr/produit/accelerometre/393b05/][PCB 393B05]] |
|
|
| (4x) 3D accelerometer [[https://www.pcbpiezotronics.fr/produit/accelerometres/356b18/][PCB 356B18]] |
|
|
|
|
* Compute the 6DoF solid body motion from several inertial sensors
|
|
:PROPERTIES:
|
|
:header-args:matlab+: :tangle matlab/meas_transformation.m
|
|
:END:
|
|
<<sec:meas_transformation>>
|
|
** Introduction :ignore:
|
|
Let's consider a solid body with several accelerometers attached to it (Figure [[fig:local_to_global_coordinates]]).
|
|
|
|
#+begin_src latex :file local_to_global_coordinates.pdf :post pdf2svg(file=*this*, ext="png") :exports results
|
|
\newcommand\irregularcircle[2]{% radius, irregularity
|
|
\pgfextra {\pgfmathsetmacro\len{(#1)+rand*(#2)}}
|
|
+(0:\len pt)
|
|
\foreach \a in {10,20,...,350}{
|
|
\pgfextra {\pgfmathsetmacro\len{(#1)+rand*(#2)}}
|
|
-- +(\a:\len pt)
|
|
} -- cycle
|
|
}
|
|
\begin{tikzpicture}
|
|
\draw[rounded corners=1mm, fill=blue!30!white] (0, 0) \irregularcircle{3cm}{1mm};
|
|
|
|
\node[] (origin) at (0, 0) {$\bullet$};
|
|
\begin{scope}[shift={(origin)}]
|
|
\def\axissize{0.8cm}
|
|
\draw[->] (0, 0) -- ++(\axissize, 0) node[above left]{$x$};
|
|
\draw[->] (0, 0) -- ++(0, \axissize) node[below right]{$y$};
|
|
\draw[fill, color=black] (0, 0) circle (0.05*\axissize);
|
|
\node[draw, circle, inner sep=0pt, minimum size=0.4*\axissize, label=left:$z$] (yaxis) at (0, 0){};
|
|
\node[below right] at (0, 0){$\{O\}$};
|
|
\end{scope}
|
|
|
|
\coordinate[] (p1) at (-1.5, -1.5);
|
|
\coordinate[] (p2) at (-1.5, -1.5);
|
|
\coordinate[] (p3) at ( 1.5, -1.5);
|
|
\coordinate[] (p4) at ( 1.5, -1.5);
|
|
\coordinate[] (p5) at ( 1.5, 1.5);
|
|
\coordinate[] (p6) at ( 1.5, 1.5);
|
|
|
|
\draw[->] (p1)node[]{$\bullet$} -- ++(0, 1)node[right]{$a_1$};
|
|
\node[draw, circle, inner sep=0pt, minimum size=0.3cm, label=left:$a_2$] at (p2){};
|
|
\draw[fill, color=black] (p2) circle (0.05);
|
|
\draw[->] (p3)node[]{$\bullet$} -- ++(1, 0)node[right]{$a_3$};
|
|
\node[draw, circle, inner sep=0pt, minimum size=0.3cm, label=left:$a_4$] at (p4){};
|
|
\draw[fill, color=black] (p4) circle (0.05);
|
|
\draw[->] (p5)node[]{$\bullet$} -- ++(1, 0)node[right]{$a_5$};
|
|
\node[draw, circle, inner sep=0pt, minimum size=0.3cm, label=left:$a_6$] at (p6){};
|
|
\draw[fill, color=black] (p6) circle (0.04);
|
|
\end{tikzpicture}
|
|
#+end_src
|
|
|
|
#+name: fig:local_to_global_coordinates
|
|
#+caption: Schematic of the measured motions of a solid body
|
|
#+RESULTS:
|
|
[[file:figs/local_to_global_coordinates.png]]
|
|
|
|
The goal of this section is to see how to compute the acceleration/angular acceleration of the solid body from the accelerations $\vec{a}_i$ measured by the accelerometers.
|
|
|
|
The acceleration/angular acceleration of the solid body is defined as the vector ${}^O\vec{x}$:
|
|
\begin{equation}
|
|
{}^O\vec{x} = \begin{bmatrix}
|
|
\dot{v}_x \\ \dot{v}_y \\ \dot{v}_z \\ \dot{\omega}_x \\ \dot{\omega}_y \\ \dot{\omega}_z
|
|
\end{bmatrix}
|
|
\end{equation}
|
|
|
|
As we want to measure 6dof, we suppose that we have 6 uniaxial acceleremoters (we could use more, but 6 is enough).
|
|
The measurement of the individual vectors is defined as the vector $\vec{a}$:
|
|
\begin{equation}
|
|
\vec{a} = \begin{bmatrix}
|
|
a_1 \\ a_2 \\ a_3 \\ a_4 \\ a_5 \\ a_6
|
|
\end{bmatrix}
|
|
\end{equation}
|
|
|
|
From the positions and orientations of the acceleremoters (defined in Section [[sec:accelerometer_pos]]), it is quite straightforward to compute the accelerations measured by the sensors from the acceleration/angular acceleration of the solid body (Section [[sec:transformation_motion_to_acc]]).
|
|
From this, we can easily build a transformation matrix $M$, such that:
|
|
\begin{equation}
|
|
\vec{a} = M \cdot {}^O\vec{x}
|
|
\end{equation}
|
|
|
|
If the matrix is invertible, we can just take the inverse in order to obtain the transformation matrix giving the 6dof acceleration of the solid body from the accelerometer measurements (Section [[sec:transformation_acc_to_motion]]):
|
|
\begin{equation}
|
|
{}^O\vec{x} = M^{-1} \cdot \vec{a}
|
|
\end{equation}
|
|
|
|
If it is not invertible, then it means that it is not possible to compute all 6dof of the solid body from the measurements.
|
|
The solution is then to change the location/orientation of some of the accelerometers.
|
|
|
|
** Matlab Init :noexport:ignore:
|
|
#+begin_src matlab :tangle no :exports none :results silent :noweb yes :var current_dir=(file-name-directory buffer-file-name)
|
|
<<matlab-dir>>
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none :results silent :noweb yes
|
|
<<matlab-init>>
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no
|
|
addpath('matlab/')
|
|
#+end_src
|
|
|
|
** Define accelerometers positions/orientations
|
|
<<sec:accelerometer_pos>>
|
|
Let's first define the position and orientation of all measured accelerations with respect to a defined frame $\{O\}$.
|
|
|
|
#+begin_src matlab
|
|
Opm = [-0.1875, -0.1875, -0.245;
|
|
-0.1875, -0.1875, -0.245;
|
|
0.1875, -0.1875, -0.245;
|
|
0.1875, -0.1875, -0.245;
|
|
0.1875, 0.1875, -0.245;
|
|
0.1875, 0.1875, -0.245]';
|
|
#+end_src
|
|
|
|
There are summarized in Table [[tab:accelerometers_table_positions]].
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
|
|
data2orgtable(Opm, {'x', 'y', 'z'}, {'$a_1$', '$a_2$', '$a_3$', '$a_4$', '$a_5$', '$a_6$'}, ' %.3f ');
|
|
#+end_src
|
|
|
|
#+name: tab:accelerometers_table_positions
|
|
#+caption: Positions of the accelerometers fixed to the vibration table with respect to $\{O\}$
|
|
#+attr_latex: :environment tabularx :width 0.55\linewidth :align Xcccccc
|
|
#+attr_latex: :center t :booktabs t :float t
|
|
#+RESULTS:
|
|
| | $a_1$ | $a_2$ | $a_3$ | $a_4$ | $a_5$ | $a_6$ |
|
|
|---+--------+--------+--------+--------+--------+--------|
|
|
| x | -0.188 | -0.188 | 0.188 | 0.188 | 0.188 | 0.188 |
|
|
| y | -0.188 | -0.188 | -0.188 | -0.188 | 0.188 | 0.188 |
|
|
| z | -0.245 | -0.245 | -0.245 | -0.245 | -0.245 | -0.245 |
|
|
|
|
We then define the direction of the measured accelerations (unit vectors):
|
|
#+begin_src matlab
|
|
Osm = [0, 1, 0;
|
|
0, 0, 1;
|
|
1, 0, 0;
|
|
0, 0, 1;
|
|
1, 0, 0;
|
|
0, 0, 1;]';
|
|
#+end_src
|
|
|
|
They are summarized in Table [[tab:accelerometers_table_orientations]].
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
|
|
data2orgtable(Osm, {'x', 'y', 'z'}, {'$\hat{s}_1$', '$\hat{s}_2$', '$\hat{s}_3$', '$\hat{s}_4$', '$\hat{s}_5$', '$\hat{s}_6$'}, ' %.0f ');
|
|
#+end_src
|
|
|
|
#+name: tab:accelerometers_table_orientations
|
|
#+caption: Orientations of the accelerometers fixed to the vibration table expressed in $\{O\}$
|
|
#+attr_latex: :environment tabularx :width 0.35\linewidth :align Xcccccc
|
|
#+attr_latex: :center t :booktabs t :float t
|
|
#+RESULTS:
|
|
| | $\hat{s}_1$ | $\hat{s}_2$ | $\hat{s}_3$ | $\hat{s}_4$ | $\hat{s}_5$ | $\hat{s}_6$ |
|
|
|---+-------------+-------------+-------------+-------------+-------------+-------------|
|
|
| x | 0 | 0 | 1 | 0 | 1 | 0 |
|
|
| y | 1 | 0 | 0 | 0 | 0 | 0 |
|
|
| z | 0 | 1 | 0 | 1 | 0 | 1 |
|
|
|
|
** Transformation matrix from motion of the solid body to accelerometer measurements
|
|
<<sec:transformation_motion_to_acc>>
|
|
|
|
Let's try to estimate the x-y-z acceleration of any point of the solid body from the acceleration/angular acceleration of the solid body expressed in $\{O\}$.
|
|
For any point $p_i$ of the solid body (corresponding to an accelerometer), we can write:
|
|
\begin{equation}
|
|
\begin{bmatrix}
|
|
a_{i,x} \\ a_{i,y} \\ a_{i,z}
|
|
\end{bmatrix} = \begin{bmatrix}
|
|
\dot{v}_x \\ \dot{v}_y \\ \dot{v}_z
|
|
\end{bmatrix} + p_i \times \begin{bmatrix}
|
|
\dot{\omega}_x \\ \dot{\omega}_y \\ \dot{\omega}_z
|
|
\end{bmatrix}
|
|
\end{equation}
|
|
|
|
We can write the cross product as a matrix product using the skew-symmetric transformation:
|
|
\begin{equation}
|
|
\begin{bmatrix}
|
|
a_{i,x} \\ a_{i,y} \\ a_{i,z}
|
|
\end{bmatrix} = \begin{bmatrix}
|
|
\dot{v}_x \\ \dot{v}_y \\ \dot{v}_z
|
|
\end{bmatrix} + \underbrace{\begin{bmatrix}
|
|
0 & p_{i,z} & -p_{i,y} \\
|
|
-p_{i,z} & 0 & p_{i,x} \\
|
|
p_{i,y} & -p_{i,x} & 0
|
|
\end{bmatrix}}_{P_{i,[\times]}} \cdot \begin{bmatrix}
|
|
\dot{\omega}_x \\ \dot{\omega}_y \\ \dot{\omega}_z
|
|
\end{bmatrix}
|
|
\end{equation}
|
|
|
|
If we now want to know the (scalar) acceleration $a_i$ of the point $p_i$ in the direction of the accelerometer direction $\hat{s}_i$, we can just project the 3d acceleration on $\hat{s}_i$:
|
|
\begin{equation}
|
|
a_i = \hat{s}_i^T \cdot \begin{bmatrix}
|
|
a_{i,x} \\ a_{i,y} \\ a_{i,z}
|
|
\end{bmatrix} = \hat{s}_i^T \cdot \begin{bmatrix}
|
|
\dot{v}_x \\ \dot{v}_y \\ \dot{v}_z
|
|
\end{bmatrix} + \left( \hat{s}_i^T \cdot P_{i,[\times]} \right) \cdot \begin{bmatrix}
|
|
\dot{\omega}_x \\ \dot{\omega}_y \\ \dot{\omega}_z
|
|
\end{bmatrix}
|
|
\end{equation}
|
|
|
|
Which is equivalent as a simple vector multiplication:
|
|
\begin{equation}
|
|
a_i = \begin{bmatrix}
|
|
\hat{s}_i^T & \hat{s}_i^T \cdot P_{i,[\times]}
|
|
\end{bmatrix}
|
|
\begin{bmatrix}
|
|
\dot{v}_x \\ \dot{v}_y \\ \dot{v}_z \\ \dot{\omega}_x \\ \dot{\omega}_y \\ \dot{\omega}_z
|
|
\end{bmatrix} = \begin{bmatrix}
|
|
\hat{s}_i^T & \hat{s}_i^T \cdot P_{i,[\times]}
|
|
\end{bmatrix} {}^O\vec{x}
|
|
\end{equation}
|
|
|
|
And finally we can combine the 6 (line) vectors for the 6 accelerometers to write that in a matrix form.
|
|
We obtain Eq. eqref:eq:M_matrix.
|
|
#+begin_important
|
|
The transformation from solid body acceleration ${}^O\vec{x}$ from sensor measured acceleration $\vec{a}$ is:
|
|
\begin{equation} \label{eq:M_matrix}
|
|
\vec{a} = \underbrace{\begin{bmatrix}
|
|
\hat{s}_1^T & \hat{s}_1^T \cdot P_{1,[\times]} \\
|
|
\vdots & \vdots \\
|
|
\hat{s}_6^T & \hat{s}_6^T \cdot P_{6,[\times]}
|
|
\end{bmatrix}}_{M} {}^O\vec{x}
|
|
\end{equation}
|
|
|
|
with $\hat{s}_i$ the unit vector representing the measured direction of the i'th accelerometer expressed in frame $\{O\}$ and $P_{i,[\times]}$ the skew-symmetric matrix representing the cross product of the position of the i'th accelerometer expressed in frame $\{O\}$.
|
|
#+end_important
|
|
|
|
Let's define such matrix using matlab:
|
|
#+begin_src matlab
|
|
M = zeros(length(Opm), 6);
|
|
|
|
for i = 1:length(Opm)
|
|
Ri = [0, Opm(3,i), -Opm(2,i);
|
|
-Opm(3,i), 0, Opm(1,i);
|
|
Opm(2,i), -Opm(1,i), 0];
|
|
M(i, 1:3) = Osm(:,i)';
|
|
M(i, 4:6) = Osm(:,i)'*Ri;
|
|
end
|
|
#+end_src
|
|
|
|
The obtained matrix is shown in Table [[tab:effect_motion_on_meas]].
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
|
|
data2orgtable(M, {'$a_1$', '$a_2$', '$a_3$', '$a_4$', '$a_5$', '$a_6$'}, {'$\dot{x}_x$', '$\dot{x}_y$', '$\dot{x}_z$', '$\dot{\omega}_x$', '$\dot{\omega}_y$', '$\dot{\omega}_z$'}, ' %.2f ');
|
|
#+end_src
|
|
|
|
#+name: tab:effect_motion_on_meas
|
|
#+caption: Effect of a displacement/rotation on the 6 measurements
|
|
#+attr_latex: :environment tabularx :width 0.45\linewidth :align Xcccccc
|
|
#+attr_latex: :center t :booktabs t :float t
|
|
#+RESULTS:
|
|
| | $\dot{x}_x$ | $\dot{x}_y$ | $\dot{x}_z$ | $\dot{\omega}_x$ | $\dot{\omega}_y$ | $\dot{\omega}_z$ |
|
|
|-------+-------------+-------------+-------------+------------------+------------------+------------------|
|
|
| $a_1$ | 0.0 | 1.0 | 0.0 | 0.24 | 0.0 | -0.19 |
|
|
| $a_2$ | 0.0 | 0.0 | 1.0 | -0.19 | 0.19 | 0.0 |
|
|
| $a_3$ | 1.0 | 0.0 | 0.0 | 0.0 | -0.24 | 0.19 |
|
|
| $a_4$ | 0.0 | 0.0 | 1.0 | -0.19 | -0.19 | 0.0 |
|
|
| $a_5$ | 1.0 | 0.0 | 0.0 | 0.0 | -0.24 | -0.19 |
|
|
| $a_6$ | 0.0 | 0.0 | 1.0 | 0.19 | -0.19 | 0.0 |
|
|
|
|
** Compute the transformation matrix from accelerometer measurement to motion of the solid body
|
|
<<sec:transformation_acc_to_motion>>
|
|
|
|
In order to compute the motion of the solid body ${}^O\vec{x}$ with respect to frame $\{O\}$ from the accelerometer measurements $\vec{a}$, we have to inverse the transformation matrix $M$.
|
|
\begin{equation}
|
|
{}^O\vec{x} = M^{-1} \vec{a}
|
|
\end{equation}
|
|
|
|
We therefore need the determinant of $M$ to be non zero:
|
|
#+begin_src matlab :results value replace :exports both :tangle no
|
|
det(M)
|
|
#+end_src
|
|
|
|
The obtained inverse of the matrix is shown in Table [[tab:compute_motion_from_meas]].
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
|
|
data2orgtable(inv(M), {'$\dot{x}_x$', '$\dot{x}_y$', '$\dot{x}_z$', '$\dot{\omega}_x$', '$\dot{\omega}_y$', '$\dot{\omega}_z$'}, {'$a_1$', '$a_2$', '$a_3$', '$a_4$', '$a_5$', '$a_6$'}, ' %.1f ');
|
|
#+end_src
|
|
|
|
#+name: tab:compute_motion_from_meas
|
|
#+caption: Compute the displacement/rotation from the 6 measurements
|
|
#+attr_latex: :environment tabularx :width 0.45\linewidth :align Xcccccc
|
|
#+attr_latex: :center t :booktabs t :float t
|
|
#+RESULTS:
|
|
| | $a_1$ | $a_2$ | $a_3$ | $a_4$ | $a_5$ | $a_6$ |
|
|
|------------------+-------+-------+-------+-------+-------+-------|
|
|
| $\dot{x}_x$ | 0.0 | 0.7 | 0.5 | -0.7 | 0.5 | 0.0 |
|
|
| $\dot{x}_y$ | 1.0 | 0.0 | 0.5 | 0.7 | -0.5 | -0.7 |
|
|
| $\dot{x}_z$ | 0.0 | 0.5 | 0.0 | 0.0 | 0.0 | 0.5 |
|
|
| $\dot{\omega}_x$ | 0.0 | 0.0 | 0.0 | -2.7 | 0.0 | 2.7 |
|
|
| $\dot{\omega}_y$ | 0.0 | 2.7 | 0.0 | -2.7 | 0.0 | 0.0 |
|
|
| $\dot{\omega}_z$ | 0.0 | 0.0 | 2.7 | 0.0 | -2.7 | 0.0 |
|
|
|
|
* Simscape Model
|
|
:PROPERTIES:
|
|
:header-args:matlab+: :tangle matlab/simscape_model.m
|
|
:END:
|
|
<<sec:simscape_model>>
|
|
** Introduction :ignore:
|
|
In this section, the Simscape model of the vibration table is described.
|
|
|
|
#+name: fig:simscape_vibration_table
|
|
#+caption: 3D representation of the simscape model
|
|
#+attr_latex: :width 0.8\linewidth
|
|
[[file:figs/simscape_vibration_table.png]]
|
|
|
|
** Matlab Init :noexport:ignore:
|
|
#+begin_src matlab :tangle no :exports none :results silent :noweb yes :var current_dir=(file-name-directory buffer-file-name)
|
|
<<matlab-dir>>
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none :results silent :noweb yes
|
|
<<matlab-init>>
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no
|
|
addpath('matlab/')
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
% Run meas_transformation.m script in order to get the tranformation (jacobian) matrix
|
|
run('meas_transformation.m')
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :eval no
|
|
%% Used to export to another version
|
|
open('vibration_table')
|
|
sim_versions = {'R2019B'}
|
|
for sim_version = sim_versions
|
|
save_system('vibration_table', ['matlab/' sim_version{1} '/vibration_table.slx'], 'ExportToVersion', sim_version{1});
|
|
save_system('accelerometer_3d', ['matlab/' sim_version{1} '/accelerometer_3d.slx'], 'ExportToVersion', sim_version{1});
|
|
save_system('accelerometer_3d_perfect', ['matlab/' sim_version{1} '/accelerometer_3d_perfect.slx'], 'ExportToVersion', sim_version{1});
|
|
end
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
% Open the Simulink File
|
|
open('vibration_table')
|
|
#+end_src
|
|
|
|
** Simscape Sub-systems
|
|
<<sec:simscape_parameters>>
|
|
|
|
Parameters for sub-components of the simscape model are defined below.
|
|
|
|
*** Springs
|
|
<<sec:simscape_springs>>
|
|
|
|
The 4 springs supporting the suspended optical table are modelled with "bushing joints" having stiffness and damping in the x-y-z directions:
|
|
|
|
#+begin_src matlab
|
|
spring.kx = 1e4; % X- Stiffness [N/m]
|
|
spring.cx = 1e1; % X- Damping [N/(m/s)]
|
|
|
|
spring.ky = 1e4; % Y- Stiffness [N/m]
|
|
spring.cy = 1e1; % Y- Damping [N/(m/s)]
|
|
|
|
spring.kz = 1e4; % Z- Stiffness [N/m]
|
|
spring.cz = 1e1; % Z- Damping [N/(m/s)]
|
|
|
|
spring.z0 = 32e-3; % Equilibrium z-length [m]
|
|
#+end_src
|
|
|
|
*** Inertial Shaker (IS20)
|
|
<<sec:simscape_inertial_shaker>>
|
|
|
|
The inertial shaker is defined as two solid bodies:
|
|
- the "housing" that is fixed to the element that we want to excite
|
|
- the "inertial mass" that is suspended inside the housing
|
|
|
|
The inertial mass is guided inside the housing and an actuator (coil and magnet) can be used to apply a force between the inertial mass and the support.
|
|
The "reacting" force on the support is then used as an excitation.
|
|
|
|
#+name: tab:is20_characteristics
|
|
#+caption: Summary of the IS20 datasheet
|
|
#+attr_latex: :environment tabularx :width 0.4\linewidth :align lX
|
|
#+attr_latex: :center t :booktabs t :float t
|
|
| Characteristic | Value |
|
|
|-----------------+------------|
|
|
| Output Force | 20 N |
|
|
| Frequency Range | 10-3000 Hz |
|
|
| Moving Mass | 0.1 kg |
|
|
| Total Mass | 0.3 kg |
|
|
|
|
From the datasheet in Table [[tab:is20_characteristics]], we can estimate the parameters of the physical shaker.
|
|
|
|
These parameters are defined below
|
|
#+begin_src matlab
|
|
shaker.w0 = 2*pi*10; % Resonance frequency of moving mass [rad/s]
|
|
shaker.m = 0.1; % Moving mass [m]
|
|
shaker.m_tot = 0.3; % Total mass [m]
|
|
shaker.k = shaker.m*shaker.w0^2; % Spring constant [N/m]
|
|
shaker.c = 0.2*sqrt(shaker.k*shaker.m); % Damping [N/(m/s)]
|
|
#+end_src
|
|
|
|
*** 3D accelerometer (356B18)
|
|
<<sec:simscape_accelerometers>>
|
|
|
|
An accelerometer consists of 2 solids:
|
|
- a "housing" rigidly fixed to the measured body
|
|
- an "inertial mass" suspended inside the housing by springs and guided in the measured direction
|
|
|
|
The relative motion between the housing and the inertial mass gives a measurement of the acceleration of the measured body (up to the suspension mode of the inertial mass).
|
|
|
|
#+name: tab:356b18_characteristics
|
|
#+caption: Summary of the 356B18 datasheet
|
|
#+attr_latex: :environment tabularx :width 0.5\linewidth :align lX
|
|
#+attr_latex: :center t :booktabs t :float t
|
|
| Characteristic | Value |
|
|
|---------------------+---------------------|
|
|
| Sensitivity | 0.102 V/(m/s2) |
|
|
| Frequency Range | 0.5 to 3000 Hz |
|
|
| Resonance Frequency | > 20 kHz |
|
|
| Resolution | 0.0005 m/s2 rms |
|
|
| Weight | 0.025 kg |
|
|
| Size | 20.3x26.1x20.3 [mm] |
|
|
|
|
Here are defined the parameters for the triaxial accelerometer:
|
|
#+begin_src matlab
|
|
acc_3d.m = 0.005; % Inertial mass [kg]
|
|
acc_3d.m_tot = 0.025; % Total mass [m]
|
|
|
|
acc_3d.w0 = 2*pi*20e3; % Resonance frequency [rad/s]
|
|
|
|
acc_3d.kx = acc_3d.m*acc_3d.w0^2; % Spring constant [N/m]
|
|
acc_3d.ky = acc_3d.m*acc_3d.w0^2; % Spring constant [N/m]
|
|
acc_3d.kz = acc_3d.m*acc_3d.w0^2; % Spring constant [N/m]
|
|
|
|
acc_3d.cx = 1e2; % Damping [N/(m/s)]
|
|
acc_3d.cy = 1e2; % Damping [N/(m/s)]
|
|
acc_3d.cz = 1e2; % Damping [N/(m/s)]
|
|
#+end_src
|
|
|
|
DC gain between support acceleration and inertial mass displacement is $-m/k$:
|
|
#+begin_src matlab
|
|
acc_3d.g_x = 1/(-acc_3d.m/acc_3d.kx); % [m/s^2/m]
|
|
acc_3d.g_y = 1/(-acc_3d.m/acc_3d.ky); % [m/s^2/m]
|
|
acc_3d.g_z = 1/(-acc_3d.m/acc_3d.kz); % [m/s^2/m]
|
|
#+end_src
|
|
|
|
We also define the sensitivity in order to have the outputs in volts.
|
|
#+begin_src matlab
|
|
acc_3d.gV_x = 0.102; % [V/(m/s^2)]
|
|
acc_3d.gV_y = 0.102; % [V/(m/s^2)]
|
|
acc_3d.gV_z = 0.102; % [V/(m/s^2)]
|
|
#+end_src
|
|
|
|
The problem with using such model for accelerometers is that this adds states to the identified models (2x3 states for each triaxial accelerometer).
|
|
These states represents the dynamics of the suspended inertial mass.
|
|
In the frequency band of interest (few Hz up to ~1 kHz), the dynamics of the inertial mass can be ignore (its resonance is way above 1kHz).
|
|
Therefore, we might as well use idealized "transform sensors" blocks as they will give the same result up to ~20kHz while allowing to reduce the number of identified states.
|
|
|
|
The accelerometer model can be chosen by setting the =type= property:
|
|
#+begin_src matlab
|
|
acc_3d.type = 2; % 1: inertial mass, 2: perfect
|
|
#+end_src
|
|
|
|
** Identification
|
|
<<sec:simscape_parameters>>
|
|
|
|
*** Number of states
|
|
Let's first use perfect 3d accelerometers:
|
|
#+begin_src matlab
|
|
acc_3d.type = 2; % 1: inertial mass, 2: perfect
|
|
#+end_src
|
|
|
|
And identify the dynamics from the shaker force to the measured accelerations:
|
|
#+begin_src matlab
|
|
%% Name of the Simulink File
|
|
mdl = 'vibration_table';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/F'], 1, 'openinput'); io_i = io_i + 1;
|
|
io(io_i) = linio([mdl, '/acc'], 1, 'openoutput'); io_i = io_i + 1;
|
|
|
|
%% Run the linearization
|
|
Gp = linearize(mdl, io);
|
|
Gp.InputName = {'F'};
|
|
Gp.OutputName = {'a1', 'a2', 'a3', 'a4', 'a5', 'a6'};
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results output replace :exports results :tangle no
|
|
size(Gp)
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: size(Gp)
|
|
: State-space model with 6 outputs, 1 inputs, and 12 states.
|
|
|
|
We indeed have the 12 states corresponding to the 6 DoF of the suspended optical table.
|
|
|
|
Let's now consider the inertial masses for the triaxial accelerometers:
|
|
#+begin_src matlab
|
|
acc_3d.type = 1; % 1: inertial mass, 2: perfect
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
%% Name of the Simulink File
|
|
mdl = 'vibration_table';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/F'], 1, 'openinput'); io_i = io_i + 1;
|
|
io(io_i) = linio([mdl, '/acc'], 1, 'openoutput'); io_i = io_i + 1;
|
|
|
|
%% Run the linearization
|
|
Ga = linearize(mdl, io);
|
|
Ga.InputName = {'F'};
|
|
Ga.OutputName = {'a1', 'a2', 'a3', 'a4', 'a5', 'a6'};
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results output replace :exports results :tangle no
|
|
size(Ga)
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: size(Ga)
|
|
: State-space model with 6 outputs, 1 inputs, and 30 states.
|
|
|
|
And we can see that 18 states have been added.
|
|
This corresponds to 6 states for each triaxial accelerometers.
|
|
|
|
*** Resonance frequencies and mode shapes
|
|
Let's now identify the resonance frequency and mode shapes associated with the suspension modes of the optical table.
|
|
|
|
#+begin_src matlab
|
|
acc_3d.type = 2; % 1: inertial mass, 2: perfect
|
|
|
|
%% Name of the Simulink File
|
|
mdl = 'vibration_table';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/F'], 1, 'openinput'); io_i = io_i + 1;
|
|
io(io_i) = linio([mdl, '/acc_O'], 1, 'openoutput'); io_i = io_i + 1;
|
|
|
|
%% Run the linearization
|
|
G = linearize(mdl, io);
|
|
G.InputName = {'F'};
|
|
G.OutputName = {'ax', 'ay', 'az', 'wx', 'wy', 'wz'};
|
|
#+end_src
|
|
|
|
Compute the resonance frequencies
|
|
#+begin_src matlab
|
|
ws = eig(G.A);
|
|
ws = ws(imag(ws) > 0);
|
|
#+end_src
|
|
|
|
And the associated response of the optical table
|
|
#+begin_src matlab
|
|
x_mod = zeros(6, 6);
|
|
|
|
for i = 1:length(ws)
|
|
xi = evalfr(G, ws(i));
|
|
x_mod(:,i) = xi./norm(xi);
|
|
end
|
|
#+end_src
|
|
|
|
The results are shown in Table [[tab:mode_shapes]].
|
|
The motion associated to the mode shapes are just indicative.
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
|
|
data2orgtable([flip(imag(ws)/2/pi)'; flip(abs(x_mod),2)], {'$\omega_0$ [Hz]', 'x', 'y', 'z', 'Rx', 'Ry', 'Rz'}, {}, ' %.1f ');
|
|
#+end_src
|
|
|
|
#+name: tab:mode_shapes
|
|
#+caption: Resonance frequency and approximation of the mode shapes
|
|
#+attr_latex: :environment tabularx :width 0.4\linewidth :align Xcccccc
|
|
#+attr_latex: :center t :booktabs t :float t
|
|
#+RESULTS:
|
|
| $\omega_0$ [Hz] | 5.6 | 5.6 | 5.7 | 8.2 | 8.2 | 8.2 |
|
|
|-----------------+-----+-----+-----+-----+-----+-----|
|
|
| x | 0.1 | 0.7 | 0.0 | 0.0 | 0.2 | 0.0 |
|
|
| y | 0.7 | 0.1 | 0.0 | 0.0 | 0.0 | 0.2 |
|
|
| z | 0.0 | 0.0 | 1.0 | 0.0 | 0.0 | 0.0 |
|
|
| Rx | 0.7 | 0.1 | 0.0 | 0.0 | 0.1 | 1.0 |
|
|
| Ry | 0.1 | 0.7 | 0.0 | 0.0 | 1.0 | 0.1 |
|
|
| Rz | 0.0 | 0.0 | 0.0 | 1.0 | 0.0 | 0.0 |
|
|
|
|
** Verify transformation
|
|
#+begin_src matlab
|
|
%% Options for Linearized
|
|
options = linearizeOptions;
|
|
options.SampleTime = 0;
|
|
|
|
%% Name of the Simulink File
|
|
mdl = 'vibration_table';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/F'], 1, 'openinput'); io_i = io_i + 1;
|
|
io(io_i) = linio([mdl, '/acc'], 1, 'openoutput'); io_i = io_i + 1;
|
|
io(io_i) = linio([mdl, '/acc_O'], 1, 'openoutput'); io_i = io_i + 1;
|
|
|
|
%% Run the linearization
|
|
G = linearize(mdl, io, 0.0, options);
|
|
G.InputName = {'F'};
|
|
G.OutputName = {'a1', 'a2', 'a3', 'a4', 'a5', 'a6', ...
|
|
'Dx', 'Dy', 'Dz', 'Rx', 'Ry', 'Rz'};
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
G_acc = inv(M)*G(1:6, 1);
|
|
G_id = G(7:12, 1);
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
bodeFig({G_acc(1), G_id(1)})
|
|
bodeFig({G_acc(2), G_id(2)})
|
|
bodeFig({G_acc(3), G_id(3)})
|
|
bodeFig({G_acc(4), G_id(4)})
|
|
bodeFig({G_acc(5), G_id(5)})
|
|
bodeFig({G_acc(6), G_id(6)})
|
|
#+end_src
|
|
|
|
* Identification of the table's dynamics
|
|
<<sec:table_dynamics>>
|
|
** Matlab Init :noexport:ignore:
|
|
#+begin_src matlab :tangle no :exports none :results silent :noweb yes :var current_dir=(file-name-directory buffer-file-name)
|
|
<<matlab-dir>>
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none :results silent :noweb yes
|
|
<<matlab-init>>
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no
|
|
addpath('matlab/')
|
|
#+end_src
|
|
|
|
* Bibliography :ignore:
|
|
#+latex: \printbibliography
|