[37] | 1 | function output = callmosek(interfacedata) |
---|
| 2 | |
---|
| 3 | % Author Johan Löfberg |
---|
| 4 | % $Id: callmosek.m,v 1.15 2006/05/22 09:58:52 joloef Exp $ |
---|
| 5 | |
---|
| 6 | % Retrieve needed data |
---|
| 7 | options = interfacedata.options; |
---|
| 8 | F_struc = interfacedata.F_struc; |
---|
| 9 | c = interfacedata.c; |
---|
| 10 | Q = interfacedata.Q; |
---|
| 11 | K = interfacedata.K; |
---|
| 12 | x0 = interfacedata.x0; |
---|
| 13 | integer_variables = interfacedata.integer_variables; |
---|
| 14 | binary_variables = interfacedata.binary_variables; |
---|
| 15 | extended_variables = interfacedata.extended_variables; |
---|
| 16 | ub = interfacedata.ub; |
---|
| 17 | lb = interfacedata.lb; |
---|
| 18 | mt = interfacedata.monomtable; |
---|
| 19 | |
---|
| 20 | % ********************************* |
---|
| 21 | % What type of variables do we have |
---|
| 22 | % ********************************* |
---|
| 23 | linear_variables = find((sum(abs(mt),2)==1) & (any(mt==1,2))); |
---|
| 24 | nonlinear_variables = setdiff((1:size(mt,1))',linear_variables); |
---|
| 25 | sigmonial_variables = find(any(0>mt,2) | any(mt-fix(mt),2)); |
---|
| 26 | |
---|
| 27 | if ~isempty(sigmonial_variables) | isequal(interfacedata.solver.version,'GEOMETRIC') |
---|
| 28 | [x,D_struc,problem,res,solvertime,prob] = call_mosek_geometric(options,F_struc,c,Q,K,ub,lb,mt,linear_variables,extended_variables); |
---|
| 29 | else |
---|
| 30 | [x,D_struc,problem,res,solvertime,prob] = call_mosek_lpqp(options,F_struc,c,Q,K,ub,lb,mt,linear_variables,integer_variables); |
---|
| 31 | end |
---|
| 32 | |
---|
| 33 | infostr = yalmiperror(problem,'MOSEK'); |
---|
| 34 | |
---|
| 35 | % Save all data sent to solver? |
---|
| 36 | if options.savesolverinput |
---|
| 37 | solverinput.prob = prob; |
---|
| 38 | else |
---|
| 39 | solverinput = []; |
---|
| 40 | end |
---|
| 41 | |
---|
| 42 | % Save all data from the solver? |
---|
| 43 | if options.savesolveroutput |
---|
| 44 | solveroutput.res = res; |
---|
| 45 | else |
---|
| 46 | solveroutput = []; |
---|
| 47 | end |
---|
| 48 | |
---|
| 49 | % Standard interface |
---|
| 50 | output.Primal = x; |
---|
| 51 | output.Dual = D_struc; |
---|
| 52 | output.Slack = []; |
---|
| 53 | output.problem = problem; |
---|
| 54 | output.infostr = infostr; |
---|
| 55 | output.solverinput = solverinput; |
---|
| 56 | output.solveroutput= solveroutput; |
---|
| 57 | output.solvertime = solvertime; |
---|
| 58 | |
---|
| 59 | |
---|
| 60 | function [x,D_struc,problem,res,solvertime,prob] = call_mosek_lpqp(options,F_struc,c,Q,K,ub,lb,mt,linear_variables,integer_variables); |
---|
| 61 | |
---|
| 62 | prob.c = c; |
---|
| 63 | if ~isempty(F_struc) |
---|
| 64 | prob.a = -F_struc(:,2:end); |
---|
| 65 | prob.buc = full(F_struc(:,1)); |
---|
| 66 | prob.blc = repmat(-inf,length(prob.buc),1); |
---|
| 67 | else |
---|
| 68 | prob.a = sparse(ones(1,length(c))); % Dummy constraint |
---|
| 69 | prob.buc = inf; |
---|
| 70 | prob.blc = -inf; |
---|
| 71 | end |
---|
| 72 | if isempty(lb) |
---|
| 73 | prob.blx = repmat(-inf,1,length(c)); |
---|
| 74 | else |
---|
| 75 | prob.blx = lb; |
---|
| 76 | end |
---|
| 77 | if isempty(ub) |
---|
| 78 | prob.bux = repmat(inf,1,length(c)); |
---|
| 79 | else |
---|
| 80 | prob.bux = ub; |
---|
| 81 | end |
---|
| 82 | |
---|
| 83 | |
---|
| 84 | if K.f>0 |
---|
| 85 | prob.blc(1:K.f) = prob.buc(1:K.f); |
---|
| 86 | end |
---|
| 87 | |
---|
| 88 | [prob.qosubi,prob.qosubj,prob.qoval] = find(tril(sparse(2*Q))); |
---|
| 89 | |
---|
| 90 | if K.q(1)>0 |
---|
| 91 | nof_new = sum(K.q); |
---|
| 92 | prob.a = [prob.a [zeros(K.f,nof_new);zeros(K.l,nof_new);eye(nof_new)]]; |
---|
| 93 | prob.a(1+K.f+K.l:end,1:length(c)) = prob.a(1+K.f+K.l:end,1:length(c)); |
---|
| 94 | prob.blc(1+K.f+K.l:end) = prob.buc(1+K.f+K.l:end); |
---|
| 95 | prob.buc(1+K.f+K.l:end) = prob.buc(1+K.f+K.l:end); |
---|
| 96 | prob.c = [prob.c;zeros(nof_new,1)]; |
---|
| 97 | top = size(F_struc,2)-1; |
---|
| 98 | for i = 1:length(K.q) |
---|
| 99 | prob.cones{i}.type = 'MSK_CT_QUAD'; |
---|
| 100 | prob.cones{i}.sub = top+1:top+K.q(i); |
---|
| 101 | prob.blx(top+1:top+K.q(i)) = -inf; |
---|
| 102 | prob.bux(top+1:top+K.q(i)) = inf; |
---|
| 103 | top = top + K.q(i); |
---|
| 104 | end |
---|
| 105 | end |
---|
| 106 | |
---|
| 107 | if ~isempty(integer_variables) |
---|
| 108 | prob.ints.sub = integer_variables; |
---|
| 109 | end |
---|
| 110 | |
---|
| 111 | prob.param = options.mosek.param; |
---|
| 112 | |
---|
| 113 | % Debug? |
---|
| 114 | if options.savedebug |
---|
| 115 | save mosekdebug prob |
---|
| 116 | end |
---|
| 117 | |
---|
| 118 | % Call MOSEK |
---|
| 119 | solvertime = clock; |
---|
| 120 | showprogress('Calling MOSEK',options.showprogress); |
---|
| 121 | if options.verbose == 0 |
---|
| 122 | [r,res] = mosekopt('minimize echo(0)',prob); |
---|
| 123 | else |
---|
| 124 | [r,res] = mosekopt('minimize',prob); |
---|
| 125 | end |
---|
| 126 | solvertime = etime(clock,solvertime); |
---|
| 127 | |
---|
| 128 | if (r == 1010) | (r == 1011) |
---|
| 129 | problem = -5; |
---|
| 130 | x = []; |
---|
| 131 | D_struc = []; |
---|
| 132 | elseif r == 1200 |
---|
| 133 | problem = 7; |
---|
| 134 | x = []; |
---|
| 135 | D_struc = []; |
---|
| 136 | else |
---|
| 137 | % Recover solutions |
---|
| 138 | sol = res.sol; |
---|
| 139 | if isempty(integer_variables) |
---|
| 140 | x = sol.itr.xx(1:length(c)); % Might have added new ones |
---|
| 141 | D_struc = (sol.itr.suc-sol.itr.slc); |
---|
| 142 | error_message = sol.itr.prosta; |
---|
| 143 | else |
---|
| 144 | x = sol.int.xx(1:length(c)); % Might have added new ones |
---|
| 145 | D_struc = []; |
---|
| 146 | error_message = sol.int.prosta; |
---|
| 147 | end |
---|
| 148 | |
---|
| 149 | switch error_message |
---|
| 150 | case {'PRIMAL_AND_DUAL_FEASIBLE','PRIMAL_FEASIBLE'} |
---|
| 151 | problem = 0; |
---|
| 152 | case 'PRIMAL_INFEASIBLE' |
---|
| 153 | problem = 1; |
---|
| 154 | case 'DUAL_INFEASIBLE' |
---|
| 155 | problem = 2; |
---|
| 156 | case 'PRIMAL_INFEASIBLE_OR_UNBOUNDED' |
---|
| 157 | problem = 12; |
---|
| 158 | otherwise |
---|
| 159 | problem = -1; |
---|
| 160 | end |
---|
| 161 | end |
---|
| 162 | |
---|
| 163 | |
---|
| 164 | |
---|
| 165 | function [x,D_struc,problem,res,solvertime,prob] = call_mosek_qclp(options,F_struc,c,Q,K,ub,lb,mt,linear_variables); |
---|
| 166 | |
---|
| 167 | nonlinear_variables =setdiff(1:length(c),linear_variables); |
---|
| 168 | var_prod = []; |
---|
| 169 | for i = 1:length(nonlinear_variables) |
---|
| 170 | vars = find(mt(nonlinear_variables(i),:)); |
---|
| 171 | if length(vars)==1 |
---|
| 172 | var_prod = [var_prod;vars vars]; |
---|
| 173 | else |
---|
| 174 | var_prod = [var_prod;sort(vars)]; |
---|
| 175 | end |
---|
| 176 | end |
---|
| 177 | |
---|
| 178 | A = -F_struc(:,1+linear_variables); |
---|
| 179 | u_c = full(F_struc(:,1)); |
---|
| 180 | l_c = repmat(-inf,length(u_c),1); |
---|
| 181 | l_x = lb; |
---|
| 182 | u_x = ub; |
---|
| 183 | |
---|
| 184 | if K.f>0 |
---|
| 185 | l_c(1:K.f) = u_c(1:K.f); |
---|
| 186 | end |
---|
| 187 | |
---|
| 188 | prob.c = c(linear_variables);; |
---|
| 189 | qcsubk = []; |
---|
| 190 | qcsubi = []; |
---|
| 191 | qcsubj = []; |
---|
| 192 | qcval = []; |
---|
| 193 | for k = 1:size(F_struc,1) |
---|
| 194 | nonlins = find(F_struc(k,nonlinear_variables+1)); |
---|
| 195 | for i = 1:length(nonlins) |
---|
| 196 | qcsubk = [qcsubk;k]; |
---|
| 197 | qcsubi = [qcsubi;var_prod(nonlins(i),1)]; |
---|
| 198 | qcsubj = [qcsubj;var_prod(nonlins(i),2)]; |
---|
| 199 | if var_prod(nonlins(i),2)==var_prod(nonlins(i),1) |
---|
| 200 | qcval = [qcval;-2*F_struc(k,1+nonlinear_variables(nonlins(i)))]; |
---|
| 201 | else |
---|
| 202 | qcval = [qcval;-F_struc(k,1+nonlinear_variables(nonlins(i)))]; |
---|
| 203 | end |
---|
| 204 | end |
---|
| 205 | end |
---|
| 206 | prob.qcsubk = qcsubk; |
---|
| 207 | prob.qcsubi = qcsubi; |
---|
| 208 | prob.qcsubj = qcsubj; |
---|
| 209 | prob.qcval = qcval; |
---|
| 210 | prob.a = A; |
---|
| 211 | prob.buc = u_c; |
---|
| 212 | |
---|
| 213 | |
---|
| 214 | % For debugging |
---|
| 215 | if options.savedebug |
---|
| 216 | save mosekdebug prob |
---|
| 217 | end |
---|
| 218 | |
---|
| 219 | % Call MOSEK |
---|
| 220 | solvertime = clock; |
---|
| 221 | showprogress('Calling MOSEK',options.showprogress); |
---|
| 222 | solvertime = clock; |
---|
| 223 | if options.verbose == 0 |
---|
| 224 | [r,res] = mosekopt('minimize echo(0)',prob); |
---|
| 225 | else |
---|
| 226 | [r,res] = mosekopt('minimize',prob); |
---|
| 227 | end |
---|
| 228 | solvertime = etime(clock,solvertime); |
---|
| 229 | |
---|
| 230 | x = zeros(length(c),1); |
---|
| 231 | x(linear_variables)=res.sol.itr.xx; |
---|
| 232 | D_struc = res.sol.itr.suc; |
---|
| 233 | if K.f>0 |
---|
| 234 | D_struc(1:K.f) = (sol.itr.suc(1:K.f)-sol.itr.slc(1:K.f)); |
---|
| 235 | end |
---|
| 236 | % Check, currently not exhaustive... |
---|
| 237 | switch res.sol.itr.prosta |
---|
| 238 | case 'PRIMAL_AND_DUAL_FEASIBLE' |
---|
| 239 | problem = 0; |
---|
| 240 | case 'PRIMAL_INFEASIBLE' |
---|
| 241 | problem = 1; |
---|
| 242 | case 'DUAL_INFEASIBLE' |
---|
| 243 | problem = 2; |
---|
| 244 | otherwise |
---|
| 245 | problem = -1; |
---|
| 246 | end |
---|
| 247 | |
---|
| 248 | function [x,D_struc,problem,res,solvertime,prob] = call_mosek_geometric(options,F_struc,c,Q,K,ub,lb,mt,linear_variables,extended_variables); |
---|
| 249 | |
---|
| 250 | x = []; |
---|
| 251 | D_struc = []; |
---|
| 252 | res = []; |
---|
| 253 | solvertime = 0; |
---|
| 254 | [prob,problem] = yalmip2geometric(options,F_struc,c,Q,K,ub,lb,mt,linear_variables,extended_variables); |
---|
| 255 | if problem == 0 |
---|
| 256 | |
---|
| 257 | % Mosek does not support equalities |
---|
| 258 | if ~isempty(prob.G) |
---|
| 259 | prob.A = [prob.A;prob.G;-prob.G]; |
---|
| 260 | prob.b = [prob.b;prob.h;1./prob.h]; |
---|
| 261 | prob.map = [prob.map;max(prob.map) + (1:2*length(prob.h))']; |
---|
| 262 | end |
---|
| 263 | if options.savedebug |
---|
| 264 | save mosekdebug prob |
---|
| 265 | end |
---|
| 266 | |
---|
| 267 | prob.param = options.mosek.param; |
---|
| 268 | |
---|
| 269 | % Call MOSEK |
---|
| 270 | solvertime = clock; |
---|
| 271 | showprogress('Calling MOSEK',options.showprogress); |
---|
| 272 | solvertime = clock; |
---|
| 273 | if options.verbose == 0 |
---|
| 274 | res = mskgpopt(prob.b,prob.A,prob.map,[],'minimize echo(0)'); |
---|
| 275 | else |
---|
| 276 | res = mskgpopt(prob.b,prob.A,prob.map,[],'minimize'); |
---|
| 277 | end |
---|
| 278 | sol = res.sol; |
---|
| 279 | solvertime = etime(clock,solvertime); |
---|
| 280 | |
---|
| 281 | x = zeros(length(c),1); |
---|
| 282 | x(linear_variables)=exp(res.sol.itr.xx); |
---|
| 283 | D_struc = [];%exp(res.sol.itr.suc); |
---|
| 284 | |
---|
| 285 | % Check, currently not exhaustive... |
---|
| 286 | switch res.sol.itr.prosta |
---|
| 287 | case 'PRIMAL_AND_DUAL_FEASIBLE' |
---|
| 288 | problem = 0; |
---|
| 289 | case 'PRIMAL_INFEASIBLE' |
---|
| 290 | problem = 1; |
---|
| 291 | case 'DUAL_INFEASIBLE' |
---|
| 292 | problem = 2; |
---|
| 293 | case 'UNKNOWN' |
---|
| 294 | problem = 9; |
---|
| 295 | otherwise |
---|
| 296 | problem = -1; |
---|
| 297 | end |
---|
| 298 | end |
---|