[37] | 1 | function varargout = dilate(F,w) |
---|
| 2 | % DILATE Derives a matrix dilation |
---|
| 3 | % |
---|
| 4 | % [G,H,M] = DILATE(X,w) where X is a symmetric variable derives the |
---|
| 5 | % decomposition and orthogonal complement used in a matrix dilation. |
---|
| 6 | % |
---|
| 7 | % X is decomposed as M(w)ŽG(x)M(w), and H(w) is an orthogonal complement to |
---|
| 8 | % M(w), with affine dependence in w. These matrices can be used to apply |
---|
| 9 | % the matrix dilation lemma to obtain a constraint affine in w |
---|
| 10 | % |
---|
| 11 | % M(w)ŽG(x)M(w) > 0 <==> existence of W s.t with G + W*H' + H*W' > 0 |
---|
| 12 | % |
---|
| 13 | % F = DILATE(F,w) where F is a SET object is used to construct the dilated |
---|
| 14 | % uncertain constraint, i.e an SDP constraint where polynomial dependence |
---|
| 15 | % w.r.t uncertain variables in all SDP constraints in F are converted |
---|
| 16 | % (conservatively) to affine dependence using the matrix dilation approach. |
---|
| 17 | % |
---|
| 18 | % See Yasuaki OISHI, A Region-Dividing Approach to Robust Semidefinite |
---|
| 19 | % Programming and Its Error Bound,DEPARTMENT OF MATHEMATICAL INFORMATICS |
---|
| 20 | % GRADUATE SCHOOL OF INFORMATION SCIENCE AND TECHNOLOGY THE UNIVERSITY OF |
---|
| 21 | % TOKYO BUNKYO-KU, TOKYO 113-8656, JAPAN , February 2006 |
---|
| 22 | % |
---|
| 23 | % See also ROBUSTIFY, SOLVEROBUST, UNCERTAIN |
---|
| 24 | |
---|
| 25 | % Author Johan Löfberg |
---|
| 26 | % $Id: dilate.m,v 1.2 2006/10/25 09:18:47 joloef Exp $ |
---|
| 27 | |
---|
| 28 | if isa(F,'sdpvar') |
---|
| 29 | [G,H,M] = matrix_dilate(F,w) |
---|
| 30 | varargout{1} = G; |
---|
| 31 | varargout{2} = H; |
---|
| 32 | varargout{3} = M; |
---|
| 33 | elseif isa(F,'lmi') |
---|
| 34 | Fnew = []; |
---|
| 35 | if nargin == 1 |
---|
| 36 | w = []; |
---|
| 37 | else |
---|
| 38 | w = w(:); |
---|
| 39 | end |
---|
| 40 | unc_declarations = is(F,'uncertain'); |
---|
| 41 | if any(unc_declarations) |
---|
| 42 | w = [w;recover(getvariables(sdpvar(F(find(unc_declarations)))))]; |
---|
| 43 | end |
---|
| 44 | for i = 1:length(F) |
---|
| 45 | if (is(F(i),'lmi') | ((length(sdpvar(F(i))) == 1) & is(F(i),'elementwise'))) & max(degree(sdpvar(F(i)),w))>0 |
---|
| 46 | [G,H,M] = matrix_dilate(sdpvar(F(i)),w); |
---|
| 47 | W = sdpvar(size(G,1),size(H,2)); |
---|
| 48 | Fnew = Fnew + set(G + W*H' + H*W' >= 0); |
---|
| 49 | else |
---|
| 50 | Fnew = Fnew + F(i); |
---|
| 51 | end |
---|
| 52 | end |
---|
| 53 | varargout{1} = Fnew; |
---|
| 54 | end |
---|
| 55 | |
---|
| 56 | |
---|
| 57 | function [G,H,M] = matrix_dilate(F,w) |
---|
| 58 | % Given an SDP F(x,w) with F polynomial in w, DILATE rewrites the problem |
---|
| 59 | % to F(x,w)=M(w)'G(x)M(w), and derives the orhogonal complements H(w) to |
---|
| 60 | % M(w), to be used in the dilated constraint G + W*H' + H*W' |
---|
| 61 | |
---|
| 62 | allvariables = getvariables(F); |
---|
| 63 | wvariables = getvariables(w); |
---|
| 64 | Fbasis = getbase(F); |
---|
| 65 | |
---|
| 66 | % Make sure it is ordered according to internal index |
---|
| 67 | w = recover(wvariables); |
---|
| 68 | |
---|
| 69 | % Degrees w.r.t the uncertain variables |
---|
| 70 | d = degree(recover(allvariables),w); |
---|
| 71 | |
---|
| 72 | % robustify code by removing unused variables |
---|
| 73 | w = w(find(d)); |
---|
| 74 | wvariables = wvariables(find(d)); |
---|
| 75 | |
---|
| 76 | % Degrees w.r.t the uncertain variables |
---|
| 77 | monomtable = yalmip('monomtable'); |
---|
| 78 | d = max(sum(monomtable(allvariables,wvariables),2)); |
---|
| 79 | n = size(F,1); |
---|
| 80 | |
---|
| 81 | % Sufficiently many monomials |
---|
| 82 | v = monolist(w,max(d)); |
---|
| 83 | % Some numerical format on these variables |
---|
| 84 | for i = 2:length(v) |
---|
| 85 | vvariables(i,1) = getvariables(v(i)); |
---|
| 86 | end |
---|
| 87 | monomtable = yalmip('monomtable'); |
---|
| 88 | wmonoms = [zeros(1,length(w));monomtable(vvariables(2:end),wvariables)]; |
---|
| 89 | |
---|
| 90 | % Find the linear certain term in the matrix |
---|
| 91 | linear_indicies = find(sum(monomtable(allvariables,wvariables),2) ==0); |
---|
| 92 | F0 = reshape(Fbasis(:,1) + Fbasis(:,1+linear_indicies)*recover(allvariables(linear_indicies)),n,n); |
---|
| 93 | |
---|
| 94 | % now find the matrix that multiplies with each monomial in M |
---|
| 95 | Fi = []; |
---|
| 96 | Fmonoms = monomtable(allvariables,wvariables); |
---|
| 97 | for i = 2:length(v) |
---|
| 98 | vmonoms = monomtable(vvariables(i),wvariables); |
---|
| 99 | index = findrows(Fmonoms,vmonoms); |
---|
| 100 | if isempty(index) |
---|
| 101 | Fi = [Fi zeros(n)]; |
---|
| 102 | else |
---|
| 103 | temp = monomtable(allvariables(index),:); |
---|
| 104 | temp(wvariables) = 0; |
---|
| 105 | base = reshape(Fbasis(:,1+index),n,n); |
---|
| 106 | if nnz(temp) == 0 |
---|
| 107 | Fi = [Fi base]; |
---|
| 108 | else |
---|
| 109 | Fi = [Fi base*recover(allvariables(find(temp)))]; |
---|
| 110 | end |
---|
| 111 | end |
---|
| 112 | end |
---|
| 113 | G = [F0 Fi/2;Fi'/2 zeros(n*(length(v)-1))]; |
---|
| 114 | |
---|
| 115 | % The outer factor |
---|
| 116 | M = kron(v,eye(n)); |
---|
| 117 | |
---|
| 118 | % Now create an orthogonal complement to M |
---|
| 119 | ii = []; |
---|
| 120 | jj = []; |
---|
| 121 | ss = []; |
---|
| 122 | for i = 2:length(v) |
---|
| 123 | monom = wmonoms(i,:); |
---|
| 124 | [dummy,index] = max(monom); |
---|
| 125 | monomnew = monom; |
---|
| 126 | monomnew(index) = monomnew(index) - 1; |
---|
| 127 | ii = [ii i]; |
---|
| 128 | jj = [jj i-1]; |
---|
| 129 | ss = [ss 1]; |
---|
| 130 | ii = [ii findrows(wmonoms,monomnew)]; |
---|
| 131 | jj = [jj i-1]; |
---|
| 132 | ss = [ss -recover(wvariables(index))]; |
---|
| 133 | end |
---|
| 134 | H = kron(sparse(ii,jj,ss),eye(n)); |
---|