2020-10-08 10:52:20 +02:00
<?xml version="1.0" encoding="utf-8"?>
< !DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-strict.dtd">
< html xmlns = "http://www.w3.org/1999/xhtml" lang = "en" xml:lang = "en" >
< head >
2020-11-12 10:44:04 +01:00
<!-- 2020 - 11 - 12 jeu. 10:43 -->
2020-10-08 10:52:20 +02:00
< meta http-equiv = "Content-Type" content = "text/html;charset=utf-8" / >
2020-11-12 10:44:04 +01:00
< title > Sensor Fusion Paper - Matlab Computation< / title >
2020-10-08 10:52:20 +02:00
< meta name = "generator" content = "Org mode" / >
< meta name = "author" content = "Thomas Dehaeze" / >
2020-11-12 10:44:04 +01:00
< link rel = "stylesheet" type = "text/css" href = "https://research.tdehaeze.xyz/css/style.css" / >
< script type = "text/javascript" src = "https://research.tdehaeze.xyz/js/script.js" > < / script >
< script > M a t h J a x = {
tex: {
tags: 'ams',
macros: {bm: ["\\boldsymbol{#1}",1],}
}
};
< / script >
< script type = "text/javascript" src = "https://cdn.jsdelivr.net/npm/mathjax@3/es5/tex-mml-chtml.js" > < / script >
2020-10-08 10:52:20 +02:00
< / head >
< body >
2020-11-12 10:44:04 +01:00
< div id = "org-div-home-and-up" >
< a accesskey = "h" href = "../index.html" > UP < / a >
|
< a accesskey = "H" href = "../index.html" > HOME < / a >
< / div > < div id = "content" >
< h1 class = "title" > Sensor Fusion Paper - Matlab Computation< / h1 >
2020-10-08 10:52:20 +02:00
< div id = "table-of-contents" >
< h2 > Table of Contents< / h2 >
< div id = "text-table-of-contents" >
< ul >
2020-11-12 10:44:04 +01:00
< li > < a href = "#org917402b" > 1. Definition of the plant< / a > < / li >
< li > < a href = "#orgde878e2" > 2. Multiplicative input uncertainty< / a > < / li >
< li > < a href = "#orgc908813" > 3. Specifications and performance weights< / a > < / li >
< li > < a href = "#orgce7235f" > 4. Upper bounds on the norm of the complementary filters for NP, RS and RP< / a > < / li >
< li > < a href = "#orgaf8f3b1" > 5. H-Infinity synthesis of complementary filters< / a > < / li >
< li > < a href = "#orga721ae7" > 6. Complementary filters using analytical formula< / a > < / li >
< li > < a href = "#orgaf265db" > 7. Comparison of complementary filters< / a > < / li >
< li > < a href = "#orga1bd48e" > 8. Controller Analysis< / a > < / li >
< li > < a href = "#org1eb8c07" > 9. Nominal Stability and Nominal Performance< / a > < / li >
< li > < a href = "#org8c33736" > 10. Robust Stability and Robust Performance< / a > < / li >
< li > < a href = "#org367d76d" > 11. Pre-filter< / a > < / li >
< li > < a href = "#org846e5c4" > 12. Matlab Functions< / a >
< ul >
< li > < a href = "#org2f50732" > 12.1. < code > createWeight< / code > < / a > < / li >
< li > < a href = "#org8ef5bfb" > 12.2. < code > plotMagUncertainty< / code > < / a > < / li >
< li > < a href = "#org3b855e9" > 12.3. < code > plotPhaseUncertainty< / code > < / a > < / li >
< / ul >
< / li >
2020-10-08 10:52:20 +02:00
< / ul >
< / div >
< / div >
< p >
2020-11-12 10:44:04 +01:00
The control architecture studied here is shown on figure < a href = "#org16aa33e" > 1< / a > where:
2020-10-08 10:52:20 +02:00
< / p >
< ul class = "org-ul" >
< li > \(G^{\prime}\) is the plant to control< / li >
< li > \(K = G^{-1} H_L^{-1}\) is the controller used with \(G\) a model of the plant< / li >
< li > \(H_L\) and \(H_H\) are complementary filters (\(H_L + H_H = 1\))< / li >
< li > \(K_r\) is a pre-filter that can be added< / li >
< / ul >
2020-11-12 10:44:04 +01:00
< div id = "org16aa33e" class = "figure" >
< p > < img src = "figs-paper/sf_arch_class_prefilter.png" alt = "sf_arch_class_prefilter.png" / >
2020-10-08 10:52:20 +02:00
< / p >
< p > < span class = "figure-number" > Figure 1: < / span > Control Architecture< / p >
< / div >
< p >
Here is the outline of the < code > matlab< / code > analysis for this control architecture:
< / p >
< ul class = "org-ul" >
2020-11-12 10:44:04 +01:00
< li > Section < a href = "#orgfbfc3f4" > 1< / a > :
The plant model \(G\) is defined< / li >
< li > Section < a href = "#orge119055" > 2< / a > :
The plant uncertainty set \(\Pi_I\) is defined using the multiplicative input uncertainty: \(\Pi_I: \ G^\prime = G (1 + w_I \Delta)\).
Thus the weight \(w_I\) is defined such that the true system dynamics is included in the set \(\Pi_I\)< / li >
< li > Section < a href = "#org3dfa664" > 3< / a > :
From the specifications on performance that are expressed in terms of upper bounds of \(S\) and \(T\), performance weights \(w_S\) and \(w_T\) are derived such that the goal is to obtain \(|S| < \frac{1}{|w_S|}\) and \(|T| < \frac{1}{|w_T|}, \ \forall \omega\)< / li >
< li > Section < a href = "#orgacb23fe" > 4< / a > :
Upper bounds on the magnitude of the complementary filters \(|H_L|\) and \(|H_H|\) are defined in order to ensure Nominal Performance (NP), Robust Stability (RS) and Robust Performance (RP)< / li >
< li > Section < a href = "#org4f2901c" > 5< / a > :
\(H_L\) and \(H_H\) are synthesize using the \(\mathcal{H}_\infty\) synthesis such that \(|H_L|\) and \(|H_H|\) are within the specified bounds and such that \(H_L + H_H = 1\) (complementary property)< / li >
< li > Section < a href = "#org1e02524" > 6< / a > :
Analytical formulas of complementary filters are used< / li >
< li > Section < a href = "#orgcb66bdb" > 7< / a > :
The obtained complementary filters using the \(\mathcal{H}_\infty\) synthesis and the analytical formula are compared< / li >
< li > Section < a href = "#org06b9c47" > 8< / a > :
The obtain controller \(K = G^{-1} H_H^{-1}\) is analyzed< / li >
< li > Section < a href = "#org55c5c37" > 9< / a > :
The Nominal Stability (NS) and Nominal Performance conditions are verified< / li >
< li > Section < a href = "#org4d2c1c6" > 10< / a > :
Robust Stability and Robust Performance conditions are studied< / li >
< li > Section < a href = "#org7ea3388" > 11< / a > :
A pre-filter that is used to limit the input usage due to the change of the reference is added< / li >
2020-10-08 10:52:20 +02:00
< / ul >
2020-11-12 10:44:04 +01:00
< div class = "note" id = "orgbb67345" >
2020-10-08 10:52:20 +02:00
< p >
All the files (data and Matlab scripts) are accessible < a href = "data/sensor_fusion.zip" > here< / a > .
< / p >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-org917402b" class = "outline-2" >
< h2 id = "org917402b" > < span class = "section-number-2" > 1< / span > Definition of the plant< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-1" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "orgfbfc3f4" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
The studied system consists of a solid positioned on top of a motorized uni-axial soft suspension.
< / p >
< p >
The absolute position \(x\) of the solid is measured using an inertial sensor and a force \(F\) can be applied to the mass using a voice coil actuator.
< / p >
< p >
2020-11-12 10:44:04 +01:00
The model of the system is represented on figure < a href = "#orgf41ac95" > 2< / a > where the mass of the solid is \(m = 10\ [kg]\), the stiffness of the suspension is \(k = 10^4\ [N/m]\) and the damping of the system is \(c = 10^2\ [N/(m/s)]\).
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "orgf41ac95" class = "figure" >
< p > < img src = "figs-paper/mech_sys_alone.png" alt = "mech_sys_alone.png" / >
2020-10-08 10:52:20 +02:00
< / p >
< p > < span class = "figure-number" > Figure 2: < / span > One degree of freedom system< / p >
< / div >
< p >
2020-11-12 10:44:04 +01:00
The plant \(G\) is defined on matlab and its bode plot is shown on figure < a href = "#org8998aea" > 3< / a > .
2020-10-08 10:52:20 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > m = 10; < span class = "org-comment" > % mass [kg]< / span >
k = 1e4; < span class = "org-comment" > % stiffness [N/m]< / span >
c = 1e2; < span class = "org-comment" > % damping [N/(m/s)]< / span >
2020-10-08 10:52:20 +02:00
2020-11-12 10:44:04 +01:00
G = 1< span class = "org-type" > /< / span > (m< span class = "org-type" > *< / span > s< span class = "org-type" > ^< / span > 2 < span class = "org-type" > +< / span > c< span class = "org-type" > *< / span > s < span class = "org-type" > +< / span > k);
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "org8998aea" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/bode_plot_mech_sys.png" alt = "bode_plot_mech_sys.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 3: < / span > Bode plot of \(G\)< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-orgde878e2" class = "outline-2" >
< h2 id = "orgde878e2" > < span class = "section-number-2" > 2< / span > Multiplicative input uncertainty< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-2" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "orge119055" > < / a >
2020-10-08 10:52:20 +02:00
We choose to use the multiplicative input uncertainty to model the plant uncertainty:
< / p >
2020-11-12 10:44:04 +01:00
\begin{equation}
\Pi_I: \ G^\prime(s) = G(s) (1 + w_I(s) \Delta(s)),\text{ with } |\Delta(j\omega)| < 1 \ \ forall \ omega
\end{equation}
2020-10-08 10:52:20 +02:00
< p >
2020-11-12 10:44:04 +01:00
The uncertainty weight is shown in figure < a href = "#org4aeff4d" > 4< / a > .
2020-10-08 10:52:20 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > wI = createWeight(< span class = "org-string" > 'n'< / span > , 2, < span class = "org-string" > 'w0'< / span > , 2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > 80, < span class = "org-string" > 'G0'< / span > , 0.1, < span class = "org-string" > 'G1'< / span > , 10, < span class = "org-string" > 'Gc'< / span > , 1);
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "org4aeff4d" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/bode_wi.png" alt = "bode_wi.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 4: < / span > Bode plot of \(w_I\)< / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
2020-11-12 10:44:04 +01:00
Elements in the uncertainty set \(\Pi_I\) are computed and their bode plot is shown on figure < a href = "#org51c57f2" > 5< / a > .
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org51c57f2" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/plant_uncertainty_bode_plot.png" alt = "plant_uncertainty_bode_plot.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 5: < / span > Some elements in the uncertainty set \(\Pi_I\)< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-orgc908813" class = "outline-2" >
< h2 id = "orgc908813" > < span class = "section-number-2" > 3< / span > Specifications and performance weights< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-3" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org3dfa664" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
The control objective is to isolate the displacement \(x\) of the mass from the ground motion \(w\).
< / p >
< p >
The specifications are described below:
< / p >
< ul class = "org-ul" >
< li > at least a factor \(10\) of disturbance rejection at \(2\ \text{Hz}\) and with a slope of \(2\) below \(2\ \text{Hz}\) until a rejection of \(10^3\)< / li >
< li > the noise attenuation should be at least \(10\) above \(100\ \text{Hz}\) and with a slope of \(-2\) above< / li >
< / ul >
< p >
2020-11-12 10:44:04 +01:00
These specifications can be represented as upper bounds on the closed loop transfer functions \(S\) and \(T\) (see figure < a href = "#org8dcc955" > 6< / a > ).
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org8dcc955" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/bode_requirements.png" alt = "bode_requirements.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 6: < / span > Upper bounds on \(S\) and \(T\)< / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
We now define two weights, \(w_S(s)\) and \(w_T(s)\) such that \(1/|w_S|\) and \(1/|w_T|\) are lower than the previously defined upper bounds.
Then, the performance specifications are satisfied if the following condition is valid:
\[ \big|S(j\omega)\big| < \frac{1}{|w_S(j\omega)|} ; \quad \big|T(j\omega)\big| < \frac{1}{|w_T(j\omega)|}, \quad \forall \omega \]
< / p >
< p >
2020-11-12 10:44:04 +01:00
The weights are defined as follow. They magnitude is compared with the upper bounds on \(S\) and \(T\) on figure < a href = "#org702c691" > 7< / a > .
2020-10-08 10:52:20 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > wS = 1600< span class = "org-type" > /< / span > (s< span class = "org-type" > +< / span > 0.13)< span class = "org-type" > ^< / span > 2;
wT = 1000< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > (2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > 1000)))< span class = "org-type" > ^< / span > 2;
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "org702c691" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/compare_weights_upper_bounds_S_T.png" alt = "compare_weights_upper_bounds_S_T.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 7: < / span > Weights \(w_S\) and \(w_T\) with the upper bounds on \(S\) and \(T\) obtained from the specifications< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-orgce7235f" class = "outline-2" >
< h2 id = "orgce7235f" > < span class = "section-number-2" > 4< / span > Upper bounds on the norm of the complementary filters for NP, RS and RP< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-4" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "orgacb23fe" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
Now that we have defined \(w_I\), \(w_S\) and \(w_T\), we can derive conditions for Nominal Performance, Robust Stability and Robust Performance (\(j\omega\) is omitted here for readability):
< / p >
\begin{align*}
\text{NP} & \Leftrightarrow |H_H| < \frac{1}{|w_S|} \text{ and } |H_L| < \frac{1}{|w_T|} \quad \forall \omega \\
\text{RS} & \Leftrightarrow |H_L| < \frac{1}{|w_I| (2 + |w_I|)} \quad \forall \omega \\
\text{RP for } S & \Leftarrow |H_H| < \frac{1 + |w_I|}{|w_S| (2 + |w_I|)} \quad \forall \omega \\
\text{RP for } T & \Leftrightarrow |H_L| < \frac{1}{|w_T| (1 + |w_I|) + |w_I|} \quad \forall \omega
\end{align*}
< p >
These conditions are upper bounds on the complementary filters used for control.
< / p >
< p >
2020-11-12 10:44:04 +01:00
We plot these conditions on figure < a href = "#org3739909" > 8< / a > .
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org3739909" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/weights_NP_RS_RP.png" alt = "weights_NP_RS_RP.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 8: < / span > Upper bounds on the norm of the complementary filters for NP, RS and RP< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-orgaf8f3b1" class = "outline-2" >
< h2 id = "orgaf8f3b1" > < span class = "section-number-2" > 5< / span > H-Infinity synthesis of complementary filters< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-5" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org4f2901c" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
We here synthesize the complementary filters using the \(\mathcal{H}_\infty\) synthesis.
The goal is to specify upper bounds on the norms of \(H_L\) and \(H_H\) while ensuring their complementary property (\(H_L + H_H = 1\)).
< / p >
< p >
2020-11-12 10:44:04 +01:00
In order to do so, we use the generalized plant shown on figure < a href = "#org3fbed5b" > 9< / a > where \(w_L\) and \(w_H\) weighting transfer functions that will be used to shape \(H_L\) and \(H_H\) respectively.
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org3fbed5b" class = "figure" >
< p > < img src = "figs-paper/sf_hinf_filters_plant_b.png" alt = "sf_hinf_filters_plant_b.png" / >
2020-10-08 10:52:20 +02:00
< / p >
< p > < span class = "figure-number" > Figure 9: < / span > Generalized plant used for the \(\mathcal{H}_\infty\) synthesis of the complementary filters< / p >
< / div >
< p >
2020-11-12 10:44:04 +01:00
The \(\mathcal{H}_\infty\) synthesis applied on this generalized plant will give a transfer function \(H_L\) (figure < a href = "#org3939c2e" > 10< / a > ) such that the \(\mathcal{H}_\infty\) norm of the transfer function from \(w\) to \([z_H,\ z_L]\) is less than one:
2020-10-08 10:52:20 +02:00
\[ \left\| \begin{array}{c} H_L w_L \\ (1 - H_L) w_H \end{array} \right\|_\infty < 1 \ ]
< / p >
< p >
Thus, if the above condition is verified, we can define \(H_H = 1 - H_L\) and we have that:
\[ \left\| \begin{array}{c} H_L w_L \\ H_H w_H \end{array} \right\|_\infty < 1 \ ]
Which is almost (with an maximum error of \(\sqrt{2}\)) equivalent to:
< / p >
\begin{align*}
|H_L| & < \frac{1}{|w_L|}, \quad \forall \omega \\
|H_H| & < \frac{1}{|w_H|}, \quad \forall \omega
\end{align*}
< p >
We then see that \(w_L\) and \(w_H\) can be used to shape both \(H_L\) and \(H_H\) while ensuring (by definition of \(H_H = 1 - H_L\)) their complementary property.
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org3939c2e" class = "figure" >
< p > < img src = "figs-paper/sf_hinf_filters_b.png" alt = "sf_hinf_filters_b.png" / >
2020-10-08 10:52:20 +02:00
< / p >
< p > < span class = "figure-number" > Figure 10: < / span > \(\mathcal{H}_\infty\) synthesis of the complementary filters< / p >
< / div >
< p >
2020-11-12 10:44:04 +01:00
Thus, if we choose \(w_L\) and \(w_H\) such that \(1/|w_L|\) and \(1/|w_H|\) lie below the upper bounds of figure < a href = "#org3739909" > 8< / a > , we will ensure the NP, RS and RP of the controlled system.
2020-10-08 10:52:20 +02:00
< / p >
< p >
Depending if we are interested only in NP, RS or RP, we can adjust the weights \(w_L\) and \(w_H\).
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > omegab = 2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > 9;
wH = (omegab)< span class = "org-type" > ^< / span > 2< span class = "org-type" > /< / span > (s < span class = "org-type" > +< / span > omegab< span class = "org-type" > *< / span > sqrt(1e< span class = "org-type" > -< / span > 5))< span class = "org-type" > ^< / span > 2;
omegab = 2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > 28;
wL = (s < span class = "org-type" > +< / span > omegab< span class = "org-type" > /< / span > (4.5)< span class = "org-type" > ^< / span > (1< span class = "org-type" > /< / span > 3))< span class = "org-type" > ^< / span > 3< span class = "org-type" > /< / span > (s< span class = "org-type" > *< / span > (1e< span class = "org-type" > -< / span > 4)< span class = "org-type" > ^< / span > (1< span class = "org-type" > /< / span > 3) < span class = "org-type" > +< / span > omegab)< span class = "org-type" > ^< / span > 3;
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "org90a3139" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/weights_wl_wh.png" alt = "weights_wl_wh.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 11: < / span > Weights on the complementary filters \(w_L\) and \(w_H\) and the associated performance weights< / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
We define the generalized plant \(P\) on matlab.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > P = [0 wL;
2020-10-08 10:52:20 +02:00
wH < span class = "org-type" > -< / span > wH;
2020-11-12 10:44:04 +01:00
1 0];
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
< p >
And we do the \(\mathcal{H}_\infty\) synthesis using the < code > hinfsyn< / code > command.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > [Hl_hinf, < span class = "org-type" > ~< / span > , gamma, < span class = "org-type" > ~< / span > ] = hinfsyn(P, 1, 1,< span class = "org-string" > 'TOLGAM'< / span > , 0.001, < span class = "org-string" > 'METHOD'< / span > , < span class = "org-string" > 'ric'< / span > , < span class = "org-string" > 'DISPLAY'< / span > , < span class = "org-string" > 'on'< / span > );
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< pre class = "example" id = "org67390b9" >
2020-10-08 10:52:20 +02:00
[Hl_hinf, ~, gamma, ~] = hinfsyn(P, 1, 1,'TOLGAM', 0.001, 'METHOD', 'ric', 'DISPLAY', 'on');
Test bounds: 0.0000 < gamma < = 1.7285
gamma hamx_eig xinf_eig hamy_eig yinf_eig nrho_xy p/f
1.729 4.1e+01 8.4e-12 1.8e-01 0.0e+00 0.0000 p
0.864 3.9e+01 -5.8e-02# 1.8e-01 0.0e+00 0.0000 f
1.296 4.0e+01 8.4e-12 1.8e-01 0.0e+00 0.0000 p
1.080 4.0e+01 8.5e-12 1.8e-01 0.0e+00 0.0000 p
0.972 3.9e+01 -4.2e-01# 1.8e-01 0.0e+00 0.0000 f
1.026 4.0e+01 8.5e-12 1.8e-01 0.0e+00 0.0000 p
0.999 3.9e+01 8.5e-12 1.8e-01 0.0e+00 0.0000 p
0.986 3.9e+01 -1.2e+00# 1.8e-01 0.0e+00 0.0000 f
0.993 3.9e+01 -8.2e+00# 1.8e-01 0.0e+00 0.0000 f
0.996 3.9e+01 8.5e-12 1.8e-01 0.0e+00 0.0000 p
0.994 3.9e+01 8.5e-12 1.8e-01 0.0e+00 0.0000 p
0.993 3.9e+01 -3.2e+01# 1.8e-01 0.0e+00 0.0000 f
Gamma value achieved: 0.9942
< / pre >
< p >
2020-11-12 10:44:04 +01:00
We then define the high pass filter \(H_H = 1 - H_L\). The bode plot of both \(H_L\) and \(H_H\) is shown on figure < a href = "#orgebb3cd3" > 12< / a > .
2020-10-08 10:52:20 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > Hh_hinf = 1 < span class = "org-type" > -< / span > Hl_hinf;
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "orgebb3cd3" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/hinf_filters_results.png" alt = "hinf_filters_results.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 12: < / span > Obtained complementary filters using \(\mathcal{H}_\infty\) synthesis< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-orga721ae7" class = "outline-2" >
< h2 id = "orga721ae7" > < span class = "section-number-2" > 6< / span > Complementary filters using analytical formula< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-6" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org1e02524" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< ul class = "org-ul" >
< li class = "off" > < code > [  ]< / code > The problem here is that the High pass filter gain goes to zero at high frequency, and when inverted, it will go to infinity.< / li >
< / ul >
2020-10-08 10:52:20 +02:00
< p >
We here use analytical formula for the complementary filters \(H_L\) and \(H_H\).
< / p >
< p >
The first two formulas that are used to generate complementary filters are:
< / p >
\begin{align*}
H_L(s) & = \frac{(1+\alpha) (\frac{s}{\omega_0})+1}{\left((\frac{s}{\omega_0})+1\right) \left((\frac{s}{\omega_0})^2 + \alpha (\frac{s}{\omega_0}) + 1\right)}\\
H_H(s) & = \frac{(\frac{s}{\omega_0})^2 \left((\frac{s}{\omega_0})+1+\alpha\right)}{\left((\frac{s}{\omega_0})+1\right) \left((\frac{s}{\omega_0})^2 + \alpha (\frac{s}{\omega_0}) + 1\right)}
\end{align*}
< p >
where:
< / p >
< ul class = "org-ul" >
< li > \(\omega_0\) is the blending frequency in rad/s.< / li >
< li > \(\alpha\) is used to change the shape of the filters:
< ul class = "org-ul" >
< li > Small values for \(\alpha\) will produce high magnitude of the filters \(|H_L(j\omega)|\) and \(|H_H(j\omega)|\) near \(\omega_0\) but smaller value for \(|H_L(j\omega)|\) above \(\approx 1.5 \omega_0\) and for \(|H_H(j\omega)|\) below \(\approx 0.7 \omega_0\)< / li >
< li > A large \(\alpha\) will do the opposite< / li >
< / ul > < / li >
< / ul >
< p >
2020-11-12 10:44:04 +01:00
This is illustrated on figure < a href = "#org72b4c79" > 13< / a > .
2020-10-08 10:52:20 +02:00
As it is usually wanted to have the \(\| S \|_\infty < 2 \ ) , \ ( \ alpha \ ) between \ ( 0 . 5 \ ) and \ ( 1 \ ) gives a good trade-off between the performance and the robustness .
The slope of those filters at high and low frequencies is \(-2\) and \(2\) respectively for \(H_L\) and \(H_H\).
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org72b4c79" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/comp_filters_param_alpha.png" alt = "comp_filters_param_alpha.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 13: < / span > Effect of the parameter \(\alpha\) on the shape of the generated second order complementary filters< / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
The parameters \(\alpha\) and \(\omega_0\) are chosen in order to have that the complementary filters stay below the defined upper bounds.
< / p >
< p >
2020-11-12 10:44:04 +01:00
The obtained complementary filters are shown on figure < a href = "#org8e27bc0" > 14< / a > .
2020-10-08 10:52:20 +02:00
The Robust Performance is not fulfilled for \(T\), and we see that the RP condition as a slop of \(-3\). We thus have to use different formula for the complementary filters here.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > w0 = 2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > 13;
alpha = 0.8;
2020-10-08 10:52:20 +02:00
2020-11-12 10:44:04 +01:00
Hh2_ana = (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1< span class = "org-type" > +< / span > alpha)< span class = "org-type" > /< / span > (((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2 < span class = "org-type" > +< / span > alpha< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0) < span class = "org-type" > +< / span > 1));
Hl2_ana = ((1< span class = "org-type" > +< / span > alpha)< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1)< span class = "org-type" > /< / span > (((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2 < span class = "org-type" > +< / span > alpha< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0) < span class = "org-type" > +< / span > 1));
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "org8e27bc0" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/complementary_filters_second_order.png" alt = "complementary_filters_second_order.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 14: < / span > Second order complementary filters using the analytical formula< / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
The following formula gives complementary filters with slopes of \(-3\) and \(3\):
< / p >
\begin{align*}
H_L(s) & = \frac{\left(1+(\alpha+1)(\beta+1)\right) (\frac{s}{\omega_0})^2 + (1+\alpha+\beta)(\frac{s}{\omega_0}) + 1}{\left(\frac{s}{\omega_0} + 1\right) \left( (\frac{s}{\omega_0})^2 + \alpha (\frac{s}{\omega_0}) + 1 \right) \left( (\frac{s}{\omega_0})^2 + \beta (\frac{s}{\omega_0}) + 1 \right)}\\
H_H(s) & = \frac{(\frac{s}{\omega_0})^3 \left( (\frac{s}{\omega_0})^2 + (1+\alpha+\beta) (\frac{s}{\omega_0}) + (1+(\alpha+1)(\beta+1)) \right)}{\left(\frac{s}{\omega_0} + 1\right) \left( (\frac{s}{\omega_0})^2 + \alpha (\frac{s}{\omega_0}) + 1 \right) \left( (\frac{s}{\omega_0})^2 + \beta (\frac{s}{\omega_0}) + 1 \right)}
\end{align*}
< p >
The parameters are:
< / p >
< ul class = "org-ul" >
< li > \(\omega_0\) is the blending frequency in rad/s< / li >
< li > \(\alpha\) and \(\beta\) that are used to change the shape of the filters similarly to the parameter \(\alpha\) for the second order complementary filters< / li >
< / ul >
< p >
2020-11-12 10:44:04 +01:00
The filters are defined below and the result is shown on figure < a href = "#org9e9b7f8" > 15< / a > where we can see that the complementary filters are below the defined upper bounds.
2020-10-08 10:52:20 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > alpha = 1;
beta = 10;
w0 = 2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > 14;
2020-10-08 10:52:20 +02:00
2020-11-12 10:44:04 +01:00
Hh3_ana = (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 3 < span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2 < span class = "org-type" > +< / span > (1< span class = "org-type" > +< / span > alpha< span class = "org-type" > +< / span > beta)< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0) < span class = "org-type" > +< / span > (1< span class = "org-type" > +< / span > (alpha< span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > (beta< span class = "org-type" > +< / span > 1)))< span class = "org-type" > /< / span > ((s< span class = "org-type" > /< / span > w0 < span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2< span class = "org-type" > +< / span > alpha< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2< span class = "org-type" > +< / span > beta< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1));
Hl3_ana = ((1< span class = "org-type" > +< / span > (alpha< span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > (beta< span class = "org-type" > +< / span > 1))< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2 < span class = "org-type" > +< / span > (1< span class = "org-type" > +< / span > alpha< span class = "org-type" > +< / span > beta)< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0) < span class = "org-type" > +< / span > 1)< span class = "org-type" > /< / span > ((s< span class = "org-type" > /< / span > w0 < span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2< span class = "org-type" > +< / span > alpha< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1)< span class = "org-type" > *< / span > ((s< span class = "org-type" > /< / span > w0)< span class = "org-type" > ^< / span > 2< span class = "org-type" > +< / span > beta< span class = "org-type" > *< / span > (s< span class = "org-type" > /< / span > w0)< span class = "org-type" > +< / span > 1));
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "org9e9b7f8" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/complementary_filters_third_order.png" alt = "complementary_filters_third_order.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 15: < / span > Third order complementary filters using the analytical formula< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-orgaf265db" class = "outline-2" >
< h2 id = "orgaf265db" > < span class = "section-number-2" > 7< / span > Comparison of complementary filters< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-7" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "orgcb66bdb" > < / a >
The generated complementary filters using \(\mathcal{H}_\infty\) and the analytical formulas are compared on figure < a href = "#org0e60c80" > 16< / a > .
2020-10-08 10:52:20 +02:00
< / p >
< p >
Although they are very close to each other, there is some difference to note here:
< / p >
< ul class = "org-ul" >
< li > the analytical formula provides a very simple way to generate the complementary filters (and thus the controller), they could even be used to tune the controller online using the parameters \(\alpha\) and \(\omega_0\). However, these formula have the property that \(|H_H|\) and \(|H_L|\) are symmetrical with the frequency \(\omega_0\) which may not be desirable.< / li >
< li > while the \(\mathcal{H}_\infty\) synthesis of the complementary filters is not as straightforward as using the analytical formula, it provides a more optimized procedure to obtain the complementary filters< / li >
< / ul >
< p >
The complementary filters obtained with the \(\mathcal{H}_\infty\) will be used for further analysis.
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org0e60c80" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/comp_hinf_analytical.png" alt = "comp_hinf_analytical.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 16: < / span > Comparison of the complementary filters obtained with \(\mathcal{H}_\infty\) synthesis and with the analytical formula< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-orga1bd48e" class = "outline-2" >
< h2 id = "orga1bd48e" > < span class = "section-number-2" > 8< / span > Controller Analysis< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-8" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org06b9c47" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
The controller \(K\) is computed from the plant model \(G\) and the low pass filter \(H_H\):
\[ K = G^{-1} H_H^{-1} \]
< / p >
< p >
As this is not proper and thus realizable, a second order low pass filter is added with a crossover frequency much larger than the control bandwidth.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > omega = 2< span class = "org-type" > *< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > 500;
K = 1< span class = "org-type" > /< / span > (Hh_hinf< span class = "org-type" > *< / span > G) < span class = "org-type" > *< / span > 1< span class = "org-type" > /< / span > ((1< span class = "org-type" > +< / span > s< span class = "org-type" > /< / span > omega)< span class = "org-type" > *< / span > (1< span class = "org-type" > +< / span > s< span class = "org-type" > /< / span > omega< span class = "org-type" > +< / span > (s< span class = "org-type" > /< / span > omega)< span class = "org-type" > ^< / span > 2));
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< pre class = "example" id = "org4380854" >
2020-10-08 10:52:20 +02:00
zpk(K)
ans =
4.961e12 (s+9.915e04) (s^2 + 5s + 500) (s^2 + 284.6s + 2.135e04) (s^2 + 130.5s + 9887)
--------------------------------------------------------------------------------------------------
(s+9.914e04) (s+6283) (s^2 + 0.3576s + 0.03198) (s^2 + 413.8s + 6.398e04) (s^2 + 6283s + 3.948e07)
Continuous-time zero/pole/gain model.
< / pre >
< p >
2020-11-12 10:44:04 +01:00
The bode plot of the controller is shown on figure < a href = "#orga448a8a" > 17< / a > :
2020-10-08 10:52:20 +02:00
< / p >
< ul class = "org-ul" >
< li > two integrator are present at low frequency< / li >
< li > the resonance of the plant at \(3.5\ \text{Hz}\) is inverted (notched)< / li >
< li > a lead is added at \(10\ \text{Hz}\)< / li >
< / ul >
2020-11-12 10:44:04 +01:00
< div id = "orga448a8a" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/bode_plot_controller.png" alt = "bode_plot_controller.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 17: < / span > Bode plot of the obtained controller \(K\)< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-org1eb8c07" class = "outline-2" >
< h2 id = "org1eb8c07" > < span class = "section-number-2" > 9< / span > Nominal Stability and Nominal Performance< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-9" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org55c5c37" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
The nominal stability of the system is first checked with the < code > allmargin< / code > matlab command.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > allmargin(K< span class = "org-type" > *< / span > G< span class = "org-type" > *< / span > Hl_hinf)
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< pre class = "example" id = "orgab687e2" >
2020-10-08 10:52:20 +02:00
allmargin(K*G*Hl_hinf)
ans =
struct with fields:
2020-11-12 10:44:04 +01:00
GainMargin: 3.86196691340257
GMFrequency: 225.992831144248
PhaseMargin: 34.0890006269983
PMFrequency: 88.3632737129296
DelayMargin: 0.0067331740287082
DMFrequency: 88.3632737129296
2020-10-08 10:52:20 +02:00
Stable: 1
< / pre >
< p >
The system is stable and the stability margins are good.
< / p >
< p >
2020-11-12 10:44:04 +01:00
The bode plot of the loop gain \(L = K*G*H_L\) is shown on figure < a href = "#orgcce877a" > 18< / a > .
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "orgcce877a" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/bode_plot_loop_gain.png" alt = "bode_plot_loop_gain.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 18: < / span > Bode Plot of the Loop Gain \(L = K G H_L\)< / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
In order to check the Nominal Performance of the system, we compute the sensibility and the complementary sensibility transfer functions.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > S = 1< span class = "org-type" > /< / span > (K< span class = "org-type" > *< / span > G< span class = "org-type" > *< / span > Hl_hinf < span class = "org-type" > +< / span > 1);
T = K< span class = "org-type" > *< / span > G< span class = "org-type" > *< / span > Hl_hinf< span class = "org-type" > /< / span > (K< span class = "org-type" > *< / span > G< span class = "org-type" > *< / span > Hl_hinf < span class = "org-type" > +< / span > 1);
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
< p >
2020-11-12 10:44:04 +01:00
We then compare their norms with the upper bounds on the performance of the system (figure < a href = "#orgbbc226a" > 19< / a > ).
2020-10-08 10:52:20 +02:00
As expected, we guarantee the Nominal Performance of the system.
< / p >
2020-11-12 10:44:04 +01:00
< div id = "orgbbc226a" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/verification_NP.png" alt = "verification_NP.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 19: < / span > Bode plot of \(S\) and \(T\) in order to verify the nominal performance of the system< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-org8c33736" class = "outline-2" >
< h2 id = "org8c33736" > < span class = "section-number-2" > 10< / span > Robust Stability and Robust Performance< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-10" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org4d2c1c6" > < / a >
2020-10-08 10:52:20 +02:00
In order to verify the Robust stability of the system, we can use the following equivalence:
\[ \text{RS} \Leftrightarrow \left| w_I T \right| < 1 \ quad \ forall \ omega \ ]
< / p >
< p >
2020-11-12 10:44:04 +01:00
This is shown on figure < a href = "#orga0177a7" > 20< / a > .
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "orga0177a7" class = "figure" >
< p > < img src = "figs/robust_stability.png" alt = "robust_stability.png" / >
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 20: < / span > Robust Stability Check: \(|w_I T| < 1 , \ quad \ forall \ omega \ ) < / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
2020-11-12 10:44:04 +01:00
To check Robust Stability, we can also look at the loop gain of the uncertain system (figure < a href = "#orgfb25f29" > 21< / a > ) or the Nyquist plot (figure < a href = "#orgf81e346" > 22< / a > ).
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "orgfb25f29" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/loop_gain_robustness.png" alt = "loop_gain_robustness.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 21: < / span > Loop Gain of the uncertain system< / p >
2020-10-08 10:52:20 +02:00
< / div >
2020-11-12 10:44:04 +01:00
< div id = "orgf81e346" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/nyquist_robustness.png" alt = "nyquist_robustness.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 22: < / span > Nyquist plot of the uncertain system< / p >
2020-10-08 10:52:20 +02:00
< / div >
< p >
The Robust Performance is verified by plotting \(|S|\) and \(|T|\) for the uncertain system along side the upper bounds defined for performance.
2020-11-12 10:44:04 +01:00
This is shown on figure < a href = "#org58bb0a6" > 23< / a > and we can indeed confirmed that the robust performance property of the system is valid.
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div id = "org58bb0a6" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/robust_performance_result.png" alt = "robust_performance_result.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 23: < / span > Verification of the Robust Performance< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-org367d76d" class = "outline-2" >
< h2 id = "org367d76d" > < span class = "section-number-2" > 11< / span > Pre-filter< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-11" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org7ea3388" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
For now, we have not specified any performance requirements on the input usage due to change of reference.
2020-11-12 10:44:04 +01:00
Do limit the input usage due to change of reference, we can use a pre-filter \(K_r\) as shown on figure < a href = "#org16aa33e" > 1< / a > .
2020-10-08 10:52:20 +02:00
< / p >
< p >
If we want a response time of 0.01 second, we can choose a first order low pass filter with a crossover frequency at \(1/0.01 = 100\ \text{Hz}\) as defined below.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > Kr = 1< span class = "org-type" > /< / span > (1< span class = "org-type" > +< / span > s< span class = "org-type" > /< / span > 2< span class = "org-type" > /< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > /< / span > 100);
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
< p >
2020-11-12 10:44:04 +01:00
We then run a simulation for a step of \(1\mu m\) for the system without and with the pre-filter \(K_r\) (figure < a href = "#orgb70a4d3" > 24< / a > ).
2020-10-08 10:52:20 +02:00
This confirms that a pre-filter can be used to limit the input usage due to change of reference using this architecture.
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > t = linspace(0, 0.01, 500);
2020-10-08 10:52:20 +02:00
opts = stepDataOptions;
2020-11-12 10:44:04 +01:00
opts.StepAmplitude = 1e< span class = "org-type" > -< / span > 6;
2020-10-08 10:52:20 +02:00
2020-11-12 10:44:04 +01:00
u = step((K)< span class = "org-type" > /< / span > (1< span class = "org-type" > +< / span > G< span class = "org-type" > *< / span > K< span class = "org-type" > *< / span > Hl_hinf), t, opts);
uf = step((Kr< span class = "org-type" > *< / span > K)< span class = "org-type" > /< / span > (1< span class = "org-type" > +< / span > G< span class = "org-type" > *< / span > K< span class = "org-type" > *< / span > Hl_hinf), t, opts);
y = step((K< span class = "org-type" > *< / span > G)< span class = "org-type" > /< / span > (1< span class = "org-type" > +< / span > G< span class = "org-type" > *< / span > K< span class = "org-type" > *< / span > Hl_hinf), t, opts);
yf = step((Kr< span class = "org-type" > *< / span > G< span class = "org-type" > *< / span > K)< span class = "org-type" > /< / span > (1< span class = "org-type" > +< / span > G< span class = "org-type" > *< / span > K< span class = "org-type" > *< / span > Hl_hinf), t, opts);
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "orgb70a4d3" class = "figure" >
2020-10-08 10:52:20 +02:00
< p > < img src = "figs/u_and_y_with_Kr.png" alt = "u_and_y_with_Kr.png" / >
< / p >
2020-11-12 10:44:04 +01:00
< p > < span class = "figure-number" > Figure 24: < / span > Input usage and response due to a step change of reference when using a pre-filter \(K_r\)< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-org846e5c4" class = "outline-2" >
< h2 id = "org846e5c4" > < span class = "section-number-2" > 12< / span > Matlab Functions< / h2 >
2020-10-08 10:52:20 +02:00
< div class = "outline-text-2" id = "text-12" >
< p >
2020-11-12 10:44:04 +01:00
< a id = "org4fb0718" > < / a >
< / p >
< / div >
< div id = "outline-container-org2f50732" class = "outline-3" >
< h3 id = "org2f50732" > < span class = "section-number-3" > 12.1< / span > < code > createWeight< / code > < / h3 >
< div class = "outline-text-3" id = "text-12-1" >
< p >
< a id = "org01d0043" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
2020-11-12 10:44:04 +01:00
This Matlab function is accessible < a href = "src/createWeight.m" > here< / a > .
2020-10-08 10:52:20 +02:00
< / p >
< div class = "org-src-container" >
2020-11-12 10:44:04 +01:00
< pre class = "src src-matlab" > < span class = "org-keyword" > function< / span > < span class = "org-variable-name" > [W]< / span > = < span class = "org-function-name" > createWeight< / span > (< span class = "org-variable-name" > args< / span > )
< span class = "org-comment" > % createWeight -< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Syntax: [in_data] = createWeight(in_data)< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Inputs:< / span >
< span class = "org-comment" > % - n - Weight Order< / span >
< span class = "org-comment" > % - G0 - Low frequency Gain< / span >
< span class = "org-comment" > % - G1 - High frequency Gain< / span >
< span class = "org-comment" > % - Gc - Gain of W at frequency w0< / span >
< span class = "org-comment" > % - w0 - Frequency at which |W(j w0)| = Gc< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Outputs:< / span >
< span class = "org-comment" > % - W - Generated Weight< / span >
arguments
args.n (1,1) double {mustBeInteger, mustBePositive} = 1
args.G0 (1,1) double {mustBeNumeric, mustBePositive} = 0.1
args.G1 (1,1) double {mustBeNumeric, mustBePositive} = 10
args.Gc (1,1) double {mustBeNumeric, mustBePositive} = 1
args.w0 (1,1) double {mustBeNumeric, mustBePositive} = 1
< span class = "org-keyword" > end< / span >
mustBeBetween(args.G0, args.Gc, args.G1);
s = tf(< span class = "org-string" > 's'< / span > );
W = (((1< span class = "org-type" > /< / span > args.w0)< span class = "org-type" > *< / span > sqrt((1< span class = "org-type" > -< / span > (args.G0< span class = "org-type" > /< / span > args.Gc)< span class = "org-type" > ^< / span > (2< span class = "org-type" > /< / span > args.n))< span class = "org-type" > /< / span > (1< span class = "org-type" > -< / span > (args.Gc< span class = "org-type" > /< / span > args.G1)< span class = "org-type" > ^< / span > (2< span class = "org-type" > /< / span > args.n)))< span class = "org-type" > *< / span > s < span class = "org-type" > +< / span > (args.G0< span class = "org-type" > /< / span > args.Gc)< span class = "org-type" > ^< / span > (1< span class = "org-type" > /< / span > args.n))< span class = "org-type" > /< / span > ((1< span class = "org-type" > /< / span > args.G1)< span class = "org-type" > ^< / span > (1< span class = "org-type" > /< / span > args.n)< span class = "org-type" > *< / span > (1< span class = "org-type" > /< / span > args.w0)< span class = "org-type" > *< / span > sqrt((1< span class = "org-type" > -< / span > (args.G0< span class = "org-type" > /< / span > args.Gc)< span class = "org-type" > ^< / span > (2< span class = "org-type" > /< / span > args.n))< span class = "org-type" > /< / span > (1< span class = "org-type" > -< / span > (args.Gc< span class = "org-type" > /< / span > args.G1)< span class = "org-type" > ^< / span > (2< span class = "org-type" > /< / span > args.n)))< span class = "org-type" > *< / span > s < span class = "org-type" > +< / span > (1< span class = "org-type" > /< / span > args.Gc)< span class = "org-type" > ^< / span > (1< span class = "org-type" > /< / span > args.n)))< span class = "org-type" > ^< / span > args.n;
< span class = "org-keyword" > end< / span >
< span class = "org-comment" > % Custom validation function< / span >
< span class = "org-keyword" > function< / span > < span class = "org-function-name" > mustBeBetween< / span > (< span class = "org-variable-name" > a< / span > ,< span class = "org-variable-name" > b< / span > ,< span class = "org-variable-name" > c< / span > )
< span class = "org-keyword" > if< / span > < span class = "org-type" > ~< / span > ((a < span class = "org-type" > > < / span > b < span class = "org-type" > & & < / span > b < span class = "org-type" > > < / span > c) < span class = "org-type" > ||< / span > (c < span class = "org-type" > > < / span > b < span class = "org-type" > & & < / span > b < span class = "org-type" > > < / span > a))
eid = < span class = "org-string" > 'createWeight:inputError'< / span > ;
msg = < span class = "org-string" > 'Gc should be between G0 and G1.'< / span > ;
throwAsCaller(MException(eid,msg))
< span class = "org-keyword" > end< / span >
< span class = "org-keyword" > end< / span >
2020-10-08 10:52:20 +02:00
< / pre >
< / div >
2020-11-12 10:44:04 +01:00
< / div >
< / div >
2020-10-08 10:52:20 +02:00
2020-11-12 10:44:04 +01:00
< div id = "outline-container-org8ef5bfb" class = "outline-3" >
< h3 id = "org8ef5bfb" > < span class = "section-number-3" > 12.2< / span > < code > plotMagUncertainty< / code > < / h3 >
< div class = "outline-text-3" id = "text-12-2" >
2020-10-08 10:52:20 +02:00
< p >
2020-11-12 10:44:04 +01:00
< a id = "orgc93a4bb" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
< p >
2020-11-12 10:44:04 +01:00
This Matlab function is accessible < a href = "src/plotMagUncertainty.m" > here< / a > .
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div class = "org-src-container" >
< pre class = "src src-matlab" > < span class = "org-keyword" > function< / span > < span class = "org-variable-name" > [p]< / span > = < span class = "org-function-name" > plotMagUncertainty< / span > (< span class = "org-variable-name" > W< / span > , < span class = "org-variable-name" > freqs< / span > , < span class = "org-variable-name" > args< / span > )
< span class = "org-comment" > % plotMagUncertainty -< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Syntax: [p] = plotMagUncertainty(W, freqs, args)< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Inputs:< / span >
< span class = "org-comment" > % - W - Multiplicative Uncertainty Weight< / span >
< span class = "org-comment" > % - freqs - Frequency Vector [Hz]< / span >
< span class = "org-comment" > % - args - Optional Arguments:< / span >
< span class = "org-comment" > % - G< / span >
< span class = "org-comment" > % - color_i< / span >
< span class = "org-comment" > % - opacity< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Outputs:< / span >
< span class = "org-comment" > % - p - Plot Handle< / span >
arguments
W
freqs double {mustBeNumeric, mustBeNonnegative}
args.G = tf(1)
args.color_i (1,1) double {mustBeInteger, mustBeNonnegative} = 0
args.opacity (1,1) double {mustBeNumeric, mustBeNonnegative} = 0.3
args.DisplayName char = < span class = "org-string" > ''< / span >
< span class = "org-keyword" > end< / span >
< span class = "org-comment" > % Get defaults colors< / span >
colors = < span class = "org-type" > get< / span > (< span class = "org-variable-name" > groot< / span > , < span class = "org-string" > 'defaultAxesColorOrder'< / span > );
p = < span class = "org-type" > patch< / span > ([freqs flip(freqs)], ...
[abs(squeeze(freqresp(args.G, freqs, < span class = "org-string" > 'Hz'< / span > )))< span class = "org-type" > .*< / span > (1 < span class = "org-type" > +< / span > abs(squeeze(freqresp(W, freqs, < span class = "org-string" > 'Hz'< / span > )))); ...
flip(abs(squeeze(freqresp(args.G, freqs, < span class = "org-string" > 'Hz'< / span > )))< span class = "org-type" > .*< / span > max(1 < span class = "org-type" > -< / span > abs(squeeze(freqresp(W, freqs, < span class = "org-string" > 'Hz'< / span > ))), 1e< span class = "org-type" > -< / span > 6))], < span class = "org-string" > 'w'< / span > , ...
< span class = "org-string" > 'DisplayName'< / span > , args.DisplayName);
< span class = "org-keyword" > if< / span > args.color_i < span class = "org-type" > ==< / span > 0
p.FaceColor = [0; 0; 0];
< span class = "org-keyword" > else< / span >
p.FaceColor = colors(args.color_i, < span class = "org-type" > :< / span > );
< span class = "org-keyword" > end< / span >
p.EdgeColor = < span class = "org-string" > 'none'< / span > ;
p.FaceAlpha = args.opacity;
< span class = "org-keyword" > end< / span >
< / pre >
< / div >
< / div >
2020-10-08 10:52:20 +02:00
< / div >
2020-11-12 10:44:04 +01:00
< div id = "outline-container-org3b855e9" class = "outline-3" >
< h3 id = "org3b855e9" > < span class = "section-number-3" > 12.3< / span > < code > plotPhaseUncertainty< / code > < / h3 >
< div class = "outline-text-3" id = "text-12-3" >
2020-10-08 10:52:20 +02:00
< p >
2020-11-12 10:44:04 +01:00
< a id = "org9d5c371" > < / a >
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< p >
This Matlab function is accessible < a href = "src/plotPhaseUncertainty.m" > here< / a > .
2020-10-08 10:52:20 +02:00
< / p >
2020-11-12 10:44:04 +01:00
< div class = "org-src-container" >
< pre class = "src src-matlab" > < span class = "org-keyword" > function< / span > < span class = "org-variable-name" > [p]< / span > = < span class = "org-function-name" > plotPhaseUncertainty< / span > (< span class = "org-variable-name" > W< / span > , < span class = "org-variable-name" > freqs< / span > , < span class = "org-variable-name" > args< / span > )
< span class = "org-comment" > % plotPhaseUncertainty -< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Syntax: [p] = plotPhaseUncertainty(W, freqs, args)< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Inputs:< / span >
< span class = "org-comment" > % - W - Multiplicative Uncertainty Weight< / span >
< span class = "org-comment" > % - freqs - Frequency Vector [Hz]< / span >
< span class = "org-comment" > % - args - Optional Arguments:< / span >
< span class = "org-comment" > % - G< / span >
< span class = "org-comment" > % - color_i< / span >
< span class = "org-comment" > % - opacity< / span >
< span class = "org-comment" > %< / span >
< span class = "org-comment" > % Outputs:< / span >
< span class = "org-comment" > % - p - Plot Handle< / span >
arguments
W
freqs double {mustBeNumeric, mustBeNonnegative}
args.G = tf(1)
args.color_i (1,1) double {mustBeInteger, mustBeNonnegative} = 0
args.opacity (1,1) double {mustBeNumeric, mustBePositive} = 0.3
args.DisplayName char = < span class = "org-string" > ''< / span >
< span class = "org-keyword" > end< / span >
< span class = "org-comment" > % Get defaults colors< / span >
colors = < span class = "org-type" > get< / span > (< span class = "org-variable-name" > groot< / span > , < span class = "org-string" > 'defaultAxesColorOrder'< / span > );
< span class = "org-comment" > % Compute Phase Uncertainty< / span >
Dphi = 180< span class = "org-type" > /< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > asin(abs(squeeze(freqresp(W, freqs, < span class = "org-string" > 'Hz'< / span > ))));
Dphi(abs(squeeze(freqresp(W, freqs, < span class = "org-string" > 'Hz'< / span > ))) < span class = "org-type" > > < / span > 1) = 360;
< span class = "org-comment" > % Compute Plant Phase< / span >
G_ang = 180< span class = "org-type" > /< / span > < span class = "org-constant" > pi< / span > < span class = "org-type" > *< / span > angle(squeeze(freqresp(args.G, freqs, < span class = "org-string" > 'Hz'< / span > )));
p = < span class = "org-type" > patch< / span > ([freqs flip(freqs)], [G_ang< span class = "org-type" > +< / span > Dphi; flip(G_ang< span class = "org-type" > -< / span > Dphi)], < span class = "org-string" > 'w'< / span > , ...
< span class = "org-string" > 'DisplayName'< / span > , args.DisplayName);
< span class = "org-keyword" > if< / span > args.color_i < span class = "org-type" > ==< / span > 0
p.FaceColor = [0; 0; 0];
< span class = "org-keyword" > else< / span >
p.FaceColor = colors(args.color_i, < span class = "org-type" > :< / span > );
< span class = "org-keyword" > end< / span >
p.EdgeColor = < span class = "org-string" > 'none'< / span > ;
p.FaceAlpha = args.opacity;
< span class = "org-keyword" > end< / span >
< / pre >
< / div >
2020-10-08 10:52:20 +02:00
< / div >
< / div >
< / div >
< / div >
< div id = "postamble" class = "status" >
< p class = "author" > Author: Thomas Dehaeze< / p >
2020-11-12 10:44:04 +01:00
< p class = "date" > Created: 2020-11-12 jeu. 10:43< / p >
2020-10-08 10:52:20 +02:00
< / div >
< / body >
< / html >