59 #include "scip/cons_linear.h"
60 #include "scip/scip.h"
69 #define CONSHDLR_NAME "SDP"
70 #define CONSHDLR_DESC "SDP constraints of the form \\sum_{j} A_j y_j - A_0 psd"
71 #define CONSHDLR_SEPAPRIORITY +1000000
72 #define CONSHDLR_ENFOPRIORITY -2000000
73 #define CONSHDLR_CHECKPRIORITY -2000000
74 #define CONSHDLR_SEPAFREQ 1
75 #define CONSHDLR_EAGERFREQ 100
77 #define CONSHDLR_MAXPREROUNDS -1
78 #define CONSHDLR_DELAYSEPA FALSE
79 #define CONSHDLR_NEEDSCONS TRUE
81 #define CONSHDLR_PRESOLTIMING SCIP_PRESOLTIMING_FAST
82 #define PARSE_STARTSIZE 1
83 #define PARSE_SIZEFACTOR 10
84 #define DEFAULT_DIAGGEZEROCUTS TRUE
85 #define DEFAULT_DIAGZEROIMPLCUTS TRUE
87 #define DEFAULT_NTHREADS 1
105 SCIP_Real maxrhsentry;
109 struct SCIP_ConshdlrData
112 SCIP_Bool diaggezerocuts;
115 SCIP_Bool diagzeroimplcuts;
134 return i*(i+1)/2 + j;
137 #define compLowerTriangPos(i, j) (i*(i+1)/2 + j)
153 assert( symMat != NULL );
154 assert( fullMat != NULL );
157 for (i = 0; i < size; i++)
159 for (j = 0; j <= i; j++)
162 fullMat[i*size + j] = symMat[ind];
163 fullMat[j*size + i] = symMat[ind];
182 SCIP_CONSDATA* consdata;
189 assert( cons != NULL );
190 assert( matrix != NULL );
192 consdata = SCIPconsGetData(cons);
193 nvars = consdata->nvars;
194 blocksize = consdata->blocksize;
197 for (i = 0; i < (blocksize * (blocksize + 1))/2; i++)
201 for (i = 0; i < nvars; i++)
203 yval = SCIPgetSolVal(scip, y, consdata->vars[i]);
204 for (ind = 0; ind < consdata->nvarnonz[i]; ind++)
205 matrix[
compLowerTriangPos(consdata->row[i][ind], consdata->col[i][ind])] += yval * consdata->val[i][ind];
209 for (ind = 0; ind < consdata->constnnonz; ind++)
210 matrix[
compLowerTriangPos(consdata->constrow[ind], consdata->constcol[ind])] -= consdata->constval[ind];
224 SCIP_CONSDATA* consdata;
227 assert( cons != NULL );
229 assert( vAv != NULL );
231 consdata = SCIPconsGetData(cons);
232 assert( consdata != NULL );
234 assert( j < consdata->nvars );
239 for (i = 0; i < consdata->nvarnonz[j]; i++)
241 if (consdata->col[j][i] == consdata->row[j][i])
242 *vAv += v[consdata->col[j][i]] * consdata->val[j][i] * v[consdata->row[j][i]];
246 *vAv += 2.0 * v[consdata->col[j][i]] * consdata->val[j][i] * v[consdata->row[j][i]];
261 SCIP_CONSDATA* consdata;
265 consdata = SCIPconsGetData(cons);
266 assert( consdata != NULL );
272 for (i = 0; i < consdata->constnnonz; i++)
274 if ( REALABS(consdata->constval[i]) > max )
275 max = REALABS(consdata->constval[i]);
278 consdata->maxrhsentry = max;
298 SCIP_CONSDATA* consdata;
299 SCIP_Real* eigenvalues = NULL;
300 SCIP_Real* matrix = NULL;
301 SCIP_Real* fullmatrix = NULL;
302 SCIP_Real* fullconstmatrix = NULL;
303 SCIP_Real* eigenvector = NULL;
304 SCIP_Real* output_vector = NULL;
308 assert( cons != NULL );
309 assert( lhs != NULL );
313 consdata = SCIPconsGetData(cons);
314 assert( consdata != NULL );
316 blocksize = consdata->blocksize;
318 SCIP_CALL( SCIPallocBufferArray(scip, &eigenvalues, blocksize) );
319 SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize+1))/2 ) );
320 SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, blocksize * blocksize ) );
321 SCIP_CALL( SCIPallocBufferArray(scip, &fullconstmatrix, blocksize * blocksize) );
322 SCIP_CALL( SCIPallocBufferArray(scip, &eigenvector, blocksize) );
323 SCIP_CALL( SCIPallocBufferArray(scip, &output_vector, blocksize) );
339 for (j = 0; j < blocksize; ++j)
340 *lhs += eigenvector[j] * output_vector[j];
343 for (j = 0; j < consdata->nvars; ++j)
348 SCIPfreeBufferArray(scip, &output_vector);
349 SCIPfreeBufferArray(scip, &eigenvector);
350 SCIPfreeBufferArray(scip, &fullconstmatrix);
351 SCIPfreeBufferArray(scip, &fullmatrix);
352 SCIPfreeBufferArray(scip, &matrix);
353 SCIPfreeBufferArray(scip, &eigenvalues);
363 SCIP_Bool checkintegrality,
364 SCIP_Bool checklprows,
365 SCIP_Bool printreason,
369 SCIP_CONSDATA* consdata;
371 SCIP_Real check_value;
372 SCIP_Real eigenvalue;
373 SCIP_Real* matrix = NULL;
374 SCIP_Real* fullmatrix = NULL;
376 assert( scip != NULL );
377 assert( cons != NULL );
378 assert( result != NULL );
379 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(cons)),
"SDP") == 0);
381 consdata = SCIPconsGetData(cons);
382 assert( consdata != NULL );
383 blocksize = consdata->blocksize;
385 SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize+1)) / 2) );
386 SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, blocksize * blocksize) );
395 check_value = (-eigenvalue) / (1.0 + consdata->maxrhsentry);
397 check_value = -eigenvalue;
400 if ( SCIPisFeasLE(scip, check_value, 0.0) )
401 *result = SCIP_FEASIBLE;
404 *result = SCIP_INFEASIBLE;
408 SCIP_CALL( SCIPprintCons(scip, cons, NULL) );
409 SCIP_CALL( SCIPprintSol(scip, sol, NULL, FALSE) );
410 SCIPinfoMessage(scip, NULL,
"non psd matrix (eigenvalue %f, dimacs error norm = %f).\n", eigenvalue, check_value);
415 SCIPfreeBufferArray(scip, &fullmatrix);
416 SCIPfreeBufferArray(scip, &matrix);
425 SCIP_CONSHDLR* conshdlr,
431 char cutname[SCIP_MAXSTRLEN];
432 SCIP_CONSDATA* consdata;
433 SCIP_CONSHDLRDATA* conshdlrdata;
435 SCIP_Real* coeff = NULL;
446 assert( cons != NULL );
448 consdata = SCIPconsGetData(cons);
449 assert( consdata != NULL );
451 nvars = consdata->nvars;
452 SCIP_CALL( SCIPallocBufferArray(scip, &coeff, nvars ) );
456 SCIP_CALL( SCIPallocBufferArray(scip, &cols, nvars ) );
457 SCIP_CALL( SCIPallocBufferArray(scip, &vals, nvars ) );
460 for (j = 0; j < nvars; ++j)
462 if ( SCIPisZero(scip, coeff[j]) )
465 cols[len] = SCIPvarGetCol(SCIPvarGetProbvar(consdata->vars[j]));
466 vals[len] = coeff[j];
470 conshdlrdata = SCIPconshdlrGetData(conshdlr);
471 assert( conshdlrdata != NULL );
474 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
475 assert( snprintfreturn < SCIP_MAXSTRLEN );
477 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
479 SCIP_CALL( SCIPcreateRowCons(scip, &row, conshdlr, cutname , len, cols, vals, lhs, SCIPinfinity(scip), FALSE, FALSE, TRUE) );
481 if ( SCIPisCutEfficacious(scip, sol, row) )
483 SCIP_Bool infeasible;
484 #ifdef SCIP_MORE_DEBUG
485 SCIPinfoMessage(scip, NULL,
"Added cut %s: ", cutname);
486 SCIPinfoMessage(scip, NULL,
"%f <= ", lhs);
487 for (j = 0; j < len; j++)
488 SCIPinfoMessage(scip, NULL,
"+ (%f)*%s", vals[j], SCIPvarGetName(SCIPcolGetVar(cols[j])));
489 SCIPinfoMessage(scip, NULL,
"\n");
492 SCIP_CALL( SCIPaddCut(scip, sol, row, FALSE, &infeasible) );
494 *result = SCIP_CUTOFF;
496 *result = SCIP_SEPARATED;
497 SCIP_CALL( SCIPresetConsAge(scip, cons) );
500 SCIP_CALL( SCIPreleaseRow(scip, &row) );
502 SCIPfreeBufferArray(scip, &vals);
503 SCIPfreeBufferArray(scip, &cols);
504 SCIPfreeBufferArray(scip, &coeff);
520 char cutname[SCIP_MAXSTRLEN];
521 SCIP_CONSDATA* consdata;
523 SCIP_Real* cons_array;
524 SCIP_Real* lhs_array;
536 for (c = 0; c < nconss; ++c)
538 SCIP_CONSHDLR* conshdlr;
540 conshdlr = SCIPconsGetHdlr(conss[c]);
541 assert( conshdlr != NULL );
543 assert( strcmp(SCIPconshdlrGetName(conshdlr),
"SDP") == 0);
546 consdata = SCIPconsGetData(conss[c]);
547 assert( consdata != NULL );
549 blocksize = consdata->blocksize;
550 nvars = consdata->nvars;
551 rhs = SCIPinfinity(scip);
553 SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize + 1)) / 2) );
556 SCIP_CALL( SCIPallocBufferArray(scip, &cons_array, blocksize * nvars) );
557 SCIP_CALL( SCIPallocBufferArray(scip, &lhs_array, blocksize) );
560 for (k = 0; k < blocksize; ++k)
564 for (j = 0; j < nvars; ++j)
566 for (k = 0; k < blocksize; ++k)
569 cons_array[k * nvars + j] = 0.0;
573 for (i = 0; i < consdata->nvarnonz[j]; i++)
575 if ( consdata->col[j][i] == consdata->row[j][i] )
576 cons_array[consdata->col[j][i] * nvars + j] = consdata->val[j][i];
581 for (k = 0; k < blocksize; ++k)
584 SCIP_CONSHDLRDATA* conshdlrdata;
586 conshdlrdata = SCIPconshdlrGetData(conshdlr);
588 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_ge_zero_%d", ++(conshdlrdata->ndiaggezerocuts));
589 assert( snprintfreturn < SCIP_MAXSTRLEN );
591 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_ge_zero_%d", ++(conshdlrdata->ndiaggezerocuts));
594 #ifdef SCIP_MORE_DEBUG
595 SCIPinfoMessage(scip, NULL,
"Added lp-constraint %s: ", cutname);
596 SCIPinfoMessage(scip, NULL,
"%f <= ", lhs_array[k]);
597 for (i = 0; i < consdata->nvars; i++)
599 if ( ! SCIPisZero(scip, cons_array[k * consdata->nvars + i]) )
600 SCIPinfoMessage(scip, NULL,
"+ (%f)*%s", cons_array[k * consdata->nvars + i], SCIPvarGetName(consdata->vars[i]));
602 SCIPinfoMessage(scip, NULL,
"\n");
605 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, cutname, consdata->nvars, consdata->vars, cons_array + k * consdata->nvars, lhs_array[k], rhs,
606 TRUE, TRUE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
608 SCIP_CALL( SCIPaddCons(scip, cons) );
609 SCIP_CALL( SCIPreleaseCons(scip, &cons) );
613 SCIPfreeBufferArray(scip, &lhs_array);
614 SCIPfreeBufferArray(scip, &cons_array);
615 SCIPfreeBufferArray(scip, &matrix);
636 char cutname[SCIP_MAXSTRLEN];
641 int* nconstnonzeroentries;
642 int** constnonzeroentries;
643 SCIP_CONSDATA* consdata;
656 SCIP_Bool anycutvalid;
661 assert( scip != NULL );
662 assert( conss != NULL );
663 assert( nconss >= 0 );
664 assert( naddconss != NULL );
666 for (i = 0; i < nconss; ++i)
668 assert( conss[i] != NULL );
669 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[i])),
"SDP") == 0 );
671 consdata = SCIPconsGetData(conss[i]);
672 assert( consdata != NULL );
674 blocksize = consdata->blocksize;
675 nvars = consdata->nvars;
676 SCIP_CALL( SCIPallocBufferArray(scip, &nconstnonzeroentries, blocksize) );
677 SCIP_CALL( SCIPallocBufferArray(scip, &constnonzeroentries, blocksize) );
678 for (j = 0; j < blocksize; j++)
680 nconstnonzeroentries[j] = 0;
681 SCIP_CALL( SCIPallocBufferArray(scip, &constnonzeroentries[j], blocksize) );
685 for (j = 0; j < consdata->constnnonz; j++)
689 if ( (consdata->constcol[j] != consdata->constrow[j]) )
691 assert( ! SCIPisZero(scip, consdata->constval[j]) );
692 if ( nconstnonzeroentries[consdata->constcol[j]] >= 0 )
694 constnonzeroentries[consdata->constcol[j]][nconstnonzeroentries[consdata->constcol[j]]] = consdata->constrow[j];
695 nconstnonzeroentries[consdata->constcol[j]]++;
697 if ( nconstnonzeroentries[consdata->constrow[j]] >= 0 )
699 constnonzeroentries[consdata->constrow[j]][nconstnonzeroentries[consdata->constrow[j]]] = consdata->constcol[j];
700 nconstnonzeroentries[consdata->constrow[j]]++;
706 assert( ! SCIPisZero(scip, consdata->constval[j]) );
707 nconstnonzeroentries[consdata->constcol[j]] = -1;
713 SCIP_CALL( SCIPallocBufferArray(scip, &diagvars, blocksize) );
714 SCIP_CALL( SCIPallocBufferArray(scip, &ndiagvars, blocksize) );
716 for (j = 0; j < blocksize; ++j)
719 if ( nconstnonzeroentries[j] > 0 )
721 SCIP_CALL( SCIPallocBufferArray(scip, &(diagvars[j]), nvars) );
729 SCIPfreeBufferArray(scip, &ndiagvars);
730 SCIPfreeBufferArray(scip, &diagvars);
731 for (j = blocksize - 1; j >= 0; j--)
733 SCIPfreeBufferArray(scip, &constnonzeroentries[j]);
735 SCIPfreeBufferArray(scip, &constnonzeroentries);
736 SCIPfreeBufferArray(scip, &nconstnonzeroentries);
742 for (v = 0; v < nvars; v++)
744 for (j = 0; j < consdata->nvarnonz[v]; j++)
748 if ( (consdata->col[v][j] == consdata->row[v][j]) && (nconstnonzeroentries[consdata->col[v][j]] > 0) )
750 if ( SCIPvarIsIntegral(consdata->vars[v]) )
752 assert( ! SCIPisEQ(scip, consdata->val[v][j], 0.0) );
753 diagvars[consdata->col[v][j]][ndiagvars[consdata->col[v][j]]] = v;
754 ndiagvars[consdata->col[v][j]]++;
758 nconstnonzeroentries[consdata->col[v][j]] = -2;
763 else if ( consdata->col[v][j] != consdata->row[v][j] )
765 if ( nconstnonzeroentries[consdata->col[v][j]] > 0 )
768 for (k = 0; k < nconstnonzeroentries[consdata->col[v][j]]; k++)
770 if ( constnonzeroentries[consdata->col[v][j]][k] == consdata->row[v][j] )
773 if ( nconstnonzeroentries[consdata->col[v][j]] > k + 1 )
775 for (l = k + 1; l < nconstnonzeroentries[consdata->col[v][j]]; l++)
776 constnonzeroentries[consdata->col[v][j]][l - 1] = constnonzeroentries[consdata->col[v][j]][l];
778 nconstnonzeroentries[consdata->col[v][j]]--;
780 if ( nconstnonzeroentries[consdata->col[v][j]] == 0 )
781 nconstnonzeroentries[consdata->col[v][j]] = -2;
787 if ( nconstnonzeroentries[consdata->row[v][j]] > 0 )
790 for (k = 0; k < nconstnonzeroentries[consdata->row[v][j]]; k++)
792 if ( constnonzeroentries[consdata->row[v][j]][k] == consdata->col[v][j] )
795 if ( nconstnonzeroentries[consdata->row[v][j]] > k + 1 )
797 for (l = k + 1; l < nconstnonzeroentries[consdata->row[v][j]]; l++)
798 constnonzeroentries[consdata->row[v][j]][l - 1] = constnonzeroentries[consdata->row[v][j]][l];
800 nconstnonzeroentries[consdata->row[v][j]]--;
802 if ( nconstnonzeroentries[consdata->row[v][j]] == 0 )
803 nconstnonzeroentries[consdata->row[v][j]] = -2;
812 for (j = 0; j < blocksize; ++j)
814 if ( nconstnonzeroentries[j] > 0 )
816 SCIP_CALL( SCIPallocBufferArray(scip, &vals, ndiagvars[j]) );
817 SCIP_CALL( SCIPallocBufferArray(scip, &vars, ndiagvars[j]) );
820 for (v = 0; v < ndiagvars[j]; ++v)
822 vars[v] = consdata->vars[diagvars[j][v]];
826 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_zero_impl_block_%d_row_%d", i, j);
827 assert( snprintfreturn < SCIP_MAXSTRLEN );
829 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_zero_impl_block_%d_row_%d", i, j);
832 #ifdef SCIP_MORE_DEBUG
833 SCIPinfoMessage(scip, NULL,
"Added lp-constraint %s: ", cutname);
834 SCIPinfoMessage(scip, NULL,
"1 <=");
835 for (l = 0; l < ndiagvars[j]; l++)
836 SCIPinfoMessage(scip, NULL,
" + (%f)*%s", vals[l], SCIPvarGetName(vars[l]));
837 SCIPinfoMessage(scip, NULL,
"\n");
841 SCIP_CALL(SCIPcreateConsLinear(scip, &cons, cutname , ndiagvars[j], vars, vals, 1.0, SCIPinfinity(scip), TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE));
842 SCIP_CALL(SCIPaddCons(scip, cons));
843 SCIP_CALL(SCIPreleaseCons(scip, &cons));
846 SCIPfreeBufferArray(scip, &vars);
847 SCIPfreeBufferArray(scip, &vals);
849 if ( nconstnonzeroentries[j] == -2 || nconstnonzeroentries[j] > 0 )
851 SCIPfreeBufferArray(scip, &diagvars[j]);
855 SCIPfreeBufferArray(scip, &ndiagvars);
856 SCIPfreeBufferArray(scip, &diagvars);
857 for (j = blocksize - 1; j >= 0; j--)
859 SCIPfreeBufferArray(scip, &constnonzeroentries[j]);
861 SCIPfreeBufferArray(scip, &constnonzeroentries);
862 SCIPfreeBufferArray(scip, &nconstnonzeroentries);
879 char cutname[SCIP_MAXSTRLEN];
881 SCIP_CONSDATA* consdata;
882 SCIP_CONSHDLRDATA* conshdlrdata;
897 assert( result != NULL );
898 *result = SCIP_SUCCESS;
900 for (i = 0; i < nconss; ++i)
902 hdlr = SCIPconsGetHdlr(conss[i]);
903 assert(hdlr != NULL);
906 assert( strcmp(SCIPconshdlrGetName(hdlr),
"SDP") == 0);
909 consdata = SCIPconsGetData(conss[i]);
910 assert( consdata != NULL );
913 if ( consdata->blocksize == 1 )
915 nvars = consdata->nvars;
916 nnonz = consdata->nnonz;
917 SCIP_CALL( SCIPallocBufferArray(scip, &vars, nvars) );
918 SCIP_CALL( SCIPallocBufferArray(scip, &coeffs, nnonz) );
923 for (var = 0; var < nvars; var++)
925 assert( consdata->nvarnonz[var] <= 1 );
927 for (j = 0; j < consdata->nvarnonz[var]; j++)
929 assert( consdata->col[var][j] == 0 && consdata->row[var][j] == 0 );
930 coeffs[count] = consdata->val[var][j];
931 vars[count] = consdata->vars[var];
937 assert( consdata->constnnonz <= 1 );
939 rhs = (consdata->constnnonz == 1) ? consdata->constval[0] : 0.0;
945 conshdlrdata = SCIPconshdlrGetData(hdlr);
947 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"1x1block_%d", ++(conshdlrdata->n1x1blocks));
948 assert( snprintfreturn < SCIP_MAXSTRLEN );
950 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"1x1block_%d", ++(conshdlrdata->n1x1blocks));
953 #ifdef SCIP_MORE_DEBUG
954 SCIPinfoMessage(scip, NULL,
"Added lp-constraint %s: ", cutname);
955 for (i = 0; i < consdata->nvars; i++)
956 SCIPinfoMessage(scip, NULL,
"+ (%f)*%s", coeffs[i], SCIPvarGetName(vars[i]));
957 SCIPinfoMessage(scip, NULL,
" <= %f \n", rhs);
960 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, cutname, consdata->nvars, vars, coeffs, rhs, SCIPinfinity(scip),
961 TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
963 SCIP_CALL( SCIPaddCons(scip, cons) );
964 SCIP_CALL( SCIPreleaseCons(scip, &cons) );
972 if ( coeffs[0] > 0.0 )
975 if ( SCIPisFeasGT(scip, rhs / coeffs[0], SCIPvarGetLbLocal(vars[0])) )
978 if( SCIPisFeasGT(scip, rhs / coeffs[0], SCIPvarGetUbLocal(vars[0])) )
980 SCIPdebugMessage(
"Problem detected to be infeasible during presolving, 1x1-SDP-constraint %s caused change"
981 "of lower bound for variable %s from %f to %f, which is bigger than upper bound of %f\n",
982 SCIPconsGetName(conss[i]), SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]), rhs / coeffs[0],
983 SCIPvarGetUbLocal(vars[0]));
985 *result = SCIP_CUTOFF;
988 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
991 SCIPfreeBufferArray(scip, &coeffs);
992 SCIPfreeBufferArray(scip, &vars);
997 SCIPdebugMessage(
"Changing lower bound of variable %s from %f to %f because of 1x1-SDP-constraint %s!\n",
998 SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]), rhs / coeffs[0], SCIPconsGetName(conss[i]));
999 SCIP_CALL( SCIPchgVarLb(scip, vars[0], rhs / coeffs[0]) );
1003 SCIPdebugMessage(
"Deleting 1x1-SDP-constraint %s, new lower bound %f for variable %s no improvement over old bound %f!\n",
1004 SCIPconsGetName(conss[i]), rhs / coeffs[0], SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]));
1007 else if ( coeffs[0] < 0.0 )
1010 if (SCIPisFeasLT(scip, rhs / coeffs[0], SCIPvarGetUbLocal(vars[0])))
1013 if( SCIPisFeasLT(scip, rhs / coeffs[0], SCIPvarGetLbLocal(vars[0])) )
1015 SCIPdebugMessage(
"Problem detected to be infeasible during presolving, 1x1-SDP-constraint %s caused change"
1016 "of upper bound for variable %s from %f to %f, which is less than lower bound of %f\n",
1017 SCIPconsGetName(conss[i]), SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]), rhs / coeffs[0],
1018 SCIPvarGetLbLocal(vars[0]));
1020 *result = SCIP_CUTOFF;
1023 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
1026 SCIPfreeBufferArray(scip, &coeffs);
1027 SCIPfreeBufferArray(scip, &vars);
1032 SCIPdebugMessage(
"Changing upper bound of variable %s from %f to %f because of 1x1-SDP-constraint %s!\n",
1033 SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]), -rhs / coeffs[0], SCIPconsGetName(conss[i]));
1034 SCIP_CALL( SCIPchgVarUb(scip, vars[0], rhs / coeffs[0]) );
1038 SCIPdebugMessage(
"Deleting 1x1-SDP-constraint %s, new upper bound %f for variable %s no improvement over old bound %f!\n",
1039 SCIPconsGetName(conss[i]), rhs / coeffs[0], SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]));
1044 SCIPdebugMessage(
"Detected 1x1 SDP-block without any nonzero coefficients \n");
1045 if ( SCIPisFeasGT(scip, rhs, 0.0) )
1047 SCIPdebugMessage(
"Detected infeasibility in 1x1 SDP-block without any nonzero coefficients but with strictly positive rhs\n");
1048 *result = SCIP_CUTOFF;
1051 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
1054 SCIPfreeBufferArray(scip, &coeffs);
1055 SCIPfreeBufferArray(scip, &vars);
1063 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
1066 SCIPfreeBufferArray(scip, &coeffs);
1067 SCIPfreeBufferArray(scip, &vars);
1079 SCIP_VAR** aggrvars,
1085 SCIP_Real* savedval,
1091 SCIP_CONSDATA* consdata;
1094 int aggrtargetlength;
1099 assert( scip != NULL );
1100 assert( cons != NULL );
1101 assert( scalars != NULL );
1102 assert( naggrvars > 0 );
1103 assert( savedcol != NULL );
1104 assert( savedrow != NULL );
1105 assert( savedval != NULL );
1106 assert( nfixednonz != NULL );
1108 consdata = SCIPconsGetData(cons);
1109 assert( consdata != NULL );
1111 SCIP_CALL( SCIPgetRealParam(scip,
"numerics/epsilon", &epsilon) );
1114 startind = *nfixednonz;
1116 if ( SCIPisEQ(scip, constant, 0.0) )
1122 for (i = 0; i < consdata->nvarnonz[*v]; i++)
1124 savedcol[*nfixednonz] = consdata->col[*v][i];
1125 savedrow[*nfixednonz] = consdata->row[*v][i];
1126 savedval[*nfixednonz] = consdata->val[*v][i];
1132 for (i = 0; i < consdata->nvarnonz[*v]; i++)
1134 savedcol[*nfixednonz] = consdata->col[*v][i];
1135 savedrow[*nfixednonz] = consdata->row[*v][i];
1136 savedval[*nfixednonz] = consdata->val[*v][i] * constant;
1140 assert( *nfixednonz - startind == consdata->nvarnonz[*v] );
1148 SCIPfreeBlockMemoryArray(scip, &(consdata->val[*v]), consdata->nvarnonz[*v]);
1149 SCIPfreeBlockMemoryArray(scip, &(consdata->row[*v]), consdata->nvarnonz[*v]);
1150 SCIPfreeBlockMemoryArray(scip, &(consdata->col[*v]), consdata->nvarnonz[*v]);
1153 SCIP_CALL( SCIPreleaseVar(scip, &consdata->vars[*v]) );
1154 consdata->col[*v] = consdata->col[consdata->nvars - 1];
1155 consdata->row[*v] = consdata->row[consdata->nvars - 1];
1156 consdata->val[*v] = consdata->val[consdata->nvars - 1];
1157 consdata->nvarnonz[*v] = consdata->nvarnonz[consdata->nvars - 1];
1158 consdata->vars[*v] = consdata->vars[consdata->nvars - 1];
1159 (consdata->nvars)--;
1163 for (aggrind = 0; aggrind < naggrvars; aggrind++)
1167 for (i = 0; i < consdata->nvars; i++)
1169 if ( consdata->vars[i] == aggrvars[aggrind] )
1176 if ( aggrconsind > -1 )
1181 aggrtargetlength = consdata->nvarnonz[aggrconsind] + *nfixednonz - startind;
1182 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->row[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1183 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->col[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1184 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->val[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1186 if ( SCIPisEQ(scip, constant, 0.0) )
1190 SCIP_CALL(
SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), epsilon, savedrow + startind, savedcol + startind, savedval + startind,
1191 *nfixednonz - startind, TRUE, scalars[aggrind], consdata->row[aggrconsind], consdata->col[aggrconsind],
1192 consdata->val[aggrconsind], &(consdata->nvarnonz[aggrconsind]), aggrtargetlength) );
1198 SCIP_CALL(
SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), epsilon, savedrow + startind, savedcol + startind, savedval + startind,
1199 *nfixednonz - startind, TRUE, scalars[aggrind] / constant, consdata->row[aggrconsind], consdata->col[aggrconsind],
1200 consdata->val[aggrconsind], &(consdata->nvarnonz[aggrconsind]), aggrtargetlength) );
1204 assert( consdata->nvarnonz[aggrconsind] <= aggrtargetlength );
1205 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->row[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1206 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->col[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1207 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->val[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1213 SCIPdebugMessage(
"adding variable %s to SDP constraint %s because of (multi-)aggregation\n", SCIPvarGetName(aggrvars[aggrind]), SCIPconsGetName(cons));
1216 if ( consdata->nvars == *vararraylength )
1218 globalnvars = SCIPgetNVars(scip);
1221 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col, *vararraylength, globalnvars) );
1222 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row, *vararraylength, globalnvars) );
1223 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val, *vararraylength, globalnvars) );
1224 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->nvarnonz, *vararraylength, globalnvars) );
1225 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->vars, *vararraylength, globalnvars) );
1226 *vararraylength = globalnvars;
1230 SCIP_CALL( SCIPcaptureVar(scip, aggrvars[aggrind]) );
1231 consdata->vars[consdata->nvars] = aggrvars[aggrind];
1234 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->col[consdata->nvars]), savedcol + startind, *nfixednonz - startind) );
1235 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->row[consdata->nvars]), savedrow + startind, *nfixednonz - startind) );
1241 if ( SCIPisEQ(scip, scalars[aggrind], constant) || SCIPisEQ(scip, constant, 0.0) )
1243 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->val[consdata->nvars]), savedval + startind, *nfixednonz - startind) );
1244 consdata->nvarnonz[consdata->nvars] = *nfixednonz - startind;
1249 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->val[consdata->nvars]), *nfixednonz - startind) );
1251 consdata->nvarnonz[consdata->nvars] = 0;
1253 for (i = 0; i < *nfixednonz - startind; i++)
1256 if ( (scalars[i] / constant) * savedval[startind + i] >= epsilon )
1258 consdata->val[consdata->nvars][consdata->nvarnonz[consdata->nvars]] = (scalars[i] / constant) * savedval[startind + i];
1259 consdata->nvarnonz[consdata->nvars]++;
1264 if ( consdata->nvarnonz[consdata->nvars] > 0 )
1271 if ( SCIPisEQ(scip, constant, 0.0) )
1272 *nfixednonz = startind;
1287 SCIP_CONSDATA* consdata;
1292 SCIP_Real* savedval;
1297 SCIP_VAR** aggrvars;
1314 assert( scip != NULL );
1315 assert( conss != NULL );
1316 assert( nconss >= 0 );
1318 SCIPdebugMessage(
"Calling fixAndAggrVars with aggregate = %u\n", aggregate);
1320 SCIP_CALL( SCIPgetRealParam(scip,
"numerics/epsilon", &epsilon) );
1322 for (c = 0; c < nconss; ++c)
1324 assert( conss[c] != NULL );
1325 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])),
"SDP") == 0);
1327 consdata = SCIPconsGetData(conss[c]);
1328 assert( consdata != NULL );
1331 SCIP_CALL( SCIPallocBufferArray(scip, &savedcol, consdata->nnonz) );
1332 SCIP_CALL( SCIPallocBufferArray(scip, &savedrow, consdata->nnonz) );
1333 SCIP_CALL( SCIPallocBufferArray(scip, &savedval, consdata->nnonz) );
1338 vararraylength = consdata->nvars;
1339 globalnvars = SCIPgetNVars(scip);
1341 for (v = 0; v < consdata->nvars; v++)
1345 if ( SCIPvarIsBinary(consdata->vars[v]) && SCIPvarIsNegated(consdata->vars[v]) )
1348 var = SCIPvarGetNegationVar(consdata->vars[v]);
1351 var = consdata->vars[v];
1354 if ( SCIPvarGetStatus(var) == SCIP_VARSTATUS_FIXED || SCIPisEQ(scip, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var)))
1356 assert( SCIPisEQ(scip, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var)) );
1358 SCIPdebugMessage(
"Globally fixing Variable %s to value %f !\n", SCIPvarGetName(var), SCIPvarGetLbGlobal(var));
1360 if ( ((! negated) && (! SCIPisEQ(scip, SCIPvarGetLbGlobal(var), 0.0))) || (negated && SCIPisEQ(scip, SCIPvarGetLbGlobal(var), 0.0)) )
1365 for (i = 0; i < consdata->nvarnonz[v]; i++)
1367 savedcol[nfixednonz] = consdata->col[v][i];
1368 savedrow[nfixednonz] = consdata->row[v][i];
1371 savedval[nfixednonz] = consdata->val[v][i] * SCIPvarGetLbGlobal(var);
1373 savedval[nfixednonz] = consdata->val[v][i];
1382 consdata->nnonz -= consdata->nvarnonz[v];
1385 SCIPfreeBlockMemoryArrayNull(scip, &(consdata->val[v]), consdata->nvarnonz[v]);
1386 SCIPfreeBlockMemoryArrayNull(scip, &(consdata->row[v]), consdata->nvarnonz[v]);
1387 SCIPfreeBlockMemoryArrayNull(scip, &(consdata->col[v]), consdata->nvarnonz[v]);
1390 SCIP_CALL( SCIPreleaseVar(scip, &(consdata->vars[v])) );
1391 consdata->col[v] = consdata->col[consdata->nvars - 1];
1392 consdata->row[v] = consdata->row[consdata->nvars - 1];
1393 consdata->val[v] = consdata->val[consdata->nvars - 1];
1394 consdata->nvarnonz[v] = consdata->nvarnonz[consdata->nvars - 1];
1395 consdata->vars[v] = consdata->vars[consdata->nvars - 1];
1399 else if ( (SCIPvarGetStatus(var) == SCIP_VARSTATUS_AGGREGATED ||
1400 SCIPvarGetStatus(var) == SCIP_VARSTATUS_MULTAGGR)
1403 SCIP_CALL( SCIPallocBufferArray(scip, &aggrvars, globalnvars) );
1404 SCIP_CALL( SCIPallocBufferArray(scip, &scalars, globalnvars) );
1409 aggrvars[0] = consdata->vars[v];
1417 aggrvars[0] = consdata->vars[v];
1424 SCIP_CALL( SCIPgetProbvarLinearSum(scip, aggrvars, scalars, &naggrvars, globalnvars, &constant, &requiredsize, TRUE) );
1425 assert( requiredsize <= globalnvars );
1430 if ( SCIPvarGetStatus(consdata->vars[v]) == SCIP_VARSTATUS_AGGREGATED )
1431 SCIPdebugMessage(
"aggregating variable %s to ", SCIPvarGetName(var));
1433 SCIPdebugMessage(
"multiaggregating variable %s to ", SCIPvarGetName(var));
1434 for (i = 0; i < naggrvars; i++)
1435 SCIPdebugMessage(
"+ (%f2) * %s ", scalars[i], SCIPvarGetName(aggrvars[i]));
1436 SCIPdebugMessage(
"+ (%f2) \n", constant);
1441 SCIP_CALL(
multiaggrVar(scip, conss[c], &v, aggrvars, scalars, naggrvars, constant, savedcol, savedrow, savedval, &nfixednonz, &vararraylength) );
1443 SCIPfreeBufferArray(scip, &aggrvars);
1444 SCIPfreeBufferArray(scip, &scalars);
1446 else if ( negated && (SCIPvarGetStatus(var) == SCIP_VARSTATUS_LOOSE || SCIPvarGetStatus(var) == SCIP_VARSTATUS_COLUMN)
1451 SCIPdebugMessage(
"Changing variable %s to negation of variable %s !\n", SCIPvarGetName(consdata->vars[v]), SCIPvarGetName(var));
1455 SCIP_CALL(
multiaggrVar(scip, conss[c], &v, &var, &scalar, 1, 1.0, savedcol, savedrow, savedval, &nfixednonz, &vararraylength) );
1460 assert( consdata->nvars <= vararraylength );
1461 if ( consdata->nvars < vararraylength )
1463 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col, vararraylength, consdata->nvars) );
1464 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row, vararraylength, consdata->nvars) );
1465 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val, vararraylength, consdata->nvars) );
1466 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->nvarnonz, vararraylength, consdata->nvars) );
1467 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->vars, vararraylength, consdata->nvars) );
1471 arraylength = consdata->constnnonz + nfixednonz;
1472 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constcol), consdata->constnnonz, arraylength) );
1473 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constrow), consdata->constnnonz, arraylength) );
1474 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constval), consdata->constnnonz, arraylength) );
1477 SCIP_CALL(
SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), epsilon, savedrow, savedcol, savedval, nfixednonz, FALSE, -1.0, consdata->constrow,
1478 consdata->constcol, consdata->constval, &(consdata->constnnonz), arraylength) );
1480 assert( consdata->constnnonz <= arraylength );
1483 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constcol), arraylength, consdata->constnnonz) );
1484 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constrow), arraylength, consdata->constnnonz) );
1485 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constval), arraylength, consdata->constnnonz) );
1488 SCIPfreeBufferArray(scip, &savedval);
1489 SCIPfreeBufferArray(scip, &savedrow);
1490 SCIPfreeBufferArray(scip, &savedcol);
1493 consdata->nnonz = 0;
1494 for (v = 0; v < consdata->nvars; v++)
1495 consdata->nnonz += consdata->nvarnonz[v];
1505 SCIP_CONSHDLR* conshdlr,
1512 char cutname[SCIP_MAXSTRLEN];
1513 SCIP_CONSDATA* consdata;
1514 SCIP_CONSHDLRDATA* conshdlrdata;
1515 SCIP_Bool all_feasible = TRUE;
1516 SCIP_Bool separated = FALSE;
1518 SCIP_Bool infeasible;
1529 *result = SCIP_FEASIBLE;
1531 for (i = 0; i < nconss; ++i)
1533 consdata = SCIPconsGetData(conss[i]);
1535 if ( *result == SCIP_FEASIBLE )
1538 all_feasible = FALSE;
1540 nvars = consdata->nvars;
1544 SCIP_CALL( SCIPallocBufferArray(scip, &coeff, nvars) );
1547 rhs = SCIPinfinity(scip);
1548 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1551 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
1552 assert( snprintfreturn < SCIP_MAXSTRLEN );
1554 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
1557 SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, conshdlr, cutname , lhs, rhs, FALSE, FALSE, TRUE) );
1558 SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
1560 for (j = 0; j < nvars; ++j)
1562 SCIP_CALL( SCIPaddVarToRow(scip, row, consdata->vars[j], coeff[j]) );
1565 SCIP_CALL( SCIPflushRowExtensions(scip, row) );
1567 #ifdef SCIP_MORE_DEBUG
1568 SCIPinfoMessage(scip, NULL,
"Added cut %s: ", cutname);
1569 SCIPinfoMessage(scip, NULL,
"%f <= ", lhs);
1570 for (j = 0; j < nvars; j++)
1571 SCIPinfoMessage(scip, NULL,
"+ (%f)*%s", coeff[j], SCIPvarGetName(consdata->vars[j]));
1572 SCIPinfoMessage(scip, NULL,
"\n");
1575 SCIP_CALL( SCIPaddCut(scip, sol, row, FALSE, &infeasible) );
1579 *result = SCIP_CUTOFF;
1581 SCIP_CALL( SCIPreleaseRow(scip, &row) );
1582 SCIPfreeBufferArray(scip, &coeff);
1588 SCIP_CALL( SCIPaddPoolCut(scip, row) );
1590 SCIP_CALL( SCIPresetConsAge(scip, conss[i]) );
1591 *result = SCIP_SEPARATED;
1594 SCIP_CALL( SCIPreleaseRow(scip, &row) );
1595 SCIPfreeBufferArray(scip, &coeff);
1602 *result = SCIP_SEPARATED;
1611 SCIP_CONSHDLRDATA* conshdlrdata;
1613 assert( conshdlr != NULL );
1615 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1616 assert( conshdlrdata != NULL );
1618 conshdlrdata->neigveccuts = 0;
1619 conshdlrdata->ndiaggezerocuts = 0;
1620 conshdlrdata->n1x1blocks = 0;
1630 SCIP_CONSDATA* consdata;
1634 SCIP_Real eigenvalue;
1636 consdata = SCIPconsGetData(cons);
1637 assert( consdata != NULL );
1638 blocksize = consdata->blocksize;
1639 nvars = consdata->nvars;
1641 SCIP_CALL( SCIPallocBufferArray(scip, &Aj, blocksize * blocksize) );
1643 for (var = 0; var < nvars; var++)
1649 if ( SCIPisNegative(scip, eigenvalue) )
1653 SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlocksneg, nlockspos) );
1657 if ( SCIPisPositive(scip, eigenvalue) )
1661 SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlockspos, nlocksneg) );
1667 if ( SCIPisPositive(scip, eigenvalue) )
1671 SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlockspos, nlocksneg) );
1676 SCIPfreeBufferArray(scip, &Aj);
1685 assert( scip != NULL );
1687 if ( conss == NULL )
1699 SCIP_CONSHDLRDATA* conshdlrdata;
1701 assert( conshdlr != NULL );
1702 assert( result != 0 );
1704 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1705 assert( conshdlrdata != NULL );
1710 if ( conshdlrdata->diaggezerocuts )
1712 SCIP_CALL(
diagGEzero(scip, conss, nconss, naddconss) );
1714 if ( conshdlrdata->diagzeroimplcuts )
1716 SCIP_CALL(
diagZeroImpl(scip, conss, nconss, naddconss) );
1727 SCIP_CONSDATA* sourcedata;
1728 SCIP_CONSDATA* targetdata;
1730 SCIP_CONSHDLRDATA* conshdlrdata;
1736 char transname[SCIP_MAXSTRLEN];
1738 sourcedata = SCIPconsGetData(sourcecons);
1739 assert( sourcedata != NULL );
1741 SCIPdebugMessage(
"Transforming constraint <%s>\n", SCIPconsGetName(sourcecons));
1744 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1745 SCIPdebugMessage(
"Setting number of threads to %d via OpenMP in Openblas.\n", conshdlrdata->nthreads);
1746 omp_set_num_threads(conshdlrdata->nthreads);
1749 SCIP_CALL( SCIPallocBlockMemory(scip, &targetdata) );
1752 targetdata->nvars = sourcedata->nvars;
1753 targetdata->nnonz = sourcedata->nnonz;
1754 targetdata->blocksize = sourcedata->blocksize;
1756 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->nvarnonz), sourcedata->nvarnonz, sourcedata->nvars) );
1759 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->col), sourcedata->nvars) );
1760 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->row), sourcedata->nvars) );
1761 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->val), sourcedata->nvars) );
1763 for (i = 0; i < sourcedata->nvars; i++)
1765 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->col[i]), sourcedata->col[i], sourcedata->nvarnonz[i]) );
1766 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->row[i]), sourcedata->row[i], sourcedata->nvarnonz[i]) );
1767 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->val[i]), sourcedata->val[i], sourcedata->nvarnonz[i]) );
1769 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->vars), sourcedata->nvars));
1772 for (i = 0; i < sourcedata->nvars; i++)
1774 targetdata->vars[i] = SCIPvarGetTransVar(sourcedata->vars[i]);
1775 SCIP_CALL( SCIPcaptureVar(scip, targetdata->vars[i]) );
1779 targetdata->constnnonz = sourcedata->constnnonz;
1781 if ( sourcedata->constnnonz > 0 )
1783 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constcol), sourcedata->constcol, sourcedata->constnnonz));
1784 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constrow), sourcedata->constrow, sourcedata->constnnonz));
1785 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constval), sourcedata->constval, sourcedata->constnnonz));
1789 targetdata->constcol = NULL;
1790 targetdata->constrow = NULL;
1791 targetdata->constval = NULL;
1795 targetdata->maxrhsentry = sourcedata->maxrhsentry;
1799 snprintfreturn = SCIPsnprintf(transname, SCIP_MAXSTRLEN,
"t_%s", SCIPconsGetName(sourcecons));
1800 assert( snprintfreturn < SCIP_MAXSTRLEN );
1802 (void) SCIPsnprintf(transname, SCIP_MAXSTRLEN,
"t_%s", SCIPconsGetName(sourcecons));
1806 SCIP_CALL( SCIPcreateCons(scip, targetcons, transname, conshdlr, targetdata,
1807 SCIPconsIsInitial(sourcecons), SCIPconsIsSeparated(sourcecons), SCIPconsIsEnforced(sourcecons),
1808 SCIPconsIsChecked(sourcecons), SCIPconsIsPropagated(sourcecons), SCIPconsIsLocal(sourcecons),
1809 SCIPconsIsModifiable(sourcecons), SCIPconsIsDynamic(sourcecons), SCIPconsIsRemovable(sourcecons),
1810 SCIPconsIsStickingAtNode(sourcecons)) );
1821 assert( scip != NULL );
1822 assert( result != NULL );
1824 *result = SCIP_FEASIBLE;
1826 for (i = 0; i < nconss; ++i)
1828 SCIP_CALL(
SCIPconsSdpCheckSdpCons(scip, conss[i], sol, checkintegrality, checklprows, printreason, result) );
1829 if ( *result == SCIP_INFEASIBLE )
1845 assert( scip != NULL );
1846 assert( result != NULL );
1847 assert( conss != NULL );
1849 *result = SCIP_DIDNOTRUN;
1851 if ( objinfeasible )
1853 SCIPdebugMessage(
"-> pseudo solution is objective infeasible, return.\n");
1857 for (i = 0; i < nconss; ++i)
1861 if (*result == SCIP_INFEASIBLE)
1864 SCIPdebugMessage(
"-> pseudo solution infeasible for SDP-constraint %s, return.\n", SCIPconsGetName(conss[i]));
1869 SCIPdebugMessage(
"-> pseudo solution feasible for all SDP-constraints.\n");
1880 assert( scip != NULL );
1881 assert( conshdlr != NULL );
1882 assert( conss != NULL );
1883 assert( result != NULL );
1893 assert( scip != NULL );
1894 assert( conshdlr != NULL );
1895 assert( conss != NULL );
1896 assert( result != NULL );
1907 assert( result != NULL );
1908 *result = SCIP_DIDNOTFIND;
1910 for (i = 0; i < nusefulconss; ++i)
1912 SCIP_CALL(
separateSol(scip, conshdlr, conss[i], sol, result) );
1924 assert( result != NULL );
1925 *result = SCIP_DIDNOTFIND;
1927 for (i = 0; i < nusefulconss; ++i)
1929 SCIP_CALL(
separateSol(scip, conshdlr, conss[i], NULL, result) );
1941 assert( cons != NULL );
1942 assert( consdata != NULL );
1944 SCIPdebugMessage(
"deleting SDP constraint <%s>.\n", SCIPconsGetName(cons));
1946 for (i = 0; i < (*consdata)->nvars; i++)
1948 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->val[i], (*consdata)->nvarnonz[i]);
1949 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->row[i], (*consdata)->nvarnonz[i]);
1950 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->col[i], (*consdata)->nvarnonz[i]);
1954 for (i = 0; i < (*consdata)->nvars; i++)
1956 SCIP_CALL( SCIPreleaseVar(scip, &((*consdata)->vars[i])) );
1959 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->vars, (*consdata)->nvars);
1960 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constval, (*consdata)->constnnonz);
1961 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constrow, (*consdata)->constnnonz);
1962 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constcol, (*consdata)->constnnonz);
1963 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->val, (*consdata)->nvars);
1964 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->row, (*consdata)->nvars);
1965 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->col, (*consdata)->nvars);
1966 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->nvarnonz, (*consdata)->nvars);
1967 SCIPfreeBlockMemory(scip, consdata);
1976 SCIP_CONSHDLRDATA* conshdlrdata;
1978 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1979 assert( conshdlrdata != NULL );
1981 SCIPfreeMemory(scip, &conshdlrdata);
1982 SCIPconshdlrSetData(conshdlr, NULL);
1991 assert(scip != NULL);
1992 assert(conshdlr != NULL);
1993 assert(strcmp(SCIPconshdlrGetName(conshdlr),
CONSHDLR_NAME) == 0);
2006 SCIP_CONSDATA* sourcedata;
2008 SCIP_VAR** targetvars;
2012 assert( scip != NULL );
2013 assert( sourcescip != NULL );
2014 assert( sourcecons != NULL );
2015 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(sourcecons)),
CONSHDLR_NAME) == 0 );
2016 assert( valid != NULL );
2018 SCIPdebugMessage(
"Copying SDP constraint <%s>\n", SCIPconsGetName(sourcecons));
2024 if ( SCIPgetStage(sourcescip) <= SCIP_STAGE_EXITPRESOLVE )
2029 sourcedata = SCIPconsGetData(sourcecons);
2030 assert( sourcedata != NULL );
2032 SCIP_CALL( SCIPallocBufferArray(scip, &targetvars, sourcedata->nvars) );
2035 for (i = 0; i < sourcedata->nvars; i++)
2037 SCIP_CALL( SCIPgetVarCopy(sourcescip, scip, sourcedata->vars[i], &var, varmap, consmap, global, &success) );
2039 targetvars[i] = var;
2050 char copyname[SCIP_MAXSTRLEN];
2054 snprintfreturn = SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s",
name);
2055 assert( snprintfreturn < SCIP_MAXSTRLEN );
2057 (void) SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s",
name);
2059 SCIP_CALL(
SCIPcreateConsSdp( scip, cons, copyname, sourcedata->nvars, sourcedata->nnonz, sourcedata->blocksize, sourcedata->nvarnonz,
2060 sourcedata->col, sourcedata->row, sourcedata->val, targetvars, sourcedata->constnnonz,
2061 sourcedata->constcol, sourcedata->constrow, sourcedata->constval) );
2068 char copyname[SCIP_MAXSTRLEN];
2072 snprintfreturn = SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s", SCIPconsGetName(sourcecons));
2073 assert( snprintfreturn < SCIP_MAXSTRLEN );
2075 (void) SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s", SCIPconsGetName(sourcecons));
2077 SCIP_CALL(
SCIPcreateConsSdp( scip, cons, SCIPconsGetName(sourcecons), sourcedata->nvars, sourcedata->nnonz, sourcedata->blocksize,
2078 sourcedata->nvarnonz, sourcedata->col, sourcedata->row, sourcedata->val, targetvars, sourcedata->constnnonz,
2079 sourcedata->constcol, sourcedata->constrow, sourcedata->constval) );
2082 SCIPfreeBufferArray(scip, &targetvars);
2091 #ifdef PRINT_HUMAN_READABLE
2092 SCIP_CONSDATA* consdata;
2093 SCIP_Real* fullmatrix;
2098 assert( scip != NULL );
2099 assert( cons != NULL );
2101 consdata = SCIPconsGetData(cons);
2102 assert( consdata != NULL );
2104 SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, consdata->blocksize * consdata->blocksize) );
2107 for (v = 0; v < consdata->nvars; v++)
2112 for (i = 0; i < consdata->blocksize; i++)
2114 for (j = 0; j < consdata->blocksize; j++)
2115 fullmatrix[i * consdata->blocksize + j] = 0.0;
2119 for (i = 0; i < consdata->nvarnonz[v]; i++)
2121 fullmatrix[consdata->row[v][i] * consdata->blocksize + consdata->col[v][i]] = consdata->val[v][i];
2122 fullmatrix[consdata->col[v][i] * consdata->blocksize + consdata->row[v][i]] = consdata->val[v][i];
2126 SCIPinfoMessage(scip, file,
"+\n");
2127 for (i = 0; i < consdata->blocksize; i++)
2129 SCIPinfoMessage(scip, file,
"( ");
2130 for (j = 0; j < consdata->blocksize; j++)
2131 SCIPinfoMessage(scip, file,
"%f ", fullmatrix[i * consdata->blocksize + j]);
2132 SCIPinfoMessage(scip, file,
")\n");
2134 SCIPinfoMessage(scip, file,
"* %s\n", SCIPvarGetName(consdata->vars[v]));
2142 for (i = 0; i < consdata->blocksize; i++)
2144 for (j = 0; j < consdata->blocksize; j++)
2145 fullmatrix[i * consdata->blocksize + j] = 0.0;
2149 for (i = 0; i < consdata->constnnonz; i++)
2151 fullmatrix[consdata->constrow[i] * consdata->blocksize + consdata->constcol[i]] = consdata->constval[i];
2152 fullmatrix[consdata->constcol[i] * consdata->blocksize + consdata->constrow[i]] = consdata->constval[i];
2156 SCIPinfoMessage(scip, file,
"-\n");
2157 for (i = 0; i < consdata->blocksize; i++)
2159 SCIPinfoMessage(scip, file,
"( ");
2160 for (j = 0; j < consdata->blocksize; j++)
2161 SCIPinfoMessage(scip, file,
"%f ", fullmatrix[i * consdata->blocksize + j]);
2162 SCIPinfoMessage(scip, file,
")\n");
2164 SCIPinfoMessage(scip, file,
">=0\n");
2166 SCIPfreeBufferArray(scip, &fullmatrix);
2170 SCIP_CONSDATA* consdata;
2174 assert( scip != NULL );
2175 assert( cons != NULL );
2177 consdata = SCIPconsGetData(cons);
2180 SCIPinfoMessage(scip, file,
"%d\n", consdata->blocksize);
2183 if ( consdata->constnnonz > 0 )
2185 SCIPinfoMessage(scip, file,
"A_0: ");
2187 for (i = 0; i < consdata->constnnonz; i++)
2189 SCIPinfoMessage(scip, file,
"(%d,%d):%f, ", consdata->constrow[i], consdata->constcol[i], consdata->constval[i]);
2191 SCIPinfoMessage(scip, file,
"\n");
2195 for (v = 0; v < consdata->nvars; v++)
2197 SCIPinfoMessage(scip, file,
"<%s>: ", SCIPvarGetName(consdata->vars[v]));
2198 for (i = 0; i < consdata->nvarnonz[v]; i++)
2200 SCIPinfoMessage(scip, file,
"(%d,%d):%f, ", consdata->row[v][i], consdata->col[v][i], consdata->val[v][i]);
2203 if (v < consdata->nvars - 1)
2205 SCIPinfoMessage(scip, file,
"\n");
2217 SCIP_Bool parsesuccess;
2218 SCIP_CONSDATA* consdata = NULL;
2224 assert( scip != NULL );
2225 assert( str != NULL );
2227 nvars = SCIPgetNVars(scip);
2229 assert( success != NULL );
2233 SCIP_CALL( SCIPallocBlockMemory(scip, &consdata) );
2234 consdata->nvars = 0;
2235 consdata->nnonz = 0;
2236 consdata->constnnonz = 0;
2237 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->nvarnonz, nvars) );
2238 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col, nvars) );
2239 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row, nvars) );
2240 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val, nvars) );
2241 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->vars, nvars));
2242 consdata->constcol = NULL;
2243 consdata->constrow = NULL;
2244 consdata->constval = NULL;
2247 parsesuccess = SCIPstrToIntValue(str, &(consdata->blocksize), &pos);
2248 *success = *success && parsesuccess;
2251 while( isspace((
unsigned char)*pos) )
2255 if ( pos[0] ==
'A' && pos[1] ==
'_' && pos[2] ==
'0' )
2259 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constcol,
PARSE_STARTSIZE) );
2260 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constrow,
PARSE_STARTSIZE) );
2261 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constval,
PARSE_STARTSIZE) );
2266 while (pos[0] ==
'(')
2271 if ( consdata->constnnonz == currentsize )
2273 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constcol, currentsize,
PARSE_SIZEFACTOR * currentsize) );
2274 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constrow, currentsize,
PARSE_SIZEFACTOR * currentsize) );
2275 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constval, currentsize,
PARSE_SIZEFACTOR * currentsize) );
2279 parsesuccess = SCIPstrToIntValue(pos, &(consdata->constrow[consdata->constnnonz]), &pos);
2280 *success = *success && parsesuccess;
2281 assert( consdata->constrow[consdata->constnnonz] < consdata->blocksize );
2283 parsesuccess = SCIPstrToIntValue(pos, &(consdata->constcol[consdata->constnnonz]), &pos);
2284 *success = *success && parsesuccess;
2285 assert( consdata->constcol[consdata->constnnonz] < consdata->blocksize );
2287 parsesuccess = SCIPstrToRealValue(pos, &(consdata->constval[consdata->constnnonz]), &pos);
2288 *success = *success && parsesuccess;
2292 if ( consdata->constcol[consdata->constnnonz] > consdata->constrow[consdata->constnnonz] )
2294 i = consdata->constcol[consdata->constnnonz];
2295 consdata->constcol[consdata->constnnonz] = consdata->constrow[consdata->constnnonz];
2296 consdata->constrow[consdata->constnnonz] = i;
2299 consdata->constnnonz++;
2302 while( isspace((
unsigned char)*pos) )
2307 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constcol, currentsize, consdata->constnnonz) );
2308 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constrow, currentsize, consdata->constnnonz) );
2309 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constval, currentsize, consdata->constnnonz) );
2313 while( isspace((
unsigned char)*pos) )
2319 while (pos[0] ==
'<')
2322 SCIP_CALL( SCIPparseVarName(scip, pos, &(consdata->vars[consdata->nvars]), &pos) );
2323 SCIP_CALL( SCIPcaptureVar(scip, consdata->vars[consdata->nvars]) );
2325 consdata->nvarnonz[consdata->nvars] = 0;
2326 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->col[consdata->nvars]),
PARSE_STARTSIZE));
2327 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->row[consdata->nvars]),
PARSE_STARTSIZE));
2328 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->val[consdata->nvars]),
PARSE_STARTSIZE));
2335 while (pos[0] ==
'(')
2340 if ( consdata->nvarnonz[consdata->nvars - 1] == currentsize )
2342 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col[consdata->nvars - 1], currentsize,
PARSE_SIZEFACTOR * currentsize) );
2343 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row[consdata->nvars - 1], currentsize,
PARSE_SIZEFACTOR * currentsize) );
2344 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val[consdata->nvars - 1], currentsize,
PARSE_SIZEFACTOR * currentsize) );
2348 parsesuccess = SCIPstrToIntValue(pos, &(consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]]), &pos);
2349 *success = *success && parsesuccess;
2350 assert( consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] < consdata->blocksize );
2352 parsesuccess = SCIPstrToIntValue(pos, &(consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]]), &pos);
2353 *success = *success && parsesuccess;
2354 assert( consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] < consdata->blocksize );
2356 parsesuccess = SCIPstrToRealValue(pos, &(consdata->val[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]]), &pos);
2357 *success = *success && parsesuccess;
2361 if ( consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] >
2362 consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] )
2364 i = consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]];
2365 consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] =
2366 consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]];
2367 consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] = i;
2370 consdata->nvarnonz[consdata->nvars - 1]++;
2373 while( isspace((
unsigned char)*pos) )
2378 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col[consdata->nvars - 1], currentsize, consdata->nvarnonz[consdata->nvars - 1]) );
2379 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row[consdata->nvars - 1], currentsize, consdata->nvarnonz[consdata->nvars - 1]) );
2380 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val[consdata->nvars - 1], currentsize, consdata->nvarnonz[consdata->nvars - 1]) );
2383 while( isspace((
unsigned char)*pos) )
2388 SCIP_CALL( SCIPcreateCons(scip, cons,
name, conshdlr, consdata, initial, separate, enforce, check, propagate, local, modifiable,
2389 dynamic, removable, stickingatnode) );
2394 #ifdef SCIP_MORE_DEBUG
2395 SCIP_CALL( SCIPprintCons(scip, *cons, NULL) );
2405 SCIP_CONSDATA* consdata;
2409 assert( scip != NULL );
2410 assert( cons != NULL );
2411 assert( vars != NULL );
2412 assert( success != NULL );
2413 assert( varssize >= 0 );
2415 consdata = SCIPconsGetData(cons);
2416 assert( consdata != NULL );
2418 nvars = consdata->nvars;
2420 if ( nvars > varssize )
2422 SCIPdebugMessage(
"consGetVarsIndicator called for array of size %d, needed size %d.\n", varssize, nvars);
2427 for (i = 0; i < nvars; i++)
2428 vars[i] = consdata->vars[i];
2439 SCIP_CONSDATA* consdata;
2441 assert( scip != NULL );
2442 assert( cons != NULL );
2443 assert( nvars != NULL );
2444 assert( success != NULL );
2446 consdata = SCIPconsGetData(cons);
2447 assert( consdata != NULL );
2449 *nvars = consdata->nvars;
2460 SCIP_CONSHDLR* conshdlr = NULL;
2461 SCIP_CONSHDLRDATA* conshdlrdata = NULL;
2463 assert( scip != 0 );
2466 SCIP_CALL( SCIPallocMemory(scip, &conshdlrdata) );
2471 consEnfolpSdp, consEnfopsSdp, consCheckSdp, consLockSdp, conshdlrdata) );
2473 assert( conshdlr != NULL );
2476 SCIP_CALL( SCIPsetConshdlrDelete(scip, conshdlr, consDeleteSdp) );
2477 SCIP_CALL( SCIPsetConshdlrFree(scip, conshdlr, consFreeSdp) );
2478 SCIP_CALL( SCIPsetConshdlrCopy(scip, conshdlr, conshdlrCopySdp, consCopySdp) );
2479 SCIP_CALL( SCIPsetConshdlrInitpre(scip, conshdlr,consInitpreSdp) );
2480 SCIP_CALL( SCIPsetConshdlrExitpre(scip, conshdlr, consExitpreSdp) );
2482 SCIP_CALL( SCIPsetConshdlrSepa(scip, conshdlr, consSepalpSdp, consSepasolSdp,
CONSHDLR_SEPAFREQ,
2484 SCIP_CALL( SCIPsetConshdlrEnforelax(scip, conshdlr, consEnforelaxSdp) );
2485 SCIP_CALL( SCIPsetConshdlrTrans(scip, conshdlr, consTransSdp) );
2486 SCIP_CALL( SCIPsetConshdlrPrint(scip, conshdlr, consPrintSdp) );
2487 SCIP_CALL( SCIPsetConshdlrParse(scip, conshdlr, consParseSdp) );
2488 SCIP_CALL( SCIPsetConshdlrGetVars(scip, conshdlr, consGetVarsSdp) );
2489 SCIP_CALL( SCIPsetConshdlrGetNVars(scip, conshdlr, consGetNVarsSdp) );
2493 SCIP_CALL( SCIPaddIntParam(scip,
"constraints/SDP/threads",
"number of threads used for OpenBLAS",
2494 &(conshdlrdata->nthreads), TRUE, DEFAULT_NTHREADS, 1, INT_MAX, NULL, NULL) );
2496 SCIP_CALL( SCIPaddBoolParam(scip,
"constraints/SDP/diaggezerocuts",
2497 "Should linear cuts enforcing the non-negativity of diagonal entries of SDP-matrices be added?",
2499 SCIP_CALL( SCIPaddBoolParam(scip,
"constraints/SDP/diagzeroimplcuts",
2500 "Should linear cuts enforcing the implications of diagonal entries of zero in SDP-matrices be added?",
2532 SCIP_CONSDATA* consdata;
2536 assert( scip != NULL );
2537 assert( cons != NULL );
2538 assert( nvars != NULL );
2539 assert( nnonz != NULL );
2540 assert( blocksize != NULL );
2541 assert( arraylength != NULL );
2542 assert( nvarnonz != NULL );
2543 assert( col != NULL );
2544 assert( row != NULL );
2545 assert( val != NULL );
2546 assert( vars != NULL );
2547 assert( constnnonz != NULL );
2549 consdata = SCIPconsGetData(cons);
2550 name = SCIPconsGetName(cons);
2552 assert( consdata->constnnonz == 0 || ( constcol != NULL && constrow != NULL && constval != NULL ) );
2554 *nvars = consdata->nvars;
2555 *nnonz = consdata->nnonz;
2556 *blocksize = consdata->blocksize;
2558 for (i = 0; i < consdata->nvars; i++)
2559 vars[i] = consdata->vars[i];
2562 if ( *arraylength < consdata->nvars )
2564 SCIPdebugMessage(
"nvarnonz, col, row and val arrays were not long enough to store the information for cons %s, they need to be at least"
2565 "size %d, given was only length %d! \n", name, consdata->nvars, *arraylength);
2566 *arraylength = consdata->nvars;
2570 for (i = 0; i < consdata->nvars; i++)
2572 nvarnonz[i] = consdata->nvarnonz[i];
2574 col[i] = consdata->col[i];
2575 row[i] = consdata->row[i];
2576 val[i] = consdata->val[i];
2581 if ( consdata->constnnonz > 0 )
2583 if ( consdata->constnnonz > *constnnonz )
2585 SCIPdebugMessage(
"The constant nonzeros arrays were not long enough to store the information for cons %s, they need to be at least"
2586 "size %d, given was only length %d! \n", name, consdata->constnnonz, *constnnonz);
2590 for (i = 0; i < consdata->constnnonz; i++)
2592 constcol[i] = consdata->constcol[i];
2593 constrow[i] = consdata->constrow[i];
2594 constval[i] = consdata->constval[i];
2599 *constnnonz = consdata->constnnonz;
2615 SCIP_CONSDATA* consdata;
2617 assert( scip != NULL );
2618 assert( cons != NULL );
2620 consdata = SCIPconsGetData(cons);
2621 assert( consdata != NULL );
2623 if ( nnonz != NULL )
2624 *nnonz = consdata->nnonz;
2626 if ( constnnonz != NULL )
2627 *constnnonz = consdata->constnnonz;
2638 SCIP_CONSDATA* consdata;
2640 assert( scip != NULL );
2641 assert( cons != NULL );
2643 consdata = SCIPconsGetData(cons);
2644 assert( consdata != NULL );
2646 return consdata->blocksize;
2657 SCIP_CONSDATA* consdata;
2661 assert( scip != NULL );
2662 assert( cons != NULL );
2664 assert( Aj != NULL );
2666 consdata = SCIPconsGetData(cons);
2667 assert( consdata != NULL );
2668 blocksize = consdata->blocksize;
2670 assert( j < consdata->nvars );
2672 for (i = 0; i < blocksize * blocksize; i++)
2675 for (i = 0; i < consdata->nvarnonz[j]; i++)
2677 Aj[consdata->col[j][i] * blocksize + consdata->row[j][i]] = consdata->val[j][i];
2678 Aj[consdata->row[j][i] * blocksize + consdata->col[j][i]] = consdata->val[j][i];
2691 SCIP_CONSDATA* consdata;
2696 assert( scip != NULL );
2697 assert( cons != NULL );
2698 assert( mat != NULL );
2700 consdata = SCIPconsGetData(cons);
2701 blocksize = consdata->blocksize;
2703 for (i = 0; i < blocksize; i++)
2705 for (j = 0; j < blocksize; j++)
2706 mat[i * blocksize + j] = 0.0;
2709 for (i = 0; i < consdata->constnnonz; i++)
2711 mat[consdata->constcol[i] * blocksize + consdata->constrow[i]] = consdata->constval[i];
2712 mat[consdata->constrow[i] * blocksize + consdata->constcol[i]] = consdata->constval[i];
2725 SCIP_CONSDATA* consdata;
2729 assert( scip != NULL );
2730 assert( cons != NULL );
2731 assert( mat != NULL );
2733 consdata = SCIPconsGetData(cons);
2734 assert( consdata != NULL );
2736 blocksize = consdata->blocksize;
2739 for (i = 0; i < (blocksize * (blocksize + 1)) / 2; i++)
2742 for (i = 0; i < consdata->constnnonz; i++)
2743 mat[
compLowerTriangPos(consdata->constrow[i], consdata->constcol[i])] = consdata->constval[i];
2761 SCIP_Real* lambdastar
2764 SCIP_CONSDATA* consdata;
2766 SCIP_Real maxinfnorm;
2768 SCIP_Real mininfnorm;
2771 SCIP_Real primalguess;
2772 SCIP_Real dualguess;
2778 assert( scip != NULL );
2779 assert( cons != NULL );
2780 assert( lambdastar != NULL );
2781 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(cons)),
CONSHDLR_NAME) == 0 );
2783 consdata = SCIPconsGetData(cons);
2784 assert( consdata != NULL );
2790 if ( consdata->nnonz == 0 )
2792 *lambdastar = 100.0;
2796 blocksize = consdata->blocksize;
2798 sparsity = consdata->nnonz / (0.5 * blocksize * (blocksize + 1));
2802 mininfnorm = SCIPinfinity(scip);
2803 for (v = 0; v < consdata->nvars; v++)
2805 for (i = 0; i < consdata->nvarnonz[v]; i++)
2807 if ( SCIPisGT(scip, REALABS(consdata->val[v][i]), maxinfnorm ) )
2808 maxinfnorm = REALABS(consdata->val[v][i]);
2809 if ( SCIPisLT(scip, REALABS(consdata->val[v][i]), mininfnorm) )
2810 mininfnorm = REALABS(consdata->val[v][i]);
2814 for (i = 0; i < consdata->constnnonz; i++)
2816 if ( SCIPisGT(scip, REALABS(consdata->constval[i]), maxconst ) )
2817 maxconst = REALABS(consdata->constval[i]);
2820 assert( SCIPisGT(scip, mininfnorm, 0.0) );
2825 for (v = 0; v < consdata->nvars; v++)
2827 if ( SCIPisGT(scip, REALABS(SCIPvarGetObj(consdata->vars[v])), maxobj) )
2828 maxobj = REALABS(SCIPvarGetObj(consdata->vars[v]));
2829 compval = SCIPisInfinity(scip, REALABS(SCIPvarGetUbGlobal(consdata->vars[v]))) ? 1e+6 : REALABS(SCIPvarGetUbGlobal(consdata->vars[v]));
2830 if ( SCIPisGT(scip, compval, maxbound) )
2832 compval = SCIPisInfinity(scip, REALABS(SCIPvarGetLbGlobal(consdata->vars[v]))) ? 1e+6 : REALABS(SCIPvarGetUbGlobal(consdata->vars[v]));
2833 if ( SCIPisGT(scip, compval, maxbound) )
2838 if ( SCIPisEQ(scip, maxbound, 0.0) )
2842 primalguess = maxobj / (sparsity * mininfnorm);
2843 dualguess = sparsity * maxinfnorm * maxbound + maxconst;
2845 if ( SCIPisGT(scip, primalguess, dualguess) )
2846 *lambdastar = primalguess;
2848 *lambdastar = dualguess;
2872 SCIP_CONSHDLR* conshdlr;
2873 SCIP_CONSDATA* consdata = NULL;
2877 assert( scip != NULL );
2878 assert( cons != NULL );
2879 assert( name != NULL );
2880 assert( nvars >= 0 );
2881 assert( nnonz >= 0 );
2882 assert( blocksize >= 0 );
2883 assert( constnnonz >= 0 );
2884 assert( nvars == 0 || vars != NULL );
2885 assert( nnonz == 0 || (nvarnonz != NULL && col != NULL && row != NULL && val != NULL ));
2886 assert( constnnonz == 0 || (constcol != NULL && constrow != NULL && constval != NULL ));
2888 conshdlr = SCIPfindConshdlr(scip,
"SDP");
2889 if ( conshdlr == NULL )
2891 SCIPerrorMessage(
"SDP constraint handler not found\n");
2892 return SCIP_PLUGINNOTFOUND;
2896 SCIP_CALL( SCIPallocBlockMemory(scip, &consdata) );
2897 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->nvarnonz, nvars) );
2898 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col, nvars) );
2899 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row, nvars) );
2900 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val, nvars) );
2901 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constcol, constnnonz) );
2902 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constrow, constnnonz) );
2903 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constval, constnnonz) );
2904 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->vars, nvars));
2906 for (i = 0; i < nvars; i++)
2908 assert( nvarnonz[i] >= 0 );
2910 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col[i], nvarnonz[i]));
2911 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row[i], nvarnonz[i]));
2912 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val[i], nvarnonz[i]));
2915 consdata->nvars = nvars;
2916 consdata->nnonz = nnonz;
2917 consdata->constnnonz = constnnonz;
2918 consdata->blocksize = blocksize;
2920 for (i = 0; i < nvars; i++)
2922 consdata->nvarnonz[i] = nvarnonz[i];
2924 for (j = 0; j < nvarnonz[i]; j++)
2926 assert( col[i][j] >= 0 );
2927 assert( col[i][j] < blocksize );
2928 assert( row[i][j] >= 0 );
2929 assert( row[i][j] < blocksize );
2931 consdata->col[i][j] = col[i][j];
2932 consdata->row[i][j] = row[i][j];
2933 consdata->val[i][j] = val[i][j];
2937 for (i = 0; i < constnnonz; i++)
2939 consdata->constcol[i] = constcol[i];
2940 consdata->constrow[i] = constrow[i];
2941 consdata->constval[i] = constval[i];
2944 for (i = 0; i < nvars; i++)
2946 consdata->vars[i] = vars[i];
2947 SCIP_CALL( SCIPcaptureVar(scip, consdata->vars[i]) );
2949 SCIPdebugMessage(
"creating cons %s\n", name);
2951 SCIP_CALL( SCIPcreateCons(scip, cons, name, conshdlr, consdata, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE) );
static SCIP_DECL_CONSTRANS(consTransSdp)
EXTERN SCIP_RETCODE SCIPlapackMatrixVectorMult(int nrows, int ncols, SCIP_Real *matrix, SCIP_Real *vector, SCIP_Real *result)
EXTERN SCIP_RETCODE SCIPlapackComputeIthEigenvalue(BMS_BUFMEM *bufmem, SCIP_Bool geteigenvectors, int n, SCIP_Real *A, int i, SCIP_Real *eigenvalue, SCIP_Real *eigenvector)
static SCIP_RETCODE separateSol(SCIP *scip, SCIP_CONSHDLR *conshdlr, SCIP_CONS *cons, SCIP_SOL *sol, SCIP_RESULT *result)
SCIP_RETCODE SCIPconsSdpGetFullAj(SCIP *scip, SCIP_CONS *cons, int j, SCIP_Real *Aj)
void SCIPsdpVarfixerSortRowCol(int *row, int *col, SCIP_Real *val, int length)
#define CONSHDLR_PRESOLTIMING
SCIP_RETCODE SCIPincludeConshdlrSdp(SCIP *scip)
static SCIP_RETCODE multiplyConstraintMatrix(SCIP_CONS *cons, int j, SCIP_Real *v, SCIP_Real *vAv)
static SCIP_DECL_CONSENFORELAX(consEnforelaxSdp)
#define CONSHDLR_NEEDSCONS
static SCIP_RETCODE diagGEzero(SCIP *scip, SCIP_CONS **conss, int nconss, int *naddconss)
static SCIP_DECL_CONSDELETE(consDeleteSdp)
#define CONSHDLR_SEPAFREQ
#define CONSHDLR_DELAYSEPA
static SCIP_RETCODE computeSdpMatrix(SCIP *scip, SCIP_CONS *cons, SCIP_SOL *y, SCIP_Real *matrix)
#define CONSHDLR_EAGERFREQ
static SCIP_RETCODE move_1x1_blocks_to_lp(SCIP *scip, SCIP_CONS **conss, int nconss, int *naddconss, int *ndelconss, SCIP_RESULT *result)
#define DEFAULT_DIAGGEZEROCUTS
SCIP_RETCODE SCIPconsSdpGetLowerTriangConstMatrix(SCIP *scip, SCIP_CONS *cons, SCIP_Real *mat)
static SCIP_DECL_CONSSEPALP(consSepalpSdp)
static SCIP_RETCODE setMaxRhsEntry(SCIP_CONS *cons)
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)
static SCIP_DECL_CONSCHECK(consCheckSdp)
Constraint handler for SDP-constraints.
static SCIP_RETCODE expandSymMatrix(int size, SCIP_Real *symMat, SCIP_Real *fullMat)
class that maps SCIP variables to SDP indices (the SCIP variables are given SDP indices in the order ...
SCIP_RETCODE SCIPconsSdpGetFullConstMatrix(SCIP *scip, SCIP_CONS *cons, SCIP_Real *mat)
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)
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)
#define CONSHDLR_CHECKPRIORITY
static SCIP_DECL_CONSENFOPS(consEnfopsSdp)
#define DEFAULT_DIAGZEROIMPLCUTS
static int compLowerTriangPos(int i, int j)
adds the main functionality to fix/unfix/(multi-)aggregate variables by merging two three-tuple-array...
#define CONSHDLR_ENFOPRIORITY
static SCIP_DECL_CONSFREE(consFreeSdp)
static SCIP_DECL_CONSPRESOL(consPresolSdp)
static SCIP_DECL_CONSEXITPRE(consExitpreSdp)
static SCIP_DECL_CONSPRINT(consPrintSdp)
static SCIP_DECL_CONSHDLRCOPY(conshdlrCopySdp)
#define CONSHDLR_MAXPREROUNDS
static SCIP_DECL_CONSINITPRE(consInitpreSdp)
static SCIP_DECL_CONSPARSE(consParseSdp)
SCIP_RETCODE SCIPconsSdpCheckSdpCons(SCIP *scip, SCIP_CONS *cons, SCIP_SOL *sol, SCIP_Bool checkintegrality, SCIP_Bool checklprows, SCIP_Bool printreason, SCIP_RESULT *result)
static SCIP_RETCODE diagZeroImpl(SCIP *scip, SCIP_CONS **conss, int nconss, int *naddconss)
int SCIPconsSdpGetBlocksize(SCIP *scip, SCIP_CONS *cons)
static SCIP_DECL_CONSSEPASOL(consSepasolSdp)
SCIP_RETCODE SCIPsdpVarfixerMergeArrays(BMS_BLKMEM *blkmem, SCIP_Real epsilon, 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)
static SCIP_RETCODE EnforceConstraint(SCIP *scip, SCIP_CONSHDLR *conshdlr, SCIP_CONS **conss, int nconss, SCIP_SOL *sol, SCIP_RESULT *result)
static SCIP_DECL_CONSGETVARS(consGetVarsSdp)
#define CONSHDLR_SEPAPRIORITY
SCIP_RETCODE SCIPconsSdpGetNNonz(SCIP *scip, SCIP_CONS *cons, int *nnonz, int *constnnonz)
static SCIP_RETCODE fixAndAggrVars(SCIP *scip, SCIP_CONS **conss, int nconss, SCIP_Bool aggregate)
static SCIP_DECL_CONSLOCK(consLockSdp)
char name[SCIP_MAXSTRLEN]
static SCIP_DECL_CONSGETNVARS(consGetNVarsSdp)
static SCIP_DECL_CONSCOPY(consCopySdp)
interface methods for eigenvector computation and matrix multiplication using different versions of L...
SCIP_RETCODE SCIPconsSdpGuessInitialPoint(SCIP *scip, SCIP_CONS *cons, SCIP_Real *lambdastar)
static SCIP_RETCODE cutUsingEigenvector(SCIP *scip, SCIP_CONS *cons, SCIP_SOL *sol, SCIP_Real *coeff, SCIP_Real *lhs)
static SCIP_DECL_CONSENFOLP(consEnfolpSdp)