All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
cons_quadratic.c
Go to the documentation of this file.
17 * @brief constraint handler for quadratic constraints \f$\textrm{lhs} \leq \sum_{i,j=1}^n a_{i,j} x_ix_j + \sum_{i=1}^n b_i x_i \leq \textrm{rhs}\f$
20 * @todo SCIP might fix linear variables on +/- infinity; remove them in presolve and take care later
23 * @todo check if some quadratic terms appear in several constraints and try to simplify (e.g., nous1)
26 * @todo sort order in bilinvar1/bilinvar2 such that the var which is involved in more terms is in bilinvar1, and use this info propagate and AddLinearReform
27 * @todo catch/drop events in consEnable/consDisable, do initsol/exitsol stuff also when a constraint is enabled/disabled during solve
31 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
55 #define CONSHDLR_ENFOPRIORITY -50 /**< priority of the constraint handler for constraint enforcing */
56 #define CONSHDLR_CHECKPRIORITY -4000000 /**< priority of the constraint handler for checking feasibility */
57 #define CONSHDLR_SEPAFREQ 1 /**< frequency for separating cuts; zero means to separate only in the root node */
58 #define CONSHDLR_PROPFREQ 1 /**< frequency for propagating domains; zero means only preprocessing propagation */
59 #define CONSHDLR_EAGERFREQ 100 /**< frequency for using all instead of only the useful constraints in separation,
61 #define CONSHDLR_MAXPREROUNDS -1 /**< maximal number of presolving rounds the constraint handler participates in (-1: no limit) */
62 #define CONSHDLR_DELAYSEPA FALSE /**< should separation method be delayed, if other separators found cuts? */
63 #define CONSHDLR_DELAYPROP FALSE /**< should propagation method be delayed, if other propagators found reductions? */
64 #define CONSHDLR_DELAYPRESOL FALSE /**< should presolving method be delayed, if other presolvers found reductions? */
65 #define CONSHDLR_NEEDSCONS TRUE /**< should the constraint handler be skipped, if no constraints are available? */
71 #define INITLPMAXVARVAL 1000.0 /**< maximal absolute value of variable for still generating a linearization cut at that point in initlp */
73 /* Activating this define enables reformulation of bilinear terms x*y with implications from x to y into linear terms.
74 * However, implications are not enforced by SCIP. Thus, if, e.g., the used implication was derived from this constraint and we then reformulate the constraint,
91 int varidx; /**< the index of the variable which bound change is caught, positive for linear variables, negative for quadratic variables */
128 unsigned int ispropagated:1; /**< was the constraint propagated with respect to the current bounds ? */
129 unsigned int ispresolved:1; /**< did we checked for possibilities of upgrading or implicit integer variables ? */
130 unsigned int initialmerge:1; /**< did we perform an initial merge and clean in presolving yet ? */
132 unsigned int isimpladded:1; /**< has there been an implication added for a binary variable in a bilinear term? */
135 SCIP_Real minlinactivity; /**< sum of minimal activities of all linear terms with finite minimal activity */
136 SCIP_Real maxlinactivity; /**< sum of maximal activities of all linear terms with finite maximal activity */
139 SCIP_INTERVAL quadactivitybounds; /**< bounds on the activity of the quadratic term, if up to date, otherwise empty interval */
141 SCIP_Real lhsviol; /**< violation of lower bound by current solution (used temporarily inside constraint handler) */
142 SCIP_Real rhsviol; /**< violation of lower bound by current solution (used temporarily inside constraint handler) */
144 int linvar_maydecrease; /**< index of a variable in linvars that may be decreased without making any other constraint infeasible, or -1 if none */
145 int linvar_mayincrease; /**< index of a variable in linvars that may be increased without making any other constraint infeasible, or -1 if none */
147 SCIP_VAR** sepaquadvars; /**< variables corresponding to quadvarterms to use in separation, only available in solving stage */
148 int* sepabilinvar2pos; /**< position of second variable in bilinear terms to use in separation, only available in solving stage */
149 SCIP_Real lincoefsmin; /**< minimal absolute value of coefficients in linear part, only available in solving stage */
150 SCIP_Real lincoefsmax; /**< maximal absolute value of coefficients in linear part, only available in solving stage */
153 SCIP_Real* factorright; /**< coefficients of right factor if constraint function is factorable */
159 SCIP_DECL_QUADCONSUPGD((*quadconsupgd)); /**< method to call for upgrading quadratic constraint */
163 typedef struct SCIP_QuadConsUpgrade SCIP_QUADCONSUPGRADE; /**< quadratic constraint update method */
168 int replacebinaryprodlength; /**< length of linear term which when multiplied with a binary variable is replaced by an auxiliary variable and an equivalent linear formulation */
169 int empathy4and; /**< how much empathy we have for using the AND constraint handler: 0 avoid always; 1 use sometimes; 2 use as often as possible */
170 SCIP_Bool binreforminitial; /**< whether to make constraints added due to replacing products with binary variables initial */
171 SCIP_Real binreformmaxcoef; /**< factor on 1/feastol to limit coefficients and coef range in linear constraints created by binary reformulation */
172 SCIP_Real mincutefficacysepa; /**< minimal efficacy of a cut in order to add it to relaxation during separation */
173 SCIP_Real mincutefficacyenfofac; /**< minimal target efficacy of a cut in order to add it to relaxation during enforcement as factor of feasibility tolerance (may be ignored) */
175 SCIP_Real cutmaxrange; /**< maximal range (maximal coef / minimal coef) of a cut in order to be added to LP */
176 SCIP_Bool linearizeheursol; /**< whether linearizations of convex quadratic constraints should be added to cutpool when some heuristics finds a new solution */
181 int maxproprounds; /**< limit on number of propagation rounds for a single constraint within one round of SCIP propagation during solve */
182 int maxproproundspresolve; /**< limit on number of propagation rounds for a single constraint within one presolving round */
183 SCIP_Real sepanlpmincont; /**< minimal required fraction of continuous variables in problem to use solution of NLP relaxation in root for separation */
184 SCIP_Bool enfocutsremovable; /**< are cuts added during enforcement removable from the LP in the same node? */
195 SCIP_NODE* lastenfolpnode; /**< the node for which enforcement was called the last time (and some constraint was violated) */
197 SCIP_QUADCONSUPGRADE** quadconsupgrades; /**< quadratic constraint upgrade methods for specializing quadratic constraints */
213 SCIP_DECL_QUADCONSUPGD((*quadconsupgd)), /**< method to call for upgrading quadratic constraint */
228 SCIPwarningMessage(scip, "Try to add already known upgrade message for constraint handler <%s>.\n", conshdlrname);
278 /* if right hand side is finite, then a tightening in the lower bound of coef*linvar is of interest
279 * since we also want to keep activities in consdata up-to-date, we also need to know when the corresponding bound is relaxed */
287 /* if left hand side is finite, then a tightening in the upper bound of coef*linvar is of interest
288 * since we also want to keep activities in consdata up-to-date, we also need to know when the corresponding bound is relaxed */
295 SCIP_CALL( SCIPcatchVarEvent(scip, consdata->linvars[linvarpos], eventtype, eventhdlr, (SCIP_EVENTDATA*)eventdata, &eventdata->filterpos) );
332 /* if right hand side is finite, then a tightening in the lower bound of coef*linvar is of interest
333 * since we also want to keep activities in consdata up-to-date, we also need to know when the corresponding bound is relaxed */
341 /* if left hand side is finite, then a tightening in the upper bound of coef*linvar is of interest
342 * since we also want to keep activities in consdata up-to-date, we also need to know when the corresponding bound is relaxed */
349 SCIP_CALL( SCIPdropVarEvent(scip, consdata->linvars[linvarpos], eventtype, eventhdlr, (SCIP_EVENTDATA*)consdata->lineventdata[linvarpos], consdata->lineventdata[linvarpos]->filterpos) );
388 SCIP_CALL( SCIPcatchVarEvent(scip, consdata->quadvarterms[quadvarpos].var, eventtype, eventhdlr, (SCIP_EVENTDATA*)eventdata, &eventdata->filterpos) );
426 SCIP_CALL( SCIPdropVarEvent(scip, consdata->quadvarterms[quadvarpos].var, eventtype, eventhdlr, (SCIP_EVENTDATA*)consdata->quadvarterms[quadvarpos].eventdata, consdata->quadvarterms[quadvarpos].eventdata->filterpos) );
460 consdata->isremovedfixings = consdata->isremovedfixings && SCIPvarIsActive(consdata->linvars[i]);
469 consdata->isremovedfixings = consdata->isremovedfixings && SCIPvarIsActive(consdata->quadvarterms[i].var);
539 SCIP_CALL( SCIPlockVarCons(scip, var, cons, !SCIPisInfinity(scip, -consdata->lhs), !SCIPisInfinity(scip, consdata->rhs)) );
543 SCIP_CALL( SCIPlockVarCons(scip, var, cons, !SCIPisInfinity(scip, consdata->rhs), !SCIPisInfinity(scip, -consdata->lhs)) );
570 SCIP_CALL( SCIPunlockVarCons(scip, var, cons, !SCIPisInfinity(scip, -consdata->lhs), !SCIPisInfinity(scip, consdata->rhs)) );
574 SCIP_CALL( SCIPunlockVarCons(scip, var, cons, !SCIPisInfinity(scip, consdata->rhs), !SCIPisInfinity(scip, -consdata->lhs)) );
625 /* if variable bounds are not strictly consistent, then the activity update methods may yield inconsistent activities
628 if( consdata->minlinactivity != SCIP_INVALID && consdata->maxlinactivity != SCIP_INVALID && /*lint !e777 */
629 (consdata->minlinactivityinf > 0 || consdata->maxlinactivityinf > 0 || consdata->minlinactivity <= consdata->maxlinactivity) )
721 assert(consdata->minlinactivityinf > 0 || consdata->maxlinactivityinf > 0 || consdata->minlinactivity <= consdata->maxlinactivity);
742 /* @todo since we check the linear activity for consistency later anyway, we may skip changing the rounding mode here */
845 /* @todo since we check the linear activity for consistency later anyway, we may skip changing the rounding mode here */
963 consdataUpdateLinearActivityLbChange(scip, consdata, consdata->lincoefs[varidx], SCIPeventGetOldbound(event), SCIPeventGetNewbound(event));
965 consdataUpdateLinearActivityUbChange(scip, consdata, consdata->lincoefs[varidx], SCIPeventGetOldbound(event), SCIPeventGetNewbound(event));
981 /* if variable is binary (quite likely if an implication has been added) and occurs in a bilinear term, then mark that we should check implications */
982 if( SCIPvarIsBinary(SCIPeventGetVar(event)) && consdata->quadvarterms[-varidx-1].nadjbilin > 0 )
1007 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->linvars, consdata->linvarssize, newsize) );
1008 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->lincoefs, consdata->linvarssize, newsize) );
1011 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->lineventdata, consdata->linvarssize, newsize) );
1037 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->quadvarterms, consdata->quadvarssize, newsize) );
1062 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &quadvarterm->adjbilin, quadvarterm->adjbilinsize, newsize) );
1087 SCIP_CALL( SCIPreallocBlockMemoryArray(scip, &consdata->bilinterms, consdata->bilintermssize, newsize) );
1196 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->quadvarterms, quadvarterms, nquadvars) );
1203 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->quadvarterms[i].adjbilin, quadvarterms[i].adjbilin, quadvarterms[i].nadjbilin) );
1231 SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->bilinterms, bilinterms, nbilinterms) );
1273 assert((*consdata)->lineventdata == NULL || (*consdata)->lineventdata[i] == NULL); /* variable events should have been dropped earlier */
1287 assert((*consdata)->quadvarterms[i].eventdata == NULL); /* variable events should have been dropped earlier */
1288 SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->quadvarterms[i].adjbilin, (*consdata)->quadvarterms[i].adjbilinsize);
1327 SCIPsortPtrReal((void**)consdata->linvars, consdata->lincoefs, SCIPvarComp, consdata->nlinvars);
1333 SCIPsortPtrPtrReal((void**)consdata->linvars, (void**)consdata->lineventdata, consdata->lincoefs, SCIPvarComp, consdata->nlinvars);
1345 /** returns the position of variable in the linear coefficients array of a constraint, or -1 if not found */
1362 if( !SCIPsortedvecFindPtr((void**)consdata->linvars, SCIPvarComp, (void*)var, consdata->nlinvars, &pos) )
1369 /** index comparison method for quadratic variable terms: compares two indices of the quadratic variable set in the quadratic constraint */
1451 /** returns the position of variable in the quadratic variable terms array of a constraint, or -1 if not found */
1457 int* pos /**< buffer where to store position of var in quadvarterms array, or -1 if not found */
1504 /** index comparison method for bilinear terms: compares two index pairs of the bilinear term set in the quadratic constraint */
1745 consdata->linvarssorted = consdata->linvarssorted && (SCIPvarCompare(consdata->linvars[consdata->nlinvars-2], consdata->linvars[consdata->nlinvars-1]) == -1);
1746 /* always set too FALSE since the new linear variable should be checked if already existing as quad var term */
1999 consdata->quadvarssorted = consdata->quadvarssorted && (SCIPvarCompare(consdata->quadvarterms[consdata->nquadvars-2].var, consdata->quadvarterms[consdata->nquadvars-1].var) == -1);
2000 /* also set to FALSE if nquadvars == 1, since the new variable should be checked for linearity and other stuff in mergeAndClean ... */
2053 SCIPfreeBlockMemoryArrayNull(scip, &consdata->quadvarterms[pos].adjbilin, consdata->quadvarterms[pos].adjbilinsize);
2078 * allows to replace x by coef*y+offset, thereby maintaining linear and square coefficients and bilinear terms */
2169 * thus, for now we just set the coefficient to 0.0 and clear in later when the bilinear terms are merged */
2277 if( SCIPvarCompare(consdata->quadvarterms[var1pos].var, consdata->quadvarterms[var2pos].var) < 0 )
2291 SCIPerrorMessage("tried to add bilinear term where both variables are the same, but appear at different positions in quadvarterms array\n");
2296 SCIP_CALL( consdataEnsureAdjBilinSize(scip, &consdata->quadvarterms[var1pos], consdata->quadvarterms[var1pos].nadjbilin + 1) );
2297 SCIP_CALL( consdataEnsureAdjBilinSize(scip, &consdata->quadvarterms[var2pos], consdata->quadvarterms[var2pos].nadjbilin + 1) );
2299 consdata->quadvarterms[var1pos].adjbilin[consdata->quadvarterms[var1pos].nadjbilin] = consdata->nbilinterms;
2300 consdata->quadvarterms[var2pos].adjbilin[consdata->quadvarterms[var2pos].nadjbilin] = consdata->nbilinterms;
2433 * also replaces squares of binary variables by the binary variables, i.e., adds sqrcoef to lincoef
2463 /* make sure quad var terms are sorted (do this in every round, since we may move variables around) */
2468 for( j = i+1; j < consdata->nquadvars && consdata->quadvarterms[j].var == quadvarterm->var; ++j )
2476 SCIP_CALL( consdataEnsureAdjBilinSize(scip, quadvarterm, quadvarterm->nadjbilin + consdata->quadvarterms[j].nadjbilin) );
2477 BMScopyMemoryArray(&quadvarterm->adjbilin[quadvarterm->nadjbilin], consdata->quadvarterms[j].adjbilin, consdata->quadvarterms[j].nadjbilin); /*lint !e866*/
2495 if( quadvarterm->sqrcoef != 0.0 && SCIPvarIsBinary(quadvarterm->var) && quadvarterm->nadjbilin == 0 )
2497 SCIPdebugMessage("replace square of binary variable by itself: <%s>^2 --> <%s>\n", SCIPvarGetName(quadvarterm->var), SCIPvarGetName(quadvarterm->var));
2531 /* merges entries with same linear variable into one entry and cleans up entries with coefficient 0.0 */
2561 /* make sure linear variables are sorted (do this in every round, since we may move variables around) */
2578 /* add newcoef to linear coefficient of quadratic variable and mark linear variable as to delete */
2603 /* merges bilinear terms with same variables into a single term, removes bilinear terms with coefficient 0.0 */
2644 for( j = i+1; j < consdata->nbilinterms && bilinterm->var1 == consdata->bilinterms[j].var1 && bilinterm->var2 == consdata->bilinterms[j].var2; ++j )
2714 SCIPdebugMessage(" linear term %g*<%s> is replaced by %g * <%s> + %g\n", consdata->lincoefs[i], SCIPvarGetName(consdata->linvars[i]), coef, SCIPvarGetName(var), offset);
2786 SCIPdebugMessage(" quadratic variable <%s> with status %d is replaced by %g * <%s> + %g\n", SCIPvarGetName(consdata->quadvarterms[i].var), SCIPvarGetStatus(consdata->quadvarterms[i].var), coef, SCIPvarGetName(var), offset);
2791 /* if not fixed to 0.0, add to linear coefs of vars in bilinear terms, and deal with linear and square term as constant */
2810 offset = consdata->quadvarterms[i].lincoef * offset + consdata->quadvarterms[i].sqrcoef * offset * offset;
2818 SCIP_CALL( removeBilinearTermsPos(scip, cons, consdata->quadvarterms[i].nadjbilin, consdata->quadvarterms[i].adjbilin) );
2828 /* if GetProbvar gave an active variable, replace the quad var term so that it uses the new variable */
2838 /* if GetProbVar gave a multi-aggregated variable, add new quad var terms and new bilinear terms
2891 SCIP_CALL( consdataEnsureBilinSize(scip, consdata, consdata->nquadvars + (scoef != 0.0 ? (naggrs * (naggrs-1))/2 : 0) + consdata->quadvarterms[j].nadjbilin * naggrs) );
2912 /* this is not efficient, but we cannot sort the quadratic terms here, since we currently iterate over them */
2937 SCIP_CALL( removeBilinearTermsPos(scip, cons, consdata->quadvarterms[i].nadjbilin, consdata->quadvarterms[i].adjbilin) );
2961 * in this case, we want the linear variable to be removed, which happens in mergeAndCleanLinearVars
2994 SCIP_VAR** quadlinvars; /* variables of linear terms using variables that are in quadratic terms */
2995 SCIP_Real* quadlincoefs; /* coefficients of linear terms using variables that are in quadratic terms */
3057 * thus, it makes sense to remember the index of the previous first variable for the case a series of bilinear terms with the same first var appears */
3089 SCIP_CALL( SCIPaddLinearCoefsToNlRow(scip, consdata->nlrow, nquadlinterms, quadlinvars, quadlincoefs) );
3123 if( SCIPisEQ(scip, consdata->lhs, consdata->rhs) && consdata->nlinvars == 0 && consdata->nquadvars == 2 &&
3124 ((SCIPvarIsBinary(consdata->quadvarterms[0].var) && consdata->quadvarterms[1].sqrcoef == 0.0) ||
3125 (SCIPvarIsBinary(consdata->quadvarterms[1].var) && consdata->quadvarterms[0].sqrcoef == 0.0)) )
3140 * if b != 0 and b+c != 0, then y = (rhs-a)/(b+c) * x + rhs/b * (1-x) = ((rhs-a)/(b+c) - rhs/b) * x + rhs/b
3143 binvaridx = (SCIPvarIsBinary(consdata->quadvarterms[0].var) && consdata->quadvarterms[1].sqrcoef == 0.0) ? 0 : 1;
3151 assert(consdata->nbilinterms <= 1); /* should actually be 1, since constraint is otherwise linear */
3156 SCIPdebugMessage("<%s> = 0 -> %g*<%s> = %g and <%s> = 1 -> %g*<%s> = %g\n", SCIPvarGetName(x), b, SCIPvarGetName(y), consdata->rhs, SCIPvarGetName(x), b+c, SCIPvarGetName(y), consdata->rhs - a);
3157 SCIPdebugMessage("=> attempt aggregation <%s> = %g*<%s> + %g\n", SCIPvarGetName(y), (consdata->rhs-a)/(b+c) - consdata->rhs/b, SCIPvarGetName(x), consdata->rhs/b);
3159 SCIP_CALL( SCIPaggregateVars(scip, x, y, (consdata->rhs-a)/(b+c) - consdata->rhs/b, -1.0, -consdata->rhs/b, &infeasible, redundant, &aggregated) );
3181 * For a product x*y, with x and y binary variables, the product is replaced by a new auxiliary variable z and the constraint z = {x and y} is added.
3234 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "prod%s_%s_%s", SCIPvarGetName(vars[0]), SCIPvarGetName(vars[1]), SCIPconsGetName(cons));
3236 SCIPvarIsInitial(vars[0]) || SCIPvarIsInitial(vars[1]), SCIPvarIsRemovable(vars[0]) && SCIPvarIsRemovable(vars[1]), NULL, NULL, NULL, NULL, NULL) );
3250 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "%sAND%s", SCIPvarGetName(vars[0]), SCIPvarGetName(vars[1]));
3300 SCIPintervalSetBounds(resultant, MIN(SCIPvarGetLbGlobal(y), SCIPvarGetUbGlobal(y)), MAX(SCIPvarGetLbGlobal(y), SCIPvarGetUbGlobal(y))); /*lint !e666 */
3326 if( !SCIPsortedvecFindPtr((void**)&implvars[nbinimpls], SCIPvarComp, (void*)y, nimpls - nbinimpls, &pos) )
3351 /** Reformulates products of binary times bounded continuous variables as system of linear inequalities (plus auxiliary variable).
3354 * an auxiliary variable z and the inequalities \f$ x^L * y \leq z \leq x^U * y \f$ and \f$ x - (1-y)*x^U \leq z \leq x - (1-y)*x^L \f$ are added.
3356 * If x is a linear term consisting of more than one variable, it is split up in groups of linear terms of length at most maxnrvar.
3357 * For each product of linear term of length at most maxnrvar with y, an auxiliary z and linear inequalities are added.
3359 * If y is a binary variable, the AND constraint \f$ z = x \wedge y \f$ may be added instead of linear constraints.
3429 SCIP_CALL( SCIPreallocBufferArray(scip, &xvars, MIN(maxnrvar, nbilinterms)+2) ); /* add 2 for later use when creating linear constraints */
3439 /* setup a list of bounded variables x_i with coefficients a_i that are multiplied with binary y: y*(sum_i a_i*x_i)
3467 if( SCIPisInfinity(scip, -SCIPvarGetLbGlobal(bvar)) || SCIPisInfinity(scip, SCIPvarGetUbGlobal(bvar)) )
3486 if( SCIPfeastol(scip) * REALABS(act0.inf) >= conshdlrdata->binreformmaxcoef || SCIPfeastol(scip) * REALABS(act0.sup) >= conshdlrdata->binreformmaxcoef )
3489 bilincoef, SCIPvarGetName(y), SCIPvarGetName(bvar), SCIPintervalGetInf(act0), SCIPintervalGetSup(act0), SCIPvarGetName(y));
3492 if( SCIPfeastol(scip) * REALABS(act1.inf) >= conshdlrdata->binreformmaxcoef || SCIPfeastol(scip) * REALABS(act1.sup) >= conshdlrdata->binreformmaxcoef )
3495 bilincoef, SCIPvarGetName(y), SCIPvarGetName(bvar), SCIPintervalGetInf(act1), SCIPintervalGetSup(act1), SCIPvarGetName(y));
3499 SCIPfeastol(scip) * MAX(REALABS(act0.inf), REALABS(act0.sup)) / MIN(REALABS(act0.inf), REALABS(act0.sup)) >= conshdlrdata->binreformmaxcoef )
3502 bilincoef, SCIPvarGetName(y), SCIPvarGetName(bvar), MAX(REALABS(act0.inf), REALABS(act0.sup)) / MIN(REALABS(act0.inf), REALABS(act0.sup)), SCIPvarGetName(y));
3506 SCIPfeastol(scip) * MAX(REALABS(act1.inf), REALABS(act1.sup)) / MIN(REALABS(act1.inf), REALABS(act1.sup)) >= conshdlrdata->binreformmaxcoef )
3509 bilincoef, SCIPvarGetName(y), SCIPvarGetName(bvar), MAX(REALABS(act1.inf), REALABS(act1.sup)) / MIN(REALABS(act1.inf), REALABS(act1.sup)), SCIPvarGetName(y));
3528 integral &= (SCIPvarGetType(bvar) < SCIP_VARTYPE_CONTINUOUS) && SCIPisIntegral(scip, bilincoef); /*lint !e514 */
3560 SCIPdebugMessage("got different bounds for y = 0: [%g, %g] and y = 1: [%g, %g]\n", xbndszero.inf, xbndszero.sup, xbndsone.inf, xbndsone.sup);
3568 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "prod%s_%s_%s", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3586 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "%sAND%s_%s", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3606 /* product of binary variable with more than one binary or with continuous variables or with binary and user did not like AND -> replace by auxvar and linear constraints */
3610 * if all coefficients are integral, take a value that preserves integrality (-> gcd), so we can make the auxiliary variable impl. integer
3619 /* scaling by the only coefficient gives auxiliary variable = x * y, which thus will be implicit integral provided y is not continuous */
3637 /* if x-term is always negative for y = 1, negate scale so we get a positive auxiliary variable; maybe this is better sometimes? */
3641 SCIPdebugMessage("binary reformulation using scale %g, nxvars = %d, integral = %u\n", scale, nxvars, integral);
3652 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "prod%s_%s_%s", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3654 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "prod%s_%s_more_%s", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3655 SCIP_CALL( SCIPcreateVar(scip, &auxvar, name, MIN(0., SCIPintervalGetInf(xbndsone)), MAX(0., SCIPintervalGetSup(xbndsone)),
3688 * it seems to be advantageous to make the varbound constraints initial and the linear constraints not initial
3689 * maybe because it is more likely that a binary variable takes value 0 instead of 1, and thus the varbound constraints
3697 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s_%s_1", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3699 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s*more_%s_1", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3700 SCIP_CALL( SCIPcreateConsVarbound(scip, &auxcons, name, auxvar, y, -SCIPintervalGetInf(xbndsone), 0.0, SCIPinfinity(scip),
3715 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s_%s_2", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3717 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s*more_%s_2", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3718 SCIP_CALL( SCIPcreateConsVarbound(scip, &auxcons, name, auxvar, y, -SCIPintervalGetSup(xbndsone), -SCIPinfinity(scip), 0.0,
3737 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s_%s_3", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3739 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s*more_%s_3", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3740 SCIP_CALL( SCIPcreateConsLinear(scip, &auxcons, name, nxvars+2, xvars, xcoef, SCIPintervalGetInf(xbndszero), SCIPinfinity(scip),
3755 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s_%s_4", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3757 (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "linreform%s*%s*more_%s_4", SCIPvarGetName(y), SCIPvarGetName(xvars[0]), SCIPconsGetName(cons));
3758 SCIP_CALL( SCIPcreateConsLinear(scip, &auxcons, name, nxvars+2, xvars, xcoef, -SCIPinfinity(scip), SCIPintervalGetSup(xbndszero),
3791 /** tries to automatically convert a quadratic constraint (or a part of it) into a more specific and more specialized constraint */
3799 int* naddconss /**< buffer to increase with number of additional constraints created during upgrade */
3920 integral = integral && SCIPisRelEQ(scip, lb, ub) && SCIPisIntegral(scip, lincoef * lb + quadcoef * lb * lb);
3933 if( SCIPvarGetType(consdata->bilinterms[i].var1) < SCIP_VARTYPE_CONTINUOUS && SCIPvarGetType(consdata->bilinterms[i].var2) < SCIP_VARTYPE_CONTINUOUS )
3944 SCIPdebugMessage(" binlin=%d binquad=%d intlin=%d intquad=%d impllin=%d implquad=%d contlin=%d contquad=%d integral=%u\n",
3988 /* count the first upgrade constraint as constraint upgrade and the remaining ones as added constraints */
4044 consdata->bilinterms[quadvarterm->adjbilin[i]].var2 : consdata->bilinterms[quadvarterm->adjbilin[i]].var1;
4047 SCIP_CALL( presolveDisaggregateMarkComponent(scip, consdata, othervaridx, var2component, componentnr) );
4053 /** for quadratic constraints that consists of a sum of quadratic terms, disaggregates the sum into a set of constraints by introducing auxiliary variables */
4087 /* sort quadratic variable terms here, so we can later search in it without reordering the array */
4091 * in other words, the number of components in the sparsity graph of the quadratic term matrix */
4093 SCIP_CALL( SCIPhashmapCreate(&var2component, SCIPblkmem(scip), SCIPcalcHashtableSize(consdata->nquadvars)) );
4106 * @todo we could still split the constraint into several while keeping the number of variables sharing several constraints as small as possible
4123 SCIP_CALL( SCIPcreateVar(scip, &auxvars[comp], name, -SCIPinfinity(scip), SCIPinfinity(scip), 0.0,
4126 SCIP_CALL( SCIPcreateConsQuadratic2(scip, &auxconss[comp], name, 0, NULL, NULL, 0, NULL, 0, NULL,
4130 SCIPconsIsChecked(cons), SCIPconsIsPropagated(cons), SCIPconsIsLocal(cons), SCIPconsIsModifiable(cons),
4145 SCIP_CALL( SCIPaddQuadVarQuadratic(scip, auxconss[comp], consdata->quadvarterms[i].var, consdata->quadvarterms[i].lincoef, consdata->quadvarterms[i].sqrcoef) );
4148 if( !SCIPisZero(scip, consdata->quadvarterms[i].lincoef) && ABS(consdata->quadvarterms[i].lincoef) < auxcoefs[comp] )
4150 if( !SCIPisZero(scip, consdata->quadvarterms[i].sqrcoef) && ABS(consdata->quadvarterms[i].sqrcoef) < auxcoefs[comp] )
4153 SCIPfreeBlockMemoryArray(scip, &consdata->quadvarterms[i].adjbilin, consdata->quadvarterms[i].adjbilinsize);
4162 assert(comp == (int)(size_t) SCIPhashmapGetImage(var2component, consdata->bilinterms[i].var2));
4189 SCIPdebugMessage("add %d constraints for disaggregation of quadratic constraint <%s>\n", ncomponents, SCIPconsGetName(cons));
4218 /** checks if there are bilinear terms x*y with a binary variable x and an implication x = {0,1} -> y = 0
4276 y = consdata->bilinterms[k].var1 == x ? consdata->bilinterms[k].var2 : consdata->bilinterms[k].var1;
4284 SCIPdebugMessage("remove bilinear term %g<%s><%s> from <%s> due to implication\n", consdata->bilinterms[k].coef, SCIPvarGetName(x), SCIPvarGetName(y), SCIPconsGetName(cons));
4295 * we only move the coefficient to the linear coef of y here and mark the bilinterms as not merged */
4296 SCIPdebugMessage("replace bilinear term %g<%s><%s> by %g<%s> in <%s> due to implication\n", consdata->bilinterms[k].coef, SCIPvarGetName(x), SCIPvarGetName(y), consdata->bilinterms[k].coef, SCIPvarGetName(y), SCIPconsGetName(cons));
4327 /** checks a quadratic constraint for convexity and/or concavity without checking multivariate functions */
4333 SCIP_Bool checkmultivariate /**< whether curvature will be checked later on for multivariate functions */
4352 SCIPdebugMessage("Checking curvature of constraint <%s> without multivariate functions\n", SCIPconsGetName(cons));
4376 consdata->isconvex = consdata->isconvex && !SCIPisNegative(scip, consdata->quadvarterms[v].sqrcoef);
4377 consdata->isconcave = consdata->isconcave && !SCIPisPositive(scip, consdata->quadvarterms[v].sqrcoef);
4397 SCIP_Bool checkmultivariate /**< whether curvature should also be checked for multivariate functions */
4429 SCIPdebugMessage("Checking curvature of constraint <%s> with multivariate functions\n", SCIPconsGetName(cons));
4438 4 * consdata->quadvarterms[0].sqrcoef * consdata->quadvarterms[1].sqrcoef >= consdata->bilinterms[0].coef * consdata->bilinterms[0].coef;
4442 4 * consdata->quadvarterms[0].sqrcoef * consdata->quadvarterms[1].sqrcoef >= consdata->bilinterms[0].coef * consdata->bilinterms[0].coef;
4499 SCIPwarningMessage(scip, "Failed to compute eigenvalues of quadratic coefficient matrix of constraint %s. Assuming matrix is indefinite.\n", SCIPconsGetName(cons));
4505 /* deconvexification reformulates a stricly convex quadratic function in binaries such that it becomes not-strictly convex
4507 * the result is still a convex form "but less so" (ref. papers by Guignard et.al.), but with hopefully tighter value for the continuous relaxation
4511 printf("cons <%s>[%g,%g] spectrum = [%g,%g]\n", SCIPconsGetName(cons), consdata->lhs, consdata->rhs, alleigval[0], alleigval[n-1]);
4524 printf("deconvexify cons <%s> by shifting hessian by %g\n", SCIPconsGetName(cons), alleigval[0]);
4534 printf("deconcavify cons <%s> by shifting hessian by %g\n", SCIPconsGetName(cons), alleigval[n-1]);
4550 consdata->iscurvchecked = TRUE; /* set to TRUE since it does not help to repeat this procedure again and again (that will not bring Ipopt in) */
4559 /** check whether indefinite constraint function is factorable and store corresponding coefficients */
4594 * compute an eigenvalue factorization of A and check if there are one positive and one negative eigenvalue
4595 * if so, then let sigma1^2 and -sigma2^2 be these eigenvalues and v1 and v2 be the first two rows of the inverse eigenvector matrix
4598 * we then store sigma1 v1^T - sigma2 v2^T as left factor coef, and sigma1 v1^T + sigma2 v2^T as right factor coef
4601 /* if we already know that there are only positive or only negative eigenvalues, then don't try */
4646 SCIPdebugMessage("Failed to compute eigenvalues and eigenvectors of augmented quadratic form matrix for constraint <%s>.\n", SCIPconsGetName(cons));
4672 SCIPdebugMessage("Augmented quadratic form of constraint <%s> is not factorable.\n", SCIPconsGetName(cons));
4697 SCIPdebugMessage("constraint <%s> has factorable quadratic form: (%g", SCIPconsGetName(cons), consdata->factorleft[n-1]);
4701 SCIPdebugPrintf(" %+g<%s>", consdata->factorleft[i], SCIPvarGetName(consdata->quadvarterms[i].var));
4707 SCIPdebugPrintf(" %+g<%s>", consdata->factorright[i], SCIPvarGetName(consdata->quadvarterms[i].var));
4743 assert(SCIPisRelEQ(scip, consdata->factorleft[i] * consdata->factorright[j] + consdata->factorleft[j] * consdata->factorright[i], 2*a[j*n+i]));
4841 /* @todo Take better care of variables at +/- infinity: e.g., run instance waste in debug mode with a short timelimit (30s). */
4857 /* project onto local box, in case the LP solution is slightly outside the bounds (which is not our job to enforce) */
4860 #if 0 /* with non-initial columns, this might fail because variables can shortly be a column variable before entering the LP and have value 0.0 in this case */
4884 /* project onto local box, in case the LP solution is slightly outside the bounds (which is not our job to enforce) */
4887 #if 0 /* with non-initial columns, this might fail because variables can shortly be a column variable before entering the LP and have value 0.0 in this case */
4894 consdata->activity += (consdata->quadvarterms[j].lincoef + consdata->quadvarterms[j].sqrcoef * varval) * varval;
4904 /* project onto local box, in case the LP solution is slightly outside the bounds (which is not our job to enforce) */
4907 #if 0 /* with non-initial columns, this might fail because variables can shortly be a column variable before entering the LP and have value 0.0 in this case */
4913 #if 0 /* with non-initial columns, this might fail because variables can shortly be a column variable before entering the LP and have value 0.0 in this case */
4982 SCIP_CONS** maxviolcon /**< buffer to store constraint with largest violation, or NULL if solution is feasible */
5019 /** tries to compute cut for multleft * <coefleft, x'> * multright <= rhs / (multright * <coefright, x'>) where x'=(x,1) */
5054 * 1 / multright*<coefright, x'> <= 1/minact + 1/maxact - 1/(minact * maxact) multright*<coefright, x'>
5058 * multright/rightminactivity + multright/rightmaxactivity - multright/(rightminactivity * rightmaxactivity) *<coefright, x'>
5061 * rhs / (multright * <coefright, x'>) <= rhs * multright * (1/rightminactivity + 1/rightmaxactivity - 1/(rightminactivity * rightmaxactivity) * <coefright, x'>)
5075 constant += rhs * multright * coefright[consdata->nquadvars] / (rightminactivity * rightmaxactivity);
5083 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_factorablesecant_%d", SCIPconsGetName(cons), SCIPgetNLPs(scip));
5091 * <= rhs / (multright * <coefright, ref'>) - rhs / (multright * <coefright, ref'>)^2 * (multright * <coefright, x'> - multright * <coefright, ref'>)
5092 * = 2*rhs / (multright * <coefright, ref'>) - rhs / (multright * <coefright, ref'>)^2 * (multright * <coefright, x'>)
5102 /* should not happen, since we checked activity of <coefright,x> before, and assume ref within bounds */
5115 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_factorablelinearization_%d", SCIPconsGetName(cons), SCIPgetNLPs(scip));
5125 /** tries to generate a cut if constraint quadratic function is factorable and there are no linear variables
5126 * (ax+b)(cx+d) <= rhs and cx+d >= 0 -> (ax+b) <= rhs / (cx+d), where the right hand side is concave and can be linearized
5219 rightminactivity += consdata->factorright[i] * SCIPvarGetLbLocal(consdata->quadvarterms[i].var);
5226 rightminactivity += consdata->factorright[i] * SCIPvarGetUbLocal(consdata->quadvarterms[i].var);
5236 rightmaxactivity += consdata->factorright[i] * SCIPvarGetUbLocal(consdata->quadvarterms[i].var);
5243 rightmaxactivity += consdata->factorright[i] * SCIPvarGetLbLocal(consdata->quadvarterms[i].var);
5269 if( !SCIPisFeasPositive(scip, rightminactivity) && !SCIPisFeasNegative(scip, rightmaxactivity) )
5277 /* write violated constraint as multleft * factorleft * multright * (multright * factorright) <= rhs
5286 generateCutFactorableDo(scip, cons, ref, multleft, consdata->factorleft, multright, consdata->factorright, rightminactivity, rightmaxactivity, rhs, cutcoef, cutrhs, islocal, success, name);
5288 else if( !SCIPisFeasPositive(scip, rightminactivity) && !SCIPisFeasNegative(scip, rightmaxactivity) )
5295 /* write violated constraint as multleft * factorright * multright * (multright * factorleft) <= rhs
5304 generateCutFactorableDo(scip, cons, ref, multleft, consdata->factorright, multright, consdata->factorleft, leftminactivity, leftmaxactivity, rhs, cutcoef, cutrhs, islocal, success, name);
5307 (!SCIPisInfinity(scip, -rightminactivity) && !SCIPisInfinity(scip, rightmaxactivity) && rightmaxactivity - rightminactivity < leftmaxactivity - leftminactivity) )
5309 /* both factors are bounded away from 0, but the right one has a smaller activity range, so divide by that one */
5311 /* write violated constraint as multleft * factorleft * multright * (multright * factorright) <= rhs
5320 generateCutFactorableDo(scip, cons, ref, multleft, consdata->factorleft, multright, consdata->factorright, rightminactivity, rightmaxactivity, rhs, cutcoef, cutrhs, islocal, success, name);
5324 /* both factors are bounded away from 0, but the left one has a smaller activity range, so divide by that one */
5326 /* write violated constraint as multleft * factorright * multright * (multright * factorleft) <= rhs
5335 generateCutFactorableDo(scip, cons, ref, multleft, consdata->factorright, multright, consdata->factorleft, leftminactivity, leftmaxactivity, rhs, cutcoef, cutrhs, islocal, success, name);
5341 /* finds intersections of a parametric line (x,y) = (x0,y0) + t [(x1,y1) - (x0,y0)] on curves x*y = wl and x*y = wu
5461 * such that intersecting it with one of them (the first if whichuse is FALSE, the second otherwise)
5520 * The code is an adaptation of the methods in exprMul-upperHull.cpp in Couenne/stable/0.4 rev773,
5597 /* both variables have mixed sign (xl < 0 && xu > 0 && yl < 0 && yu > 0) and both xl*yl and xu*yu are feasible
5617 /* both variables have mixed sign (xl < 0 && xu > 0 && yl < 0 && yu > 0) and both xl*yu and xu*yl are feasible
5662 /* project refpoint into box not only for numerical reasons, but also due to preliminary bound tightening above */
5686 if( generateCutLTIfindIntersection(scip, 0.0, 0.0, x0, y0_, wl, wu, &xlow, &ylow, &xupp, &yupp) )
5693 /* Case 1: If both are outside of bounding box, either NW or SE, then McCormick is sufficient, so return */
5697 /* There will be at least one cut. Define coefficients and rhs ---will have to change them back if (flipX || flipY) */
5727 if( generateCutLTIfindIntersection(scip, xupp, yupp, x0, y0_, wl, SCIP_INVALID, &xlow, &ylow, NULL, NULL) )
5761 if( generateCutLTIfindIntersection(scip, xlow, ylow, x0, y0_, SCIP_INVALID, wu, NULL, NULL, &xupp, &yupp) )
5812 /* Nothing to find. Just separate two inequalities at the same point, just using different support */
5830 /* find the intersection on the lower (upper) curve on the line through xLP and the upper (lower) point
5833 if( generateCutLTIfindIntersection(scip, xlow, ylow, x0, y0_, SCIP_INVALID, wu, NULL, NULL, &xupp2, &yupp2) ||
5836 if( generateCutLTIfindIntersection(scip, xupp, yupp, x0, y0_, wl, SCIP_INVALID, &xlow2, &ylow2, NULL, NULL) ||
5877 /** tries to generate a cut if constraint quadratic function is factorable and there are linear variables
5879 * Belotti, Miller, Namazifar, Lifted inequalities for bounded products of variables, SIAG/OPT Views-and-News 22:1, 2011
5888 SCIP_Real** cutcoeflin, /**< buffer to store pointer to array with coefficients for linear variables */
5922 /* currently only separate LP solution or solutions given as SCIP_SOL, i.e., no cutgeneration during initlp */
6047 rightminactivity += consdata->factorright[i] * SCIPvarGetLbLocal(consdata->quadvarterms[i].var);
6054 rightminactivity += consdata->factorright[i] * SCIPvarGetUbLocal(consdata->quadvarterms[i].var);
6064 rightmaxactivity += consdata->factorright[i] * SCIPvarGetUbLocal(consdata->quadvarterms[i].var);
6071 rightmaxactivity += consdata->factorright[i] * SCIPvarGetLbLocal(consdata->quadvarterms[i].var);
6077 /* if any of the factors is essentially fixed, give up and do usual method (numerically less sensitive, I hope) */
6078 if( SCIPisRelEQ(scip, leftminactivity, leftmaxactivity) || SCIPisRelEQ(scip, rightminactivity, rightmaxactivity) )
6096 leftminactivity, leftmaxactivity, rightminactivity, rightmaxactivity, rhsminactivity, rhsmaxactivity,
6101 if( coefleft * leftrefactivity + coefright * rightrefactivity + coefrhs * rhsrefactivity >= *cutlhs )
6109 * coefleft * leftfactor + coefright * rightfactor + coefrhs * w >= cutlhs, where conslhs - lincoefs <= w <= consrhs - lincoefs
6135 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_lti_%d", SCIPconsGetName(cons), SCIPgetNLPs(scip));
6148 SCIP_Bool isint, /**< whether corresponding variable is a discrete variable, and thus linearization could be moved */
6151 SCIP_Bool* success /**< buffer to set to FALSE if linearization has failed due to large numbers */
6188 /* sqrcoef * x^2 -> secant between f=floor(refpoint) and f+1 = sqrcoef * (f^2 + ((f+1)^2 - f^2) * (x-f)) = sqrcoef * (-f*(f+1) + (2*f+1)*x) */
6219 SCIP_Bool* success /**< buffer to set to FALSE if secant has failed due to large numbers or unboundedness */
6245 /* sqrcoef * x^2 -> sqrcoef * (lb * lb + (ub*ub - lb*lb)/(ub-lb) * (x-lb)) = sqrcoef * (lb*lb + (ub+lb)*(x-lb)) = sqrcoef * ((lb+ub)*x - lb*ub) */
6269 SCIP_Bool* success /**< buffer to set to FALSE if linearization has failed due to large numbers */
6289 /* bilincoef * x * y -> bilincoef * (refpointx * refpointy + refpointy * (x - refpointx) + refpointx * (y - refpointy))
6295 if( SCIPisInfinity(scip, REALABS(bilincoef * refpointx)) || SCIPisInfinity(scip, REALABS(bilincoef * refpointy)) || SCIPisInfinity(scip, REALABS(constant)) )
6317 SCIP_Bool overestimate, /**< whether to compute an overestimator instead of an underestimator */
6321 SCIP_Bool* success /**< buffer to set to FALSE if linearization has failed due to large numbers */
6376 (SCIPisInfinity(scip, ubx) || SCIPisInfinity(scip, uby) || (uby - refpointy) * (ubx - refpointx) >= (refpointy - lby) * (refpointx - lbx)) )
6380 /* x*y = lbx * y + (x-lbx) * y >= lbx * y + (x-lbx) * lby >= lbx * y + min{(ubx-lbx) * lby, 0 * lby} */
6387 /* x*y = lby * x + (y-lby) * x >= lby * x + (y-lby) * lbx >= lby * x + min{(uby-lby) * lbx, 0 * lbx} */
6403 /* x*y = ubx * y + (x-ubx) * y >= ubx * y + (x-ubx) * uby >= ubx * y + min{(lbx-ubx) * uby, 0 * uby} */
6410 /* x*y = uby * x + (y-uby) * x >= uby * x + (y-uby) * ubx >= uby * x + min{(lby-uby) * ubx, 0 * ubx} */
6432 (SCIPisInfinity(scip, -lbx) || SCIPisInfinity(scip, uby) || (ubx - lbx) * (refpointy - lby) <= (uby - lby) * (refpointx - lbx)) )
6436 /* x*y = ubx * y + (x-ubx) * y <= ubx * y + (x-ubx) * lby <= ubx * y + max{(lbx-ubx) * lby, 0 * lby} */
6443 /* x*y = lby * x + (y-lby) * x <= lby * x + (y-lby) * ubx <= lby * x + max{(uby-lby) * ubx, 0 * ubx} */
6459 /* x*y = lbx * y + (x-lbx) * y <= lbx * y + (x-lbx) * uby <= lbx * y + max{(ubx-lbx) * uby, 0 * uby} */
6466 /* x*y = uby * x + (y-uby) * x <= uby * x + (y-uby) * lbx <= uby * x + max{(lby-uby) * lbx, 0 * lbx} */
6485 if( SCIPisInfinity(scip, REALABS(coefx)) || SCIPisInfinity(scip, REALABS(coefy)) || SCIPisInfinity(scip, REALABS(constant)) )
6498 /* SCIPdebugMessage("%.20g * x[%.20g,%.20g] * y[%.20g,%.20g] %c= %.20g * x %+.20g * y %+.20g\n", bilincoef, lbx, ubx, lby, uby, overestimate ? '<' : '>', coefx, coefy, constant); */
6552 addSquareLinearization(scip, consdata->quadvarterms[j].sqrcoef, ref[j], consdata->quadvarterms[j].nadjbilin == 0 && SCIPvarGetType(var) < SCIP_VARTYPE_CONTINUOUS, &coef[j], &constant, success);
6568 addBilinLinearization(scip, bilinterm->coef, ref[j], ref[var2pos], &coef[j], &coef[var2pos], &constant, success);
6574 SCIPdebugMessage("no success in linearization of <%s> in reference point\n", SCIPconsGetName(cons));
6589 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_side%d_linearization_%d", SCIPconsGetName(cons), violside, SCIPgetNLPs(scip));
6646 if( (violside == SCIP_SIDETYPE_LEFT && sqrcoef <= 0.0) || (violside == SCIP_SIDETYPE_RIGHT && sqrcoef > 0.0) )
6649 addSquareLinearization(scip, sqrcoef, ref[j], SCIPvarGetType(var) < SCIP_VARTYPE_CONTINUOUS, &coef[j], &constant, success);
6654 addSquareSecant(scip, sqrcoef, SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), ref[j], &coef[j], &constant, success);
6695 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_side%d_estimation_%d", SCIPconsGetName(cons), violside, SCIPgetNLPs(scip));
6700 /** generates a cut based on linearization (if convex) or McCormick (if nonconvex) in a given reference point */
6710 SCIP_Real* efficacy, /**< buffer to store efficacy of row in reference solution, or NULL if not of interest */
6711 SCIP_Bool checkcurvmultivar, /**< are we allowed to check the curvature of a multivariate quadratic function, if not done yet */
6712 SCIP_Real minefficacy /**< minimal required efficacy (violation scaled by maximal absolute coefficient) */
6763 SCIP_CALL( generateCutFactorable(scip, cons, violside, ref, coef, &lhs, &rhs, &islocal, &success, cutname) );
6769 /* generateCutLTI needs reference values also for the linear variables, which we only have if sol is given or LP has been solved */
6770 SCIP_CALL( generateCutLTI(scip, cons, violside, ref, sol, &lincoefs, coef, &lhs, &rhs, &islocal, &success, cutname) );
6772 /* in case of LTI cuts, we have to recompute the min and max of lincoefs, since they may have been modified */
6788 if( (violside == SCIP_SIDETYPE_LEFT && consdata->isconcave) || (violside == SCIP_SIDETYPE_RIGHT && consdata->isconvex) )
6790 SCIP_CALL( generateCutConvex(scip, cons, violside, ref, coef, &lhs, &rhs, &islocal, &success, cutname) );
6794 SCIP_CALL( generateCutNonConvex(scip, cons, violside, ref, coef, &lhs, &rhs, &islocal, &success, cutname) );
6817 * it is required if we need to check or return cut efficacy, for some debug output below, and some assert
6818 * round almost integral coefficients in integers, since this will happen when adding coefs to row (see comments below)
6826 * HOWEVER, they become column variables when they are added to a row (via SCIPaddVarsToRow below).
6828 * So when calculating the row activity, we treat loose variable as if they were already column variables.
6831 refactivitylinpart += (SCIPisIntegral(scip, lincoefs[j]) ? SCIPround(scip, lincoefs[j]) : lincoefs[j]) * SCIPgetSolVal(scip, sol, consdata->linvars[j]);
6848 * thus, we anticipate by rounding coef here, but also modify constant so that cut is still valid (if possible)
6871 SCIPdebugMessage("var <%s> [%g,%g] has almost integral coef %.20g, round coefficient to %g and add constant %g\n",
6872 SCIPvarGetName(var), SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var), coef[j], roundcoef, (coef[j]-roundcoef) * xbnd);
6881 /* As above: When calculating the row activity, we treat loose variable as if they were already column variables. */
6905 SCIPdebugMessage("skip cut for constraint <%s> since all coefficients are zero and it's always satisfied\n", SCIPconsGetName(cons));
6918 SCIPdebugMessage("cut coefficients for constraint <%s> have very large range: mincoef = %g maxcoef = %g\n", SCIPconsGetName(cons), mincoef, maxcoef);
6923 * since we use local bounds, we need to make the row local if they are different from their global counterpart
6925 if( ((coef[mincoefidx] > 0.0 && !SCIPisInfinity(scip, rhs)) || (coef[mincoefidx] < 0.0 && !SCIPisInfinity(scip, -lhs))) &&
6928 SCIPdebugMessage("eliminate coefficient %g for <%s> [%g, %g]\n", coef[mincoefidx], SCIPvarGetName(var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var));
6934 else if( ((coef[mincoefidx] < 0.0 && !SCIPisInfinity(scip, rhs)) || (coef[mincoefidx] > 0.0 && !SCIPisInfinity(scip, -lhs))) &&
6937 SCIPdebugMessage("eliminate coefficient %g for <%s> [%g, %g]\n", coef[mincoefidx], SCIPvarGetName(var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var));
6966 SCIPdebugMessage("skip cut for constraint <%s> because both sides are not finite: lhs = %g, rhs = %g\n", SCIPconsGetName(cons), lhs, rhs);
6980 /* in difference to SCIPgetCutEfficacy, we scale by norm only if the norm is > 1.0 this avoid finding cuts
6981 * efficient which are only very slightly violated CPLEX does not seem to scale row coefficients up too also we
6982 * use infinity norm, since that seem to be the usual scaling strategy in LP solvers (equilibrium scaling) */
7009 if( success && !SCIPisInfinity(scip, -minefficacy) && rowefficacy < minefficacy ) /*lint !e644*/
7011 SCIPdebugMessage("skip cut for constraint <%s> because efficacy %g too low (< %g)\n", SCIPconsGetName(cons), rowefficacy, minefficacy);
7018 SCIP_CALL( SCIPcreateEmptyRowCons(scip, row, SCIPconsGetHdlr(cons), cutname, lhs, rhs, islocal && (SCIPgetDepth(scip) > 0), FALSE, TRUE) );
7027 SCIPdebugMessage("found cut <%s>, lhs=%g, rhs=%g, mincoef=%g, maxcoef=%g, range=%g, nnz=%d, violation=%g, efficacy=%g\n",
7037 * computing efficacy w.r.t. the LP solution makes only sense if the LP was solved to optimality (see bug 612)
7039 * disabled these asserts as they can fail due to numerical reasons (cancelation when substracting big numbers),
7040 * as the order in which we add up the activity for the single terms can be different than the one that lp.c uses
7044 assert((conshdlrdata->scaling != 'g') || (SCIPisFeasPositive(scip, rowefficacy) == SCIPisFeasPositive(scip, -SCIPgetRowSolFeasibility(scip, *row, sol)/MAX(1.0,SCIPgetRowMaxCoef(scip, *row)))));
7045 assert((conshdlrdata->scaling != 's') || (SCIPisFeasPositive(scip, rowefficacy) == SCIPisFeasPositive(scip, -SCIPgetRowSolFeasibility(scip, *row, sol)/MAX(1.0,MIN(REALABS(lhs),REALABS(rhs))))));
7046 assert((conshdlrdata->scaling != 'o') || (SCIPisFeasPositive(scip, rowefficacy) == SCIPisFeasPositive(scip, -SCIPgetRowSolFeasibility(scip, *row, sol))));
7053 /* if coefficients for linear variables are different than those in constraint, then free array */
7062 /** generates a cut based on linearization (if convex) or McCormick (if nonconvex) in a solution
7073 SCIP_Real* efficacy, /**< buffer to store efficacy of row in reference solution, or NULL if not of interest */
7074 SCIP_Bool checkcurvmultivar, /**< are we allowed to check the curvature of a multivariate quadratic function, if not done yet */
7075 SCIP_Real minefficacy /**< minimal required efficacy (violation scaled by maximal absolute coefficient) */
7110 SCIP_CALL( generateCut(scip, conshdlr, cons, ref, sol, violside, row, efficacy, checkcurvmultivar, minefficacy) );
7119 * for nonconvex functions, we just call generateCutSol with the unbounded solution as reference point */
7127 SCIP_Real* rowrayprod, /**< buffer to store product of ray with row coefficients, or NULL if not of interest */
7128 SCIP_Bool checkcurvmultivar /**< are we allowed to check the curvature of a multivariate quadratic function, if not done yet */
7164 SCIP_CALL( generateCutSol(scip, conshdlr, cons, NULL, NULL, violside, row, NULL, FALSE, -SCIPinfinity(scip)) );
7183 /* we seek for a linearization of the quadratic function such that it intersects with the unbounded ray
7184 * that is, we need a reference point ref such that for the gradient g of xAx+bx in ref, we have
7187 * initially, for finite rhs, we set ref_i = 1.0 if (A*ray)_i > 0.0 and ref_i = -1.0 if (A*ray)_i < 0.0 (for finite lhs analog)
7188 * <ref, 2*A*ray> + <b,ray> is sufficiently larger 0.0, we call generateCut for this point, otherwise, we scale up ref
7204 matrixrayprod += bilinterm->coef * SCIPgetPrimalRayVal(scip, bilinterm->var1 == var ? bilinterm->var2 : bilinterm->var1);
7217 assert((violside == SCIP_SIDETYPE_RIGHT && quadrayprod >= 0.0) || (violside == SCIP_SIDETYPE_LEFT && quadrayprod <= 0.0));
7230 SCIPdebugMessage("initially have <b,ray> = %g and <ref, 2*A*ref> = %g\n", linrayprod, quadrayprod);
7232 /* we scale the refpoint up, such that <ref, 2*A*ray> >= -2*<b, ray> (rhs finite) or <ref, 2*A*ray> <= -2*<b, ray> (lhs finite), if <b,ray> is not zero
7235 if( (!SCIPisZero(scip, linrayprod) && violside == SCIP_SIDETYPE_RIGHT && quadrayprod < -2*linrayprod) ||
7236 ( !SCIPisZero(scip, linrayprod) && violside == SCIP_SIDETYPE_LEFT && quadrayprod > -2*linrayprod) ||
7256 SCIP_CALL( generateCut(scip, conshdlr, cons, ref, NULL, violside, row, NULL, FALSE, -SCIPinfinity(scip)) );
7278 SCIP_Real* bestefficacy /**< buffer to store best efficacy of a cut that was added to the LP, if found; or NULL if not of interest */
7309 if( SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) || SCIPisGT(scip, consdata->rhsviol, SCIPfeastol(scip)) )
7315 violside = SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) ? SCIP_SIDETYPE_LEFT : SCIP_SIDETYPE_RIGHT;
7318 actminefficacy = inenforcement && ((violside == SCIP_SIDETYPE_RIGHT && consdata->isconvex ) || (violside == SCIP_SIDETYPE_LEFT && consdata->isconcave))
7325 /* if the LP is unbounded, then we need a cut that cuts into the direction of a hopefully existing primal ray
7326 * that is, assume a ray r is given such that p + t*r is feasible for the LP for all t >= t_0 and some p
7327 * given a cut lhs <= <c,x> <= rhs, we check whether it imposes an upper bound on t and thus bounds the ray
7328 * this is given if rhs < infinity and <c,r> > 0, since we then enforce <c,p+t*r> = <c,p> + t<c,r> <= rhs, i.e., t <= (rhs - <c,p>)/<c,r>
7335 SCIP_CALL( generateCutUnboundedLP(scip, conshdlr, conss[c], violside, &row, &rayprod, conshdlrdata->checkcurvature) );
7352 /* in difference to SCIPgetCutEfficacy, we scale by norm only if the norm is > 1.0 this avoid finding
7354 * coefficients up too also we use infinity norm, since that seem to be the usual scaling strategy
7379 /* @todo If convex, can we easily move the refpoint closer to the feasible region to get a stronger cut?
7382 SCIP_CALL( generateCutSol(scip, conshdlr, conss[c], sol, NULL, violside, &row, &efficacy, conshdlrdata->checkcurvature, actminefficacy) );
7383 /* @todo If generation failed not because of low efficacy, then probably because of numerical issues;
7384 * if the constraint is convex and we are desperate to get a cut, then we may try again with a better chosen reference point
7391 if( SCIPisGT(scip, efficacy, actminefficacy) && SCIPisCutApplicable(scip, row) ) /*lint !e644 */
7399 SCIPdebugMessage("cut for constraint <%s> is infeasible -> cutoff.\n", SCIPconsGetName(conss[c]));
7424 * others are only checked and enforced if we are still feasible or have not found a separating cut yet
7433 /** adds linearizations cuts for convex constraints w.r.t. a given reference point to cutpool and sepastore
7434 * if separatedlpsol is not NULL, then a cut that separates the LP solution is added to the sepastore and is forced to enter the LP
7435 * if separatedlpsol is not NULL, but cut does not separate the LP solution, then it is added to the cutpool only
7445 SCIP_Bool* separatedlpsol, /**< buffer to store whether a cut that separates the current LP solution was found and added to LP, or NULL if adding to cutpool only */
7446 SCIP_Real minefficacy /**< minimal efficacy of a cut when checking for separation of LP solution */
7479 SCIP_CALL( generateCutSol(scip, conshdlr, conss[c], NULL, ref, SCIP_SIDETYPE_RIGHT, &row, NULL, conshdlrdata->checkcurvature, -SCIPinfinity(scip)) ); /*lint !e613 */
7483 SCIP_CALL( generateCutSol(scip, conshdlr, conss[c], NULL, ref, SCIP_SIDETYPE_LEFT, &row, NULL, conshdlrdata->checkcurvature, -SCIPinfinity(scip)) ); /*lint !e613 */
7506 /* in difference to SCIPgetCutEfficacy, we scale by norm only if the norm is > 1.0 this avoid finding cuts
7507 * efficient which are only very slightly violated CPLEX does not seem to scale row coefficients up too
7508 * also we use infinity norm, since that seem to be the usual scaling strategy in LP solvers (equilibrium
7538 SCIPdebugMessage("added linearization cut <%s> to LP, efficacy = %g\n", SCIProwGetName(row), efficacy);
7584 /* we are only interested in solution coming from some heuristic other than trysol, but not from the tree
7585 * the reason for ignoring trysol solutions is that they may come from an NLP solve in sepalp, where we already added linearizations,
7594 SCIPdebugMessage("caught new sol event %x from heur <%s>; have %d conss\n", SCIPeventGetType(event), SCIPheurGetName(SCIPsolGetHeur(sol)), nconss);
7601 /** computes the infeasibilities of variables from the convexification gaps in the constraints and notifies the branching rule about them
7650 SCIPdebugMessage("cons %s violation: %g %g convex: %u %u\n", SCIPconsGetName(conss[c]), consdata->lhsviol, consdata->rhsviol, consdata->isconvex, consdata->isconcave);
7656 if( (SCIPisGT(scip, consdata->rhsviol, SCIPfeastol(scip)) && consdata->quadvarterms[j].sqrcoef < 0) ||
7657 ( SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) && consdata->quadvarterms[j].sqrcoef > 0) )
7663 SCIPdebugMessage("ignore fixed variable <%s>[%g, %g], diff %g\n", SCIPvarGetName(x), xlb, xub, xub-xlb);
7686 /* if any of the variables if fixed, then it actually behaves like a linear term, so we don't need to branch on it */
7709 /* if both variables are at one of its bounds, then no need to branch, since McCormick is exact there */
7717 coef_ = SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) ? -consdata->bilinterms[j].coef : consdata->bilinterms[j].coef;
7738 /* if one of the variables is binary or integral with domain width 1, then branching on this makes the term linear, so prefer this */
7768 /* if both variables are integral, prefer the one with the smaller domain, so variable gets fixed soon
7769 * does not seem to work well on tln instances, so disable for now and may look at it later again
7790 * this is, because after branching both variables will be one the bounds, and McCormick will be exact then */
7809 /** registers a quadratic variable from a violated constraint as branching candidate that has a large absolute value in the LP relaxation */
7836 if( !SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) && !SCIPisGT(scip, consdata->rhsviol, SCIPfeastol(scip)) )
7842 if( SCIPisRelEQ(scip, SCIPvarGetLbLocal(consdata->quadvarterms[i].var), SCIPvarGetUbLocal(consdata->quadvarterms[i].var)) )
7861 /** replaces violated quadratic constraints where all quadratic variables are fixed by linear constraints */
7897 if( !SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) && !SCIPisGT(scip, consdata->rhsviol, SCIPfeastol(scip)) )
7905 assert(SCIPisRelEQ(scip, SCIPvarGetLbLocal(consdata->quadvarterms[i].var), SCIPvarGetUbLocal(consdata->quadvarterms[i].var)));
7907 val1 = (SCIPvarGetUbLocal(consdata->quadvarterms[i].var) + SCIPvarGetLbLocal(consdata->quadvarterms[i].var)) / 2.0;
7908 constant += (consdata->quadvarterms[i].lincoef + consdata->quadvarterms[i].sqrcoef * val1) * val1;
7913 val1 = (SCIPvarGetUbLocal(consdata->bilinterms[i].var1) + SCIPvarGetLbLocal(consdata->bilinterms[i].var1)) / 2.0;
7914 val2 = (SCIPvarGetUbLocal(consdata->bilinterms[i].var2) + SCIPvarGetLbLocal(consdata->bilinterms[i].var2)) / 2.0;
7939 SCIPdebugMessage("Linear constraint with one variable: %g <= %g <%s> <= %g\n", lhs, coef, SCIPvarGetName(*consdata->linvars), rhs);
7964 SCIPdebugMessage("Linear constraint is a bound: %g <= <%s> <= %g\n", lhs, SCIPvarGetName(*consdata->linvars), rhs);
8009 SCIPdebugMessage("replace quadratic constraint <%s> by linear constraint after all quadratic vars have been fixed\n", SCIPconsGetName(conss[c]) );
8074 SCIPdebugMessage("%s found constraint <%s> infeasible due to tightened lower bound %g for variable <%s>\n", SCIPinProbing(scip) ? "in probing" : "", SCIPconsGetName(cons), bnd, SCIPvarGetName(var));
8081 SCIPdebugMessage("%s tightened lower bound of variable <%s> in constraint <%s> to %g\n", SCIPinProbing(scip) ? "in probing" : "", SCIPvarGetName(var), SCIPconsGetName(cons), bnd);
8133 SCIPdebugMessage("%s found constraint <%s> infeasible due to tightened upper bound %g for variable <%s>\n", SCIPinProbing(scip) ? "in probing" : "", SCIPconsGetName(cons), bnd, SCIPvarGetName(var));
8140 SCIPdebugMessage("%s tightened upper bound of variable <%s> in constraint <%s> to %g\n", SCIPinProbing(scip) ? "in probing" : "", SCIPvarGetName(var), SCIPconsGetName(cons), bnd);
8149 /** solves a quadratic equation \f$ a x^2 + b x \in rhs \f$ (with b an interval) and reduces bounds on x or deduces infeasibility if possible
8177 SCIPdebugMessage("found <%s> infeasible due to domain propagation for quadratic variable <%s>\n", SCIPconsGetName(cons), SCIPvarGetName(var));
8201 SCIPdebugMessage("found <%s> infeasible due to domain propagation for quadratic variable <%s>\n", SCIPconsGetName(cons), SCIPvarGetName(var));
8216 /* SCIPdebugMessage("%g x^2 + [%g, %g] x in [%g, %g] -> [%g, %g]\n", a, b.inf, b.sup, rhs.inf, rhs.sup, newrange.inf, newrange.sup); */
8218 if( SCIPisInfinity(scip, SCIPintervalGetInf(newrange)) || SCIPisInfinity(scip, -SCIPintervalGetSup(newrange)) )
8220 SCIPdebugMessage("found <%s> infeasible because propagated domain of quadratic variable <%s> is outside of (-infty, +infty)\n", SCIPconsGetName(cons), SCIPvarGetName(var));
8228 SCIPdebugMessage("found <%s> infeasible due to domain propagation for quadratic variable <%s>\n", SCIPconsGetName(cons), SCIPvarGetName(var));
8235 SCIP_CALL( propagateBoundsTightenVarLb(scip, cons, intervalinfty, var, SCIPintervalGetInf(newrange), result, nchgbds) );
8242 SCIP_CALL( propagateBoundsTightenVarUb(scip, cons, intervalinfty, var, SCIPintervalGetSup(newrange), result, nchgbds) );
8250 /* the new version below computes potentially tighter bounds, but also always adds a small safety area since it is not implemented roundingsafe,
8255 /** tries to deduce domain reductions for x in xsqrcoef x^2 + xlincoef x + ysqrcoef y^2 + ylincoef y + bilincoef x y \\in rhs
8293 SCIPintervalSetBounds(&varbnds, MIN(SCIPvarGetLbLocal(y), SCIPvarGetUbLocal(y)), MAX(SCIPvarGetLbLocal(y), SCIPvarGetUbLocal(y))); /*lint !e666 */
8298 /* if rhs is unbounded by above, it is sufficient to get an upper bound on ysqrcoef*y^2 + ylincoef * y */
8311 /* if rhs is unbounded by below, it is sufficient to get a lower bound on ysqrcoef*y^2 + ylincoef * y */
8337 SCIP_CALL( propagateBoundsQuadVar(scip, cons, intervalinfty, x, xsqrcoef, lincoef, myrhs, result, nchgbds) );
8342 /** tries to deduce domain reductions for x in xsqrcoef x^2 + xlincoef x + ysqrcoef y^2 + ylincoef y + bilincoef x y \\in rhs
8379 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(x), SCIPvarGetUbLocal(x))), /*lint !e666*/
8380 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(x), SCIPvarGetUbLocal(x)))); /*lint !e666*/
8382 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(y), SCIPvarGetUbLocal(y))), /*lint !e666*/
8383 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(y), SCIPvarGetUbLocal(y)))); /*lint !e666*/
8386 SCIPintervalSolveBivariateQuadExpressionAllScalar(intervalinfty, &xbnds, xsqrcoef, ysqrcoef, bilincoef, xlincoef, ylincoef, rhs, xbnds, ybnds);
8390 SCIPdebugMessage("found <%s> infeasible due to domain propagation for quadratic variable <%s>\n", SCIPconsGetName(cons), SCIPvarGetName(x));
8397 SCIP_CALL( propagateBoundsTightenVarLb(scip, cons, intervalinfty, x, SCIPintervalGetInf(xbnds), result, nchgbds) );
8404 SCIP_CALL( propagateBoundsTightenVarUb(scip, cons, intervalinfty, x, SCIPintervalGetSup(xbnds), result, nchgbds) );
8423 SCIP_Real* minquadactivity, /**< minimal activity of quadratic variable terms where only terms with finite minimal activity contribute */
8424 SCIP_Real* maxquadactivity, /**< maximal activity of quadratic variable terms where only terms with finite maximal activity contribute */
8425 int* minactivityinf, /**< number of quadratic variables that contribute -infinity to minimal activity */
8426 int* maxactivityinf, /**< number of quadratic variables that contribute +infinity to maximal activity */
8471 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(consdata->quadvarterms[i].var), SCIPvarGetUbLocal(consdata->quadvarterms[i].var))),
8472 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(consdata->quadvarterms[i].var), SCIPvarGetUbLocal(consdata->quadvarterms[i].var))));
8482 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))),
8483 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))));
8491 bnd = SCIPintervalQuadUpperBound(intervalinfty, consdata->quadvarterms[i].sqrcoef, lincoef, xrng);
8498 /* if maximal activity is below value for -infinity, let's take -1e10 as upper bound on maximal activity
8519 bnd = -SCIPintervalQuadUpperBound(intervalinfty, -consdata->quadvarterms[i].sqrcoef, lincoef, xrng);
8527 /* if minimal activity is above value for infinity, let's take 1e10 as lower bound on minimal activity
8560 SCIP_Bool* redundant /**< buffer where to store whether constraint has been found to be redundant */
8569 int quadminactinf; /* number of quadratic variables that contribute -infinity to minimal activity of quadratic term */
8570 int quadmaxactinf; /* number of quadratic variables that contribute +infinity to maximal activity of quadratic term */
8571 SCIP_INTERVAL* quadactcontr; /* contribution of each quadratic variable term to quadactivity */
8615 /* sort quadratic variable terms, in case we need to search for terms occuring in bilinear terms later
8628 propagateBoundsGetQuadActivity(scip, consdata, intervalinfty, &minquadactivity, &maxquadactivity, &quadminactinf, &quadmaxactinf, quadactcontr);
8643 SCIPintervalSetBounds(&consactivity, consdata->minlinactivityinf > 0 ? -intervalinfty : consdata->minlinactivity, consdata->maxlinactivityinf > 0 ? intervalinfty : consdata->maxlinactivity);
8647 SCIPdebugMessage("found constraint <%s> to be redundant: sides: [%g, %g], activity: [%g, %g]\n",
8648 SCIPconsGetName(cons), consdata->lhs, consdata->rhs, SCIPintervalGetInf(consactivity), SCIPintervalGetSup(consactivity));
8653 /* was SCIPintervalAreDisjoint(consbounds, consactivity), but that would allow violations up to eps only
8656 if( SCIPisFeasGT(scip, consbounds.inf, consactivity.sup) || SCIPisFeasLT(scip, consbounds.sup, consactivity.inf) )
8658 SCIPdebugMessage("found constraint <%s> to be infeasible; sides: [%g, %g], activity: [%g, %g], infeas: %g\n",
8659 SCIPconsGetName(cons), consdata->lhs, consdata->rhs, SCIPintervalGetInf(consactivity), SCIPintervalGetSup(consactivity),
8660 MAX(consdata->lhs - SCIPintervalGetSup(consactivity), SCIPintervalGetInf(consactivity) - consdata->rhs));
8665 /* propagate linear part \in rhs = consbounds - quadactivity (use the one from consdata, since that includes infinities) */
8678 * we can't expect much more tightening, but may detect infeasiblity, but shouldn't the check on the constraints activity detect that?
8700 SCIP_CALL( propagateBoundsTightenVarUb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8714 SCIP_CALL( propagateBoundsTightenVarUb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8736 SCIP_CALL( propagateBoundsTightenVarLb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8750 SCIP_CALL( propagateBoundsTightenVarLb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8775 SCIP_CALL( propagateBoundsTightenVarUb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8789 SCIP_CALL( propagateBoundsTightenVarUb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8811 SCIP_CALL( propagateBoundsTightenVarLb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8825 SCIP_CALL( propagateBoundsTightenVarLb(scip, cons, intervalinfty, var, bnd, result, nchgbds) );
8840 consdataUpdateLinearActivity(scip, consdata, intervalinfty); /* make sure, activities of linear part did not become invalid by above bound changes, if any */
8841 assert(consdata->minlinactivityinf > 0 || consdata->maxlinactivityinf > 0 || consdata->minlinactivity <= consdata->maxlinactivity);
8859 SCIP_CALL( propagateBoundsQuadVar(scip, cons, intervalinfty, var, consdata->quadvarterms[0].sqrcoef, lincoef, rhs, result, nchgbds) );
8863 /* quadratic part is just ax*x^2+bx*x + ay*y^2+by*y + c*xy -> a common case that we treat directly */
8864 assert(consdata->bilinterms[0].var1 == consdata->quadvarterms[0].var || consdata->bilinterms[0].var1 == consdata->quadvarterms[1].var);
8865 assert(consdata->bilinterms[0].var2 == consdata->quadvarterms[0].var || consdata->bilinterms[0].var2 == consdata->quadvarterms[1].var);
8869 consdata->quadvarterms[0].var, consdata->quadvarterms[0].sqrcoef, consdata->quadvarterms[0].lincoef,
8870 consdata->quadvarterms[1].var, consdata->quadvarterms[1].sqrcoef, consdata->quadvarterms[1].lincoef,
8877 consdata->quadvarterms[1].var, consdata->quadvarterms[1].sqrcoef, consdata->quadvarterms[1].lincoef,
8878 consdata->quadvarterms[0].var, consdata->quadvarterms[0].sqrcoef, consdata->quadvarterms[0].lincoef,
8892 propagateBoundsGetQuadActivity(scip, consdata, intervalinfty, &minquadactivity, &maxquadactivity, &quadminactinf, &quadmaxactinf, quadactcontr);
8898 /* if the quad activities are not hopelessly unbounded on useful sides, try to deduce domain reductions on quad vars */
8913 * we can't expect much more tightening, but may detect infeasiblity, but shouldn't the check on the constraints activity detect that?
8920 /* setup rhs2.sup = rhs.sup - (quadactivity.inf - quadactcontr[i].inf), if everything were finite
8921 * if only quadactcontr[i].inf is infinite (i.e., the other i are all finite), we just get rhs2.sup = rhs.sup
8925 if( quadminactinf == 0 || (quadminactinf == 1 && SCIPintervalGetInf(quadactcontr[i]) <= -intervalinfty) )
8934 /* the residual quad min activity w.r.t. quad var term i is finite and nonzero, so add it to right hand side */
8942 /* there are either >= 2 quad var terms contributing -infinity, or there is one which is not i */
8954 if( quadmaxactinf == 0 || (quadmaxactinf == 1 && SCIPintervalGetSup(quadactcontr[i]) >= intervalinfty) )
8956 /* the residual quad max activity w.r.t. quad var term i is finite and nonzero, so add it to right hand side */
8971 /* there are either >= 2 quad var terms contributing infinity, or there is one which is not i */
8993 /* bilinear term k contributes to the activity of quad var term i, so just add bounds to linear coef */
8995 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))),
8996 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))));
9006 * if the bounds on bilinear term k as added to rhs2 are old due to recent bound tightening, we may not do best possible, but still correct
9007 * HOWEVER: when computing rhs2, we may not just have added the bounds for the bilinear term, but for the associated quadratic term
9009 * since we do not want to repeat a call to SCIPintervalQuad for that quadratic term with bilinear term k removed,
9010 * we only remove the bounds for the bilinear term k from rhs2 if the associated quadratic term consists only of this bilinear term,
9011 * i.e., the quadratic term corresponding to var1 should be only var1*var2, but have no square or linear coefs or other bilinear terms
9012 * (for efficiency reasons, we check here only if there are any other bilinear terms than var1*var2 associated with var1, even if they are not associated with the quad var term for var1)
9025 if( (consdata->quadvarterms[otherpos].sqrcoef != 0.0) || consdata->quadvarterms[otherpos].lincoef != 0.0 || consdata->quadvarterms[otherpos].nadjbilin > 1 )
9030 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(consdata->bilinterms[k].var1), SCIPvarGetUbLocal(consdata->bilinterms[k].var1))),
9031 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(consdata->bilinterms[k].var1), SCIPvarGetUbLocal(consdata->bilinterms[k].var1))));
9036 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))),
9037 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))));
9073 -infty2infty(SCIPinfinity(scip), intervalinfty, -MIN(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))),
9074 +infty2infty(SCIPinfinity(scip), intervalinfty, MAX(SCIPvarGetLbLocal(consdata->bilinterms[k].var2), SCIPvarGetUbLocal(consdata->bilinterms[k].var2))));
9081 SCIP_CALL( propagateBoundsQuadVar(scip, cons, intervalinfty, var, consdata->quadvarterms[i].sqrcoef, lincoef, rhs2, result, nchgbds) );
9120 assert(SCIPgetStage(scip) == SCIP_STAGE_PRESOLVING || SCIPgetStage(scip) == SCIP_STAGE_SOLVING);
9137 SCIPdebugMessage("starting domain propagation round %d of %d for %d constraints\n", roundnr, maxproprounds, nconss);
9164 /* checks for a linear variable that can be increase or decreased without harming feasibility */
9178 /* check for a linear variable that can be increase or decreased without harming feasibility */
9197 /* if we have already one candidate, then take the one where the loss in the objective function is less */
9199 (SCIPvarGetObj(consdata->linvars[consdata->linvar_maydecrease]) / consdata->lincoefs[consdata->linvar_maydecrease] > SCIPvarGetObj(consdata->linvars[i]) / consdata->lincoefs[i]) )
9206 /* if we have already one candidate, then take the one where the loss in the objective function is less */
9208 (SCIPvarGetObj(consdata->linvars[consdata->linvar_mayincrease]) / consdata->lincoefs[consdata->linvar_mayincrease] > SCIPvarGetObj(consdata->linvars[i]) / consdata->lincoefs[i]) )
9216 SCIPdebugMessage("may increase <%s> to become feasible\n", SCIPvarGetName(consdata->linvars[consdata->linvar_mayincrease]));
9220 SCIPdebugMessage("may decrease <%s> to become feasible\n", SCIPvarGetName(consdata->linvars[consdata->linvar_maydecrease]));
9225 /** Given a solution where every quadratic constraint is either feasible or can be made feasible by
9226 * moving a linear variable, construct the corresponding feasible solution and pass it to the trysol heuristic.
9227 * The method assumes that this is always possible and that not all constraints are feasible already.
9236 SCIP_Bool* success /**< buffer to store whether we succeeded to construct a solution that satisfies all provided constraints */
9269 sol != NULL ? (SCIPsolGetHeur(sol) != NULL ? SCIPheurGetName(SCIPsolGetHeur(sol)) : "tree") : "LP");
9297 ((viol > 0.0 && consdata->lincoefs[consdata->linvar_mayincrease] > 0.0) || (viol < 0.0 && consdata->lincoefs[consdata->linvar_mayincrease] < 0.0)) )
9317 SCIPdebugMessage("increase <%s> by %g to %g\n", SCIPvarGetName(var), delta, SCIPgetSolVal(scip, newsol, var));
9328 ((viol > 0.0 && consdata->lincoefs[consdata->linvar_maydecrease] < 0.0) || (viol < 0.0 && consdata->lincoefs[consdata->linvar_maydecrease] > 0.0)) )
9347 SCIPdebugMessage("increase <%s> by %g to %g\n", SCIPvarGetName(var), delta, SCIPgetSolVal(scip, newsol, var));
9356 /* still here... so probably we could not make constraint feasible due to variable bounds, thus give up */
9361 /* if we have a solution that should satisfy all quadratic constraints and has a better objective than the current upper bound,
9364 if( c == nconss && (SCIPisInfinity(scip, SCIPgetUpperbound(scip)) || SCIPisSumLT(scip, SCIPgetSolTransObj(scip, newsol), SCIPgetUpperbound(scip))) )
9366 SCIPdebugMessage("pass solution with objective val %g to trysol heuristic\n", SCIPgetSolTransObj(scip, newsol));
9398 /* if a quadratic expression has been simplified, then all children of the node should be variables */
9431 /* these do not look like an quadratic expression (assuming the expression graph simplifier did run) */
9443 SCIPwarningMessage(scip, "unexpected expression operator %d in nonlinear constraint <%s>\n", SCIPexprgraphGetNodeOperator(node), SCIPconsGetName(cons));
9458 SCIPgetNLinearVarsNonlinear(scip, cons), SCIPgetLinearVarsNonlinear(scip, cons), SCIPgetLinearCoefsNonlinear(scip, cons),
9472 SCIP_CALL( SCIPaddQuadVarQuadratic(scip, upgdconss[0], (SCIP_VAR*)SCIPexprgraphGetNodeVar(exprgraph, SCIPexprgraphGetNodeChildren(node)[i]), 0.0, 0.0) );
9530 (SCIP_VAR*)SCIPexprgraphGetNodeVar(exprgraph, SCIPexprgraphGetNodeChildren(node)[quadelems[i].idx1]),
9537 (SCIP_VAR*)SCIPexprgraphGetNodeVar(exprgraph, SCIPexprgraphGetNodeChildren(node)[quadelems[i].idx1]),
9538 (SCIP_VAR*)SCIPexprgraphGetNodeVar(exprgraph, SCIPexprgraphGetNodeChildren(node)[quadelems[i].idx2]),
9574 /** destructor of constraint handler to free constraint handler data (called when SCIP is exiting) */
9625 /** deinitialization method of constraint handler (called before transformed problem is freed) */
9650 /** presolving initialization method of constraint handler (called when presolving is about to begin) */
9670 /** presolving deinitialization method of constraint handler (called after presolving has been finished) */
9721 /** solving process initialization method of constraint handler (called when branch and bound process is about to begin) */
9743 /* check for a linear variable that can be increase or decreased without harming feasibility */
9751 consdata->lincoefsmin = MIN(consdata->lincoefsmin, REALABS(consdata->lincoefs[i])); /*lint !e666 */
9752 consdata->lincoefsmax = MAX(consdata->lincoefsmax, REALABS(consdata->lincoefs[i])); /*lint !e666 */
9772 SCIP_CALL( SCIPallocBlockMemoryArray(scip, &consdata->sepabilinvar2pos, consdata->nbilinterms) );
9782 SCIP_CALL( consdataFindQuadVarTerm(scip, consdata, consdata->bilinterms[i].var2, &consdata->sepabilinvar2pos[i]) );
9788 /* check if constraint function is factorable, i.e., can be written as product of two linear functions */
9802 SCIP_CALL( SCIPcatchEvent(scip, SCIP_EVENTTYPE_SOLFOUND, eventhdlr, (SCIP_EVENTDATA*)conshdlr, &conshdlrdata->newsoleventfilterpos) );
9807 SCIPverbMessage(scip, SCIP_VERBLEVEL_HIGH, NULL, "Quadratic constraint handler does not have LAPACK for eigenvalue computation. Will assume that matrices (with size > 2x2) are indefinite.\n");
9818 /** solving process deinitialization method of constraint handler (called before branch and bound process data is freed) */
9840 SCIP_CALL( SCIPdropEvent(scip, SCIP_EVENTTYPE_SOLFOUND, eventhdlr, (SCIP_EVENTDATA*)conshdlr, conshdlrdata->newsoleventfilterpos) );
9920 SCIP_CALL( SCIPgetTransformedVar(scip, targetdata->quadvarterms[i].var, &targetdata->quadvarterms[i].var) );
9926 SCIP_CALL( SCIPgetTransformedVar(scip, targetdata->bilinterms[i].var1, &targetdata->bilinterms[i].var1) );
9927 SCIP_CALL( SCIPgetTransformedVar(scip, targetdata->bilinterms[i].var2, &targetdata->bilinterms[i].var2) );
9940 SCIPconsIsInitial(sourcecons), SCIPconsIsSeparated(sourcecons), SCIPconsIsEnforced(sourcecons),
9942 SCIPconsIsModifiable(sourcecons), SCIPconsIsDynamic(sourcecons), SCIPconsIsRemovable(sourcecons),
9951 /** LP initialization method of constraint handler (called before the initial LP relaxation at a node is solved) */
9986 SCIP_CALL( SCIPcreateEmptyRowCons(scip, &row, SCIPconsGetHdlr(conss[c]), SCIPconsGetName(conss[c]), consdata->lhs, consdata->rhs,
9988 SCIP_CALL( SCIPaddVarsToRow(scip, row, consdata->nlinvars, consdata->linvars, consdata->lincoefs) );
10033 SCIP_CALL( generateCut(scip, conshdlr, conss[c], x, NULL, consdata->isconvex ? SCIP_SIDETYPE_RIGHT : SCIP_SIDETYPE_LEFT, &row, NULL, FALSE, -SCIPinfinity(scip)) ); /*lint !e613 */
10038 SCIPdebugMessage("initlp adds row <%s> for lambda = %g of conss <%s>\n", SCIProwGetName(row), lambda, SCIPconsGetName(conss[c])); /*lint !e613 */
10067 /* set reference point to 0 projected on bounds for unbounded variables or in between lower and upper bound for bounded variables
10068 * in the first round, we set it closer to the best bound, in the second closer to the worst bound
10069 * the reason is, that for a bilinear term with bounded variables, there are always two linear underestimators
10070 * if the reference point is set to the middle, then rounding and luck decides which underestimator is chosen
10071 * we thus choose the reference point not to be the middle, so both McCormick terms are definitely chosen one time
10072 * of course, the possible number of cuts is something in the order of 2^nquadvars, and we choose two of them here
10096 x[i] = lambda * SCIPvarGetBestBoundLocal(var) + (1.0-lambda) * SCIPvarGetWorstBoundLocal(var);
10105 SCIP_CALL( generateCut(scip, conshdlr, conss[c], x, NULL, SCIP_SIDETYPE_RIGHT, &row, NULL, conshdlrdata->checkcurvature, -SCIPinfinity(scip)) ); /*lint !e613 */
10110 SCIPdebugMessage("initlp adds row <%s> for rhs of conss <%s>, round %d\n", SCIProwGetName(row), SCIPconsGetName(conss[c]), k); /*lint !e613 */
10121 SCIP_CALL( generateCut(scip, conshdlr, conss[c], x, NULL, SCIP_SIDETYPE_LEFT, &row, NULL, conshdlrdata->checkcurvature, -SCIPinfinity(scip)) ); /*lint !e613 */
10126 SCIPdebugMessage("initlp adds row <%s> for lhs of conss <%s>, round %d\n", SCIProwGetName(row), SCIPconsGetName(conss[c]), k); /*lint !e613 */
10136 /* if there are unbounded variables, then there is typically only at most one possible underestimator, so don't try another round
10137 * similar, if there are no bilinear terms and no linearizations of square terms, then the reference point does not matter, so don't do another round */
10175 /* at root, check if we want to solve the NLP relaxation and use its solutions as reference point
10176 * if there is something convex, then linearizing in the solution of the NLP relaxation can be very useful
10179 (SCIPgetNContVars(scip) >= conshdlrdata->sepanlpmincont * SCIPgetNVars(scip) || (SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_UNBOUNDEDRAY && conshdlrdata->sepanlpmincont <= 1.0)) &&
10192 * but first check whether there is a violated constraint side which corresponds to a convex function
10202 if( !SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) && !SCIPisGT(scip, consdata->rhsviol, SCIPfeastol(scip)) )
10270 SCIPdebugMessage("pass solution with obj. value %g to trysol\n", SCIPgetSolOrigObj(scip, nlpsol));
10275 SCIP_CALL( addLinearizationCuts(scip, conshdlr, conss, nconss, nlpsol, &lpsolseparated, conshdlrdata->mincutefficacysepa) );
10279 /* if a cut that separated the LP solution was added, then return, otherwise continue with usual separation in LP solution */
10289 /* if we do not want to try solving the NLP, or have no NLP, or have no NLP solver, or solving the NLP failed,
10290 * or separating with NLP solution as reference point failed, then try (again) with LP solution as reference point
10293 SCIP_CALL( separatePoint(scip, conshdlr, conss, nconss, nusefulconss, NULL, conshdlrdata->mincutefficacysepa, FALSE, result, NULL) );
10320 SCIP_CALL( separatePoint(scip, conshdlr, conss, nconss, nusefulconss, sol, conshdlrdata->mincutefficacysepa, FALSE, result, NULL) );
10362 SCIPdebugMessage("enfolp with max violation %g in cons <%s>\n", maxviol, SCIPconsGetName(maxviolcon));
10365 * (maybe the LP does not think that the cuts we add are violated, or we do ECP on a high-dimensional convex function)
10366 * in this case, check if some limit is hit or SCIP should stop for some other reason and terminate enforcement by creating a dummy node
10367 * (in optimized more, returning SCIP_INFEASIBLE in *result would be sufficient, but in debug mode this would give an assert in scip.c)
10368 * the reason to wait for 100 rounds is to avoid calls to SCIPisStopped in normal runs, which may be expensive
10379 SCIP_CALL( SCIPcreateChild(scip, &child, 1.0, SCIPnodeGetEstimate(SCIPgetCurrentNode(scip))) );
10388 /* cut off the current subtree, if a limit on the enforcement rounds should be applied. At this point, feasible
10389 * solutions might get cut off; the enfolplimit parameter should therefore only be set if SCIP is used as a
10392 if( conshdlrdata->enfolplimit != -1 && conshdlrdata->nenfolprounds > conshdlrdata->enfolplimit )
10395 "cut off subtree because enforcement limit was reached; this might lead to incorrect results\n");
10411 SCIPdebugMessage("propagation succeeded (%s)\n", propresult == SCIP_CUTOFF ? "cutoff" : "reduceddom");
10418 * thus, in the latter case, we are also happy if the efficacy is at least, say, 75% of the maximal violation
10421 minefficacy = MIN(0.75*maxviol, conshdlrdata->mincutefficacyenfofac * SCIPfeastol(scip)); /*lint !e666 */
10423 SCIP_CALL( separatePoint(scip, conshdlr, conss, nconss, nusefulconss, NULL, minefficacy, TRUE, &separateresult, &sepaefficacy) );
10432 SCIPdebugMessage("separation succeeded (bestefficacy = %g, minefficacy = %g)\n", sepaefficacy, minefficacy);
10441 SCIPdebugMessage("separation failed (bestefficacy = %g < %g = minefficacy ); max viol: %g\n", sepaefficacy, minefficacy, maxviol);
10446 /* if sepastore can decrease LP feasibility tolerance, we can add cuts with efficacy in [eps, feastol] */
10447 leastpossibleefficacy = SCIPgetRelaxFeastolFactor(scip) > 0.0 ? SCIPepsilon(scip) : SCIPfeastol(scip);
10451 SCIP_CALL( separatePoint(scip, conshdlr, conss, nconss, nusefulconss, NULL, leastpossibleefficacy, TRUE, &separateresult, &sepaefficacy) );
10468 /* fallback 2: separation probably failed because of numerical difficulties with a convex constraint;
10469 * if noone declared solution infeasible yet and we had not even found a weak cut, try to resolve by branching
10480 SCIP_CALL( replaceByLinearConstraints(scip, conss, nconss, &addedcons, &reduceddom, &infeasible) );
10481 /* if the linear constraints are actually feasible, then adding them and returning SCIP_CONSADDED confuses SCIP
10482 * when it enforces the new constraints again and nothing resolves the infeasiblity that we declare here thus,
10483 * we only add them if considered violated, and otherwise claim the solution is feasible (but print a
10494 SCIPwarningMessage(scip, "could not enforce feasibility by separating or branching; declaring solution with viol %g as feasible\n", maxviol);
10501 SCIPdebugMessage("Could not find any usual branching variable candidate. Proposed variable <%s> with LP value %g for branching.\n",
10558 if( !SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) && !SCIPisGT(scip, consdata->rhsviol, SCIPfeastol(scip)) )
10566 SCIP_CALL( SCIPaddExternBranchCand(scip, var, MAX(consdata->lhsviol, consdata->rhsviol), SCIP_INVALID) );
10576 SCIP_CALL( SCIPaddExternBranchCand(scip, var, MAX(consdata->lhsviol, consdata->rhsviol), SCIP_INVALID) );
10584 SCIPdebugMessage("All variables in violated constraints fixed (up to epsilon). Cannot find branching candidate. Forcing solution of LP.\n");
10630 * then try the reformulations (replacing products with binaries, disaggregation, setting default variable bounds)
10632 * @todo first do all usual presolving steps, then check SCIPisPresolveFinished(scip), and if true then do reformulations (and usual steps again)
10634 doreformulations = (nrounds > 0 || SCIPconshdlrWasPresolvingDelayed(conshdlr)) && SCIPisPresolveFinished(scip);
10635 SCIPdebugMessage("presolving will %swait with reformulation\n", doreformulations ? "not " : "");
10672 /* call upgrade methods if the constraint has not been presolved yet or there has been a bound tightening or possibly be a change in variable type
10673 * we want to do this before (multi)aggregated variables are replaced, since that may change structure, e.g., introduce bilinear terms
10698 SCIPdebugMessage("solving constraint <%s> says problem is infeasible in presolve\n", SCIPconsGetName(conss[c]));
10729 /* user not so empathic about AND, or we don't have products of two binaries, so try this more general reformulation */
10754 SCIP_CALL( dropVarEvents(scip, conshdlrdata->eventhdlr, conss[c]) ); /* well, there shouldn't be any variables left anyway */
10766 SCIPdebugMessage("constraint <%s> is constant and feasible, deleting\n", SCIPconsGetName(conss[c]));
10775 /* try domain propagation if there were bound changes or constraint has changed (in which case, processVarEvents may have set ispropagated to false) */
10784 SCIPdebugMessage("starting domain propagation round %d of %d\n", roundnr, conshdlrdata->maxproproundspresolve);
10790 SCIPdebugMessage("propagation on constraint <%s> says problem is infeasible in presolve\n", SCIPconsGetName(conss[c]));
10856 SCIPdebugMessage("make variable <%s> implicit integer due to constraint <%s>\n", SCIPvarGetName(candidate), SCIPconsGetName(conss[c]));
10861 SCIPdebugMessage("infeasible upgrade of variable <%s> to integral type, domain is empty\n", SCIPvarGetName(candidate));
10889 /* if we did not try reformulations, ensure that presolving is called again even if there were only a few changes (< abortfac) */
10942 /* @todo try to be more clever, but variable locks that depend on the bounds of other variables are not trival to maintain */
10943 SCIP_CALL( SCIPaddVarLocks(scip, consdata->quadvarterms[i].var, nlockspos+nlocksneg, nlockspos+nlocksneg) );
11046 SCIP_CALL( SCIPwriteVarsPolynomial(scip, file, monomialvars, monomialexps, monomialcoefs, monomialnvars, nmonomials, TRUE) );
11112 if( SCIPisGT(scip, consdata->lhsviol, SCIPfeastol(scip)) || SCIPisGT(scip, consdata->rhsviol, SCIPfeastol(scip)) )
11121 SCIPinfoMessage(scip, NULL, "violation: left hand side is violated by %.15g (scaled: %.15g)\n", consdata->lhs - consdata->activity, consdata->lhsviol);
11125 SCIPinfoMessage(scip, NULL, "violation: right hand side is violated by %.15g (scaled: %.15g)\n", consdata->activity - consdata->rhs, consdata->rhsviol);
11133 /* do not try to shift linear variables if activity is at infinity (leads to setting variable to infinity in solution, which is not allowed) */
11147 if( !(consdata->linvar_mayincrease >= 0 && consdata->lincoefs[consdata->linvar_mayincrease] > 0.0) &&
11148 ! (consdata->linvar_maydecrease >= 0 && consdata->lincoefs[consdata->linvar_maydecrease] < 0.0) )
11156 if( !(consdata->linvar_mayincrease >= 0 && consdata->lincoefs[consdata->linvar_mayincrease] < 0.0) &&
11157 ! (consdata->linvar_maydecrease >= 0 && consdata->lincoefs[consdata->linvar_maydecrease] > 0.0) )
11218 SCIP_CALL( SCIPgetVarCopy(sourcescip, scip, consdata->linvars[i], &linvars[i], varmap, consmap, global, valid) );
11237 SCIP_CALL( SCIPgetVarCopy(sourcescip, scip, consdata->quadvarterms[i].var, &quadvarterms[i].var, varmap, consmap, global, valid) );
11331 if( isdigit((unsigned char)str[0]) || ((str[0] == '-' || str[0] == '+') && isdigit((unsigned char)str[1])) )
11360 SCIP_CALL( SCIPparseVarsPolynomial(scip, str, &monomialvars, &monomialexps, &monomialcoefs, &monomialnvars, &nmonomials, &endptr, success) );
11484 SCIPfreeParseVarsPolynomialData(scip, &monomialvars, &monomialexps, &monomialcoefs, &monomialnvars, nmonomials);
11520 /** constraint method of constraint handler which returns the number of variables (if possible) */
11576 SCIP_CALL( SCIPsetConshdlrPresol(scip, conshdlr, consPresolQuadratic, CONSHDLR_MAXPREROUNDS, CONSHDLR_DELAYPRESOL) );
11578 SCIP_CALL( SCIPsetConshdlrProp(scip, conshdlr, consPropQuadratic, CONSHDLR_PROPFREQ, CONSHDLR_DELAYPROP,
11580 SCIP_CALL( SCIPsetConshdlrSepa(scip, conshdlr, consSepalpQuadratic, consSepasolQuadratic, CONSHDLR_SEPAFREQ,
11586 "max. length of linear term which when multiplied with a binary variables is replaced by an auxiliary variable and a linear reformulation (0 to turn off)",
11590 "empathy level for using the AND constraint handler: 0 always avoid using AND; 1 use AND sometimes; 2 use AND as often as possible",
11594 "whether to make non-varbound linear constraints added due to replacing products with binary variables initial",
11598 "limit (as factor on 1/feastol) on coefficients and coef. range in linear constraints created when replacing products with binary variables",
11602 "minimal efficacy for a cut to be added to the LP during separation; overwrites separating/efficacy",
11606 "minimal target efficacy of a cut in order to add it to relaxation during enforcement as a factor of the feasibility tolerance (may be ignored)",
11610 "whether scaling of infeasibility is 'o'ff, by sup-norm of function 'g'radient, or by left/right hand 's'ide",
11614 "maximal coef range of a cut (maximal coefficient divided by minimal coefficient) in order to be added to LP relaxation",
11618 "whether linearizations of convex quadratic constraints should be added to cutpool in a solution found by some heuristic",
11630 "whether to try to make solutions in check function feasible by shifting a linear variable (esp. useful if constraint was actually objective function)",
11634 "whether to disaggregate quadratic parts that decompose into a sum of non-overlapping quadratic terms",
11638 "limit on number of propagation rounds for a single constraint within one round of SCIP propagation during solve",
11642 "limit on number of propagation rounds for a single constraint within one round of SCIP presolve",
11646 "maximum number of enforcement rounds before declaring the LP relaxation infeasible (-1: no limit); WARNING: changing this parameter might lead to incorrect results!",
11650 "minimal required fraction of continuous variables in problem to use solution of NLP relaxation in root for separation",
11658 SCIP_CALL( SCIPincludeEventhdlrBasic(scip, &(conshdlrdata->eventhdlr),CONSHDLR_NAME"_boundchange", "signals a bound change to a quadratic constraint",
11662 SCIP_CALL( SCIPincludeEventhdlrBasic(scip, NULL, CONSHDLR_NAME"_newsolution", "handles the event that a new primal solution has been found",
11666 SCIP_CALL( SCIPincludeNonlinconsUpgrade(scip, nonlinconsUpgdQuadratic, NULL, NONLINCONSUPGD_PRIORITY, TRUE, CONSHDLR_NAME) );
11674 SCIP_DECL_QUADCONSUPGD((*quadconsupgd)), /**< method to call for upgrading quadratic constraint */
11721 for( i = conshdlrdata->nquadconsupgrades; i > 0 && conshdlrdata->quadconsupgrades[i-1]->priority < quadconsupgrade->priority; --i )
11728 (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "constraints/"CONSHDLR_NAME"/upgrade/%s", conshdlrname);
11729 (void) SCIPsnprintf(paramdesc, SCIP_MAXSTRLEN, "enable quadratic upgrading for constraint handler <%s>", conshdlrname);
11746 * @note the constraint gets captured, hence at one point you have to release it using the method SCIPreleaseCons()
11779 SCIP_Bool removable /**< should the relaxation be removed from the LP due to aging or cleanup?
11809 SCIP_CALL( SCIPcreateCons(scip, cons, name, conshdlr, consdata, initial, separate, enforce, check, propagate,
11813 SCIP_CALL( SCIPhashmapCreate(&quadvaridxs, SCIPblkmem(scip), SCIPcalcHashtableSize(5 * nquadterms)) );
11833 SCIP_CALL( SCIPhashmapInsert(quadvaridxs, quadvars1[i], (void*)(size_t)(consdata->nquadvars-1)) );
11854 SCIP_CALL( SCIPhashmapInsert(quadvaridxs, quadvars2[i], (void*)(size_t)(consdata->nquadvars-1)) );
11890 /* if it's a linear coefficient for a quadratic variable, add it there, otherwise add as linear variable */
11934 * @note the constraint gets captured, hence at one point you have to release it using the method SCIPreleaseCons()
11962 * \ell \leq \sum_{i=1}^n b_i x_i + \sum_{j=1}^m (a_j y_j^2 + b_j y_j) + \sum_{k=1}^p c_kv_kw_k \leq u.
11965 * @note the constraint gets captured, hence at one point you have to release it using the method SCIPreleaseCons()
11988 SCIP_Bool removable /**< should the constraint be removed from the LP due to aging or cleanup? */
12013 SCIP_CALL( SCIPcreateCons(scip, cons, name, conshdlr, consdata, initial, separate, enforce, check, propagate,
12033 * \ell \leq \sum_{i=1}^n b_i x_i + \sum_{j=1}^m (a_j y_j^2 + b_j y_j) + \sum_{k=1}^p c_kv_kw_k \leq u.
12036 * @note the constraint gets captured, hence at one point you have to release it using the method SCIPreleaseCons()
12061 /** Adds a constant to the constraint function, that is, subtracts a constant from both sides */
12168 consdata->ispresolved = consdata->ispresolved && !SCIPisZero(scip, consdata->quadvarterms[pos].lincoef);
12217 consdata->ispresolved = consdata->ispresolved && !SCIPisZero(scip, consdata->quadvarterms[pos].sqrcoef);
12398 * @note If the quadratic variable terms have not been sorted before, then a search may reorder the current order of the terms.
12470 /** Check the quadratic function of a quadratic constraint for its semi-definiteness, if not done yet.
12484 /** Indicates whether the quadratic function of a quadratic constraint is (known to be) convex.
12502 /** Indicates whether the quadratic function of a quadratic constraint is (known to be) concave.
12550 * That is, checks whether each variable with a square term is fixed and for each bilinear term at least one variable is fixed.
12603 SCIP_HASHMAP* scipvar2nlpivar, /**< mapping from SCIP variables to variable indices in NLPI */
|