tdehaeze/nass-simscape
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity
48b2d36418
nass-simscape/org/optimal_stiffness_disturbances.org

794 lines
28 KiB
Org Mode
Raw Blame History

#+TITLE: Determination of the optimal nano-hexapod's stiffness for reducing the effect of disturbances
:DRAWER:
#+STARTUP: overview
#+LANGUAGE: en
#+EMAIL: dehaeze.thomas@gmail.com
#+AUTHOR: Dehaeze Thomas
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="./css/htmlize.css"/>
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="./css/readtheorg.css"/>
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="./css/zenburn.css"/>
#+HTML_HEAD: <script type="text/javascript" src="./js/jquery.min.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="./js/bootstrap.min.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="./js/jquery.stickytableheaders.min.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="./js/readtheorg.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 ../matlab/optimal_stiffness.m
#+PROPERTY: header-args:matlab+ :mkdirp yes
#+PROPERTY: header-args:shell :eval no-export
#+PROPERTY: header-args:latex :headers '("\\usepackage{tikz}" "\\usepackage{import}" "\\import{$HOME/Cloud/thesis/latex/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+ :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:
* Introduction :ignore:
In this document is studied how the stiffness of the nano-hexapod will impact the effect of disturbances on the position error of the sample.
It is divided in the following sections:
- Section [[sec:psd_disturbances]]: the disturbances are listed and their Power Spectral Densities (PSD) are shown
- Section [[sec:effect_disturbances]]: the transfer functions from disturbances to the position error of the sample are computed for a wide range of nano-hexapod stiffnesses
- Section [[sec:granite_stiffness]]:
- Section [[sec:open_loop_budget_error]]: from both the PSD of the disturbances and the transfer function from disturbances to sample's position errors, we compute the resulting PSD and Cumulative Amplitude Spectrum (CAS)
- Section [[sec:closed_loop_budget_error]]: from a simplistic model is computed the required control bandwidth to reduce the position error to acceptable values
* Disturbances
<<sec:psd_disturbances>>
** Introduction :ignore:
The main disturbances considered here are:
- $D_w$: Ground displacement in the $x$, $y$ and $z$ directions
- $F_{ty}$: Forces applied by the Translation stage in the $x$ and $z$ directions
- $F_{rz}$: Forces applied by the Spindle in the $z$ direction
- $F_d$: Direct forces applied at the center of mass of the Payload
The level of these disturbances has been identified form experiments which are detailed in [[file:disturbances.org][this]] document.
** 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
simulinkproject('../');
#+end_src
#+begin_src matlab
load('mat/conf_simulink.mat');
open('nass_model.slx')
#+end_src
** Plots :ignore:
The measured Amplitude Spectral Densities (ASD) of these forces are shown in Figures [[fig:opt_stiff_dist_gm]] and [[fig:opt_stiff_dist_fty_frz]].
In this study, the expected frequency content of the direct forces applied to the payload is not considered.
#+begin_src matlab :exports none
load('./mat/dist_psd.mat', 'dist_f');
#+end_src
#+begin_src matlab :exports none
figure;
hold on;
plot(dist_f.f, sqrt(dist_f.psd_gm));
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('$\Gamma_{D_w}$ $\left[\frac{m}{\sqrt{Hz}}\right]$')
xlim([1, 500]);
#+end_src
#+header: :tangle no :exports results :results none :noweb yes
#+begin_src matlab :var filepath="figs/opt_stiff_dist_gm.pdf" :var figsize="wide-normal" :post pdf2svg(file=*this*, ext="png")
<<plt-matlab>>
#+end_src
#+name: fig:opt_stiff_dist_gm
#+caption: Amplitude Spectral Density of the Ground Displacement ([[./figs/opt_stiff_dist_gm.png][png]], [[./figs/opt_stiff_dist_gm.pdf][pdf]])
[[file:figs/opt_stiff_dist_gm.png]]
#+begin_src matlab :exports none
figure;
hold on;
plot(dist_f.f, sqrt(dist_f.psd_ty), 'DisplayName', '$F_{T_y}$');
plot(dist_f.f, sqrt(dist_f.psd_rz), 'DisplayName', '$F_{R_z}$');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('ASD $\left[\frac{F}{\sqrt{Hz}}\right]$')
legend('Location', 'southwest');
xlim([1, 500]);
#+end_src
#+header: :tangle no :exports results :results none :noweb yes
#+begin_src matlab :var filepath="figs/opt_stiff_dist_fty_frz.pdf" :var figsize="wide-normal" :post pdf2svg(file=*this*, ext="png")
<<plt-matlab>>
#+end_src
#+name: fig:opt_stiff_dist_fty_frz
#+caption: Amplitude Spectral Density of the "parasitic" forces comming from the Translation stage and the spindle ([[./figs/opt_stiff_dist_fty_frz.png][png]], [[./figs/opt_stiff_dist_fty_frz.pdf][pdf]])
[[file:figs/opt_stiff_dist_fty_frz.png]]
* Effect of disturbances on the position error
<<sec:effect_disturbances>>
** Introduction :ignore:
In this section, we use the Simscape model to identify the transfer function from disturbances to the position error of the sample.
We do that for a wide range of nano-hexapod stiffnesses and we compare the obtained results.
** 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
simulinkproject('../');
#+end_src
#+begin_src matlab
load('mat/conf_simulink.mat');
open('nass_model.slx')
#+end_src
** Initialization
We initialize all the stages with the default parameters.
#+begin_src matlab
initializeGround();
initializeGranite();
initializeTy();
initializeRy();
initializeRz();
initializeMicroHexapod();
initializeAxisc();
initializeMirror();
#+end_src
We use a sample mass of 10kg.
#+begin_src matlab
initializeSample('mass', 10);
#+end_src
We include gravity, and we use no controller.
#+begin_src matlab
initializeSimscapeConfiguration('gravity', true);
initializeController();
initializeDisturbances('enable', false);
initializeLoggingConfiguration('log', 'none');
#+end_src
** Identification
The considered inputs are:
- =Dwx=: Ground displacement in the $x$ direction
- =Dwy=: Ground displacement in the $y$ direction
- =Dwz=: Ground displacement in the $z$ direction
- =Fty_x=: Forces applied by the Translation stage in the $x$ direction
- =Fty_z=: Forces applied by the Translation stage in the $z$ direction
- =Frz_z=: Forces applied by the Spindle in the $z$ direction
- =Fd=: Direct forces applied at the center of mass of the Payload
The outputs are =Ex=, =Ey=, =Ez=, =Erx=, =Ery=, =Erz= which are the 3 positions and 3 orientations errors of the sample.
We initialize the set of the nano-hexapod stiffnesses, and for each of them, we identify the dynamics from defined inputs to defined outputs.
#+begin_src matlab
Ks = logspace(3,9,7); % [N/m]
#+end_src
#+begin_src matlab :exports none
%% Name of the Simulink File
mdl = 'nass_model';
%% Micro-Hexapod
clear io; io_i = 1;
io(io_i) = linio([mdl, '/Disturbances'], 1, 'openinput', [], 'Dwx'); io_i = io_i + 1; % X Ground motion
io(io_i) = linio([mdl, '/Disturbances'], 1, 'openinput', [], 'Dwy'); io_i = io_i + 1; % Y Ground motion
io(io_i) = linio([mdl, '/Disturbances'], 1, 'openinput', [], 'Dwz'); io_i = io_i + 1; % Z Ground motion
io(io_i) = linio([mdl, '/Disturbances'], 1, 'openinput', [], 'Fty_x'); io_i = io_i + 1; % Parasitic force Ty - X
io(io_i) = linio([mdl, '/Disturbances'], 1, 'openinput', [], 'Fty_z'); io_i = io_i + 1; % Parasitic force Ty - Z
io(io_i) = linio([mdl, '/Disturbances'], 1, 'openinput', [], 'Frz_z'); io_i = io_i + 1; % Parasitic force Rz - Z
io(io_i) = linio([mdl, '/Disturbances'], 1, 'openinput', [], 'Fd'); io_i = io_i + 1; % Direct forces
io(io_i) = linio([mdl, '/Tracking Error'], 1, 'openoutput', [], 'En'); io_i = io_i + 1; % Position Error
#+end_src
#+begin_src matlab :exports none
Gd = {zeros(length(Ks), 1)};
for i = 1:length(Ks)
initializeNanoHexapod('k', Ks(i));
% Run the linearization
G = linearize(mdl, io);
G.InputName = {'Dwx', 'Dwy', 'Dwz', 'Fty_x', 'Fty_z', 'Frz_z', 'Fdx', 'Fdy', 'Fdz', 'Mdx', 'Mdy', 'Mdz'};
G.OutputName = {'Ex', 'Ey', 'Ez', 'Erx', 'Ery', 'Erz'};
Gd(i) = {minreal(G)};
end
#+end_src
** Sensitivity to Stages vibration (Filtering)
The sensitivity the stage vibrations are displayed:
- Figure [[fig:opt_stiff_sensitivity_Frz]]: sensitivity to vertical spindle vibrations
- Figure [[fig:opt_stiff_sensitivity_Fty_z]]: sensitivity to vertical translation stage vibrations
- Figure [[fig:opt_stiff_sensitivity_Fty_x]]: sensitivity to horizontal (x) translation stage vibrations
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ez', 'Frz_z'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('Effect of $F_{rz}$ on $E_z$ [m/N]'); xlabel('Frequency [Hz]');
legend('location', 'southwest');
#+end_src
#+header: :tangle no :exports results :results none :noweb yes
#+begin_src matlab :var filepath="figs/opt_stiff_sensitivity_Frz.pdf" :var figsize="full-tall" :post pdf2svg(file=*this*, ext="png")
<<plt-matlab>>
#+end_src
#+name: fig:opt_stiff_sensitivity_Frz
#+caption: Sensitivity to Spindle vertical motion error ($F_{rz}$) to the vertical error position of the sample ($E_z$) ([[./figs/opt_stiff_sensitivity_Frz.png][png]], [[./figs/opt_stiff_sensitivity_Frz.pdf][pdf]])
[[file:figs/opt_stiff_sensitivity_Frz.png]]
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ez', 'Fty_z'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('Effect of $F_{ty}$ on $E_z$ [m/N]'); xlabel('Frequency [Hz]');
legend('location', 'southwest');
#+end_src
#+header: :tangle no :exports results :results none :noweb yes
#+begin_src matlab :var filepath="figs/opt_stiff_sensitivity_Fty_z.pdf" :var figsize="full-tall" :post pdf2svg(file=*this*, ext="png")
<<plt-matlab>>
#+end_src
#+name: fig:opt_stiff_sensitivity_Fty_z
#+caption: Sensitivity to Translation stage vertical motion error ($F_{ty,z}$) to the vertical error position of the sample ($E_z$) ([[./figs/opt_stiff_sensitivity_Fty_z.png][png]], [[./figs/opt_stiff_sensitivity_Fty_z.pdf][pdf]])
[[file:figs/opt_stiff_sensitivity_Fty_z.png]]
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ex', 'Fty_x'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('Effect of $F_{ty}$ on $E_x$ [m/N]'); xlabel('Frequency [Hz]');
legend('location', 'northeast');
#+end_src
#+header: :tangle no :exports results :results none :noweb yes
#+begin_src matlab :var filepath="figs/opt_stiff_sensitivity_Fty_x.pdf" :var figsize="full-tall" :post pdf2svg(file=*this*, ext="png")
<<plt-matlab>>
#+end_src
#+name: fig:opt_stiff_sensitivity_Fty_x
#+caption: Sensitivity to Translation stage $x$ motion error ($F_{ty,x}$) to the error position of the sample in the $x$ direction ($E_x$) ([[./figs/opt_stiff_sensitivity_Fty_x.png][png]], [[./figs/opt_stiff_sensitivity_Fty_x.pdf][pdf]])
[[file:figs/opt_stiff_sensitivity_Fty_x.png]]
** Effect of Ground motion (Transmissibility).
The effect of Ground motion on the position error of the sample is shown in Figure [[fig:opt_stiff_sensitivity_Dw]].
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
ax1 = subplot(1, 2, 1);
hold on;
for i = 1:length(Ks)
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ey', 'Dwy'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('$E_y/D_{wy}$ [m/m]'); xlabel('Frequency [Hz]');
ax2 = subplot(1, 2, 2);
hold on;
for i = 1:length(Ks)
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ez', 'Dwz'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('$E_z/D_{wz}$ [m/m]'); xlabel('Frequency [Hz]');
#+end_src
#+header: :tangle no :exports results :results none :noweb yes
#+begin_src matlab :var filepath="figs/opt_stiff_sensitivity_Dw.pdf" :var figsize="full-tall" :post pdf2svg(file=*this*, ext="png")
<<plt-matlab>>
#+end_src
#+name: fig:opt_stiff_sensitivity_Dw
#+caption: Sensitivity to Ground motion ($D_{w}$) to the position error of the sample ($E_y$ and $E_z$) ([[./figs/opt_stiff_sensitivity_Dw.png][png]], [[./figs/opt_stiff_sensitivity_Dw.pdf][pdf]])
[[file:figs/opt_stiff_sensitivity_Dw.png]]
** Direct Forces (Compliance).
The effect of direct forces/torques applied on the sample (cable forces for instance) on the position error of the sample is shown in Figure [[fig:opt_stiff_sensitivity_Fd]].
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
ax1 = subplot(1, 2, 1);
hold on;
for i = 1:length(Ks)
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ery', 'Mdy'), freqs, 'Hz'))), '-');
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('$E_{ry}/M_{d,y}\ \left[\frac{rad}{N m}\right]$'); xlabel('Frequency [Hz]');
ax2 = subplot(1, 2, 2);
hold on;
for i = 1:length(Ks)
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ez', 'Fdz'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('$E_{z}/F_{d,z}$ [m/N]'); xlabel('Frequency [Hz]');
legend('location', 'northeast');
linkaxes([ax1 ax2], 'xy')
#+end_src
#+header: :tangle no :exports results :results none :noweb yes
#+begin_src matlab :var filepath="figs/opt_stiff_sensitivity_Fd.pdf" :var figsize="full-tall" :post pdf2svg(file=*this*, ext="png")
<<plt-matlab>>
#+end_src
#+name: fig:opt_stiff_sensitivity_Fd
#+caption: Sensitivity to Direct forces and torques applied to the sample ($F_d$, $M_d$) to the position error of the sample ([[./figs/opt_stiff_sensitivity_Fd.png][png]], [[./figs/opt_stiff_sensitivity_Fd.pdf][pdf]])
[[file:figs/opt_stiff_sensitivity_Fd.png]]
** Save :noexport:
#+begin_src matlab
save('./mat/opt_stiffness_disturbances.mat', 'Ks', 'Gd')
#+end_src
** Conclusion
#+begin_important
#+end_important
* Effect of granite stiffness
<<sec:granite_stiffness>>
** Analytical Analysis
#+begin_src latex :file 2dof_system_granite_stiffness.pdf
\begin{tikzpicture}
% ====================
% Parameters
% ====================
\def\massw{2.2} % Width of the masses
\def\massh{0.8} % Height of the masses
\def\spaceh{1.2} % Height of the springs/dampers
\def\dispw{0.3} % Width of the dashed line for the displacement
\def\disph{0.5} % Height of the arrow for the displacements
\def\bracs{0.05} % Brace spacing vertically
\def\brach{-10pt} % Brace shift horizontaly
% ====================
% ====================
% Ground
% ====================
\draw (-0.5*\massw, 0) -- (0.5*\massw, 0);
\draw[dashed] (0.5*\massw, 0) -- ++(\dispw, 0) coordinate(dlow);
\draw[->] (0.5*\massw+0.5*\dispw, 0) -- ++(0, \disph) node[right]{$x_{w}$};
% ====================
% Micro Station
% ====================
\begin{scope}[shift={(0, 0)}]
% Mass
\draw[fill=white] (-0.5*\massw, \spaceh) rectangle (0.5*\massw, \spaceh+\massh) node[pos=0.5]{$m$};
% Spring, Damper, and Actuator
\draw[spring] (-0.4*\massw, 0) -- (-0.4*\massw, \spaceh) node[midway, left=0.1]{$k$};
\draw[damper] (0, 0) -- ( 0, \spaceh) node[midway, left=0.2]{$c$};
% Displacements
\draw[dashed] (0.5*\massw, \spaceh) -- ++(\dispw, 0);
\draw[->] (0.5*\massw+0.5*\dispw, \spaceh) -- ++(0, \disph) node[right]{$x$};
% Legend
\draw[decorate, decoration={brace, amplitude=8pt}, xshift=\brach] %
(-0.5*\massw, \bracs) -- (-0.5*\massw, \spaceh+\massh-\bracs) %
node[midway,rotate=90,anchor=south,yshift=10pt,align=center]{Granite};
\end{scope}
% ====================
% Nano Station
% ====================
\begin{scope}[shift={(0, \spaceh+\massh)}]
% Mass
\draw[fill=white] (-0.5*\massw, \spaceh) rectangle (0.5*\massw, \spaceh+\massh) node[pos=0.5]{$m^\prime$};
% Spring, Damper, and Actuator
\draw[spring] (-0.4*\massw, 0) -- (-0.4*\massw, \spaceh) node[midway, left=0.1]{$k^\prime$};
\draw[damper] (0, 0) -- ( 0, \spaceh) node[midway, left=0.2]{$c^\prime$};
% Displacements
\draw[dashed] (0.5*\massw, \spaceh) -- ++(\dispw, 0) coordinate(dhigh);
\draw[->] (0.5*\massw+0.5*\dispw, \spaceh) -- ++(0, \disph) node[right]{$x^\prime$};
% Legend
\draw[decorate, decoration={brace, amplitude=8pt}, xshift=\brach] %
(-0.5*\massw, \bracs) -- (-0.5*\massw, \spaceh+\massh-\bracs) %
node[midway,rotate=90,anchor=south,yshift=10pt,align=center]{Positioning\\Stages};
\end{scope}
\end{tikzpicture}
#+end_src
#+name: fig:2dof_system_granite_stiffness
#+caption: Figure caption
#+RESULTS:
[[file:figs/2dof_system_granite_stiffness.png]]
If we write the equation of motion of the system in Figure [[fig:2dof_system_granite_stiffness]], we obtain:
\begin{align}
m^\prime s^2 x^\prime &= (c^\prime s + k^\prime) (x - x^\prime) \\
ms^2 x &= (c^\prime s + k^\prime) (x^\prime - x) + (cs + k) (x_w - x)
\end{align}
If we note $d = x^\prime - x$, we obtain:
#+name: eq:plant_ground_transmissibility
\begin{equation}
\frac{d}{x_w} = \frac{-m^\prime s^2 (cs + k)}{ (m^\prime s^2 + c^\prime s + k^\prime) (ms^2 + cs + k) + m^\prime s^2(c^\prime s + k^\prime)}
\end{equation}
** Soft Granite
Let's initialize a soft granite that will act as an isolation stage from ground motion.
#+begin_src matlab
initializeGranite('K', 5e5*ones(3,1), 'C', 5e3*ones(3,1));
#+end_src
#+begin_src matlab
Ks = logspace(3,9,7); % [N/m]
#+end_src
#+begin_src matlab :exports none
Gdr = {zeros(length(Ks), 1)};
#+end_src
#+begin_src matlab
for i = 1:length(Ks)
initializeNanoHexapod('k', Ks(i));
G = linearize(mdl, io);
G.InputName = {'Dwx', 'Dwy', 'Dwz', 'Fty_x', 'Fty_z', 'Frz_z', 'Fdx', 'Fdy', 'Fdz', 'Mdx', 'Mdy', 'Mdz'};
G.OutputName = {'Ex', 'Ey', 'Ez', 'Erx', 'Ery', 'Erz'};
Gdr(i) = {minreal(G)};
end
#+end_src
** Effect of the Granite transfer function
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
hold on;
for i = 1:length(Ks)
set(gca,'ColorOrderIndex',i);
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ez', 'Dwz'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
set(gca,'ColorOrderIndex',i);
plot(freqs, abs(squeeze(freqresp(Gdr{i}('Ez', 'Dwz'), freqs, 'Hz'))), '--', ...
'HandleVisibility', 'off');
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('Amplitude [m/m]'); xlabel('Frequency [Hz]');
legend('location', 'southeast');
#+end_src
#+begin_src matlab :exports none
freqs = logspace(0, 3, 1000);
figure;
hold on;
for i = 1:length(Ks)
set(gca,'ColorOrderIndex',i);
plot(freqs, abs(squeeze(freqresp(Gd{i}('Ez', 'Frz_z'), freqs, 'Hz'))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
set(gca,'ColorOrderIndex',i);
plot(freqs, abs(squeeze(freqresp(Gdr{i}('Ez', 'Frz_z'), freqs, 'Hz'))), '--', ...
'HandleVisibility', 'off');
end
hold off;
set(gca, 'XScale', 'log'); set(gca, 'YScale', 'log');
ylabel('Amplitude [m/N]'); xlabel('Frequency [Hz]');
legend('location', 'southeast');
#+end_src
* Open Loop Budget Error
<<sec:open_loop_budget_error>>
** 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
simulinkproject('../');
#+end_src
#+begin_src matlab
load('mat/conf_simulink.mat');
open('nass_model.slx')
#+end_src
** Load of the identified disturbances and transfer functions
#+begin_src matlab
load('./mat/dist_psd.mat', 'dist_f');
load('./mat/opt_stiffness_disturbances.mat', 'Gd')
#+end_src
** Equations
** Results
Effect of all disturbances
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, sqrt(dist_f.psd_rz).*abs(squeeze(freqresp(Gd{i}('Ez', 'Frz_z'), freqs, 'Hz'))));
end
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('ASD $\left[\frac{m}{\sqrt{Hz}}\right]$')
legend('Location', 'southwest');
xlim([2, 500]);
#+end_src
** Cumulative Amplitude Spectrum
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(dist_f.psd_ty.*abs(squeeze(freqresp(Gd{i}('Ez', 'Fty_z'), freqs, 'Hz'))).^2)))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
plot([freqs(1) freqs(end)], [10e-9 10e-9], 'k--', 'HandleVisibility', 'off');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('CAS $[m]$')
legend('Location', 'southwest');
xlim([2, 500]); ylim([1e-10 1e-6]);
#+end_src
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(dist_f.psd_rz.*abs(squeeze(freqresp(Gd{i}('Ez', 'Frz_z'), freqs, 'Hz'))).^2)))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
plot([freqs(1) freqs(end)], [10e-9 10e-9], 'k--', 'HandleVisibility', 'off');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('CAS $[m]$')
legend('Location', 'southwest');
xlim([2, 500]); ylim([1e-10 1e-6]);
#+end_src
Ground motion
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(dist_f.psd_gm.*abs(squeeze(freqresp(Gd{i}('Ez', 'Dwz'), freqs, 'Hz'))).^2)))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
plot([freqs(1) freqs(end)], [10e-9 10e-9], 'k--', 'HandleVisibility', 'off');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('CAS $E_y$ $[m]$')
legend('Location', 'northeast');
xlim([2, 500]); ylim([1e-10 1e-6]);
#+end_src
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(dist_f.psd_gm.*abs(squeeze(freqresp(Gd{i}('Ex', 'Dwx'), freqs, 'Hz'))).^2)))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
plot([freqs(1) freqs(end)], [10e-9 10e-9], 'k--', 'HandleVisibility', 'off');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'lin');
xlabel('Frequency [Hz]'); ylabel('CAS $E_y$ $[m]$')
legend('Location', 'northeast');
xlim([2, 500]);
#+end_src
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(dist_f.psd_gm.*abs(squeeze(freqresp(Gd{i}('Ey', 'Dwy'), freqs, 'Hz'))).^2)))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
plot([freqs(1) freqs(end)], [10e-9 10e-9], 'k--', 'HandleVisibility', 'off');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'lin');
xlabel('Frequency [Hz]'); ylabel('CAS $E_y$ $[m]$')
legend('Location', 'northeast');
xlim([2, 500]);
#+end_src
Sum of all perturbations
#+begin_src matlab
psd_tot = zeros(length(freqs), length(Ks));
for i = 1:length(Ks)
psd_tot(:,i) = dist_f.psd_gm.*abs(squeeze(freqresp(Gd{i}('Ez', 'Dwz' ), freqs, 'Hz'))).^2 + ...
dist_f.psd_ty.*abs(squeeze(freqresp(Gd{i}('Ez', 'Fty_z'), freqs, 'Hz'))).^2 + ...
dist_f.psd_rz.*abs(squeeze(freqresp(Gd{i}('Ez', 'Frz_z'), freqs, 'Hz'))).^2;
end
#+end_src
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(psd_tot(:,i))))), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)));
end
plot([freqs(1) freqs(end)], [10e-9 10e-9], 'k--', 'HandleVisibility', 'off');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('CAS $E_z$ $[m]$')
legend('Location', 'northeast');
xlim([1, 500]); ylim([1e-10 1e-6]);
#+end_src
* Closed Loop Budget Error
<<sec:closed_loop_budget_error>>
** Introduction :ignore:
** Reduction thanks to feedback - Required bandwidth
#+begin_src matlab
wc = 1*2*pi; % [rad/s]
xic = 0.5;
S = (s/wc)/(1 + s/wc);
bodeFig({S}, logspace(-1,2,1000))
#+end_src
#+begin_src matlab
wc = [1, 5, 10, 20, 50, 100, 200];
S1 = {zeros(length(wc), 1)};
S2 = {zeros(length(wc), 1)};
for j = 1:length(wc)
L = (2*pi*wc(j))/s; % Simple integrator
S1{j} = 1/(1 + L);
L = ((2*pi*wc(j))/s)^2*(1 + s/(2*pi*wc(j)/2))/(1 + s/(2*pi*wc(j)*2));
S2{j} = 1/(1 + L);
end
#+end_src
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
i = 6;
for j = 1:length(wc)
set(gca,'ColorOrderIndex',j);
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(abs(squeeze(freqresp(S1{j}, freqs, 'Hz'))).^2.*psd_tot(:,i))))), '-', ...
'DisplayName', sprintf('$\\omega_c = %.0f$ [Hz]', wc(j)));
end
plot(freqs, sqrt(flip(-cumtrapz(flip(freqs), flip(psd_tot(:,i))))), 'k-', ...
'DisplayName', 'Open-Loop');
plot([freqs(1) freqs(end)], [10e-9 10e-9], 'k--', 'HandleVisibility', 'off');
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Frequency [Hz]'); ylabel('CAS $E_y$ $[m]$')
legend('Location', 'northeast');
xlim([0.5, 500]); ylim([1e-10 1e-6]);
#+end_src
#+begin_src matlab
wc = logspace(0, 3, 100);
Dz1_rms = zeros(length(Ks), length(wc));
Dz2_rms = zeros(length(Ks), length(wc));
for i = 1:length(Ks)
for j = 1:length(wc)
L = (2*pi*wc(j))/s;
Dz1_rms(i, j) = sqrt(trapz(freqs, abs(squeeze(freqresp(1/(1 + L), freqs, 'Hz'))).^2.*psd_tot(:,i)));
L = ((2*pi*wc(j))/s)^2*(1 + s/(2*pi*wc(j)/2))/(1 + s/(2*pi*wc(j)*2));
Dz2_rms(i, j) = sqrt(trapz(freqs, abs(squeeze(freqresp(1/(1 + L), freqs, 'Hz'))).^2.*psd_tot(:,i)));
end
end
#+end_src
#+begin_src matlab
freqs = dist_f.f;
figure;
hold on;
for i = 1:length(Ks)
set(gca,'ColorOrderIndex',i);
plot(wc, Dz1_rms(i, :), '-', ...
'DisplayName', sprintf('$k = %.0g$ [N/m]', Ks(i)))
set(gca,'ColorOrderIndex',i);
plot(wc, Dz2_rms(i, :), '--', ...
'HandleVisibility', 'off')
end
hold off;
set(gca, 'xscale', 'log'); set(gca, 'yscale', 'log');
xlabel('Control Bandwidth [Hz]'); ylabel('$E_z\ [m, rms]$')
legend('Location', 'southwest');
xlim([1, 500]);
#+end_src
* Conclusion
Reference in New Issue View Git Blame Copy Permalink