All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
sepa_cgmip.c
Go to the documentation of this file.
21 * Separate Chvátal-Gomory cuts using a sub-MIP. The approach is based on the following papers.
25 * in: M. Jünger and V. Kaibel (eds.) Integer Programming and Combinatorial Optimization IPCO 2005,@n
39 * - The CMIR-routines of SCIP can be used (if @p usecmir is true). One can determine which bound is
44 * false). The cut is not take from the solution of the MIP, but is recomputed, and some care (but
50 * - If paramter @a earlyterm is true, the separation is run until the first cut that is violated is
51 * found. (Note that these cuts are not necessarily added to the LP, because here also the norm of
52 * the cuts are taken into account - which cannot easily be included into the separation subscip.)
61 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
79 #define SEPA_DELAY FALSE /**< should separation method be delayed, if other separators found cuts? */
82 #define DEFAULT_MAXROUNDSROOT 50 /**< maximal number of separation rounds in the root node (-1: unlimited) */
85 #define DEFAULT_TIMELIMIT 1e20 /**< time limit for sub-MIP (set to infinity in order to be deterministic) */
88 #define DEFAULT_MINNODELIMIT 500LL /**< minimum number of nodes considered for sub-MIP (-1: unlimited) */
89 #define DEFAULT_MAXNODELIMIT 5000LL /**< maximum number of nodes considered for sub-MIP (-1: unlimited) */
93 #define DEFAULT_ONLYINTVARS FALSE /**< Generate cuts for problems with only integer variables? */
95 #define DEFAULT_CONTCONVERT FALSE /**< Convert some integral variables to be continuous to reduce the size of the sub-MIP? */
96 #define DEFAULT_CONTCONVFRAC 0.1 /**< fraction of integral variables converted to be continuous (if contconvert) */
97 #define DEFAULT_CONTCONVMIN 100 /**< minimum number of integral variables before some are converted to be continuous */
98 #define DEFAULT_INTCONVERT FALSE /**< Convert some integral variables attaining fractional values to have integral value? */
99 #define DEFAULT_INTCONVFRAC 0.1 /**< fraction of fractional integral variables converted to have integral value (if intconvert) */
100 #define DEFAULT_INTCONVMIN 100 /**< minimum number of integral variables before some are converted to have integral value */
101 #define DEFAULT_SKIPMULTBOUNDS TRUE /**< Skip the upper bounds on the multipliers in the sub-MIP? */
102 #define DEFAULT_OBJLONE FALSE /**< Should the objective of the sub-MIP only minimize the l1-norm of the multipliers? */
105 #define DEFAULT_DYNAMICCUTS TRUE /**< Should generated cuts be removed from the LP if they are no longer tight? */
108 #define DEFAULT_CMIROWNBOUNDS FALSE /**< Tell CMIR-generator which bounds to used in rounding? */
109 #define DEFAULT_USECUTPOOL TRUE /**< Use cutpool to store CG-cuts even if the are not efficient? */
110 #define DEFAULT_PRIMALSEPARATION TRUE /**< Only separate cuts that are tight for the best feasible solution? */
111 #define DEFAULT_EARLYTERM TRUE /**< Terminate separation if a violated (but possibly sub-optimal) cut has been found? */
112 #define DEFAULT_ADDVIOLATIONCONS FALSE /**< Add constraint to subscip that only allows violated cuts (otherwise add obj. limit)?*/
113 #define DEFAULT_ADDVIOLCONSHDLR FALSE /**< Add constraint handler to filter out violated cuts? */
114 #define DEFAULT_CONSHDLRUSENORM TRUE /**< Should the violation constraint handler use the norm of a cut to check for feasibility? */
115 #define DEFAULT_USEOBJUB FALSE /**< Use upper bound on objective function (via primal solution)? */
117 #define DEFAULT_SUBSCIPFAST TRUE /**< Should the settings for the sub-MIP be optimized for speed? */
118 #define DEFAULT_OUTPUT FALSE /**< Should information about the sub-MIP and cuts be displayed? */
121 #define NCOLSTOOSMALL 5 /**< only separate if the number of columns is larger than this number */
124 #define BETAEPSILONVALUE 1e-02 /**< epsilon value for fracbeta - is larger than EPSILONVALUE for numerical stability */
126 #define CONSHDLRFULLNORM FALSE /**< compute real cut and compute norm for this (if addviolconshdlr and conshdlrusenorm are true) */
138 #define MAXWEIGHTRANGE 1e+05 /**< maximal valid range max(|weights|)/min(|weights|) of row weights */
155 SCIP_Longint minnodelimit; /**< minimum number of nodes considered for sub-MIP (-1: unlimited) */
156 SCIP_Longint maxnodelimit; /**< maximum number of nodes considered for sub-MIP (-1: unlimited) */
163 SCIP_Bool contconvert; /**< Convert some integral variables to be continuous to reduce the size of the sub-MIP? */
164 SCIP_Real contconvfrac; /**< fraction of integral variables converted to be continuous (if contconvert) */
165 int contconvmin; /**< minimum number of integral variables before some are converted to be continuous */
166 SCIP_Bool intconvert; /**< Convert some integral variables attaining fractional values to have integral value? */
167 SCIP_Real intconvfrac; /**< fraction of frac. integral variables converted to have integral value (if intconvert) */
168 int intconvmin; /**< minimum number of integral variables before some are converted to have integral value */
170 SCIP_Bool objlone; /**< Should the objective of the sub-MIP only minimize the l1-norm of the multipliers? */
173 SCIP_Bool dynamiccuts; /**< Should generated cuts be removed from the LP if they are no longer tight? */
178 SCIP_Bool primalseparation; /**< Only separate cuts that are tight for the best feasible solution? */
179 SCIP_Bool earlyterm; /**< Terminate separation if a violated (but possibly sub-optimal) cut has been found? */
182 SCIP_Bool conshdlrusenorm; /**< Should the violation constraint handler use the cut-norm to check for feasibility? */
196 colAtUb = 3, /**< variable corresponding to column was at it's upper bound and was complemented */
197 colAtLb = 4 /**< variable corresponding to column was at it's lower bound (possibly complemented) */
229 SCIP_Bool conshdlrfullnorm; /**< compute real cut and compute norm for this (if addviolconshdlr and conshdlrusenorm are true) */
262 SCIP_Bool* localboundsused, /**< pointer to store whether local bounds were used in summation */
317 SCIP_CALL( computeCut(mipdata->scip, mipdata->sepa, mipdata, mipdata->sepadata, sol, cutcoefs, &rhs, &localrowsused, &localboundsused, &cutrank, &success) );
447 SCIPdebugMessage("Violated cut from solution - act: %f, rhs: %f, norm: %f, eff.: %f\n", act, rhs, norm, (act-rhs)/norm);
451 SCIPdebugMessage("Rejected cut from solution - act: %f, rhs: %f, norm: %f, eff.: %f\n", act, rhs, norm, (act-rhs)/norm);
461 /** destructor of constraint handler to free constraint handler data (called when SCIP is exiting) */
499 *result = SCIP_CUTOFF; /* cutoff, since all integer variables are integer, but the solution is not feasible */
587 /** stores nonzero elements of dense coefficient vector as sparse vector and calculates activity and norm
712 * \gamma + \sigma\, a_j\,s \leq \sum_{i \neq j} a_i\, x_i' + \sigma a_j\, x_j' \leq \delta + \sigma\, a_j\, s
716 * \frac{1}{\sigma} \ell_j + s \leq x_j' \leq \frac{1}{\sigma} u_j + s, \qquad \mbox{ if }\sigma > 0
720 * \frac{1}{\sigma} u_j + s \leq x_j' \leq \frac{1}{\sigma} \ell_j + s, \qquad \mbox{ if }\sigma < 0.
725 * - If \f$x_j \geq \ell_j\f$ has a (nonzero) lower bound, one can use \f$s = -\ell_j\f$, \f$\sigma = 1\f$,
726 * and obtain \f$\gamma - a_j\,\ell_j \leq a^T x' \leq \delta - a_j\,\ell_j\f$, \f$0 \leq x_j' \leq u_j - \ell_j\f$.
728 * - If \f$x_j \leq u_j\f$ has a (nonzero) upper bound, one can use \f$s = u_j\f$, \f$\sigma = -1\f$,
729 * and obtain \f$\gamma - a_j\,u_j \leq \sum_{i \neq j} a_i\, x_i' - a_j\, x_j' \leq \delta - a_j\, u_j\f$,
838 * where \f$r\f$ is the size of the current row, \f$a \in [0,1]\f$ is a parameter, and \f$r_{\mbox{max}}\f$ and
922 * Note that we use zero multipliers for the bounds on the continuous variables \f$r\f$. Moreover,
942 * Following Fischetti and Lodi [2005], let \f$(x^*,r^*)\f$ be a fractional solution of the above
956 * Here, \f$w\f$ is a weight vector; it's idea is to make the sum over all components of \f$y\f$ as
961 * - If the lower bounds on \f$x_i\f$ or \f$r_j\f$ are finite, we shift the variable to have a zero
962 * lower bound, i.e., we replace it by \f$x_i - \ell_i\f$ (or \f$r_j - u_j\f$). This is helpful in
964 * allows to drop a variable if \f$x^*_i = 0\f$, see the next comment. If the lower bounds are not
968 * practice), we require \f$f_j = 0\f$, i.e., we force that \f$(\tilde{A}^T y + \overline{z})_j =
971 * - If \f$x^*_j = 0 = \ell_j\f$ (after the above preprocessing), we drop variable \f$\alpha_j\f$
976 * - If \f$x^*_i = u_i\f$, we complement the variable and drop it from the formulation, since the
982 * - If a continuous variable \f$r_j\f$ is free, we have to force equality for the corresponding components in
994 * - If required, i.e., parameter @p primalseparation is true, we force a primal separation step. For
995 * this we require that the cut is tight at the currently best solution. To get reliable solutions
998 * - If required (via parameters @p useobjub or @p useobjlb), we add a row corresponding to the objective function with
1076 (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "sepa_cgmip separating MIP (%s)", SCIPgetProbName(scip));
1097 /* store lhs/rhs for complementing (see below) and compute maximal nonzeros of candidate rows */
1121 if ( ! sepadata->onlyactiverows && sepadata->maxrowage > 0 && SCIProwGetAge(row) > sepadata->maxrowage )
1228 /* integral variables taking integral values are not interesting - will be substituted out below */
1273 /* integral variables taking integral values are not interesting - will be substituted out below */
1295 /* if integer variable is at its upper bound -> complementing (this also generates a 0 lower bound) */
1299 SCIP_CALL( transformColumn(scip, sepadata, mipdata, col, ub[j], -1.0, lhs, rhs, &(lb[j]), &(ub[j]), &(primsol[j])) );
1311 SCIP_CALL( transformColumn(scip, sepadata, mipdata, col, -lb[j], 1.0, lhs, rhs, &(lb[j]), &(ub[j]), &(primsol[j])) );
1330 SCIP_CALL( transformColumn(scip, sepadata, mipdata, col, ub[j], -1.0, lhs, rhs, &(lb[j]), &(ub[j]), &(primsol[j])) );
1352 SCIPdebugMessage("Converted %u fractional integral variables to have integral value.\n", nintconverted);
1359 SCIPdebugMessage("original variables: %d integral, %d continuous, %u shifted, %u complemented, %u at lb, %u at ub\n",
1360 SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip) + SCIPgetNImplVars(scip), SCIPgetNContVars(scip),
1386 if ( ! sepadata->onlyactiverows && sepadata->maxrowage > 0 && SCIProwGetAge(row) > sepadata->maxrowage )
1402 assert( SCIPisFeasEQ(scip, SCIPgetRowLPActivity(scip, row), SCIProwGetLhs(row)) ); /* equations should always be active */
1406 weight = - sepadata->objweight * computeObjWeightSize(SCIProwGetNLPNonz(row), minrowsize, maxrowsize);
1442 weight = -sepadata->objweight * computeObjWeightSize(SCIProwGetNLPNonz(row), minrowsize, maxrowsize);
1457 SCIPdebugMessage("Created variable <%s> for >= inequality <%s> (weight: %f).\n", name, SCIProwGetName(row), weight);
1473 weight = -sepadata->objweight * computeObjWeightSize(SCIProwGetNLPNonz(row), minrowsize, maxrowsize);
1487 SCIPdebugMessage("Created variable <%s> for <= inequality <%s> (weight: %f).\n", name, SCIProwGetName(row), weight);
1507 weight = -sepadata->objweight * computeObjWeightSize(SCIPgetNObjVars(scip), minrowsize, maxrowsize);
1521 SCIPdebugMessage("Created variable <%s> for upper bound on objective (weight: %f).\n", name, weight);
1535 SCIPdebugMessage("Created variable <%s> for lower bound on objective (weight: %f).\n", name, weight);
1563 SCIP_CALL( SCIPcreateVar(subscip, &(mipdata->alpha[j]), name, -sepadata->cutcoefbnd, sepadata->cutcoefbnd, obj,
1602 SCIP_CALL( SCIPcreateVar(subscip, &(mipdata->beta), "beta", -sepadata->cutcoefbnd, sepadata->cutcoefbnd, 0.0,
1607 SCIP_CALL( SCIPcreateVar(subscip, &(mipdata->beta), "beta", -sepadata->cutcoefbnd, sepadata->cutcoefbnd, -1.0,
1613 SCIP_CALL( SCIPcreateVar(subscip, &(mipdata->fracbeta), "fracbeta", 0.0, 1.0-BETAEPSILONVALUE, 0.0,
1697 if ( (sepadata->useobjub || sepadata->useobjlb) && ! SCIPisZero(scip, SCIPcolGetObj(cols[j])) )
1774 if ( (sepadata->useobjub || sepadata->useobjlb) && ! SCIPisZero(scip, SCIPcolGetObj(cols[j])) )
1896 SCIP_CALL( SCIPcreateConsLinear(subscip, &cons, "beta", nconsvars, consvars, consvals, 0.0, 0.0,
1929 SCIP_CALL( SCIPcreateConsLinear(subscip, &cons, "primalseparation", nconsvars, consvars, consvals, -EPSILONVALUE, EPSILONVALUE,
1956 SCIP_CALL( SCIPcreateConsLinear(subscip, &cons, "violationConstraint", nconsvars, consvars, consvals, MINEFFICACY, SCIPinfinity(subscip),
2023 /* disable memory saving mode: this is likely to result in the maximal depth being reached. This is because DFS
2024 * results in a repeated branching on the alpha-variables, which often have large bounds resulting in deep levels of
2173 /* substract the memory already used by the main SCIP and the estimated memory usage of external software */
2204 SCIPdebugMessage("Solving sub-SCIP (time limit: %f mem limit: %f node limit: %"SCIP_LONGINT_FORMAT") ...\n", timelimit, memorylimit, nodelimit);
2212 * hence, the return code is caught and a warning is printed, only in debug mode, SCIP will stop. */
2218 SCIPwarningMessage(scip, "Error while solving subproblem in CGMIP separator; sub-SCIP terminated with code <%d>\n", retcode);
2234 /* if the solution process was terminated or the problem is infeasible (can happen because of violation constraint) */
2235 if ( status == SCIP_STATUS_TIMELIMIT || status == SCIP_STATUS_USERINTERRUPT || status == SCIP_STATUS_INFEASIBLE || status == SCIP_STATUS_INFORUNBD )
2244 SCIPerrorMessage("Solution of subscip for CG-separation returned with invalid status %d.\n", status);
2258 * hence, the return code is caught and a warning is printed, only in debug mode, SCIP will stop. */
2264 SCIPwarningMessage(scip, "Error while solving subproblem in CGMIP separator; sub-SCIP terminated with code <%d>\n", retcode);
2271 /* if the solution process was terminated or the problem is infeasible (can happen because of violation constraint) */
2272 if ( status == SCIP_STATUS_TIMELIMIT || status == SCIP_STATUS_USERINTERRUPT || status == SCIP_STATUS_NODELIMIT ||
2273 status == SCIP_STATUS_INFEASIBLE || status == SCIP_STATUS_INFORUNBD || status == SCIP_STATUS_MEMLIMIT )
2288 /* all other statuses except optimal or bestsollimit are invalid - (problem cannot be infeasible) */
2291 SCIPerrorMessage("Solution of subscip for CG-separation returned with invalid status %d.\n", status);
2305 * hence, the return code is caught and a warning is printed, only in debug mode, SCIP will stop. */
2311 SCIPwarningMessage(scip, "Error while solving subproblem in CGMIP separator; sub-SCIP terminated with code <%d>\n", retcode);
2329 if ( status == SCIP_STATUS_TIMELIMIT || status == SCIP_STATUS_USERINTERRUPT || status == SCIP_STATUS_MEMLIMIT )
2336 if ( status != SCIP_STATUS_OPTIMAL && status != SCIP_STATUS_STALLNODELIMIT && status != SCIP_STATUS_NODELIMIT )
2338 SCIPerrorMessage("Solution of subscip for CG-separation returned with invalid status %d.\n", status);
2352 * multiplier is positive w.r.t. the feasibility tolerance. In the sub-MIP, however, the rows are
2353 * combined in any case. This makes a difference, if the coefficients in the matrix are large and
2358 * If variable \f$x_j\f$ was complemented, we have \f$x'_j = u_j - x_j\f$. If in the transformed
2359 * system the lower bound is used, its corresponding multiplier is \f$y^T A'_j - \lfloor y^T A'_j
2362 * y^T A'_j - \lfloor y^T A'_j \rfloor = - y^T A_j - \lfloor - y^T A_j \rfloor = - y^T A_j + \lceil y^T A_j \rceil
2366 * If such a variable was at its upper bound before the transformation, it is at its lower bound
2372 * Furthermore, note that if continuous variables have been shifted, the computed violated may be
2387 SCIP_Bool* localboundsused, /**< pointer to store whether local bounds were used in summation */
2468 /* in a suboptimal solution both values may be positive - take the one with larger absolute value */
2500 /* in a suboptimal solution both values may be positive - take the one with larger absolute value */
2552 /* in a suboptimal solution both values may be positive - take the one with larger absolute value */
2641 /* in a suboptimal solution both values may be positive - take the one with larger absolute value */
2765 /* a variable may have status COLUMN, but the corresponding column may not (yet) be in the LP */
2766 if ( pos >= 0 && mipdata->coltype[pos] != colContinuous && mipdata->coltype[pos] != colConverted )
2846 assert( pos == -1 || mipdata->coltype[pos] == colContinuous || mipdata->coltype[pos] == colConverted );
2924 SCIP_CALL( computeCut(scip, sepa, mipdata, sepadata, sol, cutcoefs, &cutrhs, &localrowsused, &localboundsused, &cutrank, &success) );
2939 /* the following test should be treated with care because of numerical differences - see computeCut() */
2971 if ( (mipdata->coltype[j] == colContinuous || mipdata->coltype[j] == colConverted) && mipdata->isshifted[j] )
2993 SCIPdebugMessage("act=%f, rhs=%f, norm=%f, eff=%f\n", cutact, cutrhs, cutnorm, (cutact - cutrhs)/cutnorm);
3006 SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &cut, sepa, name, -SCIPinfinity(scip), cutrhs, cutislocal, FALSE, sepadata->dynamiccuts) );
3042 SCIPdebugMessage(" -> CG-cut <%s>: act=%f, rhs=%f, norm=%f, eff=%f, min=%f, max=%f (range=%f)\n",
3172 /* in a suboptimal solution both values may be positive - take the one with larger absolute value */
3217 if ( mipdata->z[pos] != NULL && SCIPisSumPositive(subscip, SCIPgetSolVal(subscip, sol, mipdata->z[pos])) )
3237 SCIP_CALL( SCIPcalcMIR(scip, NULL, BOUNDSWITCH, USEVBDS, sepadata->allowlocal, FIXINTEGRALRHS, boundsfortrans, boundtypesfortrans,
3239 weights, NULL, 1.0, NULL, NULL, cutcoefs, &cutrhs, &cutact, &success, &cutislocal, &cutrank) );
3261 /* only if the cut if violated - if it is not violated we might store non-local cuts in the pool */
3268 SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &cut, sepa, name, -SCIPinfinity(scip), cutrhs, cutislocal, FALSE, sepadata->dynamiccuts) );
3327 SCIPdebugMessage(" -> CG-cut <%s>: act=%f, rhs=%f, norm=%f, eff=%f, rank=%d, min=%f, max=%f (range=%f)\n",
3353 SCIPdebugMessage(" -> CG-cut <%s> could not be scaled to integral coefficients: act=%f, rhs=%f, norm=%f, eff=%f\n",
3459 /* in a suboptimal solution both values may be positive - take the one with larger absolute value */
3468 SCIP_CALL( SCIPcalcStrongCG(scip, BOUNDSWITCH, USEVBDS, sepadata->allowlocal, (int) MAXAGGRLEN(nvars), MAXWEIGHTRANGE, MINFRAC, MAXFRAC,
3471 SCIPdebugMessage("Strong-CG: success = %u, cut is%sviolated (cutact: %g, cutrhs: %g)\n", success,
3491 /* only if the cut if violated - if it is not violated we might store non-local cuts in the pool */
3498 SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &cut, sepa, name, -SCIPinfinity(scip), cutrhs, cutislocal, FALSE, sepadata->dynamiccuts) );
3529 SCIPdebugMessage(" -> CG-cut <%s> no longer efficacious: act=%f, rhs=%f, norm=%f, eff=%f rank=%d\n",
3556 SCIPdebugMessage(" -> CG-cut <%s>: act=%f, rhs=%f, norm=%f, eff=%f, min=%f, max=%f (range=%f)\n",
3584 SCIPdebugMessage(" -> CG-cut <%s> could not be scaled to integral coefficients: act=%f, rhs=%f, norm=%f, eff=%f\n",
3651 SCIPdebugMessage("Generating CG-cuts from %d sols (cmir: %u, strong-cg: %u) ...\n", nsols, sepadata->usecmir, sepadata->usestrongcg);
3700 SCIP_CALL( createCGCutCMIR(scip, sepa, sepadata, mipdata, sol, cutcoefs, cutvars, cutvals, varsolvals, weights,
3706 SCIP_CALL( createCGCutStrongCG(scip, sepa, sepadata, mipdata, sol, cutcoefs, cutvars, cutvals, varsolvals, weights,
3712 SCIP_CALL( createCGCutDirect(scip, sepa, sepadata, mipdata, sol, cutcoefs, cutvars, cutvals, varsolvals, weights,
3924 if ( ( sepadata->useobjub || sepadata->useobjlb ) && ( sepadata->usecmir || sepadata->usestrongcg ) )
3927 "WARNING - sepa_cgmip: Using objective function bounds and CMIR or Strong-CG functions is useless. Turning off usage of objective function bounds.\n");
3955 else if ( 333 < nrows && nrows <= 874 && 0.15 < firstlptime && firstlptime <= 0.25 && 2614 < ncols && ncols <= 5141 )
3963 else if ( 875 < nrows && nrows <= 1676 && 0.15 < firstlptime && firstlptime <= 0.25 && 1291 < ncols && ncols <= 2613 )
3965 else if ( nrows > 8146 && 0.75 < firstlptime && firstlptime <= 6.25 && 655 < ncols && ncols <= 1290 )
3967 else if ( nrows > 8146 && 0.75 < firstlptime && firstlptime <= 6.25 && 1291 < ncols && ncols <= 2613 )
4061 SCIP_CALL( SCIPincludeSepaBasic(scip, &sepa, SEPA_NAME, SEPA_DESC, SEPA_PRIORITY, SEPA_FREQ, SEPA_MAXBOUNDDIST,
|