Teljes Matlab script kiegészítő függvényekkel.
Check controllability and observability of an LTI system using LMI conditions
File: LTI_ctrb_obsv.m Directory: 2_demonstrations/workspace/ccs/ccs_2018 Author: Peter Polcz (ppolcz@gmail.com)
Created on 2018. November 06.
a = 4; % Nr. of contr. and obs.
b = 2; % Nr. of contr. and unobs.
c = 1; % Nr. of uncontr. and obs.
d = 2; % Nr. of uncontr. and unobs.
n = a+b+c+d; % Nr. of states
m = 2; % Nr. of inputs
p = 3; % Nr. of outputs
is_stabilisable = 0;
[A,B,C,D] = LTI_generate_sys(a,b,c,d,m,p,is_stabilisable);
pcz_display(A,B,C,D)
A [9x9] =
Columns 1 through 7
-0.3854 -0.1115 -0.8659 0.9982 0.1392 -0.1227 -0.5507
0.3527 -0.3931 -0.2876 0.0364 -0.3432 0.2381 -0.2856
0.3185 1.5245 1.0902 0.1159 1.3377 0.2412 -0.4362
-0.7679 0.6399 0.7834 -0.6353 -0.3211 -0.9986 -0.4801
-0.5263 -1.8969 -0.3223 -1.4173 -1.7810 0.2332 0.2312
0.7449 -0.6363 0.0445 -0.1514 -0.3093 -0.7250 0.7354
-0.3823 -1.2174 0.8619 0.3559 0.6969 0.2617 1.4206
0.8841 1.6209 -1.1647 1.1125 0.5730 -1.0474 -1.1842
0.5016 -0.4082 -0.5664 1.0788 0.0259 -1.0732 0.2587
Columns 8 through 9
0.5190 0.4326
-0.2757 -0.2659
0.3040 0.4320
-0.6928 -0.0070
-0.9877 -0.6107
-0.4421 0.0664
0.6699 0.3284
0.2541 0.7315
0.1711 -0.5239
B [9x2] =
-0.6677 1.0997
-0.0440 0.1737
-0.6871 -0.6192
-0.9191 0.2573
0.5206 0.9395
-0.4330 1.8049
0.6751 -0.0365
-0.3556 -1.2189
-0.1271 0.6115
C [3x9] =
Columns 1 through 7
-1.5725 -0.0285 -1.0454 -0.7619 0.2302 -1.1239 1.3237
1.0396 0.4461 0.4953 0.5334 -0.1662 0.1683 -0.3672
-0.0657 0.7354 0.2316 0.0232 0.1102 -0.4217 0.1611
Columns 8 through 9
1.3796 0.4617
-0.8704 0.6409
-0.0781 0.5073
D [3x2] =
0 0
0 0
0 0
Check wether the system is stabilisable.
Cn = ctrb(A,B);
[U,Sigma,V] = svd(Cn);
T = U';
A_Ctrb_Staircase = round(T*A*T',10)
B_Ctrb_Staircase = round(T*B,10)
k = c+d-1;
A_Unctrb = A_Ctrb_Staircase(end-k:end,end-k:end)
Uncontrollable_Eigenvalues = eig(A_Unctrb)
Check wether the system is detectable.
On = obsv(A,C);
[U,Sigma,V] = svd(On);
T = V';
A_Obsv_Staircase = round(T*A*T',10)
C_Obsv_Staircase = round(C*T',10)
k = b+d-1;
A_Unobsv = A_Obsv_Staircase(end-k:end,end-k:end)
Unobservable_Eigenvalues = eig(A_Unobsv)
A_Ctrb_Staircase =
Columns 1 through 7
-1.9598 1.1465 -0.2252 0.7011 0.2094 0.4224 1.6623
-0.3071 0.6537 -1.8117 -0.6719 0.0106 -0.3603 2.3914
-0.0673 0.8104 0.0956 -0.6325 -1.6739 -1.0993 0.0284
-0.0318 0.0655 -0.0627 -0.6385 0.9473 1.5515 1.4845
0.0094 -0.0190 0.0170 0.1275 -0.7077 0.5763 0.5926
0.0028 -0.0064 0.0062 0.0336 -0.2870 0.4305 0.1154
0 0 0 0 0 0 0.6165
0 0 0 0 0 0 0.4374
0 0 0 0 0 0 0.1456
Columns 8 through 9
-1.1410 -1.2883
-1.2424 -0.2010
0.1954 0.2247
-0.7585 -0.0685
-0.5822 -0.2886
0.0241 0.0826
0.3062 2.5902
-0.3872 0.1614
-0.3808 0.2182
B_Ctrb_Staircase =
-0.1374 -1.6309
-1.0618 0.0152
0.6935 0.6369
-0.5930 -0.9557
0.4561 -1.2795
-0.7966 1.4403
0 0
0 0
0 0
A_Unctrb =
0.6165 0.3062 2.5902
0.4374 -0.3872 0.1614
0.1456 -0.3808 0.2182
Uncontrollable_Eigenvalues =
-0.8166
0.7065
0.5575
A_Obsv_Staircase =
Columns 1 through 7
0.1586 -0.9420 0.0117 -0.0009 -0.0149 0 0
1.4624 0.4363 0.0625 0.0313 0.0037 0 0
-0.0182 0.1673 -1.1065 -0.0118 0.0434 0 0
-0.3312 -1.7904 0.3777 0.2239 -0.0169 0 0
1.0048 1.2020 -0.4829 0.3304 -0.7166 0 0
0.2582 0.0235 0.2845 0.6016 0.4536 -1.7981 0.0095
-0.0123 -0.2227 -0.5979 -0.2324 0.1160 0.6361 0.9252
-0.1025 0.5023 1.2264 0.1758 -0.5110 -0.7362 -1.9214
-0.5494 1.0979 2.2412 -1.2708 -1.9260 -0.9490 1.7835
Columns 8 through 9
0 0
0 0
0 0
0 0
0 0
0.0028 0.0012
0.0140 -0.0108
0.2743 0.0040
-2.3174 -0.0757
C_Obsv_Staircase =
Columns 1 through 7
-2.0664 0.2258 -0.1648 2.2170 -0.2618 0 0
0.7404 -0.1379 0.0075 -1.2983 -0.9570 0 0
-0.2656 0.5642 0.1135 -0.0946 -0.8173 0 0
Columns 8 through 9
0 0
0 0
0 0
A_Unobsv =
-1.7981 0.0095 0.0028 0.0012
0.6361 0.9252 0.0140 -0.0108
-0.7362 -1.9214 0.2743 0.0040
-0.9490 1.7835 -2.3174 -0.0757
Unobservable_Eigenvalues =
-1.7987
0.7065
0.5575
-0.1396
Q = sdpvar(n,n,'symmetric');
N = sdpvar(m,n,'full');
CONS = [
Q - 0.01*eye(n) >= 0
Q*A' + A*Q - B*N - N'*B' + 0.01*eye(n) <= 0
];
sol = optimize(CONS)
Q = double(Q);
N = double(N);
P = inv(Q);
K = N/Q;
Check solution
Eigenvalues_of_the_closed_loop = eig(A - B*K)
Problem
Name :
Objective sense : min
Type : CONIC (conic optimization problem)
Constraints : 63
Cones : 0
Scalar variables : 0
Matrix variables : 2
Integer variables : 0
Optimizer started.
Presolve started.
Linear dependency checker started.
Linear dependency checker terminated.
Eliminator - tries : 0 time : 0.00
Lin. dep. - tries : 1 time : 0.00
Lin. dep. - number : 0
Presolve terminated. Time: 0.00
Problem
Name :
Objective sense : min
Type : CONIC (conic optimization problem)
Constraints : 63
Cones : 0
Scalar variables : 0
Matrix variables : 2
Integer variables : 0
Optimizer - threads : 4
Optimizer - solved problem : the primal
Optimizer - Constraints : 63
Optimizer - Cones : 0
Optimizer - Scalar variables : 0 conic : 0
Optimizer - Semi-definite variables: 2 scalarized : 90
Factor - setup time : 0.00 dense det. time : 0.00
Factor - ML order time : 0.00 GP order time : 0.00
Factor - nonzeros before factor : 2016 after factor : 2016
Factor - dense dim. : 0 flops : 2.16e+05
ITE PFEAS DFEAS GFEAS PRSTATUS POBJ DOBJ MU TIME
0 1.8e+00 1.0e+00 8.2e-01 0.00e+00 -1.800000000e-01 0.000000000e+00 1.0e+00 0.00
1 3.2e-01 1.9e-01 3.0e-01 9.00e-01 -8.174551064e-02 0.000000000e+00 1.8e-01 0.01
2 2.9e-02 1.7e-02 4.2e-02 5.85e-01 -8.305168886e-02 0.000000000e+00 1.7e-02 0.01
3 2.6e-04 1.5e-04 5.0e-05 -6.13e-01 -5.818354522e+00 0.000000000e+00 1.5e-04 0.01
4 1.4e-06 7.7e-07 1.8e-08 -1.00e+00 -1.259788337e+03 0.000000000e+00 7.6e-07 0.01
5 4.1e-08 2.3e-08 8.6e-11 -1.02e+00 -4.801030346e+04 0.000000000e+00 2.3e-08 0.01
6 2.0e-09 1.3e-09 9.0e-13 -1.02e+00 -1.043649970e+06 0.000000000e+00 1.1e-09 0.01
7 6.6e-11 2.8e-09 5.9e-15 -1.01e+00 -2.954177094e+07 0.000000000e+00 3.8e-11 0.01
8 4.7e-12 3.2e-07 1.1e-16 -9.81e-01 -3.953550118e+08 0.000000000e+00 2.7e-12 0.01
Optimizer terminated. Time: 0.02
MOSEK DUAL INFEASIBILITY REPORT.
Problem status: The problem is dual infeasible
Interior-point solution summary
Problem status : DUAL_INFEASIBLE
Solution status : DUAL_INFEASIBLE_CER
Primal. obj: -5.2587724064e-02 nrm: 1e+00 Viol. con: 2e-11 barvar: 0e+00
Optimizer summary
Optimizer - time: 0.02
Interior-point - iterations : 8 time: 0.02
Basis identification - time: 0.00
Primal - iterations : 0 time: 0.00
Dual - iterations : 0 time: 0.00
Clean primal - iterations : 0 time: 0.00
Clean dual - iterations : 0 time: 0.00
Simplex - time: 0.00
Primal simplex - iterations : 0 time: 0.00
Dual simplex - iterations : 0 time: 0.00
Mixed integer - relaxations: 0 time: 0.00
sol =
struct with fields:
yalmiptime: 0.2210
solvertime: 0.0275
info: 'Infeasible problem (MOSEK)'
problem: 1
Eigenvalues_of_the_closed_loop =
-0.5157 + 2.5562i
-0.5157 - 2.5562i
0.7068 + 0.0000i
0.5573 + 0.0000i
-0.6098 + 0.9554i
-0.6098 - 0.9554i
-1.4739 + 0.0000i
-0.6405 + 0.0000i
-0.8166 + 0.0000i