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 assert( symMat != NULL );
135 assert( fullMat != NULL );
138 for (i = 0; i < size; i++)
140 for (j = 0; j <= i; j++)
143 fullMat[i*size + j] = symMat[ind];
144 fullMat[j*size + i] = symMat[ind];
163 SCIP_CONSDATA* consdata;
170 assert( cons != NULL );
171 assert( matrix != NULL );
173 consdata = SCIPconsGetData(cons);
174 nvars = consdata->nvars;
175 blocksize = consdata->blocksize;
178 for (i = 0; i < (blocksize * (blocksize + 1))/2; i++)
182 for (i = 0; i < nvars; i++)
184 yval = SCIPgetSolVal(scip, y, consdata->vars[i]);
185 for (ind = 0; ind < consdata->nvarnonz[i]; ind++)
190 for (ind = 0; ind < consdata->constnnonz; ind++)
205 SCIP_CONSDATA* consdata;
208 assert( cons != NULL );
210 assert( vAv != NULL );
212 consdata = SCIPconsGetData(cons);
213 assert( consdata != NULL );
215 assert( j < consdata->nvars );
220 for (i = 0; i < consdata->nvarnonz[j]; i++)
222 if (consdata->col[j][i] == consdata->row[j][i])
223 *vAv += v[consdata->col[j][i]] * consdata->val[j][i] * v[consdata->row[j][i]];
227 *vAv += 2.0 * v[consdata->col[j][i]] * consdata->val[j][i] * v[consdata->row[j][i]];
242 SCIP_CONSDATA* consdata;
246 consdata = SCIPconsGetData(cons);
247 assert( consdata != NULL );
253 for (i = 0; i < consdata->constnnonz; i++)
255 if ( REALABS(consdata->constval[i]) > max )
256 max = REALABS(consdata->constval[i]);
259 consdata->maxrhsentry = max;
279 SCIP_CONSDATA* consdata;
280 SCIP_Real* eigenvalues = NULL;
281 SCIP_Real* matrix = NULL;
282 SCIP_Real* fullmatrix = NULL;
283 SCIP_Real* fullconstmatrix = NULL;
284 SCIP_Real* eigenvector = NULL;
285 SCIP_Real* output_vector = NULL;
289 assert( cons != NULL );
290 assert( lhs != NULL );
294 consdata = SCIPconsGetData(cons);
295 assert( consdata != NULL );
297 blocksize = consdata->blocksize;
299 SCIP_CALL( SCIPallocBufferArray(scip, &eigenvalues, blocksize) );
300 SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize+1))/2 ) );
301 SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, blocksize * blocksize ) );
302 SCIP_CALL( SCIPallocBufferArray(scip, &fullconstmatrix, blocksize * blocksize) );
303 SCIP_CALL( SCIPallocBufferArray(scip, &eigenvector, blocksize) );
304 SCIP_CALL( SCIPallocBufferArray(scip, &output_vector, blocksize) );
320 for (j = 0; j < blocksize; ++j)
321 *lhs += eigenvector[j] * output_vector[j];
324 for (j = 0; j < consdata->nvars; ++j)
329 SCIPfreeBufferArray(scip, &output_vector);
330 SCIPfreeBufferArray(scip, &eigenvector);
331 SCIPfreeBufferArray(scip, &fullconstmatrix);
332 SCIPfreeBufferArray(scip, &fullmatrix);
333 SCIPfreeBufferArray(scip, &matrix);
334 SCIPfreeBufferArray(scip, &eigenvalues);
344 SCIP_Bool checkintegrality,
345 SCIP_Bool checklprows,
346 SCIP_Bool printreason,
350 SCIP_CONSDATA* consdata;
352 SCIP_Real check_value;
353 SCIP_Real eigenvalue;
354 SCIP_Real* matrix = NULL;
355 SCIP_Real* fullmatrix = NULL;
357 assert( scip != NULL );
358 assert( cons != NULL );
359 assert( result != NULL );
360 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(cons)),
"SDP") == 0);
362 consdata = SCIPconsGetData(cons);
363 assert( consdata != NULL );
364 blocksize = consdata->blocksize;
366 SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize+1)) / 2) );
367 SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, blocksize * blocksize) );
376 check_value = (-eigenvalue) / (1.0 + consdata->maxrhsentry);
378 check_value = -eigenvalue;
381 if ( SCIPisFeasLE(scip, check_value, 0.0) )
382 *result = SCIP_FEASIBLE;
385 *result = SCIP_INFEASIBLE;
389 SCIP_CALL( SCIPprintCons(scip, cons, NULL) );
390 SCIP_CALL( SCIPprintSol(scip, sol, NULL, FALSE) );
391 SCIPinfoMessage(scip, NULL,
"non psd matrix (eigenvalue %f, dimacs error norm = %f).\n", eigenvalue, check_value);
397 SCIPupdateSolConsViolation(scip, sol, -eigenvalue, (-eigenvalue) / (1.0 + consdata->maxrhsentry));
399 SCIPfreeBufferArray(scip, &fullmatrix);
400 SCIPfreeBufferArray(scip, &matrix);
409 SCIP_CONSHDLR* conshdlr,
415 char cutname[SCIP_MAXSTRLEN];
416 SCIP_CONSDATA* consdata;
417 SCIP_CONSHDLRDATA* conshdlrdata;
419 SCIP_Real* coeff = NULL;
430 assert( cons != NULL );
432 consdata = SCIPconsGetData(cons);
433 assert( consdata != NULL );
435 nvars = consdata->nvars;
436 SCIP_CALL( SCIPallocBufferArray(scip, &coeff, nvars ) );
440 SCIP_CALL( SCIPallocBufferArray(scip, &cols, nvars ) );
441 SCIP_CALL( SCIPallocBufferArray(scip, &vals, nvars ) );
444 for (j = 0; j < nvars; ++j)
446 if ( SCIPisZero(scip, coeff[j]) )
449 cols[len] = SCIPvarGetCol(SCIPvarGetProbvar(consdata->vars[j]));
450 vals[len] = coeff[j];
454 conshdlrdata = SCIPconshdlrGetData(conshdlr);
455 assert( conshdlrdata != NULL );
458 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
459 assert( snprintfreturn < SCIP_MAXSTRLEN );
461 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
463 SCIP_CALL( SCIPcreateRowCons(scip, &row, conshdlr, cutname , len, cols, vals, lhs, SCIPinfinity(scip), FALSE, FALSE, TRUE) );
465 if ( SCIPisCutEfficacious(scip, sol, row) )
467 SCIP_Bool infeasible;
468 #ifdef SCIP_MORE_DEBUG
469 SCIPinfoMessage(scip, NULL,
"Added cut %s: ", cutname);
470 SCIPinfoMessage(scip, NULL,
"%f <= ", lhs);
471 for (j = 0; j < len; j++)
472 SCIPinfoMessage(scip, NULL,
"+ (%f)*%s", vals[j], SCIPvarGetName(SCIPcolGetVar(cols[j])));
473 SCIPinfoMessage(scip, NULL,
"\n");
476 SCIP_CALL( SCIPaddRow(scip, row, FALSE, &infeasible) );
478 *result = SCIP_CUTOFF;
480 *result = SCIP_SEPARATED;
481 SCIP_CALL( SCIPresetConsAge(scip, cons) );
484 SCIP_CALL( SCIPreleaseRow(scip, &row) );
486 SCIPfreeBufferArray(scip, &vals);
487 SCIPfreeBufferArray(scip, &cols);
488 SCIPfreeBufferArray(scip, &coeff);
504 char cutname[SCIP_MAXSTRLEN];
505 SCIP_CONSDATA* consdata;
507 SCIP_Real* cons_array;
508 SCIP_Real* lhs_array;
520 for (c = 0; c < nconss; ++c)
522 SCIP_CONSHDLR* conshdlr;
524 conshdlr = SCIPconsGetHdlr(conss[c]);
525 assert( conshdlr != NULL );
527 assert( strcmp(SCIPconshdlrGetName(conshdlr),
"SDP") == 0);
530 consdata = SCIPconsGetData(conss[c]);
531 assert( consdata != NULL );
533 blocksize = consdata->blocksize;
534 nvars = consdata->nvars;
535 rhs = SCIPinfinity(scip);
537 SCIP_CALL( SCIPallocBufferArray(scip, &matrix, (blocksize * (blocksize + 1)) / 2) );
540 SCIP_CALL( SCIPallocBufferArray(scip, &cons_array, blocksize * nvars) );
541 SCIP_CALL( SCIPallocBufferArray(scip, &lhs_array, blocksize) );
544 for (k = 0; k < blocksize; ++k)
548 for (j = 0; j < nvars; ++j)
550 for (k = 0; k < blocksize; ++k)
553 cons_array[k * nvars + j] = 0.0;
557 for (i = 0; i < consdata->nvarnonz[j]; i++)
559 if ( consdata->col[j][i] == consdata->row[j][i] )
560 cons_array[consdata->col[j][i] * nvars + j] = consdata->val[j][i];
565 for (k = 0; k < blocksize; ++k)
568 SCIP_CONSHDLRDATA* conshdlrdata;
570 conshdlrdata = SCIPconshdlrGetData(conshdlr);
572 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_ge_zero_%d", ++(conshdlrdata->ndiaggezerocuts));
573 assert( snprintfreturn < SCIP_MAXSTRLEN );
575 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_ge_zero_%d", ++(conshdlrdata->ndiaggezerocuts));
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++)
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]));
586 SCIPinfoMessage(scip, NULL,
"\n");
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) );
592 SCIP_CALL( SCIPaddCons(scip, cons) );
593 SCIP_CALL( SCIPreleaseCons(scip, &cons) );
597 SCIPfreeBufferArray(scip, &lhs_array);
598 SCIPfreeBufferArray(scip, &cons_array);
599 SCIPfreeBufferArray(scip, &matrix);
620 char cutname[SCIP_MAXSTRLEN];
625 int* nconstnonzeroentries;
626 int** constnonzeroentries;
627 SCIP_CONSDATA* consdata;
640 SCIP_Bool anycutvalid;
645 assert( scip != NULL );
646 assert( conss != NULL );
647 assert( nconss >= 0 );
648 assert( naddconss != NULL );
650 for (i = 0; i < nconss; ++i)
652 assert( conss[i] != NULL );
653 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[i])),
"SDP") == 0 );
655 consdata = SCIPconsGetData(conss[i]);
656 assert( consdata != NULL );
658 blocksize = consdata->blocksize;
659 nvars = consdata->nvars;
660 SCIP_CALL( SCIPallocBufferArray(scip, &nconstnonzeroentries, blocksize) );
661 SCIP_CALL( SCIPallocBufferArray(scip, &constnonzeroentries, blocksize) );
662 for (j = 0; j < blocksize; j++)
664 nconstnonzeroentries[j] = 0;
665 SCIP_CALL( SCIPallocBufferArray(scip, &constnonzeroentries[j], blocksize) );
669 for (j = 0; j < consdata->constnnonz; j++)
673 if ( (consdata->constcol[j] != consdata->constrow[j]) )
675 assert( ! SCIPisZero(scip, consdata->constval[j]) );
676 if ( nconstnonzeroentries[consdata->constcol[j]] >= 0 )
678 constnonzeroentries[consdata->constcol[j]][nconstnonzeroentries[consdata->constcol[j]]] = consdata->constrow[j];
679 nconstnonzeroentries[consdata->constcol[j]]++;
681 if ( nconstnonzeroentries[consdata->constrow[j]] >= 0 )
683 constnonzeroentries[consdata->constrow[j]][nconstnonzeroentries[consdata->constrow[j]]] = consdata->constcol[j];
684 nconstnonzeroentries[consdata->constrow[j]]++;
690 assert( ! SCIPisZero(scip, consdata->constval[j]) );
691 nconstnonzeroentries[consdata->constcol[j]] = -1;
697 SCIP_CALL( SCIPallocBufferArray(scip, &diagvars, blocksize) );
698 SCIP_CALL( SCIPallocBufferArray(scip, &ndiagvars, blocksize) );
700 for (j = 0; j < blocksize; ++j)
703 if ( nconstnonzeroentries[j] > 0 )
705 SCIP_CALL( SCIPallocBufferArray(scip, &(diagvars[j]), nvars) );
713 SCIPfreeBufferArray(scip, &ndiagvars);
714 SCIPfreeBufferArray(scip, &diagvars);
715 for (j = blocksize - 1; j >= 0; j--)
717 SCIPfreeBufferArray(scip, &constnonzeroentries[j]);
719 SCIPfreeBufferArray(scip, &constnonzeroentries);
720 SCIPfreeBufferArray(scip, &nconstnonzeroentries);
726 for (v = 0; v < nvars; v++)
728 for (j = 0; j < consdata->nvarnonz[v]; j++)
732 if ( (consdata->col[v][j] == consdata->row[v][j]) && (nconstnonzeroentries[consdata->col[v][j]] > 0) )
734 if ( SCIPvarIsIntegral(consdata->vars[v]) )
736 assert( ! SCIPisEQ(scip, consdata->val[v][j], 0.0) );
737 diagvars[consdata->col[v][j]][ndiagvars[consdata->col[v][j]]] = v;
738 ndiagvars[consdata->col[v][j]]++;
742 nconstnonzeroentries[consdata->col[v][j]] = -2;
747 else if ( consdata->col[v][j] != consdata->row[v][j] )
749 if ( nconstnonzeroentries[consdata->col[v][j]] > 0 )
752 for (k = 0; k < nconstnonzeroentries[consdata->col[v][j]]; k++)
754 if ( constnonzeroentries[consdata->col[v][j]][k] == consdata->row[v][j] )
757 if ( nconstnonzeroentries[consdata->col[v][j]] > k + 1 )
759 for (l = k + 1; l < nconstnonzeroentries[consdata->col[v][j]]; l++)
760 constnonzeroentries[consdata->col[v][j]][l - 1] = constnonzeroentries[consdata->col[v][j]][l];
762 nconstnonzeroentries[consdata->col[v][j]]--;
764 if ( nconstnonzeroentries[consdata->col[v][j]] == 0 )
765 nconstnonzeroentries[consdata->col[v][j]] = -2;
771 if ( nconstnonzeroentries[consdata->row[v][j]] > 0 )
774 for (k = 0; k < nconstnonzeroentries[consdata->row[v][j]]; k++)
776 if ( constnonzeroentries[consdata->row[v][j]][k] == consdata->col[v][j] )
779 if ( nconstnonzeroentries[consdata->row[v][j]] > k + 1 )
781 for (l = k + 1; l < nconstnonzeroentries[consdata->row[v][j]]; l++)
782 constnonzeroentries[consdata->row[v][j]][l - 1] = constnonzeroentries[consdata->row[v][j]][l];
784 nconstnonzeroentries[consdata->row[v][j]]--;
786 if ( nconstnonzeroentries[consdata->row[v][j]] == 0 )
787 nconstnonzeroentries[consdata->row[v][j]] = -2;
796 for (j = 0; j < blocksize; ++j)
798 if ( nconstnonzeroentries[j] > 0 )
800 SCIP_CALL( SCIPallocBufferArray(scip, &vals, ndiagvars[j]) );
801 SCIP_CALL( SCIPallocBufferArray(scip, &vars, ndiagvars[j]) );
804 for (v = 0; v < ndiagvars[j]; ++v)
806 vars[v] = consdata->vars[diagvars[j][v]];
810 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_zero_impl_block_%d_row_%d", i, j);
811 assert( snprintfreturn < SCIP_MAXSTRLEN );
813 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"diag_zero_impl_block_%d_row_%d", i, j);
816 #ifdef SCIP_MORE_DEBUG
817 SCIPinfoMessage(scip, NULL,
"Added lp-constraint %s: ", cutname);
818 SCIPinfoMessage(scip, NULL,
"1 <=");
819 for (l = 0; l < ndiagvars[j]; l++)
820 SCIPinfoMessage(scip, NULL,
" + (%f)*%s", vals[l], SCIPvarGetName(vars[l]));
821 SCIPinfoMessage(scip, NULL,
"\n");
825 SCIP_CALL(SCIPcreateConsLinear(scip, &cons, cutname , ndiagvars[j], vars, vals, 1.0, SCIPinfinity(scip), TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE));
826 SCIP_CALL(SCIPaddCons(scip, cons));
827 SCIP_CALL(SCIPreleaseCons(scip, &cons));
830 SCIPfreeBufferArray(scip, &vars);
831 SCIPfreeBufferArray(scip, &vals);
833 if ( nconstnonzeroentries[j] == -2 || nconstnonzeroentries[j] > 0 )
835 SCIPfreeBufferArray(scip, &diagvars[j]);
839 SCIPfreeBufferArray(scip, &ndiagvars);
840 SCIPfreeBufferArray(scip, &diagvars);
841 for (j = blocksize - 1; j >= 0; j--)
843 SCIPfreeBufferArray(scip, &constnonzeroentries[j]);
845 SCIPfreeBufferArray(scip, &constnonzeroentries);
846 SCIPfreeBufferArray(scip, &nconstnonzeroentries);
863 char cutname[SCIP_MAXSTRLEN];
865 SCIP_CONSDATA* consdata;
866 SCIP_CONSHDLRDATA* conshdlrdata;
881 assert( result != NULL );
882 *result = SCIP_SUCCESS;
884 for (i = 0; i < nconss; ++i)
886 hdlr = SCIPconsGetHdlr(conss[i]);
887 assert(hdlr != NULL);
890 assert( strcmp(SCIPconshdlrGetName(hdlr),
"SDP") == 0);
893 consdata = SCIPconsGetData(conss[i]);
894 assert( consdata != NULL );
897 if ( consdata->blocksize == 1 )
899 nvars = consdata->nvars;
900 nnonz = consdata->nnonz;
901 SCIP_CALL( SCIPallocBufferArray(scip, &vars, nvars) );
902 SCIP_CALL( SCIPallocBufferArray(scip, &coeffs, nnonz) );
907 for (var = 0; var < nvars; var++)
909 assert( consdata->nvarnonz[var] <= 1 );
911 for (j = 0; j < consdata->nvarnonz[var]; j++)
913 assert( consdata->col[var][j] == 0 && consdata->row[var][j] == 0 );
914 coeffs[count] = consdata->val[var][j];
915 vars[count] = consdata->vars[var];
921 assert( consdata->constnnonz <= 1 );
923 rhs = (consdata->constnnonz == 1) ? consdata->constval[0] : 0.0;
929 conshdlrdata = SCIPconshdlrGetData(hdlr);
931 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"1x1block_%d", ++(conshdlrdata->n1x1blocks));
932 assert( snprintfreturn < SCIP_MAXSTRLEN );
934 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"1x1block_%d", ++(conshdlrdata->n1x1blocks));
937 #ifdef SCIP_MORE_DEBUG
938 SCIPinfoMessage(scip, NULL,
"Added lp-constraint %s: ", cutname);
939 for (i = 0; i < consdata->nvars; i++)
940 SCIPinfoMessage(scip, NULL,
"+ (%f)*%s", coeffs[i], SCIPvarGetName(vars[i]));
941 SCIPinfoMessage(scip, NULL,
" <= %f \n", rhs);
944 SCIP_CALL( SCIPcreateConsLinear(scip, &cons, cutname, consdata->nvars, vars, coeffs, rhs, SCIPinfinity(scip),
945 TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
947 SCIP_CALL( SCIPaddCons(scip, cons) );
948 SCIP_CALL( SCIPreleaseCons(scip, &cons) );
956 if ( coeffs[0] > 0.0 )
959 if ( SCIPisFeasGT(scip, rhs / coeffs[0], SCIPvarGetLbLocal(vars[0])) )
962 if( SCIPisFeasGT(scip, rhs / coeffs[0], SCIPvarGetUbLocal(vars[0])) )
964 SCIPdebugMessage(
"Problem detected to be infeasible during presolving, 1x1-SDP-constraint %s caused change"
965 "of lower bound for variable %s from %f to %f, which is bigger than upper bound of %f\n",
966 SCIPconsGetName(conss[i]), SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]), rhs / coeffs[0],
967 SCIPvarGetUbLocal(vars[0]));
969 *result = SCIP_CUTOFF;
972 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
975 SCIPfreeBufferArray(scip, &coeffs);
976 SCIPfreeBufferArray(scip, &vars);
981 SCIPdebugMessage(
"Changing lower bound of variable %s from %f to %f because of 1x1-SDP-constraint %s!\n",
982 SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]), rhs / coeffs[0], SCIPconsGetName(conss[i]));
983 SCIP_CALL( SCIPchgVarLb(scip, vars[0], rhs / coeffs[0]) );
987 SCIPdebugMessage(
"Deleting 1x1-SDP-constraint %s, new lower bound %f for variable %s no improvement over old bound %f!\n",
988 SCIPconsGetName(conss[i]), rhs / coeffs[0], SCIPvarGetName(vars[0]), SCIPvarGetLbLocal(vars[0]));
991 else if ( coeffs[0] < 0.0 )
994 if (SCIPisFeasLT(scip, rhs / coeffs[0], SCIPvarGetUbLocal(vars[0])))
997 if( SCIPisFeasLT(scip, rhs / coeffs[0], SCIPvarGetLbLocal(vars[0])) )
999 SCIPdebugMessage(
"Problem detected to be infeasible during presolving, 1x1-SDP-constraint %s caused change"
1000 "of upper bound for variable %s from %f to %f, which is less than lower bound of %f\n",
1001 SCIPconsGetName(conss[i]), SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]), rhs / coeffs[0],
1002 SCIPvarGetLbLocal(vars[0]));
1004 *result = SCIP_CUTOFF;
1007 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
1010 SCIPfreeBufferArray(scip, &coeffs);
1011 SCIPfreeBufferArray(scip, &vars);
1016 SCIPdebugMessage(
"Changing upper bound of variable %s from %f to %f because of 1x1-SDP-constraint %s!\n",
1017 SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]), -rhs / coeffs[0], SCIPconsGetName(conss[i]));
1018 SCIP_CALL( SCIPchgVarUb(scip, vars[0], rhs / coeffs[0]) );
1022 SCIPdebugMessage(
"Deleting 1x1-SDP-constraint %s, new upper bound %f for variable %s no improvement over old bound %f!\n",
1023 SCIPconsGetName(conss[i]), rhs / coeffs[0], SCIPvarGetName(vars[0]), SCIPvarGetUbLocal(vars[0]));
1028 SCIPdebugMessage(
"Detected 1x1 SDP-block without any nonzero coefficients \n");
1029 if ( SCIPisFeasGT(scip, rhs, 0.0) )
1031 SCIPdebugMessage(
"Detected infeasibility in 1x1 SDP-block without any nonzero coefficients but with strictly positive rhs\n");
1032 *result = SCIP_CUTOFF;
1035 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
1038 SCIPfreeBufferArray(scip, &coeffs);
1039 SCIPfreeBufferArray(scip, &vars);
1047 SCIP_CALL( SCIPdelCons(scip, conss[i]) );
1050 SCIPfreeBufferArray(scip, &coeffs);
1051 SCIPfreeBufferArray(scip, &vars);
1063 SCIP_VAR** aggrvars,
1069 SCIP_Real* savedval,
1075 SCIP_CONSDATA* consdata;
1078 int aggrtargetlength;
1083 assert( scip != NULL );
1084 assert( cons != NULL );
1085 assert( scalars != NULL );
1086 assert( naggrvars > 0 );
1087 assert( savedcol != NULL );
1088 assert( savedrow != NULL );
1089 assert( savedval != NULL );
1090 assert( nfixednonz != NULL );
1092 consdata = SCIPconsGetData(cons);
1093 assert( consdata != NULL );
1095 SCIP_CALL( SCIPgetRealParam(scip,
"numerics/epsilon", &epsilon) );
1098 startind = *nfixednonz;
1100 if ( SCIPisEQ(scip, constant, 0.0) )
1106 for (i = 0; i < consdata->nvarnonz[*v]; i++)
1108 savedcol[*nfixednonz] = consdata->col[*v][i];
1109 savedrow[*nfixednonz] = consdata->row[*v][i];
1110 savedval[*nfixednonz] = consdata->val[*v][i];
1116 for (i = 0; i < consdata->nvarnonz[*v]; i++)
1118 savedcol[*nfixednonz] = consdata->col[*v][i];
1119 savedrow[*nfixednonz] = consdata->row[*v][i];
1120 savedval[*nfixednonz] = consdata->val[*v][i] * constant;
1124 assert( *nfixednonz - startind == consdata->nvarnonz[*v] );
1132 SCIPfreeBlockMemoryArray(scip, &(consdata->val[*v]), consdata->nvarnonz[*v]);
1133 SCIPfreeBlockMemoryArray(scip, &(consdata->row[*v]), consdata->nvarnonz[*v]);
1134 SCIPfreeBlockMemoryArray(scip, &(consdata->col[*v]), consdata->nvarnonz[*v]);
1137 SCIP_CALL( SCIPreleaseVar(scip, &consdata->vars[*v]) );
1138 consdata->col[*v] = consdata->col[consdata->nvars - 1];
1139 consdata->row[*v] = consdata->row[consdata->nvars - 1];
1140 consdata->val[*v] = consdata->val[consdata->nvars - 1];
1141 consdata->nvarnonz[*v] = consdata->nvarnonz[consdata->nvars - 1];
1142 consdata->vars[*v] = consdata->vars[consdata->nvars - 1];
1143 (consdata->nvars)--;
1147 for (aggrind = 0; aggrind < naggrvars; aggrind++)
1151 for (i = 0; i < consdata->nvars; i++)
1153 if ( consdata->vars[i] == aggrvars[aggrind] )
1160 if ( aggrconsind > -1 )
1165 aggrtargetlength = consdata->nvarnonz[aggrconsind] + *nfixednonz - startind;
1166 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->row[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1167 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->col[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1168 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->val[aggrconsind]), consdata->nvarnonz[aggrconsind], aggrtargetlength) );
1170 if ( SCIPisEQ(scip, constant, 0.0) )
1174 SCIP_CALL(
SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), epsilon, savedrow + startind, savedcol + startind, savedval + startind,
1175 *nfixednonz - startind, TRUE, scalars[aggrind], consdata->row[aggrconsind], consdata->col[aggrconsind],
1176 consdata->val[aggrconsind], &(consdata->nvarnonz[aggrconsind]), aggrtargetlength) );
1182 SCIP_CALL(
SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), epsilon, savedrow + startind, savedcol + startind, savedval + startind,
1183 *nfixednonz - startind, TRUE, scalars[aggrind] / constant, consdata->row[aggrconsind], consdata->col[aggrconsind],
1184 consdata->val[aggrconsind], &(consdata->nvarnonz[aggrconsind]), aggrtargetlength) );
1188 assert( consdata->nvarnonz[aggrconsind] <= aggrtargetlength );
1189 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->row[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1190 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->col[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1191 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->val[aggrconsind]), aggrtargetlength, consdata->nvarnonz[aggrconsind]) );
1197 SCIPdebugMessage(
"adding variable %s to SDP constraint %s because of (multi-)aggregation\n", SCIPvarGetName(aggrvars[aggrind]), SCIPconsGetName(cons));
1200 if ( consdata->nvars == *vararraylength )
1202 globalnvars = SCIPgetNVars(scip);
1205 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col, *vararraylength, globalnvars) );
1206 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row, *vararraylength, globalnvars) );
1207 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val, *vararraylength, globalnvars) );
1208 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->nvarnonz, *vararraylength, globalnvars) );
1209 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->vars, *vararraylength, globalnvars) );
1210 *vararraylength = globalnvars;
1214 SCIP_CALL( SCIPcaptureVar(scip, aggrvars[aggrind]) );
1215 consdata->vars[consdata->nvars] = aggrvars[aggrind];
1218 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->col[consdata->nvars]), savedcol + startind, *nfixednonz - startind) );
1219 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->row[consdata->nvars]), savedrow + startind, *nfixednonz - startind) );
1225 if ( SCIPisEQ(scip, scalars[aggrind], constant) || SCIPisEQ(scip, constant, 0.0) )
1227 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(consdata->val[consdata->nvars]), savedval + startind, *nfixednonz - startind) );
1228 consdata->nvarnonz[consdata->nvars] = *nfixednonz - startind;
1233 SCIP_CALL( SCIPgetRealParam(scip,
"numerics/epsilon", &epsilon) );
1235 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->val[consdata->nvars]), *nfixednonz - startind) );
1237 consdata->nvarnonz[consdata->nvars] = 0;
1239 for (i = 0; i < *nfixednonz - startind; i++)
1242 if ( (scalars[i] / constant) * savedval[startind + i] >= epsilon )
1244 consdata->val[consdata->nvars][consdata->nvarnonz[consdata->nvars]] = (scalars[i] / constant) * savedval[startind + i];
1245 consdata->nvarnonz[consdata->nvars]++;
1250 if ( consdata->nvarnonz[consdata->nvars] > 0 )
1257 if ( SCIPisEQ(scip, constant, 0.0) )
1258 *nfixednonz = startind;
1273 SCIP_CONSDATA* consdata;
1278 SCIP_Real* savedval;
1283 SCIP_VAR** aggrvars;
1300 assert( scip != NULL );
1301 assert( conss != NULL );
1302 assert( nconss >= 0 );
1304 SCIPdebugMessage(
"Calling fixAndAggrVars with aggregate = %u\n", aggregate);
1306 SCIP_CALL( SCIPgetRealParam(scip,
"numerics/epsilon", &epsilon) );
1308 for (c = 0; c < nconss; ++c)
1310 assert( conss[c] != NULL );
1311 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(conss[c])),
"SDP") == 0);
1313 consdata = SCIPconsGetData(conss[c]);
1314 assert( consdata != NULL );
1317 SCIP_CALL( SCIPallocBufferArray(scip, &savedcol, consdata->nnonz) );
1318 SCIP_CALL( SCIPallocBufferArray(scip, &savedrow, consdata->nnonz) );
1319 SCIP_CALL( SCIPallocBufferArray(scip, &savedval, consdata->nnonz) );
1324 vararraylength = consdata->nvars;
1325 globalnvars = SCIPgetNVars(scip);
1327 for (v = 0; v < consdata->nvars; v++)
1331 if ( SCIPvarIsBinary(consdata->vars[v]) && SCIPvarIsNegated(consdata->vars[v]) )
1334 var = SCIPvarGetNegationVar(consdata->vars[v]);
1337 var = consdata->vars[v];
1340 if ( SCIPvarGetStatus(var) == SCIP_VARSTATUS_FIXED || SCIPisEQ(scip, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var)))
1342 assert( SCIPisEQ(scip, SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var)) );
1344 SCIPdebugMessage(
"Globally fixing Variable %s to value %f !\n", SCIPvarGetName(var), SCIPvarGetLbGlobal(var));
1346 if ( ((! negated) && (! SCIPisEQ(scip, SCIPvarGetLbGlobal(var), 0.0))) || (negated && SCIPisEQ(scip, SCIPvarGetLbGlobal(var), 0.0)) )
1351 for (i = 0; i < consdata->nvarnonz[v]; i++)
1353 savedcol[nfixednonz] = consdata->col[v][i];
1354 savedrow[nfixednonz] = consdata->row[v][i];
1357 savedval[nfixednonz] = consdata->val[v][i] * SCIPvarGetLbGlobal(var);
1359 savedval[nfixednonz] = consdata->val[v][i];
1368 consdata->nnonz -= consdata->nvarnonz[v];
1371 SCIPfreeBlockMemoryArrayNull(scip, &(consdata->val[v]), consdata->nvarnonz[v]);
1372 SCIPfreeBlockMemoryArrayNull(scip, &(consdata->row[v]), consdata->nvarnonz[v]);
1373 SCIPfreeBlockMemoryArrayNull(scip, &(consdata->col[v]), consdata->nvarnonz[v]);
1376 SCIP_CALL( SCIPreleaseVar(scip, &(consdata->vars[v])) );
1377 consdata->col[v] = consdata->col[consdata->nvars - 1];
1378 consdata->row[v] = consdata->row[consdata->nvars - 1];
1379 consdata->val[v] = consdata->val[consdata->nvars - 1];
1380 consdata->nvarnonz[v] = consdata->nvarnonz[consdata->nvars - 1];
1381 consdata->vars[v] = consdata->vars[consdata->nvars - 1];
1385 else if ( (SCIPvarGetStatus(var) == SCIP_VARSTATUS_AGGREGATED ||
1386 SCIPvarGetStatus(var) == SCIP_VARSTATUS_MULTAGGR)
1389 SCIP_CALL( SCIPallocBufferArray(scip, &aggrvars, globalnvars) );
1390 SCIP_CALL( SCIPallocBufferArray(scip, &scalars, globalnvars) );
1395 aggrvars[0] = consdata->vars[v];
1403 aggrvars[0] = consdata->vars[v];
1410 SCIP_CALL( SCIPgetProbvarLinearSum(scip, aggrvars, scalars, &naggrvars, globalnvars, &constant, &requiredsize, TRUE) );
1411 assert( requiredsize <= globalnvars );
1416 if ( SCIPvarGetStatus(consdata->vars[v]) == SCIP_VARSTATUS_AGGREGATED )
1417 SCIPdebugMessage(
"aggregating variable %s to ", SCIPvarGetName(var));
1419 SCIPdebugMessage(
"multiaggregating variable %s to ", SCIPvarGetName(var));
1420 for (i = 0; i < naggrvars; i++)
1421 SCIPdebugMessage(
"+ (%f2) * %s ", scalars[i], SCIPvarGetName(aggrvars[i]));
1422 SCIPdebugMessage(
"+ (%f2) \n", constant);
1427 SCIP_CALL(
multiaggrVar(scip, conss[c], &v, aggrvars, scalars, naggrvars, constant, savedcol, savedrow, savedval, &nfixednonz, &vararraylength) );
1429 SCIPfreeBufferArray(scip, &aggrvars);
1430 SCIPfreeBufferArray(scip, &scalars);
1432 else if ( negated && (SCIPvarGetStatus(var) == SCIP_VARSTATUS_LOOSE || SCIPvarGetStatus(var) == SCIP_VARSTATUS_COLUMN)
1437 SCIPdebugMessage(
"Changing variable %s to negation of variable %s !\n", SCIPvarGetName(consdata->vars[v]), SCIPvarGetName(var));
1441 SCIP_CALL(
multiaggrVar(scip, conss[c], &v, &var, &scalar, 1, 1.0, savedcol, savedrow, savedval, &nfixednonz, &vararraylength) );
1446 assert( consdata->nvars <= vararraylength );
1447 if ( consdata->nvars < vararraylength )
1449 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col, vararraylength, consdata->nvars) );
1450 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row, vararraylength, consdata->nvars) );
1451 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val, vararraylength, consdata->nvars) );
1452 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->nvarnonz, vararraylength, consdata->nvars) );
1453 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->vars, vararraylength, consdata->nvars) );
1457 arraylength = consdata->constnnonz + nfixednonz;
1458 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constcol), consdata->constnnonz, arraylength) );
1459 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constrow), consdata->constnnonz, arraylength) );
1460 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constval), consdata->constnnonz, arraylength) );
1463 SCIP_CALL(
SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), epsilon, savedrow, savedcol, savedval, nfixednonz, FALSE, -1.0, consdata->constrow,
1464 consdata->constcol, consdata->constval, &(consdata->constnnonz), arraylength) );
1466 assert( consdata->constnnonz <= arraylength );
1469 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constcol), arraylength, consdata->constnnonz) );
1470 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constrow), arraylength, consdata->constnnonz) );
1471 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &(consdata->constval), arraylength, consdata->constnnonz) );
1474 SCIPfreeBufferArray(scip, &savedval);
1475 SCIPfreeBufferArray(scip, &savedrow);
1476 SCIPfreeBufferArray(scip, &savedcol);
1479 consdata->nnonz = 0;
1480 for (v = 0; v < consdata->nvars; v++)
1481 consdata->nnonz += consdata->nvarnonz[v];
1491 SCIP_CONSHDLR* conshdlr,
1498 char cutname[SCIP_MAXSTRLEN];
1499 SCIP_CONSDATA* consdata;
1500 SCIP_CONSHDLRDATA* conshdlrdata;
1501 SCIP_Bool all_feasible = TRUE;
1502 SCIP_Bool separated = FALSE;
1504 SCIP_Bool infeasible;
1515 *result = SCIP_FEASIBLE;
1517 for (i = 0; i < nconss; ++i)
1519 consdata = SCIPconsGetData(conss[i]);
1521 if ( *result == SCIP_FEASIBLE )
1524 all_feasible = FALSE;
1526 nvars = consdata->nvars;
1530 SCIP_CALL( SCIPallocBufferArray(scip, &coeff, nvars) );
1533 rhs = SCIPinfinity(scip);
1534 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1537 snprintfreturn = SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
1538 assert( snprintfreturn < SCIP_MAXSTRLEN );
1540 (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN,
"sepa_eig_sdp_%d", ++(conshdlrdata->neigveccuts));
1543 SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, conshdlr, cutname , lhs, rhs, FALSE, FALSE, TRUE) );
1544 SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
1546 for (j = 0; j < nvars; ++j)
1548 SCIP_CALL( SCIPaddVarToRow(scip, row, consdata->vars[j], coeff[j]) );
1551 SCIP_CALL( SCIPflushRowExtensions(scip, row) );
1553 #ifdef SCIP_MORE_DEBUG
1554 SCIPinfoMessage(scip, NULL,
"Added cut %s: ", cutname);
1555 SCIPinfoMessage(scip, NULL,
"%f <= ", lhs);
1556 for (j = 0; j < nvars; j++)
1557 SCIPinfoMessage(scip, NULL,
"+ (%f)*%s", coeff[j], SCIPvarGetName(consdata->vars[j]));
1558 SCIPinfoMessage(scip, NULL,
"\n");
1561 SCIP_CALL( SCIPaddRow(scip, row, FALSE, &infeasible) );
1565 *result = SCIP_CUTOFF;
1567 SCIP_CALL( SCIPreleaseRow(scip, &row) );
1568 SCIPfreeBufferArray(scip, &coeff);
1574 SCIP_CALL( SCIPaddPoolCut(scip, row) );
1576 SCIP_CALL( SCIPresetConsAge(scip, conss[i]) );
1577 *result = SCIP_SEPARATED;
1580 SCIP_CALL( SCIPreleaseRow(scip, &row) );
1581 SCIPfreeBufferArray(scip, &coeff);
1588 *result = SCIP_SEPARATED;
1597 SCIP_CONSHDLRDATA* conshdlrdata;
1599 assert( conshdlr != NULL );
1601 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1602 assert( conshdlrdata != NULL );
1604 conshdlrdata->neigveccuts = 0;
1605 conshdlrdata->ndiaggezerocuts = 0;
1606 conshdlrdata->n1x1blocks = 0;
1616 SCIP_CONSDATA* consdata;
1620 SCIP_Real eigenvalue;
1622 consdata = SCIPconsGetData(cons);
1623 assert( consdata != NULL );
1624 blocksize = consdata->blocksize;
1625 nvars = consdata->nvars;
1627 SCIP_CALL( SCIPallocBufferArray(scip, &Aj, blocksize * blocksize) );
1629 for (var = 0; var < nvars; var++)
1635 if ( SCIPisNegative(scip, eigenvalue) )
1639 SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlocksneg, nlockspos) );
1643 if ( SCIPisPositive(scip, eigenvalue) )
1647 SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlockspos, nlocksneg) );
1653 if ( SCIPisPositive(scip, eigenvalue) )
1657 SCIP_CALL( SCIPaddVarLocks(scip, consdata->vars[var], nlockspos, nlocksneg) );
1662 SCIPfreeBufferArray(scip, &Aj);
1671 assert( scip != NULL );
1673 if ( conss == NULL )
1685 SCIP_CONSHDLRDATA* conshdlrdata;
1687 assert( conshdlr != NULL );
1688 assert( result != 0 );
1690 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1691 assert( conshdlrdata != NULL );
1696 if ( conshdlrdata->diaggezerocuts )
1698 SCIP_CALL(
diagGEzero(scip, conss, nconss, naddconss) );
1700 if ( conshdlrdata->diagzeroimplcuts )
1702 SCIP_CALL(
diagZeroImpl(scip, conss, nconss, naddconss) );
1713 SCIP_CONSDATA* sourcedata;
1714 SCIP_CONSDATA* targetdata;
1716 SCIP_CONSHDLRDATA* conshdlrdata;
1722 char transname[SCIP_MAXSTRLEN];
1724 sourcedata = SCIPconsGetData(sourcecons);
1725 assert( sourcedata != NULL );
1727 SCIPdebugMessage(
"Transforming constraint <%s>\n", SCIPconsGetName(sourcecons));
1730 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1731 SCIPdebugMessage(
"Setting number of threads to %d via OpenMP in Openblas.\n", conshdlrdata->nthreads);
1732 omp_set_num_threads(conshdlrdata->nthreads);
1735 SCIP_CALL( SCIPallocBlockMemory(scip, &targetdata) );
1738 targetdata->nvars = sourcedata->nvars;
1739 targetdata->nnonz = sourcedata->nnonz;
1740 targetdata->blocksize = sourcedata->blocksize;
1742 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->nvarnonz), sourcedata->nvarnonz, sourcedata->nvars) );
1745 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->col), sourcedata->nvars) );
1746 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->row), sourcedata->nvars) );
1747 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->val), sourcedata->nvars) );
1749 for (i = 0; i < sourcedata->nvars; i++)
1751 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->col[i]), sourcedata->col[i], sourcedata->nvarnonz[i]) );
1752 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->row[i]), sourcedata->row[i], sourcedata->nvarnonz[i]) );
1753 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->val[i]), sourcedata->val[i], sourcedata->nvarnonz[i]) );
1755 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(targetdata->vars), sourcedata->nvars));
1758 for (i = 0; i < sourcedata->nvars; i++)
1760 targetdata->vars[i] = SCIPvarGetTransVar(sourcedata->vars[i]);
1761 SCIP_CALL( SCIPcaptureVar(scip, targetdata->vars[i]) );
1765 targetdata->constnnonz = sourcedata->constnnonz;
1767 if ( sourcedata->constnnonz > 0 )
1769 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constcol), sourcedata->constcol, sourcedata->constnnonz));
1770 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constrow), sourcedata->constrow, sourcedata->constnnonz));
1771 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(targetdata->constval), sourcedata->constval, sourcedata->constnnonz));
1775 targetdata->constcol = NULL;
1776 targetdata->constrow = NULL;
1777 targetdata->constval = NULL;
1781 targetdata->maxrhsentry = sourcedata->maxrhsentry;
1785 snprintfreturn = SCIPsnprintf(transname, SCIP_MAXSTRLEN,
"t_%s", SCIPconsGetName(sourcecons));
1786 assert( snprintfreturn < SCIP_MAXSTRLEN );
1788 (void) SCIPsnprintf(transname, SCIP_MAXSTRLEN,
"t_%s", SCIPconsGetName(sourcecons));
1792 SCIP_CALL( SCIPcreateCons(scip, targetcons, transname, conshdlr, targetdata,
1793 SCIPconsIsInitial(sourcecons), SCIPconsIsSeparated(sourcecons), SCIPconsIsEnforced(sourcecons),
1794 SCIPconsIsChecked(sourcecons), SCIPconsIsPropagated(sourcecons), SCIPconsIsLocal(sourcecons),
1795 SCIPconsIsModifiable(sourcecons), SCIPconsIsDynamic(sourcecons), SCIPconsIsRemovable(sourcecons),
1796 SCIPconsIsStickingAtNode(sourcecons)) );
1807 assert( scip != NULL );
1808 assert( result != NULL );
1810 *result = SCIP_FEASIBLE;
1812 for (i = 0; i < nconss; ++i)
1814 SCIP_CALL(
SCIPconsSdpCheckSdpCons(scip, conss[i], sol, checkintegrality, checklprows, printreason, result) );
1815 if ( *result == SCIP_INFEASIBLE )
1831 assert( scip != NULL );
1832 assert( result != NULL );
1833 assert( conss != NULL );
1835 *result = SCIP_DIDNOTRUN;
1837 if ( objinfeasible )
1839 SCIPdebugMessage(
"-> pseudo solution is objective infeasible, return.\n");
1843 for (i = 0; i < nconss; ++i)
1847 if (*result == SCIP_INFEASIBLE)
1850 SCIPdebugMessage(
"-> pseudo solution infeasible for SDP-constraint %s, return.\n", SCIPconsGetName(conss[i]));
1855 SCIPdebugMessage(
"-> pseudo solution feasible for all SDP-constraints.\n");
1866 assert( scip != NULL );
1867 assert( conshdlr != NULL );
1868 assert( conss != NULL );
1869 assert( result != NULL );
1879 assert( scip != NULL );
1880 assert( conshdlr != NULL );
1881 assert( conss != NULL );
1882 assert( result != NULL );
1893 assert( result != NULL );
1894 *result = SCIP_DIDNOTFIND;
1896 for (i = 0; i < nusefulconss; ++i)
1898 SCIP_CALL(
separateSol(scip, conshdlr, conss[i], sol, result) );
1910 assert( result != NULL );
1911 *result = SCIP_DIDNOTFIND;
1913 for (i = 0; i < nusefulconss; ++i)
1915 SCIP_CALL(
separateSol(scip, conshdlr, conss[i], NULL, result) );
1927 assert( cons != NULL );
1928 assert( consdata != NULL );
1930 SCIPdebugMessage(
"deleting SDP constraint <%s>.\n", SCIPconsGetName(cons));
1932 for (i = 0; i < (*consdata)->nvars; i++)
1934 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->val[i], (*consdata)->nvarnonz[i]);
1935 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->row[i], (*consdata)->nvarnonz[i]);
1936 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->col[i], (*consdata)->nvarnonz[i]);
1940 for (i = 0; i < (*consdata)->nvars; i++)
1942 SCIP_CALL( SCIPreleaseVar(scip, &((*consdata)->vars[i])) );
1945 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->vars, (*consdata)->nvars);
1946 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constval, (*consdata)->constnnonz);
1947 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constrow, (*consdata)->constnnonz);
1948 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->constcol, (*consdata)->constnnonz);
1949 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->val, (*consdata)->nvars);
1950 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->row, (*consdata)->nvars);
1951 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->col, (*consdata)->nvars);
1952 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->nvarnonz, (*consdata)->nvars);
1953 SCIPfreeBlockMemory(scip, consdata);
1962 SCIP_CONSHDLRDATA* conshdlrdata;
1964 conshdlrdata = SCIPconshdlrGetData(conshdlr);
1965 assert( conshdlrdata != NULL );
1967 SCIPfreeMemory(scip, &conshdlrdata);
1968 SCIPconshdlrSetData(conshdlr, NULL);
1977 assert(scip != NULL);
1978 assert(conshdlr != NULL);
1979 assert(strcmp(SCIPconshdlrGetName(conshdlr),
CONSHDLR_NAME) == 0);
1992 SCIP_CONSDATA* sourcedata;
1994 SCIP_VAR** targetvars;
1998 assert( scip != NULL );
1999 assert( sourcescip != NULL );
2000 assert( sourcecons != NULL );
2001 assert( strcmp(SCIPconshdlrGetName(SCIPconsGetHdlr(sourcecons)),
CONSHDLR_NAME) == 0 );
2002 assert( valid != NULL );
2004 SCIPdebugMessage(
"Copying SDP constraint <%s>\n", SCIPconsGetName(sourcecons));
2010 if ( SCIPgetStage(sourcescip) <= SCIP_STAGE_EXITPRESOLVE )
2015 sourcedata = SCIPconsGetData(sourcecons);
2016 assert( sourcedata != NULL );
2018 SCIP_CALL( SCIPallocBufferArray(scip, &targetvars, sourcedata->nvars) );
2021 for (i = 0; i < sourcedata->nvars; i++)
2023 SCIP_CALL( SCIPgetVarCopy(sourcescip, scip, sourcedata->vars[i], &var, varmap, consmap, global, &success) );
2025 targetvars[i] = var;
2036 char copyname[SCIP_MAXSTRLEN];
2040 snprintfreturn = SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s",
name);
2041 assert( snprintfreturn < SCIP_MAXSTRLEN );
2043 (void) SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s",
name);
2045 SCIP_CALL(
SCIPcreateConsSdp( scip, cons, copyname, sourcedata->nvars, sourcedata->nnonz, sourcedata->blocksize, sourcedata->nvarnonz,
2046 sourcedata->col, sourcedata->row, sourcedata->val, targetvars, sourcedata->constnnonz,
2047 sourcedata->constcol, sourcedata->constrow, sourcedata->constval) );
2054 char copyname[SCIP_MAXSTRLEN];
2058 snprintfreturn = SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s", SCIPconsGetName(sourcecons));
2059 assert( snprintfreturn < SCIP_MAXSTRLEN );
2061 (void) SCIPsnprintf(copyname, SCIP_MAXSTRLEN,
"c_%s", SCIPconsGetName(sourcecons));
2063 SCIP_CALL(
SCIPcreateConsSdp( scip, cons, SCIPconsGetName(sourcecons), sourcedata->nvars, sourcedata->nnonz, sourcedata->blocksize,
2064 sourcedata->nvarnonz, sourcedata->col, sourcedata->row, sourcedata->val, targetvars, sourcedata->constnnonz,
2065 sourcedata->constcol, sourcedata->constrow, sourcedata->constval) );
2068 SCIPfreeBufferArray(scip, &targetvars);
2077 #ifdef PRINT_HUMAN_READABLE
2078 SCIP_CONSDATA* consdata;
2079 SCIP_Real* fullmatrix;
2084 assert( scip != NULL );
2085 assert( cons != NULL );
2087 consdata = SCIPconsGetData(cons);
2088 assert( consdata != NULL );
2090 SCIP_CALL( SCIPallocBufferArray(scip, &fullmatrix, consdata->blocksize * consdata->blocksize) );
2093 for (v = 0; v < consdata->nvars; v++)
2098 for (i = 0; i < consdata->blocksize; i++)
2100 for (j = 0; j < consdata->blocksize; j++)
2101 fullmatrix[i * consdata->blocksize + j] = 0.0;
2105 for (i = 0; i < consdata->nvarnonz[v]; i++)
2107 fullmatrix[consdata->row[v][i] * consdata->blocksize + consdata->col[v][i]] = consdata->val[v][i];
2108 fullmatrix[consdata->col[v][i] * consdata->blocksize + consdata->row[v][i]] = consdata->val[v][i];
2112 SCIPinfoMessage(scip, file,
"+\n");
2113 for (i = 0; i < consdata->blocksize; i++)
2115 SCIPinfoMessage(scip, file,
"( ");
2116 for (j = 0; j < consdata->blocksize; j++)
2117 SCIPinfoMessage(scip, file,
"%f ", fullmatrix[i * consdata->blocksize + j]);
2118 SCIPinfoMessage(scip, file,
")\n");
2120 SCIPinfoMessage(scip, file,
"* %s\n", SCIPvarGetName(consdata->vars[v]));
2128 for (i = 0; i < consdata->blocksize; i++)
2130 for (j = 0; j < consdata->blocksize; j++)
2131 fullmatrix[i * consdata->blocksize + j] = 0.0;
2135 for (i = 0; i < consdata->constnnonz; i++)
2137 fullmatrix[consdata->constrow[i] * consdata->blocksize + consdata->constcol[i]] = consdata->constval[i];
2138 fullmatrix[consdata->constcol[i] * consdata->blocksize + consdata->constrow[i]] = consdata->constval[i];
2142 SCIPinfoMessage(scip, file,
"-\n");
2143 for (i = 0; i < consdata->blocksize; i++)
2145 SCIPinfoMessage(scip, file,
"( ");
2146 for (j = 0; j < consdata->blocksize; j++)
2147 SCIPinfoMessage(scip, file,
"%f ", fullmatrix[i * consdata->blocksize + j]);
2148 SCIPinfoMessage(scip, file,
")\n");
2150 SCIPinfoMessage(scip, file,
">=0\n");
2152 SCIPfreeBufferArray(scip, &fullmatrix);
2156 SCIP_CONSDATA* consdata;
2160 assert( scip != NULL );
2161 assert( cons != NULL );
2163 consdata = SCIPconsGetData(cons);
2166 SCIPinfoMessage(scip, file,
"%d\n", consdata->blocksize);
2169 if ( consdata->constnnonz > 0 )
2171 SCIPinfoMessage(scip, file,
"A_0: ");
2173 for (i = 0; i < consdata->constnnonz; i++)
2175 SCIPinfoMessage(scip, file,
"(%d,%d):%f, ", consdata->constrow[i], consdata->constcol[i], consdata->constval[i]);
2177 SCIPinfoMessage(scip, file,
"\n");
2181 for (v = 0; v < consdata->nvars; v++)
2183 SCIPinfoMessage(scip, file,
"<%s>: ", SCIPvarGetName(consdata->vars[v]));
2184 for (i = 0; i < consdata->nvarnonz[v]; i++)
2186 SCIPinfoMessage(scip, file,
"(%d,%d):%f, ", consdata->row[v][i], consdata->col[v][i], consdata->val[v][i]);
2189 if (v < consdata->nvars - 1)
2191 SCIPinfoMessage(scip, file,
"\n");
2203 SCIP_Bool parsesuccess;
2204 SCIP_CONSDATA* consdata = NULL;
2210 assert( scip != NULL );
2211 assert( str != NULL );
2213 nvars = SCIPgetNVars(scip);
2215 assert( success != NULL );
2219 SCIP_CALL( SCIPallocBlockMemory(scip, &consdata) );
2220 consdata->nvars = 0;
2221 consdata->nnonz = 0;
2222 consdata->constnnonz = 0;
2223 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->nvarnonz, nvars) );
2224 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col, nvars) );
2225 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row, nvars) );
2226 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val, nvars) );
2227 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->vars, nvars));
2228 consdata->constcol = NULL;
2229 consdata->constrow = NULL;
2230 consdata->constval = NULL;
2233 parsesuccess = SCIPstrToIntValue(str, &(consdata->blocksize), &pos);
2234 *success = *success && parsesuccess;
2237 while( isspace((
unsigned char)*pos) )
2241 if ( pos[0] ==
'A' && pos[1] ==
'_' && pos[2] ==
'0' )
2245 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constcol,
PARSE_STARTSIZE) );
2246 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constrow,
PARSE_STARTSIZE) );
2247 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constval,
PARSE_STARTSIZE) );
2252 while (pos[0] ==
'(')
2257 if ( consdata->constnnonz == currentsize )
2259 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constcol, currentsize,
PARSE_SIZEFACTOR * currentsize) );
2260 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constrow, currentsize,
PARSE_SIZEFACTOR * currentsize) );
2261 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constval, currentsize,
PARSE_SIZEFACTOR * currentsize) );
2265 parsesuccess = SCIPstrToIntValue(pos, &(consdata->constrow[consdata->constnnonz]), &pos);
2266 *success = *success && parsesuccess;
2267 assert( consdata->constrow[consdata->constnnonz] < consdata->blocksize );
2269 parsesuccess = SCIPstrToIntValue(pos, &(consdata->constcol[consdata->constnnonz]), &pos);
2270 *success = *success && parsesuccess;
2271 assert( consdata->constcol[consdata->constnnonz] < consdata->blocksize );
2273 parsesuccess = SCIPstrToRealValue(pos, &(consdata->constval[consdata->constnnonz]), &pos);
2274 *success = *success && parsesuccess;
2278 if ( consdata->constcol[consdata->constnnonz] > consdata->constrow[consdata->constnnonz] )
2280 i = consdata->constcol[consdata->constnnonz];
2281 consdata->constcol[consdata->constnnonz] = consdata->constrow[consdata->constnnonz];
2282 consdata->constrow[consdata->constnnonz] = i;
2285 consdata->constnnonz++;
2288 while( isspace((
unsigned char)*pos) )
2293 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constcol, currentsize, consdata->constnnonz) );
2294 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constrow, currentsize, consdata->constnnonz) );
2295 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->constval, currentsize, consdata->constnnonz) );
2299 while( isspace((
unsigned char)*pos) )
2305 while (pos[0] ==
'<')
2308 SCIP_CALL( SCIPparseVarName(scip, pos, &(consdata->vars[consdata->nvars]), &pos) );
2309 SCIP_CALL( SCIPcaptureVar(scip, consdata->vars[consdata->nvars]) );
2311 consdata->nvarnonz[consdata->nvars] = 0;
2312 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->col[consdata->nvars]),
PARSE_STARTSIZE));
2313 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->row[consdata->nvars]),
PARSE_STARTSIZE));
2314 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(consdata->val[consdata->nvars]),
PARSE_STARTSIZE));
2321 while (pos[0] ==
'(')
2326 if ( consdata->nvarnonz[consdata->nvars - 1] == currentsize )
2328 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col[consdata->nvars - 1], currentsize,
PARSE_SIZEFACTOR * currentsize) );
2329 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row[consdata->nvars - 1], currentsize,
PARSE_SIZEFACTOR * currentsize) );
2330 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val[consdata->nvars - 1], currentsize,
PARSE_SIZEFACTOR * currentsize) );
2334 parsesuccess = SCIPstrToIntValue(pos, &(consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]]), &pos);
2335 *success = *success && parsesuccess;
2336 assert( consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] < consdata->blocksize );
2338 parsesuccess = SCIPstrToIntValue(pos, &(consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]]), &pos);
2339 *success = *success && parsesuccess;
2340 assert( consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] < consdata->blocksize );
2342 parsesuccess = SCIPstrToRealValue(pos, &(consdata->val[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]]), &pos);
2343 *success = *success && parsesuccess;
2347 if ( consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] >
2348 consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] )
2350 i = consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]];
2351 consdata->col[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] =
2352 consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]];
2353 consdata->row[consdata->nvars - 1][consdata->nvarnonz[consdata->nvars - 1]] = i;
2356 consdata->nvarnonz[consdata->nvars - 1]++;
2359 while( isspace((
unsigned char)*pos) )
2364 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->col[consdata->nvars - 1], currentsize, consdata->nvarnonz[consdata->nvars - 1]) );
2365 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->row[consdata->nvars - 1], currentsize, consdata->nvarnonz[consdata->nvars - 1]) );
2366 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->val[consdata->nvars - 1], currentsize, consdata->nvarnonz[consdata->nvars - 1]) );
2369 while( isspace((
unsigned char)*pos) )
2374 SCIP_CALL( SCIPcreateCons(scip, cons,
name, conshdlr, consdata, initial, separate, enforce, check, propagate, local, modifiable,
2375 dynamic, removable, stickingatnode) );
2380 #ifdef SCIP_MORE_DEBUG
2381 SCIP_CALL( SCIPprintCons(scip, *cons, NULL) );
2391 SCIP_CONSDATA* consdata;
2395 assert( scip != NULL );
2396 assert( cons != NULL );
2397 assert( vars != NULL );
2398 assert( success != NULL );
2399 assert( varssize >= 0 );
2401 consdata = SCIPconsGetData(cons);
2402 assert( consdata != NULL );
2404 nvars = consdata->nvars;
2406 if ( nvars > varssize )
2408 SCIPdebugMessage(
"consGetVarsIndicator called for array of size %d, needed size %d.\n", varssize, nvars);
2413 for (i = 0; i < nvars; i++)
2414 vars[i] = consdata->vars[i];
2425 SCIP_CONSDATA* consdata;
2427 assert( scip != NULL );
2428 assert( cons != NULL );
2429 assert( nvars != NULL );
2430 assert( success != NULL );
2432 consdata = SCIPconsGetData(cons);
2433 assert( consdata != NULL );
2435 *nvars = consdata->nvars;
2446 SCIP_CONSHDLR* conshdlr = NULL;
2447 SCIP_CONSHDLRDATA* conshdlrdata = NULL;
2449 assert( scip != 0 );
2452 SCIP_CALL( SCIPallocMemory(scip, &conshdlrdata) );
2457 consEnfolpSdp, consEnfopsSdp, consCheckSdp, consLockSdp, conshdlrdata) );
2459 assert( conshdlr != NULL );
2462 SCIP_CALL( SCIPsetConshdlrDelete(scip, conshdlr, consDeleteSdp) );
2463 SCIP_CALL( SCIPsetConshdlrFree(scip, conshdlr, consFreeSdp) );
2464 SCIP_CALL( SCIPsetConshdlrCopy(scip, conshdlr, conshdlrCopySdp, consCopySdp) );
2465 SCIP_CALL( SCIPsetConshdlrInitpre(scip, conshdlr,consInitpreSdp) );
2466 SCIP_CALL( SCIPsetConshdlrExitpre(scip, conshdlr, consExitpreSdp) );
2468 SCIP_CALL( SCIPsetConshdlrSepa(scip, conshdlr, consSepalpSdp, consSepasolSdp,
CONSHDLR_SEPAFREQ,
2470 SCIP_CALL( SCIPsetConshdlrEnforelax(scip, conshdlr, consEnforelaxSdp) );
2471 SCIP_CALL( SCIPsetConshdlrTrans(scip, conshdlr, consTransSdp) );
2472 SCIP_CALL( SCIPsetConshdlrPrint(scip, conshdlr, consPrintSdp) );
2473 SCIP_CALL( SCIPsetConshdlrParse(scip, conshdlr, consParseSdp) );
2474 SCIP_CALL( SCIPsetConshdlrGetVars(scip, conshdlr, consGetVarsSdp) );
2475 SCIP_CALL( SCIPsetConshdlrGetNVars(scip, conshdlr, consGetNVarsSdp) );
2479 SCIP_CALL( SCIPaddIntParam(scip,
"constraints/SDP/threads",
"number of threads used for OpenBLAS",
2480 &(conshdlrdata->nthreads), TRUE, DEFAULT_NTHREADS, 1, INT_MAX, NULL, NULL) );
2482 SCIP_CALL( SCIPaddBoolParam(scip,
"constraints/SDP/diaggezerocuts",
2483 "Should linear cuts enforcing the non-negativity of diagonal entries of SDP-matrices be added?",
2485 SCIP_CALL( SCIPaddBoolParam(scip,
"constraints/SDP/diagzeroimplcuts",
2486 "Should linear cuts enforcing the implications of diagonal entries of zero in SDP-matrices be added?",
2503 return i*(i+1)/2 + j;
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++)
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;
2859 SCIP_CONSDATA* consdata;
2861 assert( scip != NULL );
2862 assert( cons != NULL );
2864 consdata = SCIPconsGetData(cons);
2866 return consdata->maxrhsentry;
2875 SCIP_CONSDATA* consdata;
2880 assert( scip != NULL );
2881 assert( cons != NULL );
2883 consdata = SCIPconsGetData(cons);
2887 for (v = 0; v < consdata->nvars; v++)
2889 for (i = 0; i < consdata->nvarnonz[v]; i++)
2891 if ( SCIPisGT(scip, REALABS(consdata->val[v][i]), maxcoef) )
2892 maxcoef = REALABS(consdata->val[v][i]);
2909 SCIP_CONSDATA* consdata;
2914 assert( cons != NULL );
2916 consdata = SCIPconsGetData(cons);
2917 assert( consdata != NULL );
2919 ub = consdata->constnnonz;
2921 for (v = 0; v < consdata->nvars; v++)
2922 ub += consdata->nvarnonz[v];
2924 denselength = consdata->blocksize * (consdata->blocksize + 1) / 2;
2926 return (ub <= denselength ? ub : denselength);
2944 SCIP_CONSDATA* consdata;
2950 assert( scip != NULL );
2951 assert( cons != NULL );
2952 assert( sol != NULL );
2953 assert( length != NULL );
2954 assert( row != NULL );
2955 assert( col != NULL );
2956 assert( val != NULL );
2958 consdata = SCIPconsGetData(cons);
2959 assert( consdata != NULL );
2962 nnonz = consdata->constnnonz;
2963 if ( *length < nnonz )
2966 SCIPerrorMessage(
"Arrays not long enough in SCIPconsSdpComputeSparseSdpMatrix, length %d given, need at least %d (probably more)\n",
2970 for (i = 0; i < consdata->constnnonz; i++)
2972 row[i] = consdata->constrow[i];
2973 col[i] = consdata->constcol[i];
2974 val[i] = -1 * consdata->constval[i];
2978 SCIP_CALL( SCIPgetRealParam(scip,
"numerics/epsilon", &epsilon) );
2980 for (v = 0; v < consdata->nvars; v++)
2982 SCIP_CALL(
SCIPsdpVarfixerMergeArrays(SCIPblkmem(scip), epsilon, consdata->row[v], consdata->col[v], consdata->val[v], consdata->nvarnonz[v],
2983 FALSE, SCIPgetSolVal(scip, sol, consdata->vars[v]), row, col, val, &nnonz, *length) );
2984 if ( nnonz > *length )
2987 SCIPerrorMessage(
"Arrays not long enough in SCIPconsSdpComputeSparseSdpMatrix, length %d given, need at least %d (probably more)\n",
3018 SCIP_CONSHDLR* conshdlr;
3019 SCIP_CONSDATA* consdata = NULL;
3023 assert( scip != NULL );
3024 assert( cons != NULL );
3025 assert( name != NULL );
3026 assert( nvars >= 0 );
3027 assert( nnonz >= 0 );
3028 assert( blocksize >= 0 );
3029 assert( constnnonz >= 0 );
3030 assert( nvars == 0 || vars != NULL );
3031 assert( nnonz == 0 || (nvarnonz != NULL && col != NULL && row != NULL && val != NULL ));
3032 assert( constnnonz == 0 || (constcol != NULL && constrow != NULL && constval != NULL ));
3034 conshdlr = SCIPfindConshdlr(scip,
"SDP");
3035 if ( conshdlr == NULL )
3037 SCIPerrorMessage(
"SDP constraint handler not found\n");
3038 return SCIP_PLUGINNOTFOUND;
3042 SCIP_CALL( SCIPallocBlockMemory(scip, &consdata) );
3043 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->nvarnonz, nvars) );
3044 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col, nvars) );
3045 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row, nvars) );
3046 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val, nvars) );
3047 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constcol, constnnonz) );
3048 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constrow, constnnonz) );
3049 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->constval, constnnonz) );
3050 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->vars, nvars));
3052 for (i = 0; i < nvars; i++)
3054 assert( nvarnonz[i] >= 0 );
3056 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->col[i], nvarnonz[i]));
3057 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->row[i], nvarnonz[i]));
3058 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->val[i], nvarnonz[i]));
3061 consdata->nvars = nvars;
3062 consdata->nnonz = nnonz;
3063 consdata->constnnonz = constnnonz;
3064 consdata->blocksize = blocksize;
3066 for (i = 0; i < nvars; i++)
3068 consdata->nvarnonz[i] = nvarnonz[i];
3070 for (j = 0; j < nvarnonz[i]; j++)
3072 assert( col[i][j] >= 0 );
3073 assert( col[i][j] < blocksize );
3074 assert( row[i][j] >= 0 );
3075 assert( row[i][j] < blocksize );
3077 consdata->col[i][j] = col[i][j];
3078 consdata->row[i][j] = row[i][j];
3079 consdata->val[i][j] = val[i][j];
3083 for (i = 0; i < constnnonz; i++)
3085 consdata->constcol[i] = constcol[i];
3086 consdata->constrow[i] = constrow[i];
3087 consdata->constval[i] = constval[i];
3090 for (i = 0; i < nvars; i++)
3092 consdata->vars[i] = vars[i];
3093 SCIP_CALL( SCIPcaptureVar(scip, consdata->vars[i]) );
3095 SCIPdebugMessage(
"creating cons %s\n", name);
3097 SCIP_CALL( SCIPcreateCons(scip, cons, name, conshdlr, consdata, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE) );
SCIP_Real SCIPconsSdpGetMaxConstEntry(SCIP *scip, SCIP_CONS *cons)
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 SCIPconsSdpComputeSparseSdpMatrix(SCIP *scip, SCIP_CONS *cons, SCIP_SOL *sol, int *length, int *row, int *col, SCIP_Real *val)
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)
maps SCIP variables to SDP indices (the SCIP variables are given SDP indices in the order in which th...
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_Real SCIPconsSdpGetMaxSdpCoef(SCIP *scip, SCIP_CONS *cons)
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)
int SCIPconsSdpCompLowerTriangPos(int i, int j)
#define CONSHDLR_CHECKPRIORITY
static SCIP_DECL_CONSENFOPS(consEnfopsSdp)
#define DEFAULT_DIAGZEROIMPLCUTS
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)
int SCIPconsSdpComputeUbSparseSdpMatrixLength(SCIP *scip, SCIP_CONS *cons)
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)