1 | function varargout = iff_internal(varargin) |
---|
2 | X = varargin{1}; |
---|
3 | Y = varargin{2}; |
---|
4 | switch class(varargin{1}) |
---|
5 | |
---|
6 | case 'sdpvar' |
---|
7 | |
---|
8 | if length(X)>1 |
---|
9 | error('IMPLIES not implemented for this case'); |
---|
10 | end |
---|
11 | |
---|
12 | switch class(Y) |
---|
13 | case 'sdpvar' % X <--> Y |
---|
14 | varargout{1} = set(Y == X); |
---|
15 | |
---|
16 | case {'lmi','constraint'} |
---|
17 | Y=set(Y,[],[],1); |
---|
18 | switch settype(Y) |
---|
19 | case 'elementwise' % X <--> Y(:)>=0 |
---|
20 | varargout{1} = binary_iff_lp(X,-sdpvar(Y)); |
---|
21 | case 'equality' % X <--> Y(:)==0 |
---|
22 | varargout{1} = binary_iff_eq(X,sdpvar(Y)); |
---|
23 | otherwise |
---|
24 | error('IFF not implemented for this case'); |
---|
25 | end |
---|
26 | |
---|
27 | otherwise |
---|
28 | error('IFF not implemented for this case'); |
---|
29 | end |
---|
30 | |
---|
31 | case {'lmi','constraint'} |
---|
32 | |
---|
33 | if isa(X,'constraint') |
---|
34 | X = set(X,[],[],1); % FIX: passes one to avoid pruning infeasible constraints |
---|
35 | end |
---|
36 | switch class(Y) |
---|
37 | case 'sdpvar' |
---|
38 | switch settype(X) |
---|
39 | case 'elementwise' |
---|
40 | varargout{1} = binary_iff_lp(Y,-sdpvar(X)); |
---|
41 | case 'equality' |
---|
42 | binvar Z W |
---|
43 | |
---|
44 | % varargout{1} = set(implies(Z&W,Y))+binary_iff_lp(Z,sdpvar(X)+eps)+binary_iff_lp(W,-sdpvar(X)+eps); |
---|
45 | X = [sdpvar(X)-eps;eps-sdpvar(X)]; |
---|
46 | varargout{1} = binary_iff_lp(Y,X);%,sdpvar(X)+eps)+binary_iff_lp(W,-sdpvar(X)+eps); |
---|
47 | % varargout{1} = binary_iff_eq(Y,sdpvar(X)); |
---|
48 | otherwise |
---|
49 | error('IFF not implemented for this case'); |
---|
50 | end |
---|
51 | |
---|
52 | case {'lmi','constraint'} % F(X) <--> F(Y) |
---|
53 | d = binvar(1,1); |
---|
54 | varargout{1} = iff_internal(X,d)+iff_internal(Y,d); |
---|
55 | |
---|
56 | otherwise |
---|
57 | error('IFF not implemented for this case'); |
---|
58 | end |
---|
59 | |
---|
60 | otherwise |
---|
61 | error('IFF not implemented for this case'); |
---|
62 | end |
---|
63 | |
---|
64 | |
---|
65 | function F = binary_iff_eq(X,f) |
---|
66 | [M,m,infbound] = derivebounds(f); |
---|
67 | if infbound |
---|
68 | warning('You have unbounded variables in IFF leading to a lousy big-M relaxation.'); |
---|
69 | end |
---|
70 | eps = 1e-2; |
---|
71 | |
---|
72 | [nf,mf]=size(f); |
---|
73 | if mf>1 |
---|
74 | f = reshape(f,nf*mf,1); |
---|
75 | end |
---|
76 | |
---|
77 | if 1%nf*mf>1 |
---|
78 | di1 = binvar(nf*mf,1); |
---|
79 | di2 = binvar(nf*mf,1); |
---|
80 | F = set(M*(1-X) >= f >= m*(1-X)); |
---|
81 | % F = F + set(-eps + (M+eps).*(1-di) >= f >= eps + (m-eps).*(1-di))+set(X>=sum(di)-length(di)+1); |
---|
82 | % F = F + set(-eps + (M+eps).*di >= f >= eps + (m-eps).*di)+set(X>=sum(di)-length(di)+1); |
---|
83 | |
---|
84 | F = F + set(f>=eps+(m-eps).*di1)+set(-f>=-eps+(-M+eps).*di2)+set(X>=sum(di1)-length(di1)+1)+set(X>=sum(di2)-length(di2)+1); |
---|
85 | |
---|
86 | else |
---|
87 | x1 = binvar(1,1); |
---|
88 | x2 = binvar(1,1); |
---|
89 | F = set(M*(1-x1) + eps >= f >= -eps + m*(1-x2)); |
---|
90 | F = F + set(f >= eps + m.*x1); |
---|
91 | F = F + set(f <= -eps + M.*x2); |
---|
92 | F = F + set(x1+x2-1 <=X); |
---|
93 | |
---|
94 | % F = F + set(iff(~X,~W | ~Z));% >= 1-Z) + set(X >= 1-W); |
---|
95 | % F = F + set(f >= eps + (m-eps).*X)+set(-f >= eps + (-M-eps).*X); |
---|
96 | % F = F + set(f >= eps + (m-eps)*Z)+set(-f >= eps + (-M-eps)*W); |
---|
97 | % F = F + set(X == (Z | W)); |
---|
98 | end |
---|
99 | |
---|
100 | function F = binary_iff_lp(X,f) |
---|
101 | [M,m,infbound] = derivebounds(f); |
---|
102 | if infbound |
---|
103 | warning('You have unbounded variables in IFF leading to a lousy big-M relaxation.'); |
---|
104 | end |
---|
105 | eps = 1e-5; |
---|
106 | |
---|
107 | [nf,mf]=size(f); |
---|
108 | if nf*mf==1 |
---|
109 | F = set(f <= M*(1-X)) + set(f>=eps+(m-eps)*X); |
---|
110 | else |
---|
111 | if mf>1 |
---|
112 | f = reshape(f,nf*mf,1); |
---|
113 | end |
---|
114 | if nf*mf>1 |
---|
115 | di = binvar(nf*mf,1); |
---|
116 | % di=0 means the ith hypeplane is violated |
---|
117 | % X=1 means we are in the polytope |
---|
118 | F = set(f <= M*(1-X)) + set(f>=eps+(m-eps).*di)+set(X>=sum(di)-length(di)+1) + set(X <= di); |
---|
119 | |
---|
120 | % Add some cuts for c < a'x+b < d |
---|
121 | [bA] = getbase(f); |
---|
122 | b = bA(:,1); |
---|
123 | A = bA(:,2:end); |
---|
124 | S = zeros(0,length(di)); |
---|
125 | for i = 1:length(b) |
---|
126 | j = findrows(abs(A),abs(A(i,:))); |
---|
127 | j = j(j > i); |
---|
128 | if length(j)==1 |
---|
129 | S(end+1,[i j]) = 1; |
---|
130 | end |
---|
131 | end |
---|
132 | if size(S,1) > 0 |
---|
133 | % Add cut cannot be outside both constraints |
---|
134 | F = F + set(S*di >= 1); |
---|
135 | end |
---|
136 | else |
---|
137 | F = set(f <= M*(1-X)) + set(f>=eps+(m-eps)*X); |
---|
138 | end |
---|
139 | end |
---|