1453 lines
46 KiB
Org Mode
1453 lines
46 KiB
Org Mode
#+TITLE: Test Bench - Amplified Piezoelectric Actuator
|
|
:DRAWER:
|
|
#+STARTUP: overview
|
|
|
|
#+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>
|
|
|
|
#+HTML_MATHJAX: align: center tagside: right font: TeX
|
|
|
|
#+PROPERTY: header-args:matlab :session *MATLAB*
|
|
#+PROPERTY: header-args:matlab+ :comments org
|
|
#+PROPERTY: header-args:matlab+ :results none
|
|
#+PROPERTY: header-args:matlab+ :exports both
|
|
#+PROPERTY: header-args:matlab+ :eval no-export
|
|
#+PROPERTY: header-args:matlab+ :output-dir figs
|
|
#+PROPERTY: header-args:matlab+ :tangle no
|
|
#+PROPERTY: header-args:matlab+ :mkdirp yes
|
|
|
|
#+PROPERTY: header-args:shell :eval no-export
|
|
:END:
|
|
|
|
* Introduction :ignore:
|
|
|
|
- Section [[sec:experimental_setup]]:
|
|
- Section [[sec:simscape_model]]:
|
|
- Section [[]]:
|
|
- Section [[]]:
|
|
- Section [[]]:
|
|
|
|
* Experimental Setup
|
|
<<sec:experimental_setup>>
|
|
|
|
A schematic of the test-bench is shown in Figure [[fig:test_bench_apa_schematic]].
|
|
|
|
A mass can be vertically moved using the amplified piezoelectric actuator (APA95ML).
|
|
The displacement of the mass (relative to the mechanical frame) is measured by the interferometer.
|
|
|
|
The APA95ML has three stacks that can be used as actuator or as sensors.
|
|
|
|
Pictures of the test bench are shown in Figure [[fig:setup_picture]] and [[fig:setup_zoom]].
|
|
|
|
#+name: fig:test_bench_apa_schematic
|
|
#+caption: Schematic of the Setup
|
|
[[file:figs/test_bench_apa_schematic.png]]
|
|
|
|
#+name: fig:setup_picture
|
|
#+caption: Picture of the Setup
|
|
[[file:figs/setup_picture.png]]
|
|
|
|
#+name: fig:setup_zoom
|
|
#+caption: Zoom on the APA
|
|
[[file:figs/setup_zoom.png]]
|
|
|
|
#+begin_note
|
|
Here are the equipment used in the test bench:
|
|
- Attocube interferometer ([[file:doc/IDS3010.pdf][doc]])
|
|
- Cedrat Amplified Piezoelectric Actuator APA95ML ([[file:doc/APA95ML.pdf][doc]])
|
|
- Voltage Amplifier LA75B ([[file:doc/LA75B.pdf][doc]])
|
|
- Speedgoat IO131 with 16bits ADC and DAC ([[file:doc/IO130 IO131 OEM Datasheet.pdf][doc]])
|
|
- Low Noise Voltage Preamplifier from Ametek ([[file:doc/model_5113.pdf][doc]])
|
|
#+end_note
|
|
|
|
* Estimation of electrical/mechanical relationships
|
|
** Introduction :ignore:
|
|
In order to correctly model the piezoelectric actuator with Simscape, we need to determine:
|
|
1. $g_a$: the ratio of the generated force $F_a$ to the supply voltage $V_a$ across the piezoelectric stack
|
|
2. $g_s$: the ratio of the generated voltage $V_s$ across the piezoelectric stack when subject to a strain $\Delta h$
|
|
|
|
We estimate $g_a$ and $g_s$ using different approaches:
|
|
1. Section [[sec:estimation_datasheet]]: $g_a$ is estimated from the datasheet of the piezoelectric stack
|
|
2. Section [[sec:estimation_piezo_params]]: $g_a$ and $g_s$ are estimated using the piezoelectric constants
|
|
3. Section [[sec:estimation_identification]]: $g_a$ and $g_s$ are estimated experimentally
|
|
|
|
** 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/mat/');
|
|
addpath('./matlab/');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :eval no
|
|
addpath('./mat/');
|
|
#+end_src
|
|
|
|
** Estimation from Data-sheet
|
|
<<sec:estimation_datasheet>>
|
|
|
|
The stack parameters taken from the data-sheet are shown in Table [[tab:stack_parameters]].
|
|
|
|
#+name: tab:stack_parameters
|
|
#+caption: Stack Parameters
|
|
| Parameter | Unit | Value |
|
|
|----------------+-----------+----------|
|
|
| Nominal Stroke | $\mu m$ | 20 |
|
|
| Blocked force | $N$ | 4700 |
|
|
| Stiffness | $N/\mu m$ | 235 |
|
|
| Voltage Range | $V$ | -20..150 |
|
|
| Capacitance | $\mu F$ | 4.4 |
|
|
| Length | $mm$ | 20 |
|
|
| Stack Area | $mm^2$ | 10x10 |
|
|
|
|
Let's compute the generated force
|
|
|
|
The stroke is $L_{\max} = 20\mu m$ for a voltage range of $V_{\max} = 170 V$.
|
|
Furthermore, the stiffness is $k_a = 235 \cdot 10^6 N/m$.
|
|
|
|
The relation between the applied voltage and the generated force can be estimated as follows:
|
|
\begin{equation}
|
|
g_a = k_a \frac{L_{\max}}{V_{\max}}
|
|
\end{equation}
|
|
|
|
#+begin_src matlab
|
|
ka = 235e6; % [N/m]
|
|
Lmax = 20e-6; % [m]
|
|
Vmax = 170; % [V]
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results value replace
|
|
ka*Lmax/Vmax % [N/V]
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: 27.647
|
|
|
|
From the parameters of the stack, it seems not possible to estimate the relation between the strain and the generated voltage.
|
|
|
|
** Estimation from Piezoelectric parameters
|
|
<<sec:estimation_piezo_params>>
|
|
|
|
In order to make the conversion from applied voltage to generated force or from the strain to the generated voltage, we need to defined some parameters corresponding to the piezoelectric material:
|
|
#+begin_src matlab
|
|
d33 = 300e-12; % Strain constant [m/V]
|
|
n = 80; % Number of layers per stack
|
|
ka = 235e6; % Stack stiffness [N/m]
|
|
#+end_src
|
|
|
|
The ratio of the developed force to applied voltage is:
|
|
#+name: eq:piezo_voltage_to_force
|
|
\begin{equation}
|
|
F_a = g_a V_a, \quad g_a = d_{33} n k_a
|
|
\end{equation}
|
|
where:
|
|
- $F_a$: developed force in [N]
|
|
- $n$: number of layers of the actuator stack
|
|
- $d_{33}$: strain constant in [m/V]
|
|
- $k_a$: actuator stack stiffness in [N/m]
|
|
- $V_a$: applied voltage in [V]
|
|
|
|
If we take the numerical values, we obtain:
|
|
#+begin_src matlab
|
|
ga = d33*n*ka; % [N/V]
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results value replace :exports results :tangle no
|
|
sprintf('ga = %.1f [N/V]', ga)
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: ga = 5.6 [N/V]
|
|
|
|
|
|
From cite:fleming14_desig_model_contr_nanop_system (page 123), the relation between relative displacement of the sensor stack and generated voltage is:
|
|
#+name: eq:piezo_strain_to_voltage
|
|
\begin{equation}
|
|
V_s = \frac{d_{33}}{\epsilon^T s^D n} \Delta h
|
|
\end{equation}
|
|
where:
|
|
- $V_s$: measured voltage in [V]
|
|
- $d_{33}$: strain constant in [m/V]
|
|
- $\epsilon^T$: permittivity under constant stress in [F/m]
|
|
- $s^D$: elastic compliance under constant electric displacement in [m^2/N]
|
|
- $n$: number of layers of the sensor stack
|
|
- $\Delta h$: relative displacement in [m]
|
|
|
|
If we take the numerical values, we obtain:
|
|
#+begin_src matlab
|
|
d33 = 300e-12; % Strain constant [m/V]
|
|
n = 80; % Number of layers per stack
|
|
eT = 5.3e-9; % Permittivity under constant stress [F/m]
|
|
sD = 2e-11; % Compliance under constant electric displacement [m2/N]
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
gs = d33/(eT*sD*n); % [V/m]
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results value replace :exports results :tangle no
|
|
sprintf('gs = %.1f [V/um]', 1e-6*gs)
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: gs = 35.4 [V/um]
|
|
|
|
** Estimation from Experiment
|
|
<<sec:estimation_identification>>
|
|
*** Introduction :ignore:
|
|
The idea here is to obtain the parameters $g_a$ and $g_s$ from the comparison of an experimental identification and the identification using Simscape.
|
|
|
|
Using the experimental identification, we can easily obtain the gain from the applied voltage to the generated displacement, but not to the generated force.
|
|
However, from the Simscape model, we can easily have the link from the generated force to the displacement, them we can computed $g_a$.
|
|
|
|
Similarly, it is fairly easy to experimentally obtain the gain from the stack displacement to the generated voltage across the stack.
|
|
To link that to the strain of the sensor stack, the simscape model is used.
|
|
|
|
*** From actuator voltage $V_a$ to actuator force $F_a$
|
|
The data from the identification test is loaded.
|
|
#+begin_src matlab
|
|
load('apa95ml_5kg_Amp_E505.mat', 't', 'um', 'y');
|
|
|
|
% Any offset value is removed:
|
|
um = detrend(um, 0); % Amplifier Input Voltage [V]
|
|
y = detrend(y , 0); % Mass displacement [m]
|
|
#+end_src
|
|
|
|
Now we add a factor 10 to take into account the gain of the voltage amplifier and thus obtain the voltage across the piezoelectric stack.
|
|
#+begin_src matlab
|
|
um = 10*um; % Stack Actuator Input Voltage [V]
|
|
#+end_src
|
|
|
|
Then, the transfer function from the stack voltage to the vertical displacement is computed.
|
|
#+begin_src matlab
|
|
Ts = t(end)/(length(t)-1);
|
|
Fs = 1/Ts;
|
|
|
|
win = hanning(ceil(1*Fs));
|
|
|
|
[tf_est, f] = tfestimate(um, y, win, [], [], 1/Ts);
|
|
#+end_src
|
|
|
|
The gain from input voltage of the stack to the vertical displacement is determined:
|
|
#+begin_src matlab
|
|
g_d_Va = 4e-7; % [m/V]
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(0, 4, 1000);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(f, abs(tf_est), 'DisplayName', 'Identification')
|
|
plot([f(2) f(end)], [g_d_Va g_d_Va], '--', 'DisplayName', sprintf('$d/V_a \\approx %.0g\\ m/V$', g_d_Va))
|
|
hold off;
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
xlabel('Frequency [Hz]'); ylabel('Amplitude $d/V_a$ [m/V]');
|
|
xlim([10, 5e3]); ylim([1e-8, 1e-5]);
|
|
legend('location', 'southwest');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/gain_Va_to_d.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:gain_Va_to_d
|
|
#+caption: Transfer function from actuator stack voltage $V_a$ to vertical displacement of the mass $d$
|
|
#+RESULTS:
|
|
[[file:figs/gain_Va_to_d.png]]
|
|
|
|
Then, the transfer function from forces applied by the stack actuator to the vertical displacement of the mass is identified from the Simscape model.
|
|
#+begin_src matlab :exports none
|
|
% Load Parameters
|
|
K = readmatrix('APA95ML_K.CSV');
|
|
M = readmatrix('APA95ML_M.CSV');
|
|
[int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes('APA95ML_out_nodes_3D.txt');
|
|
|
|
% Define parameters just for the simulation
|
|
ga = 1; % [N/V]
|
|
gs = 1; % [V/m]
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
m = 5.5;
|
|
|
|
%% Name of the Simulink File
|
|
mdl = 'piezo_amplified_3d';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/Fa'], 1, 'openinput'); io_i = io_i + 1; % Actuator Force [N]
|
|
io(io_i) = linio([mdl, '/y'], 1, 'openoutput'); io_i = io_i + 1; % Vertical Displacement [m]
|
|
|
|
Gd = linearize(mdl, io);
|
|
#+end_src
|
|
|
|
The DC gain the the identified dynamics
|
|
#+begin_src matlab
|
|
g_d_Fa = abs(dcgain(Gd)); % [m/N]
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results value replace :exports results :tangle no
|
|
sprintf('G_d_Fa = %.2g [m/N]', g_d_Fa)
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: G_d_Fa = 1.2e-08 [m/N]
|
|
|
|
|
|
And finally, the gain $g_a$ from the the actuator voltage $V_a$ to the generated force $F_a$ can be computed:
|
|
\begin{equation}
|
|
g_a = \frac{F_a}{V_a} = \frac{F_a}{d} \frac{d}{V_a}
|
|
\end{equation}
|
|
|
|
#+begin_src matlab
|
|
ga = g_d_Va/g_d_Fa;
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results value replace :exports results :tangle no
|
|
sprintf('ga = %.1f [N/V]', ga)
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: ga = 33.7 [N/V]
|
|
|
|
The obtained comparison between the Simscape model and the identified dynamics is shown in Figure [[fig:compare_Gd_id_simscape]].
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(1, 4, 1000);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(f, abs(tf_est), 'DisplayName', 'Identification')
|
|
plot(f, ga*abs(squeeze(freqresp(Gd, f, 'Hz'))), 'DisplayName', 'Simscape Model')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude $d/V_a$ [m/V]'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
xlim([1, 5e3]); ylim([1e-9, 1e-5]);
|
|
legend('location', 'southwest');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/compare_Gd_id_simscape.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:compare_Gd_id_simscape
|
|
#+caption: Comparison of the identified transfer function between actuator voltage $V_a$ and vertical mass displacement $d$
|
|
#+RESULTS:
|
|
[[file:figs/compare_Gd_id_simscape.png]]
|
|
|
|
*** From stack strain $\Delta h$ to generated voltage $V_s$
|
|
Now, the gain from the stack strain $\Delta h$ to the generated voltage $V_s$ is estimated.
|
|
|
|
We can determine the gain from actuator voltage $V_a$ to sensor voltage $V_s$ thanks to the identification.
|
|
Using the simscape model, we can have the transfer function from the actuator voltage $V_a$ (using the previously estimated gain $g_a$) to the sensor stack strain $\Delta h$.
|
|
|
|
Finally, using these two values, we can compute the gain $g_s$ from the stack strain $\Delta h$ to the generated Voltage $V_s$.
|
|
|
|
Identification data is loaded.
|
|
#+begin_src matlab
|
|
load('apa95ml_5kg_2a_1s.mat', 't', 'u', 'v');
|
|
|
|
u = detrend(u, 0); % Input Voltage of the Amplifier [V]
|
|
v = detrend(v, 0); % Voltage accross the stack sensor [V]
|
|
#+end_src
|
|
|
|
Here, an amplifier with a gain of 20 is used.
|
|
#+begin_src matlab
|
|
u = 20*u; % Input Voltage of the Amplifier [V]
|
|
#+end_src
|
|
|
|
Then, the transfer function from $V_a$ to $V_s$ is identified and its DC gain is estimated (Figure [[fig:gain_Va_to_Vs]]).
|
|
#+begin_src matlab
|
|
Ts = t(end)/(length(t)-1);
|
|
Fs = 1/Ts;
|
|
|
|
win = hann(ceil(10/Ts));
|
|
|
|
[tf_est, f] = tfestimate(u, v, win, [], [], 1/Ts);
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
g_Vs_Va = 0.022; % [V/V]
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(1, 4, 1000);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(f, abs(tf_est), 'DisplayName', 'Identification')
|
|
plot([f(2); f(end)], [g_Vs_Va; g_Vs_Va], '--', 'DisplayName', sprintf('$V_s/V_a \\approx %.0g\\ m/V$', g_Vs_Va))
|
|
hold off;
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude $V_s/V_a$ [V/V]'); xlabel('Frequency [Hz]');
|
|
xlim([5e-1 5e3]); ylim([1e-3, 1e1]);
|
|
legend('location', 'northwest');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/gain_Va_to_Vs.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:gain_Va_to_Vs
|
|
#+caption: Transfer function from actuator stack voltage $V_a$ to sensor stack voltage $V_s$
|
|
#+RESULTS:
|
|
[[file:figs/gain_Va_to_Vs.png]]
|
|
|
|
|
|
Now the transfer function from the actuator stack voltage to the sensor stack strain is estimated using the Simscape model.
|
|
#+begin_src matlab
|
|
m = 5.5;
|
|
|
|
%% Name of the Simulink File
|
|
mdl = 'piezo_amplified_3d';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/Va'], 1, 'openinput'); io_i = io_i + 1; % Actuator Voltage [V]
|
|
io(io_i) = linio([mdl, '/dL'], 1, 'openoutput'); io_i = io_i + 1; % Sensor Stack displacement [m]
|
|
|
|
Gf = linearize(mdl, io);
|
|
#+end_src
|
|
|
|
The gain from the actuator stack voltage to the sensor stack strain is estimated below.
|
|
#+begin_src matlab
|
|
G_dh_Va = abs(dcgain(Gf));
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results value replace :exports results :tangle no
|
|
sprintf('G_dh_Va = %.2g [m/V]', G_dh_Va);
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: G_dh_Va = 6.2e-09 [m/V]
|
|
|
|
And finally, the gain $g_s$ from the sensor stack strain to the generated voltage can be estimated:
|
|
\begin{equation}
|
|
g_s = \frac{V_s}{\Delta h} = \frac{V_s}{V_a} \frac{V_a}{\Delta h}
|
|
\end{equation}
|
|
|
|
#+begin_src matlab
|
|
gs = g_Vs_Va/G_dh_Va; % [V/m]
|
|
#+end_src
|
|
|
|
#+begin_src matlab :results value replace :exports results :tangle no
|
|
sprintf('gs = %.1f [V/um]', 1e-6*gs)
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
: gs = 3.5 [V/um]
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(1, 4, 1000);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(f, abs(tf_est), 'DisplayName', 'Identification')
|
|
plot(f, gs*abs(squeeze(freqresp(Gf, f, 'Hz'))), 'DisplayName', 'Simscape Model')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude $V_s/V_a$ [V/V]'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
xlim([1, 5e3]); ylim([1e-3, 1e1]);
|
|
legend('location', 'southwest');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/compare_Gf_id_simscape.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:compare_Gf_id_simscape
|
|
#+caption: Comparison of the identified transfer function between actuator voltage $V_a$ and sensor stack voltage
|
|
#+RESULTS:
|
|
[[file:figs/compare_Gf_id_simscape.png]]
|
|
|
|
** Conclusion
|
|
The obtained parameters $g_a$ and $g_s$ are not consistent between the different methods.
|
|
|
|
The one using the experimental data are saved and further used.
|
|
#+begin_src matlab :tangle no
|
|
save('./matlab/mat/apa95ml_params.mat', 'ga', 'gs');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none :eval no
|
|
save('./mat/apa95ml_params.mat', 'ga', 'gs');
|
|
#+end_src
|
|
|
|
* Simscape model of the test-bench
|
|
:PROPERTIES:
|
|
:header-args:matlab+: :tangle matlab/simscape_model.m
|
|
:header-args:matlab+: :comments org :mkdirp yes
|
|
:END:
|
|
<<sec:simscape_model>>
|
|
|
|
** Introduction :ignore:
|
|
The idea here is to model the test-bench using Simscape.
|
|
|
|
Whereas the suspended mass and metrology frame can be considered as rigid bodies in the frequency range of interest, the Amplified Piezoelectric Actuator (APA) is flexible.
|
|
|
|
To model the APA, a Finite Element Model (FEM) is used (Figure [[fig:APA95ML_FEM]]) and imported into Simscape.
|
|
|
|
#+name: fig:APA95ML_FEM
|
|
#+caption: Finite Element Model of the APA95ML
|
|
[[file:figs/APA95ML_FEM.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/mat/');
|
|
addpath('./matlab/');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :eval no
|
|
addpath('./mat/');
|
|
#+end_src
|
|
|
|
** Import Mass Matrix, Stiffness Matrix, and Interface Nodes Coordinates
|
|
We first extract the stiffness and mass matrices.
|
|
#+begin_src matlab
|
|
K = readmatrix('APA95ML_K.CSV');
|
|
M = readmatrix('APA95ML_M.CSV');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no
|
|
data2orgtable(K(1:10, 1:10), {}, {}, ' %.1g ');
|
|
#+end_src
|
|
|
|
#+caption: First 10x10 elements of the Stiffness matrix
|
|
#+RESULTS:
|
|
| 300000000.0 | -1000.0 | -30000.0 | -40.0 | 70000.0 | 300.0 | 20000000.0 | -30.0 | -5000.0 | 5 |
|
|
| -1000.0 | 50000000.0 | -7000.0 | 800000.0 | -20.0 | 300.0 | 3000.0 | 5000000.0 | 400.0 | -40000.0 |
|
|
| -30000.0 | -7000.0 | 100000000.0 | -200.0 | -60.0 | 70.0 | 3000.0 | 3000.0 | -8000000.0 | -30.0 |
|
|
| -40.0 | 800000.0 | -200.0 | 20000.0 | -0.4 | 4 | 30.0 | 40000.0 | 7 | -300.0 |
|
|
| 70000.0 | -20.0 | -60.0 | -0.4 | 3000.0 | 1 | -6000.0 | 10.0 | 8 | -0.1 |
|
|
| 300.0 | 300.0 | 70.0 | 4 | 1 | 40000.0 | -10.0 | -10.0 | 30.0 | 0.1 |
|
|
| 20000000.0 | 3000.0 | 3000.0 | 30.0 | -6000.0 | -10.0 | 300000000.0 | 2000.0 | 9000.0 | -100.0 |
|
|
| -30.0 | 5000000.0 | 3000.0 | 40000.0 | 10.0 | -10.0 | 2000.0 | 50000000.0 | -3000.0 | -800000.0 |
|
|
| -5000.0 | 400.0 | -8000000.0 | 7 | 8 | 30.0 | 9000.0 | -3000.0 | 100000000.0 | 100.0 |
|
|
| 5 | -40000.0 | -30.0 | -300.0 | -0.1 | 0.1 | -100.0 | -800000.0 | 100.0 | 20000.0 |
|
|
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no
|
|
data2orgtable(M(1:10, 1:10), {}, {}, ' %.1g ');
|
|
#+end_src
|
|
|
|
#+caption: First 10x10 elements of the Mass matrix
|
|
#+RESULTS:
|
|
| 0.03 | 7e-08 | 2e-06 | -3e-09 | -0.0002 | -6e-08 | -0.001 | 8e-07 | 6e-07 | -8e-09 |
|
|
| 7e-08 | 0.02 | -1e-06 | 9e-05 | -3e-09 | -4e-09 | -1e-06 | -0.0006 | -4e-08 | 5e-06 |
|
|
| 2e-06 | -1e-06 | 0.02 | -3e-08 | -4e-08 | 1e-08 | 1e-07 | -2e-07 | 0.0003 | 1e-09 |
|
|
| -3e-09 | 9e-05 | -3e-08 | 1e-06 | -3e-11 | -3e-13 | -7e-09 | -5e-06 | -3e-10 | 3e-08 |
|
|
| -0.0002 | -3e-09 | -4e-08 | -3e-11 | 2e-06 | 6e-10 | 2e-06 | -7e-09 | -2e-09 | 7e-11 |
|
|
| -6e-08 | -4e-09 | 1e-08 | -3e-13 | 6e-10 | 1e-06 | 1e-08 | 3e-09 | -2e-09 | 2e-13 |
|
|
| -0.001 | -1e-06 | 1e-07 | -7e-09 | 2e-06 | 1e-08 | 0.03 | 4e-08 | -2e-06 | 8e-09 |
|
|
| 8e-07 | -0.0006 | -2e-07 | -5e-06 | -7e-09 | 3e-09 | 4e-08 | 0.02 | -9e-07 | -9e-05 |
|
|
| 6e-07 | -4e-08 | 0.0003 | -3e-10 | -2e-09 | -2e-09 | -2e-06 | -9e-07 | 0.02 | 2e-08 |
|
|
| -8e-09 | 5e-06 | 1e-09 | 3e-08 | 7e-11 | 2e-13 | 8e-09 | -9e-05 | 2e-08 | 1e-06 |
|
|
|
|
|
|
Then, we extract the coordinates of the interface nodes.
|
|
#+begin_src matlab
|
|
[int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes('APA95ML_out_nodes_3D.txt');
|
|
#+end_src
|
|
|
|
The interface nodes are shown in Figure [[fig:APA95ML_nodes_1]] and their coordinates are listed in Table [[tab:apa95ml_nodes_coordinates]].
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no
|
|
data2orgtable([length(n_i); length(int_i); size(M,1) - 6*length(int_i); size(M,1)], {'Total number of Nodes', 'Number of interface Nodes', 'Number of Modes', 'Size of M and K matrices'}, {}, ' %.0f ');
|
|
#+end_src
|
|
|
|
#+RESULTS:
|
|
| Total number of Nodes | 7 |
|
|
| Number of interface Nodes | 7 |
|
|
| Number of Modes | 6 |
|
|
| Size of M and K matrices | 48 |
|
|
|
|
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
|
|
data2orgtable([[1:length(int_i)]', int_i, int_xyz], {}, {'Node i', 'Node Number', 'x [m]', 'y [m]', 'z [m]'}, ' %f ');
|
|
#+end_src
|
|
|
|
#+name: tab:apa95ml_nodes_coordinates
|
|
#+caption: Coordinates of the interface nodes
|
|
#+RESULTS:
|
|
| Node i | Node Number | x [m] | y [m] | z [m] |
|
|
|--------+-------------+--------+-------+-----------|
|
|
| 1.0 | 40467.0 | 0.0 | 0.0 | 0.029997 |
|
|
| 2.0 | 40469.0 | 0.0 | 0.0 | -0.029997 |
|
|
| 3.0 | 40470.0 | -0.035 | 0.0 | 0.0 |
|
|
| 4.0 | 40475.0 | -0.015 | 0.0 | 0.0 |
|
|
| 5.0 | 40476.0 | -0.005 | 0.0 | 0.0 |
|
|
| 6.0 | 40477.0 | 0.015 | 0.0 | 0.0 |
|
|
| 7.0 | 40478.0 | 0.035 | 0.0 | 0.0 |
|
|
|
|
#+name: fig:APA95ML_nodes_1
|
|
#+caption: Interface Nodes chosen for the APA95ML
|
|
[[file:figs/APA95ML_nodes_1.png]]
|
|
|
|
Using =K=, =M= and =int_xyz=, we can use the =Reduced Order Flexible Solid= simscape block.
|
|
|
|
** Simscape Model
|
|
The flexible element is imported using the =Reduced Order Flexible Solid= Simscape block.
|
|
|
|
To model the actuator, an =Internal Force= block is added between the nodes 3 and 12.
|
|
A =Relative Motion Sensor= block is added between the nodes 1 and 2 to measure the displacement and the amplified piezo.
|
|
|
|
One mass is fixed at one end of the piezo-electric stack actuator, the other end is fixed to the world frame.
|
|
#+begin_src matlab
|
|
m = 5;
|
|
#+end_src
|
|
|
|
** Dynamics from Actuator Voltage to Vertical Mass Displacement
|
|
The identified dynamics is shown in Figure [[fig:dynamics_act_disp_comp_mass]].
|
|
|
|
#+begin_src matlab
|
|
%% Name of the Simulink File
|
|
mdl = 'piezo_amplified_3d';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/Va'], 1, 'openinput'); io_i = io_i + 1; % Actuator Voltage [V]
|
|
io(io_i) = linio([mdl, '/y'], 1, 'openoutput'); io_i = io_i + 1; % Vertical Displacement [m]
|
|
|
|
Ghm = -linearize(mdl, io);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(0, 4, 5000);
|
|
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(freqs, abs(squeeze(freqresp(Ghm, freqs, 'Hz'))));
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(Ghm, freqs, 'Hz')))));
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
yticks(-360:90:360);
|
|
ylim([-360 0]);
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([freqs(1), freqs(end)]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/dynamics_act_disp_comp_mass.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:dynamics_act_disp_comp_mass
|
|
#+caption: Dynamics from $F$ to $d$ without a payload and with a 5kg payload
|
|
#+RESULTS:
|
|
[[file:figs/dynamics_act_disp_comp_mass.png]]
|
|
|
|
** Dynamics from Actuator Voltage to Force Sensor Voltage
|
|
The obtained dynamics is shown in Figure [[fig:dynamics_force_force_sensor_comp_mass]].
|
|
|
|
#+begin_src matlab
|
|
%% Name of the Simulink File
|
|
mdl = 'piezo_amplified_3d';
|
|
|
|
%% Input/Output definition
|
|
clear io; io_i = 1;
|
|
io(io_i) = linio([mdl, '/Va'], 1, 'openinput'); io_i = io_i + 1; % Voltage Actuator [V]
|
|
io(io_i) = linio([mdl, '/Vs'], 1, 'openoutput'); io_i = io_i + 1; % Sensor Voltage [V]
|
|
|
|
Gfm = linearize(mdl, io);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(1, 5, 1000);
|
|
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(freqs, abs(squeeze(freqresp(Gfm, freqs, 'Hz'))));
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
ylabel('Amplitude'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(Gfm, freqs, 'Hz')))));
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
|
|
yticks(-360:90:360);
|
|
ylim([-390 30]);
|
|
xlabel('Frequency [Hz]'); ylabel('Phase [deg]');
|
|
hold off;
|
|
linkaxes([ax1,ax2],'x');
|
|
xlim([freqs(1), freqs(end)]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/dynamics_force_force_sensor_comp_mass.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:dynamics_force_force_sensor_comp_mass
|
|
#+caption: Dynamics from $F$ to $F_m$ for $m=0$ and $m = 10kg$
|
|
#+RESULTS:
|
|
[[file:figs/dynamics_force_force_sensor_comp_mass.png]]
|
|
|
|
** Save Data for further use
|
|
#+begin_src matlab :tangle no
|
|
save('matlab/mat/fem_simscape_models.mat', 'Ghm', 'Gfm')
|
|
#+end_src
|
|
|
|
#+begin_src matlab :eval no
|
|
save('mat/fem_simscape_models.mat', 'Ghm', 'Gfm')
|
|
#+end_src
|
|
|
|
* Huddle Test
|
|
:PROPERTIES:
|
|
:header-args:matlab+: :tangle matlab/huddle_test.m
|
|
:header-args:matlab+: :comments org :mkdirp yes
|
|
:END:
|
|
<<sec:huddle_test>>
|
|
** Introduction :ignore:
|
|
** 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/mat/');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :eval no
|
|
addpath('./mat/');
|
|
#+end_src
|
|
|
|
** Load Data :noexport:
|
|
#+begin_src matlab
|
|
load('huddle_test.mat', 't', 'y');
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
y = y - mean(y(1000:end));
|
|
#+end_src
|
|
|
|
** Time Domain Data
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
plot(t(1000:end), y(1000:end))
|
|
ylabel('Output Displacement [m]'); xlabel('Time [s]');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/huddle_test_time_domain.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:huddle_test_time_domain
|
|
#+caption: Measurement of the Mass displacement during Huddle Test
|
|
#+RESULTS:
|
|
[[file:figs/huddle_test_time_domain.png]]
|
|
|
|
** PSD of Measurement Noise
|
|
#+begin_src matlab
|
|
Ts = t(end)/(length(t)-1);
|
|
Fs = 1/Ts;
|
|
|
|
win = hanning(ceil(1*Fs));
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
[pxx, f] = pwelch(y(1000:end), win, [], [], Fs);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
plot(f, sqrt(pxx));
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
xlabel('Frequency [Hz]'); ylabel('ASD [$m/\sqrt{Hz}$]');
|
|
xlim([1, Fs/2]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/huddle_test_pdf.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:huddle_test_pdf
|
|
#+caption: Amplitude Spectral Density of the Displacement during Huddle Test
|
|
#+RESULTS:
|
|
[[file:figs/huddle_test_pdf.png]]
|
|
|
|
* TODO Identification of the dynamics from actuator to displacement
|
|
:PROPERTIES:
|
|
:header-args:matlab+: :tangle matlab/motion_identification.m
|
|
:header-args:matlab+: :comments org :mkdirp yes
|
|
:END:
|
|
<<sec:motion_identification>>
|
|
** Introduction :ignore:
|
|
|
|
- [ ] List of equipment
|
|
- [ ] Schematic
|
|
- [ ] Problem of matching between the models? (there is a factor 10)
|
|
|
|
E505 with gain of 10.
|
|
|
|
** 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/mat/');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :eval no
|
|
addpath('./mat/');
|
|
#+end_src
|
|
|
|
** Load Data
|
|
The data from the "noise test" and the identification test are loaded.
|
|
#+begin_src matlab
|
|
ht = load('huddle_test.mat', 't', 'u', 'y');
|
|
load('apa95ml_5kg_Amp_E505.mat', 't', 'um', 'y');
|
|
#+end_src
|
|
|
|
Any offset value is removed:
|
|
#+begin_src matlab
|
|
um = detrend(um, 0); % Input Voltage [V]
|
|
y = detrend(y , 0); % Mass displacement [m]
|
|
|
|
ht.u = detrend(ht.u, 0);
|
|
ht.y = detrend(ht.y, 0);
|
|
#+end_src
|
|
|
|
Now we add a factor 10 to take into account the gain of the voltage amplifier.
|
|
#+begin_src matlab
|
|
um = 10*um;
|
|
ht.u = 10*ht.u;
|
|
#+end_src
|
|
|
|
** Comparison of the PSD with Huddle Test
|
|
#+begin_src matlab
|
|
Ts = t(end)/(length(t)-1);
|
|
Fs = 1/Ts;
|
|
|
|
win = hanning(ceil(1*Fs));
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
[pxx, f] = pwelch(y, win, [], [], Fs);
|
|
[pht, ~] = pwelch(ht.y, win, [], [], Fs);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
hold on;
|
|
plot(f, sqrt(pxx), 'DisplayName', '5kg');
|
|
plot(f, sqrt(pht), 'DisplayName', 'Huddle Test');
|
|
hold off;
|
|
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
|
|
xlabel('Frequency [Hz]'); ylabel('ASD [$m/\sqrt{Hz}$]');
|
|
legend('location', 'northeast');
|
|
xlim([1, Fs/2]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/apa95ml_5kg_PI_pdf_comp_huddle.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:apa95ml_5kg_PI_pdf_comp_huddle
|
|
#+caption: Comparison of the ASD for the identification test and the huddle test
|
|
#+RESULTS:
|
|
[[file:figs/apa95ml_5kg_PI_pdf_comp_huddle.png]]
|
|
|
|
** Compute TF estimate and Coherence
|
|
#+begin_src matlab
|
|
[tf_est, f] = tfestimate(um, -y, win, [], [], 1/Ts);
|
|
[co_est, ~] = mscohere( um, -y, win, [], [], 1/Ts);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
|
|
hold on;
|
|
plot(f, co_est, 'k-')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Coherence'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
xlim([10, 5e3]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/apa95ml_5kg_PI_coh.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:apa95ml_5kg_PI_coh
|
|
#+caption: Coherence
|
|
#+RESULTS:
|
|
[[file:figs/apa95ml_5kg_PI_coh.png]]
|
|
|
|
Comparison with the FEM model
|
|
#+begin_src matlab
|
|
load('mat/fem_simscape_models.mat', 'Ghm');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(0, 4, 1000);
|
|
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(tf_est), 'DisplayName', 'Identification')
|
|
plot(freqs, abs(squeeze(freqresp(Ghm, freqs, 'Hz'))), 'DisplayName', 'FEM')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude [m/V]'); set(gca, 'XTickLabel', []);
|
|
legend('location', 'northeast')
|
|
hold off;
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(tf_est))
|
|
plot(freqs, 180/pi*angle(squeeze(freqresp(Ghm, freqs, 'Hz'))))
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
ylim([-180, 180]);
|
|
yticks(-180:90:180);
|
|
|
|
linkaxes([ax1,ax2], 'x');
|
|
xlim([10, 5e3]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/apa95ml_5kg_pi_comp_fem.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:apa95ml_5kg_pi_comp_fem
|
|
#+caption: Comparison of the identified transfer function and the one estimated from the FE model
|
|
#+RESULTS:
|
|
[[file:figs/apa95ml_5kg_pi_comp_fem.png]]
|
|
|
|
* Identification of the dynamics from actuator to force sensor
|
|
:PROPERTIES:
|
|
:header-args:matlab+: :tangle matlab/force_sensor_identification.m
|
|
:header-args:matlab+: :comments org :mkdirp yes
|
|
:END:
|
|
<<sec:force_sensor_identification>>
|
|
** Introduction :ignore:
|
|
Two measurements are performed:
|
|
- Speedgoat DAC => Voltage Amplifier (x20) => 1 Piezo Stack => ... => 2 Stacks as Force Sensor (parallel) => Speedgoat ADC
|
|
- Speedgoat DAC => Voltage Amplifier (x20) => 2 Piezo Stacks (parallel) => ... => 1 Stack as Force Sensor => Speedgoat ADC
|
|
|
|
The obtained dynamics from force actuator to force sensor are compare with the FEM model.
|
|
|
|
** 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/mat/');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :eval no
|
|
addpath('./mat/');
|
|
#+end_src
|
|
|
|
** Load Data :ignore:
|
|
The data are loaded:
|
|
#+begin_src matlab
|
|
load('apa95ml_5kg_2a_1s.mat', 't', 'u', 'v');
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
u = detrend(u, 0);
|
|
v = detrend(v, 0);
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
u = 20*u;
|
|
#+end_src
|
|
|
|
** Adjust gain :ignore:
|
|
Let's use the amplifier gain to obtain the true voltage applied to the actuator stack(s)
|
|
|
|
The parameters of the piezoelectric stacks are defined below:
|
|
#+begin_src matlab
|
|
d33 = 3e-10; % Strain constant [m/V]
|
|
n = 80; % Number of layers per stack
|
|
eT = 1.6e-8; % Permittivity under constant stress [F/m]
|
|
sD = 2e-11; % Elastic compliance under constant electric displacement [m2/N]
|
|
ka = 235e6; % Stack stiffness [N/m]
|
|
#+end_src
|
|
|
|
From the FEM, we construct the transfer function from DAC voltage to ADC voltage.
|
|
#+begin_src matlab
|
|
Gfem_aa_s = exp(-s/1e4)*20*(2*d33*n*ka)*(G(3,1)+G(3,2))*d33/(eT*sD*n);
|
|
Gfem_a_ss = exp(-s/1e4)*20*( d33*n*ka)*(G(3,1)+G(2,1))*d33/(eT*sD*n);
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
Gfem_aa_s = exp(-s/1e4)*20*(2*d33*n*ka)*Gfm*d33/(eT*sD*n);
|
|
Gfem_a_ss = exp(-s/1e4)*20*( d33*n*ka)*Gfm*d33/(eT*sD*n);
|
|
#+end_src
|
|
|
|
** Compute TF estimate and Coherence :ignore:
|
|
The transfer function from input voltage to output voltage are computed and shown in Figure [[fig:bode_plot_force_sensor_voltage_comp_fem]].
|
|
|
|
#+begin_src matlab
|
|
Ts = t(end)/(length(t)-1);
|
|
Fs = 1/Ts;
|
|
|
|
win = hann(ceil(10/Ts));
|
|
|
|
[tf_est, f] = tfestimate(u, v, win, [], [], 1/Ts);
|
|
[coh, ~] = mscohere( u, v, win, [], [], 1/Ts);
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
load('mat/fem_simscape_models.mat', 'Gfm');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(1, 4, 1000);
|
|
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(tf_est))
|
|
plot(freqs, abs(squeeze(freqresp(Gfm, freqs, 'Hz'))))
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
ylim([1e-3, 1e1]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(tf_est), ...
|
|
'DisplayName', 'Identification')
|
|
plot(freqs, 180/pi*angle(squeeze(freqresp(Gfm, freqs, 'Hz'))), ...
|
|
'DisplayName', 'FEM')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
ylim([-180, 180]);
|
|
yticks(-180:90:180);
|
|
legend('location', 'northeast')
|
|
|
|
linkaxes([ax1,ax2], 'x');
|
|
xlim([10, 5e3]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/bode_plot_force_sensor_voltage_comp_fem.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:bode_plot_force_sensor_voltage_comp_fem
|
|
#+caption: Comparison of the identified dynamics from voltage output to voltage input and the FEM
|
|
#+RESULTS:
|
|
[[file:figs/bode_plot_force_sensor_voltage_comp_fem.png]]
|
|
|
|
** System Identification
|
|
#+begin_src matlab
|
|
w_z = 2*pi*111; % Zeros frequency [rad/s]
|
|
w_p = 2*pi*255; % Pole frequency [rad/s]
|
|
xi_z = 0.05;
|
|
xi_p = 0.015;
|
|
G_inf = 0.1;
|
|
|
|
Gi = G_inf*(s^2 - 2*xi_z*w_z*s + w_z^2)/(s^2 + 2*xi_p*w_p*s + w_p^2);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
freqs = logspace(1, 4, 1000);
|
|
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2,1]);
|
|
hold on;
|
|
plot(f, abs(tf_est), '-')
|
|
plot(freqs, abs(squeeze(freqresp(Gi, freqs, 'Hz'))), '--')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude'); set(gca, 'XTickLabel', []);
|
|
hold off;
|
|
ylim([1e-3, 1e1]);
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(tf_est), '-', 'DisplayName', '2 Act - 1 Sen')
|
|
plot(freqs, 180/pi*angle(squeeze(freqresp(Gi, freqs, 'Hz'))), '--', 'DisplayName', '2 Act - 1 Sen, - FEM')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
ylim([-180, 180]);
|
|
yticks(-180:90:180);
|
|
legend('location', 'northeast')
|
|
|
|
linkaxes([ax1,ax2], 'x');
|
|
xlim([10, 5e3]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/iff_plant_identification_apa95ml.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:iff_plant_identification_apa95ml
|
|
#+caption: Identification of the IFF plant
|
|
#+RESULTS:
|
|
[[file:figs/iff_plant_identification_apa95ml.png]]
|
|
|
|
|
|
** Integral Force Feedback
|
|
#+begin_src matlab :exports none
|
|
gains = logspace(0, 5, 1000);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(real(pole(Gi)), imag(pole(Gi)), 'kx');
|
|
plot(real(tzero(Gi)), imag(tzero(Gi)), 'ko');
|
|
for i = 1:length(gains)
|
|
cl_poles = pole(feedback(Gi, (gains(i)/(s + 2*2*pi)*s/(s + 0.5*2*pi))));
|
|
plot(real(cl_poles), imag(cl_poles), 'k.');
|
|
end
|
|
ylim([0, 1800]);
|
|
xlim([-1600,200]);
|
|
xlabel('Real Part')
|
|
ylabel('Imaginary Part')
|
|
axis square
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/root_locus_iff_apa95ml_identification.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:root_locus_iff_apa95ml_identification
|
|
#+caption: Root Locus for IFF
|
|
#+RESULTS:
|
|
[[file:figs/root_locus_iff_apa95ml_identification.png]]
|
|
|
|
* Integral Force Feedback
|
|
:PROPERTIES:
|
|
:header-args:matlab+: :tangle matlab/integral_force_feedback.m
|
|
:header-args:matlab+: :comments org :mkdirp yes
|
|
:END:
|
|
<<sec:integral_force_feedback>>
|
|
|
|
** Introduction :ignore:
|
|
|
|
#+name: fig:test_bench_apa_schematic_iff
|
|
#+caption: Schematic of the test bench using IFF
|
|
[[file:figs/test_bench_apa_schematic_iff.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/mat/');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :eval no
|
|
addpath('./mat/');
|
|
#+end_src
|
|
|
|
** First tests with few gains
|
|
#+begin_src matlab
|
|
iff_g10 = load('apa95ml_iff_g10_res.mat', 'u', 't', 'y', 'v');
|
|
iff_g100 = load('apa95ml_iff_g100_res.mat', 'u', 't', 'y', 'v');
|
|
iff_of = load('apa95ml_iff_off_res.mat', 'u', 't', 'y', 'v');
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
Ts = 1e-4;
|
|
win = hann(ceil(10/Ts));
|
|
|
|
[tf_iff_g10, f] = tfestimate(iff_g10.u, iff_g10.y, win, [], [], 1/Ts);
|
|
[co_iff_g10, ~] = mscohere(iff_g10.u, iff_g10.y, win, [], [], 1/Ts);
|
|
|
|
[tf_iff_g100, ~] = tfestimate(iff_g100.u, iff_g100.y, win, [], [], 1/Ts);
|
|
[co_iff_g100, ~] = mscohere(iff_g100.u, iff_g100.y, win, [], [], 1/Ts);
|
|
|
|
[tf_iff_of, ~] = tfestimate(iff_of.u, iff_of.y, win, [], [], 1/Ts);
|
|
[co_iff_of, ~] = mscohere(iff_of.u, iff_of.y, win, [], [], 1/Ts);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
|
|
hold on;
|
|
plot(f, co_iff_of, '-', 'DisplayName', 'g=0')
|
|
plot(f, co_iff_g10, '-', 'DisplayName', 'g=10')
|
|
plot(f, co_iff_g100, '-', 'DisplayName', 'g=100')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Coherence'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
legend();
|
|
xlim([60, 600])
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/iff_first_test_coherence.pdf', 'width', 'wide', 'height', 'normal');
|
|
#+end_src
|
|
|
|
#+name: fig:iff_first_test_coherence
|
|
#+caption: Coherence
|
|
#+RESULTS:
|
|
[[file:figs/iff_first_test_coherence.png]]
|
|
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2, 1]);
|
|
hold on;
|
|
plot(f, abs(tf_iff_of), '-', 'DisplayName', 'g=0')
|
|
plot(f, abs(tf_iff_g10), '-', 'DisplayName', 'g=10')
|
|
plot(f, abs(tf_iff_g100), '-', 'DisplayName', 'g=100')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend();
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
plot(f, 180/pi*angle(-tf_iff_of), '-')
|
|
plot(f, 180/pi*angle(-tf_iff_g10), '-')
|
|
plot(f, 180/pi*angle(-tf_iff_g100), '-')
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
yticks(-360:90:360);
|
|
|
|
linkaxes([ax1,ax2], 'x');
|
|
xlim([60, 600]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/iff_first_test_bode_plot.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:iff_first_test_bode_plot
|
|
#+caption: Bode plot for different values of IFF gain
|
|
#+RESULTS:
|
|
[[file:figs/iff_first_test_bode_plot.png]]
|
|
|
|
** Second test with many Gains
|
|
#+begin_src matlab
|
|
load('apa95ml_iff_test.mat', 'results');
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
Ts = 1e-4;
|
|
win = hann(ceil(10/Ts));
|
|
#+end_src
|
|
|
|
#+begin_src matlab
|
|
tf_iff = {zeros(1, length(results))};
|
|
co_iff = {zeros(1, length(results))};
|
|
g_iff = [0, 1, 5, 10, 50, 100];
|
|
|
|
for i=1:length(results)
|
|
[tf_est, f] = tfestimate(results{i}.u, results{i}.y, win, [], [], 1/Ts);
|
|
[co_est, ~] = mscohere(results{i}.u, results{i}.y, win, [], [], 1/Ts);
|
|
|
|
tf_iff(i) = {tf_est};
|
|
co_iff(i) = {co_est};
|
|
end
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
|
|
hold on;
|
|
for i = 1:length(results)
|
|
plot(f, co_iff{i}, '-', 'DisplayName', sprintf('g = %0.f', g_iff(i)))
|
|
end
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Coherence'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
legend();
|
|
xlim([60, 600])
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2, 1]);
|
|
hold on;
|
|
for i = 1:length(results)
|
|
plot(f, abs(tf_iff{i}), '-', 'DisplayName', sprintf('g = %0.f', g_iff(i)))
|
|
end
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend();
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
for i = 1:length(results)
|
|
plot(f, 180/pi*angle(-tf_iff{i}), '-')
|
|
end
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
|
yticks(-360:90:360);
|
|
axis padded 'auto x'
|
|
|
|
linkaxes([ax1,ax2], 'x');
|
|
xlim([60, 600]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/iff_results_bode_plots.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:iff_results_bode_plots
|
|
#+caption:
|
|
#+RESULTS:
|
|
[[file:figs/iff_results_bode_plots.png]]
|
|
|
|
#+begin_src matlab
|
|
G_id = {zeros(1,length(results))};
|
|
|
|
f_start = 70; % [Hz]
|
|
f_end = 500; % [Hz]
|
|
|
|
for i = 1:length(results)
|
|
tf_id = tf_iff{i}(sum(f<f_start):length(f)-sum(f>f_end));
|
|
f_id = f(sum(f<f_start):length(f)-sum(f>f_end));
|
|
|
|
gfr = idfrd(tf_id, 2*pi*f_id, Ts);
|
|
G_id(i) = {procest(gfr,'P2UDZ')};
|
|
end
|
|
#+end_src
|
|
|
|
#+begin_src matlab :exports none
|
|
figure;
|
|
tiledlayout(3, 1, 'TileSpacing', 'None', 'Padding', 'None');
|
|
|
|
ax1 = nexttile([2, 1]);
|
|
hold on;
|
|
for i = 1:length(results)
|
|
set(gca,'ColorOrderIndex',i)
|
|
plot(f, abs(tf_iff{i}), '-', 'DisplayName', sprintf('g = %0.f', g_iff(i)))
|
|
set(gca,'ColorOrderIndex',i)
|
|
plot(f, abs(squeeze(freqresp(G_id{i}, f, 'Hz'))), '--', 'HandleVisibility', 'off')
|
|
end
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'log');
|
|
ylabel('Amplitude'); set(gca, 'XTickLabel',[]);
|
|
hold off;
|
|
legend();
|
|
|
|
ax2 = nexttile;
|
|
hold on;
|
|
for i = 1:length(results)
|
|
set(gca,'ColorOrderIndex',i)
|
|
plot(f, 180/pi*angle(tf_iff{i}), '-')
|
|
set(gca,'ColorOrderIndex',i)
|
|
plot(f, 180/pi*angle(squeeze(freqresp(G_id{i}, f, 'Hz'))), '--')
|
|
end
|
|
set(gca, 'Xscale', 'log'); set(gca, 'Yscale', 'lin');
|
|
ylabel('Phase'); xlabel('Frequency [Hz]');
|
|
hold off;
|
|
yticks(-360:90:360);
|
|
axis padded 'auto x'
|
|
|
|
linkaxes([ax1,ax2], 'x');
|
|
xlim([60, 600]);
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/iff_results_bode_plots_identification.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:iff_results_bode_plots_identification
|
|
#+caption:
|
|
#+RESULTS:
|
|
[[file:figs/iff_results_bode_plots_identification.png]]
|
|
|
|
#+begin_src matlab :exports none
|
|
w_z = 2*pi*111; % Zeros frequency [rad/s]
|
|
w_p = 2*pi*255; % Pole frequency [rad/s]
|
|
xi_z = 0.05;
|
|
xi_p = 0.015;
|
|
G_inf = 2;
|
|
|
|
Gi = G_inf*(s^2 - 2*xi_z*w_z*s + w_z^2)/(s^2 + 2*xi_p*w_p*s + w_p^2);
|
|
|
|
|
|
gains = logspace(0, 5, 1000);
|
|
|
|
figure;
|
|
hold on;
|
|
plot(real(pole(Gi)), imag(pole(Gi)), 'kx', 'HandleVisibility', 'off');
|
|
plot(real(tzero(Gi)), imag(tzero(Gi)), 'ko', 'HandleVisibility', 'off');
|
|
for i = 1:length(results)
|
|
set(gca,'ColorOrderIndex',i)
|
|
plot(real(pole(G_id{i})), imag(pole(G_id{i})), 'o', 'DisplayName', sprintf('g = %0.f', g_iff(i)));
|
|
end
|
|
for i = 1:length(gains)
|
|
cl_poles = pole(feedback(Gi, (gains(i)/(s + 2*pi*2))));
|
|
plot(real(cl_poles), imag(cl_poles), 'k.', 'HandleVisibility', 'off');
|
|
end
|
|
ylim([0, 1800]);
|
|
xlim([-1600,200]);
|
|
xlabel('Real Part')
|
|
ylabel('Imaginary Part')
|
|
axis square
|
|
legend('location', 'northwest');
|
|
#+end_src
|
|
|
|
#+begin_src matlab :tangle no :exports results :results file replace
|
|
exportFig('figs/iff_results_root_locus.pdf', 'width', 'wide', 'height', 'tall');
|
|
#+end_src
|
|
|
|
#+name: fig:iff_results_root_locus
|
|
#+caption:
|
|
#+RESULTS:
|
|
[[file:figs/iff_results_root_locus.png]]
|