Tartalomjegyzék

LPV decomposition (version 2)

Teljes Matlab script kiegészítő függvényekkel.

File: d2018_01_13_Kalman_decomp_LPV.m
Directory: projects/3_outsel/2018_01_10_LPV_inversion
Author: Peter Polcz (ppolcz@gmail.com)
Created on 2018. January 13.

Inherited from

Kalman decomposition
File: Kalman_decomposition.m
Directory: demonstrations/oktatas/ccs/2017fall
Author: Peter Polcz (ppolcz@gmail.com)
Created on 2018. January 13.

Inhereted from:

CCS 2017 fall. Homework 2. Solutions
File: d2017_10_26_hf2_mo.m
Directory: demonstrations/oktatas/ccs/2017fall
Author: Peter Polcz (ppolcz@gmail.com)
Created on 2017. October 04.
Reviewed on 2017. October 30.

Nice matrices

Given the following LPV system:

$$ \begin{aligned} &\dot x = A(\rho) x + B(\rho) u,~~~ \rho \in \mathcal P \\ &y = C x \\ &\begin{aligned} \text{where: } & A(\rho) = A_0 + A_1 \rho \in \mathbb{R}^{n\times n} \\ & B(\rho) = B_0 + B_1 \rho \in \mathbb{R}^{n\times r} \\ & C \in \mathbb{R}^{m\times n} \\ & D = 0_{m\times r} \end{aligned} \end{aligned} $$

load LPV_18_01_13_Time211859
A_fh = @(rho) A0 + A1*rho;
B_fh = @(rho) B0 + B1*rho;
C_fh = @(rho) C0 + C1*rho;
D_fh = @(rho) D0 + D1*rho;

A = {A0, A1};
B = {B0, B1};
C = {C0, C1};
D = {D0, D1};

AA = [A{:}];
BB = [B{:}];
\begin{equation} A(\rho) = \left(\begin{array}{ccccccc} 1.0 & -1.1 & 2.7 & -0.068 & 1.0 & -0.39 & -0.68 \\ 0 & -1.210^{-3} & -0.67 & 0 & 1.6 & 0 & -0.82 \\ 0 & 0.17 & -0.7 & 0 & 1.1 & 0 & -0.55 \\ -0.73 & 2.1 & -3.4 & 0.3 & -3.4 & 2.2 & 0.93 \\ -0.73 & 0.82\rho +0.93 & -1.6 & -1.1 & -2.0 & 0.33 & 0.93 \\ 0 & 1.0 & -1.6 & 1.0 & -1.0 & -1.7 & 0 \\ -1.5 & 2.2 & -4.6 & -2.2 & -6.5 & 0.67 & 3.1 \\ \end{array}\right) \end{equation} \begin{equation} B(\rho) = \left(\begin{array}{cc} 0.98 & 0.88 \\ 0 & 0 \\ 0 & 0 \\ 0.85\rho -0.5 & -0.44 \\ 0.15 & 0.13 \\ -0.23 & -0.21 \\ 0.29 & 0.26 \\ \end{array}\right) \end{equation} \begin{equation} C(\rho) = C = \left(\begin{array}{ccccccc} 0 & -1.7 & 2.6 & -0.23 & 0.23 & -2.2 & 0 \\ 0 & -0.4 & 0.61 & -1.3 & 1.3 & 1.0 & 0 \\ \end{array}\right) \end{equation}
tol = 1e-10;

Im_B = pcz_vecalg_orth(BB,tol);
Ker_C = pcz_vecalg_ints(pcz_vecalg_null(C0,tol), pcz_vecalg_null(C1,tol));

[R,V] = CSA(AA,Im_B,Ker_C);

% Check the strong invertibility condition, the followings has to be empty
pcz_vecalg_ints(V,Im_B,tol);

