SCIP-SDP  3.2.0
sdpisolver_mosek.c
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of SCIPSDP - a solving framework for mixed-integer */
4 /* semidefinite programs based on SCIP. */
5 /* */
6 /* Copyright (C) 2011-2013 Discrete Optimization, TU Darmstadt */
7 /* EDOM, FAU Erlangen-Nürnberg */
8 /* 2014-2020 Discrete Optimization, TU Darmstadt */
9 /* */
10 /* */
11 /* This program is free software; you can redistribute it and/or */
12 /* modify it under the terms of the GNU Lesser General Public License */
13 /* as published by the Free Software Foundation; either version 3 */
14 /* of the License, or (at your option) any later version. */
15 /* */
16 /* This program is distributed in the hope that it will be useful, */
17 /* but WITHOUT ANY WARRANTY; without even the implied warranty of */
18 /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
19 /* GNU Lesser General Public License for more details. */
20 /* */
21 /* You should have received a copy of the GNU Lesser General Public License */
22 /* along with this program; if not, write to the Free Software */
23 /* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.*/
24 /* */
25 /* */
26 /* Based on SCIP - Solving Constraint Integer Programs */
27 /* Copyright (C) 2002-2020 Zuse Institute Berlin */
28 /* SCIP is distributed under the terms of the SCIP Academic Licence, */
29 /* see file COPYING in the SCIP distribution. */
30 /* */
31 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
32 
33 /*#define SCIP_DEBUG*/
34 /*#define SCIP_MORE_DEBUG*/
35 /*#define SCIP_DEBUG_PRINTTOFILE *//* prints each problem inserted into MOSEK to the file mosek.task */
36 
64 #include <assert.h>
65 #include <sys/time.h>
66 
67 #include "sdpi/sdpisolver.h"
68 
69 #include "blockmemshell/memory.h" /* for memory allocation */
70 #include "scip/def.h" /* for SCIP_Real, _Bool, ... */
71 #include "scip/pub_misc.h" /* for sorting */
72 #include "mosek.h" /* for MOSEK routines */
73 #include "sdpi/sdpsolchecker.h" /* to check solution with regards to feasibility tolerance */
74 #include "scip/pub_message.h" /* for debug and error message */
75 
76 /* @todo Use MSK_putexitfunc to catch errors.
77  * @todo Think about what to do with near optimality etc. (If MOSEK cannot compute a solution that has the prescribed accuracy, then it will
78  * multiply the termination tolerances with MSK_DPAR_INTPNT_CO_TOL_NEAR_REL. If the solution then satisfies the termination criteria, then
79  * the solution is denoted near optimal, near feasible and so forth.)
80  */
81 
82 #define MIN_PENALTYPARAM 1e5
83 #define MAX_PENALTYPARAM 1e10
84 #define PENALTYPARAM_FACTOR 1e6
85 #define MAX_MAXPENALTYPARAM 1e15
86 #define MAXPENALTYPARAM_FACTOR 1e6
87 #define PENALTYBOUNDTOL 1E-3
89 #define INFEASFEASTOLCHANGE 0.1
90 #define INFEASMINFEASTOL 1E-9
92 #define CONVERT_ABSOLUTE_TOLERANCES TRUE
93 #if MSK_VERSION_MAJOR >= 9
94 #define NEAR_REL_TOLERANCE 1.0
95 #endif
96 
98 struct SCIP_SDPiSolver
99 {
100  SCIP_MESSAGEHDLR* messagehdlr;
101  BMS_BLKMEM* blkmem;
102  BMS_BUFMEM* bufmem;
103  MSKenv_t mskenv;
104  MSKtask_t msktask;
105  int nvars;
106  int nactivevars;
107  int* inputtomosekmapper;
110  int* mosektoinputmapper;
111  SCIP_Real* fixedvarsval;
112  SCIP_Real fixedvarsobjcontr;
113  SCIP_Real* objcoefs;
114  int nvarbounds;
115  int* varboundpos;
119  SCIP_Bool solved;
120  int sdpcounter;
121  SCIP_Real epsilon;
122  SCIP_Real gaptol;
123  SCIP_Real feastol;
124  SCIP_Real sdpsolverfeastol;
125  SCIP_Real objlimit;
126  SCIP_Bool sdpinfo;
127  SCIP_Bool penalty;
128  SCIP_Bool feasorig;
130  SCIP_Bool rbound;
131  MSKrescodee terminationcode;
132  SCIP_Bool timelimit;
133  SCIP_Bool timelimitinitial;
134  int nthreads;
135  int niterations;
136  int nsdpcalls;
137 };
138 
139 
140 /*
141  * Local Methods
142  */
143 
145 #define MOSEK_CALL(x) do \
146  { \
147  MSKrescodee _mosekerrorcode_; \
148  if ( (_mosekerrorcode_ = (x)) != MSK_RES_OK ) \
149  { \
150  SCIPerrorMessage("MOSEK-Error <%d> in function call.\n", (int)_mosekerrorcode_); \
151  return SCIP_LPERROR; \
152  } \
153  } \
154  while( FALSE )
155 
157 #define MOSEK_CALL_BOOL(x) do \
158  { \
159  MSKrescodee _mosekerrorcode_; \
160  if ( (_mosekerrorcode_ = (x)) != MSK_RES_OK ) \
161  { \
162  SCIPerrorMessage("MOSEK-Error <%d> in function call.\n", (int)_mosekerrorcode_); \
163  return FALSE; \
164  } \
165  } \
166  while( FALSE )
167 
169 #define MOSEK_CALLM(x) do \
170  { \
171  MSKrescodee _mosekerrorcode_; \
172  if ( (_mosekerrorcode_ = (x)) != MSK_RES_OK ) \
173  { \
174  SCIPerrorMessage("MOSEK-Error <%d> in function call.\n", (int)_mosekerrorcode_); \
175  return SCIP_NOMEMORY; \
176  } \
177  } \
178  while( FALSE )
179 
181 #define BMS_CALL(x) do \
182  { \
183  if( NULL == (x) ) \
184  { \
185  SCIPerrorMessage("No memory in function call.\n"); \
186  return SCIP_NOMEMORY; \
187  } \
188  } \
189  while( FALSE )
190 
192 #define TIMEOFDAY_CALL(x) do \
193  { \
194  if ( (x) != 0 ) \
195  { \
196  SCIPerrorMessage("Error in gettimeofday! \n"); \
197  return SCIP_ERROR; \
198  } \
199  } \
200  while( FALSE )
201 
203 #define CHECK_IF_SOLVED(sdpisolver) do \
204  { \
205  if (!(sdpisolver->solved)) \
206  { \
207  SCIPerrorMessage("Tried to access solution information for SDP %d ahead of solving!\n", sdpisolver->sdpcounter); \
208  return SCIP_LPERROR; \
209  } \
210  } \
211  while( FALSE )
212 
214 #define CHECK_IF_SOLVED_BOOL(sdpisolver) do \
215  { \
216  if (!(sdpisolver->solved)) \
217  { \
218  SCIPerrorMessage("Tried to access solution information for SDP %d ahead of solving!\n", sdpisolver->sdpcounter); \
219  assert( 0 ); \
220  return FALSE; \
221  } \
222  } \
223  while( FALSE )
224 
226 static
227 void MSKAPI printstr(
228  void* handle,
229 #if MSK_VERSION_MAJOR >= 9
230  const char str[]
231 #else
232  MSKCONST char str[]
233 #endif
234  )
235 { /*lint --e{715}*/
236  printf("%s",str);
237 }
238 
239 #ifndef NDEBUG
240 
241 static
242 SCIP_Bool isFixed(
243  SCIP_SDPISOLVER* sdpisolver,
244  SCIP_Real lb,
245  SCIP_Real ub
246  )
247 {
248  assert( sdpisolver != NULL );
249  assert( lb < ub + sdpisolver->feastol );
250 
251  return (ub-lb <= sdpisolver->epsilon);
252 }
253 #else
254 #define isFixed(sdpisolver,lb,ub) (ub-lb <= sdpisolver->epsilon)
255 #endif
256 
257 
258 /*
259  * Miscellaneous Methods
260  */
261 
265 char name[SCIP_MAXSTRLEN];
266 
268 const char* SCIPsdpiSolverGetSolverName(
269  void
270  )
271 {
272  MSKrescodee rescodee;
273  int majorver = 0;
274  int minorver = 0;
275 #if MSK_VERSION_MAJOR < 9
276  int build = 0;
277 #endif
278  int revision = 0;
279 #ifndef NDEBUG
280  int snprintfreturn; /* used to check the return code of snprintf */
281 #endif
282 
283 #if MSK_VERSION_MAJOR < 9
284  rescodee = MSK_getversion(&majorver, &minorver, &build, &revision);/*lint !e123*/
285 #else
286  rescodee = MSK_getversion(&majorver, &minorver, &revision);/*lint !e123*/
287 #endif
288 
289  if ( rescodee != MSK_RES_OK )
290  return "MOSEK";
291 
292 #ifndef NDEBUG
293  #if MSK_VERSION_MAJOR < 9
294  snprintfreturn = SCIPsnprintf(name, SCIP_MAXSTRLEN, "MOSEK %d.%d.%d.%d", majorver, minorver, build, revision);/*lint !e123*/
295  #else
296  snprintfreturn = SCIPsnprintf(name, SCIP_MAXSTRLEN, "MOSEK %d.%d.%d", majorver, minorver, revision);/*lint !e123*/
297  #endif
298  assert( snprintfreturn < SCIP_MAXSTRLEN ); /* check whether the name fits into the string */
299 #else
300  #if MSK_VERSION_MAJOR < 9
301  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "MOSEK %d.%d.%d.%d", majorver, minorver, build, revision);
302  #else
303  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "MOSEK %d.%d.%d", majorver, minorver, revision);
304  #endif
305 #endif
306 
307  return name;
308 }
309 
311 const char* SCIPsdpiSolverGetSolverDesc(
312  void
313  )
314 {
315  return "Homogeneous and self-dual interior-point solver for semidefinite programming developed by MOSEK ApS"
316  "(http://www.mosek.com)";
317 }
318 
326  SCIP_SDPISOLVER* sdpisolver
327  )
328 {
329  assert( sdpisolver != NULL );
330  return (void*) NULL;
331 }
332 
335  void
336  )
337 {
338  return 1E-7;
339 }
340 
343  void
344  )
345 {
346  return 1E-5;
347 }
348 
351  void
352  )
353 {
354  return 8;
355 }
356 
359  void
360  )
361 {
362  return FALSE;
363 }
364 
368 /*
369  * SDPI Creation and Destruction Methods
370  */
371 
376 SCIP_RETCODE SCIPsdpiSolverCreate(
377  SCIP_SDPISOLVER** sdpisolver,
378  SCIP_MESSAGEHDLR* messagehdlr,
379  BMS_BLKMEM* blkmem,
380  BMS_BUFMEM* bufmem
381  )
382 {
383  assert( sdpisolver != NULL );
384  assert( blkmem != NULL );
385  assert( bufmem != NULL );
386 
387  SCIPdebugMessage("Calling SCIPsdpiCreate \n");
388 
389  BMS_CALL( BMSallocBlockMemory(blkmem, sdpisolver) );
390 
391  (*sdpisolver)->messagehdlr = messagehdlr;
392  (*sdpisolver)->blkmem = blkmem;
393  (*sdpisolver)->bufmem = bufmem;
394 
395  MOSEK_CALLM( MSK_makeenv(&((*sdpisolver)->mskenv), NULL) );/*lint !e641*/ /* the NULL-argument is a debug file, but setting this will spam the whole folder */
396 
397  /* this will be properly initialized then calling solve */
398  (*sdpisolver)->msktask = NULL;
399 
400  (*sdpisolver)->nvars = 0;
401  (*sdpisolver)->nactivevars = 0;
402  (*sdpisolver)->inputtomosekmapper = NULL;
403  (*sdpisolver)->mosektoinputmapper = NULL;
404  (*sdpisolver)->fixedvarsval = NULL;
405  (*sdpisolver)->fixedvarsobjcontr = 0.0;
406  (*sdpisolver)->objcoefs = NULL;
407  (*sdpisolver)->nvarbounds = 0;
408  (*sdpisolver)->varboundpos = NULL;
409  (*sdpisolver)->solved = FALSE;
410  (*sdpisolver)->sdpcounter = 0;
411 
412  (*sdpisolver)->epsilon = 1e-9;
413  (*sdpisolver)->gaptol = 1e-4;
414  (*sdpisolver)->feastol = 1e-6;
415  (*sdpisolver)->sdpsolverfeastol = 1e-6;
416  (*sdpisolver)->objlimit = SCIPsdpiSolverInfinity(*sdpisolver);
417  (*sdpisolver)->sdpinfo = FALSE;
418  (*sdpisolver)->nthreads = -1;
419  (*sdpisolver)->timelimit = FALSE;
420  (*sdpisolver)->timelimitinitial = FALSE;
421 
422  return SCIP_OKAY;
423 }
424 
426 SCIP_RETCODE SCIPsdpiSolverFree(
427  SCIP_SDPISOLVER** sdpisolver
428  )
429 {
430  assert( sdpisolver != NULL );
431  assert( *sdpisolver != NULL );
432 
433  SCIPdebugMessage("Freeing SDPISolver\n");
434 
435  if ( (*sdpisolver)->msktask != NULL )
436  {
437  MOSEK_CALL( MSK_deletetask(&((*sdpisolver)->msktask)) );/*lint !e641*/
438  }
439 
440  if ( (*sdpisolver)->mskenv != NULL )
441  {
442  MOSEK_CALL( MSK_deleteenv(&((*sdpisolver)->mskenv)) );/*lint !e641*/
443  }
444 
445  BMSfreeBlockMemoryArrayNull((*sdpisolver)->blkmem, &(*sdpisolver)->varboundpos, 2 * (*sdpisolver)->nactivevars);
446  BMSfreeBlockMemoryArrayNull((*sdpisolver)->blkmem, &(*sdpisolver)->inputtomosekmapper, (*sdpisolver)->nvars);
447  BMSfreeBlockMemoryArrayNull((*sdpisolver)->blkmem, &(*sdpisolver)->mosektoinputmapper, (*sdpisolver)->nactivevars);
448  BMSfreeBlockMemoryArrayNull((*sdpisolver)->blkmem, &(*sdpisolver)->fixedvarsval, (*sdpisolver)->nvars - (*sdpisolver)->nactivevars);
449  BMSfreeBlockMemoryArrayNull((*sdpisolver)->blkmem, &(*sdpisolver)->objcoefs, (*sdpisolver)->nactivevars);
450 
451  BMSfreeBlockMemory((*sdpisolver)->blkmem, sdpisolver);
452 
453  return SCIP_OKAY;
454 }
455 
457 SCIP_RETCODE SCIPsdpiSolverIncreaseCounter(
458  SCIP_SDPISOLVER* sdpisolver
459  )
460 {
461  assert( sdpisolver != NULL );
462 
463  sdpisolver->sdpcounter++;
464 
465  return SCIP_OKAY;
466 }
467 
469 SCIP_RETCODE SCIPsdpiSolverResetCounter(
470  SCIP_SDPISOLVER* sdpisolver
471  )
472 {
473  assert( sdpisolver != NULL );
474 
475  SCIPdebugMessage("Resetting counter of SDP-Interface from %d to 0.\n", sdpisolver->sdpcounter);
476  sdpisolver->sdpcounter = 0;
477 
478  return SCIP_OKAY;
479 }
480 
484 /*
485  * Solving Methods
486  */
487 
506 SCIP_RETCODE SCIPsdpiSolverLoadAndSolve(
507  SCIP_SDPISOLVER* sdpisolver,
508  int nvars,
509  SCIP_Real* obj,
510  SCIP_Real* lb,
511  SCIP_Real* ub,
512  int nsdpblocks,
513  int* sdpblocksizes,
514  int* sdpnblockvars,
515  int sdpconstnnonz,
516  int* sdpconstnblocknonz,
518  int** sdpconstrow,
519  int** sdpconstcol,
520  SCIP_Real** sdpconstval,
521  int sdpnnonz,
522  int** sdpnblockvarnonz,
524  int** sdpvar,
526  int*** sdprow,
527  int*** sdpcol,
528  SCIP_Real*** sdpval,
529  int** indchanges,
531  int* nremovedinds,
532  int* blockindchanges,
533  int nremovedblocks,
534  int nlpcons,
535  int noldlpcons,
536  SCIP_Real* lplhs,
537  SCIP_Real* lprhs,
538  int* rownactivevars,
539  int lpnnonz,
540  int* lprow,
541  int* lpcol,
542  SCIP_Real* lpval,
543  SCIP_Real* starty,
544  int* startZnblocknonz,
546  int** startZrow,
548  int** startZcol,
550  SCIP_Real** startZval,
552  int* startXnblocknonz,
554  int** startXrow,
556  int** startXcol,
558  SCIP_Real** startXval,
560  SCIP_SDPSOLVERSETTING startsettings,
562  SCIP_Real timelimit
563  )
564 {
565  return SCIPsdpiSolverLoadAndSolveWithPenalty(sdpisolver, 0.0, TRUE, FALSE, nvars, obj, lb, ub, nsdpblocks, sdpblocksizes, sdpnblockvars, sdpconstnnonz,
566  sdpconstnblocknonz, sdpconstrow, sdpconstcol, sdpconstval, sdpnnonz, sdpnblockvarnonz, sdpvar, sdprow, sdpcol, sdpval, indchanges,
567  nremovedinds, blockindchanges, nremovedblocks, nlpcons, noldlpcons, lplhs, lprhs, rownactivevars, lpnnonz, lprow, lpcol, lpval,
568  starty, startZnblocknonz, startZrow, startZcol, startZval, startXnblocknonz, startXrow, startXcol, startXval, startsettings,
569  timelimit, NULL, NULL);
570 }
571  /*lint --e{715}*/
595  SCIP_SDPISOLVER* sdpisolver,
596  SCIP_Real penaltyparam,
597  SCIP_Bool withobj,
598  SCIP_Bool rbound,
599  int nvars,
600  SCIP_Real* obj,
601  SCIP_Real* lb,
602  SCIP_Real* ub,
603  int nsdpblocks,
604  int* sdpblocksizes,
605  int* sdpnblockvars,
606  int sdpconstnnonz,
607  int* sdpconstnblocknonz,
609  int** sdpconstrow,
610  int** sdpconstcol,
611  SCIP_Real** sdpconstval,
612  int sdpnnonz,
613  int** sdpnblockvarnonz,
615  int** sdpvar,
617  int*** sdprow,
618  int*** sdpcol,
619  SCIP_Real*** sdpval,
620  int** indchanges,
622  int* nremovedinds,
623  int* blockindchanges,
624  int nremovedblocks,
625  int nlpcons,
626  int noldlpcons,
627  SCIP_Real* lplhs,
628  SCIP_Real* lprhs,
629  int* rownactivevars,
630  int lpnnonz,
631  int* lprow,
632  int* lpcol,
633  SCIP_Real* lpval,
634  SCIP_Real* starty,
635  int* startZnblocknonz,
637  int** startZrow,
639  int** startZcol,
641  SCIP_Real** startZval,
643  int* startXnblocknonz,
645  int** startXrow,
647  int** startXcol,
649  SCIP_Real** startXval,
651  SCIP_SDPSOLVERSETTING startsettings,
653  SCIP_Real timelimit,
654  SCIP_Bool* feasorig,
656  SCIP_Bool* penaltybound
658 )
659 {/*lint --e{715}*/
660  int b;
661  int i;
662  int j;
663  int blockvar;
664  int v;
665  int k;
666  long long ind;
667  int mosekind = 0;
668  SCIP_Real* mosekvarbounds;
669  int nfixedvars;
670  int oldnactivevars;
671  int* vartorowmapper; /* maps the lpvars to the corresponding left- and right-hand-sides of the LP constraints */
672  int* vartolhsrhsmapper; /* maps the lpvars to the corresponding entries in lplhs and lprhs */
673  int nlpvars;
674  int pos;
675  int newpos;
676  int* mosekblocksizes;
677  SCIP_Real one; /* MOSEK always wants a pointer to factors for a sum of matrices, we always use a single matrix with factor one */
678  int* mosekrow;
679  int* mosekcol;
680  SCIP_Real* mosekval;
681  int row;
682  SCIP_Real val;
683  struct timeval starttime;
684  struct timeval currenttime;
685  SCIP_Real startseconds;
686  SCIP_Real currentseconds;
687  SCIP_Real elapsedtime;
688  MSKsolstae solstat;
689 #ifdef SCIP_MORE_DEBUG
690  int nmosekconss;
691  int nmosekvars;
692  int nmosekcones;
693  char varname[SCIP_MAXSTRLEN];
694 #endif
695 #if CONVERT_ABSOLUTE_TOLERANCES
696  SCIP_Real maxrhscoef; /* MOSEK uses a relative feasibility tolerance, the largest rhs-coefficient is needed for converting the absolute tolerance */
697 #endif
698 
699  assert( sdpisolver != NULL );
700  assert( sdpisolver->mskenv != NULL );
701  assert( penaltyparam > -1 * sdpisolver->epsilon );
702  assert( penaltyparam < sdpisolver->epsilon || ( feasorig != NULL ) );
703  assert( nvars > 0 );
704  assert( obj != NULL );
705  assert( lb != NULL );
706  assert( ub != NULL );
707  assert( nsdpblocks >= 0 );
708  assert( nsdpblocks == 0 || sdpblocksizes != NULL );
709  assert( nsdpblocks == 0 || sdpnblockvars != NULL );
710  assert( sdpconstnnonz >= 0 );
711  assert( nsdpblocks == 0 || sdpconstnnonz == 0 || sdpconstnblocknonz != NULL );
712  assert( nsdpblocks == 0 || sdpconstnnonz == 0 || sdpconstrow != NULL );
713  assert( nsdpblocks == 0 || sdpconstnnonz == 0 || sdpconstcol != NULL );
714  assert( nsdpblocks == 0 || sdpconstnnonz == 0 || sdpconstval != NULL );
715  assert( sdpnnonz >= 0 );
716  assert( nsdpblocks == 0 || sdpnblockvarnonz != NULL );
717  assert( nsdpblocks == 0 || sdpvar != NULL );
718  assert( nsdpblocks == 0 || sdprow != NULL );
719  assert( nsdpblocks == 0 || sdpcol != NULL );
720  assert( nsdpblocks == 0 || sdpval != NULL );
721  assert( nsdpblocks == 0 || indchanges != NULL );
722  assert( nsdpblocks == 0 || nremovedinds != NULL );
723  assert( nsdpblocks == 0 || blockindchanges != NULL );
724  assert( 0 <= nremovedblocks && nremovedblocks <= nsdpblocks );
725  assert( nlpcons >= 0 );
726  assert( noldlpcons >= nlpcons );
727  assert( nlpcons == 0 || lplhs != NULL );
728  assert( nlpcons == 0 || lprhs != NULL );
729  assert( nlpcons == 0 || rownactivevars != NULL );
730  assert( lpnnonz >= 0 );
731  assert( nlpcons == 0 || lprow != NULL );
732  assert( nlpcons == 0 || lpcol != NULL );
733  assert( nlpcons == 0 || lpval != NULL );
734 
735  /* check the timelimit */
736  if ( timelimit <= 0.0 )
737  {
738  sdpisolver->timelimit = TRUE;
739  sdpisolver->timelimitinitial = TRUE;
740  sdpisolver->solved = FALSE;
741  return SCIP_OKAY;
742  }
743  sdpisolver->timelimit = FALSE;
744  sdpisolver->timelimitinitial = FALSE;
745  sdpisolver->feasorig = FALSE;
746 
747  /* start the timing */
748  TIMEOFDAY_CALL( gettimeofday(&starttime, NULL) );/*lint !e438, !e550, !e641 */
749 
750  one = 1.0;
751 #if CONVERT_ABSOLUTE_TOLERANCES
752  maxrhscoef = 0.0;
753 #endif
754 
755  /* create an empty task (second and third argument are guesses for maximum number of constraints and variables), if there already is one, delete it */
756  if ((sdpisolver->msktask) != NULL)
757  {
758  MOSEK_CALL( MSK_deletetask(&(sdpisolver->msktask)) );/*lint !e641*/
759  }
760  if ( penaltyparam < sdpisolver->epsilon )
761  {
762  MOSEK_CALLM( MSK_maketask(sdpisolver->mskenv, nvars, nsdpblocks - nremovedblocks + nlpcons + 2 * nvars, &(sdpisolver->msktask)) );/*lint !e641*/
763  }
764  else
765  {
766  MOSEK_CALLM( MSK_maketask(sdpisolver->mskenv, nvars + 1, nsdpblocks - nremovedblocks + nlpcons + 2 * nvars, &(sdpisolver->msktask)) );/*lint !e641*/
767  }
768 
769 #if MSK_VERSION_MAJOR >= 9
770  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_NEAR_REL, NEAR_REL_TOLERANCE) );
771 #endif
772 
773 #ifdef SCIP_MORE_DEBUG
774  MOSEK_CALL( MSK_linkfunctotaskstream (sdpisolver->msktask, MSK_STREAM_LOG, NULL, printstr) ); /* output to console */
775 #else
776  /* if sdpinfo is true, redirect output to console */
777  if ( sdpisolver->sdpinfo )
778  {
779  MOSEK_CALL( MSK_linkfunctotaskstream (sdpisolver->msktask, MSK_STREAM_LOG, NULL, printstr) );/*lint !e641*/
780  }
781 #endif
782 
783  /* set number of threads */
784  if ( sdpisolver->nthreads > 0 )
785  {
786  MOSEK_CALL( MSK_putintparam(sdpisolver->msktask, MSK_IPAR_NUM_THREADS, sdpisolver->nthreads) );/*lint !e641*/
787  }
788 
789  /* only increase the counter if we don't use the penalty formulation to stay in line with the numbers in the general interface (where this is still the
790  * same SDP) */
791  if ( penaltyparam < sdpisolver->epsilon )
792  SCIPdebugMessage("Inserting Data into MOSEK for SDP (%d) \n", ++sdpisolver->sdpcounter);
793  else
794  SCIPdebugMessage("Inserting Data again into MOSEK for SDP (%d) \n", sdpisolver->sdpcounter);
795 
796  /* set the penalty and rbound flags accordingly */
797  sdpisolver->penalty = (penaltyparam < sdpisolver->epsilon) ? FALSE : TRUE;
798  sdpisolver->rbound = rbound;
799 
800  /* allocate memory for inputtomosekmapper, mosektoinputmapper and the fixed and active variable information, for the latter this will
801  * later be shrinked if the needed size is known */
802  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->inputtomosekmapper), sdpisolver->nvars, nvars) );
803  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->mosektoinputmapper), sdpisolver->nactivevars, nvars) );
804  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->fixedvarsval), sdpisolver->nvars - sdpisolver->nactivevars, nvars) ); /*lint !e776*/
805  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->objcoefs), sdpisolver->nactivevars, nvars) ); /*lint !e776*/
806 
807  oldnactivevars = sdpisolver->nactivevars; /* we need to save this to realloc the varboundpos-array if needed */
808  sdpisolver->nvars = nvars;
809  sdpisolver->nactivevars = 0;
810  nfixedvars = 0;
811 
812  /* find the fixed variables */
813  sdpisolver->fixedvarsobjcontr = 0.0;
814  for (i = 0; i < nvars; i++)
815  {
816  if ( isFixed(sdpisolver, lb[i], ub[i]) )
817  {
818  sdpisolver->fixedvarsobjcontr += obj[i] * lb[i]; /* this is the value this fixed variable contributes to the objective */
819  sdpisolver->fixedvarsval[nfixedvars] = lb[i]; /* if lb=ub, then this is the value the variable will have in every solution */
820  nfixedvars++;
821  sdpisolver->inputtomosekmapper[i] = -nfixedvars;
822  SCIPdebugMessage("Fixing variable %d locally to %f for SDP %d in MOSEK\n", i, lb[i], sdpisolver->sdpcounter);
823  }
824  else
825  {
826  sdpisolver->mosektoinputmapper[sdpisolver->nactivevars] = i;
827  sdpisolver->inputtomosekmapper[i] = sdpisolver->nactivevars;
828  sdpisolver->objcoefs[sdpisolver->nactivevars] = obj[i];
829  sdpisolver->nactivevars++;
830 #ifdef SCIP_MORE_DEBUG
831  SCIPdebugMessage("Variable %d becomes variable %d for SDP %d in MOSEK\n", i, sdpisolver->inputtomosekmapper[i], sdpisolver->sdpcounter);
832 #endif
833  }
834  }
835  assert( sdpisolver->nactivevars + nfixedvars == sdpisolver->nvars );
836 
837  /* if we want to solve without objective, we reset fixedvarsobjcontr */
838  if ( ! withobj )
839  sdpisolver->fixedvarsobjcontr = 0.0;
840 
841  /* shrink the fixedvars, objcoefs and mosektoinputmapper arrays to the right size */
842  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->objcoefs), nvars, sdpisolver->nactivevars) );
843  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->fixedvarsval), nvars, nfixedvars) );
844  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->mosektoinputmapper), nvars, sdpisolver->nactivevars) );
845 
846  /* compute number of variable bounds and save them in mosekvarbounds */
847  sdpisolver->nvarbounds = 0;
848  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekvarbounds, 2 * sdpisolver->nactivevars) ); /*lint !e647*/
849 
850  if ( sdpisolver->nactivevars != oldnactivevars )
851  {
852  if ( sdpisolver->varboundpos == NULL )
853  {
854  BMS_CALL( BMSallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->varboundpos), 2 * sdpisolver->nactivevars) ); /*lint !e647*/
855  }
856  else
857  {
858  BMS_CALL( BMSreallocBlockMemoryArray(sdpisolver->blkmem, &(sdpisolver->varboundpos), 2 * oldnactivevars, 2 * sdpisolver->nactivevars) ); /*lint !e647*/
859  }
860  }
861  assert( sdpisolver->varboundpos != NULL );
862 
863  for (i = 0; i < sdpisolver->nactivevars; i++)
864  {
865  assert( 0 <= sdpisolver->mosektoinputmapper[i] && sdpisolver->mosektoinputmapper[i] < nvars );
866  if ( ! SCIPsdpiSolverIsInfinity(sdpisolver, lb[sdpisolver->mosektoinputmapper[i]]) )
867  {
868  mosekvarbounds[sdpisolver->nvarbounds] = lb[sdpisolver->mosektoinputmapper[i]];
869  sdpisolver->varboundpos[sdpisolver->nvarbounds] = -(i + 1); /* negative sign means lower bound */
870  (sdpisolver->nvarbounds)++;
871  }
872  if ( ! SCIPsdpiSolverIsInfinity(sdpisolver, ub[sdpisolver->mosektoinputmapper[i]]) )
873  {
874  mosekvarbounds[sdpisolver->nvarbounds] = -1 * ub[sdpisolver->mosektoinputmapper[i]]; /* we give the upper bounds a negative sign for the objective */
875  sdpisolver->varboundpos[sdpisolver->nvarbounds] = +(i + 1); /* positive sign means upper bound */
876  (sdpisolver->nvarbounds)++;
877  }
878  }
879 
880  if ( nlpcons > 0 )
881  {
882  /* allocate memory to save which lpvariable corresponds to which lp constraint, negative signs correspond to left-hand-sides of lp constraints,
883  * entry i or -i corresponds to the constraint in position |i|-1, as we have to add +1 to make the entries strictly positive or strictly negative */
884  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &vartorowmapper, 2*noldlpcons) ); /*lint !e647*/ /*lint !e530*/
885  /* allocate memory to save at which indices the corresponding lhss and rhss of the lpvars are saved */
886  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &vartolhsrhsmapper, 2*noldlpcons) ); /*lint !e647*/ /*lint !e530*/
887 
888  /* compute the number of LP constraints after splitting the ranged rows and compute the rowtovarmapper */
889  pos = 0;
890  newpos = 0; /* the position in the lhs and rhs arrays */
891  for (i = 0; i < noldlpcons; i++)
892  {
893  assert( newpos <= nlpcons );
894  if ( rownactivevars[i] >= 2 )
895  {
896  if ( lplhs[newpos] > - SCIPsdpiSolverInfinity(sdpisolver) )
897  {
898  vartorowmapper[pos] = -(i+1);
899  vartolhsrhsmapper[pos] = newpos;
900  pos++;
901 
902 #if CONVERT_ABSOLUTE_TOLERANCES
903  /* update largest rhs-entry */
904  if ( REALABS(lplhs[newpos]) > maxrhscoef )
905  maxrhscoef = REALABS(lplhs[newpos]);
906 #endif
907 
908  }
909  if ( lprhs[newpos] < SCIPsdpiSolverInfinity(sdpisolver) )
910  {
911  vartorowmapper[pos] = i+1;
912  vartolhsrhsmapper[pos] = newpos;
913  pos++;
914 
915 #if CONVERT_ABSOLUTE_TOLERANCES
916  /* update largest rhs-entry */
917  if ( REALABS(lprhs[newpos]) > maxrhscoef )
918  maxrhscoef = REALABS(lprhs[newpos]);
919 #endif
920  }
921  newpos++;
922  }
923  }
924  nlpvars = pos;
925  assert( nlpvars <= 2*nlpcons ); /* *2 comes from left- and right-hand-sides */
926  }
927  else
928  nlpvars = 0;
929 
930  /* create matrix variables */
931  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekblocksizes, nsdpblocks - nremovedblocks) ); /*lint !e679 !e776*/
932 
933  for (b = 0; b < nsdpblocks; b++)
934  {
935  if ( blockindchanges[b] > -1 )
936  {
937  assert( 0 <= blockindchanges[b] && blockindchanges[b] <= b && (b - blockindchanges[b]) <= (nsdpblocks - nremovedblocks) );
938  mosekblocksizes[b - blockindchanges[b]] = sdpblocksizes[b] - nremovedinds[b];/*lint !e679*/
939  }
940  }
941  MOSEK_CALLM( MSK_appendbarvars(sdpisolver->msktask, nsdpblocks - nremovedblocks, mosekblocksizes) );/*lint !e641*/
942 
943  /* create scalar variables (since we solve the primal problem, these are not the active variables but the dual variables for the
944  * lp constraints and variable bounds) */
945  MOSEK_CALLM( MSK_appendvars(sdpisolver->msktask, nlpvars + sdpisolver->nvarbounds) );/*lint !e641*/
946 
947  /* all of these are non-negative */
948  for (v = 0; v < nlpvars + sdpisolver->nvarbounds; v++)
949  {
950  MOSEK_CALL( MSK_putvarbound(sdpisolver->msktask, v, MSK_BK_LO, 0.0, (double) MSK_DPAR_DATA_TOL_BOUND_INF) );/*lint !e641*/
951  }
952 
953  /* append empty constraints (since we solve the primal problem, we get one constraint for each active variable) */
954  MOSEK_CALLM( MSK_appendcons(sdpisolver->msktask, (penaltyparam < sdpisolver->epsilon) ? sdpisolver->nactivevars : sdpisolver->nactivevars + 1) );/*lint !e641*/
955 
956  /* set objective values for the matrix variables */
957  i = 0;
958 
959  if ( sdpconstnnonz > 0 )
960  {
961  for (b = 0; b < nsdpblocks; b++)
962  {
963  if ( blockindchanges[b] > -1 )
964  {
965  /* if some indices in the block were removed, we need to change indices accordingly */
966  if ( nremovedinds[b] > 0 )
967  {
968  int* moseksdpconstrow;
969  int* moseksdpconstcol;
970 
971  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &moseksdpconstrow, sdpconstnblocknonz[b]) );
972  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &moseksdpconstcol, sdpconstnblocknonz[b]) );
973 
974  for (k = 0; k < sdpconstnblocknonz[b]; k++)
975  {
976  /* rows and cols with active nonzeros should not be removed */
977  assert( -1 < indchanges[b][sdpconstrow[b][k]] && indchanges[b][sdpconstrow[b][k]] <= sdpconstrow[b][k] );
978  assert( -1 < indchanges[b][sdpconstcol[b][k]] && indchanges[b][sdpconstcol[b][k]] <= sdpconstcol[b][k] );
979 
980  assert( 0 <= sdpconstrow[b][k] && sdpconstrow[b][k] <= sdpblocksizes[b] );
981  assert( 0 <= sdpconstcol[b][k] && sdpconstcol[b][k] <= sdpblocksizes[b] );
982 
983  moseksdpconstrow[k] = sdpconstrow[b][k] - indchanges[b][sdpconstrow[b][k]];
984  moseksdpconstcol[k] = sdpconstcol[b][k] - indchanges[b][sdpconstcol[b][k]];
985 
986 #if CONVERT_ABSOLUTE_TOLERANCES
987  /* update largest rhs-entry */
988  if ( REALABS(sdpconstval[b][k]) > maxrhscoef )
989  maxrhscoef = REALABS(sdpconstval[b][k]);
990 #endif
991  }
992 
993  MOSEK_CALL( MSK_appendsparsesymmat(sdpisolver->msktask, mosekblocksizes[b - blockindchanges[b]], sdpconstnblocknonz[b],
994  moseksdpconstrow, moseksdpconstcol, sdpconstval[b], &ind) );/*lint !e641, !e679, !e747*/
995 
996  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &moseksdpconstcol);
997  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &moseksdpconstrow);
998  }
999  else
1000  {
1001  MOSEK_CALL( MSK_appendsparsesymmat(sdpisolver->msktask, mosekblocksizes[b - blockindchanges[b]], sdpconstnblocknonz[b],
1002  sdpconstrow[b], sdpconstcol[b], sdpconstval[b], &ind) );/*lint !e641, !e679, !e747*/
1003  }
1004  MOSEK_CALL( MSK_putbarcj(sdpisolver->msktask, i, 1, &ind, &one) );/*lint !e641, !e747*/
1005  i++;
1006  }
1007  }
1008  }
1009 
1010  /* set objective values for the scalar variables */
1011  /* TODO: check for equality constraints */
1012  /* first for those corresponding to LP constraints in the dual */
1013  for (i = 0; i < nlpvars; i++)
1014  {
1015  assert( vartolhsrhsmapper != NULL ); /* for lint */
1016  if ( vartorowmapper[i] > 0 )/*lint !e644*/ /* right-hand side */
1017  {
1018  MOSEK_CALL( MSK_putcj(sdpisolver->msktask, i, -1 * lprhs[vartolhsrhsmapper[i]]) );/*lint !e641, !e644*/
1019 #ifdef SCIP_MORE_DEBUG
1020  /* give the variable a meaningful name for debug output */
1021  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "rhs-%d", vartorowmapper[i] - 1);
1022  MOSEK_CALL( MSK_putvarname ( sdpisolver->msktask, i, varname) );
1023 #endif
1024  }
1025  else /* left-hand side */
1026  {
1027  assert( vartorowmapper[i] < 0 ); /* we should not have value 0 so that we can clearly differentiate between positive and negative */
1028  MOSEK_CALL( MSK_putcj(sdpisolver->msktask, i, lplhs[vartolhsrhsmapper[i]]) );/*lint !e641*/
1029 #ifdef SCIP_MORE_DEBUG
1030  /* give the variable a meaningful name for debug output */
1031  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "lhs-%d", -1 * vartorowmapper[i] - 1);
1032  MOSEK_CALL( MSK_putvarname ( sdpisolver->msktask, i, varname) );
1033 #endif
1034  }
1035  }
1036 
1037  /* finally for those corresponding to variable bounds in the dual */
1038  for (i = 0; i < sdpisolver->nvarbounds; i++)
1039  {
1040  MOSEK_CALL( MSK_putcj(sdpisolver->msktask, nlpvars + i, mosekvarbounds[i]) );/*lint !e641*/ /* for the ub's we already added a negative sign in mosekvarbounds*/
1041 #ifdef SCIP_MORE_DEBUG
1042  if ( sdpisolver->varboundpos[i] < 0 ) /* lower bound */
1043  {
1044  /* give the variable a meaningful name for debug output */
1045  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "lb-%d", sdpisolver->mosektoinputmapper[-1 * sdpisolver->varboundpos[i] - 1]);
1046  MOSEK_CALL( MSK_putvarname ( sdpisolver->msktask, nlpvars + i, varname) );
1047  }
1048  else /* upper bound */
1049  {
1050  assert( sdpisolver->varboundpos[i] > 0 ); /* we should not have value 0 so that we can clearly differentiate between positive and negative */
1051  /* give the variable a meaningful name for debug output */
1052  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "ub-%d", sdpisolver->mosektoinputmapper[sdpisolver->varboundpos[i] - 1]);
1053  MOSEK_CALL( MSK_putvarname ( sdpisolver->msktask, nlpvars + i, varname) );
1054  }
1055 #endif
1056  }
1057 
1058  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &mosekvarbounds);
1059 
1060  /* set objective sense (since we want to minimize in the dual, we maximize in the primal) */
1061  MOSEK_CALL( MSK_putobjsense(sdpisolver->msktask, MSK_OBJECTIVE_SENSE_MAXIMIZE) );/*lint !e641*/
1062 
1063  /* start inserting the constraints */
1064 
1065  /* first add the matrices A_i */
1066  for (b = 0; b < nsdpblocks; b++)
1067  {
1068  if ( blockindchanges[b] > -1 )
1069  {
1070  for (blockvar = 0; blockvar < sdpnblockvars[b]; blockvar++)
1071  {
1072  v = sdpisolver->inputtomosekmapper[sdpvar[b][blockvar]];
1073 
1074  /* check if the variable is active */
1075  if ( v > -1 )
1076  {
1077  assert( v < sdpisolver->nactivevars );
1078  /* if there are removed indices, we have to adjust the column and row indices accordingly */
1079  if ( nremovedinds[b] > 0 )
1080  {
1081  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekrow, sdpnblockvarnonz[b][blockvar]) );
1082  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekcol, sdpnblockvarnonz[b][blockvar]) );
1083 
1084  for (k = 0; k < sdpnblockvarnonz[b][blockvar]; k++)
1085  {
1086  /* rows and cols with active nonzeros should not be removed */
1087  assert( -1 < indchanges[b][sdprow[b][blockvar][k]] && indchanges[b][sdprow[b][blockvar][k]] <= sdprow[b][blockvar][k] );
1088  assert( -1 < indchanges[b][sdpcol[b][blockvar][k]] && indchanges[b][sdpcol[b][blockvar][k]] <= sdpcol[b][blockvar][k] );
1089 
1090  assert( 0 <= sdprow[b][blockvar][k] && sdprow[b][blockvar][k] < sdpblocksizes[b] );
1091  assert( 0 <= sdpcol[b][blockvar][k] && sdpcol[b][blockvar][k] < sdpblocksizes[b] );
1092 
1093  mosekrow[k] = sdprow[b][blockvar][k] - indchanges[b][sdprow[b][blockvar][k]];
1094  mosekcol[k] = sdpcol[b][blockvar][k] - indchanges[b][sdpcol[b][blockvar][k]];
1095  }
1096 
1097  MOSEK_CALL( MSK_appendsparsesymmat(sdpisolver->msktask, mosekblocksizes[b - blockindchanges[b]], (long long) sdpnblockvarnonz[b][blockvar],
1098  mosekrow, mosekcol, sdpval[b][blockvar], &ind) );/*lint !e641, !e679*/
1099 
1100  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &mosekcol);
1101  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &mosekrow);
1102  }
1103  else
1104  {
1105  MOSEK_CALL( MSK_appendsparsesymmat(sdpisolver->msktask, mosekblocksizes[b - blockindchanges[b]], (long long) sdpnblockvarnonz[b][blockvar],
1106  sdprow[b][blockvar], sdpcol[b][blockvar], sdpval[b][blockvar], &ind) );/*lint !e641, !e679*/
1107  }
1108 
1109  MOSEK_CALL( MSK_putbaraij(sdpisolver->msktask, v, b - blockindchanges[b], (long long) 1, &ind, &one) );/*lint !e641*/
1110  }
1111  }
1112  }
1113  }
1114 
1115  /* add the identity matrix corresponding to the penalty variable */
1116  if ( penaltyparam >= sdpisolver->epsilon )
1117  {
1118  int* identityindices;
1119  SCIP_Real* identityvalues;
1120 
1121  for (b = 0; b < nsdpblocks; b++)
1122  {
1123  if ( blockindchanges[b] > -1 )
1124  {
1125  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &identityindices, mosekblocksizes[b - blockindchanges[b]]) );
1126  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &identityvalues, mosekblocksizes[b - blockindchanges[b]]) );
1127 
1128  for (i = 0; i < mosekblocksizes[b - blockindchanges[b]]; i++)
1129  {
1130  identityindices[i] = i;
1131  identityvalues[i] = 1.0;
1132  }
1133  MOSEK_CALL( MSK_appendsparsesymmat(sdpisolver->msktask, mosekblocksizes[b - blockindchanges[b]], (long long) mosekblocksizes[b - blockindchanges[b]],
1134  identityindices, identityindices, identityvalues, &ind) );/*lint !e641, !e679*/
1135  MOSEK_CALL( MSK_putbaraij(sdpisolver->msktask, sdpisolver->nactivevars, b - blockindchanges[b], (long long) 1, &ind, &one) );/*lint !e641, !e679*/
1136 
1137  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &identityvalues);
1138  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &identityindices);
1139  }
1140  }
1141  }
1142 
1143  /* now add the entries corresponding to the lp-constraints in the dual problem */
1144  if ( penaltyparam < sdpisolver->epsilon )
1145  {
1146  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekrow, lpnnonz) );
1147  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekval, lpnnonz) );
1148  }
1149  else
1150  {
1151  /* one extra entry for the penalty-constraint Trace = Gamma */
1152  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekrow, lpnnonz + 1) );/*lint !e776*/
1153  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &mosekval, lpnnonz + 1) );/*lint !e776*/
1154  }
1155 
1156  ind = 0;
1157  for (i = 0; i < nlpvars; i++)
1158  {
1159  assert( vartorowmapper != NULL ); /* for lint */
1160  if ( vartorowmapper[i] < 0 ) /* left-hand side */
1161  {
1162  mosekind = 0;
1163  /* find the first lp-entry belonging to this variable (those in between have to belong to constraints with less than two active variables and
1164  * will therefore not be used) */
1165  while (ind < lpnnonz && lprow[ind] < -1 * vartorowmapper[i] - 1)
1166  ind++;
1167  /* iterate over all nonzeros to input them to the array given to MOSEK */
1168  while (ind < lpnnonz && lprow[ind] == -1 * vartorowmapper[i] - 1) /* they should already be sorted by rows in the sdpi */
1169  {
1170  v = sdpisolver->inputtomosekmapper[lpcol[ind]];
1171  if ( v > -1 )
1172  {
1173  mosekrow[mosekind] = v;
1174  mosekval[mosekind] = lpval[ind];
1175  mosekind++;
1176  }
1177  ind++;
1178  }
1179  assert( mosekind <= lpnnonz );
1180  }
1181  else /* right-hand side */
1182  {
1183  assert( vartorowmapper[i] > 0 ); /* we should not have value 0 so that we can clearly differentiate between positive and negative */
1184 
1185  if ( i > 0 && vartorowmapper[i] == -1 * vartorowmapper[i - 1] )
1186  {
1187  /* we already iterated over this constraint in the lp-arrays, so we keep the current ind position and as we currenlty have
1188  * the corresponding entries in the mosek-arrays, we iterate over them again just changing the signs (except for the penalty-entry) */
1189  for (j = 0; j < (penaltyparam < sdpisolver->epsilon ? mosekind : mosekind - 1); j++)/*lint !e644*/
1190  mosekval[j] *= -1;
1191  }
1192  else
1193  {
1194  /* no left hand side for this constraint exists, so we didnot yet iterate over this constraint in the lp arrays */
1195  mosekind = 0;
1196  /* find the first lp-entry belonging to this variable (those in between have to belong to constraints with less than two active variables and
1197  * will therefore not be used) */
1198  while (lprow[ind] < vartorowmapper[i] - 1)
1199  ind++;
1200  while (ind < lpnnonz && lprow[ind] == vartorowmapper[i] - 1)
1201  {
1202  v = sdpisolver->inputtomosekmapper[lpcol[ind]];
1203  if ( v > -1 )
1204  {
1205  mosekrow[mosekind] = v;
1206  mosekval[mosekind] = -1 * lpval[ind]; /* because we need to change the <= to a >= constraint */
1207  mosekind++;
1208  }
1209  ind++;
1210  }
1211  assert( mosekind <= lpnnonz );
1212  }
1213  }
1214 
1215  /* add the additional entry for the penalty constraint Trace = Gamma */
1216  if ( penaltyparam >= sdpisolver->epsilon )
1217  {
1218  /* check if we reset mosekind, in case we did not and just "copied" the lhs-entries for the rhs, we do not need to reset the extra entry,
1219  * since it is already there */
1220  if ( ! (i > 0 && vartorowmapper[i] == -1 * vartorowmapper[i - 1] ))
1221  {
1222  mosekrow[mosekind] = sdpisolver->nactivevars;
1223  mosekval[mosekind] = 1.0;
1224  mosekind++;
1225  }
1226  assert( mosekind <= lpnnonz + 1 );
1227  }
1228 
1229  MOSEK_CALL( MSK_putacol(sdpisolver->msktask, i, mosekind, mosekrow, mosekval) );/*lint !e641*/
1230  }
1231 
1232  BMSfreeBufferMemoryArrayNull(sdpisolver->bufmem, &mosekval);
1233  BMSfreeBufferMemoryArrayNull(sdpisolver->bufmem, &mosekrow);
1234 
1235  /* finally add the entries corresponding to varbounds in the dual problem, we get exactly one entry per variable */
1236  for (i = 0; i < sdpisolver->nvarbounds; i++)
1237  {
1238  if ( sdpisolver->varboundpos[i] < 0 )
1239  {
1240  /* lower bound */
1241  row =-1 * sdpisolver->varboundpos[i] - 1; /* minus one because we added one to get strictly positive/negative values */
1242  val = 1.0;
1243  }
1244  else
1245  {
1246  /* upper bound */
1247  assert( sdpisolver->varboundpos[i] > 0 ); /* we should not have a zero as we wanted a clear differentiation between positive and negative */
1248  row = sdpisolver->varboundpos[i] - 1; /* minus one because we added one to get strictly positive/negative values */
1249  val = -1.0;
1250  }
1251  MOSEK_CALL( MSK_putacol(sdpisolver->msktask, nlpvars + i, 1, &row, &val) );/*lint !e641*/
1252  }
1253 
1254  /* make all constraints equality constraints with right-hand side b_i (or 0 if we solve without objective) */
1255  for (i = 0; i < sdpisolver->nactivevars; i++)
1256  {
1257  if ( withobj )
1258  {
1259  MOSEK_CALL( MSK_putconbound(sdpisolver->msktask, i, MSK_BK_FX, obj[sdpisolver->mosektoinputmapper[i]], obj[sdpisolver->mosektoinputmapper[i]]) );/*lint !e641*/
1260  }
1261  else
1262  {
1263  MOSEK_CALL( MSK_putconbound(sdpisolver->msktask, i, MSK_BK_FX, 0.0, 0.0) );/*lint !e641*/
1264  }
1265 #ifdef SCIP_MORE_DEBUG
1266  /* give the constraint a meaningful name for debug output */
1267  (void) SCIPsnprintf(varname, SCIP_MAXSTRLEN, "var-%d", sdpisolver->mosektoinputmapper[i]);
1268  MOSEK_CALL( MSK_putconname ( sdpisolver->msktask, i, varname) );
1269 #endif
1270  }
1271 
1272  /* the penalty constraint has right-hand side Gamma, it is a <=-inequality if r was bounded and an equality constraint otherwise */
1273  if ( penaltyparam >= sdpisolver->epsilon )
1274  {
1275  if ( rbound )
1276  {
1277  MOSEK_CALL( MSK_putconbound(sdpisolver->msktask, sdpisolver->nactivevars, MSK_BK_UP, (double) -1 * MSK_DPAR_DATA_TOL_BOUND_INF, penaltyparam) );/*lint !e641*/
1278  }
1279  else
1280  {
1281  MOSEK_CALL( MSK_putconbound(sdpisolver->msktask, sdpisolver->nactivevars, MSK_BK_FX, penaltyparam, penaltyparam) );/*lint !e641*/
1282  }
1283  }
1284 
1285  /* write to file if asked to */
1286 #ifdef SCIP_DEBUG_PRINTTOFILE
1287  SCIP_CALL( SCIPsdpiSolverWriteSDP(sdpisolver, "mosek.task") );
1288 #endif
1289 
1290  /* print whole problem and parameters if asked to */
1291 #ifdef SCIP_MORE_DEBUG
1292 #if MSK_VERSION_MAJOR < 9
1293  MOSEK_CALL( MSK_getnumcon (sdpisolver->msktask, &nmosekconss) );
1294  MOSEK_CALL( MSK_getnumvar (sdpisolver->msktask, &nmosekvars) );
1295  MOSEK_CALL( MSK_getnumcone (sdpisolver->msktask, &nmosekcones) );
1296 
1297  MOSEK_CALL( MSK_printdata (sdpisolver->msktask, MSK_STREAM_LOG, 0, nmosekconss, 0, nmosekvars, 0, nmosekcones, 1, 1, 1, 1, 1, 1, 1, 1) );
1298 #ifdef SCIP_PRINT_PARAMETERS
1299  MOSEK_CALL( MSK_printparam (sdpisolver->msktask) );
1300 #endif
1301 #endif
1302 #endif
1303 
1304  /* set the time limit */
1305  startseconds = (SCIP_Real) starttime.tv_sec + (SCIP_Real) starttime.tv_usec / 1e6;
1306 
1307  TIMEOFDAY_CALL( gettimeofday(&currenttime, NULL) );/*lint !e438, !e550, !e641 */
1308  currentseconds = (SCIP_Real) currenttime.tv_sec + (SCIP_Real) currenttime.tv_usec / 1e6;
1309 
1310  elapsedtime = currentseconds - startseconds;
1311 
1312  if ( timelimit <= elapsedtime )
1313  {
1314  sdpisolver->timelimit = TRUE;
1315  sdpisolver->solved = FALSE;
1316  }
1317  else
1318  {
1319  SCIP_Real feastol;
1320 
1321  /* set epsilon and feasibility tolerance (note that the dual in MOSEK is the problem we are interested in, so this is where we use feastol,
1322  * since MOSEK works with relative tolerance, we adjust our absolute tolerance accordingly, so that any solution satisfying the relative
1323  * tolerance in MOSEK satisfies our absolute tolerance) */
1324  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_PFEAS, sdpisolver->gaptol) );/*lint !e641*/
1325 #if CONVERT_ABSOLUTE_TOLERANCES
1326  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_DFEAS, sdpisolver->sdpsolverfeastol / (1 + maxrhscoef)) );/*lint !e641*/
1327  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_INFEAS, sdpisolver->sdpsolverfeastol / (1 + maxrhscoef)) );/*lint !e641*/
1328  SCIPdebugMessage("Setting relative feasibility tolerance for MOSEK to %.10f / %f = %.12f\n", sdpisolver->sdpsolverfeastol,
1329  1+maxrhscoef, sdpisolver->sdpsolverfeastol / (1 + maxrhscoef));
1330 #else
1331  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_DFEAS, sdpisolver->sdpsolverfeastol) );/*lint !e641*/
1332  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_INFEAS, sdpisolver->sdpsolverfeastol) );/*lint !e641*/
1333 #endif
1334  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_MU_RED, sdpisolver->gaptol) );/*lint !e641*/
1335  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_REL_GAP, sdpisolver->gaptol) );/*lint !e641*/
1336 
1337  if ( ! SCIPsdpiSolverIsInfinity(sdpisolver, timelimit - elapsedtime) )
1338  {
1339  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_OPTIMIZER_MAX_TIME, timelimit - elapsedtime) );/*lint !e641*/
1340  }
1341 
1342  /* set objective cutoff */
1343  if ( ! SCIPsdpiSolverIsInfinity(sdpisolver, sdpisolver->objlimit) )
1344  {
1345  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_UPPER_OBJ_CUT, sdpisolver->objlimit) );/*lint !e641*/
1346  }
1347 
1348  /* solve the problem */
1349  MOSEK_CALL( MSK_optimizetrm(sdpisolver->msktask, &(sdpisolver->terminationcode)) );/*lint !e641*/
1350 
1351  if ( sdpisolver->sdpinfo )
1352  {
1353  MOSEK_CALL( MSK_optimizersummary(sdpisolver->msktask, MSK_STREAM_LOG) );/*lint !e641*/
1354  MOSEK_CALL( MSK_analyzesolution(sdpisolver->msktask, MSK_STREAM_LOG, MSK_SOL_ITR) );/*lint !e641*/
1355  }
1356 
1357  SCIPdebugMessage("Solving problem using MOSEK, return code %d\n", sdpisolver->terminationcode);
1358 
1359  sdpisolver->solved = TRUE;
1360 
1361  sdpisolver->nsdpcalls = 1;
1362  MOSEK_CALL( MSK_getnaintinf(sdpisolver->msktask, "MSK_IINF_INTPNT_ITER", &(sdpisolver->niterations)) );/*lint !e641*/
1363 
1364  /* if the problem has been stably solved but did not reach the required feasibility tolerance, even though the solver
1365  * reports feasibility, resolve it with adjusted tolerance */
1366 #if CONVERT_ABSOLUTE_TOLERANCES
1367  feastol = sdpisolver->sdpsolverfeastol / (1 + maxrhscoef);
1368 #else
1369  feastol = sdpisolver->sdpsolverfeastol;
1370 #endif
1371 
1372  while ( SCIPsdpiSolverIsAcceptable(sdpisolver) && SCIPsdpiSolverIsDualFeasible(sdpisolver) && penaltyparam < sdpisolver->epsilon && feastol >= INFEASMINFEASTOL )
1373  {
1374  SCIP_Real* solvector;
1375  int nvarspointer;
1376  SCIP_Bool infeasible;
1377  int newiterations;
1378 
1379  /* get current solution */
1380  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &solvector, nvars) );
1381  nvarspointer = nvars;
1382  SCIP_CALL( SCIPsdpiSolverGetSol(sdpisolver, NULL, solvector, &nvarspointer) );
1383  assert( nvarspointer == nvars );
1384 
1385  /* check the solution for feasibility with regards to our tolerance */
1386  SCIP_CALL( SCIPsdpSolcheckerCheck(sdpisolver->bufmem, nvars, lb, ub, nsdpblocks, sdpblocksizes, sdpnblockvars, sdpconstnnonz,
1387  sdpconstnblocknonz, sdpconstrow, sdpconstcol, sdpconstval, sdpnnonz, sdpnblockvarnonz, sdpvar, sdprow, sdpcol, sdpval,
1388  indchanges, nremovedinds, blockindchanges, nlpcons, noldlpcons, lplhs, lprhs, rownactivevars, lpnnonz, lprow, lpcol, lpval,
1389  solvector, sdpisolver->feastol, sdpisolver->epsilon, &infeasible) );
1390 
1391  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &solvector);
1392 
1393  if ( infeasible )
1394  {
1395  SCIPdebugMessage("Solution feasible for MOSEK but outside feasibility tolerance, changing MOSEK feasibility tolerance from %f to %f\n",
1396  feastol, feastol * INFEASFEASTOLCHANGE);
1397  feastol *= INFEASFEASTOLCHANGE;
1398 
1399  if ( feastol >= INFEASMINFEASTOL )
1400  {
1401  /* update settings */
1402  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_DFEAS, feastol) );/*lint !e641*/
1403  MOSEK_CALL( MSK_putdouparam(sdpisolver->msktask, MSK_DPAR_INTPNT_CO_TOL_INFEAS, feastol) );/*lint !e641*/
1404 
1405  /* set the time limit */
1406  startseconds = (SCIP_Real) starttime.tv_sec + (SCIP_Real) starttime.tv_usec / 1e6;
1407 
1408  TIMEOFDAY_CALL( gettimeofday(&currenttime, NULL) );/*lint !e438, !e550, !e641 */
1409  currentseconds = (SCIP_Real) currenttime.tv_sec + (SCIP_Real) currenttime.tv_usec / 1e6;
1410 
1411  elapsedtime = currentseconds - startseconds;
1412 
1413  if ( timelimit <= elapsedtime )
1414  {
1415  sdpisolver->timelimit = TRUE;
1416  sdpisolver->solved = FALSE;
1417  }
1418 
1419  /* solve the problem */
1420  MOSEK_CALL( MSK_optimizetrm(sdpisolver->msktask, &(sdpisolver->terminationcode)) );/*lint !e641*/
1421 
1422  if ( sdpisolver->sdpinfo )
1423  {
1424  MOSEK_CALL( MSK_optimizersummary(sdpisolver->msktask, MSK_STREAM_LOG) );/*lint !e641*/
1425  MOSEK_CALL( MSK_analyzesolution(sdpisolver->msktask, MSK_STREAM_LOG, MSK_SOL_ITR) );/*lint !e641*/
1426  }
1427 
1428  /* update number of SDP-iterations and -calls */
1429  sdpisolver->nsdpcalls++;
1430  MOSEK_CALL( MSK_getnaintinf(sdpisolver->msktask, "MSK_IINF_INTPNT_ITER", &newiterations) );/*lint !e641*/
1431  sdpisolver->niterations += newiterations;
1432  }
1433  else
1434  {
1435  sdpisolver->solved = FALSE;
1436  SCIPmessagePrintInfo(sdpisolver->messagehdlr, "MOSEK failed to reach required feasibility tolerance! \n");
1437  }
1438  }
1439  else
1440  break;
1441  }
1442  }
1443 
1444 
1445  /* if using a penalty formulation, check if the solution is feasible for the original problem
1446  * we should always count it as infeasible if the penalty problem was unbounded */
1447  MOSEK_CALL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1448  if ( penaltyparam >= sdpisolver->epsilon && (solstat == MSK_SOL_STA_PRIM_INFEAS_CER) )
1449  {
1450  assert( feasorig != NULL );
1451  *feasorig = FALSE;
1452  SCIPdebugMessage("Penalty Problem unbounded!\n");
1453  }
1454  else if ( penaltyparam >= sdpisolver->epsilon && ( ! sdpisolver->timelimit ) && ( sdpisolver->terminationcode != MSK_RES_TRM_MAX_TIME ) )
1455  {
1456  SCIP_Real* moseksol;
1457  SCIP_Real trace = 0.0;
1458  SCIP_Real* x;
1459 
1460  assert( feasorig != NULL );
1461 
1462  /* get the r variable in the dual problem */
1463  BMSallocBufferMemoryArray(sdpisolver->bufmem, &moseksol, sdpisolver->nactivevars + 1);/*lint !e776*/
1464 
1465  MOSEK_CALL( MSK_gety(sdpisolver->msktask, MSK_SOL_ITR, moseksol) );/*lint !e641*/
1466 
1467  *feasorig = (moseksol[sdpisolver->nactivevars] < sdpisolver->feastol); /*lint !e413*/
1468 
1469  /* only set sdpisolver->feasorig to true if we solved with objective, because only in this case we want to compute
1470  * the objective value by hand since it is numerically more stable then the result returned by MOSEK */
1471  if ( withobj )
1472  sdpisolver->feasorig = *feasorig;
1473 
1474  /* if r > 0 also check the primal bound */
1475  if ( ! *feasorig && penaltybound != NULL )
1476  {
1477 
1478  SCIPdebugMessage("Solution not feasible in original problem, r = %f\n", moseksol[sdpisolver->nactivevars]);
1479 
1480  /* compute Tr(X) */
1481 
1482  /* start with the diagonal entries of the primal semidefinite variables */
1483  for (b = 0; b < nsdpblocks; b++)
1484  {
1485  if ( blockindchanges[b] > -1 )
1486  {
1487  SCIP_Real* X; /* the upper triangular entries of matrix X */
1488  int size;
1489 
1490  size = sdpblocksizes[b] - nremovedinds[b];
1491 
1492  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &X, 0.5 * size * (size + 1)) );
1493  MOSEK_CALL( MSK_getbarxj(sdpisolver->msktask, MSK_SOL_ITR, b - blockindchanges[b], X) );/*lint !e641*/
1494 
1495  /* iterate over all diagonal entries */
1496  for (i = 0; i < size; i++)
1497  {
1498  /* get index in the lower triangular part */
1499  ind = i * (i + 3) / 2;/*lint !e776*/ /* i*(i+1)/2 + i */
1500  assert( ind < 0.5 * size * (size + 1) );
1501  trace += X[ind];
1502  }
1503 
1504  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &X);
1505  }
1506  }
1507 
1508  /* add primal lp-variables */
1509  BMS_CALL( BMSallocBufferMemoryArray(sdpisolver->bufmem, &x, nlpvars + sdpisolver->nvarbounds) );
1510 
1511  MOSEK_CALL( MSK_getxx(sdpisolver->msktask, MSK_SOL_ITR, x) );/*lint !e641*/
1512 
1513  for (i = 0; i < nlpvars; i++)
1514  trace += x[i];
1515 
1516  BMSfreeBufferMemoryArrayNull(sdpisolver->bufmem, &x);
1517 
1518  /* if the relative gap is smaller than the tolerance, we return equality */
1519  if ( (penaltyparam - trace) / penaltyparam < PENALTYBOUNDTOL )/*lint !e414*/
1520  {
1521  assert( penaltybound != NULL );
1522  *penaltybound = TRUE;
1523  SCIPdebugMessage("Tr(X) = %f == %f = Gamma, penalty formulation not exact, Gamma should be increased or problem is infeasible\n",
1524  trace, penaltyparam);
1525  }
1526  else
1527  *penaltybound = FALSE;
1528  }
1529  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &moseksol);
1530  }
1531 
1532  /* free memory */
1533  BMSfreeBufferMemoryArrayNull(sdpisolver->bufmem, &mosekblocksizes);
1534  if ( nlpcons > 0 )
1535  {
1536  BMSfreeBufferMemoryArrayNull(sdpisolver->bufmem, &vartolhsrhsmapper);
1537  BMSfreeBufferMemoryArrayNull(sdpisolver->bufmem, &vartorowmapper);
1538  }
1539 
1540  return SCIP_OKAY;
1541 }
1547 /*
1548  * Solution Information Methods
1549  */
1550 
1555 SCIP_Bool SCIPsdpiSolverWasSolved(
1556  SCIP_SDPISOLVER* sdpisolver
1557  )
1558 {
1559  assert( sdpisolver != NULL );
1560 
1561  return sdpisolver->solved;
1562 }
1563 
1571  SCIP_SDPISOLVER* sdpisolver
1572  )
1573 {
1574  MSKsolstae solstat;
1575 
1576  assert( sdpisolver != NULL );
1577  CHECK_IF_SOLVED_BOOL( sdpisolver );
1578 
1579  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1580 
1581  switch ( solstat )
1582  {
1583  case MSK_SOL_STA_UNKNOWN:
1584  case MSK_SOL_STA_PRIM_FEAS:
1585  case MSK_SOL_STA_DUAL_FEAS:
1586  return FALSE;
1587  case MSK_SOL_STA_OPTIMAL:
1588  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1589  case MSK_SOL_STA_PRIM_INFEAS_CER:
1590  case MSK_SOL_STA_DUAL_INFEAS_CER:
1591  return TRUE;
1592  default:
1593  SCIPdebugMessage("Unknown return code in SCIPsdpiSolverFeasibilityKnown\n"); /* TODO: add illposed_cer */
1594  return FALSE;
1595  }/*lint !e788*/
1596 }
1597 
1599 SCIP_RETCODE SCIPsdpiSolverGetSolFeasibility(
1600  SCIP_SDPISOLVER* sdpisolver,
1601  SCIP_Bool* primalfeasible,
1602  SCIP_Bool* dualfeasible
1603  )
1604 {
1605  MSKsolstae solstat;
1606 
1607  assert( sdpisolver != NULL );
1608  assert( primalfeasible != NULL );
1609  assert( dualfeasible != NULL );
1610  CHECK_IF_SOLVED( sdpisolver );
1611 
1612  MOSEK_CALL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1613 
1614  switch ( solstat )
1615  {
1616  case MSK_SOL_STA_OPTIMAL:
1617  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1618  *primalfeasible = TRUE;
1619  *dualfeasible = TRUE;
1620  break;
1621  case MSK_SOL_STA_PRIM_INFEAS_CER:
1622  *primalfeasible = FALSE;
1623  *dualfeasible = TRUE;
1624  break;
1625  case MSK_SOL_STA_DUAL_INFEAS_CER:
1626  *primalfeasible = TRUE;
1627  *dualfeasible = FALSE;
1628  break;
1629  default:
1630  SCIPdebugMessage("MOSEK does not know about feasibility of solutions\n");
1631  return SCIP_LPERROR;
1632  }/*lint !e788*/
1633 
1634  return SCIP_OKAY;
1635 }
1636 
1641  SCIP_SDPISOLVER* sdpisolver
1642  )
1643 {
1644  MSKsolstae solstat;
1645 
1646  assert( sdpisolver != NULL );
1647  CHECK_IF_SOLVED_BOOL( sdpisolver );
1648 
1649  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1650 
1651  switch ( solstat )
1652  {
1653  case MSK_SOL_STA_DUAL_INFEAS_CER:
1654  return TRUE;
1655  case MSK_SOL_STA_OPTIMAL:
1656  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1657  case MSK_SOL_STA_PRIM_INFEAS_CER:
1658  break;
1659  default:
1660  SCIPdebugMessage("MOSEK does not know about feasibility of solutions\n");
1661  break;
1662  }/*lint !e788*/
1663  return FALSE;
1664 }
1665 
1670  SCIP_SDPISOLVER* sdpisolver
1671  )
1672 {
1673  MSKsolstae solstat;
1674 
1675  assert( sdpisolver != NULL );
1676  CHECK_IF_SOLVED_BOOL( sdpisolver );
1677 
1678  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1679 
1680  switch ( solstat )
1681  {
1682  case MSK_SOL_STA_PRIM_INFEAS_CER:
1683  return TRUE;
1684  case MSK_SOL_STA_OPTIMAL:
1685  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1686  case MSK_SOL_STA_DUAL_INFEAS_CER:
1687  break;
1688  default:
1689  SCIPdebugMessage("MOSEK does not know about feasibility of solutions\n");
1690  break;
1691  }/*lint !e788*/
1692  return FALSE;
1693 }
1694 
1699  SCIP_SDPISOLVER* sdpisolver
1700  )
1701 {
1702  MSKsolstae solstat;
1703 
1704  assert( sdpisolver != NULL );
1705  CHECK_IF_SOLVED_BOOL( sdpisolver );
1706 
1707  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1708 
1709  switch ( solstat )
1710  {
1711  case MSK_SOL_STA_OPTIMAL:
1712  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1713  case MSK_SOL_STA_DUAL_INFEAS_CER:
1714  return TRUE;
1715  case MSK_SOL_STA_PRIM_INFEAS_CER:
1716  break;
1717  default:
1718  SCIPdebugMessage("MOSEK does not know about feasibility of solutions\n");
1719  break;
1720  }/*lint !e788*/
1721  return FALSE;
1722 }
1723 
1728  SCIP_SDPISOLVER* sdpisolver
1729  )
1730 {
1731  MSKsolstae solstat;
1732 
1733  assert( sdpisolver != NULL );
1734  CHECK_IF_SOLVED_BOOL( sdpisolver );
1735 
1736  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1737 
1738  switch ( solstat )
1739  {
1740  case MSK_SOL_STA_PRIM_INFEAS_CER:
1741  return TRUE;
1742  case MSK_SOL_STA_OPTIMAL:
1743  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1744  case MSK_SOL_STA_DUAL_INFEAS_CER:
1745  break;
1746  default:
1747  SCIPdebugMessage("MOSEK does not know about feasibility of solutions\n");
1748  break;
1749  }/*lint !e788*/
1750  return FALSE;
1751 }
1752 
1757  SCIP_SDPISOLVER* sdpisolver
1758  )
1759 {
1760  MSKsolstae solstat;
1761 
1762  assert( sdpisolver != NULL );
1763  CHECK_IF_SOLVED_BOOL( sdpisolver );
1764 
1765  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1766 
1767  switch ( solstat )
1768  {
1769  case MSK_SOL_STA_DUAL_INFEAS_CER:
1770  return TRUE;
1771  case MSK_SOL_STA_OPTIMAL:
1772  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1773  case MSK_SOL_STA_PRIM_INFEAS_CER:
1774  break;
1775  default:
1776  SCIPdebugMessage("MOSEK does not know about feasibility of solutions\n");
1777  break;
1778  }/*lint !e788*/
1779  return FALSE;
1780 }
1781 
1786  SCIP_SDPISOLVER* sdpisolver
1787  )
1788 {
1789  MSKsolstae solstat;
1790 
1791  assert( sdpisolver != NULL );
1792  CHECK_IF_SOLVED_BOOL( sdpisolver );
1793 
1794  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1795 
1796  switch ( solstat )
1797  {
1798  case MSK_SOL_STA_OPTIMAL:
1799  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1800  case MSK_SOL_STA_PRIM_INFEAS_CER:
1801  return TRUE;
1802  case MSK_SOL_STA_DUAL_INFEAS_CER:
1803  break;
1804  default:
1805  SCIPdebugMessage("MOSEK does not know about feasibility of solutions\n");
1806  break;
1807  }/*lint !e788*/
1808  return FALSE;
1809 }
1810 
1812 SCIP_Bool SCIPsdpiSolverIsConverged(
1813  SCIP_SDPISOLVER* sdpisolver
1814  )
1815 {
1816  assert( sdpisolver != NULL );
1817 
1818  if ( sdpisolver->timelimit )
1819  return FALSE;
1820 
1821  CHECK_IF_SOLVED_BOOL( sdpisolver );
1822 
1823  /* check if Mosek stalled when it was already acceptable */
1824  if ( sdpisolver->terminationcode == MSK_RES_TRM_STALL )
1825  {
1826  MSKsolstae solstat;
1827  SCIP_Real pobj;
1828  SCIP_Real dobj;
1829  SCIP_Real gapnormalization;
1830 
1831  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );
1832 
1833  /* check the solution status */
1834  switch ( solstat )
1835  {
1836  case MSK_SOL_STA_UNKNOWN:
1837  case MSK_SOL_STA_PRIM_FEAS:
1838  case MSK_SOL_STA_DUAL_FEAS:
1839  return FALSE;
1840  case MSK_SOL_STA_OPTIMAL:
1841  case MSK_SOL_STA_PRIM_AND_DUAL_FEAS:
1842  case MSK_SOL_STA_PRIM_INFEAS_CER:
1843  case MSK_SOL_STA_DUAL_INFEAS_CER:
1844  /* check duality gap */
1845  MOSEK_CALL_BOOL( MSK_getdualobj(sdpisolver->msktask, MSK_SOL_ITR, &dobj) );
1846  MOSEK_CALL_BOOL( MSK_getprimalobj(sdpisolver->msktask, MSK_SOL_ITR, &pobj) );
1847  /* for the relative gap we divide by max(1.0, min(pobj, dobj)), as this is also done in Mosek */
1848  gapnormalization = dobj > pobj ? (pobj > 1.0 ? pobj : 1.0) : (dobj > 1.0 ? dobj : 1.0);
1849  if ( REALABS((pobj-dobj) / gapnormalization) < sdpisolver->gaptol )
1850  return TRUE;
1851  else
1852  return FALSE;
1853  default:
1854  return FALSE;
1855  }
1856  }
1857 
1858  return sdpisolver->terminationcode == MSK_RES_OK;
1859 }
1860 
1862 SCIP_Bool SCIPsdpiSolverIsObjlimExc(
1863  SCIP_SDPISOLVER* sdpisolver
1864  )
1865 {
1866  assert( sdpisolver != NULL );
1867  CHECK_IF_SOLVED_BOOL( sdpisolver );
1868 
1869  return sdpisolver->terminationcode == MSK_RES_TRM_OBJECTIVE_RANGE;
1870 }
1871 
1873 SCIP_Bool SCIPsdpiSolverIsIterlimExc(
1874  SCIP_SDPISOLVER* sdpisolver
1875  )
1876 {
1877  assert( sdpisolver != NULL );
1878  CHECK_IF_SOLVED_BOOL( sdpisolver );
1879 
1880  return sdpisolver->terminationcode == MSK_RES_TRM_MAX_ITERATIONS;
1881 }
1882 
1884 SCIP_Bool SCIPsdpiSolverIsTimelimExc(
1885  SCIP_SDPISOLVER* sdpisolver
1886  )
1887 {
1888  assert( sdpisolver != NULL );
1889 
1890  if ( sdpisolver->timelimit )
1891  return TRUE;
1892 
1893  if ( ! sdpisolver->solved )
1894  return FALSE;
1895 
1896  return sdpisolver->terminationcode == MSK_RES_TRM_MAX_TIME;
1897 }
1898 
1911  SCIP_SDPISOLVER* sdpisolver
1912  )
1913 {
1914  assert( sdpisolver != NULL );
1915 
1916  if ( ! sdpisolver->solved )
1917  return -1;
1918 
1919  if ( sdpisolver->timelimit )
1920  return 5;
1921 
1922  switch ( sdpisolver->terminationcode )
1923  {
1924  case MSK_RES_OK:
1925  return 0;
1926  case MSK_RES_TRM_MAX_NUM_SETBACKS:
1927  case MSK_RES_TRM_NUMERICAL_PROBLEM:
1928  case MSK_RES_TRM_STALL:
1929  return 2;
1930  case MSK_RES_TRM_OBJECTIVE_RANGE:
1931  return 3;
1932  case MSK_RES_TRM_MAX_ITERATIONS:
1933  return 4;
1934  case MSK_RES_TRM_MAX_TIME:
1935  return 5;
1936  default:
1937  return 7;
1938  }/*lint !e788*/
1939 }
1940 
1942 SCIP_Bool SCIPsdpiSolverIsOptimal(
1943  SCIP_SDPISOLVER* sdpisolver
1944  )
1945 {
1946  MSKsolstae solstat;
1947 
1948  assert( sdpisolver != NULL );
1949 
1950  if ( sdpisolver->timelimit )
1951  return FALSE;
1952 
1953  CHECK_IF_SOLVED_BOOL( sdpisolver );
1954 
1955  if ( sdpisolver->terminationcode != MSK_RES_OK )
1956  return FALSE;
1957 
1958  MOSEK_CALL_BOOL( MSK_getsolsta(sdpisolver->msktask, MSK_SOL_ITR, &solstat) );/*lint !e641*/
1959 
1960  if ( solstat != MSK_SOL_STA_OPTIMAL )
1961  return FALSE;
1962 
1963  return TRUE;
1964 }
1965 
1969 SCIP_Bool SCIPsdpiSolverIsAcceptable(
1970  SCIP_SDPISOLVER* sdpisolver
1971  )
1972 {
1973  assert( sdpisolver != NULL );
1974 
1975  if ( sdpisolver->timelimit )
1976  return FALSE;
1977 
1978  if ( ! sdpisolver->solved )
1979  return FALSE;
1980 
1981  return SCIPsdpiSolverIsConverged(sdpisolver) && SCIPsdpiSolverFeasibilityKnown(sdpisolver);
1982 }
1983 
1985 SCIP_RETCODE SCIPsdpiSolverIgnoreInstability(
1986  SCIP_SDPISOLVER* sdpisolver,
1987  SCIP_Bool* success
1988  )
1989 {/*lint --e{715}*/
1990  SCIPdebugMessage("Not implemented yet\n");
1991  return SCIP_LPERROR;
1992 }
1993 
1995 SCIP_RETCODE SCIPsdpiSolverGetObjval(
1996  SCIP_SDPISOLVER* sdpisolver,
1997  SCIP_Real* objval
1998  )
1999 {
2000  SCIP_Real* moseksol;
2001 
2002  assert( sdpisolver != NULL );
2003  CHECK_IF_SOLVED( sdpisolver );
2004  assert( objval != NULL );
2005 
2006  /* check for unboundedness */
2007  if ( SCIPsdpiSolverIsDualUnbounded(sdpisolver) )
2008  {
2009  *objval = -SCIPsdpiSolverInfinity(sdpisolver);
2010  return SCIP_OKAY;
2011  }
2012 
2013  if ( sdpisolver->penalty && ( ! sdpisolver->feasorig ) )
2014  {
2015  /* in this case we cannot really trust the solution given by MOSEK, since changes in the value of r much less than epsilon can
2016  * cause huge changes in the objective, so using the objective value given by MOSEK is numerically more stable */
2017  MOSEK_CALL( MSK_getdualobj(sdpisolver->msktask, MSK_SOL_ITR, objval) );
2018  }
2019  else
2020  {
2021  int v;
2022 
2023  /* since the objective value given by MOSEK sometimes differs slightly from the correct value for the given solution,
2024  * we get the solution from MOSEK and compute the correct objective value */
2025  BMSallocBufferMemoryArray(sdpisolver->bufmem, &moseksol, sdpisolver->penalty ? sdpisolver->nactivevars + 1 : sdpisolver->nactivevars);
2026  MOSEK_CALL( MSK_gety(sdpisolver->msktask, MSK_SOL_ITR, moseksol) );/*lint !e641*/
2027 
2028  *objval = 0.0;
2029  for (v = 0; v < sdpisolver->nactivevars; v++)
2030  *objval += moseksol[v] * sdpisolver->objcoefs[v];
2031  }
2032 
2033  /* as we didn't add the fixed (lb = ub) variables to MOSEK, we have to add their contributions to the objective as well */
2034  *objval += sdpisolver->fixedvarsobjcontr;
2035 
2036  if ( ( ! sdpisolver->penalty ) || sdpisolver->feasorig)
2037  {
2038  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &moseksol);
2039  }
2040 
2041  return SCIP_OKAY;
2042 }
2043 
2048 SCIP_RETCODE SCIPsdpiSolverGetSol(
2049  SCIP_SDPISOLVER* sdpisolver,
2050  SCIP_Real* objval,
2051  SCIP_Real* dualsol,
2052  int* dualsollength
2054  )
2055 {
2056  int v;
2057  SCIP_Real* moseksol;
2058 
2059  assert( sdpisolver != NULL );
2060  CHECK_IF_SOLVED( sdpisolver );
2061  assert( dualsollength != NULL );
2062 
2063  if ( *dualsollength > 0 )
2064  {
2065  assert( dualsol != NULL );
2066  if ( *dualsollength < sdpisolver->nvars )
2067  {
2068  SCIPdebugMessage("The given array in SCIPsdpiSolverGetSol only had length %d, but %d was needed", *dualsollength, sdpisolver->nvars);
2069  *dualsollength = sdpisolver->nvars;
2070 
2071  return SCIP_OKAY;
2072  }
2073 
2074  BMSallocBufferMemoryArray(sdpisolver->bufmem, &moseksol, sdpisolver->penalty ? sdpisolver->nactivevars + 1 : sdpisolver->nactivevars);
2075 
2076  MOSEK_CALL( MSK_gety(sdpisolver->msktask, MSK_SOL_ITR, moseksol) );/*lint !e641*/
2077 
2078  /* insert the entries into dualsol, for non-fixed vars we copy those from MOSEK, the rest are the saved entries from inserting (they equal lb=ub) */
2079  for (v = 0; v < sdpisolver->nvars; v++)
2080  {
2081  if ( sdpisolver->inputtomosekmapper[v] >= 0 )
2082  dualsol[v] = moseksol[sdpisolver->inputtomosekmapper[v]];
2083  else
2084  {
2085  /* this is the value that was saved when inserting, as this variable has lb=ub */
2086  assert( -sdpisolver->inputtomosekmapper[v] <= sdpisolver->nvars - sdpisolver->nactivevars );
2087  dualsol[v] = sdpisolver->fixedvarsval[(-1 * sdpisolver->inputtomosekmapper[v]) - 1]; /*lint !e679*/ /* -1 because we wanted strictly negative vals */
2088  }
2089  }
2090 
2091  /* if both solution and objective should be printed, we can use the solution to compute the objective */
2092  if ( objval != NULL )
2093  {
2094  if ( sdpisolver->penalty && ( ! sdpisolver->feasorig ))
2095  {
2096  /* in this case we cannot really trust the solution given by MOSEK, since changes in the value of r much less than epsilon can
2097  * cause huge changes in the objective, so using the objective value given by MOSEK is numerically more stable */
2098  MOSEK_CALL( MSK_getdualobj(sdpisolver->msktask, MSK_SOL_ITR, objval) );
2099  }
2100  else
2101  {
2102  /* since the objective value given by MOSEK sometimes differs slightly from the correct value for the given solution,
2103  * we get the solution from MOSEK and compute the correct objective value */
2104  *objval = 0.0;
2105  for (v = 0; v < sdpisolver->nactivevars; v++)
2106  *objval += moseksol[v] * sdpisolver->objcoefs[v];
2107  }
2108 
2109  /* as we didn't add the fixed (lb = ub) variables to MOSEK, we have to add their contributions to the objective as well */
2110  *objval += sdpisolver->fixedvarsobjcontr;
2111  }
2112 
2113  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &moseksol);
2114  }
2115  else if ( objval != NULL )
2116  {
2117  SCIP_CALL( SCIPsdpiSolverGetObjval(sdpisolver, objval) );
2118  }
2119 
2120  return SCIP_OKAY;
2121 }
2122 
2125  SCIP_SDPISOLVER* sdpisolver,
2126  int nblocks,
2127  int* startXnblocknonz
2129  )
2130 { /*lint --e{715}*/
2131  SCIPdebugMessage("Not implemented yet\n");
2132 
2133  return SCIP_PLUGINNOTFOUND;
2134 }
2135 
2143 SCIP_RETCODE SCIPsdpiSolverGetPreoptimalSol(
2144  SCIP_SDPISOLVER* sdpisolver,
2145  SCIP_Bool* success,
2146  SCIP_Real* dualsol,
2147  int* dualsollength,
2149  int nblocks,
2150  int* startXnblocknonz,
2152  int** startXrow,
2153  int** startXcol,
2154  SCIP_Real** startXval
2155  )
2156 {/*lint !e1784*/
2157  SCIPdebugMessage("Not implemented yet\n");
2158  return SCIP_LPERROR;
2159 }/*lint !e715*/
2160 
2169  SCIP_SDPISOLVER* sdpisolver,
2170  SCIP_Real* lbvars,
2171  SCIP_Real* ubvars,
2172  int* arraylength
2174  )
2175 {
2176  SCIP_Real* primalvars;
2177  int nprimalvars;
2178  int i;
2179 
2180  assert( sdpisolver != NULL );
2181  CHECK_IF_SOLVED( sdpisolver );
2182  assert( arraylength != NULL );
2183  assert( lbvars != NULL );
2184  assert( ubvars != NULL );
2185 
2186  /* check if the arrays are long enough */
2187  if ( *arraylength < sdpisolver->nvars )
2188  {
2189  *arraylength = sdpisolver->nvars;
2190  SCIPdebugMessage("Insufficient length of array in SCIPsdpiSolverGetPrimalBoundVars (gave %d, needed %d)\n", *arraylength, sdpisolver->nvars);
2191  return SCIP_OKAY;
2192  }
2193 
2194  /* initialize the return-arrays with zero */
2195  for (i = 0; i < sdpisolver->nvars; i++)
2196  {
2197  lbvars[i] = 0.0;
2198  ubvars[i] = 0.0;
2199  }
2200 
2201  /* get number of primal variables in MOSEK */
2202  MOSEK_CALL( MSK_getnumvar(sdpisolver->msktask, &nprimalvars) );/*lint !e641*/
2203 
2204  BMSallocBufferMemoryArray(sdpisolver->bufmem, &primalvars, nprimalvars);
2205 
2206  MOSEK_CALL( MSK_getxx(sdpisolver->msktask, MSK_SOL_ITR, primalvars) );/*lint !e641*/
2207 
2208  /* iterate over all variable bounds and insert the corresponding primal variables in the right positions of the return-arrays */
2209  assert( sdpisolver->nvarbounds <= 2 * sdpisolver->nvars );
2210 
2211  for (i = 0; i < sdpisolver->nvarbounds; i++)
2212  {
2213  if ( sdpisolver->varboundpos[i] < 0 )
2214  {
2215  /* this is a lower bound */
2216  /* the last nvarbounds entries correspond to the varbounds */
2217  lbvars[sdpisolver->mosektoinputmapper[-1 * sdpisolver->varboundpos[i] -1]] = primalvars[nprimalvars - sdpisolver->nvarbounds + i]; /*lint !e679, !e834 */
2218  }
2219  else
2220  {
2221  /* this is an upper bound */
2222 
2223  assert( sdpisolver->varboundpos[i] > 0 );
2224 
2225  /* the last nvarbounds entries correspond to the varbounds */
2226  ubvars[sdpisolver->mosektoinputmapper[sdpisolver->varboundpos[i] - 1]] = primalvars[nprimalvars - sdpisolver->nvarbounds + i]; /*lint !e679, !e834 */
2227  }
2228  }
2229 
2230  BMSfreeBufferMemoryArray(sdpisolver->bufmem, &primalvars);
2231 
2232  return SCIP_OKAY;
2233 }
2234 
2236 SCIP_RETCODE SCIPsdpiSolverGetPrimalNonzeros(
2237  SCIP_SDPISOLVER* sdpisolver,
2238  int nblocks,
2239  int* startXnblocknonz
2240  )
2241 {/*lint --e{715}*/
2242  SCIPdebugMessage("Not implemented yet\n");
2243  return SCIP_LPERROR;
2244 }
2245 
2252 SCIP_RETCODE SCIPsdpiSolverGetPrimalMatrix(
2253  SCIP_SDPISOLVER* sdpisolver,
2254  int nblocks,
2255  int* startXnblocknonz,
2257  int** startXrow,
2258  int** startXcol,
2259  SCIP_Real** startXval
2260  )
2261 {/*lint --e{715}*/
2262  SCIPdebugMessage("Not implemented yet\n");
2263  return SCIP_LPERROR;
2264 }
2265 
2268  SCIP_SDPISOLVER* sdpisolver
2269  )
2270 {/*lint --e{715}*/
2271  SCIPdebugMessage("Not implemented yet\n");
2272  return SCIP_INVALID;
2273 }
2274 
2276 SCIP_RETCODE SCIPsdpiSolverGetIterations(
2277  SCIP_SDPISOLVER* sdpisolver,
2278  int* iterations
2279  )
2280 {
2281  if ( sdpisolver->timelimitinitial )
2282  *iterations = 0;
2283  else
2284  {
2285  *iterations = sdpisolver->niterations;
2286  }
2287 
2288  return SCIP_OKAY;
2289 }
2290 
2292 SCIP_RETCODE SCIPsdpiSolverGetSdpCalls(
2293  SCIP_SDPISOLVER* sdpisolver,
2294  int* calls
2295  )
2296 {/*lint --e{715,1784}*/
2297  assert( calls != NULL );
2298 
2299  *calls = sdpisolver->timelimitinitial ? 0 : sdpisolver->nsdpcalls;
2300 
2301  return SCIP_OKAY;
2302 }
2303 
2305 SCIP_RETCODE SCIPsdpiSolverSettingsUsed(
2306  SCIP_SDPISOLVER* sdpisolver,
2308  )
2309 {
2310  assert( sdpisolver != NULL );
2311  assert( usedsetting != NULL );
2312 
2313  if ( ! SCIPsdpiSolverIsAcceptable(sdpisolver) )
2314  *usedsetting = SCIP_SDPSOLVERSETTING_UNSOLVED;
2315  else if ( sdpisolver->penalty )
2316  *usedsetting = SCIP_SDPSOLVERSETTING_PENALTY;
2317  else
2318  *usedsetting = SCIP_SDPSOLVERSETTING_FAST;
2319 
2320  return SCIP_OKAY;
2321 }
2322 
2328 /*
2329  * Numerical Methods
2330  */
2331 
2336 SCIP_Real SCIPsdpiSolverInfinity(
2337  SCIP_SDPISOLVER* sdpisolver
2338  )
2339 {/*lint --e{715}*/
2340  return 1.0e16;
2341 }
2342 
2344 SCIP_Bool SCIPsdpiSolverIsInfinity(
2345  SCIP_SDPISOLVER* sdpisolver,
2346  SCIP_Real val
2347  )
2348 {
2349  return ((val <= -SCIPsdpiSolverInfinity(sdpisolver)) || (val >= SCIPsdpiSolverInfinity(sdpisolver)));
2350 }
2351 
2353 SCIP_RETCODE SCIPsdpiSolverGetRealpar(
2354  SCIP_SDPISOLVER* sdpisolver,
2356  SCIP_Real* dval
2357  )
2358 {
2359  assert( sdpisolver != NULL );
2360  assert( dval != NULL );
2361 
2362  switch( type )
2363  {
2364  case SCIP_SDPPAR_EPSILON:
2365  *dval = sdpisolver->epsilon;
2366  break;
2367  case SCIP_SDPPAR_GAPTOL:
2368  *dval = sdpisolver->gaptol;
2369  break;
2370  case SCIP_SDPPAR_FEASTOL:
2371  *dval = sdpisolver->feastol;
2372  break;
2374  *dval = sdpisolver->sdpsolverfeastol;
2375  break;
2376  case SCIP_SDPPAR_OBJLIMIT:
2377  *dval = sdpisolver->objlimit;
2378  break;
2379  default:
2380  return SCIP_PARAMETERUNKNOWN;
2381  }/*lint !e788*/
2382 
2383  return SCIP_OKAY;
2384 }
2385 
2387 SCIP_RETCODE SCIPsdpiSolverSetRealpar(
2388  SCIP_SDPISOLVER* sdpisolver,
2390  SCIP_Real dval
2391  )
2392 {
2393  assert( sdpisolver != NULL );
2394 
2395  switch( type )
2396  {
2397  case SCIP_SDPPAR_EPSILON:
2398  sdpisolver->epsilon = dval;
2399  SCIPdebugMessage("Setting sdpisolver epsilon to %f.\n", dval);
2400  break;
2401  case SCIP_SDPPAR_GAPTOL:
2402  sdpisolver->gaptol = dval;
2403  SCIPdebugMessage("Setting sdpisolver gaptol to %f.\n", dval);
2404  break;
2405  case SCIP_SDPPAR_FEASTOL:
2406  sdpisolver->feastol = dval;
2407  SCIPdebugMessage("Setting sdpisolver feastol to %f.\n", dval);
2408  break;
2410  sdpisolver->sdpsolverfeastol = dval;
2411  SCIPdebugMessage("Setting sdpisolver sdpsolverfeastol to %f.\n", dval);
2412  break;
2413  case SCIP_SDPPAR_OBJLIMIT:
2414  SCIPdebugMessage("Setting sdpisolver objlimit to %f.\n", dval);
2415  sdpisolver->objlimit = dval;
2416  break;
2417  default:
2418  return SCIP_PARAMETERUNKNOWN;
2419  }/*lint !e788*/
2420 
2421  return SCIP_OKAY;
2422 }
2423 
2425 SCIP_RETCODE SCIPsdpiSolverGetIntpar(
2426  SCIP_SDPISOLVER* sdpisolver,
2428  int* ival
2429  )
2430 {
2431  assert( sdpisolver != NULL );
2432 
2433  switch( type )
2434  {
2435  case SCIP_SDPPAR_SDPINFO:
2436  *ival = (int) sdpisolver->sdpinfo;
2437  SCIPdebugMessage("Getting sdpisolver information output (%d).\n", *ival);
2438  break;
2439  case SCIP_SDPPAR_NTHREADS:
2440  *ival = sdpisolver->nthreads;
2441  SCIPdebugMessage("Getting sdpisolver number of threads: %d.\n", *ival);
2442  break;
2443  default:
2444  return SCIP_PARAMETERUNKNOWN;
2445  }/*lint !e788*/
2446 
2447  return SCIP_OKAY;
2448 }
2449 
2451 SCIP_RETCODE SCIPsdpiSolverSetIntpar(
2452  SCIP_SDPISOLVER* sdpisolver,
2454  int ival
2455  )
2456 {
2457  assert( sdpisolver != NULL );
2458 
2459  switch( type )
2460  {
2461  case SCIP_SDPPAR_NTHREADS:
2462  sdpisolver->nthreads = ival;
2463  SCIPdebugMessage("Setting sdpisolver number of threads to %d.\n", ival);
2464  break;
2465  case SCIP_SDPPAR_SDPINFO:
2466  assert( 0 <= ival && ival <= 1 );
2467  sdpisolver->sdpinfo = (SCIP_Bool) ival;
2468  SCIPdebugMessage("Setting sdpisolver information output (%d).\n", ival);
2469  break;
2470  default:
2471  return SCIP_PARAMETERUNKNOWN;
2472  }/*lint !e788*/
2473 
2474  return SCIP_OKAY;
2475 }
2476 
2478 SCIP_RETCODE SCIPsdpiSolverComputeLambdastar(
2479  SCIP_SDPISOLVER* sdpisolver,
2480  SCIP_Real maxguess
2481  )
2482 {/*lint --e{715}*/
2483  SCIPdebugMessage("Lambdastar parameter not used by MOSEK"); /* this parameter is only used by SDPA */
2484 
2485  return SCIP_OKAY;
2486 }
2487 
2490  SCIP_SDPISOLVER* sdpisolver,
2491  SCIP_Real maxcoeff,
2492  SCIP_Real* penaltyparam
2493  )
2494 {/*lint --e{1784}*/
2495  SCIP_Real compval;
2496 
2497  assert( sdpisolver != NULL );
2498  assert( penaltyparam != NULL );
2499 
2500  compval = PENALTYPARAM_FACTOR * maxcoeff;
2501 
2502  if ( compval < MIN_PENALTYPARAM )
2503  {
2504  SCIPdebugMessage("Setting penaltyparameter to %f.\n", MIN_PENALTYPARAM);
2505  *penaltyparam = MIN_PENALTYPARAM;
2506  }
2507  else if ( compval > MAX_PENALTYPARAM )
2508  {
2509  SCIPdebugMessage("Setting penaltyparameter to %f.\n", MAX_PENALTYPARAM);
2510  *penaltyparam = MAX_PENALTYPARAM;
2511  }
2512  else
2513  {
2514  SCIPdebugMessage("Setting penaltyparameter to %f.\n", compval);
2515  *penaltyparam = compval;
2516  }
2517  return SCIP_OKAY;
2518 }
2519 
2522  SCIP_SDPISOLVER* sdpisolver,
2523  SCIP_Real penaltyparam,
2524  SCIP_Real* maxpenaltyparam
2525  )
2526 {/*lint --e{1784}*/
2527  SCIP_Real compval;
2528 
2529  assert( sdpisolver != NULL );
2530  assert( maxpenaltyparam != NULL );
2531 
2532  compval = penaltyparam * MAXPENALTYPARAM_FACTOR;
2533 
2534  if ( compval < MAX_MAXPENALTYPARAM )
2535  {
2536  *maxpenaltyparam = compval;
2537  SCIPdebugMessage("Setting maximum penaltyparameter to %f.\n", compval);
2538  }
2539  else
2540  {
2541  *maxpenaltyparam = MAX_MAXPENALTYPARAM;
2542  SCIPdebugMessage("Setting penaltyparameter to %f.\n", MAX_MAXPENALTYPARAM);
2543  }
2544 
2545  return SCIP_OKAY;
2546 }
2547 
2553 /*
2554  * File Interface Methods
2555  */
2556 
2561 SCIP_RETCODE SCIPsdpiSolverReadSDP(
2562  SCIP_SDPISOLVER* sdpisolver,
2563  const char* fname
2564  )
2565 {/*lint --e{715}*/
2566  SCIPdebugMessage("Not implemented yet\n");
2567  return SCIP_LPERROR;
2568 }
2569 
2571 SCIP_RETCODE SCIPsdpiSolverWriteSDP(
2572  SCIP_SDPISOLVER* sdpisolver,
2573  const char* fname
2574  )
2575 {
2576  assert( sdpisolver != NULL );
2577  assert( fname != NULL );
2578 
2579  MOSEK_CALL( MSK_writedata(sdpisolver->msktask, fname) );/*lint !e641*/
2580 
2581  return SCIP_OKAY;
2582 }
2583 
SCIP_RETCODE SCIPsdpiSolverCreate(SCIP_SDPISOLVER **sdpisolver, SCIP_MESSAGEHDLR *messagehdlr, BMS_BLKMEM *blkmem, BMS_BUFMEM *bufmem)
SCIP_Real SCIPsdpiSolverInfinity(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpiSolverIncreaseCounter(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpiSolverComputeLambdastar(SCIP_SDPISOLVER *sdpisolver, SCIP_Real maxguess)
SCIP_Bool SCIPsdpiSolverWasSolved(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpiSolverReadSDP(SCIP_SDPISOLVER *sdpisolver, const char *fname)
#define PENALTYPARAM_FACTOR
#define INFEASMINFEASTOL
SCIP_Bool SCIPsdpiSolverIsConverged(SCIP_SDPISOLVER *sdpisolver)
#define CHECK_IF_SOLVED(sdpisolver)
enum SCIP_SDPSolverSetting SCIP_SDPSOLVERSETTING
Definition: type_sdpi.h:80
SCIP_RETCODE SCIPsdpiSolverResetCounter(SCIP_SDPISOLVER *sdpisolver)
SCIP_Bool SCIPsdpiSolverIsDualFeasible(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpiSolverIgnoreInstability(SCIP_SDPISOLVER *sdpisolver, SCIP_Bool *success)
int SCIPsdpiSolverGetInternalStatus(SCIP_SDPISOLVER *sdpisolver)
void * SCIPsdpiSolverGetSolverPointer(SCIP_SDPISOLVER *sdpisolver)
#define MIN_PENALTYPARAM
SCIP_RETCODE SCIPsdpiSolverWriteSDP(SCIP_SDPISOLVER *sdpisolver, const char *fname)
interface methods for specific SDP-solvers
SCIP_RETCODE SCIPsdpiSolverComputeMaxPenaltyparam(SCIP_SDPISOLVER *sdpisolver, SCIP_Real penaltyparam, SCIP_Real *maxpenaltyparam)
SCIP_Bool SCIPsdpiSolverIsInfinity(SCIP_SDPISOLVER *sdpisolver, SCIP_Real val)
SCIP_Bool SCIPsdpiSolverIsPrimalFeasible(SCIP_SDPISOLVER *sdpisolver)
#define MAX_MAXPENALTYPARAM
SCIP_RETCODE SCIPsdpiSolverSettingsUsed(SCIP_SDPISOLVER *sdpisolver, SCIP_SDPSOLVERSETTING *usedsetting)
SCIP_RETCODE SCIPsdpiSolverFree(SCIP_SDPISOLVER **sdpisolver)
#define MAXPENALTYPARAM_FACTOR
#define MOSEK_CALL_BOOL(x)
SCIP_Bool SCIPsdpiSolverIsPrimalUnbounded(SCIP_SDPISOLVER *sdpisolver)
#define TIMEOFDAY_CALL(x)
SCIP_RETCODE SCIPsdpiSolverGetIterations(SCIP_SDPISOLVER *sdpisolver, int *iterations)
SCIP_RETCODE SCIPsdpiSolverGetPrimalMatrix(SCIP_SDPISOLVER *sdpisolver, int nblocks, int *startXnblocknonz, int **startXrow, int **startXcol, SCIP_Real **startXval)
#define PENALTYBOUNDTOL
SCIP_RETCODE SCIPsdpiSolverLoadAndSolve(SCIP_SDPISOLVER *sdpisolver, int nvars, SCIP_Real *obj, SCIP_Real *lb, SCIP_Real *ub, int nsdpblocks, int *sdpblocksizes, int *sdpnblockvars, int sdpconstnnonz, int *sdpconstnblocknonz, int **sdpconstrow, int **sdpconstcol, SCIP_Real **sdpconstval, int sdpnnonz, int **sdpnblockvarnonz, int **sdpvar, int ***sdprow, int ***sdpcol, SCIP_Real ***sdpval, int **indchanges, int *nremovedinds, int *blockindchanges, int nremovedblocks, int nlpcons, int noldlpcons, SCIP_Real *lplhs, SCIP_Real *lprhs, int *rownactivevars, int lpnnonz, int *lprow, int *lpcol, SCIP_Real *lpval, SCIP_Real *starty, int *startZnblocknonz, int **startZrow, int **startZcol, SCIP_Real **startZval, int *startXnblocknonz, int **startXrow, int **startXcol, SCIP_Real **startXval, SCIP_SDPSOLVERSETTING startsettings, SCIP_Real timelimit)
SCIP_RETCODE SCIPsdpiSolverComputePenaltyparam(SCIP_SDPISOLVER *sdpisolver, SCIP_Real maxcoeff, SCIP_Real *penaltyparam)
SCIP_RETCODE SCIPsdpiSolverSetIntpar(SCIP_SDPISOLVER *sdpisolver, SCIP_SDPPARAM type, int ival)
SCIP_Bool SCIPsdpiSolverIsPrimalInfeasible(SCIP_SDPISOLVER *sdpisolver)
int SCIPsdpiSolverGetDefaultSdpiSolverNpenaltyIncreases(void)
SCIP_Bool SCIPsdpiSolverIsOptimal(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpSolcheckerCheck(BMS_BUFMEM *bufmem, int nvars, SCIP_Real *lb, SCIP_Real *ub, int nsdpblocks, int *sdpblocksizes, int *sdpnblockvars, int sdpconstnnonz, int *sdpconstnblocknonz, int **sdpconstrow, int **sdpconstcol, SCIP_Real **sdpconstval, int sdpnnonz, int **sdpnblockvarnonz, int **sdpvar, int ***sdprow, int ***sdpcol, SCIP_Real ***sdpval, int **indchanges, int *nremovedinds, int *blockindchanges, int nlpcons, int noldlpcons, SCIP_Real *lplhs, SCIP_Real *lprhs, int *rownactivevars, int lpnnonz, int *lprow, int *lpcol, SCIP_Real *lpval, SCIP_Real *solvector, SCIP_Real feastol, SCIP_Real epsilon, SCIP_Bool *infeasible)
Definition: sdpsolchecker.c:58
#define INFEASFEASTOLCHANGE
SCIP_Real SCIPsdpiSolverGetDefaultSdpiSolverGaptol(void)
SCIP_RETCODE SCIPsdpiSolverGetSol(SCIP_SDPISOLVER *sdpisolver, SCIP_Real *objval, SCIP_Real *dualsol, int *dualsollength)
SCIP_RETCODE SCIPsdpiSolverGetObjval(SCIP_SDPISOLVER *sdpisolver, SCIP_Real *objval)
#define BMS_CALL(x)
checks a given SDP solution for feasibility
static void MSKAPI printstr(void *handle, MSKCONST char str[])
const char * SCIPsdpiSolverGetSolverDesc(void)
SCIP_RETCODE SCIPsdpiSolverGetPreoptimalPrimalNonzeros(SCIP_SDPISOLVER *sdpisolver, int nblocks, int *startXnblocknonz)
SCIP_Real SCIPsdpiSolverGetDefaultSdpiSolverFeastol(void)
SCIP_Bool SCIPsdpiSolverIsObjlimExc(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpiSolverGetPrimalNonzeros(SCIP_SDPISOLVER *sdpisolver, int nblocks, int *startXnblocknonz)
static SCIP_Bool isFixed(SCIP_SDPISOLVER *sdpisolver, SCIP_Real lb, SCIP_Real ub)
SCIP_RETCODE SCIPsdpiSolverGetRealpar(SCIP_SDPISOLVER *sdpisolver, SCIP_SDPPARAM type, SCIP_Real *dval)
SCIP_Bool SCIPsdpiSolverIsDualInfeasible(SCIP_SDPISOLVER *sdpisolver)
#define CHECK_IF_SOLVED_BOOL(sdpisolver)
SCIP_Bool SCIPsdpiSolverIsIterlimExc(SCIP_SDPISOLVER *sdpisolver)
#define MOSEK_CALLM(x)
SCIP_Real SCIPsdpiSolverGetMaxPrimalEntry(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpiSolverGetPrimalBoundVars(SCIP_SDPISOLVER *sdpisolver, SCIP_Real *lbvars, SCIP_Real *ubvars, int *arraylength)
const char * SCIPsdpiSolverGetSolverName(void)
SCIP_RETCODE SCIPsdpiSolverGetSdpCalls(SCIP_SDPISOLVER *sdpisolver, int *calls)
SCIP_RETCODE SCIPsdpiSolverGetSolFeasibility(SCIP_SDPISOLVER *sdpisolver, SCIP_Bool *primalfeasible, SCIP_Bool *dualfeasible)
SCIP_Bool SCIPsdpiSolverIsAcceptable(SCIP_SDPISOLVER *sdpisolver)
struct SCIP_SDPiSolver SCIP_SDPISOLVER
Definition: sdpisolver.h:70
SCIP_RETCODE SCIPsdpiSolverSetRealpar(SCIP_SDPISOLVER *sdpisolver, SCIP_SDPPARAM type, SCIP_Real dval)
SCIP_RETCODE SCIPsdpiSolverGetIntpar(SCIP_SDPISOLVER *sdpisolver, SCIP_SDPPARAM type, int *ival)
#define MOSEK_CALL(x)
SCIP_RETCODE SCIPsdpiSolverGetPreoptimalSol(SCIP_SDPISOLVER *sdpisolver, SCIP_Bool *success, SCIP_Real *dualsol, int *dualsollength, int nblocks, int *startXnblocknonz, int **startXrow, int **startXcol, SCIP_Real **startXval)
enum SCIP_SDPParam SCIP_SDPPARAM
Definition: type_sdpi.h:69
SCIP_Bool SCIPsdpiSolverIsDualUnbounded(SCIP_SDPISOLVER *sdpisolver)
SCIP_RETCODE SCIPsdpiSolverLoadAndSolveWithPenalty(SCIP_SDPISOLVER *sdpisolver, SCIP_Real penaltyparam, SCIP_Bool withobj, SCIP_Bool rbound, int nvars, SCIP_Real *obj, SCIP_Real *lb, SCIP_Real *ub, int nsdpblocks, int *sdpblocksizes, int *sdpnblockvars, int sdpconstnnonz, int *sdpconstnblocknonz, int **sdpconstrow, int **sdpconstcol, SCIP_Real **sdpconstval, int sdpnnonz, int **sdpnblockvarnonz, int **sdpvar, int ***sdprow, int ***sdpcol, SCIP_Real ***sdpval, int **indchanges, int *nremovedinds, int *blockindchanges, int nremovedblocks, int nlpcons, int noldlpcons, SCIP_Real *lplhs, SCIP_Real *lprhs, int *rownactivevars, int lpnnonz, int *lprow, int *lpcol, SCIP_Real *lpval, SCIP_Real *starty, int *startZnblocknonz, int **startZrow, int **startZcol, SCIP_Real **startZval, int *startXnblocknonz, int **startXrow, int **startXcol, SCIP_Real **startXval, SCIP_SDPSOLVERSETTING startsettings, SCIP_Real timelimit, SCIP_Bool *feasorig, SCIP_Bool *penaltybound)
SCIP_Bool SCIPsdpiSolverIsTimelimExc(SCIP_SDPISOLVER *sdpisolver)
char name[SCIP_MAXSTRLEN]
SCIP_Bool SCIPsdpiSolverDoesWarmstartNeedPrimal(void)
SCIP_Bool SCIPsdpiSolverFeasibilityKnown(SCIP_SDPISOLVER *sdpisolver)
#define MAX_PENALTYPARAM