Detailed Description
minor separator with intersection cuts
Let X be the matrix of auxiliary variables added for bilinear terms, X_{ij} = x_ix_j. The separator enforces quadratic constraints det(2x2 minor of X) = 0 via intersection cuts.
Definition in file sepa_interminor.c.
#include <assert.h>
#include <string.h>
#include "scip/sepa_interminor.h"
#include "scip/expr.h"
#include "scip/expr_var.h"
#include "scip/expr_pow.h"
#include "scip/expr_product.h"
#include "scip/nlpi_ipopt.h"
#include "scip/cons_nonlinear.h"
Go to the source code of this file.
Data Structures | |
struct | rowdata |
Macros | |
#define | SEPA_NAME "interminor" |
#define | SEPA_DESC "intersection cuts separator to ensure that 2x2 minors of X (= xx') have determinant 0" |
#define | SEPA_PRIORITY 0 |
#define | SEPA_FREQ -1 |
#define | SEPA_MAXBOUNDDIST 1.0 |
#define | SEPA_USESSUBSCIP FALSE |
#define | SEPA_DELAY FALSE |
#define | DEFAULT_MINCUTVIOL 1e-4 |
#define | DEFAULT_RANDSEED 157 |
#define | DEFAULT_MAXROUNDS 10 |
#define | DEFAULT_MAXROUNDSROOT -1 |
#define | BINSEARCH_MAXITERS 120 |
#define | DEFAULT_USESTRENGTHENING FALSE |
#define | DEFAULT_USEBOUNDS FALSE |
Functions | |
static SCIP_RETCODE | sepadataAddMinor (SCIP *scip, SCIP_SEPADATA *sepadata, SCIP_VAR *auxvarxik, SCIP_VAR *auxvarxil, SCIP_VAR *auxvarxjk, SCIP_VAR *auxvarxjl, SCIP_Bool isauxvarxikdiag, SCIP_Bool isauxvarxildiag, SCIP_Bool isauxvarxjkdiag, SCIP_Bool isauxvarxjldiag) |
static SCIP_RETCODE | sepadataClear (SCIP *scip, SCIP_SEPADATA *sepadata) |
static SCIP_RETCODE | getMinorVars (SCIP_SEPADATA *sepadata, int idx, SCIP_VAR **auxvarxik, SCIP_VAR **auxvarxil, SCIP_VAR **auxvarxjk, SCIP_VAR **auxvarxjl, SCIP_Bool *isauxvarxikdiag, SCIP_Bool *isauxvarxildiag, SCIP_Bool *isauxvarxjkdiag, SCIP_Bool *isauxvarxjldiag) |
static SCIP_RETCODE | insertIndex (SCIP *scip, SCIP_HASHMAP *rowmap, SCIP_VAR *row, SCIP_VAR *col, SCIP_VAR *auxvar, int *rowindices, int *nrows) |
static SCIP_RETCODE | detectMinors (SCIP *scip, SCIP_SEPADATA *sepadata) |
static SCIP_RETCODE | constructBasicVars2TableauRowMap (SCIP *scip, int *map) |
static SCIP_RETCODE | computeRestrictionToRay (SCIP *scip, SCIP_Real *ray, SCIP_VAR **vars, SCIP_Real *coefs, SCIP_Real *coefs4b, SCIP_Real *coefscondition, SCIP_Bool usebounds, SCIP_Real *ad, SCIP_Bool *success) |
static SCIP_Real | evalPhiAtRay (SCIP *scip, SCIP_Real t, SCIP_Real a, SCIP_Real b, SCIP_Real c, SCIP_Real d, SCIP_Real e) |
static void | doBinarySearch (SCIP *scip, SCIP_Real a, SCIP_Real b, SCIP_Real c, SCIP_Real d, SCIP_Real e, SCIP_Real *sol) |
static SCIP_Real | isCase4a (SCIP_Real tsol, SCIP_Real *coefs, SCIP_Real *coefscondition) |
static SCIP_Real | computeRoot (SCIP *scip, SCIP_Real *coefs) |
static SCIP_Real | computeIntersectionPoint (SCIP *scip, SCIP_Bool usebounds, SCIP_Real *coefs, SCIP_Real *coefs4b, SCIP_Real *coefscondition) |
static SCIP_RETCODE | addColToCut (SCIP *scip, SCIP_ROWPREP *rowprep, SCIP_Real cutcoef, SCIP_COL *col) |
static SCIP_RETCODE | addRowToCut (SCIP *scip, SCIP_ROWPREP *rowprep, SCIP_Real cutcoef, SCIP_ROW *row, SCIP_Bool *success) |
static SCIP_RETCODE | getTableauRows (SCIP *scip, SCIP_VAR **vars, int *basicvarpos2tableaurow, SCIP_HASHMAP *tableau, SCIP_Real **tableaurows, SCIP_Bool *success) |
static SCIP_RETCODE | addCols (SCIP *scip, SCIP_VAR **vars, SCIP_Real **tableaurows, SCIP_ROWPREP *rowprep, SCIP_Real *rays, int *nrays, int *rayslppos, SCIP_Real *interpoints, SCIP_Bool usebounds, SCIP_Real *ad, SCIP_Bool *success) |
static SCIP_RETCODE | addRows (SCIP *scip, SCIP_VAR **vars, SCIP_Real **tableaurows, SCIP_ROWPREP *rowprep, SCIP_Real *rays, int *nrays, int *rayslppos, SCIP_Real *interpoints, SCIP_Bool usebounds, SCIP_Real *ad, SCIP_Bool *success) |
static SCIP_Bool | raysAreDependent (SCIP *scip, SCIP_Real *ray1, SCIP_Real *ray2, SCIP_Real *coef) |
static SCIP_RETCODE | findRho (SCIP *scip, SCIP_Real *rays, int nrays, int idx, SCIP_Real *interpoints, SCIP_VAR **vars, SCIP_Real *rho, SCIP_Bool usebounds, SCIP_Real *ad, SCIP_Bool *success) |
static SCIP_RETCODE | computeNegCutcoefs (SCIP *scip, SCIP_VAR **vars, SCIP_Real *rays, int nrays, int *rayslppos, SCIP_Real *interpoints, SCIP_ROWPREP *rowprep, SCIP_Bool usebounds, SCIP_Real *ad, SCIP_Bool *success) |
static SCIP_RETCODE | separateDeterminant (SCIP *scip, SCIP_SEPA *sepa, SCIP_SEPADATA *sepadata, SCIP_VAR *xik, SCIP_VAR *xil, SCIP_VAR *xjk, SCIP_VAR *xjl, SCIP_Bool *isxikdiag, SCIP_Bool *isxildiag, SCIP_Bool *isxjkdiag, SCIP_Bool *isxjldiag, int *basicvarpos2tableaurow, SCIP_HASHMAP *tableau, SCIP_RESULT *result) |
static SCIP_RETCODE | separatePoint (SCIP *scip, SCIP_SEPA *sepa, SCIP_RESULT *result) |
static | SCIP_DECL_SEPACOPY (sepaCopyMinor) |
static | SCIP_DECL_SEPAFREE (sepaFreeMinor) |
static | SCIP_DECL_SEPAINIT (sepaInitMinor) |
static | SCIP_DECL_SEPAEXIT (sepaExitMinor) |
static | SCIP_DECL_SEPAINITSOL (sepaInitsolMinor) |
static | SCIP_DECL_SEPAEXITSOL (sepaExitsolMinor) |
static | SCIP_DECL_SEPAEXECLP (sepaExeclpMinor) |
SCIP_RETCODE | SCIPincludeSepaInterminor (SCIP *scip) |
Macro Definition Documentation
◆ SEPA_NAME
#define SEPA_NAME "interminor" |
Definition at line 50 of file sepa_interminor.c.
Referenced by SCIP_DECL_SEPACOPY(), and SCIPincludeSepaInterminor().
◆ SEPA_DESC
#define SEPA_DESC "intersection cuts separator to ensure that 2x2 minors of X (= xx') have determinant 0" |
Definition at line 51 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ SEPA_PRIORITY
#define SEPA_PRIORITY 0 |
Definition at line 52 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ SEPA_FREQ
#define SEPA_FREQ -1 |
Definition at line 53 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ SEPA_MAXBOUNDDIST
#define SEPA_MAXBOUNDDIST 1.0 |
Definition at line 54 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ SEPA_USESSUBSCIP
#define SEPA_USESSUBSCIP FALSE |
does the separator use a secondary SCIP instance?
Definition at line 55 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ SEPA_DELAY
#define SEPA_DELAY FALSE |
should separation method be delayed, if other separators found cuts?
Definition at line 56 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ DEFAULT_MINCUTVIOL
#define DEFAULT_MINCUTVIOL 1e-4 |
default minimum required violation of a cut
Definition at line 58 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ DEFAULT_RANDSEED
#define DEFAULT_RANDSEED 157 |
default random seed
Definition at line 59 of file sepa_interminor.c.
Referenced by SCIP_DECL_SEPAINIT().
◆ DEFAULT_MAXROUNDS
#define DEFAULT_MAXROUNDS 10 |
maximal number of separation rounds per node (-1: unlimited)
Definition at line 60 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ DEFAULT_MAXROUNDSROOT
#define DEFAULT_MAXROUNDSROOT -1 |
maximal number of separation rounds in the root node (-1: unlimited)
Definition at line 61 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ BINSEARCH_MAXITERS
#define BINSEARCH_MAXITERS 120 |
default iteration limit for binary search
Definition at line 62 of file sepa_interminor.c.
Referenced by doBinarySearch(), and findRho().
◆ DEFAULT_USESTRENGTHENING
#define DEFAULT_USESTRENGTHENING FALSE |
default for using strengthend intersection cuts to separate
Definition at line 63 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
◆ DEFAULT_USEBOUNDS
#define DEFAULT_USEBOUNDS FALSE |
default for using nonnegativity bounds when separating
Definition at line 64 of file sepa_interminor.c.
Referenced by SCIPincludeSepaInterminor().
Function Documentation
◆ sepadataAddMinor()
|
static |
helper method to store a 2x2 minor in the separation data
- Parameters
-
scip SCIP data structure sepadata separator data auxvarxik auxiliary variable X_ik = x_i * x_k auxvarxil auxiliary variable X_il = x_i * x_l auxvarxjk auxiliary variable X_jk = x_j * x_k auxvarxjl auxiliary variable X_jl = x_j * x_l isauxvarxikdiag is X_ik diagonal? (i.e. i = k) isauxvarxildiag is X_il diagonal? (i.e. i = l) isauxvarxjkdiag is X_jk diagonal? (i.e. j = k) isauxvarxjldiag is X_jl diagonal? (i.e. j = l)
Definition at line 102 of file sepa_interminor.c.
References NULL, SCIP_CALL, SCIP_OKAY, SCIPcalcMemGrowSize(), SCIPcaptureVar(), SCIPdebugMsg, SCIPreallocBlockMemoryArray, and SCIPvarGetName().
Referenced by detectMinors().
◆ sepadataClear()
|
static |
helper method to clear separation data
- Parameters
-
scip SCIP data structure sepadata separator data
Definition at line 159 of file sepa_interminor.c.
References NULL, SCIP_CALL, SCIP_OKAY, SCIPdebugMsg, SCIPfreeBlockMemoryArrayNull, and SCIPreleaseVar().
Referenced by SCIP_DECL_SEPAEXITSOL().
◆ getMinorVars()
|
static |
helper method to get the variables associated to a minor
- Parameters
-
sepadata separator data idx index of the stored minor auxvarxik auxiliary variable X_ik = x_i * x_k auxvarxil auxiliary variable X_il = x_i * x_l auxvarxjk auxiliary variable X_jk = x_j * x_k auxvarxjl auxiliary variable X_jl = x_j * x_l isauxvarxikdiag is X_ik diagonal? (i.e. i = k) isauxvarxildiag is X_il diagonal? (i.e. i = l) isauxvarxjkdiag is X_jk diagonal? (i.e. j = k) isauxvarxjldiag is X_jl diagonal? (i.e. j = l)
Definition at line 190 of file sepa_interminor.c.
References NULL, and SCIP_OKAY.
Referenced by separatePoint().
◆ insertIndex()
|
static |
adds a new entry (i.e., auxvar) of in (row, col) of matrix M.
we have a matrix, M, indexed by the variables M(xi, xk) is the auxiliary variable of xi * xk if it exists We store, for each row of the matrix, the indices of the nonzero column entries (assoc with the given row) and the auxiliary variable for xi * xk The nonzero column entries are stored as an array (struct rowdata) So we have a hashmap mapping each variable (row of the matrix) with its array representing the nonzero entries of the row.
- Parameters
-
scip SCIP data structure rowmap hashmap of the rows of the matrix row variable corresponding to row of new entry col variable corresponding to column of new entry auxvar auxvar to insert into the matrix rowindices array of indices of all variables corresponding to a row nrows number of rows
Definition at line 231 of file sepa_interminor.c.
References rowdata::auxvars, rowdata::nvals, rowdata::rowidx, SCIP_CALL, SCIP_OKAY, SCIPallocBuffer, SCIPallocBufferArray, SCIPblkmem(), SCIPcalcMemGrowSize(), SCIPdebugMsg, SCIPhashmapCreate(), SCIPhashmapExists(), SCIPhashmapGetImage(), SCIPhashmapInsert(), SCIPreallocBufferArray, SCIPvarGetName(), SCIPvarGetProbindex(), rowdata::vals, and rowdata::valssize.
Referenced by detectMinors().
◆ detectMinors()
|
static |
method to detect and store principal minors
- Parameters
-
scip SCIP data structure sepadata separator data
Definition at line 295 of file sepa_interminor.c.
References rowdata::auxvars, FALSE, insertIndex(), NULL, rowdata::nvals, rowdata::rowidx, SCIP_Bool, SCIP_CALL, SCIP_EXPRITER_DFS, SCIP_OKAY, SCIP_Real, SCIPallocBufferArray, SCIPblkmem(), SCIPcomputeArraysIntersectionInt(), SCIPconsGetName(), SCIPconshdlrGetConss(), SCIPconshdlrGetNConss(), SCIPcreateExpriter(), SCIPdebugMsg, SCIPexprGetChildren(), SCIPexprGetNChildren(), SCIPexpriterGetNext(), SCIPexpriterInit(), SCIPexpriterIsEnd(), SCIPexpriterRestartDFS(), SCIPfindConshdlr(), SCIPfreeBufferArray, SCIPfreeBufferArrayNull, SCIPfreeExpriter(), SCIPgetExponentExprPow(), SCIPgetExprAuxVarNonlinear(), SCIPgetExprNonlinear(), SCIPgetNVars(), SCIPgetProbName(), SCIPgetTotalTime(), SCIPgetVars(), SCIPhashmapCreate(), SCIPhashmapFree(), SCIPhashmapGetImage(), SCIPisExprPower(), SCIPisExprProduct(), SCIPsortInt(), SCIPstatisticMessage, sepadataAddMinor(), TRUE, and rowdata::vals.
Referenced by SCIP_DECL_SEPAEXECLP().
◆ constructBasicVars2TableauRowMap()
|
static |
constructs map between lp position of a basic variable and its row in the tableau
- Parameters
-
scip SCIP data structure map buffer to store the map
Definition at line 510 of file sepa_interminor.c.
References SCIP_CALL, SCIP_OKAY, SCIPallocBufferArray, SCIPfreeBufferArray, SCIPgetLPBasisInd(), and SCIPgetNLPRows().
Referenced by separatePoint().
◆ computeRestrictionToRay()
|
static |
The restriction of the function representing the maximal S-free set to zlp + t * ray has the form sqrt(A t^2 + B t + C) - (D t + E). This function computes the coefficients A, B, C, D, E for the given ray.
- Parameters
-
scip SCIP data structure ray coefficients of ray vars variables coefs buffer to store A, B, C, D, and E of cases 1, 2, 3, or 4a coefs4b buffer to store A, B, C, D, and E of case 4b coefscondition buffer to store coefs for checking whether we are in case 4a or 4b usebounds TRUE if we want to separate non-negative bound ad coefs a and d for the hyperplane aTx + dTy <= 0 success FALSE if we need to abort generation because of numerics
Definition at line 539 of file sepa_interminor.c.
References a, ABS, b, FALSE, NULL, SCIP_OKAY, SCIP_Real, SCIPinfinity(), SCIPinfoMessage(), SCIPisHugeValue(), SCIPvarGetLPSol(), SQR, and TRUE.
◆ evalPhiAtRay()
|
static |
returns phi(zlp + t * ray) = sqrt(A t^2 + B t + C) - (D t + E)
- Parameters
-
scip SCIP data structure t argument of phi restricted to ray a value of A b value of B c value of C d value of D e value of E
Definition at line 700 of file sepa_interminor.c.
References QUAD, QUAD_ASSIGN, QUAD_SCALE, QUAD_TO_DBL, SCIP_Real, SCIPisInfinity(), SCIPquadprecProdDD, SCIPquadprecProdQD, SCIPquadprecSqrtQ, SCIPquadprecSquareD, SCIPquadprecSumQD, and SCIPquadprecSumQQ.
Referenced by computeRoot(), and doBinarySearch().
◆ doBinarySearch()
|
static |
helper function of computeRoot: we want phi to be <= 0
- Parameters
-
scip SCIP data structure a value of A b value of B c value of C d value of D e value of E sol buffer to store solution; also gives initial point
Definition at line 757 of file sepa_interminor.c.
References BINSEARCH_MAXITERS, evalPhiAtRay(), SCIP_Real, SCIPisFeasEQ(), and SCIPisFeasZero().
Referenced by computeRoot().
◆ isCase4a()
checks if we are in case 4a, i.e., if (num(xhat_{r+1}(zlp)) / E) * sqrt(A * tsol^2 + B * tsol + C) + w(ray) * tsol + num(yhat_{s+1}(zlp)) <= 0
- Parameters
-
tsol t in the above formula coefs coefficients A, B, C, D, and E of case 4a coefscondition extra coefficients needed for the evaluation of the condition: num(xhat_{r+1}(zlp)) / E; w(ray); num(yhat_{s+1}(zlp))
Definition at line 799 of file sepa_interminor.c.
References SQR.
Referenced by computeIntersectionPoint().
◆ computeRoot()
finds smallest positive root phi by finding the smallest positive root of (A - D^2) t^2 + (B - 2 D*E) t + (C - E^2) = 0
However, we are conservative and want a solution such that phi is negative, but close to 0; thus we correct the result with a binary search
- Parameters
-
scip SCIP data structure coefs value of A
Definition at line 817 of file sepa_interminor.c.
References a, b, doBinarySearch(), evalPhiAtRay(), SCIP_INTERVAL_INFINITY, SCIP_Real, SCIPinfinity(), SCIPintervalGetInf(), SCIPintervalIsEmpty(), SCIPintervalSetBounds(), and SCIPintervalSolveUnivariateQuadExpressionPositiveAllScalar().
Referenced by computeIntersectionPoint().
◆ computeIntersectionPoint()
|
static |
The maximal S-free set is gamma(z) <= 0; we find the intersection point of the ray ray
starting from zlp with the boundary of the S-free set. That is, we find t >= 0 such that gamma(zlp + t * ray) = 0.
In cases 1,2, and 3, gamma is of the form gamma(zlp + t * ray) = sqrt(A t^2 + B t + C) - (D t + E)
In the case 4 gamma is of the form gamma(zlp + t * ray) = sqrt(A t^2 + B t + C) - (D t + E) if some condition holds sqrt(A' t^2 + B' t + C') - (D' t + E') otherwise
It can be shown (given the special properties of gamma) that the smallest positive root of each function of the form sqrt(a t^2 + b t + c) - (d t + e) is the same as the smallest positive root of the quadratic equation: (sqrt(a t^2 + b t + c) - (d t + e)) * (sqrt(a t^2 + b t + c) + (d t + e)) = 0 <==> (a - d^2) t^2 + (b - 2 d*e) t + (c - e^2) = 0
So, in cases 1, 2, and 3, this function just returns the solution of the above equation. In case 4, it first solves the equation assuming we are in the first piece. If there is no solution, then the second piece can't have a solution (first piece >= second piece for all t) Then we check if the solution satisfies the condition. If it doesn't then we solve the equation for the second piece. If it has a solution, then it has to be the solution.
- Parameters
-
scip SCIP data structure usebounds whether we are in case 4 or not coefs values of A, B, C, D, and E of cases 1, 2, 3, or 4a coefs4b values of A, B, C, D, and E of case 4b coefscondition values needed to evaluate condition of case 4
Definition at line 904 of file sepa_interminor.c.
References computeRoot(), isCase4a(), MAX, NULL, SCIP_Real, and SCIPisInfinity().
◆ addColToCut()
|
static |
adds cutcoef * (col - col*) to rowprep
- Parameters
-
scip SCIP data structure rowprep rowprep to store intersection cut cutcoef cut coefficient col column to add to rowprep
Definition at line 953 of file sepa_interminor.c.
References NULL, SCIP_BASESTAT_LOWER, SCIP_CALL, SCIP_OKAY, SCIPaddRowprepTerm(), SCIPcolGetBasisStatus(), SCIPcolGetPrimsol(), SCIPcolGetVar(), SCIPinfoMessage(), SCIProwprepAddConstant(), SCIPvarGetLbLocal(), SCIPvarGetName(), and SCIPvarGetUbLocal().
Referenced by addCols(), and computeNegCutcoefs().
◆ addRowToCut()
|
static |
adds cutcoef * (slack - slack*) to rowprep
row is lhs <= <coefs, vars> + constant <= rhs, thus slack is defined by slack + <coefs, vars> + constant = side If row (slack) is at upper, it means that <coefs,vars*> + constant = rhs, and so slack* = side - rhs –> slack - slack* = rhs - <coefs, vars> - constant. If row (slack) is at lower, then <coefs,vars*> + constant = lhs, and so slack* = side - lhs –> slack - slack* = lhs - <coefs, vars> - constant.
- Parameters
-
scip SCIP data structure rowprep rowprep to store intersection cut cutcoef cut coefficient row row, whose slack we are adding to rowprep success buffer to store whether the row is nonbasic enough
Definition at line 985 of file sepa_interminor.c.
References FALSE, NULL, SCIP_BASESTAT_LOWER, SCIP_CALL, SCIP_OKAY, SCIP_Real, SCIPaddRowprepTerm(), SCIPcolGetVar(), SCIPgetRowActivity(), SCIPinfoMessage(), SCIPisFeasEQ(), SCIPisInfinity(), SCIProwGetBasisStatus(), SCIProwGetCols(), SCIProwGetConstant(), SCIProwGetLhs(), SCIProwGetLPPos(), SCIProwGetNLPNonz(), SCIProwGetRhs(), SCIProwGetVals(), and SCIProwprepAddConstant().
Referenced by addRows(), and computeNegCutcoefs().
◆ getTableauRows()
|
static |
get the tableau rows of the variables in vars
- Parameters
-
scip SCIP data structure vars variables in the minor basicvarpos2tableaurow map between basic var and its tableau row tableau map between var an its tableau row tableaurows buffer to store tableau row success set to TRUE if no variable had basisstat = ZERO
Definition at line 1046 of file sepa_interminor.c.
References FALSE, NULL, SCIP_BASESTAT_BASIC, SCIP_BASESTAT_ZERO, SCIP_CALL, SCIP_OKAY, SCIP_Real, SCIPallocBufferArray, SCIPcolGetBasisStatus(), SCIPcolGetLPPos(), SCIPgetLPBInvARow(), SCIPgetLPBInvRow(), SCIPgetNLPCols(), SCIPgetNLPRows(), SCIPhashmapExists(), SCIPhashmapGetImage(), SCIPhashmapInsert(), SCIPvarGetCol(), and TRUE.
Referenced by separateDeterminant().
◆ addCols()
|
static |
computes the cut coefs of the non-basic (non-slack) variables (correspond to cols) and adds them to the intersection cut
- Parameters
-
scip SCIP data structure vars variables tableaurows tableau rows corresponding to the variables in vars rowprep store cut rays buffer to store rays nrays pointer to store number of nonzero rays rayslppos buffer to store lppos of nonzero rays interpoints buffer to store intersection points or NULL if not needed usebounds TRUE if we want to separate non-negative bound ad coefs a and d for the hyperplane aTx + dTy <= 0 success pointer to store whether the generation of cutcoefs was successful
Definition at line 1110 of file sepa_interminor.c.
References addColToCut(), computeIntersectionPoint(), computeRestrictionToRay(), FALSE, NULL, SCIP_BASESTAT_LOWER, SCIP_BASESTAT_UPPER, SCIP_BASESTAT_ZERO, SCIP_Bool, SCIP_CALL, SCIP_OKAY, SCIP_Real, SCIPcolGetBasisStatus(), SCIPgetLPCols(), SCIPgetNLPCols(), SCIPisInfinity(), SCIPisZero(), SCIPvarGetCol(), and TRUE.
Referenced by separateDeterminant().
◆ addRows()
|
static |
computes the cut coefs of the non-basic slack variables (correspond to rows) and adds them to the intersection cut
- Parameters
-
scip SCIP data structure vars variables tableaurows tableau rows corresponding to the variables in vars rowprep store cut rays buffer to store rays nrays pointer to store number of nonzero rays rayslppos buffer to store lppos of nonzero rays interpoints buffer to store intersection points or NULL if not needed usebounds TRUE if we want to separate non-negative bound ad coefs a and d for the hyperplane aTx + dTy <= 0 success pointer to store whether the generation of cutcoefs was successful
Definition at line 1216 of file sepa_interminor.c.
References addRowToCut(), computeIntersectionPoint(), computeRestrictionToRay(), FALSE, NULL, SCIP_BASESTAT_LOWER, SCIP_BASESTAT_UPPER, SCIP_BASESTAT_ZERO, SCIP_Bool, SCIP_CALL, SCIP_OKAY, SCIP_Real, SCIPgetLPRows(), SCIPgetNLPCols(), SCIPgetNLPRows(), SCIPisInfinity(), SCIPisZero(), SCIProwGetBasisStatus(), and TRUE.
Referenced by separateDeterminant().
◆ raysAreDependent()
|
static |
- Parameters
-
scip SCIP data structure ray1 coefficients of ray 1 ray2 coefficients of ray 2 coef pointer to store coef (s.t. r1 = coef * r2) in case rays are dependent
Definition at line 1326 of file sepa_interminor.c.
References FALSE, SCIPisFeasEQ(), SCIPisZero(), and TRUE.
Referenced by findRho().
◆ findRho()
|
static |
finds the smallest negative steplength for the current ray r_idx such that the combination of r_idx with all rays not in the recession cone is in the recession cone
- Parameters
-
scip SCIP data structure rays rays nrays number of nonzero rays idx index of current ray we want to find rho for interpoints intersection points of nonzero rays vars variables rho pointer to store the optimal rho usebounds TRUE if we want to separate non-negative bound ad coefs a and d for the hyperplane aTx + dTy <= 0 success TRUE if computation of rho was successful
Definition at line 1363 of file sepa_interminor.c.
References BINSEARCH_MAXITERS, computeIntersectionPoint(), computeRestrictionToRay(), raysAreDependent(), SCIP_CALL, SCIP_OKAY, SCIP_Real, SCIPinfinity(), SCIPisEQ(), SCIPisInfinity(), SCIPisZero(), and TRUE.
Referenced by computeNegCutcoefs().
◆ computeNegCutcoefs()
|
static |
computes negative steplengths for the rays that are in the recession cone of the S-free set, i.e., which have an infinite intersection point.
- Parameters
-
scip SCIP data structure vars variables rays rays nrays number of nonzero rays rayslppos lppos of nonzero rays interpoints intersection points rowprep rowprep for the generated cut usebounds TRUE if we want to separate non-negative bound ad coefs a and d for the hyperplane aTx + dTy <= 0 success if a cut candidate could be computed
Definition at line 1467 of file sepa_interminor.c.
References addColToCut(), addRowToCut(), findRho(), SCIP_BASESTAT_LOWER, SCIP_BASESTAT_UPPER, SCIP_CALL, SCIP_OKAY, SCIP_Real, SCIPcolGetBasisStatus(), SCIPgetLPCols(), SCIPgetLPRows(), SCIPisInfinity(), SCIProwGetBasisStatus(), and TRUE.
Referenced by separateDeterminant().
◆ separateDeterminant()
|
static |
separates cuts for stored principal minors
- Parameters
-
scip SCIP data structure sepa separator sepadata separator data xik variable X_ik = x_i * x_k xil variable X_il = x_i * x_l xjk variable X_jk = x_j * x_k xjl variable X_jl = x_j * x_l isxikdiag is X_ik diagonal? (i.e. i = k) isxildiag is X_il diagonal? (i.e. i = l) isxjkdiag is X_jk diagonal? (i.e. j = k) isxjldiag is X_jl diagonal? (i.e. j = l) basicvarpos2tableaurow map from basic var to its tableau row tableau map from var to its tableau row result pointer to store the result of the separation call
Definition at line 1538 of file sepa_interminor.c.
References addCols(), addRows(), computeNegCutcoefs(), FALSE, getTableauRows(), NULL, SCIP_Bool, SCIP_CALL, SCIP_CUTOFF, SCIP_OKAY, SCIP_Real, SCIP_SEPARATED, SCIP_SIDETYPE_LEFT, SCIPaddRow(), SCIPallocBufferArray, SCIPcleanupRowprep(), SCIPcreateRowprep(), SCIPfreeBuffer, SCIPfreeRowprep(), SCIPgetCutEfficacy(), SCIPgetNLPCols(), SCIPgetNLPRows(), SCIPgetRowprepRowSepa(), SCIPisFeasNegative(), SCIPmergeRowprepTerms(), SCIPreleaseRow(), SCIProwprepAddSide(), SCIPvarGetLPSol(), and TRUE.
Referenced by separatePoint().
◆ separatePoint()
|
static |
separates cuts for stored principal minors
- Parameters
-
scip SCIP data structure sepa separator result pointer to store the result of the separation call
Definition at line 1682 of file sepa_interminor.c.
References constructBasicVars2TableauRowMap(), getMinorVars(), NULL, SCIP_Bool, SCIP_CALL, SCIP_CUTOFF, SCIP_DIDNOTFIND, SCIP_DIDNOTRUN, SCIP_OKAY, SCIP_Real, SCIPallocBufferArray, SCIPblkmem(), SCIPfreeBufferArray, SCIPfreeBufferArrayNull, SCIPgetNLPCols(), SCIPgetNVars(), SCIPhashmapCreate(), SCIPhashmapEntryGetImage(), SCIPhashmapFree(), SCIPhashmapGetEntry(), SCIPhashmapGetNEntries(), SCIPisFeasNegative(), SCIPisFeasPositive(), SCIPisFeasZero(), SCIPsepaGetData(), SCIPvarGetLPSol(), and separateDeterminant().
Referenced by SCIP_DECL_SEPAEXECLP().
◆ SCIP_DECL_SEPACOPY()
|
static |
copy method for separator plugins (called when SCIP copies plugins)
Definition at line 1795 of file sepa_interminor.c.
References NULL, SCIP_CALL, SCIP_OKAY, SCIPincludeSepaInterminor(), SCIPsepaGetName(), and SEPA_NAME.
◆ SCIP_DECL_SEPAFREE()
|
static |
destructor of separator to free user data (called when SCIP is exiting)
Definition at line 1810 of file sepa_interminor.c.
References NULL, SCIP_OKAY, SCIPfreeBlockMemory, SCIPsepaGetData(), and SCIPsepaSetData().
◆ SCIP_DECL_SEPAINIT()
|
static |
initialization method of separator (called after problem was transformed)
Definition at line 1830 of file sepa_interminor.c.
References DEFAULT_RANDSEED, NULL, SCIP_CALL, SCIP_OKAY, SCIPcreateRandom(), SCIPsepaGetData(), and TRUE.
◆ SCIP_DECL_SEPAEXIT()
|
static |
deinitialization method of separator (called before transformed problem is freed)
Definition at line 1848 of file sepa_interminor.c.
References NULL, SCIP_OKAY, SCIPfreeRandom(), and SCIPsepaGetData().
◆ SCIP_DECL_SEPAINITSOL()
|
static |
solving process initialization method of separator (called when branch and bound process is about to begin)
Definition at line 1866 of file sepa_interminor.c.
References SCIP_OKAY.
◆ SCIP_DECL_SEPAEXITSOL()
|
static |
solving process deinitialization method of separator (called before branch and bound process data is freed)
Definition at line 1874 of file sepa_interminor.c.
References NULL, SCIP_CALL, SCIP_OKAY, SCIPsepaGetData(), and sepadataClear().
◆ SCIP_DECL_SEPAEXECLP()
|
static |
LP solution separation method of separator
Definition at line 1890 of file sepa_interminor.c.
References detectMinors(), NULL, SCIP_CALL, SCIP_OKAY, SCIPdebugMsg, SCIPgetDepth(), SCIPisIpoptAvailableIpopt(), SCIPsepaGetData(), SCIPsepaGetNCallsAtNode(), and separatePoint().