1 | yalmip('clear'); |
---|
2 | echo on |
---|
3 | clc |
---|
4 | %********************************************************* |
---|
5 | % |
---|
6 | % Simultaneuous stabilization using non-convex |
---|
7 | % semidefinite programming. |
---|
8 | % |
---|
9 | %********************************************************* |
---|
10 | % |
---|
11 | % This example solves a simultaneuous stabilization problem |
---|
12 | % |
---|
13 | % The example is taken from |
---|
14 | % "A path-following method for solving BMI problems in control" |
---|
15 | % A. Hassibi, J. How, and S. Boyd. |
---|
16 | % Proceedings of American Control Conference, 2:1385-1389, June 1999. |
---|
17 | % |
---|
18 | % |
---|
19 | % The goal is to find a state-feedback, stabilizing a set |
---|
20 | % of systems, under constraints on the feedback matrix. |
---|
21 | % |
---|
22 | % The performance measure is to maximize the smallest |
---|
23 | % decay-rate for the closed-loop systems. |
---|
24 | % |
---|
25 | % The optimization problem can be stated as: |
---|
26 | % |
---|
27 | % max t |
---|
28 | % subject to (A1+B*K)'*P1+P1*(A1+B*K) < -2*t*P1 |
---|
29 | % (A2+B*K)'*P2+P2*(A2+B*K) < -2*t*P2 |
---|
30 | % (A3+B*K)'*P3+P3*(A3+B*K) < -2*t*P3 |
---|
31 | % -50 < K |
---|
32 | % K < 50 |
---|
33 | % P1,P2,P3 > 0 |
---|
34 | % |
---|
35 | % Clearly, this is a BMI due to the products between Pi |
---|
36 | % and K, and the products t*Pi |
---|
37 | pause |
---|
38 | clc |
---|
39 | |
---|
40 | % Generate system data |
---|
41 | A1 = [1 -1 0;1 1 0;0 0 -0.5]; |
---|
42 | A2 = [1.5 -7 0;7 1.5 0;0 0 1]; |
---|
43 | A3 = [-0.5 -3 0;3 -0.5 0;0 0 2]; |
---|
44 | B = [0.2477 -0.1645;0.4070 0.8115;0.6481 0.4083]; |
---|
45 | pause |
---|
46 | |
---|
47 | % An initial guess is often important when solving BMIs |
---|
48 | % |
---|
49 | % The following set of LMIs gives us a stabilizing controller and |
---|
50 | % a reasonable initial guess |
---|
51 | % (let P1=P2=P3=P, multply BMIs with Q=P^-1, let L=KP^-1, t=0 and |
---|
52 | % neglect amplitude constraints) |
---|
53 | L = sdpvar(2,3); |
---|
54 | Q = sdpvar(3,3); |
---|
55 | F = set(Q>0) + set(trace(Q)==1); % Normalize |
---|
56 | F = F + set(Q*A1'+(B*L)'+(A1*Q+B*L) < 0); |
---|
57 | F = F + set(Q*A2'+(B*L)'+(A2*Q+B*L) < 0); |
---|
58 | F = F + set(Q*A3'+(B*L)'+(A3*Q+B*L) < 0); |
---|
59 | pause |
---|
60 | % Find a feasible solution using penbmi as convex SDP solver |
---|
61 | solvesdp(F,-trace(Q),sdpsettings('solver','penbmi')); |
---|
62 | |
---|
63 | % The solution Q and L can be used later to define initial guesses for P and K |
---|
64 | P0 = inv(double(Q)); |
---|
65 | K0 = double(L)*P0; |
---|
66 | pause |
---|
67 | clc |
---|
68 | % Define variables in actual BMI-problem |
---|
69 | t = sdpvar(1,1); |
---|
70 | K = sdpvar(2,3); |
---|
71 | P1 = sdpvar(3,3); |
---|
72 | P2 = sdpvar(3,3); |
---|
73 | P3 = sdpvar(3,3); |
---|
74 | |
---|
75 | % BMI problem (P>0 changed to P>I for numerical reasons) |
---|
76 | F = set(P1>eye(3)) + set(P2>eye(3))+set(P3>eye(3)); |
---|
77 | F = F + set(-50 < K < 50); |
---|
78 | F = F + set((A1+B*K)'*P1+P1*(A1+B*K) < -2*t*P1); |
---|
79 | F = F + set((A2+B*K)'*P2+P2*(A2+B*K) < -2*t*P2); |
---|
80 | F = F + set((A3+B*K)'*P3+P3*(A3+B*K) < -2*t*P3); |
---|
81 | pause |
---|
82 | |
---|
83 | % Solve with the calculated initial guesses (and t = 0 initially) |
---|
84 | setsdpvar(P1,P0); |
---|
85 | setsdpvar(P2,P0); |
---|
86 | setsdpvar(P3,P0); |
---|
87 | setsdpvar(K,K0); |
---|
88 | setsdpvar(t,0); |
---|
89 | solvesdp(F,-t,sdpsettings('usex0',1)); |
---|
90 | |
---|
91 | % Decay-rate obtained |
---|
92 | double(t) |
---|
93 | % Feedback matrix |
---|
94 | double(K) |
---|
95 | pause |
---|
96 | echo off |
---|