Scippy

SCIP

Solving Constraint Integer Programs

Detailed Description

minor separator with intersection cuts

Author
Felipe Serrano
Antonia Chmiela

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 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

helper method to store a 2x2 minor in the separation data

Parameters
scipSCIP data structure
sepadataseparator data
auxvarxikauxiliary variable X_ik = x_i * x_k
auxvarxilauxiliary variable X_il = x_i * x_l
auxvarxjkauxiliary variable X_jk = x_j * x_k
auxvarxjlauxiliary variable X_jl = x_j * x_l
isauxvarxikdiagis X_ik diagonal? (i.e. i = k)
isauxvarxildiagis X_il diagonal? (i.e. i = l)
isauxvarxjkdiagis X_jk diagonal? (i.e. j = k)
isauxvarxjldiagis 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 SCIP_RETCODE sepadataClear ( SCIP scip,
SCIP_SEPADATA sepadata 
)
static

helper method to clear separation data

Parameters
scipSCIP data structure
sepadataseparator 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 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

helper method to get the variables associated to a minor

Parameters
sepadataseparator data
idxindex of the stored minor
auxvarxikauxiliary variable X_ik = x_i * x_k
auxvarxilauxiliary variable X_il = x_i * x_l
auxvarxjkauxiliary variable X_jk = x_j * x_k
auxvarxjlauxiliary variable X_jl = x_j * x_l
isauxvarxikdiagis X_ik diagonal? (i.e. i = k)
isauxvarxildiagis X_il diagonal? (i.e. i = l)
isauxvarxjkdiagis X_jk diagonal? (i.e. j = k)
isauxvarxjldiagis 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 SCIP_RETCODE insertIndex ( SCIP scip,
SCIP_HASHMAP rowmap,
SCIP_VAR row,
SCIP_VAR col,
SCIP_VAR auxvar,
int *  rowindices,
int *  nrows 
)
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
scipSCIP data structure
rowmaphashmap of the rows of the matrix
rowvariable corresponding to row of new entry
colvariable corresponding to column of new entry
auxvarauxvar to insert into the matrix
rowindicesarray of indices of all variables corresponding to a row
nrowsnumber 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()

◆ constructBasicVars2TableauRowMap()

static SCIP_RETCODE constructBasicVars2TableauRowMap ( SCIP scip,
int *  map 
)
static

constructs map between lp position of a basic variable and its row in the tableau

Parameters
scipSCIP data structure
mapbuffer 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 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

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
scipSCIP data structure
raycoefficients of ray
varsvariables
coefsbuffer to store A, B, C, D, and E of cases 1, 2, 3, or 4a
coefs4bbuffer to store A, B, C, D, and E of case 4b
coefsconditionbuffer to store coefs for checking whether we are in case 4a or 4b
useboundsTRUE if we want to separate non-negative bound
adcoefs a and d for the hyperplane aTx + dTy <= 0
successFALSE 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.

Referenced by addCols(), addRows(), and findRho().

◆ evalPhiAtRay()

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

returns phi(zlp + t * ray) = sqrt(A t^2 + B t + C) - (D t + E)

Parameters
scipSCIP data structure
targument of phi restricted to ray
avalue of A
bvalue of B
cvalue of C
dvalue of D
evalue 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 void doBinarySearch ( SCIP scip,
SCIP_Real  a,
SCIP_Real  b,
SCIP_Real  c,
SCIP_Real  d,
SCIP_Real  e,
SCIP_Real sol 
)
static

helper function of computeRoot: we want phi to be <= 0

Parameters
scipSCIP data structure
avalue of A
bvalue of B
cvalue of C
dvalue of D
evalue of E
solbuffer 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()

static SCIP_Real isCase4a ( SCIP_Real  tsol,
SCIP_Real coefs,
SCIP_Real coefscondition 
)
static

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
tsolt in the above formula
coefscoefficients A, B, C, D, and E of case 4a
coefsconditionextra 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()

static SCIP_Real computeRoot ( SCIP scip,
SCIP_Real coefs 
)
static

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
scipSCIP data structure
coefsvalue 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 SCIP_Real computeIntersectionPoint ( SCIP scip,
SCIP_Bool  usebounds,
SCIP_Real coefs,
SCIP_Real coefs4b,
SCIP_Real coefscondition 
)
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
scipSCIP data structure
useboundswhether we are in case 4 or not
coefsvalues of A, B, C, D, and E of cases 1, 2, 3, or 4a
coefs4bvalues of A, B, C, D, and E of case 4b
coefsconditionvalues 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().

Referenced by addCols(), addRows(), and findRho().

◆ addColToCut()