Im_V = pcz_vecalg_orth(V,tol);
Co_V = pcz_vecalg_null(V',tol);
T = [ Im_V Co_V ]
Output:
T =
         0   -1.0000    0.0000    0.0000    0.0000   -0.0000   -0.0000
         0         0    0.4749   -0.5746    0.4671    0.4021   -0.2539
         0    0.0000   -0.1114   -0.1298    0.6797   -0.6031    0.3808
         0    0.0000   -0.7286   -0.3517    0.0646   -0.0822   -0.5784
         0   -0.0000   -0.1775   -0.6778   -0.4097    0.0822    0.5784
         0   -0.0000   -0.4468    0.2644    0.3845    0.6790    0.3486
   -1.0000         0         0         0         0         0         0
N = numel(A);
tA = cell(1,N);
tB = cell(1,N);
tC = cell(1,N);
tD = cell(1,N);

for i = 1:numel(A)
    tA{i} = round(T*A{i}/T, -log10(tol));
    tB{i} = round(T*B{i}, -log10(tol));
    tC{i} = round(C{i}/T, -log10(tol));
    tD{i} = round(D{i}, -log10(tol));
end

tA_fh = @(rho) tA{1} + tA{2}*rho;
tB_fh = @(rho) tB{1} + tB{2}*rho;
tC_fh = @(rho) tC{1} + tC{2}*rho;
tD_fh = @(rho) tD{1} + tD{2}*rho;
\begin{equation} \bar A(\rho) = \left(\begin{array}{ccccccc} -1.210^{-3} & -0.66 & -0.88 & -1.1 & 1.0 & -0.64 & 0 \\ 0.85-0.38\rho & 1.3 & 2.4 & -0.13 & -2.4 & -1.2 & -0.45 \\ -0.56\rho -0.55 & -0.67 & -1.6 & 0.47 & 4.1 & 0.35 & 0.96 \\ 2.2-0.053\rho & 3.4 & 2.2 & -2.2 & -4.2 & -1.4 & -1.1 \\ 0.34\rho +0.49 & -0.52 & 0.64 & -0.34 & 1.3 & -1.4 & 0.051 \\ -0.32\rho -2.3 & -4.6 & -2.2 & 1.5 & 4.1 & 0.14 & 0.98 \\ -1.1 & -1.8 & -0.39 & 1.4 & 1.3 & 1.3 & 1.0 \\ \end{array}\right) \end{equation} \begin{equation} \bar B(\rho) = \left(\begin{array}{cc} 0 & 0 \\ 0.18-0.49\rho & 0.17 \\ 0.42-0.11\rho & 0.37 \\ 0.034-0.3\rho & 0.03 \\ 0.43-0.58\rho & 0.38 \\ 0.23\rho -0.13 & -0.12 \\ -0.98 & -0.88 \\ \end{array}\right) \end{equation} \begin{equation} \bar C(\rho) = C = \left(\begin{array}{ccccccc} 1.7 & 0.59 & 1.2 & -1.6 & -0.58 & -2.6 & 0 \\ 0.4 & 2.1 & 0.37 & 0.021 & 0.33 & 0.6 & 0 \\ \end{array}\right) \end{equation}

Dimensions, number of parameters, dim of V and dim of Bc

np = 2;
nv = size(V,2);
nb = size(Im_B,2);
ny = size(C,1);
nu = size(B,2);
nf = size(B,2);

Kalman decomposition

1. Controllability staircase form

Ez itt egy kicsit csalás, mert a Kalman felbontást csak az $(A_0,B_0,C_0,D_0)$-ra csinálom, majd az itt kapott transzformációs mátrixot alkalmazom az $(A(\rho),B(\rho),C(\rho),D(\rho))$-n. Tudom, hogy ez így nem helyes!

Controllability matrix

Cn = ctrb(A0,B0);

Transformation matrix, which generates the controllability staircase form.

[S1,~,~] = svd(Cn)
Output:
S1 =
   -0.1355   -0.1355   -0.9709   -0.1436    0.0000    0.0000   -0.0000
    0.0000   -0.0000    0.0000   -0.0000   -0.2624    0.7129   -0.6503
    0.0000   -0.0000    0.0000   -0.0000   -0.2249    0.6102    0.7596
    0.4321    0.8806   -0.1897    0.0438   -0.0000    0.0000   -0.0000
   -0.3637    0.1968    0.0474   -0.1634   -0.8393   -0.3091   -0.0002
    0.3652   -0.1119    0.1005   -0.9187    0.0000    0.0000   -0.0000
   -0.7275    0.3937    0.0949   -0.3268    0.4197    0.1546    0.0001

The transformation.

A0_1 = S1' * A0 * S1
B0_1 = S1' * B0
C0_1 = C0 * S1

A1_1 = S1' * A1 * S1
B1_1 = S1' * B1
C1_1 = C1 * S1
Output:
A0_1 =
    1.1164    1.9637   -1.3756    0.3174   -4.0033   -0.8309    1.4621
    1.2288   -1.1523    1.8172   -1.8983    6.4017    1.1012   -5.9312
   -0.1389    0.1220    1.0448    0.1846    1.6372   -0.3675   -2.6655
    0.0646   -0.0513   -0.1960   -1.5872   -3.2284   -0.6264    3.0341
   -0.0000    0.0000    0.0000    0.0000    1.8134    0.5783   -0.2679
    0.0000   -0.0000   -0.0000   -0.0000   -1.3206   -1.1362   -0.9534
   -0.0000    0.0000   -0.0000    0.0000    0.2360    0.1240   -0.1575
B0_1 =
   -0.6978   -0.6253
   -0.3982   -0.3569
   -0.8442   -0.7565
   -0.0660   -0.0592
    0.0000    0.0000
    0.0000    0.0000
   -0.0000   -0.0000
C0_1 =
   -0.9839    0.0858   -0.1647    1.9612   -0.3243    0.2768    3.0821
   -0.6683   -1.0194    0.4181   -1.2345   -1.1380   -0.3260    0.7243
A1_1 =
   -0.0000    0.0000   -0.0000    0.0000    0.0784   -0.2130    0.1943
    0.0000   -0.0000    0.0000   -0.0000   -0.0424    0.1153   -0.1052
    0.0000   -0.0000    0.0000   -0.0000   -0.0102    0.0278   -0.0254
   -0.0000    0.0000   -0.0000    0.0000    0.0352   -0.0957    0.0873
   -0.0000    0.0000   -0.0000    0.0000    0.1809   -0.4916    0.4484
   -0.0000    0.0000   -0.0000    0.0000    0.0666   -0.1810    0.1652
   -0.0000    0.0000   -0.0000    0.0000    0.0000   -0.0001    0.0001
B1_1 =
    0.3680         0
    0.7499         0
   -0.1615         0
    0.0373         0
   -0.0000         0
    0.0000         0
   -0.0000         0
C1_1 =
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0

Order of the controllable subsystem

ab = rank(Cn);

2. Observability staircase form

2.1. Observability staircase form of the controllable subsystem

On1 = obsv(A0_1(1:ab,1:ab), C0_1(:,1:ab));
[~,~,S21] = svd(On1);

2.2. Observability staircase form of the uncontrollable subsystem

On2 = obsv(A0_1(ab+1:end,ab+1:end), C0_1(:,ab+1:end));
[~,~,S22] = svd(On2);

Second transformation matrix to generate the observability staircase form for both controllable and uncontrollable subsystems.

S2 = blkdiag(S21, S22)
Output:
S2 =
    0.1456   -0.8000    0.5611   -0.1549         0         0         0
   -0.2833   -0.5673   -0.6011    0.4864         0         0         0
    0.1574    0.1793    0.4521    0.8595         0         0         0
   -0.9347    0.0775    0.3457   -0.0268         0         0         0
         0         0         0         0    0.6113   -0.7508    0.2502
         0         0         0         0    0.0697   -0.2639   -0.9620
         0         0         0         0   -0.7883   -0.6055    0.1089

Transformation

A0_2 = round(S2' * A0_1 * S2,10);
B0_2 = round(S2' * B0_1,10);
C0_2 = round(C0_1 * S2,10);

A1_2 = round(S2' * A1_1 * S2,10);
B1_2 = round(S2' * B1_1,10);
C1_2 = round(C1_1 * S2,10);
\begin{equation} \bar A(\rho) = \left(\begin{array}{ccccccc} -2.2 & -0.075 & 0 & 0 & 0.017\rho +1.6 & 0.017\rho +0.14 & -0.035\rho -0.02 \\ -0.63 & 1.9 & 0 & 0 & 0.058\rho -1.8 & 0.058\rho -0.87 & 0.19-0.12\rho \\ -1.6 & -1.1 & -1.6 & -0.36 & -0.12\rho -7.4 & 2.7-0.12\rho & 0.24\rho +0.3 \\ 1.3 & 0.49 & 1.3 & 1.3 & 0.065\rho +7.5 & 0.065\rho -0.7 & 0.36-0.13\rho \\ 0 & 0 & 0 & 0 & 0.6-0.18\rho & -0.18\rho -0.61 & 0.36\rho +0.026 \\ 0 & 0 & 0 & 0 & 0.24\rho -1.2 & 0.24\rho +0.59 & -0.48\rho -0.028 \\ 0 & 0 & 0 & 0 & 0.029\rho +0.5 & 0.029\rho -2.1 & -0.059\rho -0.68 \\ \end{array}\right) \end{equation} \begin{equation} \bar B(\rho) = \left(\begin{array}{cc} -0.22\rho -0.06 & -0.054 \\ 0.63-0.75\rho & 0.56 \\ -0.3\rho -0.56 & -0.5 \\ 0.17\rho -0.81 & -0.73 \\ 0 & 0 \\ 0 & 0 \\ 0 & 0 \\ \end{array}\right) \end{equation} \begin{equation} \bar C(\rho) = \left(\begin{array}{ccccccc} -2.0 & 0.86 & 0 & 0 & -2.6 & -1.7 & -0.012 \\ 1.4 & 1.1 & 0 & 0 & -1.3 & 0.5 & 0.11 \\ \end{array}\right) \end{equation}