#+end_export
* Introduction :ignore:
The goal of this test bench is to extract all the important parameters of the Amplified Piezoelectric Actuator APA300ML.
This include:
- Stroke
- Stiffness
- Hysteresis
- Gain from the applied voltage $V_a$ to the generated Force $F_a$
- Gain from the sensor stack strain $\delta L$ to the generated voltage $V_s$
- Dynamical behavior
#+name: fig:apa300ML
#+caption: Picture of the APA300ML
#+attr_latex: :width 0.8\linewidth
[[file:figs/apa300ML.png]]
* Model of an Amplified Piezoelectric Actuator and Sensor
Consider a schematic of the Amplified Piezoelectric Actuator in Figure [[fig:apa_model_schematic]].
#+name: fig:apa_model_schematic
#+caption: Amplified Piezoelectric Actuator Schematic
[[file:figs/apa_model_schematic.png]]
A voltage $V_a$ applied to the actuator stacks will induce an actuator force $F_a$:
\begin{equation}
F_a = g_a \cdot V_a
\end{equation}
A change of length $dl$ of the sensor stack will induce a voltage $V_s$:
\begin{equation}
V_s = g_s \cdot dl
\end{equation}
We wish here to experimental measure $g_a$ and $g_s$.
The block-diagram model of the piezoelectric actuator is then as shown in Figure [[fig:apa-model-simscape-schematic]].
#+begin_src latex :file apa-model-simscape-schematic.pdf
\begin{tikzpicture}
\node[block={2.0cm}{2.0cm}, align=center] (model) at (0,0){Simscape\\Model};
\node[block, left=1.0 of model] (ga){$g_a(s)$};
\node[block, right=1.0 of model] (gs){$g_s(s)$};
\draw[<-] (ga.west) -- node[midway, above]{$V_a$} node[midway, below]{$[V]$} ++(-1.0, 0);
\draw[->] (ga.east) --node[midway, above]{$F_a$} node[midway, below]{$[N]$} (model.west);
\draw[->] (model.east) --node[midway, above]{$dl$} node[midway, below]{$[m]$} (gs.west);
\draw[->] (gs.east) -- node[midway, above]{$V_s$} node[midway, below]{$[V]$} ++(1.0, 0);
\end{tikzpicture}
#+end_src
#+name: fig:apa-model-simscape-schematic
#+caption: Model of the APA with Simscape/Simulink
#+RESULTS:
[[file:figs/apa-model-simscape-schematic.png]]
* Geometrical Measurements
** Introduction :ignore:
#+name: fig:received_apa
#+caption: Received APA
#+attr_latex: :width 0.6\linewidth
[[file:figs/IMG_20210224_143500.jpg]]
** Matlab Init :noexport:ignore:
#+begin_src matlab :tangle no :exports none :results silent :noweb yes :var current_dir=(file-name-directory buffer-file-name)
<>
#+end_src
#+begin_src matlab :exports none :results silent :noweb yes
<>
#+end_src
#+begin_src matlab :tangle no
addpath('./matlab/mat/');
addpath('./matlab/');
#+end_src
#+begin_src matlab :eval no
addpath('./mat/');
#+end_src
** Measurement Setup
#+name: fig:flatness_meas_setup
#+caption: Measurement Setup
#+attr_latex: :width 0.6\linewidth
[[file:figs/IMG_20210224_143809.jpg]]
** Measurement Results
Height (Z) measurements:
#+begin_src matlab
apa1 = 1e-6*[0, -0.5 , 3.5 , 3.5 , 42 , 45.5, 52.5 , 46];
apa2 = 1e-6*[0, -2.5 , -3 , 0 , -1.5 , 1 , -2 , -4];
apa3 = 1e-6*[0, -1.5 , 15 , 17.5 , 6.5 , 6.5 , 21 , 23];
apa4 = 1e-6*[0, 6.5 , 14.5 , 9 , 16 , 22 , 29.5 , 21];
apa5 = 1e-6*[0, -12.5, 16.5 , 28.5 , -43 , -52 , -22.5, -13.5];
apa6 = 1e-6*[0, -8 , -2 , 5 , -57.5, -62 , -55.5, -52.5];
apa7 = 1e-6*[0, 19.5 , -8 , -29.5, 75 , 97.5, 70 , 48];
apa7b = 1e-6*[0, 9 , -18.5, -30 , 31 , 46.5, 16.5 , 7.5];
apa = {apa1, apa2, apa3, apa4, apa5, apa6, apa7b};
#+end_src
X/Y Positions of the 8 measurement points:
#+begin_src matlab
W = 20e-3; % Width [m]
L = 61e-3; % Length [m]
d = 1e-3; % Distance from border [m]
l = 15.5e-3; % [m]
pos = [[-L/2 + d; W/2 - d], [-L/2 + l - d; W/2 - d], [-L/2 + l - d; -W/2 + d], [-L/2 + d; -W/2 + d], [L/2 - l + d; W/2 - d], [L/2 - d; W/2 - d], [L/2 - d; -W/2 + d], [L/2 - l + d; -W/2 + d]];
#+end_src
#+begin_src matlab
apa_d = zeros(1, 7);
for i = 1:7
fun = @(x)max(abs(([pos; apa{i}]-[0;0;x(1)])'*([x(2:3);1]/norm([x(2:3);1]))));
x0 = [0;0;0];
[x, min_d] = fminsearch(fun,x0);
apa_d(i) = min_d;
end
#+end_src
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
data2orgtable(1e6*apa_d', {}, {'Flatness [um]'}, ' %.1f ');
#+end_src
#+name: tab:flatness_meas
#+caption: Estimated flatness
#+attr_latex: :environment tabularx :width \linewidth :align c
#+attr_latex: :center t :booktabs t :float t
#+RESULTS:
| Flatness [um] |
|---------------|
| 8.9 |
| 3.1 |
| 9.1 |
| 3.0 |
| 1.9 |
| 7.1 |
| 18.7 |
* Electrical Measurements
#+begin_note
The capacitance of the stacks is measure with the [[https://www.gwinstek.com/en-global/products/detail/LCR-800][LCR-800 Meter]] ([[file:doc/DS_LCR-800_Series_V2_E.pdf][doc]])
#+end_note
#+name: fig:LCR_meter
#+caption: LCR Meter used for the measurements
#+attr_latex: :width 0.9\linewidth
[[file:figs/IMG_20210312_120337.jpg]]
The excitation frequency is set to be 1kHz.
#+name: tab:apa300ml_capacitance
#+caption: Capacitance measured with the LCR meter. The excitation signal is a sinus at 1kHz
#+attr_latex: :environment tabularx :width 0.6\linewidth :align lcc
#+attr_latex: :center t :booktabs t :float t
| *APA Number* | *Sensor Stack* | *Actuator Stacks* |
|--------------+----------------+-------------------|
| 1 | 5.10 | 10.03 |
| 2 | 4.99 | 9.85 |
| 3 | 1.72 | 5.18 |
| 4 | 4.94 | 9.82 |
| 5 | 4.90 | 9.66 |
| 6 | 4.99 | 9.91 |
| 7 | 4.85 | 9.85 |
#+begin_warning
There is clearly a problem with APA300ML number 3
#+end_warning
* Stiffness measurement
** APA test
#+begin_src matlab
load('meas_stiff_apa_1_x.mat', 't', 'F', 'd');
#+end_src
#+begin_src matlab
figure;
plot(t, F)
#+end_src
#+begin_src matlab
%% Automatic Zero of the force
F = F - mean(F(t > 0.1 & t < 0.3));
%% Start measurement at t = 0.2 s
d = d(t > 0.2);
F = F(t > 0.2);
t = t(t > 0.2); t = t - t(1);
#+end_src
#+begin_src matlab
i_l_start = find(F > 0.3, 1, 'first');
[~, i_l_stop] = max(F);
#+end_src
#+begin_src matlab
F_l = F(i_l_start:i_l_stop);
d_l = d(i_l_start:i_l_stop);
#+end_src
#+begin_src matlab
fit_l = polyfit(F_l, d_l, 1);
% %% Reset displacement based on fit
% d = d - fit_l(2);
% fit_s(2) = fit_s(2) - fit_l(2);
% fit_l(2) = 0;
% %% Estimated Stroke
% F_max = fit_s(2)/(fit_l(1) - fit_s(1));
% d_max = fit_l(1)*F_max;
#+end_src
#+begin_src matlab
h^2/fit_l(1)
#+end_src
#+begin_src matlab
figure;
hold on;
plot(F,d,'k')
plot(F_l, d_l)
plot(F_l, F_l*fit_l(1) + fit_l(2), '--')
#+end_src
* Test-Bench Description
#+begin_note
Here are the documentation of the equipment used for this test bench:
- Voltage Amplifier: [[file:doc/PD200-V7-R1.pdf][PD200]]
- Amplified Piezoelectric Actuator: [[file:doc/APA300ML.pdf][APA300ML]]
- DAC/ADC: Speedgoat [[file:doc/IO131-OEM-Datasheet.pdf][IO313]]
- Encoder: [[file:doc/L-9517-9678-05-A_Data_sheet_VIONiC_series_en.pdf][Renishaw Vionic]] and used [[file:doc/L-9517-9862-01-C_Data_sheet_RKLC_EN.pdf][Ruler]]
- Interferometer: [[https://www.attocube.com/en/products/laser-displacement-sensor/displacement-measuring-interferometer][Attocube IDS3010]]
#+end_note
#+name: fig:test_bench_apa_alone
#+caption: Schematic of the Test Bench
[[file:figs/test_bench_apa_alone.png]]
* Measurement Procedure
** Introduction :ignore:
** Stroke Measurement
Using the PD200 amplifier, output a voltage:
\[ V_a = 65 + 85 \sin(2\pi \cdot t) \]
To have a quasi-static excitation between -20 and 150V.
As the gain of the PD200 amplifier is 20, the DAC output voltage should be:
\[ V_{dac}(t) = 3.25 + 4.25\sin(2\pi \cdot t) \]
Verify that the voltage offset of the PD200 is zero!
Measure the output vertical displacement $d$ using the interferometer.
Then, plot $d$ as a function of $V_a$, and perform a linear regression.
Conclude on the obtained stroke.
** Stiffness Measurement
Add some (known) weight $\delta m g$ on the suspended mass and measure the deflection $\delta d$.
This can be tested when the piezoelectric stacks are open-circuit.
As the stiffness will be around $k \approx 10^6 N/m$, an added mass of $m \approx 100g$ will induce a static deflection of $\approx 1\mu m$ which should be large enough for a precise measurement using the interferometer.
Then the obtained stiffness is:
\begin{equation}
k = \frac{\delta m g}{\delta d}
\end{equation}
** Hysteresis measurement
Supply a quasi static sinusoidal excitation $V_a$ at different voltages.
The offset should be 65V, and the sin amplitude can range from 1V up to 85V.
For each excitation amplitude, the vertical displacement $d$ of the mass is measured.
Then, $d$ is plotted as a function of $V_a$ for all the amplitudes.
#+name: fig:expected_hysteresis
#+caption: Expected Hysteresis (cite:poel10_explor_activ_hard_mount_vibrat)
#+attr_latex: :width 0.8\linewidth
[[file:figs/expected_hysteresis.png]]
** Piezoelectric Actuator Constant
Using the measurement test-bench, it is rather easy the determine the static gain between the applied voltage $V_a$ to the induced displacement $d$.
Use a quasi static (1Hz) excitation signal $V_a$ on the piezoelectric stack and measure the vertical displacement $d$.
Perform a linear regression to obtain:
\begin{equation}
d = g_{d/V_a} \cdot V_a
\end{equation}
Using the Simscape model of the APA, it is possible to determine the static gain between the actuator force $F_a$ to the induced displacement $d$:
\begin{equation}
d = g_{d/F_a} \cdot F_a
\end{equation}
From the two gains, it is then easy to determine $g_a$:
\begin{equation}
g_a = \frac{F_a}{V_a} = \frac{F_a}{d} \cdot \frac{d}{V_a} = \frac{g_{d/V_a}}{g_{d/F_a}}
\end{equation}
** Piezoelectric Sensor Constant
From a quasi static excitation of the piezoelectric stack, measure the gain from $V_a$ to $V_s$:
\begin{equation}
V_s = g_{V_s/V_a} V_a
\end{equation}
Note here that there is an high pass filter formed by the piezo capacitor and parallel resistor.
The excitation frequency should then be in between the cut-off frequency of this high pass filter and the first resonance.
Alternatively, the gain can be computed from the dynamical identification and taking the gain at the wanted frequency.
Using the simscape model, compute the static gain from the actuator force $F_a$ to the strain of the sensor stack $dl$:
\begin{equation}
dl = g_{dl/F_a} F_a
\end{equation}
Then, the static gain from the sensor stack strain $dl$ to the general voltage $V_s$ is:
\begin{equation}
g_s = \frac{V_s}{dl} = \frac{V_s}{V_a} \cdot \frac{V_a}{F_a} \cdot \frac{F_a}{dl} = \frac{g_{V_s/V_a}}{g_a \cdot g_{dl/F_a}}
\end{equation}
Alternatively, we could impose an external force to add strain in the APA that should be equally present in all the 3 stacks and equal to 1/5 of the vertical strain.
This external force can be some weight added, or a piezo in parallel.
** Capacitance Measurement
Measure the capacitance of the 3 stacks individually using a precise multi-meter.
** Dynamical Behavior
Perform a system identification from $V_a$ to the measured displacement $d$ by the interferometer and by the encoder, and to the generated voltage $V_s$.
This can be performed using different excitation signals.
This can also be performed with and without the encoder fixed to the APA.
** Compare the results obtained for all 7 APA300ML
Compare all the obtained parameters for all the test APA.
* Measurement Results
* TODO Compare with the FEM/Simscape Model :noexport:
:PROPERTIES:
:header-args:matlab+: :tangle matlab/APA300ML.m
:END:
** Introduction :ignore:
In this section, the Amplified Piezoelectric Actuator APA300ML ([[file:doc/APA300ML.pdf][doc]]) is modeled using a Finite Element Software.
Then a /super element/ is exported and imported in Simscape where its dynamic is studied.
A 3D view of the Amplified Piezoelectric Actuator (APA300ML) is shown in Figure [[fig:apa300ml_ansys]].
The remote point used are also shown in this figure.
#+name: fig:apa300ml_ansys
#+caption: Ansys FEM of the APA300ML
[[file:figs/apa300ml_ansys.jpg]]
** Matlab Init :noexport:ignore:
#+begin_src matlab :tangle no :exports none :results silent :noweb yes :var current_dir=(file-name-directory buffer-file-name)
<>
#+end_src
#+begin_src matlab :exports none :results silent :noweb yes
<>
#+end_src
#+begin_src matlab :tangle no
addpath('matlab/');
addpath('matlab/APA300ML/');
#+end_src
#+begin_src matlab :eval no
addpath('APA300ML/');
#+end_src
#+begin_src matlab
open('APA300ML.slx');
#+end_src
** Import Mass Matrix, Stiffness Matrix, and Interface Nodes Coordinates
We first extract the stiffness and mass matrices.
#+begin_src matlab
K = readmatrix('APA300ML_mat_K.CSV');
M = readmatrix('APA300ML_mat_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:
| 200000000.0 | 30000.0 | -20000.0 | -70.0 | 300000.0 | 40.0 | 10000000.0 | 10000.0 | -6000.0 | 30.0 |
| 30000.0 | 30000000.0 | 2000.0 | -200000.0 | 60.0 | -10.0 | 4000.0 | 2000000.0 | -500.0 | 9000.0 |
| -20000.0 | 2000.0 | 7000000.0 | -10.0 | -30.0 | 10.0 | 6000.0 | 900.0 | -500000.0 | 3 |
| -70.0 | -200000.0 | -10.0 | 1000.0 | -0.1 | 0.08 | -20.0 | -9000.0 | 3 | -30.0 |
| 300000.0 | 60.0 | -30.0 | -0.1 | 900.0 | 0.1 | 30000.0 | 20.0 | -10.0 | 0.06 |
| 40.0 | -10.0 | 10.0 | 0.08 | 0.1 | 10000.0 | 20.0 | 9 | -5 | 0.03 |
| 10000000.0 | 4000.0 | 6000.0 | -20.0 | 30000.0 | 20.0 | 200000000.0 | 10000.0 | 9000.0 | 50.0 |
| 10000.0 | 2000000.0 | 900.0 | -9000.0 | 20.0 | 9 | 10000.0 | 30000000.0 | -500.0 | 200000.0 |
| -6000.0 | -500.0 | -500000.0 | 3 | -10.0 | -5 | 9000.0 | -500.0 | 7000000.0 | -2 |
| 30.0 | 9000.0 | 3 | -30.0 | 0.06 | 0.03 | 50.0 | 200000.0 | -2 | 1000.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.01 | -2e-06 | 1e-06 | 6e-09 | 5e-05 | -5e-09 | -0.0005 | -7e-07 | 6e-07 | -3e-09 |
| -2e-06 | 0.01 | 8e-07 | -2e-05 | -8e-09 | 2e-09 | -9e-07 | -0.0002 | 1e-08 | -9e-07 |
| 1e-06 | 8e-07 | 0.009 | 5e-10 | 1e-09 | -1e-09 | -5e-07 | 3e-08 | 6e-05 | 1e-10 |
| 6e-09 | -2e-05 | 5e-10 | 3e-07 | 2e-11 | -3e-12 | 3e-09 | 9e-07 | -4e-10 | 3e-09 |
| 5e-05 | -8e-09 | 1e-09 | 2e-11 | 6e-07 | -4e-11 | -1e-06 | -2e-09 | 1e-09 | -8e-12 |
| -5e-09 | 2e-09 | -1e-09 | -3e-12 | -4e-11 | 1e-07 | -2e-09 | -1e-09 | -4e-10 | -5e-12 |
| -0.0005 | -9e-07 | -5e-07 | 3e-09 | -1e-06 | -2e-09 | 0.01 | 1e-07 | -3e-07 | -2e-08 |
| -7e-07 | -0.0002 | 3e-08 | 9e-07 | -2e-09 | -1e-09 | 1e-07 | 0.01 | -4e-07 | 2e-05 |
| 6e-07 | 1e-08 | 6e-05 | -4e-10 | 1e-09 | -4e-10 | -3e-07 | -4e-07 | 0.009 | -2e-10 |
| -3e-09 | -9e-07 | 1e-10 | 3e-09 | -8e-12 | -5e-12 | -2e-08 | 2e-05 | -2e-10 | 3e-07 |
Then, we extract the coordinates of the interface nodes.
#+begin_src matlab
[int_xyz, int_i, n_xyz, n_i, nodes] = extractNodes('APA300ML_out_nodes_3D.txt');
#+end_src
#+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
#+caption: Coordinates of the interface nodes
#+RESULTS:
| Node i | Node Number | x [m] | y [m] | z [m] |
|--------+-------------+---------+-------+--------|
| 1.0 | 697783.0 | 0.0 | 0.0 | -0.015 |
| 2.0 | 697784.0 | 0.0 | 0.0 | 0.015 |
| 3.0 | 697785.0 | -0.0325 | 0.0 | 0.0 |
| 4.0 | 697786.0 | -0.0125 | 0.0 | 0.0 |
| 5.0 | 697787.0 | -0.0075 | 0.0 | 0.0 |
| 6.0 | 697788.0 | 0.0125 | 0.0 | 0.0 |
| 7.0 | 697789.0 | 0.0325 | 0.0 | 0.0 |
#+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
#+caption: Some extracted parameters of the FEM
#+RESULTS:
| Total number of Nodes | 7 |
| Number of interface Nodes | 7 |
| Number of Modes | 120 |
| Size of M and K matrices | 162 |
Using =K=, =M= and =int_xyz=, we can now use the =Reduced Order Flexible Solid= simscape block.
** Piezoelectric parameters
#+begin_src matlab
Ga = 1; % [N/V]
Gs = 1; % [V/m]
#+end_src
#+begin_src matlab
m = 0.1; % [kg]
#+end_src
** Simscape Model
The flexible element is imported using the =Reduced Order Flexible Solid= simscape block.
Let's say we use two stacks as a force sensor and one stack as an actuator:
- A =Relative Motion Sensor= block is added between the nodes A and C
- An =Internal Force= block is added between the remote points E and B
The interface nodes are shown in Figure [[fig:apa300ml_ansys]].
One mass is fixed at one end of the piezo-electric stack actuator (remove point F), the other end is fixed to the world frame (remote point G).
** Identification of the APA Characteristics
*** Stiffness
#+begin_src matlab :exports none
m = 0.0001;
#+end_src
The transfer function from vertical external force to the relative vertical displacement is identified.
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML';
%% Input/Output definition
clear io; io_i = 1;
io(io_i) = linio([mdl, '/Fd'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/z'], 1, 'openoutput'); io_i = io_i + 1;
G = linearize(mdl, io);
#+end_src
The inverse of its DC gain is the axial stiffness of the APA:
#+begin_src matlab :results replace value
1e-6/dcgain(G) % [N/um]
#+end_src
#+RESULTS:
: 1.753
The specified stiffness in the datasheet is $k = 1.8\, [N/\mu m]$.
*** Resonance Frequency
The resonance frequency is specified to be between 650Hz and 840Hz.
This is also the case for the FEM model (Figure [[fig:apa300ml_resonance]]).
#+begin_src matlab :exports none
freqs = logspace(2, 4, 5000);
figure;
hold on;
plot(freqs, abs(squeeze(freqresp(G, freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]'); ylabel('Amplitude');
hold off;
#+end_src
#+begin_src matlab :tangle no :exports results :results file replace
exportFig('figs/apa300ml_resonance.pdf', 'width', 'wide', 'height', 'normal');
#+end_src
#+name: fig:apa300ml_resonance
#+caption: First resonance is around 800Hz
#+RESULTS:
[[file:figs/apa300ml_resonance.png]]
*** Amplification factor
The amplification factor is the ratio of the vertical displacement to the stack displacement.
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML';
%% 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, '/z'], 1, 'openoutput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/d'], 1, 'openoutput'); io_i = io_i + 1;
G = linearize(mdl, io);
#+end_src
The ratio of the two displacement is computed from the FEM model.
#+begin_src matlab :results replace value
abs(dcgain(G(1,1))./dcgain(G(2,1)))
#+end_src
#+RESULTS:
: 5.0749
This is actually correct and approximately corresponds to the ratio of the piezo height and length:
#+begin_src matlab :results replace value
75/15
#+end_src
#+RESULTS:
: 5
*** Stroke
Estimation of the actuator stroke:
\[ \Delta H = A n \Delta L \]
with:
- $\Delta H$ Axial Stroke of the APA
- $A$ Amplification factor (5 for the APA300ML)
- $n$ Number of stack used
- $\Delta L$ Stroke of the stack (0.1% of its length)
#+begin_src matlab :results replace value
1e6 * 5 * 3 * 20e-3 * 0.1e-2
#+end_src
#+RESULTS:
: 300
This is exactly the specified stroke in the data-sheet.
*** TODO Stroke BIS
- [ ] Identified the stroke form the transfer function from V to z
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML';
%% Input/Output definition
clear io; io_i = 1;
io(io_i) = linio([mdl, '/V'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/d'], 1, 'openoutput'); io_i = io_i + 1;
G = linearize(mdl, io);
1e6*170*abs(dcgain(G))
#+end_src
** Identification of the Dynamics from actuator to replace displacement
We first set the mass to be approximately zero.
#+begin_src matlab :exports none
m = 0.01;
#+end_src
The dynamics is identified from the applied force to the measured relative displacement.
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML';
%% 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, '/z'], 1, 'openoutput'); io_i = io_i + 1;
Gh = -linearize(mdl, io);
#+end_src
The same dynamics is identified for a payload mass of 10Kg.
#+begin_src matlab
m = 10;
#+end_src
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML';
%% 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, '/z'], 1, 'openoutput'); io_i = io_i + 1;
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(Gh, freqs, 'Hz'))), '-');
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(Gh, freqs, 'Hz')))), '-', ...
'DisplayName', '$m = 0kg$');
plot(freqs, 180/pi*unwrap(angle(squeeze(freqresp(Ghm, freqs, 'Hz')))), '-', ...
'DisplayName', '$m = 10kg$');
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)]);
legend('location', 'southwest');
#+end_src
#+begin_src matlab :tangle no :exports results :results file replace
exportFig('figs/apa300ml_plant_dynamics.pdf', 'width', 'wide', 'height', 'tall');
#+end_src
#+name: fig:apa300ml_plant_dynamics
#+caption: Transfer function from forces applied by the stack to the axial displacement of the APA
#+RESULTS:
[[file:figs/apa300ml_plant_dynamics.png]]
The root locus corresponding to Direct Velocity Feedback with a mass of 10kg is shown in Figure [[fig:apa300ml_dvf_root_locus]].
#+begin_src matlab :exports none
figure;
gains = logspace(0, 5, 500);
hold on;
plot(real(pole(Ghm)), imag(pole(G)), 'kx');
plot(real(tzero(Ghm)), imag(tzero(G)), 'ko');
for k = 1:length(gains)
cl_poles = pole(feedback(Ghm, gains(k)*s));
plot(real(cl_poles), imag(cl_poles), 'k.');
end
hold off;
axis square;
xlim([-500, 10]); ylim([0, 510]);
xlabel('Real Part'); ylabel('Imaginary Part');
#+end_src
#+begin_src matlab :tangle no :exports results :results file replace
exportFig('figs/apa300ml_dvf_root_locus.pdf', 'width', 'wide', 'height', 'tall');
#+end_src
#+name: fig:apa300ml_dvf_root_locus
#+caption: Root Locus for Direct Velocity Feedback
#+RESULTS:
[[file:figs/apa300ml_dvf_root_locus.png]]
** Identification of the Dynamics from actuator to force sensor
Let's use 2 stacks as a force sensor and 1 stack as force actuator.
The transfer function from actuator voltage to sensor voltage is identified and shown in Figure [[fig:apa300ml_iff_plant]].
#+begin_src matlab :exports none
m = 10;
#+end_src
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML';
%% Input/Output definition
clear io; io_i = 1;
io(io_i) = linio([mdl, '/Va'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/Vs'], 1, 'openoutput'); io_i = io_i + 1;
Giff = -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(Giff, 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(Giff, freqs, 'Hz')))), '-');
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
yticks(-360:90:360);
ylim([-180 180]);
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/apa300ml_iff_plant.pdf', 'width', 'wide', 'height', 'tall');
#+end_src
#+name: fig:apa300ml_iff_plant
#+caption: Transfer function from actuator to force sensor
#+RESULTS:
[[file:figs/apa300ml_iff_plant.png]]
For root locus corresponding to IFF is shown in Figure [[fig:apa300ml_iff_root_locus]].
#+begin_src matlab :exports none
figure;
gains = logspace(0, 5, 500);
hold on;
plot(real(pole(Giff)), imag(pole(Giff)), 'kx');
plot(real(tzero(Giff)), imag(tzero(Giff)), 'ko');
for k = 1:length(gains)
cl_poles = pole(feedback(Giff, gains(k)/s));
plot(real(cl_poles), imag(cl_poles), 'k.');
end
hold off;
axis square;
xlim([-500, 10]); ylim([0, 510]);
xlabel('Real Part'); ylabel('Imaginary Part');
#+end_src
#+begin_src matlab :tangle no :exports results :results file replace
exportFig('figs/apa300ml_iff_root_locus.pdf', 'width', 'wide', 'height', 'tall');
#+end_src
#+name: fig:apa300ml_iff_root_locus
#+caption: Root Locus for IFF
#+RESULTS:
[[file:figs/apa300ml_iff_root_locus.png]]
** Identification for a simpler model
The goal in this section is to identify the parameters of a simple APA model from the FEM.
This can be useful is a lower order model is to be used for simulations.
The presented model is based on cite:souleille18_concep_activ_mount_space_applic.
The model represents the Amplified Piezo Actuator (APA) from Cedrat-Technologies (Figure [[fig:souleille18_model_piezo]]).
The parameters are shown in the table below.
#+name: fig:souleille18_model_piezo
#+caption: Picture of an APA100M from Cedrat Technologies. Simplified model of a one DoF payload mounted on such isolator
[[file:./figs/souleille18_model_piezo.png]]
#+caption:Parameters used for the model of the APA 100M
| | Meaning |
|-------+----------------------------------------------------------------|
| $k_e$ | Stiffness used to adjust the pole of the isolator |
| $k_1$ | Stiffness of the metallic suspension when the stack is removed |
| $k_a$ | Stiffness of the actuator |
| $c_1$ | Added viscous damping |
The goal is to determine $k_e$, $k_a$ and $k_1$ so that the simplified model fits the FEM model.
\[ \alpha = \frac{x_1}{f}(\omega=0) = \frac{\frac{k_e}{k_e + k_a}}{k_1 + \frac{k_e k_a}{k_e + k_a}} \]
\[ \beta = \frac{x_1}{F}(\omega=0) = \frac{1}{k_1 + \frac{k_e k_a}{k_e + k_a}} \]
If we can fix $k_a$, we can determine $k_e$ and $k_1$ with:
\[ k_e = \frac{k_a}{\frac{\beta}{\alpha} - 1} \]
\[ k_1 = \frac{1}{\beta} - \frac{k_e k_a}{k_e + k_a} \]
#+begin_src matlab :exports none
m = 10;
#+end_src
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML';
%% Input/Output definition
clear io; io_i = 1;
io(io_i) = linio([mdl, '/Fd'], 1, 'openinput'); io_i = io_i + 1; % External Vertical Force [N]
io(io_i) = linio([mdl, '/w'], 1, 'openinput'); io_i = io_i + 1; % Base Motion [m]
io(io_i) = linio([mdl, '/Fa'], 1, 'openinput'); io_i = io_i + 1; % Actuator Force [N]
io(io_i) = linio([mdl, '/z'], 1, 'openoutput'); io_i = io_i + 1; % Vertical Displacement [m]
io(io_i) = linio([mdl, '/Vs'], 1, 'openoutput'); io_i = io_i + 1; % Force Sensor [V]
io(io_i) = linio([mdl, '/d'], 1, 'openoutput'); io_i = io_i + 1; % Stack Displacement [m]
G = linearize(mdl, io);
G.InputName = {'Fd', 'w', 'Fa'};
G.OutputName = {'y', 'Fs', 'd'};
#+end_src
From the identified dynamics, compute $\alpha$ and $\beta$
#+begin_src matlab
alpha = abs(dcgain(G('y', 'Fa')));
beta = abs(dcgain(G('y', 'Fd')));
#+end_src
$k_a$ is estimated using the following formula:
#+begin_src matlab
ka = 0.8/abs(dcgain(G('y', 'Fa')));
#+end_src
The factor can be adjusted to better match the curves.
Then $k_e$ and $k_1$ are computed.
#+begin_src matlab
ke = ka/(beta/alpha - 1);
k1 = 1/beta - ke*ka/(ke + ka);
#+end_src
#+begin_src matlab :exports results :results value table replace :tangle no :post addhdr(*this*)
data2orgtable(1e-6*[ka; ke; k1], {'ka', 'ke', 'k1'}, {'Value [N/um]'}, ' %.1f ');
#+end_src
#+RESULTS:
| | Value [N/um] |
|----+--------------|
| ka | 40.5 |
| ke | 1.5 |
| k1 | 0.4 |
The damping in the system is adjusted to match the FEM model if necessary.
#+begin_src matlab
c1 = 1e2;
#+end_src
The analytical model of the simpler system is defined below:
#+begin_src matlab
Ga = 1/(m*s^2 + k1 + c1*s + ke*ka/(ke + ka)) * ...
[ 1 , k1 + c1*s + ke*ka/(ke + ka) , ke/(ke + ka) ;
-ke*ka/(ke + ka), ke*ka/(ke + ka)*m*s^2 , -ke/(ke + ka)*(m*s^2 + c1*s + k1)];
Ga.InputName = {'Fd', 'w', 'Fa'};
Ga.OutputName = {'y', 'Fs'};
#+end_src
And the DC gain is adjusted for the force sensor:
#+begin_src matlab
F_gain = dcgain(G('Fs', 'Fd'))/dcgain(Ga('Fs', 'Fd'));
#+end_src
The dynamics of the FEM model and the simpler model are compared in Figure [[fig:apa300ml_comp_simpler_model]].
#+begin_src matlab :exports none
freqs = logspace(0, 5, 1000);
figure;
tiledlayout(2, 3, 'TileSpacing', 'None', 'Padding', 'None');
ax1 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'y', 'w'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Ga('y', 'w'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]);
ylabel('$x_1/w$ [m/m]');
ylim([1e-6, 1e2]);
ax2 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'y', 'Fa'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Ga('y', 'Fa'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]);
ylabel('$x_1/f$ [m/N]');
ylim([1e-14, 1e-6]);
ax3 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'y', 'Fd'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Ga('y', 'Fd'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]);
ylabel('$x_1/F$ [m/N]');
ylim([1e-14, 1e-4]);
ax4 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'Fs', 'w'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(F_gain*Ga('Fs', 'w'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]');
ylabel('$F_s/w$ [m/m]');
ylim([1e2, 1e8]);
ax5 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'Fs', 'Fa'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(F_gain*Ga('Fs', 'Fa'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]');
ylabel('$F_s/f$ [m/N]');
ylim([1e-4, 1e1]);
ax6 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'Fs', 'Fd'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(F_gain*Ga('Fs', 'Fd'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]');
ylabel('$F_s/F$ [m/N]');
ylim([1e-7, 1e2]);
linkaxes([ax1,ax2,ax3,ax4,ax5,ax6],'x');
#+end_src
#+begin_src matlab :tangle no :exports results :results file replace
exportFig('figs/apa300ml_comp_simpler_model.pdf', 'width', 'full', 'height', 'full');
#+end_src
#+name: fig:apa300ml_comp_simpler_model
#+caption: Comparison of the Dynamics between the FEM model and the simplified one
#+RESULTS:
[[file:figs/apa300ml_comp_simpler_model.png]]
The simplified model has also been implemented in Simscape.
The dynamics of the Simscape simplified model is identified and compared with the FEM one in Figure [[fig:apa300ml_comp_simpler_simscape]].
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML_simplified';
%% Input/Output definition
clear io; io_i = 1;
io(io_i) = linio([mdl, '/Fd'], 1, 'openinput'); io_i = io_i + 1; % External Vertical Force [N]
io(io_i) = linio([mdl, '/w'], 1, 'openinput'); io_i = io_i + 1; % Base Motion [m]
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]
io(io_i) = linio([mdl, '/Fs'], 1, 'openoutput'); io_i = io_i + 1; % Force Sensor [V]
Gs = linearize(mdl, io);
Gs.InputName = {'Fd', 'w', 'Fa'};
Gs.OutputName = {'y', 'Fs'};
#+end_src
#+begin_src matlab :exports none
freqs = logspace(0, 5, 1000);
figure;
tiledlayout(2, 3, 'TileSpacing', 'None', 'Padding', 'None');
ax1 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'y', 'w'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Gs('y', 'w'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]);
ylabel('$x_1/w$ [m/m]');
ylim([1e-6, 1e2]);
ax2 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'y', 'Fa'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Gs('y', 'Fa'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]);
ylabel('$x_1/f$ [m/N]');
ylim([1e-14, 1e-6]);
ax3 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'y', 'Fd'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Gs('y', 'Fd'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]);
ylabel('$x_1/F$ [m/N]');
ylim([1e-14, 1e-4]);
ax4 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'Fs', 'w'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(F_gain*Gs('Fs', 'w'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]');
ylabel('$F_s/w$ [m/m]');
ylim([1e2, 1e8]);
ax5 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'Fs', 'Fa'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(F_gain*Gs('Fs', 'Fa'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]');
ylabel('$F_s/f$ [m/N]');
ylim([1e-4, 1e1]);
ax6 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G( 'Fs', 'Fd'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(F_gain*Gs('Fs', 'Fd'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]');
ylabel('$F_s/F$ [m/N]');
ylim([1e-7, 1e2]);
linkaxes([ax1,ax2,ax3,ax4,ax5,ax6],'x');
#+end_src
#+begin_src matlab :tangle no :exports results :results file replace
exportFig('figs/apa300ml_comp_simpler_simscape.pdf', 'width', 'full', 'height', 'full');
#+end_src
#+name: fig:apa300ml_comp_simpler_simscape
#+caption: Comparison of the Dynamics between the FEM model and the simplified simscape model
#+RESULTS:
[[file:figs/apa300ml_comp_simpler_simscape.png]]
** Integral Force Feedback
In this section, Integral Force Feedback control architecture is applied on the APA300ML.
First, the plant (dynamics from voltage actuator to voltage sensor is identified).
#+begin_src matlab :exports none
Kiff = tf(0);
#+end_src
The payload mass is set to 10kg.
#+begin_src matlab
m = 10;
#+end_src
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'APA300ML_IFF';
%% Input/Output definition
clear io; io_i = 1;
io(io_i) = linio([mdl, '/w'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/F'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/Fd'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/z'], 1, 'openoutput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/APA300ML'], 1, 'openoutput'); io_i = io_i + 1;
G_ol = linearize(mdl, io);
G_ol.InputName = {'w', 'f', 'F'};
G_ol.OutputName = {'x1', 'Fs'};
G = G_ol({'Fs'}, {'f'});
#+end_src
The obtained dynamics is shown in Figure [[fig:piezo_amplified_iff_plant]].
#+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(G, 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(G, 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/piezo_amplified_iff_plant.pdf', 'width', 'wide', 'height', 'tall');
#+end_src
#+name: fig:piezo_amplified_iff_plant
#+caption: IFF Plant
#+RESULTS:
[[file:figs/piezo_amplified_iff_plant.png]]
The controller is defined below and the loop gain is shown in Figure [[fig:piezo_amplified_iff_loop_gain]].
#+begin_src matlab
Kiff = -1e3/s;
#+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(G*Kiff, 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(G*Kiff, freqs, 'Hz')))));
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'lin');
yticks(-360:90:360);
ylim([-180 180]);
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/piezo_amplified_iff_loop_gain.pdf', 'width', 'wide', 'height', 'tall');
#+end_src
#+name: fig:piezo_amplified_iff_loop_gain
#+caption: IFF Loop Gain
#+RESULTS:
[[file:figs/piezo_amplified_iff_loop_gain.png]]
Now the closed-loop system is identified again and compare with the open loop system in Figure [[fig:piezo_amplified_iff_comp]].
It is the expected behavior as shown in the Figure [[fig:souleille18_results]] (from cite:souleille18_concep_activ_mount_space_applic).
#+begin_src matlab :exports none
clear io; io_i = 1;
io(io_i) = linio([mdl, '/w'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/F'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/Fd'], 1, 'openinput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/z'], 1, 'openoutput'); io_i = io_i + 1;
io(io_i) = linio([mdl, '/APA300ML'], 1, 'output'); io_i = io_i + 1;
Giff = linearize(mdl, io);
Giff.InputName = {'w', 'f', 'F'};
Giff.OutputName = {'x1', 'Fs'};
#+end_src
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
tiledlayout(2, 3, 'TileSpacing', 'None', 'Padding', 'None');
ax1 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G_ol('x1', 'w'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Giff('x1', 'w'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]); ylabel('$x_1/w$ [m/m]')
ax2 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G_ol('x1', 'f'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Giff('x1', 'f'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]); ylabel('$x_1/f$ [m/N]');
ax3 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G_ol('x1', 'F'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Giff('x1', 'F'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
set(gca, 'XTickLabel',[]); ylabel('$x_1/F$ [m/N]');
ax4 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G_ol('Fs', 'w'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Giff('Fs', 'w'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]'); ylabel('$F_s/w$ [N/m]');
ax5 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G_ol('Fs', 'f'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Giff('Fs', 'f'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]'); ylabel('$F_s/f$ [N/N]');
ax6 = nexttile;
hold on;
plot(freqs, abs(squeeze(freqresp(G_ol('Fs', 'F'), freqs, 'Hz'))));
plot(freqs, abs(squeeze(freqresp(Giff('Fs', 'F'), freqs, 'Hz'))));
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
xlabel('Frequency [Hz]'); ylabel('$F_s/F$ [N/N]');
linkaxes([ax1,ax2,ax3,ax4,ax5,ax6],'x');
#+end_src
#+begin_src matlab :tangle no :exports results :results file replace
exportFig('figs/piezo_amplified_iff_comp.pdf', 'width', 'full', 'height', 'full');
#+end_src
#+name: fig:piezo_amplified_iff_comp
#+caption: OL and CL transfer functions
#+RESULTS:
[[file:figs/piezo_amplified_iff_comp.png]]
#+name: fig:souleille18_results
#+caption: Results obtained in cite:souleille18_concep_activ_mount_space_applic
[[file:figs/souleille18_results.png]]
* Bibliography :ignore:
#+latex: \printbibliography