static SCIP_RETCODE addColToCut ( SCIP scip,
SCIP_ROWPREP rowprep,
SCIP_Real  cutcoef,
SCIP_COL col 
)
static

adds cutcoef * (col - col*) to rowprep

Parameters
scipSCIP data structure
rowpreprowprep to store intersection cut
cutcoefcut coefficient
colcolumn 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 SCIP_RETCODE addRowToCut ( SCIP scip,
SCIP_ROWPREP rowprep,
SCIP_Real  cutcoef,
SCIP_ROW row,
SCIP_Bool success 
)
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
scipSCIP data structure
rowpreprowprep to store intersection cut
cutcoefcut coefficient
rowrow, whose slack we are adding to rowprep
successbuffer 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 SCIP_RETCODE getTableauRows ( SCIP scip,
SCIP_VAR **  vars,
int *  basicvarpos2tableaurow,
SCIP_HASHMAP tableau,
SCIP_Real **  tableaurows,
SCIP_Bool success 
)
static

get the tableau rows of the variables in vars

Parameters
scipSCIP data structure
varsvariables in the minor
basicvarpos2tableaurowmap between basic var and its tableau row
tableaumap between var an its tableau row
tableaurowsbuffer to store tableau row
successset 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 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

computes the cut coefs of the non-basic (non-slack) variables (correspond to cols) and adds them to the intersection cut

Parameters
scipSCIP data structure
varsvariables
tableaurowstableau rows corresponding to the variables in vars
rowprepstore cut
raysbuffer to store rays
nrayspointer to store number of nonzero rays
rayslpposbuffer to store lppos of nonzero rays
interpointsbuffer to store intersection points or NULL if not needed
useboundsTRUE if we want to separate non-negative bound
adcoefs a and d for the hyperplane aTx + dTy <= 0
successpointer 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 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

computes the cut coefs of the non-basic slack variables (correspond to rows) and adds them to the intersection cut

Parameters
scipSCIP data structure
varsvariables
tableaurowstableau rows corresponding to the variables in vars
rowprepstore cut
raysbuffer to store rays
nrayspointer to store number of nonzero rays
rayslpposbuffer to store lppos of nonzero rays
interpointsbuffer to store intersection points or NULL if not needed
useboundsTRUE if we want to separate non-negative bound
adcoefs a and d for the hyperplane aTx + dTy <= 0
successpointer 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 SCIP_Bool raysAreDependent ( SCIP scip,
SCIP_Real ray1,
SCIP_Real ray2,
SCIP_Real coef 
)
static
Parameters
scipSCIP data structure
ray1coefficients of ray 1
ray2coefficients of ray 2
coefpointer 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 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

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
scipSCIP data structure
raysrays
nraysnumber of nonzero rays
idxindex of current ray we want to find rho for
interpointsintersection points of nonzero rays
varsvariables
rhopointer to store the optimal rho
useboundsTRUE if we want to separate non-negative bound
adcoefs a and d for the hyperplane aTx + dTy <= 0
successTRUE 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 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

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
scipSCIP data structure
varsvariables
raysrays
nraysnumber of nonzero rays
rayslpposlppos of nonzero rays
interpointsintersection points
rowpreprowprep for the generated cut
useboundsTRUE if we want to separate non-negative bound
adcoefs a and d for the hyperplane aTx + dTy <= 0
successif 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 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

separates cuts for stored principal minors

Parameters
scipSCIP data structure
sepaseparator
sepadataseparator data
xikvariable X_ik = x_i * x_k
xilvariable X_il = x_i * x_l
xjkvariable X_jk = x_j * x_k
xjlvariable X_jl = x_j * x_l
isxikdiagis X_ik diagonal? (i.e. i = k)
isxildiagis X_il diagonal? (i.e. i = l)
isxjkdiagis X_jk diagonal? (i.e. j = k)
isxjldiagis X_jl diagonal? (i.e. j = l)
basicvarpos2tableaurowmap from basic var to its tableau row
tableaumap from var to its tableau row
resultpointer 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()

◆ SCIP_DECL_SEPACOPY()

static SCIP_DECL_SEPACOPY ( sepaCopyMinor  )
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 SCIP_DECL_SEPAFREE ( sepaFreeMinor  )
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 SCIP_DECL_SEPAINIT ( sepaInitMinor  )
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 SCIP_DECL_SEPAEXIT ( sepaExitMinor  )
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 SCIP_DECL_SEPAINITSOL ( sepaInitsolMinor  )
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 SCIP_DECL_SEPAEXITSOL ( sepaExitsolMinor  )
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 SCIP_DECL_SEPAEXECLP ( sepaExeclpMinor  )
static