[37] | 1 | clc |
---|
| 2 | echo on |
---|
| 3 | %********************************************************* |
---|
| 4 | % |
---|
| 5 | % Global bilinear programming |
---|
| 6 | % |
---|
| 7 | %********************************************************* |
---|
| 8 | % |
---|
| 9 | % YALMIP comes with a simple branch-and-bound code |
---|
| 10 | % for solving problems with bilinear/quadratic constraints. |
---|
| 11 | % |
---|
| 12 | % The code is under development, and is currently only |
---|
| 13 | % applicable to small problems. |
---|
| 14 | % |
---|
| 15 | % For the examples below to work, you must have an efficient |
---|
| 16 | % LP solver (cplex,clp,nag,...) insrtalled and a solver for |
---|
| 17 | % general polynomial problems (currently only fmincon or PENBMI) |
---|
| 18 | |
---|
| 19 | |
---|
| 20 | pause |
---|
| 21 | yalmip('clear') |
---|
| 22 | clc |
---|
| 23 | %********************************************************* |
---|
| 24 | % The first example is the problem from the moment |
---|
| 25 | % relaxation demo (concave quadratic constraint). |
---|
| 26 | %********************************************************* |
---|
| 27 | x1 = sdpvar(1,1); |
---|
| 28 | x2 = sdpvar(1,1); |
---|
| 29 | x3 = sdpvar(1,1); |
---|
| 30 | |
---|
| 31 | objective = -2*x1+x2-x3; |
---|
| 32 | |
---|
| 33 | F = set(x1*(4*x1-4*x2+4*x3-20)+x2*(2*x2-2*x3+9)+x3*(2*x3-13)+24>0); |
---|
| 34 | F = F + set(4-(x1+x2+x3)>0); |
---|
| 35 | F = F + set(6-(3*x2+x3)>0); |
---|
| 36 | F = F + set(x1>0); |
---|
| 37 | F = F + set(2-x1>0); |
---|
| 38 | F = F + set(x2>0); |
---|
| 39 | F = F + set(x3>0); |
---|
| 40 | F = F + set(3-x3>0); |
---|
| 41 | pause |
---|
| 42 | |
---|
| 43 | % Explicitly specify the global solver (this is currently recommended) |
---|
| 44 | ops = sdpsettings('solver','bmibnb'); % Global solver |
---|
| 45 | pause |
---|
| 46 | |
---|
| 47 | % Solve! |
---|
| 48 | solvesdp(F,objective,ops) |
---|
| 49 | pause |
---|
| 50 | yalmip('clear') |
---|
| 51 | clc |
---|
| 52 | %*********************************************************** |
---|
| 53 | % The second example is slightly larger. It is a problem |
---|
| 54 | % with a concave quadratic objective function, and two |
---|
| 55 | % concave quadratic constraints. |
---|
| 56 | %*********************************************************** |
---|
| 57 | x1 = sdpvar(1,1); |
---|
| 58 | x2 = sdpvar(1,1); |
---|
| 59 | x3 = sdpvar(1,1); |
---|
| 60 | x4 = sdpvar(1,1); |
---|
| 61 | x5 = sdpvar(1,1); |
---|
| 62 | x6 = sdpvar(1,1); |
---|
| 63 | |
---|
| 64 | objective = -25*(x1-2)^2-(x2-2)^2-(x3-1)^2-(x4-4)^2-(x5-1)^2-(x6-4)^2; |
---|
| 65 | |
---|
| 66 | F = set((x3-3)^2+x4>4) + set((x5-3)^2+x6>4); |
---|
| 67 | F = F + set(x1-3*x2<2) + set(-x1+x2<2) + set(x1+x2>2); |
---|
| 68 | F = F + set(6>x1+x2>2); |
---|
| 69 | F = F + set(1<x3<5) + set(0<x4<6)+set(1<x5<5)+set(0<x6<10)+set(x1>0)+set(x2>0); |
---|
| 70 | |
---|
| 71 | pause |
---|
| 72 | |
---|
| 73 | % The global solver in YALMIP can currently only deal with linear objectives. |
---|
| 74 | % A simple epi-graph formulation solves this. |
---|
| 75 | t = sdpvar(1,1); |
---|
| 76 | F = F + set(objective<t); |
---|
| 77 | pause |
---|
| 78 | |
---|
| 79 | % Solve! |
---|
| 80 | solvesdp(F,t,ops); |
---|
| 81 | pause |
---|
| 82 | clc |
---|
| 83 | %*********************************************************** |
---|
| 84 | % Random bound-constrained indefinte quadratic program. |
---|
| 85 | %*********************************************************** |
---|
| 86 | |
---|
| 87 | n = 5; |
---|
| 88 | x = sdpvar(n,1); |
---|
| 89 | t = sdpvar(1,1); |
---|
| 90 | |
---|
| 91 | A = randn(n,n);A = A*A'; |
---|
| 92 | B = randn(n,n);B = B*B'; |
---|
| 93 | |
---|
| 94 | objective = x'*A*x-x'*B*x; |
---|
| 95 | |
---|
| 96 | F = set(-1 < x < 1) + set(objective<t); |
---|
| 97 | pause |
---|
| 98 | |
---|
| 99 | % Solve! |
---|
| 100 | solvesdp(F,t,ops); |
---|
| 101 | pause |
---|
| 102 | yalmip('clear') |
---|
| 103 | clc |
---|
| 104 | %*********************************************************** |
---|
| 105 | % Random boolean quadratic programming. |
---|
| 106 | % |
---|
| 107 | % By adding a nonlinear constraint x(x-1)==0, we can solve |
---|
| 108 | % boolean programming using global nonlinear programming. |
---|
| 109 | %*********************************************************** |
---|
| 110 | |
---|
| 111 | n = 10; |
---|
| 112 | x = sdpvar(n,1); |
---|
| 113 | t = sdpvar(1,1); |
---|
| 114 | |
---|
| 115 | Q = randn(n,n);Q = Q*Q'/norm(Q)^2; |
---|
| 116 | c = randn(n,1); |
---|
| 117 | |
---|
| 118 | objective = x'*Q*x+c'*x; |
---|
| 119 | |
---|
| 120 | F = set(x.*(x-1)==0) + set(objective < t); |
---|
| 121 | pause |
---|
| 122 | |
---|
| 123 | % Note, all complicating variables (variables involved in nonlinear terms) |
---|
| 124 | % must be bounded (either explicitely or implicetely via linear and relaxed |
---|
| 125 | % nonlinear constraints). |
---|
| 126 | % |
---|
| 127 | % The easible set for the relaxed problem is not bounded |
---|
| 128 | % so YALMIP will not be able to derive variable bounds. |
---|
| 129 | % We add some valid bounds to fix this |
---|
| 130 | F = F + set(0<x<1); |
---|
| 131 | pause |
---|
| 132 | |
---|
| 133 | % Solve! |
---|
| 134 | solvesdp(F,t,ops); |
---|
| 135 | double(x'*Q*x+c'*x) |
---|
| 136 | pause |
---|
| 137 | |
---|
| 138 | % Let us compare with an integer solver (use we use the internal MIQP solver) |
---|
| 139 | solvesdp(set(binary(x)),x'*Q*x+c'*x,sdpsettings(ops,'solver','bnb')); |
---|
| 140 | double(x'*Q*x+c'*x) |
---|
| 141 | pause |
---|
| 142 | |
---|
| 143 | yalmip('clear') |
---|
| 144 | clc |
---|
| 145 | %*********************************************************** |
---|
| 146 | % This example is a somewhat larger indefinite quadratic |
---|
| 147 | % programming problem, taken from the GAMS problem library. |
---|
| 148 | %*********************************************************** |
---|
| 149 | pause |
---|
| 150 | |
---|
| 151 | % The problem has 11 variables |
---|
| 152 | x1 = sdpvar(1,1); |
---|
| 153 | x2 = sdpvar(1,1); |
---|
| 154 | x3 = sdpvar(1,1); |
---|
| 155 | x4 = sdpvar(1,1); |
---|
| 156 | x5 = sdpvar(1,1); |
---|
| 157 | x6 = sdpvar(1,1); |
---|
| 158 | x7 = sdpvar(1,1); |
---|
| 159 | x8 = sdpvar(1,1); |
---|
| 160 | x9 = sdpvar(1,1); |
---|
| 161 | x10 = sdpvar(1,1); |
---|
| 162 | x11 = sdpvar(1,1); |
---|
| 163 | t = sdpvar(1,1); |
---|
| 164 | x = [x1;x2;x3;x4;x5;x6;x7;x8;x9;x10;x11]; |
---|
| 165 | |
---|
| 166 | pause |
---|
| 167 | |
---|
| 168 | % ...and 22 linear constraints |
---|
| 169 | e1 = - x1 - 2*x2 - 3*x3 - 4*x4 - 5*x5 - 6*x6 - 7*x7 - 8*x8 - 9*x9 - 10*x10 - 11*x11 - 0; |
---|
| 170 | e2 = - 2*x1 - 3*x2 - 4*x3 - 5*x4 - 6*x5 - 7*x6 - 8*x7 - 9*x8 - 10*x9 - 11*x10 - x11 - 0; |
---|
| 171 | e3 = - 3*x1 - 4*x2 - 5*x3 - 6*x4 - 7*x5 - 8*x6 - 9*x7 - 10*x8 - 11*x9 - x10 - 2*x11 - 0; |
---|
| 172 | e4 = - 4*x1 - 5*x2 - 6*x3 - 7*x4 - 8*x5 - 9*x6 - 10*x7 - 11*x8 - x9 - 2*x10 - 3*x11 - 0; |
---|
| 173 | e5 = - 5*x1 - 6*x2 - 7*x3 - 8*x4 - 9*x5 - 10*x6 - 11*x7 - x8 - 2*x9 - 3*x10 - 4*x11 - 0; |
---|
| 174 | e6 = - 6*x1 - 7*x2 - 8*x3 - 9*x4 - 10*x5 - 11*x6 - x7 - 2*x8 - 3*x9 - 4*x10 - 5*x11 - 0; |
---|
| 175 | e7 = - 7*x1 - 8*x2 - 9*x3 - 10*x4 - 11*x5 - x6 - 2*x7 - 3*x8 - 4*x9 - 5*x10 - 6*x11 - 0; |
---|
| 176 | e8 = - 8*x1 - 9*x2 - 10*x3 - 11*x4 - x5 - 2*x6 - 3*x7 - 4*x8 - 5*x9 - 6*x10 - 7*x11 - 0; |
---|
| 177 | e9 = - 9*x1 - 10*x2 - 11*x3 - x4 - 2*x5 - 3*x6 - 4*x7 - 5*x8 - 6*x9 - 7*x10 - 8*x11 - 0; |
---|
| 178 | e10 = - 10*x1 - 11*x2 - x3 - 2*x4 - 3*x5 - 4*x6 - 5*x7 - 6*x8 - 7*x9 - 8*x10 - 9*x11 - 0; |
---|
| 179 | e11 = - 11*x1 - x2 - 2*x3 - 3*x4 - 4*x5 - 5*x6 - 6*x7 - 7*x8 - 8*x9 - 9*x10 - 10*x11 - 0; |
---|
| 180 | e12 = x1 + 2*x2 + 3*x3 + 4*x4 + 5*x5 + 6*x6 + 7*x7 + 8*x8 + 9*x9 + 10*x10 + 11*x11 - 66; |
---|
| 181 | e13 = 2*x1 + 3*x2 + 4*x3 + 5*x4 + 6*x5 + 7*x6 + 8*x7 + 9*x8 + 10*x9 + 11*x10 + x11 - 66; |
---|
| 182 | e14 = 3*x1 + 4*x2 + 5*x3 + 6*x4 + 7*x5 + 8*x6 + 9*x7 + 10*x8 + 11*x9 + x10 + 2*x11 - 66; |
---|
| 183 | e15 = 4*x1 + 5*x2 + 6*x3 + 7*x4 + 8*x5 + 9*x6 + 10*x7 + 11*x8 + x9 + 2*x10 + 3*x11 - 66; |
---|
| 184 | e16 = 5*x1 + 6*x2 + 7*x3 + 8*x4 + 9*x5 + 10*x6 + 11*x7 + x8 + 2*x9 + 3*x10 + 4*x11 - 66; |
---|
| 185 | e17 = 6*x1 + 7*x2 + 8*x3 + 9*x4 + 10*x5 + 11*x6 + x7 + 2*x8 + 3*x9 + 4*x10 + 5*x11 - 66; |
---|
| 186 | e18 = 7*x1 + 8*x2 + 9*x3 + 10*x4 + 11*x5 + x6 + 2*x7 + 3*x8 + 4*x9 + 5*x10 + 6*x11 - 66; |
---|
| 187 | e19 = 8*x1 + 9*x2 + 10*x3 + 11*x4 + x5 + 2*x6 + 3*x7 + 4*x8 + 5*x9 + 6*x10 + 7*x11 - 66; |
---|
| 188 | e20 = 9*x1 + 10*x2 + 11*x3 + x4 + 2*x5 + 3*x6 + 4*x7 + 5*x8 + 6*x9 + 7*x10 + 8*x11 - 66; |
---|
| 189 | e21 = 10*x1 + 11*x2 + x3 + 2*x4 + 3*x5 + 4*x6 + 5*x7 + 6*x8 + 7*x9 + 8*x10 + 9*x11 - 66; |
---|
| 190 | e22 = 11*x1 + x2 + 2*x3 + 3*x4 + 4*x5 + 5*x6 + 6*x7 + 7*x8 + 8*x9 + 9*x10 + 10*x11 - 66; |
---|
| 191 | |
---|
| 192 | pause |
---|
| 193 | |
---|
| 194 | |
---|
| 195 | % Define the objective function and the whole program |
---|
| 196 | obj = (0.5*x1*x2 - x1*x1 + 0.5*x2*x1 - x2*x2 + 0.5*x2*x3 + 0.5*x3*x2 - x3*x3 + 0.5*x3*x4 + 0.5*x4*x3 - x4*x4 + 0.5*x4*x5 + 0.5*x5*x4 - x5*x5 + 0.5*x5 *x6 + 0.5*x6*x5 - x6*x6 + 0.5*x6*x7 + 0.5*x7*x6 - x7*x7 + 0.5*x7*x8 + 0.5 *x8*x7 - x8*x8 + 0.5*x8*x9 + 0.5*x9*x8 - x9*x9 + 0.5*x9*x10 + 0.5*x10*x9 - x10*x10 + 0.5*x10*x11 + 0.5*x11*x10 - x11*x11); |
---|
| 197 | e = -[e1;e2;e3;e4;e5;e6;e7;e8;e9;e10;e11;e12;e13;e14;e15;e16;e17;e18;e19;e20;e21;e22]; |
---|
| 198 | F = set(10>x>0); |
---|
| 199 | F = F + set(e>0); |
---|
| 200 | F = F + set(obj==t); |
---|
| 201 | |
---|
| 202 | % We will nowe try to solve this using the same strategy as before |
---|
| 203 | solvesdp(F,t,ops) |
---|
| 204 | |
---|
| 205 | pause |
---|
| 206 | clc |
---|
| 207 | |
---|
| 208 | % As an alternative, we might use additional cuts using the CUT functionality. |
---|
| 209 | % These additional constraints are only used in domain reduction and |
---|
| 210 | % solution of lower bound programs, but not in the local solver. |
---|
| 211 | |
---|
| 212 | % We square all linear constraints and add these as cuts |
---|
| 213 | % |
---|
| 214 | % Note that the triu operator leaves returns a matrix with |
---|
| 215 | % zeros below the diagonal. These are however neglected by YALMIP |
---|
| 216 | % in the call to SET. |
---|
| 217 | % |
---|
| 218 | f = sdpvar(F(find(is(F,'linear')))); |
---|
| 219 | F = F + cut(triu(f*f') > 0); |
---|
| 220 | pause |
---|
| 221 | |
---|
| 222 | % ...and solve (To speed up things, we turn off LP based domain reduction) |
---|
| 223 | % |
---|
| 224 | % NOTE : The first ~15 iterations are rather slow, but things starts to |
---|
| 225 | % speed up when YALMIP starts pruning redundnant linear constraints. |
---|
| 226 | % |
---|
| 227 | solvesdp(F,t,sdpsettings(ops,'bmibnb.lpreduce',0)) |
---|
| 228 | pause |
---|
| 229 | clc |
---|
| 230 | |
---|
| 231 | % The global solver in YALMIP can only solve bilinear programs, |
---|
| 232 | % but a pre-processor is capable of transforming higher order problems |
---|
| 233 | % to bilinear programs. As an example, the variable x^3y^2 will be replaced |
---|
| 234 | % with the the variable w and the constraints w == uv, u == zx, z == x^2, v == y^2. |
---|
| 235 | % The is done automatically if the global solver, or PENBMI, is called with |
---|
| 236 | % a higher order polynomial problem. Note that this conversion is rather |
---|
| 237 | % inefficient, and only very small problems can be addressed using this simple approach. |
---|
| 238 | % Additionally, this module has not been tested in any serious sense. |
---|
| 239 | pause |
---|
| 240 | |
---|
| 241 | % Let us solve a small trivial polynomial program |
---|
| 242 | sdpvar x y |
---|
| 243 | F = set(x^3+y^5 < 5) + set(y > 0); |
---|
| 244 | F = F + set(-100 < [x y] < 100) % Always bounded domain |
---|
| 245 | pause |
---|
| 246 | |
---|
| 247 | solvesdp(F,-x,ops) |
---|
| 248 | pause |
---|
| 249 | |
---|
| 250 | %*********************************************************** |
---|
| 251 | % The main motivation for implementing a bilinear solver in |
---|
| 252 | % YALMIP is matrix-valued problems, i.e. BMI problems. |
---|
| 253 | % |
---|
| 254 | % Of-course, BMI problems are even harder than the nonconvex |
---|
| 255 | % quadratic problems solved above. However, in some cases, |
---|
| 256 | % it is actually possible to solve BMI problems globally. |
---|
| 257 | % |
---|
| 258 | % Don't expect too much though... |
---|
| 259 | %*********************************************************** |
---|
| 260 | |
---|
| 261 | pause |
---|
| 262 | clc |
---|
| 263 | |
---|
| 264 | %*********************************************************** |
---|
| 265 | % The following is a classical problem |
---|
| 266 | %*********************************************************** |
---|
| 267 | A0 = [-10 -0.5 -2;-0.5 4.5 0;-2 0 0]; |
---|
| 268 | A1 = [9 0.5 0;0.5 0 -3 ; 0 -3 -1]; |
---|
| 269 | A2 = [-1.8 -0.1 -0.4;-0.1 1.2 -1;-0.4 -1 0]; |
---|
| 270 | K12 = [0 0 2;0 -5.5 3;2 3 0]; |
---|
| 271 | |
---|
| 272 | x = sdpvar(1,1); |
---|
| 273 | y = sdpvar(1,1); |
---|
| 274 | t = sdpvar(1,1); |
---|
| 275 | |
---|
| 276 | F = set(A0+x*A1+y*A2+x*y*K12-t*eye(3)<0); |
---|
| 277 | F = F + set(-0.5 < x < 2) + set(-3 < y < 7); |
---|
| 278 | pause |
---|
| 279 | |
---|
| 280 | % Specify solvers (this requires PENBMI!) |
---|
| 281 | ops = sdpsettings('solver','bmibnb'); % Global solver |
---|
| 282 | ops = sdpsettings(ops,'bmibnb.uppersolver','penbmi'); % Local solver |
---|
| 283 | |
---|
| 284 | pause |
---|
| 285 | |
---|
| 286 | % Solve! |
---|
| 287 | solvesdp(F,t,ops); |
---|
| 288 | pause |
---|
| 289 | clc |
---|
| 290 | %*********************************************************** |
---|
| 291 | % The next example shows how to solve a standard LQ control |
---|
| 292 | % problem with constraints on the feedback gain |
---|
| 293 | %*********************************************************** |
---|
| 294 | |
---|
| 295 | A = [-1 2;-3 -4];B = [1;1]; |
---|
| 296 | |
---|
| 297 | P = sdpvar(2,2); |
---|
| 298 | K = sdpvar(1,2); |
---|
| 299 | |
---|
| 300 | F = set(P>0)+set((A+B*K)'*P+P*(A+B*K) < -eye(2)-K'*K); |
---|
| 301 | F = F + set(K<0.1)+set(K>-0.1); |
---|
| 302 | F = F + set(1000> P(:)>-1000); |
---|
| 303 | solvesdp(F,trace(P),ops); |
---|
| 304 | pause |
---|
| 305 | clc |
---|
| 306 | %*********************************************************** |
---|
| 307 | % Decay-rate estimation as a global BMI example ... |
---|
| 308 | % |
---|
| 309 | % Note, this is a quasi-convex problem that PENBMI actually |
---|
| 310 | % solves to global optimality, so the use of a global solver |
---|
| 311 | % is pretty stupid. It is only included here to show some |
---|
| 312 | % properties of the global solver. |
---|
| 313 | % |
---|
| 314 | % may run for up to 100 iterations depending on the solvers |
---|
| 315 | %*********************************************************** |
---|
| 316 | |
---|
| 317 | A = [-1 2;-3 -4]; |
---|
| 318 | t = sdpvar(1,1); |
---|
| 319 | P = sdpvar(2,2); |
---|
| 320 | F = set(P>eye(2))+set(A'*P+P*A < -2*t*P); |
---|
| 321 | F = F + set(-1e4 < P(:) < 1e4) + set(100 > t > 0) ; |
---|
| 322 | pause |
---|
| 323 | solvesdp(F,-t,ops); |
---|
| 324 | pause |
---|
| 325 | |
---|
| 326 | % Now, that sucked, right! 100 iterations and still 10% gap... |
---|
| 327 | % |
---|
| 328 | % The way we model a problem can make a huge difference. |
---|
| 329 | % The decay-rate problem can be substantially improved by |
---|
| 330 | % using a smarter model. The original problem is homogenous |
---|
| 331 | % w.r.t P, and a better way to make the feasible set bounded |
---|
| 332 | % is to constrain the trace of P |
---|
| 333 | A = [-1 2;-3 -4]; |
---|
| 334 | t = sdpvar(1,1); |
---|
| 335 | P = sdpvar(2,2); |
---|
| 336 | F = set(P>0)+set(A'*P+P*A < -2*t*P); |
---|
| 337 | F = F + set(trace(P)==1) + set(100 > t > 0) ; |
---|
| 338 | pause |
---|
| 339 | solvesdp(F,-t,ops); |
---|
| 340 | pause |
---|
| 341 | |
---|
| 342 | % Nothing prevents us from adding an additional |
---|
| 343 | % constraint derived from trace(P)=1 |
---|
| 344 | F = F + set(trace(A'*P+P*A) < -2*t); |
---|
| 345 | pause |
---|
| 346 | |
---|
| 347 | solvesdp(F,-t,ops); |
---|
| 348 | pause |
---|
| 349 | |
---|
| 350 | % When we add redundant constraints, we improve the relaxations |
---|
| 351 | % and obtain better lower bounds. For the upper bounds however, |
---|
| 352 | % these cuts only add additional computational burden. |
---|
| 353 | % To overcome this, the command CUT can be used. |
---|
| 354 | % Constraints defined with CUT are not used when the upper |
---|
| 355 | % bound problems are solved, but are only used in relaxations |
---|
| 356 | F = set(P>0)+set(A'*P+P*A < -2*t*P); |
---|
| 357 | F = F + set(trace(P)==1) + set(100 > t > 0) ; |
---|
| 358 | F = F + cut(trace(A'*P+P*A) < -2*t); |
---|
| 359 | pause |
---|
| 360 | |
---|
| 361 | solvesdp(F,-t,ops); |
---|
| 362 | pause |
---|
| 363 | |
---|
| 364 | % Finally, it should be mentioned that the branch & bound |
---|
| 365 | % algorithm can run without any installed local BMI solver. |
---|
| 366 | % This version of the algorithms is obtained by specifying |
---|
| 367 | % the upper bound solver as 'none' |
---|
| 368 | ops = sdpsettings(ops,'bmibnb.uppersolver','none'); |
---|
| 369 | F = set(P>0)+set(A'*P+P*A < -2*t*P); |
---|
| 370 | F = F + set(trace(P)==1) + set(100 > t > 0) ; |
---|
| 371 | F = F + cut(trace(A'*P+P*A) < -2*t); |
---|
| 372 | pause |
---|
| 373 | |
---|
| 374 | solvesdp(F,-t,ops); |
---|
| 375 | pause |
---|
| 376 | |
---|
| 377 | echo off |
---|
| 378 | |
---|
| 379 | |
---|
| 380 | |
---|
| 381 | break |
---|
| 382 | |
---|
| 383 | % ******************************** |
---|
| 384 | % Sahinidis (-6.666 in 1) |
---|
| 385 | % ******************************** |
---|
| 386 | yalmip('clear'); |
---|
| 387 | x1 = sdpvar(1,1); |
---|
| 388 | x2 = sdpvar(1,1); |
---|
| 389 | F = set(0 < x1 < 6)+set(0 < x2 < 4) + set(x1*x2 < 4) |
---|
| 390 | solvesdp(F,-x1-x2,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bmibnb.maxiter',100,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','cdd','cdd.method','dual-simplex')) |
---|
| 391 | |
---|
| 392 | % **************************** |
---|
| 393 | % Decay-rate (-2.5) |
---|
| 394 | % **************************** |
---|
| 395 | yalmip('clear'); |
---|
| 396 | A = [-1 2;-3 -4]; |
---|
| 397 | t = sdpvar(1,1); |
---|
| 398 | P = sdpvar(2,2); |
---|
| 399 | F = set(P>eye(2))+set(A'*P+P*A < -2*t*P); |
---|
| 400 | F = F + set(-1e4 < P(:) < 1e4) + set(100 > t > -100) ; |
---|
| 401 | ops = sdpsettings('bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bmibnb.maxiter',100,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','cdd','cdd.method','dual-simplex'); |
---|
| 402 | solvesdp(F,-t,ops); |
---|
| 403 | |
---|
| 404 | % ******************************** |
---|
| 405 | % Classical example from Safonov etc |
---|
| 406 | % ******************************** |
---|
| 407 | yalmip('clear') |
---|
| 408 | x = sdpvar(1,1); |
---|
| 409 | y = sdpvar(1,1); |
---|
| 410 | t = sdpvar(1,1); |
---|
| 411 | A0 = [-10 -0.5 -2;-0.5 4.5 0;-2 0 0]; |
---|
| 412 | A1 = [9 0.5 0;0.5 0 -3 ; 0 -3 -1]; |
---|
| 413 | A2 = [-1.8 -0.1 -0.4;-0.1 1.2 -1;-0.4 -1 0]; |
---|
| 414 | K12 = [0 0 2;0 -5.5 3;2 3 0]; |
---|
| 415 | F = lmi(x>-0.5)+lmi(x<2) + lmi(y>-3)+lmi(y<7) + lmi(t>-1000)+lmi(t<1000);%set(x+y<2.47)+set(x+y>2.47); |
---|
| 416 | F = F + lmi(A0+x*A1+y*A2+x*y*K12-t*eye(3)<0); |
---|
| 417 | solvesdp(F,t,sdpsettings('bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.09e-1,'bmibnb.lpsolver','glpk')) |
---|
| 418 | |
---|
| 419 | % ******************************** |
---|
| 420 | % Lassere 1 |
---|
| 421 | % optimal : -4, takes 13 |
---|
| 422 | % ******************************** |
---|
| 423 | clear all |
---|
| 424 | x1 = sdpvar(1,1); |
---|
| 425 | x2 = sdpvar(1,1); |
---|
| 426 | x3 = sdpvar(1,1); |
---|
| 427 | x = [x1;x2;x3]; |
---|
| 428 | B = [0 0 1;0 -1 0;-2 1 -1]; |
---|
| 429 | b = [3;0;-4]; |
---|
| 430 | v = [0;-1;-6]; |
---|
| 431 | r = [1.5;-0.5;-5]; |
---|
| 432 | p = -2*x1+x2-x3; |
---|
| 433 | F = set(x1+x2+x3 < 4) + set(x1<2)+set(x3<3) + set(3*x2+x3 < 6); |
---|
| 434 | F = F + set(x>0) + set(x'*B'*B*x-2*r'*B*x+r'*r-0.25*(b-v)'*(b-v) >0); |
---|
| 435 | solvesdp(F,p,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.09e-1,'bmibnb.lpsolver','glpk')) |
---|
| 436 | |
---|
| 437 | % ******************************** |
---|
| 438 | % Lassere 2 |
---|
| 439 | % optimal : -310 in 3 |
---|
| 440 | % ******************************** |
---|
| 441 | clear all |
---|
| 442 | x1 = sdpvar(1,1); |
---|
| 443 | x2 = sdpvar(1,1); |
---|
| 444 | x3 = sdpvar(1,1); |
---|
| 445 | x4 = sdpvar(1,1); |
---|
| 446 | x5 = sdpvar(1,1); |
---|
| 447 | x6 = sdpvar(1,1); |
---|
| 448 | t = sdpvar(1,1); |
---|
| 449 | p = -25*(x1-2)^2-(x2-2)^2-(x3-1)^2-(x4-4)^2-(x5-1)^2-(x6-4)^2; |
---|
| 450 | F = set((x3-3)^2+x4>4)+set((x5-3)^2+x6>4); |
---|
| 451 | F = F + set(x1-3*x2<2)+set(-x1+x2<2); |
---|
| 452 | F = F + set(x1-3*x2 <2)+set(x1+x2>2); |
---|
| 453 | F = F + set(6>x1+x2>2); |
---|
| 454 | F = F + set(1<x3<5) + set(0<x4<6)+set(1<x5<5)+set(0<x6<10)+set(x1>0)+set(x2>0); |
---|
| 455 | |
---|
| 456 | solvesdp(F,p,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.09e-1,'bmibnb.lpsolver','glpk')) |
---|
| 457 | F = F + set(p<t); |
---|
| 458 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.09e-1,'bmibnb.lpsolver','glpk')) |
---|
| 459 | % ******************************** |
---|
| 460 | % Lassere 2.6 |
---|
| 461 | % optimal : -39 in 6 |
---|
| 462 | % ******************************** |
---|
| 463 | clear all |
---|
| 464 | Q = 100*eye(10); |
---|
| 465 | c = [48, 42, 48, 45, 44, 41, 47, 42, 45, 46]'; |
---|
| 466 | b = [-4, 22,-6,-23,-12]'; |
---|
| 467 | A =[-2 -6 -1 0 -3 -3 -2 -6 -2 -2; |
---|
| 468 | 6 -5 8 -3 0 1 3 8 9 -3; |
---|
| 469 | -5 6 5 3 8 -8 9 2 0 -9; |
---|
| 470 | 9 5 0 -9 1 -8 3 -9 -9 -3; |
---|
| 471 | -8 7 -4 -5 -9 1 -7 -1 3 -2]; |
---|
| 472 | x = sdpvar(10,1); |
---|
| 473 | t = sdpvar(1,1); |
---|
| 474 | p = c'*x-0.5*x'*Q*x; |
---|
| 475 | F = set(0<x<1)+set(A*x<b); |
---|
| 476 | solvesdp(F,p,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 477 | F = F + set(p<t); |
---|
| 478 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 479 | % ******************************** |
---|
| 480 | % Lassere |
---|
| 481 | % -0.37? >100 |
---|
| 482 | % ******************************** |
---|
| 483 | clear all |
---|
| 484 | x = sdpvar(10,1); |
---|
| 485 | t = sdpvar(1,1); |
---|
| 486 | x(1) = 1-sum(x(2:end)); |
---|
| 487 | obj = 0; |
---|
| 488 | for i = 1:9 |
---|
| 489 | obj = obj+x(i)*x(i+1); |
---|
| 490 | end |
---|
| 491 | for i = 1:8 |
---|
| 492 | obj = obj+x(i)*x(i+2); |
---|
| 493 | end |
---|
| 494 | obj = obj + x(1)*x(7)+ x(1)*x(9)+ x(2)*x(10)+ x(4)*x(7); |
---|
| 495 | F = set(1>x>0); |
---|
| 496 | F = F + set(-obj<t); |
---|
| 497 | solvesdp(F,-obj,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 498 | F = F + set(-obj<t); |
---|
| 499 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 500 | |
---|
| 501 | % ******************************** |
---|
| 502 | % Is J co-positive ? |
---|
| 503 | % Optimal : 0 |
---|
| 504 | % ******************************** |
---|
| 505 | clear all |
---|
| 506 | t = sdpvar(1,1); |
---|
| 507 | x = sdpvar(5,1); |
---|
| 508 | J = [1 -1 1 1 -1; |
---|
| 509 | -1 1 -1 1 1; |
---|
| 510 | 1 -1 1 -1 1; |
---|
| 511 | 1 1 -1 1 -1; |
---|
| 512 | -1 1 1 -1 1]; |
---|
| 513 | F = set(x>0) + set(sum(x)<1); |
---|
| 514 | F = F+set(-10<t<10); |
---|
| 515 | solvesdp(F,x'*J*x,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 516 | F = F + set(x'*J*x<t); |
---|
| 517 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 518 | |
---|
| 519 | % ******************************** |
---|
| 520 | % BARON example problems |
---|
| 521 | % -17 in 9 |
---|
| 522 | % ******************************** |
---|
| 523 | clear all |
---|
| 524 | x1 = sdpvar(1,1); |
---|
| 525 | x2 = sdpvar(1,1); |
---|
| 526 | x3 = sdpvar(1,1); |
---|
| 527 | x4 = sdpvar(1,1); |
---|
| 528 | x5 = sdpvar(1,1); |
---|
| 529 | t = sdpvar(1,1); |
---|
| 530 | p = 42*x1-50*x1^2+44*x2-50*x2^2+45*x3-50*x3^2+47*x4-50*x4^2+47.5*x5-50*x5^2; |
---|
| 531 | F = set([20 12 11 7 4]*[x1;x2;x3;x4;x5] < 40) + set(0<[x1;x2;x3;x4;x5]<1); |
---|
| 532 | solvesdp(F,p,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','nag','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 533 | F = F+set(p<t); |
---|
| 534 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','nag','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 535 | |
---|
| 536 | % ******************************** |
---|
| 537 | % BARON example |
---|
| 538 | % -13 in 1 |
---|
| 539 | % ******************************** |
---|
| 540 | clear all |
---|
| 541 | x1 = sdpvar(1,1); |
---|
| 542 | x2 = sdpvar(1,1); |
---|
| 543 | x3 = sdpvar(1,1); |
---|
| 544 | x4 = sdpvar(1,1); |
---|
| 545 | t = sdpvar(1,1); |
---|
| 546 | p = x1 - x2 - x3 - x1*x3 + x1*x4 + x2*x3 - x2*x4; |
---|
| 547 | F = set(x1+4*x2 <= 8) + set(4*x1+x2 <= 12) + set(3*x1+4*x2 <= 12) + set(2*x3+x4 <= 8) + set(x3+2*x4 <= 8) + set(x3+x4 <= 5)+set([x1;x2;x3;x4]>0); |
---|
| 548 | F = F+set(p<t); |
---|
| 549 | solvesdp(F,p,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.09e-1,'bmibnb.lpsolver','glpk')) |
---|
| 550 | F = F+set(p<t); |
---|
| 551 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.09e-1,'bmibnb.lpsolver','glpk')) |
---|
| 552 | |
---|
| 553 | % ******************************** |
---|
| 554 | % Control design |
---|
| 555 | % 5.46 in 12 |
---|
| 556 | % ******************************** |
---|
| 557 | A = [1 2;-3 0];B = [1;1]; |
---|
| 558 | [K0,P0] = lqr(A,B,eye(2),1); |
---|
| 559 | P = sdpvar(2,2);setsdpvar(P,2*P0);K0(K0>1)=1;K0(K0<-1)=-1; |
---|
| 560 | K = sdpvar(1,2);setsdpvar(K,-K0); |
---|
| 561 | F = set(K<1)+set(K>-1)+set(P>0)+set((A+B*K)'*P+P*(A+B*K) < -eye(2)-K'*K); |
---|
| 562 | %F = F + set([-(A+B*K)'*P-P*(A+B*K)-eye(2) K';K eye(1)] > 0) |
---|
| 563 | F = F+lmi(diag(P)>0)+lmi(P(:)>-151) + lmi(P(:)<150) + lmi(P>P0)+lmi(K>-100) + lmi(K<100); |
---|
| 564 | solvesdp(F,trace(P),sdpsettings('usex0',1,'bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk','bmibnb.lpreduce',1)) |
---|
| 565 | |
---|
| 566 | % ******************************** |
---|
| 567 | % GAMS st_qpc_m1 |
---|
| 568 | % optimal -473.77 in 3 |
---|
| 569 | % ******************************** |
---|
| 570 | clear all |
---|
| 571 | x1 = sdpvar(1,1); |
---|
| 572 | x2 = sdpvar(1,1); |
---|
| 573 | x3 = sdpvar(1,1); |
---|
| 574 | x4 = sdpvar(1,1); |
---|
| 575 | x5 = sdpvar(1,1); |
---|
| 576 | t = sdpvar(1,1); |
---|
| 577 | |
---|
| 578 | F = set([x1;x2;x3;x4;x5]>0); |
---|
| 579 | F = F + set(x1 + x2 + 2*x3 + x4 + x5 > 10); |
---|
| 580 | F = F + set(2*x1 + 3*x2 + x5 > 8); |
---|
| 581 | F = F + set(x2 + 4*x3 - x4 + 2*x5 > 12); |
---|
| 582 | F = F + set(8*x1 - x2 - x3 + 6*x4 > 20); |
---|
| 583 | F = F + set(- 2*x1 - x2 - 3*x3 - x4 - x5 > -30); |
---|
| 584 | obj = -(- (10*x1 - 0.34*x1*x1 - 0.28*x1*x2 + 10*x2 - 0.22*x1*x3 + 10*x3 - 0.24*x1*x4 + 10*x4 - 0.51*x1*x5 + 10*x5 - 0.28*x2*x1 - 0.34*x2*x2 - 0.23*x2*x3 -0.24*x2*x4 - 0.45*x2*x5 - 0.22*x3*x1 - 0.23*x3*x2 - 0.35*x3*x3 - 0.22*x3*x4 - 0.34*x3*x5 - 0.24*x4*x1 - 0.24*x4*x2 - 0.22*x4*x3 - 0.2*x4*x4 - 0.38*x4*x5 - 0.51*x5*x1 - 0.45*x5*x2 - 0.34*x5*x3 - 0.38*x5*x4 - 0.99*x5*x5)); |
---|
| 585 | |
---|
| 586 | solvesdp(F,obj,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 587 | F = F + set(obj<t); |
---|
| 588 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 589 | |
---|
| 590 | % ******************************** |
---|
| 591 | % GAMS st_qpk1 |
---|
| 592 | % optimal -3 in 7 |
---|
| 593 | % ******************************** |
---|
| 594 | clear all |
---|
| 595 | x1 = sdpvar(1,1); |
---|
| 596 | x2 = sdpvar(1,1); |
---|
| 597 | t = sdpvar(1,1); |
---|
| 598 | |
---|
| 599 | F = set([x1;x2]>0); |
---|
| 600 | F = F + set(- x1 + x2 < 1); |
---|
| 601 | F = F + set(x1 - x2 < 1); |
---|
| 602 | F = F + set(- x1 + 2*x2 < 3); |
---|
| 603 | F = F + set(2*x1 - x2 < 3); |
---|
| 604 | obj = (2*x1 - 2*x1*x1 + 2*x1*x2 + 3*x2 - 2*x2*x2); |
---|
| 605 | F = F + set(obj<t); |
---|
| 606 | solvesdp(F,t,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 607 | |
---|
| 608 | |
---|
| 609 | % ******************************** |
---|
| 610 | % GAMS st_robot |
---|
| 611 | % optimal -3.52 in 4 |
---|
| 612 | % ******************************** |
---|
| 613 | yalmip('clear') |
---|
| 614 | x1 = sdpvar(1,1); |
---|
| 615 | x2 = sdpvar(1,1); |
---|
| 616 | x3 = sdpvar(1,1); |
---|
| 617 | x4 = sdpvar(1,1); |
---|
| 618 | x5 = sdpvar(1,1); |
---|
| 619 | x6 = sdpvar(1,1); |
---|
| 620 | x7 = sdpvar(1,1); |
---|
| 621 | x8 = sdpvar(1,1); |
---|
| 622 | x = [x1;x2;x3;x4;x5;x6;x7;x8]; |
---|
| 623 | F = set(-1<x<1); |
---|
| 624 | F = F + set( 0.004731*x1*x3 - 0.1238*x1 - 0.3578*x2*x3 - 0.001637*x2 - 0.9338*x4 + x7 == 0.3571); |
---|
| 625 | F = F + set(0.2238*x1*x3 + 0.2638*x1 + 0.7623*x2*x3 - 0.07745*x2 - 0.6734*x4 - x7 == 0.6022); |
---|
| 626 | F = F + set( x6*x8 + 0.3578*x1 + 0.004731*x2 == 0); |
---|
| 627 | F = F + set( - 0.7623*x1 + 0.2238*x2 == -0.3461); |
---|
| 628 | F = F + set(x1^2 + x2^2 == 1); |
---|
| 629 | F = F + set(x3^2 + x4^2 == 1); |
---|
| 630 | F = F + set(x5^2 + x6^2 == 1); |
---|
| 631 | F = F + set(x7^2 + x8^2 == 1); |
---|
| 632 | |
---|
| 633 | solvesdp(F,sum(x),sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','pensdp','verbose',2,'solver','bmibnb','bnb.maxiter',30,'penbmi.P0',0.9e-1,'bmibnb.lpsolver','glpk')) |
---|
| 634 | |
---|
| 635 | |
---|
| 636 | % ******************************** |
---|
| 637 | % GAMS st_robot |
---|
| 638 | % optimal -382 in 5 |
---|
| 639 | % ******************************** |
---|
| 640 | yalmip('clear') |
---|
| 641 | x1 = sdpvar(1,1); |
---|
| 642 | x2 = sdpvar(1,1); |
---|
| 643 | x3 = sdpvar(1,1); |
---|
| 644 | x4 = sdpvar(1,1); |
---|
| 645 | x5 = sdpvar(1,1); |
---|
| 646 | x6 = sdpvar(1,1); |
---|
| 647 | x7 = sdpvar(1,1); |
---|
| 648 | x8 = sdpvar(1,1); |
---|
| 649 | x9 = sdpvar(1,1); |
---|
| 650 | x10 = sdpvar(1,1); |
---|
| 651 | t = sdpvar(1,1); |
---|
| 652 | x = [x1;x2;x3;x4;x5;x6;x7;x8;x9;x10]; |
---|
| 653 | |
---|
| 654 | e1= 20*x1 + 20*x2 + 60*x3 + 60*x4 + 60*x5 + 60*x6 + 5*x7 + 45*x8 + 55*x9 + 65*x10 - 600.1; |
---|
| 655 | |
---|
| 656 | e2= 5*x1 + 7*x2 + 3*x3 + 8*x4 + 13*x5 + 13*x6 + 2*x7 + 14*x8 + 14*x9 + 14*x10 - 310.5; |
---|
| 657 | |
---|
| 658 | e3= 100*x1 + 130*x2 + 50*x3 + 70*x4 + 70*x5 + 70*x6 + 20*x7 + 80*x8 + 80*x9 + 80*x10 - 1800; |
---|
| 659 | |
---|
| 660 | e4= 200*x1 + 280*x2 + 100*x3 + 200*x4 + 250*x5 + 280*x6 + 100*x7 + 180*x8 + 200*x9 + 220*x10 - 3850; |
---|
| 661 | |
---|
| 662 | e5= 2*x1 + 2*x2 + 4*x3 + 4*x4 + 4*x5 + 4*x6 + 2*x7 + 6*x8 + 6*x9 + 6*x10 - 18.6; |
---|
| 663 | |
---|
| 664 | e6= 4*x1 + 8*x2 + 2*x3 + 6*x4 + 10*x5 + 10*x6 + 5*x7 + 10*x8 + 10*x9 + 10*x10 - 198.7; |
---|
| 665 | |
---|
| 666 | e7= 60*x1 + 110*x2 + 20*x3 + 40*x4 + 60*x5 + 70*x6 + 10*x7 + 40*x8 + 50*x9 + 50*x10 - 882; |
---|
| 667 | |
---|
| 668 | e8= 150*x1 + 210*x2 + 40*x3 + 70*x4 + 90*x5 + 105*x6 + 60*x7 + 100*x8 + 140*x9 + 180*x10 - 4200; |
---|
| 669 | |
---|
| 670 | e9= 80*x1 + 100*x2 + 6*x3 + 16*x4 + 20*x5 + 22*x6 + 20*x8 + 30*x9 + 30*x10 - 40.25; |
---|
| 671 | |
---|
| 672 | e10= 40*x1 + 40*x2 + 12*x3 + 20*x4 + 24*x5 + 28*x6 + 40*x9 + 50*x10 - 327; |
---|
| 673 | |
---|
| 674 | obj = (10*x1 - 6.8*x1*x1 - 4.6*x1*x2 + 10*x2 - 7.9*x1*x3 + 10*x3 - 5.1*x1*x4 + 10*x4 - 6.9*x1*x5 + 10*x5 - 6.8*x1*x6 + 10*x6 - 4.6*x1*x7 + 10*x7 - 7.9*x1*x8 + 10*x8 - 5.1*x1*x9 + 10*x9 - 6.9*x1*x10 + 10*x10 - 4.6*x2*x1 - 5.5*x2*x2 - 5.8*x2*x3 - 4.5*x2*x4 - 6*x2*x5 - 4.6*x2*x6 - 5.5*x2*x7 - 5.8*x2*x8 - 4.5*x2*x9 - 6*x2*x10 - 7.9*x3*x1 - 5.8*x3*x2 - 13.3*x3*x3 - 6.7*x3*x4 - 8.9*x3*x5 - 7.9*x3*x6 - 5.8*x3*x7 - 13.3*x3*x8 - 6.7*x3*x9 - 8.9*x3*x10 - 5.1*x4*x1 - 4.5*x4*x2 - 6.7*x4*x3 - 6.9*x4*x4 - 5.8*x4*x5 - 5.1*x4*x6 - 4.5*x4*x7 - 6.7*x4*x8 - 6.9*x4*x9 - 5.8*x4*x10 - 6.9*x5*x1 - 6*x5*x2 - 8.9*x5*x3 - 5.8*x5*x4 - 11.9*x5*x5 - 6.9*x5*x6 - 6*x5*x7 - 8.9* x5*x8 - 5.8*x5*x9 - 11.9*x5*x10 - 6.8*x6*x1 - 4.6*x6*x2 - 7.9*x6*x3 - 5.1 *x6*x4 - 6.9*x6*x5 - 6.8*x6*x6 - 4.6*x6*x7 - 7.9*x6*x8 - 5.1*x6*x9 - 6.9* x6*x10 - 4.6*x7*x1 - 5.5*x7*x2 - 5.8*x7*x3 - 4.5*x7*x4 - 6*x7*x5 - 4.6*x7 *x6 - 5.5*x7*x7 - 5.8*x7*x8 - 4.5*x7*x9 - 6*x7*x10 - 7.9*x8*x1 - 5.8*x8* x2 - 13.3*x8*x3 - 6.7*x8*x4 - 8.9*x8*x5 - 7.9*x8*x6 - 5.8*x8*x7 - 13.3*x8 *x8 - 6.7*x8*x9 - 8.9*x8*x10 - 5.1*x9*x1 - 4.5*x9*x2 - 6.7*x9*x3 - 6.9*x9 *x4 - 5.8*x9*x5 - 5.1*x9*x6 - 4.5*x9*x7 - 6.7*x9*x8 - 6.9*x9*x9 - 5.8*x9* x10 - 6.9*x10*x1 - 6*x10*x2 - 8.9*x10*x3 - 5.8*x10*x4 - 11.9*x10*x5 - 6.9 *x10*x6 - 6*x10*x7 - 8.9*x10*x8 - 5.8*x10*x9 - 11.9*x10*x10); |
---|
| 675 | |
---|
| 676 | F = set(0<x<10000); |
---|
| 677 | F = F + set([e1;e2;e3;e4;e5;e6;e7;e8;e9;e10]<0); |
---|
| 678 | %F = F + set(obj==t); |
---|
| 679 | solvesdp(F,obj,sdpsettings('bmibnb.uppersolver','penbmi','bmibnb.lowersolver','glpk','solver','bmibnb','bmibnb.maxiter',15,'penbmi.P0',0.01,'bmibnb.lpsolver','glpk')) |
---|