Teljes Matlab script kiegészítő függvényekkel.
File: d2018_01_10_LPV_inversion.m Author: Peter Polcz (ppolcz@gmail.com)
Created on 2018. January 10.
LPV system inversion benchmark problem Code shows how to solve dynamic unknown input inversion for daLPV systems MIMO example is considered in the sequel
Unknown input direction is parameter dependent with control input action. Parameter dependent measurement map.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% last modified by B. Kulcsar 09/01/2018
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all; clc; close all;
load Penimesterfigyelmebe2
fnames = fieldnames(sys_LPV);
for i = 1:numel(fnames)
varname = fnames{i};
pcz_num2str_multiline(sys_LPV.(varname), ...
'label', sprintf('sys_%s = ', varname), ...
'format','%8.4f')
disp ' '
end
sys_A1 = [
1.0000 0.0500 0.0000 0.0000 0.0000 0.0000
0.0000 0.9828 13.4000 0.4223 0.0000 0.0000
0.0000 -0.0032 0.4345 0.0000 0.0000 0.0000
];
sys_B1 = [
0.0000 0.0000
0.0000 0.0000
0.0595 0.0000
];
sys_C1 = [
1.0000 0.0000 0.0000 0.0000 0.0000
0.0000 1.0000 0.0000 0.0000 0.0000
];
sys_D1 = [
0.0000 0.0000
0.0000 0.0000
];
sys_Ts = [
0.0500
];
sys_A = [
1.0016 0.0335 -0.1540 0.0755 0.0201 0.1662
-0.0042 0.9117 0.6542 -0.0142 -0.0990 -0.3030
0.0004 -0.0069 -0.2807 -0.1783 0.2149 0.8768
];
sys_B = [
-0.2194 0.1980
-1.4856 -1.5972
-0.1421 3.2974
];
sys_C = [
-0.9518 0.0948 0.1576 0.0000 0.0000 0.0000
-0.0394 -0.8749 0.1931 0.0000 0.0000 0.0000
];
sys_D = [
0.0329 0.0000
0.1066 0.0000
];
sys_A1 = [
1.0000 0.0500 0.0000 0.0000 0.0000 0.0000
0.0000 0.9828 13.4000 0.4223 0.0000 0.0000
0.0000 -0.0032 0.4345 0.0000 0.0000 0.0000
];
sys_B1 = [
0.0000 0.0000
0.0000 0.0000
0.0595 0.0000
];
sys_C1 = [
1.0000 0.0000 0.0000 0.0000 0.0000
0.0000 1.0000 0.0000 0.0000 0.0000
];
sys_D1 = [
0.0000 0.0000
0.0000 0.0000
];
sys_Ts = 0.05;
sys_A = [
1.0016 0.0335 -0.1540 0.0755 0.0201 0.1662
-0.0042 0.9117 0.6542 -0.0142 -0.0990 -0.3030
0.0004 -0.0069 -0.2807 -0.1783 0.2149 0.8768
];
sys_B = [
-0.2194 0.1980
-1.4856 -1.5972
-0.1421 3.2974
];
sys_C = [
-0.9518 0.0948 0.1576 0.0000 0.0000 0.0000
-0.0394 -0.8749 0.1931 0.0000 0.0000 0.0000
];
sys_D = [
0.0329 0.0000
0.1066 0.0000
];
A={};Ab={};E={};B={};Bd={};C={};
A{1}=sys_A1(:,1:3);
A{2}=sys_A1(:,4:6);
A{:}
ans =
1.0000 0.0500 0
0 0.9828 13.4000
0 -0.0032 0.4345
ans =
0 0 0
0.4223 0 0
0 0 0
E{1}=sys_B1(:,1)*0;
E{2}=[0 A{2}(2,1) 0]';
E{:}
ans =
0
0
0
ans =
0
0.4223
0
B{1}=sys_B1(:,1);
B{2}=sys_B1(:,2);
B{:}
ans =
0
0
0.0595
ans =
0
0
0
%Bd{1}=sys_B1(:,1);
%Bd{2}=sys_B1(:,2);
Bd{1}=0.1*[0.1;0;0.1];
Bd{2}=[0;.01;0];
Bd{:}
ans =
0.0100
0
0.0100
ans =
0
0.0100
0
C{1}=[1 0 0; 0 1 0];
C{2}=[0 0 0; 0 0 0];
C{:}
ans =
1 0 0
0 1 0
ans =
0 0 0
0 0 0
Dd{1}=[.05;0];
Dd{2}=[0;.05];
Dd{:}
ans =
0.0500
0
ans =
0
0.0500
D{1}=[0;0];
D{2}=[0;0];
D{:}
ans =
0
0
ans =
0
0
Nagyon sok itt mátrix és nem tudom mi hol van. Ha jól sejtem, akkor a következő az állapottérmodell:
$$ \begin{aligned} \dot x = A(\rho) x + E(\rho) u \\ y = C(\rho) x + D(\rho)u \end{aligned} $$
Ezt hogy értelmezzük?
Kiszámítjuk a következő képteret: $Im([E_1 ~ E_2 ]) = Im(E_1) + Im(E_2)$. [OK]
Ec = IMA([E{1} E{2}])
Ec =
0
-1
0
Constant kernel matrix of the output map. Kiszámítjuk a két kimeneti mátrix nullterének metszetét. [OK] $Ker(C_1) \cap Ker(C_2)$
KerC = INTS(KER(C{1}),KER(C{2}))
KerC =
0
0
1
The controllability subspace in ker(C)
Jól sejtem, hogy most ezt számoljuk: minimal $\mathcal A$-invariant subspace containg $\mathcal B$, vagyis $\mathcal R_{(\mathcal A,\mathcal)}$, ami egyben $(\mathcal A,\mathcal)$ invariáns. Ennek a térnek kell a metszetét venni a $K$ sorvektorai által kiveszített eltérrel? A tézist olvasva ennél bonyolultabb a helyzet.
"The set of all $(\mathcal A,\mathcal)$-invariant subspaces contained in a given subspace $\mathcal K$, is an upper semilattice with respect to subspace addition. This semilattice admits a maximum which will be denoted by $V^*$."
Lényegében $\mathcal R = \left< \mathcal A + B \mathcal F ~\middle|~ Im(BK) \right>$ áltércsalád maximumát számolom? Lásd [Balas et.al. 2003, Eqs. (28-30)]
[R,V]=CSA([A{1} A{2}],Ec,KerC);
%dimensions, number of parameters, dim of V and dim of Bc
np=2;
nv=size(V,2);
nb=size(Ec,2);
ny=size(C{1},1);
n=size(A{1},1);
nu=size(E{1},2);
nf=size(E{1},2);
% Check the strong inveribility condition, the followings has to be empty
INTS(V,Ec);
State transformation, invariant feedback design
Itt van egy csomó olyan transzformáció, ami ki van kommentelve. Miért is?
% ORTCO([V]);
% T=[IMA(Ec) ORTCO([V Ec]) V];
% T=[ORTCO(V) V];
% T=[ORTCO(V)'; ORTCO(Ec)']
% T=[IMA(Ec) V];
% T=[SUMS(IMA(Ec),ORTCO(V)) V];
Ezt egyáltalán nem értem. Milyen transzformációs mátrix ez?
T=[IMA(Ec) ORTCO(IMA(Ec))]
T =
0 1 0
1 0 0
0 0 1
A következőt még akár érteném is, mert a (6.4)-re hasonlít, csak fel van cserélve a két állapot helye:
T=[ORTCO(V) V];
Viszont a (9.17)-hez a következő hasonlít a legjobban hasonlít a legjobban:
$$ T = \begin{pmatrix} {V^*}^\perp \\ \Lambda \end{pmatrix}, \text{ ahol } \Lambda \subset \mathcal B^\perp $$
Ugyanakkor egy vektort el kellene hagyni ORTCO(Ec)-ből.
T=[ORTCO(V)'; ORTCO(Ec)']
%Splitting the system matrices into observable and unobservable parts,
%creating the invariance feedback F
nbb=n-nv; %not necessary n-nv =nb!
for i=1:np
Ab{i}=inv(T)*A{i}*T;
Tmp = inv(T)*A{i}*T
%Fault direction
% Eb=inv(T)*Ec
% Eb1=inv(T)*E{1};
% Eb2=inv(T)*E{2};
% E1=Eb(1:nbb,:);
% Eb11=Eb1(1:nbb,:);
% Eb12=Eb2(1:nbb,:);
% %Control input direction
% Bb1=inv(T)*B{1};
% Bb2=inv(T)*B{2};
% Bb11=Bb1(1:nbb,:);
% Bb12=Bb2(1:nbb,:);
% Bb21=Bb1(nbb+1:end,:);
% Bb22=Bb2(nbb+1:end,:);
% %Control input direction
% Bdb1=inv(T)*Bd{1};
% Bdb2=inv(T)*Bd{2};
% Bdb11=Bdb1(1:nbb,:);
% Bdb12=Bdb2(1:nbb,:);
% Bdb21=Bdb1(nbb+1:end,:);
% Bdb22=Bdb2(nbb+1:end,:);
%State matrix transformation
A12{i}=Tmp(1:nbb,nbb+1:end);
A22{i}=Tmp(nbb+1:end,nbb+1:end); %This is the LPV zero dynamics
A11{i}=Tmp(1:nbb,1:nbb);
A21{i}=Tmp(nbb+1:end,1:nbb);
%[A11{i} A12{i}; A21{i} A22{i}]==Ab{i};
%F{i}=-E1\A12{i};
end
A11_concat = [ A11{:} ]
A12_concat = [ A12{:} ]
A21_concat = [ A21{:} ]
A22_concat = [ A22{:} ]
A_concat = [ A11{:} A12{:} ; A21{:} A22{:} ]
Tmp =
0.9828 0 13.4000
0.0500 1.0000 0
-0.0032 0 0.4345
Tmp =
0 0.4223 0
0 0 0
0 0 0
A11_concat =
0.9828 0 0 0.4223
0.0500 1.0000 0 0
A12_concat =
13.4000 0
0 0
A21_concat =
-0.0032 0 0 0
A22_concat =
0.4345 0
A_concat =
0.9828 0 0 0.4223 13.4000 0
0.0500 1.0000 0 0 0 0
-0.0032 0 0 0 0.4345 0
Valóban, a második két egyenletből eltűnik az input hatása, ezt értem:
T*[E{:}]
ans =
0 0.4223
0 0
0 0