SCIP-SDP  2.0.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
cons_sdp.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 programms based on SCIP. */
5 /* */
6 /* Copyright (C) 2011-2013 Discrete Optimization, TU Darmstadt */
7 /* EDOM, FAU Erlangen-Nürnberg */
8 /* 2014-2015 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-2015 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 
40 /* #define SCIP_DEBUG*/
41 /* #define SCIP_MORE_DEBUG *//* shows all cuts added */
42 
43 #include "cons_sdp.h"
44 
45 #include <assert.h> /*lint !e451*/
46 #include <string.h> /* for NULL, strcmp */
47 
48 #include "sdpi/lapack.h"
49 
50 #include "scipsdp/SdpVarmapper.h"
51 #include "scipsdp/SdpVarfixer.h"
52 
53 #include "scip/cons_linear.h" /* for SCIPcreateConsLinear */
54 #include "scip/scip.h" /* for SCIPallocBufferArray, etc */
55 
56 #define CONSHDLR_NAME "SDP"
57 #define CONSHDLR_DESC "SDP constraints of the form \\sum_{j} A_j y_j - A_0 psd"
58 #define CONSHDLR_SEPAPRIORITY +1000000
59 #define CONSHDLR_ENFOPRIORITY -2000000
60 #define CONSHDLR_CHECKPRIORITY -2000000
61 #define CONSHDLR_SEPAFREQ 1
62 #define CONSHDLR_EAGERFREQ 100
64 #define CONSHDLR_MAXPREROUNDS -1
65 #define CONSHDLR_DELAYSEPA FALSE
66 #define CONSHDLR_NEEDSCONS TRUE
68 #define CONSHDLR_PRESOLTIMING SCIP_PRESOLTIMING_FAST
69 
70 
72 struct SCIP_ConsData
73 {
74  int nvars;
75  int nnonz;
76  int blocksize;
77  int* nvarnonz;
78  int** col;
79  int** row;
80  SCIP_Real** val;
81  SCIP_VAR** vars;
82  int constnnonz;
83  int* constcol;
84  int* constrow;
85  SCIP_Real* constval;
86  SCIP_Real maxrhsentry;
87 };
88 
90 struct SCIP_ConshdlrData
91 {
92  int neigveccuts;
93  int ndiaggezerocuts;
94  int ndiagdomcuts;
95  int n1x1blocks;
96 };
97 
98 #ifndef NDEBUG
99 
102 static
103 int compLowerTriangPos(
104  int i,
105  int j
106  )
107 {
108  assert( j >= 0 );
109  assert( i >= j );
110 
111  return i*(i+1)/2 + j;
112 }
113 #else
114 #define compLowerTriangPos(i, j) (i*(i+1)/2 + j)
115 #endif
116 
118 static
119 SCIP_RETCODE expandSymMatrix(
120  int size,
121  SCIP_Real* symMat,
122  SCIP_Real* fullMat
123  )
124 {
125  int i;
126  int j;
127  int ind = 0;
128 
129  assert( size >= 0 );
130  assert( symMat != NULL );
131  assert( fullMat != NULL );
132 
133  /* traverse the lower triangular part in the order of the indices and copy the values to both lower and upper triangular part */
134  for (i = 0; i < size; i++)
135  {
136  for (j = 0; j <= i; j++)
137  {
138  assert( ind == compLowerTriangPos(i,j) );
139  fullMat[i*size + j] = symMat[ind]; /*lint !e679*/
140  fullMat[j*size + i] = symMat[ind]; /*lint !e679*/
141  ind++;
142  }
143  }
144 
145  return SCIP_OKAY;
146 }
147 
151 static
152 SCIP_RETCODE computeSdpMatrix(
153  SCIP* scip,
154  SCIP_CONS* cons,
155  SCIP_SOL* y,
156  SCIP_Real* matrix
157  )
158 {
159  SCIP_CONSDATA* consdata;
160  int i;
161  int ind;
162  int nvars;
163  int blocksize;
164  SCIP_Real yval;
165 
166  assert( cons != NULL );
167  assert( matrix != NULL );
168 
169  consdata = SCIPconsGetData(cons);
170  nvars = consdata->nvars;
171  blocksize = consdata->blocksize;
172 
173  /* initialize the matrix with 0 */
174  for (i = 0; i < (blocksize * (blocksize + 1))/2; i++)
175  matrix[i] = 0.0;
176 
177  /* add the non-constant-part */
178  for (i = 0; i < nvars; i++)
179  {
180  yval = SCIPgetSolVal(scip, y, consdata->vars[i]);
181  for (ind = 0; ind < consdata->nvarnonz[i]; ind++)
182  matrix[compLowerTriangPos(consdata->row[i][ind], consdata->col[i][ind])] += yval * consdata->val[i][ind];
183  }
184 
185  /* substract the constant part */
186  for (ind = 0; ind < consdata->constnnonz; ind++)
187  matrix[compLowerTriangPos(consdata->constrow[ind], consdata->constcol[ind])] -= consdata->constval[ind];
188 
189  return SCIP_OKAY;
190 }
191 
193 static
194 SCIP_RETCODE multiplyConstraintMatrix(
195  SCIP_CONS* cons,
196  int j,
197  SCIP_Real* v,
198  SCIP_Real* vAv
199  )
200 {
201  SCIP_CONSDATA* consdata;
202  int i;
203 
204  assert( cons != NULL );
205  assert( j >= 0 );
206  assert( vAv != NULL );
207 
208  consdata = SCIPconsGetData(cons);
209  assert( consdata != NULL );
210 
211  assert( j < consdata->nvars );
212 
213  /* initialize the product with 0 */
214  *vAv = 0.0;
215 
216  for (i = 0; i < consdata->nvarnonz[j]; i++)
217  {
218  if (consdata->col[j][i] == consdata->row[j][i])
219  *vAv += v[consdata->col[j][i]] * consdata->val[j][i] * v[consdata->row[j][i]];
220  else
221  {
222  /* Multiply by 2, because the matrix is symmetric and there is one identical contribution each from lower and upper triangular part. */
223  *vAv += 2.0 * v[consdata->col[j][i]] * consdata->val[j][i] * v[consdata->row[j][i]];
224  }
225  }
226 
227  return SCIP_OKAY;
228 }
229 
233 static
234 SCIP_RETCODE setMaxRhsEntry(
235  SCIP_CONS* cons
236  )
237 {
238  SCIP_CONSDATA* consdata;
239  SCIP_Real max;
240  int i;
241 
242  consdata = SCIPconsGetData(cons);
243  assert( consdata != NULL );
244 
245  /* initialize max with zero (this is used if there is no constant-matrix) */
246  max = 0.0;
247 
248  /* iterate over the entries of the constant matrix, updating max if a higher absolute value is found */
249  for (i = 0; i < consdata->constnnonz; i++)
250  {
251  if ( REALABS(consdata->constval[i]) > max )
252  max = REALABS(consdata->constval[i]);
253  }
254 
255  consdata->maxrhsentry = max;
256 
257  return SCIP_OKAY;
258 }
259 
260 
266 static
267 SCIP_RETCODE cutUsingEigenvector(
268  SCIP* scip,
269  SCIP_CONS* cons,
270  SCIP_SOL* sol,
271  SCIP_Real* coeff,
272  SCIP_Real* lhs
273  )
274 {
275  SCIP_CONSDATA* consdata;
276  SCIP_Real* eigenvalues = NULL;
277  SCIP_Real* matrix = NULL;
278  SCIP_Real* fullmatrix = NULL;
279  SCIP_Real* fullconstmatrix = NULL;
280  SCIP_Real* eigenvector = NULL;
281  SCIP_Real* output_vector = NULL;
282  int blocksize;
283  int j;
284 
285  assert( cons != NULL );
286  assert( lhs != NULL );
287 
288  *lhs = 0.0;
289 
290  consdata = SCIPconsGetData(cons);
291  assert( consdata != NULL );
292 
293  blocksize = consdata->blocksize;
294 
295  SCIP_CALL( SCIPallocBufferArray(scip, &eigenvalues, blocksize) );
296  SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize+1))/2 ) );
297  SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, blocksize * blocksize ) );
298  SCIP_CALL( SCIPallocBufferArray(scip, &fullconstmatrix, blocksize * blocksize) );
299  SCIP_CALL( SCIPallocBufferArray(scip, &eigenvector, blocksize) );
300  SCIP_CALL( SCIPallocBufferArray(scip, &output_vector, blocksize) );
301 
302  /* compute the matrix \f$ \sum_j A_j y_j - A_0 \f$ */
303  SCIP_CALL( computeSdpMatrix(scip, cons, sol, matrix) );
304 
305  /* expand it because LAPACK wants the full matrix instead of the lower triangular part */
306  SCIP_CALL( expandSymMatrix(blocksize, matrix, fullmatrix) );
307 
308  SCIP_CALL( SCIPlapackComputeIthEigenvalue(SCIPblkmem(scip), TRUE, blocksize, fullmatrix, 1, eigenvalues, eigenvector) );
309 
310  /* get full constant matrix */
311  SCIP_CALL( SCIPconsSdpGetFullConstMatrix(scip, cons, fullconstmatrix) );
312 
313  /* multiply eigenvector with constant matrix to get lhs (after multiplying again with eigenvector from the left) */
314  SCIP_CALL( SCIPlapackMatrixVectorMult(blocksize, blocksize, fullconstmatrix, eigenvector, output_vector) );
315 
316  for (j = 0; j < blocksize; ++j)
317  *lhs += eigenvector[j] * output_vector[j];
318 
319  /* compute \f$ v^T A_j v \f$ for eigenvector v and each matrix \f$ A_j \f$ to get the coefficients of the LP cut */
320  for (j = 0; j < consdata->nvars; ++j)
321  {
322  SCIP_CALL( multiplyConstraintMatrix(cons, j, eigenvector, &coeff[j]) );
323  }
324 
325  SCIPfreeBufferArray(scip, &output_vector);
326  SCIPfreeBufferArray(scip, &eigenvector);
327  SCIPfreeBufferArray(scip, &fullconstmatrix);
328  SCIPfreeBufferArray(scip, &fullmatrix);
329  SCIPfreeBufferArray(scip, &matrix);
330  SCIPfreeBufferArray(scip, &eigenvalues);
331 
332  return SCIP_OKAY;
333 }
334 
336 SCIP_RETCODE SCIPconsSdpCheckSdpCons(
337  SCIP* scip,
338  SCIP_CONS* cons,
339  SCIP_SOL* sol,
340  SCIP_Bool checkintegrality,
341  SCIP_Bool checklprows,
342  SCIP_Bool printreason,
343  SCIP_RESULT* result
344  )
345 { /*lint --e{715}*/
346  SCIP_CONSDATA* consdata;
347  int blocksize;
348  double check_value;
349  SCIP_Real eigenvalue;
350  SCIP_Real* matrix = NULL;
351  SCIP_Real* fullmatrix = NULL;
352 
353  assert( scip != NULL );
354  assert( cons != NULL );
355  assert( result != NULL );
356  assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(cons)), "SDP") == 0);
357 
358  consdata = SCIPconsGetData(cons);
359  assert( consdata != NULL );
360  blocksize = consdata->blocksize;
361 
362  SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize+1)) / 2) );
363  SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, blocksize * blocksize) );
364  SCIP_CALL( computeSdpMatrix(scip, cons, sol, matrix) );
365  SCIP_CALL( expandSymMatrix(blocksize, matrix, fullmatrix) );
366 
367  SCIP_CALL( SCIPlapackComputeIthEigenvalue(SCIPblkmem(scip), FALSE, blocksize, fullmatrix, 1, &eigenvalue, NULL) );
368 
369  /* This enables checking the second DIMACS Error Norm: err=max{0, -lambda_min(x)/(1+maximumentry of rhs} */
370 #ifdef DIMACS
371  check_value = (-eigenvalue) / (1.0 + consdata->maxrhsentry);
372 #else
373  check_value = -eigenvalue;
374 #endif
375 
376  if ( SCIPisFeasLE(scip, check_value, 0.0) )
377  *result = SCIP_FEASIBLE;
378  else
379  {
380  *result = SCIP_INFEASIBLE;
381 #ifdef SCIP_DEBUG
382  if ( printreason )
383  {
384  SCIP_CALL( SCIPprintCons(scip, cons, NULL) );
385  SCIP_CALL( SCIPprintSol(scip, sol, NULL, FALSE) );
386  SCIPinfoMessage(scip, NULL, "non psd matrix (eigenvalue %f, dimacs error norm = %f).\n", eigenvalue, check_value);
387  }
388 #endif
389  }
390 
391  SCIPfreeBufferArray(scip, &fullmatrix);
392  SCIPfreeBufferArray(scip, &matrix);
393 
394  return SCIP_OKAY;
395 }
396 
398 static
399 SCIP_RETCODE separateSol(
400  SCIP* scip,
401  SCIP_CONSHDLR* conshdlr,
402  SCIP_CONS* cons,
403  SCIP_SOL* sol,
404  SCIP_RESULT* result
405  )
406 {
407  SCIP_CONSDATA* consdata;
408  SCIP_CONSHDLRDATA* conshdlrdata;
409  SCIP_Real lhs = 0.0;
410  SCIP_Real* coeff = NULL;
411  int nvars;
412  char* cutname;
413 #ifndef NDEBUG
414  int snprintfreturn; /* this is used to assert that the SCIP string concatenation works */
415 #endif
416  int j;
417  SCIP_COL** cols;
418  SCIP_Real* vals;
419  int len;
420  SCIP_ROW* row;
421 
422  assert( cons != NULL );
423 
424  consdata = SCIPconsGetData(cons);
425  assert( consdata != NULL );
426 
427  nvars = consdata->nvars;
428  SCIP_CALL( SCIPallocBufferArray(scip, &coeff, nvars ) );
429 
430  SCIP_CALL( cutUsingEigenvector(scip, cons, sol, coeff, &lhs) );
431 
432  SCIP_CALL( SCIPallocBufferArray(scip, &cols, nvars ) );
433  SCIP_CALL( SCIPallocBufferArray(scip, &vals, nvars ) );
434 
435  len = 0;
436  for (j = 0; j < nvars; ++j)
437  {
438  if ( SCIPisZero(scip, coeff[j]) )
439  continue;
440 
441  cols[len] = SCIPvarGetCol(SCIPvarGetProbvar(consdata->vars[j]));
442  vals[len] = coeff[j];
443  ++len;
444  }
445 
446  conshdlrdata = SCIPconshdlrGetData(conshdlr);
447  assert( conshdlrdata != NULL );
448 
449  SCIP_CALL( SCIPallocBufferArray(scip, &cutname, 255) );
450 #ifndef NDEBUG
451  snprintfreturn = SCIPsnprintf(cutname, 255, "sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
452  assert( snprintfreturn < 256 ); /* it returns the number of positions that would have been needed, if that is more than 255, it failed */
453 #else
454  SCIPsnprintf(cutname, 255, "sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
455 #endif
456  SCIP_CALL( SCIPcreateRowCons(scip, &row, conshdlr, cutname , len, cols, vals, lhs, SCIPinfinity(scip), FALSE, FALSE, TRUE) );
457 
458  if ( SCIPisCutEfficacious(scip, sol, row) )
459  {
460 #ifdef SCIP_MORE_DEBUG
461  SCIPinfoMessage(scip, NULL, "Added cut %s: ", cutname);
462  SCIPinfoMessage(scip, NULL, "%f <= ", lhs);
463  for (j = 0; j < len; j++)
464  SCIPinfoMessage(scip, NULL, "+ (%f)*%s", vals[j], SCIPvarGetName(SCIPcolGetVar(cols[j])));
465  SCIPinfoMessage(scip, NULL, "\n");
466 #endif
467 
468  SCIP_Bool infeasible;
469  SCIP_CALL( SCIPaddCut(scip, sol, row, FALSE, &infeasible) );
470  if ( infeasible )
471  *result = SCIP_CUTOFF;
472  else
473  *result = SCIP_SEPARATED;
474  SCIP_CALL( SCIPresetConsAge(scip, cons) );
475  }
476 
477  SCIP_CALL( SCIPreleaseRow(scip, &row) );
478  SCIPfreeBufferArray(scip, &cutname);
479 
480  SCIPfreeBufferArray(scip, &vals);
481  SCIPfreeBufferArray(scip, &cols);
482  SCIPfreeBufferArray(scip, &coeff);
483 
484  return SCIP_OKAY;
485 }
486 
490 static
491 SCIP_RETCODE diagGEzero(
492  SCIP* scip,
493  SCIP_CONS** conss,
494  int nconss,
495  int* naddconss
496  )
497 {
498  int blocksize;
499  SCIP_CONSDATA* consdata;
500  int nvars;
501  int i;
502  int j;
503  int k;
504  int c;
505  SCIP_Real rhs;
506  char* cutname;
507 #ifndef NDEBUG
508  int snprintfreturn; /* used to check if sdnprintf worked */
509 #endif
510  SCIP_Real* matrix;
511  SCIP_Real* cons_array;
512  SCIP_Real* lhs_array;
513 
514  for (c = 0; c < nconss; ++c)
515  {
516  SCIP_CONSHDLR* conshdlr;
517 #ifndef NDEBUG
518  const char* conshdlrName;
519 #endif
520 
521  conshdlr = SCIPconsGetHdlr(conss[c]);
522  assert( conshdlr != NULL );
523 #ifndef NDEBUG
524  conshdlrName = SCIPconshdlrGetName(conshdlr);
525  assert( strcmp(conshdlrName, "SDP") == 0);
526 #endif
527 
528  consdata = SCIPconsGetData(conss[c]);
529  assert( consdata != NULL );
530 
531  blocksize = consdata->blocksize;
532  nvars = consdata->nvars;
533  rhs = SCIPinfinity(scip);
534 
535  SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize + 1)) / 2) );
536  SCIP_CALL( SCIPconsSdpGetLowerTriangConstMatrix(scip, conss[c], matrix) );
537 
538  SCIP_CALL( SCIPallocBufferArray(scip, &cons_array, blocksize * nvars) );
539  SCIP_CALL( SCIPallocBufferArray(scip, &lhs_array, blocksize) );
540 
541  /* the lhs is the (k,k)-entry of the constant matrix */
542  for (k = 0; k < blocksize; ++k)
543  lhs_array[k] = matrix[compLowerTriangPos(k,k)];
544 
545  /* get the (k,k)-entry of every matrix A_j */
546  for (j = 0; j < nvars; ++j)
547  {
548  for (k = 0; k < blocksize; ++k)
549  {
550  /* initialize these with 0 */
551  cons_array[k * nvars + j] = 0.0; /*lint !e679*/
552  }
553 
554  /* go through the nonzeroes of A_j and look for diagonal entries */
555  for (i = 0; i < consdata->nvarnonz[j]; i++)
556  {
557  if ( consdata->col[j][i] == consdata->row[j][i] )
558  cons_array[consdata->col[j][i] * nvars + j] = consdata->val[j][i]; /*lint !e679*/
559  }
560  }
561 
562  SCIP_CALL( SCIPallocBufferArray(scip, &cutname, 255));
563 
564  /* add the LP-cuts to SCIP */
565  for (k = 0; k < blocksize; ++k)
566  {
567  SCIP_CONS* cons;
568  SCIP_CONSHDLRDATA* conshdlrdata;
569 
570  conshdlrdata = SCIPconshdlrGetData(conshdlr);
571 #ifndef NDEBUG
572  snprintfreturn = SCIPsnprintf(cutname, 255, "diag_ge_zero_%d", ++(conshdlrdata->ndiaggezerocuts));
573  assert( snprintfreturn < 256 ); /* this is the number of positions needed, we gave 255 */
574 #else
575  SCIPsnprintf(cutname, 255, "diag_ge_zero_%d", ++(conshdlrdata->ndiaggezerocuts));
576 #endif
577 
578 #ifdef SCIP_MORE_DEBUG
579  SCIPinfoMessage(scip, NULL, "Added lp-constraint %s: ", cutname);
580  SCIPinfoMessage(scip, NULL, "%f <= ", lhs_array[k]);
581  for (i = 0; i < consdata->nvars; i++)
582  {
583  if ( ! SCIPisZero(scip, cons_array[k * consdata->nvars + i]) )
584  SCIPinfoMessage(scip, NULL, "+ (%f)*%s", cons_array[k * consdata->nvars + i], SCIPvarGetName(consdata->vars[i]));
585  }
586  SCIPinfoMessage(scip, NULL, "\n");
587 #endif
588 
589  SCIP_CALL( SCIPcreateConsLinear(scip, &cons, cutname, consdata->nvars, consdata->vars, cons_array + k * consdata->nvars, lhs_array[k], rhs,
590  TRUE, TRUE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) ); /*lint !e679*/
591 
592  SCIP_CALL( SCIPaddCons(scip, cons) );
593  SCIP_CALL( SCIPreleaseCons(scip, &cons) );
594  ++(*naddconss);
595  }
596 
597  SCIPfreeBufferArray(scip, &cutname);
598  SCIPfreeBufferArray(scip, &lhs_array);
599  SCIPfreeBufferArray(scip, &cons_array);
600  SCIPfreeBufferArray(scip, &matrix);
601  }
602 
603  return SCIP_OKAY;
604 }
605 
606 #if 0
607 
615 static
616 SCIP_RETCODE diagDominant(
617  SCIP* scip,
618  SCIP_CONS** conss,
619  int nconss,
620  int* naddconss
621  )
622 {
623  SCIP_Bool* nonzerorows; /* entry i will be 1 if there is an entry (A_0)_ij \f$ for some \f$ j \neq i */
624  int blocksize;
625  int i;
626  int j;
627  int nvars;
628  int var;
629  SCIP_CONS* cons;
630  SCIP_CONSHDLRDATA* conshdlrdata;
631  char* cutname;
632 #ifndef NDEBUG
633  int snprintfreturn;
634 #endif
635  int** diagvars;
636  int* ndiagvars;
637 
638  assert( scip != NULL );
639  assert( conss != NULL );
640  assert( nconss >= 0 );
641  assert( naddconss != NULL );
642 
643  for (i = 0; i < nconss; ++i)
644  {
645  assert( conss[i] != NULL );
646  assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[i])), "SDP") == 0 );
647 
648  SCIP_CONSDATA* consdata = SCIPconsGetData(conss[i]);
649  assert( consdata != NULL );
650 
651  blocksize = consdata->blocksize;
652  nvars = consdata->nvars;
653  SCIP_CALL( SCIPallocBufferArray(scip, &nonzerorows, blocksize) );
654 
655  /* initialize nonzerorows with FALSE */
656  for (j = 0; j < blocksize; j++)
657  nonzerorows[j] = FALSE;
658 
659  /* iterate over all nonzeroes of the constant matrix and set all corresponding rows/cols to true */
660  for (j = 0; j < consdata->constnnonz; j++)
661  {
662  if ( consdata->constcol[j] != consdata->constrow[j] )
663  {
664  assert( ! SCIPisEQ(scip, consdata->constval[j], 0.0) );
665  nonzerorows[consdata->constcol[j]] = TRUE;
666  nonzerorows[consdata->constrow[j]] = TRUE;
667  }
668  }
669 
670  /* diagvars[i] is an array with all variables with a diagonal entry (i,i) in the corresponding matrix, if nonzerorows[i] is true or NULL otherwise
671  * the outer array goes over all rows to ease the access, but only for those that are really needed memory will be allocated */
672  SCIP_CALL( SCIPallocBufferArray(scip, &diagvars, blocksize) );
673  SCIP_CALL( SCIPallocBufferArray(scip, &ndiagvars, blocksize) );
674  for (j = 0; j < blocksize; ++j)
675  {
676  ndiagvars[j] = 0;
677  if ( nonzerorows[j] )
678  {
679  SCIP_CALL( SCIPallocBufferArray(scip, &(diagvars[j]), nvars) );
680  }
681  }
682 
683  /* find all variables with corresponding diagonal entries for a row with nonzero non-diagonal constant entry */
684  for (var = 0; var < nvars; var++)
685  {
686  for (j = 0; j < consdata->nvarnonz[var]; j++)
687  {
688  if ( consdata->col[var][j] == consdata->row[var][j] && nonzerorows[consdata->row[var][j]])
689  {
690  assert( ! SCIPisEQ(scip, consdata->val[var][j], 0.0) );
691  diagvars[consdata->col[var][j]][ndiagvars[consdata->col[var][j]]] = var;
692  ndiagvars[consdata->col[var][j]]++;
693  }
694  }
695  }
696 
697  SCIP_VAR** vars;
698  SCIP_Real* vals;
699  SCIP_CALL( SCIPallocBufferArray(scip, &cutname, 255));
700 
701  for (j = 0; j < blocksize; ++j)
702  {
703  if ( nonzerorows[j] )
704  {
705  SCIP_CALL( SCIPallocBufferArray(scip, &vals, ndiagvars[j]) );
706  SCIP_CALL( SCIPallocBufferArray(scip, &vars, ndiagvars[j]) );
707 
708  /* get the corresponding SCIP variables and set all coefficients to 1 */
709  for (var = 0; var < ndiagvars[j]; ++var)
710  {
711  vars[var] = consdata->vars[diagvars[j][var]];
712  vals[var] = 1.0;
713  }
714 
715  conshdlrdata = SCIPconshdlrGetData(SCIPconsGetHdlr(conss[i]));
716 #ifndef NDEBUG
717  snprintfreturn = SCIPsnprintf(cutname, 255, "diag_dom_%d", ++(conshdlrdata->ndiagdomcuts));
718  assert( snprintfreturn < 256 ); /* the return is the number of spots needed, we gave 255 */
719 #else
720  SCIPsnprintf(cutname, 255, "diag_dom_%d", ++(conshdlrdata->ndiagdomcuts));
721 #endif
722 
723 #ifdef SCIP_MORE_DEBUG
724  SCIPinfoMessage(scip, NULL, "Added lp-constraint %s: ", cutname);
725  SCIPinfoMessage(scip, NULL, "1 <= ");
726  for (i = 0; i < ndiagvars[j]; i++)
727  SCIPinfoMessage(scip, NULL, "+ (%f)*%s", vals[i], SCIPvarGetName(vars[i]));
728  SCIPinfoMessage(scip, NULL, "\n");
729 #endif
730 
731  /* add the linear constraint sum_j 1.0 * diagvars[j] >= 1.0 */
732  SCIP_CALL(SCIPcreateConsLinear(scip, &cons, cutname , ndiagvars[j], vars, vals, 1.0, SCIPinfinity(scip), TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE));
733  SCIP_CALL(SCIPaddCons(scip, cons));
734  SCIP_CALL(SCIPreleaseCons(scip, &cons));
735  (*naddconss)++;
736 
737  SCIPfreeBufferArray(scip, &vars);
738  SCIPfreeBufferArray(scip, &vals);
739  SCIPfreeBufferArray(scip, &diagvars[j]);
740  }
741  }
742 
743  SCIPfreeBufferArray(scip, &cutname);
744  SCIPfreeBufferArray(scip, &ndiagvars);
745  SCIPfreeBufferArray(scip, &diagvars);
746  SCIPfreeBufferArray(scip, &nonzerorows);
747  }
748 
749  return SCIP_OKAY;
750 }
751 #endif
752 
754 static
755 SCIP_RETCODE move_1x1_blocks_to_lp(
756  SCIP* scip,
757  SCIP_CONS** conss,
758  int nconss,
759  int* naddconss,
760  int* ndelconss,
761  SCIP_RESULT* result
762  )
763 {
764  SCIP_CONSHDLR* hdlr;
765  SCIP_CONS* cons;
766  SCIP_CONSHDLRDATA* conshdlrdata;
767  int nnonz;
768  SCIP_VAR** vars;
769  SCIP_Real* coeffs;
770  int nvars;
771  int i;
772  int j;
773  SCIP_Real rhs;
774  int count;
775  int var;
776  char* cutname;
777  SCIP_CONSDATA* consdata;
778 #ifndef NDEBUG
779  int snprintfreturn; /* used to assert the return code of snprintf */
780  const char* hdlrName;
781 #endif
782 
783  *result = SCIP_SUCCESS;
784 
785  for (i = 0; i < nconss; ++i)
786  {
787  hdlr = SCIPconsGetHdlr(conss[i]);
788  assert(hdlr != NULL);
789 
790 #ifndef NDEBUG
791  hdlrName = SCIPconshdlrGetName(hdlr);
792  assert( strcmp(hdlrName, "SDP") == 0);
793 #endif
794 
795  consdata = SCIPconsGetData(conss[i]);
796  assert( consdata != NULL );
797 
798  /* if there is a 1x1 SDP-Block */
799  if ( consdata->blocksize == 1 )
800  {
801  nvars = consdata->nvars;
802  nnonz = consdata->nnonz;
803  SCIP_CALL( SCIPallocBufferArray(scip, &vars, nvars) );
804  SCIP_CALL( SCIPallocBufferArray(scip, &coeffs, nnonz) );
805 
806  /* get all lhs-entries */
807  count = 0;
808 
809  for (var = 0; var < nvars; var++)
810  {
811  assert( consdata->nvarnonz[var] <= 1 ); /* in a 1x1 block there may be at most one entry per variable */
812 
813  for (j = 0; j < consdata->nvarnonz[var]; j++)
814  {
815  assert( consdata->col[var][j] == 0 && consdata->row[var][j] == 0 ); /* if the block is size one, all entries should have row and col equal to 0 */
816  coeffs[count] = consdata->val[var][j];
817  vars[count] = consdata->vars[var];
818  count++;
819  }
820  }
821 
822  /* get rhs */
823  assert( consdata->constnnonz <= 1 ); /* the 1x1 constant matrix may only have one entry */
824 
825  rhs = (consdata->constnnonz == 1) ? consdata->constval[0] : 0.0; /* if this one entry is not 0, than this is the rhs, otherwise it's 0 */
826 
827  /* if there is more than one left-hand-side-entry, add a linear constraint, otherwise update the variable bound */
828  if ( count > 1 )
829  {
830  /* add new linear cons */
831  conshdlrdata = SCIPconshdlrGetData(hdlr);
832  SCIP_CALL( SCIPallocBufferArray(scip, &cutname, 255) );
833 #ifndef NDEBUG
834  snprintfreturn = SCIPsnprintf(cutname, 255, "1x1block_%d", ++(conshdlrdata->n1x1blocks));
835  assert( snprintfreturn < 256 ); /* the return is the number of spots needed, we gave 255 */
836 #else
837  SCIPsnprintf(cutname, 255, "1x1block_%d", ++(conshdlrdata->n1x1blocks));
838 #endif
839 
840 #ifdef SCIP_MORE_DEBUG
841  SCIPinfoMessage(scip, NULL, "Added lp-constraint %s: ", cutname);
842  for (i = 0; i < consdata->nvars; i++)
843  SCIPinfoMessage(scip, NULL, "+ (%f)*%s", coeffs[i], SCIPvarGetName(vars[i]));
844  SCIPinfoMessage(scip, NULL, " <= %f \n", rhs);
845 #endif
846 
847  SCIP_CALL( SCIPcreateConsLinear(scip, &cons, cutname, consdata->nvars, vars, coeffs, rhs, SCIPinfinity(scip),
848  TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
849 
850  SCIP_CALL( SCIPaddCons(scip, cons) );
851  SCIP_CALL( SCIPreleaseCons(scip, &cons) );
852 
853  SCIPfreeBufferArray(scip, &cutname);
854 
855  (*naddconss)++;
856  }
857  else
858  {
859  /* we compare this new variable with the current (local) bounds, we don't do an epsilon check here, because 10^(-21)*x >= 10^(-19) for
860  * epsilon = 10^(-20) would be infeasible then, even though it only says x >= 100 */
861  if ( coeffs[0] > 0.0 )
862  {
863  /* this gives a lower bound, if it is bigger than the current one, we need to update it */
864  if ( SCIPisFeasGT(scip, rhs / coeffs[0], SCIPvarGetLbLocal(vars[0])) )
865  {
866  /* check if the changed bound renders the problem infeasible */
867  if( SCIPisFeasGT(scip, rhs / coeffs[0], SCIPvarGetUbLocal(vars[0])) )
868  {
869  SCIPdebugMessage("Problem found infeasible during presolving, 1x1-SDP-constraint %s caused change"
870  "of lower bound for variable %s from %f to %f, which is bigger than upper bound of %f\n",
871  SCIPconsGetName(conss[i]), SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]), rhs / coeffs[0],
872  SCIPvarGetUbLocal(vars[0]));
873 
874  *result = SCIP_CUTOFF;
875  /* delete old 1x1 sdpcone */
876  SCIP_CALL( SCIPdelCons(scip, conss[i]) );
877  (*ndelconss)++;
878 
879  SCIPfreeBufferArray(scip, &coeffs);
880  SCIPfreeBufferArray(scip, &vars);
881 
882  return SCIP_OKAY; /* the node is infeasible, we don't care for the other constraints */
883  }
884 
885  SCIPdebugMessage("Changing lower bound of variable %s from %f to %f because of 1x1-SDP-constraint %s!\n",
886  SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]), rhs / coeffs[0], SCIPconsGetName(conss[i]));
887  SCIP_CALL( SCIPchgVarLb(scip, vars[0], rhs / coeffs[0]) );
888  }
889  else
890  {
891  SCIPdebugMessage("Deleting 1x1-SDP-constraint %s, new lower bound %f for variable %s no improvement over old bound %f!\n",
892  SCIPconsGetName(conss[i]), rhs / coeffs[0], SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]));
893  }
894  }
895  else if ( coeffs[0] < 0.0 )
896  {
897  /* this gives an upper bound, if it is lower than the current one, we need to update it */
898  if (SCIPisFeasLT(scip, rhs / coeffs[0], SCIPvarGetUbLocal(vars[0])))
899  {
900  /* check if the changed bound renders the problem infeasible */
901  if( SCIPisFeasLT(scip, rhs / coeffs[0], SCIPvarGetLbLocal(vars[0])) )
902  {
903  SCIPdebugMessage("Problem found infeasible during presolving, 1x1-SDP-constraint %s caused change"
904  "of upper bound for variable %s from %f to %f, which is less than lower bound of %f\n",
905  SCIPconsGetName(conss[i]), SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]), rhs / coeffs[0],
906  SCIPvarGetLbLocal(vars[0]));
907 
908  *result = SCIP_CUTOFF;
909  /* delete old 1x1 sdpcone */
910  SCIP_CALL( SCIPdelCons(scip, conss[i]) );
911  (*ndelconss)++;
912 
913  SCIPfreeBufferArray(scip, &coeffs);
914  SCIPfreeBufferArray(scip, &vars);
915 
916  return SCIP_OKAY; /* the node is infeasible, we don't care for the other constraints */
917  }
918 
919  SCIPdebugMessage("Changing upper bound of variable %s from %f to %f because of 1x1-SDP-constraint %s!\n",
920  SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]), -rhs / coeffs[0], SCIPconsGetName(conss[i]));
921  SCIP_CALL( SCIPchgVarUb(scip, vars[0], rhs / coeffs[0]) );
922  }
923  else
924  {
925  SCIPdebugMessage("Deleting 1x1-SDP-constraint %s, new upper bound %f for variable %s no improvement over old bound %f!\n",
926  SCIPconsGetName(conss[i]), rhs / coeffs[0], SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]));
927  }
928  }
929  else
930  {
931  SCIPdebugMessage("Detected 1x1 SDP-block without any nonzero coefficients \n");
932  if ( SCIPisFeasGT(scip, rhs, 0.0) )
933  {
934  SCIPdebugMessage("Detected infeasibility in 1x1 SDP-block without any nonzero coefficients but with strictly positive rhs\n");
935  *result = SCIP_CUTOFF;
936  /* delete old 1x1 sdpcone */
937  SCIP_CALL( SCIPdelCons(scip, conss[i]) );
938  (*ndelconss)++;
939 
940  SCIPfreeBufferArray(scip, &coeffs);
941  SCIPfreeBufferArray(scip, &vars);
942 
943  return SCIP_OKAY; /* the node is infeasible, we don't care for the other constraints */
944  }
945  }
946  }
947 
948  /* delete old 1x1 sdpcone */
949  SCIP_CALL(SCIPdelCons(scip, conss[i]));
950  (*ndelconss)++;
951 
952  SCIPfreeBufferArray(scip, &coeffs);
953  SCIPfreeBufferArray(scip, &vars);
954  }
955  }
956  return SCIP_OKAY;
957 }
958 
960 static
961 SCIP_RETCODE multiaggrVar(
962  SCIP* scip,
963  SCIP_CONS* cons,
964  int* v,
965  SCIP_VAR** aggrvars,
966  SCIP_Real* scalars,
967  int naggrvars,
968  SCIP_Real constant,
969  int* savedcol,
970  int* savedrow,
971  SCIP_Real* savedval,
972  int* nfixednonz,
973  int* vararraylength
974  )
975 {
976  int i;
977  SCIP_CONSDATA* consdata;
978  int startind;
979  int aggrind;
980  int aggrtargetlength;
981  int globalnvars;
982  int aggrconsind;
983  SCIP_Real feastol;
984 
985  assert( scip != NULL );
986  assert( cons != NULL );
987  assert( scalars != NULL );
988  assert( naggrvars > 0 );
989  assert( savedcol != NULL );
990  assert( savedrow != NULL );
991  assert( savedval != NULL );
992  assert( nfixednonz != NULL );
993 
994  consdata = SCIPconsGetData(cons);
995  assert( consdata != NULL );
996 
997  SCIP_CALL( SCIPgetRealParam(scip, "numerics/feastol", &feastol) );
998 
999  /* save the current nfixednonz-index, all entries starting from here will need to be added to the variables this is aggregated to */
1000  startind = *nfixednonz;
1001 
1002  if ( SCIPisEQ(scip, constant, 0.0) )
1003  {
1004  /* if there is no constant part, we just save the nonzeros to add them to the variables they are aggregated to, we do this to be able to remove
1005  * the nonzero-arrays for this variable to be able to fill it with a newly inserted variable, as copying all variables, if we created an empty
1006  * gap in the variable arrays, will be more time consuming then copying all variables (as we will usually have more variables than nonzeros per
1007  * variable */
1008  for (i = 0; i < consdata->nvarnonz[*v]; i++)
1009  {
1010  savedcol[*nfixednonz] = consdata->col[*v][i];
1011  savedrow[*nfixednonz] = consdata->row[*v][i];
1012  savedval[*nfixednonz] = consdata->val[*v][i];
1013  (*nfixednonz)++;
1014  }
1015  }
1016  else
1017  {
1018  for (i = 0; i < consdata->nvarnonz[*v]; i++)
1019  {
1020  savedcol[*nfixednonz] = consdata->col[*v][i];
1021  savedrow[*nfixednonz] = consdata->row[*v][i];
1022  savedval[*nfixednonz] = consdata->val[*v][i] * constant;
1023  (*nfixednonz)++;
1024  }
1025  }
1026  assert( *nfixednonz - startind == consdata->nvarnonz[*v] );
1027 
1028  /* sort them by nondecreasing row and then col to make the search for already existing entries easier (this is done here, because it
1029  * only needs to be done once and not for each variable this is multiaggregated to), we add startind to the pointers to only start where we started
1030  * inserting, the number of elements added to the saved arrays for this variable is nfixednonz - startind */
1031  SCIPsdpVarfixerSortRowCol(savedrow + startind, savedcol + startind, savedval + startind, *nfixednonz - startind);
1032 
1033  /* fill the empty spot of the (multi-)aggregated variable with the last variable of this constraint (as they don't have to be sorted) */
1034  SCIP_CALL( SCIPreleaseVar(scip, &consdata->vars[*v]) );
1035  consdata->col[*v] = consdata->col[consdata->nvars - 1];
1036  consdata->row[*v] = consdata->row[consdata->nvars - 1];
1037  consdata->val[*v] = consdata->val[consdata->nvars - 1];
1038  consdata->nvarnonz[*v] = consdata->nvarnonz[consdata->nvars - 1];
1039  consdata->vars[*v] = consdata->vars[consdata->nvars - 1];
1040  (consdata->nvars)--;
1041  (*v)--; /* we need to check again if the variable we just shifted to this position also needs to be (multi-)aggregated */
1042 
1043  /* iterate over all variables this was aggregated to and insert the corresponding nonzeroes */
1044  for (aggrind = 0; aggrind < naggrvars; aggrind++)
1045  {
1046  /* check if the variable already exists in this block */
1047  aggrconsind = -1;
1048  for (i = 0; i < consdata->nvars; i++)
1049  {
1050  if ( consdata->vars[i] == aggrvars[aggrind] )
1051  {
1052  aggrconsind = i;
1053  break;
1054  }
1055  }
1056 
1057  if ( aggrconsind > -1 )
1058  {
1059  /* if the varialbe to aggregate to is already part of this sdp-constraint, just add the nonzeros of the old variable to it */
1060 
1061  /* resize the arrays to the maximally needed length */
1062  aggrtargetlength = consdata->nvarnonz[aggrconsind] + *nfixednonz - startind;
1063  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->row[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1064  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->col[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1065  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->val[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1066 
1067  if ( SCIPisEQ(scip, constant, 0.0) )
1068  {
1069  /* in this case we saved the original values in savedval, we add startind to the pointers to only add those from
1070  * the current variable, the number of entries is the current position minus the position whre we started */
1071  SCIP_CALL( SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), feastol, savedrow + startind, savedcol + startind, savedval + startind,
1072  *nfixednonz - startind, TRUE, scalars[aggrind], consdata->row[aggrconsind], consdata->col[aggrconsind],
1073  consdata->val[aggrconsind], &(consdata->nvarnonz[aggrconsind]), aggrtargetlength) );
1074  }
1075  else
1076  {
1077  /* in this case we saved the original values * constant, so we now have to divide by constant, we add startind to the pointers
1078  * to only add those from the current variable, the number of entries is the current position minus the position whre we started */
1079  SCIP_CALL( SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), feastol, savedrow + startind, savedcol + startind, savedval + startind,
1080  *nfixednonz - startind, TRUE, scalars[aggrind] / constant, consdata->row[aggrconsind], consdata->col[aggrconsind],
1081  consdata->val[aggrconsind], &(consdata->nvarnonz[aggrconsind]), aggrtargetlength) );
1082  }
1083 
1084  /* shrink them again if nonzeros could be combined */
1085  assert( consdata->nvarnonz[aggrconsind] <= aggrtargetlength );
1086  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->row[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1087  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->col[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1088  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->val[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1089  }
1090  else
1091  {
1092  /* the variable has to be added to this constraint */
1093 
1094  SCIPdebugMessage("adding variable %s to SDP constraint %s because of (multi-)aggregation\n", SCIPvarGetName(aggrvars[aggrind]), SCIPconsGetName(cons));
1095 
1096  /* check if we have to enlarge the arrays */
1097  if ( consdata->nvars == *vararraylength )
1098  {
1099  globalnvars = SCIPgetNVars(scip);
1100 
1101  /* we don't want to enlarge this by one for every variable added, so we immediately set it to the maximum possible size */
1102  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col, *vararraylength, globalnvars) );
1103  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row, *vararraylength, globalnvars) );
1104  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val, *vararraylength, globalnvars) );
1105  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->nvarnonz, *vararraylength, globalnvars) );
1106  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->vars, *vararraylength, globalnvars) );
1107  *vararraylength = globalnvars;
1108  }
1109 
1110  /* we insert this variable at the last position, as the ordering doesn't matter */
1111  SCIP_CALL( SCIPcaptureVar(scip, aggrvars[aggrind]) );
1112  consdata->vars[consdata->nvars] = aggrvars[aggrind];
1113 
1114  /* as there were no nonzeros thus far, we can just duplicate the saved arrays to get the nonzeros for the new variable */
1115  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->col[consdata->nvars]), savedcol + startind, *nfixednonz - startind) ); /*lint !e776*/
1116  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->row[consdata->nvars]), savedrow + startind, *nfixednonz - startind) ); /*lint !e776*/
1117 
1118  /* if scalars[aggrind] = constant, we would multiply with 1, if constant = 0, we didn't divide by constant, so in these cases, we can just
1119  * memcopy the array of nonzero-values */
1120  /* TODO: only checking scalar and constant for feas-equality might lead to big differences, if the nonzeros they are multiplied with are big,
1121  * e.g. scalar = 0, constant = 10^(-6), nonzero = 10^(10) leads to new nonzero of 10^4 instead of 0 */
1122  if ( SCIPisEQ(scip, scalars[aggrind], constant) || SCIPisEQ(scip, constant, 0.0) )
1123  {
1124  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->val[consdata->nvars]), savedval + startind, *nfixednonz - startind) ); /*lint !e776*/
1125  consdata->nvarnonz[consdata->nvars] = *nfixednonz - startind; /* as there were no nonzeros thus far, the number of nonzeros equals
1126  * the number of nonzeros of the aggregated variable */
1127  }
1128  else /* we have to multiply all entries by scalar before inserting them */
1129  {
1130  SCIP_Real epsilon;
1131 
1132  SCIP_CALL( SCIPgetRealParam(scip, "numerics/epsilon", &epsilon) );
1133 
1134  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->val[consdata->nvars]), *nfixednonz - startind) );
1135 
1136  consdata->nvarnonz[consdata->nvars] = 0;
1137 
1138  for (i = 0; i < *nfixednonz - startind; i++)
1139  {
1140  /* if both scalar and savedval are small this might become too small */
1141  if ( (scalars[i] / constant) * savedval[startind + i] >= epsilon ) /*lint !e679*/
1142  {
1143  consdata->val[consdata->nvars][consdata->nvarnonz[consdata->nvars]] = (scalars[i] / constant) * savedval[startind + i]; /*lint !e679*/
1144  consdata->nvarnonz[consdata->nvars]++;
1145  }
1146  }
1147  }
1148 
1149  if ( consdata->nvarnonz[consdata->nvars] > 0 ) /* if scalar and all savedvals were to small */
1150  consdata->nvars++;
1151  }
1152  }
1153 
1154  /* if the constant part is zero, we may delete the nonzero-entries from the saved arrays (by resetting the nfixednonz entry to where
1155  * it started, so that these entries will be overwritten */
1156  if ( SCIPisEQ(scip, constant, 0.0) )
1157  *nfixednonz = startind;
1158 
1159  return SCIP_OKAY;
1160 }
1161 
1162 
1164 static
1165 SCIP_RETCODE fixAndAggrVars(
1166  SCIP* scip,
1167  SCIP_CONS** conss,
1168  int nconss,
1169  SCIP_Bool aggregate
1170  )
1171 {
1172  SCIP_CONSDATA* consdata;
1173  int i;
1174  int nfixednonz;
1175  int* savedcol;
1176  int* savedrow;
1177  SCIP_Real* savedval;
1178  int c;
1179  int v;
1180  int arraylength;
1181  SCIP_VAR* var;
1182  SCIP_VAR** aggrvars;
1183  SCIP_Real scalar;
1184  SCIP_Real* scalars;
1185  int naggrvars;
1186  SCIP_Real constant;
1187  int requiredsize;
1188  int globalnvars;
1189  int vararraylength;
1190  SCIP_Bool negated;
1191  SCIP_Real feastol;
1192 
1193 
1194  /* loop over all variables once, add all fixed to savedrow/col/val, for all multiaggregated variables, if constant-scalar =!= 0, add
1195  * constant-scalar * entry to savedrow/col/val and call mergeArrays for all aggrvars for savedrow[startindex of this var] and scalar/constant-scalar,
1196  * if constant-scalar = 0, add 1*entry to savedrow/col/val, call mergeArrays for all aggrvars for savedrow[startindex of this var] and scalar and later
1197  * reduze the saved size of savedrow/col/val by the number of nonzeros of the mutliagrregated variable to not add them to the constant part later */
1198 
1199  assert( scip != NULL );
1200  assert( conss != NULL );
1201  assert( nconss >= 0 );
1202 
1203  SCIPdebugMessage("Calling fixAndAggrVars with aggregate = %u\n", aggregate);
1204 
1205  SCIP_CALL( SCIPgetRealParam(scip, "numerics/feastol", &feastol) );
1206 
1207  for (c = 0; c < nconss; ++c)
1208  {
1209  assert( conss[c] != NULL );
1210  assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])), "SDP") == 0);
1211 
1212  consdata = SCIPconsGetData(conss[c]);
1213  assert( consdata != NULL );
1214 
1215  /* allocate memory to save nonzeros that need to be fixed */
1216  SCIP_CALL( SCIPallocBufferArray(scip, &savedcol, consdata->nnonz) );
1217  SCIP_CALL( SCIPallocBufferArray(scip, &savedrow, consdata->nnonz) );
1218  SCIP_CALL( SCIPallocBufferArray(scip, &savedval, consdata->nnonz) );
1219 
1220  /* initialize this with zero for each block */
1221  nfixednonz = 0;
1222 
1223  vararraylength = consdata->nvars;
1224  globalnvars = SCIPgetNVars(scip);
1225 
1226  for (v = 0; v < consdata->nvars; v++)/*lint --e{850}*/
1227  {
1228  negated = FALSE;
1229  /* if the variable is negated, get the negation var */
1230  if ( SCIPvarIsBinary(consdata->vars[v]) && SCIPvarIsNegated(consdata->vars[v]) )
1231  {
1232  negated = TRUE;
1233  var = SCIPvarGetNegationVar(consdata->vars[v]);
1234  }
1235  else
1236  var = consdata->vars[v];
1237 
1238  /* check if the variable is fixed in SCIP */
1239  if ( SCIPvarGetStatus(var) == SCIP_VARSTATUS_FIXED || SCIPisEQ(scip, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var)))
1240  {
1241  assert( SCIPisEQ(scip, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var)) );
1242 
1243  SCIPdebugMessage("Globally fixing Variable %s to value %f !\n", SCIPvarGetName(var), SCIPvarGetLbGlobal(var));
1244 
1245  if ( ((! negated) && (! SCIPisEQ(scip, SCIPvarGetLbGlobal(var), 0.0))) || (negated && SCIPisEQ(scip, SCIPvarGetLbGlobal(var), 0.0)) )
1246  {
1247  /* the nonzeros are saved to later be inserted into the constant part (this is only done after all nonzeros of fixed variables have been
1248  * assembled, because we need to sort the constant nonzeros and loop over them, which we only want to do once and not once for each fixed
1249  * variable) */
1250  for (i = 0; i < consdata->nvarnonz[v]; i++)
1251  {
1252  savedcol[nfixednonz] = consdata->col[v][i];
1253  savedrow[nfixednonz] = consdata->row[v][i];
1254  /* this is the final value to add, we no longer have to remember from which variable this comes, minus because we have +A_i but -A_0 */
1255  if ( ! negated )
1256  savedval[nfixednonz] = consdata->val[v][i] * SCIPvarGetLbGlobal(var);
1257  else
1258  savedval[nfixednonz] = consdata->val[v][i]; /* if it is the negation of a variable fixed to zero, this variable is fixed to one */
1259 
1260  nfixednonz++;
1261  consdata->nnonz--;
1262  }
1263  }
1264  else
1265  {
1266  /* if the variable is fixed to zero, the nonzeros will just vanish, so we only reduce the number of nonzeros */
1267  consdata->nnonz -= consdata->nvarnonz[v];
1268  }
1269  /* as the variables don't need to be sorted, we just put the last variable into the empty spot and decrease sizes by one (at the end) */
1270  SCIP_CALL( SCIPreleaseVar(scip, &(consdata->vars[v])) );
1271  consdata->col[v] = consdata->col[consdata->nvars - 1];
1272  consdata->row[v] = consdata->row[consdata->nvars - 1];
1273  consdata->val[v] = consdata->val[consdata->nvars - 1];
1274  consdata->nvarnonz[v] = consdata->nvarnonz[consdata->nvars - 1];
1275  consdata->vars[v] = consdata->vars[consdata->nvars - 1];
1276  consdata->nvars--;
1277  v--; /* we need to check again if the variable we just shifted to this position also needs to be fixed */
1278  }
1279  else if ( (SCIPvarGetStatus(var) == SCIP_VARSTATUS_AGGREGATED ||
1280  SCIPvarGetStatus(var) == SCIP_VARSTATUS_MULTAGGR)
1281  && aggregate )
1282  {
1283  SCIP_CALL( SCIPallocBufferArray(scip, &aggrvars, globalnvars) );
1284  SCIP_CALL( SCIPallocBufferArray(scip, &scalars, globalnvars) );
1285 
1286  /* this is how they should be initialized before calling SCIPgetProbvarLinearSum */
1287  if (! negated)
1288  {
1289  aggrvars[0] = consdata->vars[v];
1290  naggrvars = 1;
1291  constant = 0.0;
1292  scalars[0] = 1.0;
1293  }
1294  else
1295  {
1296  /* if this variable is the negation of var, than we look for a representation of 1-var */
1297  aggrvars[0] = consdata->vars[v];
1298  naggrvars = 1;
1299  constant = 1.0;
1300  scalars[0] = -1.0;
1301  }
1302 
1303  /* get the variables this var was aggregated to */
1304  SCIP_CALL( SCIPgetProbvarLinearSum(scip, aggrvars, scalars, &naggrvars, globalnvars, &constant, &requiredsize, TRUE) );
1305  assert( requiredsize <= globalnvars ); /* requiredsize is the number of empty spots in aggrvars needed, globalnvars is the number
1306  * of spots we provided */
1307 
1308  /* Debugmessages for the (multi-)aggregation */
1309 #ifdef SCIP_DEBUG
1310  if ( SCIPvarGetStatus(consdata->vars[v]) == SCIP_VARSTATUS_AGGREGATED )
1311  SCIPdebugMessage("aggregating variable %s to ", SCIPvarGetName(var));
1312  else
1313  SCIPdebugMessage("multiaggregating variable %s to ", SCIPvarGetName(var));
1314  for (i = 0; i < naggrvars; i++)
1315  printf("+ (%f2) * %s ", scalars[i], SCIPvarGetName(aggrvars[i]));
1316  printf("+ (%f2) \n", constant);
1317 #endif
1318 
1319  /* add the nonzeros to the saved-arrays for the constant part, remove the nonzeros for the old variables and add them to the variables this variable
1320  * was (multi-)aggregated to */
1321  SCIP_CALL( multiaggrVar(scip, conss[c], &v, aggrvars, scalars, naggrvars, constant, savedcol, savedrow, savedval, &nfixednonz, &vararraylength) );
1322 
1323  SCIPfreeBufferArray(scip, &aggrvars);
1324  SCIPfreeBufferArray(scip, &scalars);
1325  }
1326  else if ( negated && (SCIPvarGetStatus(var) == SCIP_VARSTATUS_LOOSE || SCIPvarGetStatus(var) == SCIP_VARSTATUS_COLUMN)
1327  && aggregate)
1328  {
1329  /* if var1 is the negation of var2, then this is equivalent to it being aggregated to -var2 + 1 = 1 - var2 */
1330 
1331  SCIPdebugMessage("Changing variable %s to negation of variable %s !\n", SCIPvarGetName(consdata->vars[v]), SCIPvarGetName(var));
1332 
1333  scalar = -1.0;
1334 
1335  SCIP_CALL( multiaggrVar(scip, conss[c], &v, &var, &scalar, 1, 1.0, savedcol, savedrow, savedval, &nfixednonz, &vararraylength) );
1336  }
1337  }
1338 
1339  /* shrink the variable arrays if they were enlarged too much (or more vars were removed than added) */
1340  assert( consdata->nvars <= vararraylength );
1341  if ( consdata->nvars < vararraylength )
1342  {
1343  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col, vararraylength, consdata->nvars) );
1344  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row, vararraylength, consdata->nvars) );
1345  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val, vararraylength, consdata->nvars) );
1346  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->nvarnonz, vararraylength, consdata->nvars) );
1347  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->vars, vararraylength, consdata->nvars) );
1348  }
1349 
1350  /* allocate the maximally needed memory for inserting the fixed variables into the constant part */
1351  arraylength = consdata->constnnonz + nfixednonz;
1352  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constcol), consdata->constnnonz, arraylength) );
1353  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constrow), consdata->constnnonz, arraylength) );
1354  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constval), consdata->constnnonz, arraylength) );
1355 
1356  /* insert the fixed variables into the constant arrays, as we have +A_i but -A_0 we mutliply them by -1 */
1357  SCIP_CALL( SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), feastol, savedrow, savedcol, savedval, nfixednonz, FALSE, -1.0, consdata->constrow,
1358  consdata->constcol, consdata->constval, &(consdata->constnnonz), arraylength) );
1359 
1360  assert( consdata->constnnonz <= arraylength ); /* the allocated memory should always be sufficient */
1361 
1362  /* shrink the arrays if nonzeros could be combined */
1363  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constcol), arraylength, consdata->constnnonz) );
1364  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constrow), arraylength, consdata->constnnonz) );
1365  SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constval), arraylength, consdata->constnnonz) );
1366 
1367  /* free the saved arrays */
1368  SCIPfreeBufferArray(scip, &savedval);
1369  SCIPfreeBufferArray(scip, &savedrow);
1370  SCIPfreeBufferArray(scip, &savedcol);
1371 
1372  /* recompute sdpnnonz */
1373  consdata->nnonz = 0;
1374  for (v = 0; v < consdata->nvars; v++)
1375  consdata->nnonz += consdata->nvarnonz[v];
1376  }
1377 
1378  return SCIP_OKAY;
1379 }
1380 
1382 static
1383 SCIP_DECL_CONSINITPRE(consInitpreSdp)
1384 {/*lint --e{715}*/
1385  SCIP_CONSHDLRDATA* conshdlrdata;
1386 
1387  assert( conshdlr != NULL );
1388 
1389  conshdlrdata = SCIPconshdlrGetData(conshdlr);
1390  assert( conshdlrdata != NULL );
1391 
1392  conshdlrdata->neigveccuts = 0; /* this is used to give the eigenvector-cuts distinguishable names */
1393  conshdlrdata->ndiaggezerocuts = 0; /* this is used to give the diagGEzero-cuts distinguishable names */
1394  conshdlrdata->ndiagdomcuts = 0; /* this is used to give the diagDominant-cuts distinguishable names */
1395  conshdlrdata->n1x1blocks = 0; /* this is used to give the lp constraints resulting from 1x1 sdp-blocks distinguishable names */
1396 
1397  return SCIP_OKAY;
1398 }
1399 
1401 static
1402 SCIP_DECL_CONSLOCK(consLockSdp)
1403 {/*lint --e{715}*/
1404  SCIP_Real* Aj;
1405  SCIP_CONSDATA* consdata;
1406  int blocksize;
1407  int var;
1408  int nvars;
1409  SCIP_Real eigenvalue;
1410 
1411  consdata = SCIPconsGetData(cons);
1412  assert( consdata != NULL );
1413  blocksize = consdata->blocksize;
1414  nvars = consdata->nvars;
1415 
1416  SCIP_CALL( SCIPallocBufferArray(scip, &Aj, blocksize * blocksize) );
1417 
1418  for (var = 0; var < nvars; var++)
1419  {
1420  SCIP_CALL( SCIPconsSdpGetFullAj(scip, cons, var, Aj) );
1421 
1422  /* compute the smallest eigenvalue */
1423  SCIP_CALL( SCIPlapackComputeIthEigenvalue(SCIPblkmem(scip), FALSE, blocksize, Aj, 1, &eigenvalue, NULL) );
1424  if ( SCIPisNegative(scip, eigenvalue) )
1425  {
1426  /* as the lowest eigenvalue is negative, the matrix is not positive semidefinite, so adding more of it can remove positive
1427  * semidefiniteness of the SDP-matrix */
1428  SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlocksneg, nlockspos) );
1429  }
1430 
1431  /* if the smallest eigenvalue is already positive, we don't need to compute the biggest one */
1432  if ( SCIPisPositive(scip, eigenvalue) )
1433  {
1434  /* as an eigenvalue is positive, the matrix is not negative semidefinite, so substracting more of it can remove positive
1435  * semidefiniteness of the SDP-matrix */
1436  SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlockspos, nlocksneg) );
1437  }
1438  else
1439  {
1440  /* compute the biggest eigenvalue */
1441  SCIP_CALL( SCIPlapackComputeIthEigenvalue(SCIPblkmem(scip), FALSE, blocksize, Aj, blocksize, &eigenvalue, NULL) );
1442  if ( SCIPisPositive(scip, eigenvalue) )
1443  {
1444  /* as the biggest eigenvalue is positive, the matrix is not negative semidefinite, so substracting more of it can remove positive
1445  * semidefiniteness of the SDP-matrix */
1446  SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlockspos, nlocksneg) );
1447  }
1448  }
1449  }
1450 
1451  SCIPfreeBufferArray(scip, &Aj);
1452 
1453  return SCIP_OKAY;
1454 }
1455 
1457 static
1458 SCIP_DECL_CONSEXITPRE(consExitpreSdp)
1459 {/*lint --e{715}*/
1460  assert( scip != NULL );
1461  assert( conss != NULL );
1462 
1463  SCIP_CALL( fixAndAggrVars(scip, conss, nconss, TRUE) );
1464 
1465  return SCIP_OKAY;
1466 }
1467 
1469 static
1470 SCIP_DECL_CONSPRESOL(consPresolSdp)
1471 {/*lint --e{715}*/
1472  assert( result != 0 );
1473 
1474  if ( nrounds == 0 )
1475  {
1476  SCIP_CALL( diagGEzero(scip, conss, nconss, naddconss) );
1477  SCIP_CALL( move_1x1_blocks_to_lp(scip, conss, nconss, naddconss, ndelconss, result) );
1478  }
1479 
1480 #if 0
1481  if ( nrounds == 0 )
1482  {
1483  SCIP_CALL( diagDominant(scip, conss, nconss, naddconss) ); /* TODO: could be activated for some problem classes but doesn't work in the general case */
1484  }
1485 #endif
1486 
1487  return SCIP_OKAY;
1488 }
1489 
1491 static
1492 SCIP_DECL_CONSTRANS(consTransSdp)
1493 {/*lint --e{715}*/
1494  SCIP_CONSDATA* sourcedata;
1495  SCIP_CONSDATA* targetdata;
1496  int i;
1497 
1498  sourcedata = SCIPconsGetData(sourcecons);
1499  assert( sourcedata != NULL );
1500 
1501  SCIP_CALL( SCIPallocBlockMemory(scip, &targetdata) );
1502 
1503  /* copy some general data */
1504  targetdata->nvars = sourcedata->nvars;
1505  targetdata->nnonz = sourcedata->nnonz;
1506  targetdata->blocksize = sourcedata->blocksize;
1507 
1508  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->nvarnonz), sourcedata->nvarnonz, sourcedata->nvars) );
1509 
1510  /* copy the non-constant nonzeros */
1511  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->col), sourcedata->nvars) );
1512  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->row), sourcedata->nvars) );
1513  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->val), sourcedata->nvars) );
1514 
1515  for (i = 0; i < sourcedata->nvars; i++)
1516  {
1517  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->col[i]), sourcedata->col[i], sourcedata->nvarnonz[i]) );
1518  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->row[i]), sourcedata->row[i], sourcedata->nvarnonz[i]) );
1519  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->val[i]), sourcedata->val[i], sourcedata->nvarnonz[i]) );
1520  }
1521  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->vars), sourcedata->nvars));
1522 
1523  /* copy & transform the vars array */
1524  for (i = 0; i < sourcedata->nvars; i++)
1525  {
1526  targetdata->vars[i] = SCIPvarGetTransVar(sourcedata->vars[i]);
1527  SCIP_CALL( SCIPcaptureVar(scip, targetdata->vars[i]) );
1528  }
1529 
1530  /* copy the constant nonzeros */
1531  targetdata->constnnonz = sourcedata->constnnonz;
1532 
1533  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constcol), sourcedata->constcol, sourcedata->constnnonz));
1534  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constrow), sourcedata->constrow, sourcedata->constnnonz));
1535  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constval), sourcedata->constval, sourcedata->constnnonz));
1536 
1537  /* copy the maxrhsentry */
1538  targetdata->maxrhsentry = sourcedata->maxrhsentry;
1539 
1540  /* create target constraint */
1541  SCIP_CALL( SCIPcreateCons(scip, targetcons, SCIPconsGetName(sourcecons), conshdlr, targetdata,
1542  SCIPconsIsInitial(sourcecons), SCIPconsIsSeparated(sourcecons), SCIPconsIsEnforced(sourcecons),
1543  SCIPconsIsChecked(sourcecons), SCIPconsIsPropagated(sourcecons), SCIPconsIsLocal(sourcecons),
1544  SCIPconsIsModifiable(sourcecons), SCIPconsIsDynamic(sourcecons), SCIPconsIsRemovable(sourcecons),
1545  SCIPconsIsStickingAtNode(sourcecons)) );
1546 
1547  return SCIP_OKAY;
1548 }
1549 
1551 static
1552 SCIP_DECL_CONSCHECK(consCheckSdp)
1553 {/*lint --e{715}*/
1554  int i;
1555 
1556  assert( scip != NULL );
1557  assert( result != NULL );
1558 
1559  *result = SCIP_FEASIBLE;
1560 
1561  for (i = 0; i < nconss; ++i)
1562  {
1563  SCIP_CALL( SCIPconsSdpCheckSdpCons(scip, conss[i], sol, checkintegrality, checklprows, printreason, result) );
1564  if ( *result == SCIP_INFEASIBLE )
1565  return SCIP_OKAY;
1566  }
1567 
1568  return SCIP_OKAY;
1569 }
1570 
1575 static
1576 SCIP_DECL_CONSENFOPS(consEnfopsSdp)
1577 {/*lint --e{715}*/
1578  int i;
1579 
1580  assert( scip != NULL );
1581  assert( result != NULL );
1582  assert( conss != NULL );
1583 
1584  *result = SCIP_DIDNOTRUN;
1585 
1586  if ( objinfeasible )
1587  {
1588  SCIPdebugMessage("-> pseudo solution is objective infeasible, return.\n");
1589  return SCIP_OKAY;
1590  }
1591 
1592  for (i = 0; i < nconss; ++i)
1593  {
1594  SCIP_CALL( SCIPconsSdpCheckSdpCons(scip, conss[i], NULL, 0, 0, 0, result) );
1595 
1596  if (*result == SCIP_INFEASIBLE)
1597  {
1598  /* if it is infeasible for one SDP constraint, it is infeasible for the whole problem */
1599  SCIPdebugMessage("-> pseudo solution infeasible for SDP-constraint %s, return.\n", SCIPconsGetName(conss[i]));
1600  return SCIP_OKAY;
1601  }
1602  }
1603 
1604  SCIPdebugMessage("-> pseudo solution feasible for all SDP-constraints.\n");
1605 
1606  return SCIP_OKAY;
1607 }
1608 
1609 
1614 static
1615 SCIP_DECL_CONSENFOLP(consEnfolpSdp)
1616 {/*lint --e{715}*/
1617  SCIP_CONSDATA* consdata;
1618  SCIP_CONSHDLRDATA* conshdlrdata;
1619  SCIP_Bool all_feasible = TRUE;
1620  SCIP_Bool separated = FALSE;
1621  char* cutname;
1622 #ifndef NDEBUG
1623  int snprintfreturn; /* used to check the return code of snprintf */
1624 #endif
1625  int i;
1626  int j;
1627  int nvars;
1628  SCIP_ROW* row;
1629  SCIP_Bool infeasible;
1630  SCIP_Real lhs;
1631  SCIP_Real* coeff;
1632  SCIP_Real rhs;
1633 #if 0 /* TODO: see below */
1634  SCIP_VAR** vars;
1635  int count;
1636 #endif
1637 
1638  assert( result != NULL );
1639  *result = SCIP_FEASIBLE;
1640 
1641  for (i = 0; i < nconss; ++i)
1642  {
1643  consdata = SCIPconsGetData(conss[i]);
1644  SCIP_CALL( SCIPconsSdpCheckSdpCons(scip, conss[i], NULL, 0, 0, 0, result) );
1645  if ( *result == SCIP_FEASIBLE )
1646  continue;
1647 
1648  all_feasible = FALSE;
1649 
1650  nvars = consdata->nvars;
1651  lhs = 0.0;
1652  coeff = NULL;
1653 
1654  SCIP_CALL( SCIPallocBufferArray(scip, &coeff, nvars) );
1655  SCIP_CALL( cutUsingEigenvector(scip, conss[i], NULL, coeff, &lhs) );
1656 
1657  rhs = SCIPinfinity(scip);
1658  conshdlrdata = SCIPconshdlrGetData(conshdlr);
1659 
1660  SCIP_CALL( SCIPallocBufferArray(scip, &cutname, 255) );
1661 #ifndef NDEBUG
1662  snprintfreturn = SCIPsnprintf(cutname, 255, "sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
1663  assert( snprintfreturn < 256 ); /* this is the number of spots needed, we gave 255 */
1664 #else
1665  SCIPsnprintf(cutname, 255, "sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
1666 #endif
1667 
1668  SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, conshdlr, cutname , lhs, rhs, FALSE, FALSE, TRUE) );
1669  SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
1670 
1671  for (j = 0; j < nvars; ++j)
1672  {
1673  SCIP_CALL( SCIPaddVarToRow(scip, row, consdata->vars[j], coeff[j]) );
1674  }
1675 
1676  SCIP_CALL( SCIPflushRowExtensions(scip, row) );
1677 
1678 #ifdef SCIP_MORE_DEBUG
1679  SCIPinfoMessage(scip, NULL, "Added cut %s: ", cutname);
1680  SCIPinfoMessage(scip, NULL, "%f <= ", lhs);
1681  for (j = 0; j < nvars; j++)
1682  SCIPinfoMessage(scip, NULL, "+ (%f)*%s", coeff[j], SCIPvarGetName(consdata->vars[j]));
1683  SCIPinfoMessage(scip, NULL, "\n");
1684 #endif
1685 
1686  SCIP_CALL( SCIPaddCut(scip, NULL, row, FALSE, &infeasible) );
1687 
1688  if ( infeasible )
1689  *result = SCIP_CUTOFF;
1690  else
1691  {
1692  SCIP_CALL( SCIPaddPoolCut(scip, row) );
1693 
1694  SCIP_CALL( SCIPresetConsAge(scip, conss[i]) );
1695  *result = SCIP_SEPARATED;
1696  separated = TRUE;
1697  }
1698  SCIP_CALL( SCIPreleaseRow(scip, &row) );
1699  SCIPfreeBufferArray(scip, &cutname);
1700  SCIPfreeBufferArray(scip, &coeff);
1701  }
1702  if ( all_feasible )
1703  return SCIP_OKAY;
1704 
1705  if ( separated )
1706  *result = SCIP_SEPARATED;
1707 #if 0 /* TODO: should this be done here or is it the task of the conshdlr_integer ? if we do it, we should use Feastol and also check for integrality of solution and the counter should obviously be removed */
1708  vars = SCIPgetVars(scip);
1709  count = 0;
1710  for (i = 0; i < SCIPgetNVars(scip); ++i)
1711  {
1712  if ( !SCIPisRelEQ(scip, SCIPvarGetLbLocal(vars[i]), SCIPvarGetUbLocal(vars[i])) && SCIPvarIsIntegral(vars[i]))
1713  {
1714  SCIP_CALL( SCIPaddExternBranchCand(scip, vars[i], 10000.0, SCIP_INVALID) );
1715  count++;
1716  }
1717  }
1718 #endif
1719 
1720  return SCIP_OKAY;
1721 }
1722 
1724 static
1725 SCIP_DECL_CONSSEPASOL(consSepasolSdp)
1726 {/*lint --e{715}*/
1727  int i;
1728 
1729  *result = SCIP_DIDNOTFIND;
1730  for (i = 0; i < nusefulconss; ++i)
1731  {
1732  SCIP_CALL( separateSol(scip, conshdlr, conss[i], sol, result) );
1733  }
1734 
1735  return SCIP_OKAY;
1736 }
1737 
1739 static
1740 SCIP_DECL_CONSSEPALP(consSepalpSdp)
1741 {/*lint --e{715}*/
1742  int i;
1743 
1744  *result = SCIP_DIDNOTFIND;
1745  for (i = 0; i < nusefulconss; ++i)
1746  {
1747  SCIP_CALL( separateSol(scip, conshdlr, conss[i], NULL, result) );
1748  }
1749 
1750  return SCIP_OKAY;
1751 }
1752 
1754 static
1755 SCIP_DECL_CONSDELETE(consDeleteSdp)
1756 {/*lint --e{715}*/
1757  int i;
1758 
1759  assert( cons != NULL );
1760  assert( consdata != NULL );
1761 
1762  SCIPdebugMessage("deleting %s \n", SCIPconsGetName(cons));
1763 
1764  for (i = 0; i < (*consdata)->nvars; i++)
1765  {
1766  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->val[i], (*consdata)->nvarnonz[i]);
1767  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->row[i], (*consdata)->nvarnonz[i]);
1768  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->col[i], (*consdata)->nvarnonz[i]);
1769  }
1770 
1771  /* release all variables */
1772  for (i = 0; i < (*consdata)->nvars; i++)
1773  {
1774  SCIP_CALL( SCIPreleaseVar(scip, &((*consdata)->vars[i])) );
1775  }
1776 
1777  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->vars, (*consdata)->nvars);
1778  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constval, (*consdata)->constnnonz);
1779  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constrow, (*consdata)->constnnonz);
1780  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constcol, (*consdata)->constnnonz);
1781  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->val, (*consdata)->nvars);
1782  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->row, (*consdata)->nvars);
1783  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->col, (*consdata)->nvars);
1784  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->nvarnonz, (*consdata)->nvars);
1785  SCIPfreeBlockMemory(scip, consdata);
1786 
1787  return SCIP_OKAY;
1788 }
1789 
1791 static
1792 SCIP_DECL_CONSFREE(consFreeSdp)
1793 {/*lint --e{715}*/
1794  SCIP_CONSHDLRDATA* conshdlrdata;
1795 
1796  conshdlrdata = SCIPconshdlrGetData(conshdlr);
1797  assert( conshdlrdata != NULL );
1798 
1799  SCIPfreeMemory(scip, &conshdlrdata);
1800  SCIPconshdlrSetData(conshdlr, NULL);
1801 
1802  return SCIP_OKAY;
1803 }
1804 
1806 static
1807 SCIP_DECL_CONSHDLRCOPY(conshdlrCopySdp)
1808 {
1809  assert(scip != NULL);
1810  assert(conshdlr != NULL);
1811  assert(strcmp(SCIPconshdlrGetName(conshdlr), CONSHDLR_NAME) == 0);
1812 
1813  SCIP_CALL( SCIPincludeConshdlrSdp(scip) );
1814 
1815  *valid = TRUE;
1816 
1817  return SCIP_OKAY;
1818 }
1819 
1821 static
1822 SCIP_DECL_CONSCOPY(consCopySdp)
1823 {/*lint --e{715}*/
1824  SCIP_CONSDATA* sourcedata;
1825  SCIP_Bool success;
1826  SCIP_VAR** targetvars;
1827  SCIP_VAR* var;
1828  int i;
1829 
1830  assert( scip != NULL );
1831  assert( sourcescip != NULL );
1832  assert( sourcecons != NULL );
1833  assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(sourcecons)), CONSHDLR_NAME) == 0 );
1834 
1835  SCIPdebugMessage("Copying SDP constraint %s\n", SCIPconsGetName(sourcecons));
1836 
1837  *valid = TRUE;
1838 
1839  /* as we can only map active variables, we have to make sure, that the constraint contains no fixed or (multi-)aggregated vars, after
1840  * exitpresolve (stage 6) this should always be the case, earlier than that we need to call fixAndAggrVars */
1841  if ( SCIPgetStage(sourcescip) <= SCIP_STAGE_EXITPRESOLVE )
1842  {
1843  SCIP_CALL( fixAndAggrVars(sourcescip, &sourcecons, 1, TRUE) );
1844  }
1845 
1846 
1847  sourcedata = SCIPconsGetData(sourcecons);
1848  assert( sourcedata != NULL );
1849 
1850  SCIP_CALL( SCIPallocBufferArray(scip, &targetvars, sourcedata->nvars) );
1851 
1852  /* map all variables in the constraint */
1853  for (i = 0; i < sourcedata->nvars; i++)
1854  {
1855  SCIP_CALL( SCIPgetVarCopy(sourcescip, scip, sourcedata->vars[i], &var, varmap, consmap, global, &success) );
1856  if ( success )
1857  {
1858  targetvars[i] = var;
1859  SCIP_CALL( SCIPcaptureVar(scip, targetvars[i]) );
1860  }
1861  else
1862  *valid = FALSE;
1863  }
1864 
1865  /* create the new constraint */
1866  SCIP_CALL( SCIPcreateConsSdp( scip, cons, name, sourcedata->nvars, sourcedata->nnonz, sourcedata->blocksize, sourcedata->nvarnonz,
1867  sourcedata->col, sourcedata->row, sourcedata->val, targetvars, sourcedata->constnnonz,
1868  sourcedata->constcol, sourcedata->constrow, sourcedata->constval) );
1869 
1870  SCIPfreeBufferArray(scip, &targetvars);
1871 
1872  return SCIP_OKAY;
1873 }
1874 
1876 static
1877 SCIP_DECL_CONSPRINT(consPrintSdp)
1878 {/*lint --e{715}*/
1879  SCIP_CONSDATA* consdata;
1880  SCIP_Real* fullmatrix;
1881  int v;
1882  int i;
1883  int j;
1884 
1885  assert( cons != NULL );
1886 
1887  consdata = SCIPconsGetData(cons);
1888 
1889  SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, consdata->blocksize * consdata->blocksize) );
1890 
1891  /* print the non-constant matrices, for this they first have to be assembled in fullmatrix */
1892  for (v = 0; v < consdata->nvars; v++)
1893  {
1894  /* assemble the matrix */
1895 
1896  /* first initialize it with zero */
1897  for (i = 0; i < consdata->blocksize; i++)
1898  {
1899  for (j = 0; j < consdata->blocksize; j++)
1900  fullmatrix[i * consdata->blocksize + j] = 0.0; /*lint !e679*/
1901  }
1902 
1903  /* then add the nonzeros */
1904  for (i = 0; i < consdata->nvarnonz[v]; i++)
1905  {
1906  fullmatrix[consdata->row[v][i] * consdata->blocksize + consdata->col[v][i]] = consdata->val[v][i]; /* lower triangular entry */ /*lint !e679*/
1907  fullmatrix[consdata->col[v][i] * consdata->blocksize + consdata->row[v][i]] = consdata->val[v][i]; /* upper triangular entry */ /*lint !e679*/
1908  }
1909 
1910  /* print it */
1911  SCIPinfoMessage(scip, file, "+\n");
1912  for (i = 0; i < consdata->blocksize; i++)
1913  {
1914  SCIPinfoMessage(scip, file, "( ");
1915  for (j = 0; j < consdata->blocksize; j++)
1916  SCIPinfoMessage(scip, file, "%f ", fullmatrix[i * consdata->blocksize + j]); /*lint !e679*/
1917  SCIPinfoMessage(scip, file, ")\n");
1918  }
1919  SCIPinfoMessage(scip, file, "* %s\n", SCIPvarGetName(consdata->vars[v]));
1920  }
1921 
1922  /* print the constant-matrix */
1923 
1924  /* assemble the matrix */
1925 
1926  /* first initialize it with zero */
1927  for (i = 0; i < consdata->blocksize; i++)
1928  {
1929  for (j = 0; j < consdata->blocksize; j++)
1930  fullmatrix[i * consdata->blocksize + j] = 0.0; /*lint !e679*/
1931  }
1932 
1933  /* then add the nonzeros */
1934  for (i = 0; i < consdata->constnnonz; i++)
1935  {
1936  fullmatrix[consdata->constrow[i] * consdata->blocksize + consdata->constcol[i]] = consdata->constval[i]; /* lower triangular entry */ /*lint !e679*/
1937  fullmatrix[consdata->constcol[i] * consdata->blocksize + consdata->constrow[i]] = consdata->constval[i]; /* upper triangular entry */ /*lint !e679*/
1938  }
1939 
1940  /* print it */
1941  SCIPinfoMessage(scip, file, "-\n");
1942  for (i = 0; i < consdata->blocksize; i++)
1943  {
1944  SCIPinfoMessage(scip, file, "( ");
1945  for (j = 0; j < consdata->blocksize; j++)
1946  SCIPinfoMessage(scip, file, "%f ", fullmatrix[i * consdata->blocksize + j]); /*lint !e679*/
1947  SCIPinfoMessage(scip, file, ")\n");
1948  }
1949  SCIPinfoMessage(scip, file, ">=0\n");
1950 
1951  SCIPfreeBufferArray(scip, &fullmatrix);
1952 
1953  return SCIP_OKAY;
1954 }
1955 
1957 static
1958 SCIP_DECL_CONSGETVARS(consGetVarsSdp)
1959 {/*lint --e{715}*/
1960  SCIP_CONSDATA* consdata;
1961  int i;
1962  int nvars;
1963 
1964  assert( scip != NULL );
1965  assert( cons != NULL );
1966  assert( vars != NULL );
1967  assert( success != NULL );
1968  assert( varssize >= 0 );
1969 
1970  consdata = SCIPconsGetData(cons);
1971  assert( consdata != NULL );
1972 
1973  nvars = consdata->nvars;
1974 
1975  if ( nvars > varssize )
1976  {
1977  SCIPdebugMessage("consGetVarsIndicator called for array of size %d, needed size %d.\n", varssize, nvars);
1978  *success = FALSE;
1979  return SCIP_OKAY;
1980  }
1981 
1982  for (i = 0; i < nvars; i++)
1983  {
1984  vars[i] = consdata->vars[i];
1985  }
1986 
1987  *success = TRUE;
1988  return SCIP_OKAY;
1989 }
1990 
1992 static
1993 SCIP_DECL_CONSGETNVARS(consGetNVarsSdp)
1994 {/*lint --e{715}*/
1995  SCIP_CONSDATA* consdata;
1996 
1997  assert( scip != NULL );
1998  assert( cons != NULL );
1999  assert( nvars != NULL );
2000  assert( success != NULL );
2001 
2002  consdata = SCIPconsGetData(cons);
2003  assert( consdata != NULL );
2004 
2005  *nvars = consdata->nvars;
2006  *success = TRUE;
2007 
2008  return SCIP_OKAY;
2009 }
2010 
2012 SCIP_RETCODE SCIPincludeConshdlrSdp(
2013  SCIP* scip
2014  )
2015 {
2016  SCIP_CONSHDLR* conshdlr = NULL;
2017  SCIP_CONSHDLRDATA* conshdlrdata = NULL;
2018 
2019  assert( scip != 0 );
2020 
2021  /* allocate memory for the conshdlrdata */
2022  SCIP_CALL( SCIPallocMemory(scip, &conshdlrdata) );
2023 
2024  /* include constraint handler */
2025  SCIP_CALL( SCIPincludeConshdlrBasic(scip, &conshdlr, CONSHDLR_NAME, CONSHDLR_DESC,
2026  CONSHDLR_ENFOPRIORITY, CONSHDLR_CHECKPRIORITY, CONSHDLR_EAGERFREQ, CONSHDLR_NEEDSCONS,
2027  consEnfolpSdp, consEnfopsSdp, consCheckSdp, consLockSdp, conshdlrdata) );
2028 
2029  assert( conshdlr != NULL );
2030 
2031  /* set non-fundamental callbacks via specific setter functions */
2032  SCIP_CALL( SCIPsetConshdlrDelete(scip, conshdlr, consDeleteSdp) );
2033  SCIP_CALL( SCIPsetConshdlrFree(scip, conshdlr, consFreeSdp) );
2034  SCIP_CALL( SCIPsetConshdlrCopy(scip, conshdlr, conshdlrCopySdp, consCopySdp) );
2035  SCIP_CALL( SCIPsetConshdlrInitpre(scip, conshdlr,consInitpreSdp) );
2036  SCIP_CALL( SCIPsetConshdlrExitpre(scip, conshdlr, consExitpreSdp) );
2037  SCIP_CALL( SCIPsetConshdlrPresol(scip, conshdlr, consPresolSdp, CONSHDLR_MAXPREROUNDS, CONSHDLR_PRESOLTIMING) );
2038  SCIP_CALL( SCIPsetConshdlrSepa(scip, conshdlr, consSepalpSdp, consSepasolSdp, CONSHDLR_SEPAFREQ,
2039  CONSHDLR_SEPAPRIORITY, CONSHDLR_DELAYSEPA) );
2040  SCIP_CALL( SCIPsetConshdlrTrans(scip, conshdlr, consTransSdp) );
2041  SCIP_CALL( SCIPsetConshdlrPrint(scip, conshdlr, consPrintSdp) );
2042  SCIP_CALL( SCIPsetConshdlrGetVars(scip, conshdlr, consGetVarsSdp) );
2043  SCIP_CALL( SCIPsetConshdlrGetNVars(scip, conshdlr, consGetNVarsSdp) );
2044 
2045  return SCIP_OKAY;
2046 }
2047 
2054 SCIP_RETCODE SCIPconsSdpGetData(
2055  SCIP* scip,
2056  SCIP_CONS* cons,
2057  int* nvars,
2058  int* nnonz,
2059  int* blocksize,
2060  int* arraylength,
2061  int* nvarnonz,
2063  int** col,
2064  int** row,
2065  SCIP_Real** val,
2066  SCIP_VAR** vars,
2067  int* constnnonz,
2069  int* constcol,
2070  int* constrow,
2071  SCIP_Real* constval
2072  )
2073 {
2074  SCIP_CONSDATA* consdata;
2075  int i;
2076  const char* name;
2077 
2078  assert( scip != NULL );
2079  assert( cons != NULL );
2080  assert( nvars != NULL );
2081  assert( nnonz != NULL );
2082  assert( blocksize != NULL );
2083  assert( arraylength != NULL );
2084  assert( nvarnonz != NULL );
2085  assert( col != NULL );
2086  assert( row != NULL );
2087  assert( val != NULL );
2088  assert( vars != NULL );
2089  assert( constnnonz != NULL );
2090 
2091  consdata = SCIPconsGetData(cons);
2092  name = SCIPconsGetName(cons);
2093 
2094  assert( consdata->constnnonz == 0 || ( constcol != NULL && constrow != NULL && constval != NULL ) );
2095 
2096  *nvars = consdata->nvars;
2097  *nnonz = consdata->nnonz;
2098  *blocksize = consdata->blocksize;
2099 
2100  for (i = 0; i < consdata->nvars; i++)
2101  vars[i] = consdata->vars[i];
2102 
2103  /* check that the sdp-arrays are long enough to store the information */
2104  if ( *arraylength < consdata->nvars )
2105  {
2106  SCIPdebugMessage("nvarnonz, col, row and val arrays were not long enough to store the information for cons %s, they need to be at least"
2107  "size %d, given was only length %d! \n", name, consdata->nvars, *arraylength);
2108  *arraylength = consdata->nvars;
2109  }
2110  else
2111  {
2112  for (i = 0; i < consdata->nvars; i++)
2113  {
2114  nvarnonz[i] = consdata->nvarnonz[i];
2115  /* set the pointers for each variable */
2116  col[i] = consdata->col[i];
2117  row[i] = consdata->row[i];
2118  val[i] = consdata->val[i];
2119  }
2120  }
2121 
2122  /* set the constant pointers (if a constant part exists) */
2123  if ( consdata->constnnonz > 0 )
2124  {
2125  if ( consdata->constnnonz > *constnnonz )
2126  {
2127  SCIPdebugMessage("The constant nonzeros arrays were not long enough to store the information for cons %s, they need to be at least"
2128  "size %d, given was only length %d! \n", name, consdata->constnnonz, *constnnonz);
2129  }
2130  else
2131  {
2132  for (i = 0; i < consdata->constnnonz; i++)
2133  {
2134  constcol[i] = consdata->constcol[i];
2135  constrow[i] = consdata->constrow[i];
2136  constval[i] = consdata->constval[i];
2137  }
2138  }
2139  }
2140 
2141  *constnnonz = consdata->constnnonz;
2142 
2143  return SCIP_OKAY;
2144 }
2145 
2150 SCIP_RETCODE SCIPconsSdpGetNNonz(
2151  SCIP* scip,
2152  SCIP_CONS* cons,
2153  int* nnonz,
2154  int* constnnonz
2155  )
2156 {
2157  SCIP_CONSDATA* consdata;
2158 
2159  assert( scip != NULL );
2160  assert( cons != NULL );
2161 
2162  consdata = SCIPconsGetData(cons);
2163  assert( consdata != NULL );
2164 
2165  if ( nnonz != NULL )
2166  *nnonz = consdata->nnonz;
2167 
2168  if ( constnnonz != NULL )
2169  *constnnonz = consdata->constnnonz;
2170 
2171  return SCIP_OKAY;
2172 }
2173 
2175 SCIP_RETCODE SCIPconsSdpGetFullAj(
2176  SCIP* scip,
2177  SCIP_CONS* cons,
2178  int j,
2179  SCIP_Real* Aj
2180  )
2181 {
2182  SCIP_CONSDATA* consdata;
2183  int blocksize;
2184  int i;
2185 
2186  assert( scip != NULL );
2187  assert( cons != NULL );
2188  assert( j >= 0 );
2189  assert( Aj != NULL );
2190 
2191  consdata = SCIPconsGetData(cons);
2192  assert( consdata != NULL );
2193  blocksize = consdata->blocksize;
2194 
2195  assert( j < consdata->nvars );
2196 
2197  for (i = 0; i < blocksize * blocksize; i++)
2198  Aj[i] = 0;
2199 
2200  for (i = 0; i < consdata->nvarnonz[j]; i++)
2201  {
2202  Aj[consdata->col[j][i] * blocksize + consdata->row[j][i]] = consdata->val[j][i]; /*lint !e679*/
2203  Aj[consdata->row[j][i] * blocksize + consdata->col[j][i]] = consdata->val[j][i]; /*lint !e679*/
2204  }
2205 
2206  return SCIP_OKAY;
2207 }
2208 
2210 SCIP_RETCODE SCIPconsSdpGetFullConstMatrix(
2211  SCIP* scip,
2212  SCIP_CONS* cons,
2213  SCIP_Real* mat
2214  )
2215 {
2216  SCIP_CONSDATA* consdata;
2217  int blocksize;
2218  int i;
2219  int j;
2220 
2221  assert( scip != NULL );
2222  assert( cons != NULL );
2223  assert( mat != NULL );
2224 
2225  consdata = SCIPconsGetData(cons);
2226  blocksize = consdata->blocksize;
2227 
2228  for (i = 0; i < blocksize; i++)
2229  {
2230  for (j = 0; j < blocksize; j++)
2231  mat[i * blocksize + j] = 0.0; /*lint !e679*/
2232  }
2233 
2234  for (i = 0; i < consdata->constnnonz; i++)
2235  {
2236  mat[consdata->constcol[i] * blocksize + consdata->constrow[i]] = consdata->constval[i]; /*lint !e679*/
2237  mat[consdata->constrow[i] * blocksize + consdata->constcol[i]] = consdata->constval[i]; /*lint !e679*/
2238  }
2239 
2240  return SCIP_OKAY;
2241 }
2242 
2244 SCIP_RETCODE SCIPconsSdpGetLowerTriangConstMatrix(
2245  SCIP* scip,
2246  SCIP_CONS* cons,
2247  SCIP_Real* mat
2248  )
2249 {
2250  SCIP_CONSDATA* consdata;
2251  int blocksize;
2252  int i;
2253 
2254  assert( scip != NULL );
2255  assert( cons != NULL );
2256  assert( mat != NULL );
2257 
2258  consdata = SCIPconsGetData(cons);
2259  assert( consdata != NULL );
2260 
2261  blocksize = consdata->blocksize;
2262 
2263  /* initialize the matrix with 0 */
2264  for (i = 0; i < (blocksize * (blocksize + 1)) / 2; i++)
2265  mat[i] = 0.0;
2266 
2267  for (i = 0; i < consdata->constnnonz; i++)
2268  mat[compLowerTriangPos(consdata->constrow[i], consdata->constcol[i])] = consdata->constval[i];
2269 
2270  return SCIP_OKAY;
2271 }
2272 
2274 SCIP_RETCODE SCIPcreateConsSdp(
2275  SCIP* scip,
2276  SCIP_CONS** cons,
2277  const char* name,
2278  int nvars,
2279  int nnonz,
2280  int blocksize,
2281  int* nvarnonz,
2282  int** col,
2283  int** row,
2284  SCIP_Real** val,
2285  SCIP_VAR** vars,
2286  int constnnonz,
2287  int* constcol,
2288  int* constrow,
2289  SCIP_Real* constval
2290  )
2291 {
2292  SCIP_CONSHDLR* conshdlr;
2293  SCIP_CONSDATA* consdata = NULL;
2294  int i;
2295  int j;
2296 
2297  assert( scip != NULL );
2298  assert( cons != NULL );
2299  assert( name != NULL );
2300  assert( nvars >= 0 );
2301  assert( nnonz >= 0 );
2302  assert( blocksize >= 0 );
2303  assert( constnnonz >= 0 );
2304  assert( nvars == 0 || vars != NULL );
2305  assert( nnonz == 0 || (nvarnonz != NULL && col != NULL && row != NULL && val != NULL ));
2306  assert( constnnonz == 0 || (constcol != NULL && constrow != NULL && constval != NULL ));
2307 
2308  conshdlr = SCIPfindConshdlr(scip, "SDP");
2309  if ( conshdlr == NULL )
2310  {
2311  SCIPerrorMessage("SDP constraint handler not found\n");
2312  return SCIP_PLUGINNOTFOUND;
2313  }
2314 
2315  /* create constraint data */
2316  SCIP_CALL( SCIPallocBlockMemory(scip, &consdata) );
2317  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->nvarnonz, nvars) );
2318  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col, nvars) );
2319  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row, nvars) );
2320  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val, nvars) );
2321  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constcol, constnnonz) );
2322  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constrow, constnnonz) );
2323  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constval, constnnonz) );
2324  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->vars, nvars));
2325 
2326  for (i = 0; i < nvars; i++)
2327  {
2328  assert( nvarnonz[i] >= 0 );
2329 
2330  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col[i], nvarnonz[i]));
2331  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row[i], nvarnonz[i]));
2332  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val[i], nvarnonz[i]));
2333  }
2334 
2335  consdata->nvars = nvars;
2336  consdata->nnonz = nnonz;
2337  consdata->constnnonz = constnnonz;
2338  consdata->blocksize = blocksize;
2339 
2340  for (i = 0; i < nvars; i++)
2341  {
2342  consdata->nvarnonz[i] = nvarnonz[i];
2343 
2344  for (j = 0; j < nvarnonz[i]; j++)
2345  {
2346  assert( col[i][j] >= 0 );
2347  assert( col[i][j] < blocksize );
2348  assert( row[i][j] >= 0 );
2349  assert( row[i][j] < blocksize );
2350 
2351  consdata->col[i][j] = col[i][j];
2352  consdata->row[i][j] = row[i][j];
2353  consdata->val[i][j] = val[i][j];
2354  }
2355  }
2356 
2357  for (i = 0; i < constnnonz; i++)
2358  {
2359  consdata->constcol[i] = constcol[i];
2360  consdata->constrow[i] = constrow[i];
2361  consdata->constval[i] = constval[i];
2362  }
2363 
2364  for (i = 0; i < nvars; i++)
2365  {
2366  consdata->vars[i] = vars[i];
2367  SCIP_CALL( SCIPcaptureVar(scip, consdata->vars[i]) );
2368  }
2369 
2370  /* create constraint */
2371  SCIP_CALL( SCIPcreateCons(scip, cons, name, conshdlr, consdata, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE) );
2372 
2373  /* compute maximum rhs entry for later use in the DIMACS Error Norm */
2374  SCIP_CALL( setMaxRhsEntry(*cons) );
2375 
2376  return SCIP_OKAY;
2377 }
SCIP_DECL_CONSTRANS
static SCIP_DECL_CONSTRANS(consTransSdp)
Definition: cons_sdp.c:1493
SCIPlapackMatrixVectorMult
EXTERN SCIP_RETCODE SCIPlapackMatrixVectorMult(int nrows, int ncols, SCIP_Real *matrix, SCIP_Real *vector, SCIP_Real *result)
Definition: lapack_dsdp.c:213
separateSol
static SCIP_RETCODE separateSol(SCIP *scip, SCIP_CONSHDLR *conshdlr, SCIP_CONS *cons, SCIP_SOL *sol, SCIP_RESULT *result)
Definition: cons_sdp.c:400
SCIPconsSdpGetFullAj
SCIP_RETCODE SCIPconsSdpGetFullAj(SCIP *scip, SCIP_CONS *cons, int j, SCIP_Real *Aj)
Definition: cons_sdp.c:2176
SCIPsdpVarfixerSortRowCol
void SCIPsdpVarfixerSortRowCol(int *row, int *col, SCIP_Real *val, int length)
Definition: SdpVarfixer.c:46
CONSHDLR_PRESOLTIMING
#define CONSHDLR_PRESOLTIMING
Definition: cons_sdp.c:69
CONSHDLR_NAME
#define CONSHDLR_NAME
Definition: cons_sdp.c:56
SCIPincludeConshdlrSdp
SCIP_RETCODE SCIPincludeConshdlrSdp(SCIP *scip)
Definition: cons_sdp.c:2013
multiplyConstraintMatrix
static SCIP_RETCODE multiplyConstraintMatrix(SCIP_CONS *cons, int j, SCIP_Real *v, SCIP_Real *vAv)
Definition: cons_sdp.c:195
CONSHDLR_NEEDSCONS
#define CONSHDLR_NEEDSCONS
Definition: cons_sdp.c:67
diagGEzero
static SCIP_RETCODE diagGEzero(SCIP *scip, SCIP_CONS **conss, int nconss, int *naddconss)
Definition: cons_sdp.c:492
SCIP_DECL_CONSDELETE
static SCIP_DECL_CONSDELETE(consDeleteSdp)
Definition: cons_sdp.c:1756
CONSHDLR_SEPAFREQ
#define CONSHDLR_SEPAFREQ
Definition: cons_sdp.c:61
CONSHDLR_DELAYSEPA
#define CONSHDLR_DELAYSEPA
Definition: cons_sdp.c:66
computeSdpMatrix
static SCIP_RETCODE computeSdpMatrix(SCIP *scip, SCIP_CONS *cons, SCIP_SOL *y, SCIP_Real *matrix)
Definition: cons_sdp.c:153
CONSHDLR_EAGERFREQ
#define CONSHDLR_EAGERFREQ
Definition: cons_sdp.c:62
move_1x1_blocks_to_lp
static SCIP_RETCODE move_1x1_blocks_to_lp(SCIP *scip, SCIP_CONS **conss, int nconss, int *naddconss, int *ndelconss, SCIP_RESULT *result)
Definition: cons_sdp.c:756
SCIPconsSdpGetLowerTriangConstMatrix
SCIP_RETCODE SCIPconsSdpGetLowerTriangConstMatrix(SCIP *scip, SCIP_CONS *cons, SCIP_Real *mat)
Definition: cons_sdp.c:2245
SCIP_DECL_CONSSEPALP
static SCIP_DECL_CONSSEPALP(consSepalpSdp)
Definition: cons_sdp.c:1741
setMaxRhsEntry
static SCIP_RETCODE setMaxRhsEntry(SCIP_CONS *cons)
Definition: cons_sdp.c:235
multiaggrVar
static SCIP_RETCODE multiaggrVar(SCIP *scip, SCIP_CONS *cons, int *v, SCIP_VAR **aggrvars, SCIP_Real *scalars, int naggrvars, SCIP_Real constant, int *savedcol, int *savedrow, SCIP_Real *savedval, int *nfixednonz, int *vararraylength)
Definition: cons_sdp.c:962
SCIP_DECL_CONSCHECK
static SCIP_DECL_CONSCHECK(consCheckSdp)
Definition: cons_sdp.c:1553
CONSHDLR_DESC
#define CONSHDLR_DESC
Definition: cons_sdp.c:57
cons_sdp.h
Constraint handler for SDP-constraints.
expandSymMatrix
static SCIP_RETCODE expandSymMatrix(int size, SCIP_Real *symMat, SCIP_Real *fullMat)
Definition: cons_sdp.c:120
SdpVarmapper.h
class that maps SCIP variables to SDP indices (the SCIP variables are given SDP indices in the order ...
SCIPconsSdpGetFullConstMatrix
SCIP_RETCODE SCIPconsSdpGetFullConstMatrix(SCIP *scip, SCIP_CONS *cons, SCIP_Real *mat)
Definition: cons_sdp.c:2211
SCIPcreateConsSdp
SCIP_RETCODE SCIPcreateConsSdp(SCIP *scip, SCIP_CONS **cons, const char *name, int nvars, int nnonz, int blocksize, int *nvarnonz, int **col, int **row, SCIP_Real **val, SCIP_VAR **vars, int constnnonz, int *constcol, int *constrow, SCIP_Real *constval)
Definition: cons_sdp.c:2275
SCIPconsSdpGetData
SCIP_RETCODE SCIPconsSdpGetData(SCIP *scip, SCIP_CONS *cons, int *nvars, int *nnonz, int *blocksize, int *arraylength, int *nvarnonz, int **col, int **row, SCIP_Real **val, SCIP_VAR **vars, int *constnnonz, int *constcol, int *constrow, SCIP_Real *constval)
Definition: cons_sdp.c:2055
CONSHDLR_CHECKPRIORITY
#define CONSHDLR_CHECKPRIORITY
Definition: cons_sdp.c:60
SCIPsdpVarfixerMergeArrays
SCIP_RETCODE SCIPsdpVarfixerMergeArrays(BMS_BLKMEM *blkmem, SCIP_Real feastol, int *originrow, int *origincol, SCIP_Real *originval, int originlength, SCIP_Bool originsorted, SCIP_Real scalar, int *targetrow, int *targetcol, SCIP_Real *targetval, int *targetlength, int targetmemory)
Definition: SdpVarfixer.c:83
SCIP_DECL_CONSENFOPS
static SCIP_DECL_CONSENFOPS(consEnfopsSdp)
Definition: cons_sdp.c:1577
compLowerTriangPos
static int compLowerTriangPos(int i, int j)
Definition: cons_sdp.c:104
SdpVarfixer.h
adds the main functionality to fix/unfix/(multi-)aggregate variables by merging two three-tuple-array...
CONSHDLR_ENFOPRIORITY
#define CONSHDLR_ENFOPRIORITY
Definition: cons_sdp.c:59
SCIP_DECL_CONSFREE
static SCIP_DECL_CONSFREE(consFreeSdp)
Definition: cons_sdp.c:1793
SCIP_DECL_CONSPRESOL
static SCIP_DECL_CONSPRESOL(consPresolSdp)
Definition: cons_sdp.c:1471
SCIP_DECL_CONSEXITPRE
static SCIP_DECL_CONSEXITPRE(consExitpreSdp)
Definition: cons_sdp.c:1459
SCIP_DECL_CONSPRINT
static SCIP_DECL_CONSPRINT(consPrintSdp)
Definition: cons_sdp.c:1878
SCIP_DECL_CONSHDLRCOPY
static SCIP_DECL_CONSHDLRCOPY(conshdlrCopySdp)
Definition: cons_sdp.c:1808
CONSHDLR_MAXPREROUNDS
#define CONSHDLR_MAXPREROUNDS
Definition: cons_sdp.c:65
SCIP_DECL_CONSINITPRE
static SCIP_DECL_CONSINITPRE(consInitpreSdp)
Definition: cons_sdp.c:1384
SCIPconsSdpCheckSdpCons
SCIP_RETCODE SCIPconsSdpCheckSdpCons(SCIP *scip, SCIP_CONS *cons, SCIP_SOL *sol, SCIP_Bool checkintegrality, SCIP_Bool checklprows, SCIP_Bool printreason, SCIP_RESULT *result)
Definition: cons_sdp.c:337
SCIP_DECL_CONSSEPASOL
static SCIP_DECL_CONSSEPASOL(consSepasolSdp)
Definition: cons_sdp.c:1726
SCIP_DECL_CONSGETVARS
static SCIP_DECL_CONSGETVARS(consGetVarsSdp)
Definition: cons_sdp.c:1959
CONSHDLR_SEPAPRIORITY
#define CONSHDLR_SEPAPRIORITY
Definition: cons_sdp.c:58
SCIPconsSdpGetNNonz
SCIP_RETCODE SCIPconsSdpGetNNonz(SCIP *scip, SCIP_CONS *cons, int *nnonz, int *constnnonz)
Definition: cons_sdp.c:2151
fixAndAggrVars
static SCIP_RETCODE fixAndAggrVars(SCIP *scip, SCIP_CONS **conss, int nconss, SCIP_Bool aggregate)
Definition: cons_sdp.c:1166
SCIP_DECL_CONSLOCK
static SCIP_DECL_CONSLOCK(consLockSdp)
Definition: cons_sdp.c:1403
SCIP_DECL_CONSGETNVARS
static SCIP_DECL_CONSGETNVARS(consGetNVarsSdp)
Definition: cons_sdp.c:1994
SCIP_DECL_CONSCOPY
static SCIP_DECL_CONSCOPY(consCopySdp)
Definition: cons_sdp.c:1823
lapack.h
interface methods for eigenvector computation and matrix multiplication using different versions of L...
cutUsingEigenvector
static SCIP_RETCODE cutUsingEigenvector(SCIP *scip, SCIP_CONS *cons, SCIP_SOL *sol, SCIP_Real *coeff, SCIP_Real *lhs)
Definition: cons_sdp.c:268
SCIP_DECL_CONSENFOLP
static SCIP_DECL_CONSENFOLP(consEnfolpSdp)
Definition: cons_sdp.c:1616
SCIPlapackComputeIthEigenvalue
EXTERN SCIP_RETCODE SCIPlapackComputeIthEigenvalue(BMS_BLKMEM *blkmem, SCIP_Bool geteigenvectors, int n, SCIP_Real *A, int i, SCIP_Real *eigenvalue, SCIP_Real *eigenvector)
Generated on Fri Jun 26 2015 15:51:35 for SCIP-SDP by   doxygen 1.8.6