1 | function pass = akerblad |
---|
2 | |
---|
3 | G=tf(1,[1 1 1]); |
---|
4 | Gd=c2d(G,0.1); |
---|
5 | [A,B]=ssdata(Gd); |
---|
6 | [n,m]=size(B); |
---|
7 | N=100; |
---|
8 | |
---|
9 | Ci=[1 0;0 1; 0 0]; |
---|
10 | Di=[0; 0; 1]; |
---|
11 | p=3; |
---|
12 | q=4; |
---|
13 | Ei=[1 3;-1 -3;0 0; 0 0]; |
---|
14 | Fi=[0;0;1; -1]; |
---|
15 | bi=[10;10;1;1]; |
---|
16 | P=[1 0;0 1;]; |
---|
17 | W=[1 0;0 1;]; |
---|
18 | E=[];F=[];C=[];D=[]; |
---|
19 | ci=[1;1;zeros(p,1)]; |
---|
20 | x0=[2;2]; |
---|
21 | b=[]; |
---|
22 | c=[]; |
---|
23 | for k=1:N, |
---|
24 | E(:,:,k)=Ei; |
---|
25 | F(:,:,k)=Fi; |
---|
26 | C(:,:,k)=Ci; |
---|
27 | D(:,:,k)=Di; |
---|
28 | b=[b;bi]; |
---|
29 | c=[c;ci]; |
---|
30 | end |
---|
31 | c=[c;1;1;zeros(n,1);1;zeros(n,1)]; |
---|
32 | d=sparse([],[],[],(N+1)*n,1); |
---|
33 | d(1:n)=x0; |
---|
34 | H=sparse([],[],[]); |
---|
35 | H(1:n,1+1:n+1)=eye(n); |
---|
36 | M=sparse([],[],[],N*q,N*(n+m+1)+1+n); |
---|
37 | L=sparse([],[],[]); |
---|
38 | |
---|
39 | for k=1:N, |
---|
40 | M(k*q+1-q:k*q,k*(n+m+1)-n-m:k*(n+m+1))=[zeros(q,1) E(:,:,k) F(:,:,k)]; |
---|
41 | L(k*(p+2)+1-p-2:k*(p+2),k*(n+m+1)-n-m:k*(n+m+1))=[-1 zeros(1,n+m);1 zeros(1,n+m);zeros(p,1) -2*C(:,:,k) -2*D(:,:,k)]; |
---|
42 | H(k*n+1:k*n+n,k*(n+m+1)+1-n-m:k*(1+n+m)+n+1)=[-A -B zeros(n,1) ... |
---|
43 | eye(n)]; |
---|
44 | end |
---|
45 | L(k*(p+2)+1:(k)*(p+2)+2+n+1+n,k*(n+m+1)+1:k*(n+m+1)+n+1)=[-1 ... |
---|
46 | zeros(1,n);1 zeros(1,n);zeros(n,1) -2*P;zeros(1,n+ ... |
---|
47 | 1);zeros(n,1) -W]; |
---|
48 | f=sparse([],[],[],N*(n+m+1)+n+1,1); |
---|
49 | f(1:n+m+1:N*(n+m+1)+n+1)=1; |
---|
50 | [yy,ss,rr,sol]=my_socp2(f,d,c,b,L,M,H,N); |
---|
51 | |
---|
52 | norm(yy) |
---|
53 | norm(ss) |
---|
54 | |
---|
55 | |
---|
56 | |
---|
57 | |
---|
58 | |
---|
59 | |
---|
60 | |
---|
61 | |
---|
62 | |
---|
63 | function [y,s,r,sol]=my_socp2(f,b,c,d,A,B,H,N) |
---|
64 | %[t,x,s,r,gamma,lambda]=my_socp(f,c,d,A,B) |
---|
65 | %Solves the primal and dual LP+SOCP on the form |
---|
66 | %min f'y | max -d'gamma-c'lambda-b'rho |
---|
67 | %st s=c-Ay | st -A'lambda-B'gamma-H'rho=f |
---|
68 | % r=d-By | 0<gamma, ||lambda_1||<lambda_0 |
---|
69 | % Hy=b | |
---|
70 | % 0<r, ||s_1||<s_0 | |
---|
71 | [n,m]=size(A); |
---|
72 | [p,m]=size(B); |
---|
73 | [q,m]=size(H); |
---|
74 | |
---|
75 | tic |
---|
76 | yalmip('clear') |
---|
77 | |
---|
78 | y=sdpvar(m,1); |
---|
79 | s = c-A*y; |
---|
80 | r = d-B*y; |
---|
81 | F=lmi('0<r'); |
---|
82 | F=addlmi(F,'H*y=b'); |
---|
83 | for i=1:N |
---|
84 | F=addlmi(F,'||s(i*(2+3)-3:i*5)||<s(i*5-4)'); |
---|
85 | end |
---|
86 | F=addlmi(F,'||s((i+1)*(2+3)-3:(i+1)*5-1)||<s((i+1)*5-4)'); |
---|
87 | F=addlmi(F,'||s(n-1:n)||<s(n-2)'); |
---|
88 | toc |
---|
89 | sol=solvesdp(F,[],f'*y,sdpsettings('sedumi.beta',0.8)) |
---|
90 | |
---|
91 | y=double(y); s=double(s); r=double(r); |
---|
92 | |
---|
93 | |
---|
94 | |
---|