// Copyright (C) 2002, International Business Machines // Corporation and others. All Rights Reserved. #ifndef CoinPresolveMatrix_H #define CoinPresolveMatrix_H #include "CoinPragma.hpp" #include "CoinPackedMatrix.hpp" #include "CoinMessage.hpp" #include "CoinTime.hpp" #include #include #include #include /*! \file Declarations for CoinPresolveMatrix and CoinPostsolveMatrix and their common base class CoinPrePostsolveMatrix. Also declarations for CoinPresolveAction and a number of non-member utility functions. */ #if defined(_MSC_VER) // Avoid MS Compiler problem in recognizing type to delete // by casting to type. #define deleteAction(array,type) delete [] ((type) array) #else #define deleteAction(array,type) delete [] array #endif /*! \brief Zero tolerance OSL had a fixed zero tolerance; we still use that here. */ const double ZTOLDP = 1e-12; // But use a different one if we are doing doubletons etc const double ZTOLDP2 = 1e-10; //#define PRESOLVE_DEBUG 1 // Debugging macros/functions #if PRESOLVE_DEBUG || PRESOLVE_CONSISTENCY #define PRESOLVE_STMT(s) s #define PRESOLVEASSERT(x) \ ((x) ? 1 : \ ((std::cerr << "FAILED ASSERTION at line " \ << __LINE__ << ": " #x "\n"), abort(), 0)) inline void DIE(const char *s) { std::cout<i.e., an object that contains instructions for backing out the presolve transform to recover the original problem. These postsolve objects are accumulated in a linked list, with each successive presolve action adding its postsolve action to the head of the list. The end result of all this is a presolved problem object, and a list of postsolve objects. The presolved problem object is then handed to a solver for optimization, and the problem object augmented with the results. The list of postsolve objects is then traversed. Each of them (un)modifies the problem object, with the end result being the original problem, augmented with solution information. The problem object representation is CoinPrePostsolveMatrix and subclasses. Check there for details. The \c CoinPresolveAction class and subclasses represent the presolve and postsolve objects. In spite of the name, the only information held in a \c CoinPresolveAction object is the information needed to postsolve (i.e., the information needed to back out the presolve transformation). This information is not expected to change, so the fields are all \c const. A subclass of \c CoinPresolveAction, implementing a specific pre/postsolve action, is expected to declare a static function that attempts to perform a presolve transformation. This function will be handed a CoinPresolveMatrix to transform, and a pointer to the head of the list of postsolve objects. If the transform is successful, the function will create a new \c CoinPresolveAction object, link it at the head of the list of postsolve objects, and return a pointer to the postsolve object it has just created. Otherwise, it should return 0. It is expected that these static functions will be the only things that can create new \c CoinPresolveAction objects; this is expressed by making each subclass' constructor(s) private. Every subclass must also define a \c postsolve method. This function will be handed a CoinPostsolveMatrix to transform. It is the client's responsibility to implement presolve and postsolve driver routines. See OsiPresolve for examples. \note Since the only fields in a \c CoinPresolveAction are \c const, anything one can do with a variable declared \c CoinPresolveAction* can also be done with a variable declared \c const \c CoinPresolveAction* It is expected that all derived subclasses of \c CoinPresolveAction also have this property. */ class CoinPresolveAction { public: /*! \brief Stub routine to throw exceptions. Exceptions are inefficient, particularly with g++. Even with xlC, the use of exceptions adds a long prologue to a routine. Therefore, rather than use throw directly in the routine, I use it in a stub routine. */ static void throwCoinError(const char *error, const char *ps_routine) { throw CoinError(error, ps_routine, "CoinPresolve"); } /*! \brief The next presolve transformation Set at object construction. */ const CoinPresolveAction *next; /*! \brief Construct a postsolve object and add it to the transformation list. This is an `add to head' operation. This object will point to the one passed as the parameter. */ CoinPresolveAction(const CoinPresolveAction *next) : next(next) {} /*! \brief A name for debug printing. It is expected that the name is not stored in the transform itself. */ virtual const char *name() const = 0; /*! \brief Apply the postsolve transformation for this particular presolve action. */ virtual void postsolve(CoinPostsolveMatrix *prob) const = 0; /*! \brief Virtual destructor. */ virtual ~CoinPresolveAction() {} }; /* These are needed for OSI-aware constructors associated with CoinPrePostsolveMatrix, CoinPresolveMatrix, and CoinPostsolveMatrix. */ class ClpSimplex; class OsiSolverInterface; /* CoinWarmStartBasis is required for methods in CoinPrePostsolveMatrix that accept/return a CoinWarmStartBasis object. */ class CoinWarmStartBasis ; /*! \class CoinPrePostsolveMatrix \brief Collects all the information about the problem that is needed in both presolve and postsolve. In a bit more detail, a column-major representation of the constraint matrix and upper and lower bounds on variables and constraints, plus row and column solutions, reduced costs, and status. There's also a set of arrays holding the original row and column numbers. As presolve and postsolve transform the matrix, it will occasionally be necessary to expand the number of entries in a column. There are two aspects:
  • During postsolve, the constraint system is expected to grow as the smaller presolved system is transformed back to the original system.
  • During both pre- and postsolve, transforms can increase the number of coefficients in a row or column. (See the variable substitution, doubleton, and tripleton transforms.)
The first is addressed by the members #ncols0_, #nrows0_, and #nelems0_. These should be set (via constructor parameters) to values large enough for the largest size taken on by the constraint system. Typically, this will be the size of the original constraint system. The second is addressed by a generous allocation of extra (empty) space for the arrays used to hold coefficients and row indices. When columns must be expanded, they are moved into the empty space. When it is used up, the arrays are compacted. When compaction fails to produce sufficient space, presolve/postsolve will fail. CoinPrePostsolveMatrix isn't really intended to be used `bare' --- the expectation is that it'll be used through CoinPresolveMatrix or CoinPostsolveMatrix. Some of the functions needed to load a problem are defined in the derived classes. When CoinPresolve is applied when reoptimising, we need to be prepared to accept a basis and modify it in step with the presolve actions (otherwise we throw away all the advantages of warm start for reoptimization). But other solution components (#acts_, #rowduals_, #sol_, and #rcosts_) are needed only for postsolve, where they're used in places to determine the proper action(s) when restoring rows or columns. If presolve is provided with a solution, it will modify it in step with the presolve actions. Moving the solution components from CoinPrePostsolveMatrix to CoinPostsolveMatrix would break a lot of code. It's not clear that it's worth it, and it would preclude upgrades to the presolve side that might make use of any of these. -- lh, 080501 -- */ class CoinPrePostsolveMatrix { public: /*! \name Constructors & Destructors */ //@{ /*! \brief `Native' constructor This constructor creates an empty object which must then be loaded. On the other hand, it doesn't assume that the client is an OsiSolverInterface. */ CoinPrePostsolveMatrix(int ncols_alloc, int nrows_alloc, CoinBigIndex nelems_alloc) ; /*! \brief Generic OSI constructor See OSI code for the definition. */ CoinPrePostsolveMatrix(const OsiSolverInterface * si, int ncols_, int nrows_, CoinBigIndex nelems_); /*! ClpOsi constructor See Clp code for the definition. */ CoinPrePostsolveMatrix(const ClpSimplex * si, int ncols_, int nrows_, CoinBigIndex nelems_, double bulkRatio); /// Destructor ~CoinPrePostsolveMatrix(); //@} /*! \brief Enum for status of various sorts Matches CoinWarmStartBasis::Status and adds superBasic. Most code that converts between CoinPrePostsolveMatrix::Status and CoinWarmStartBasis::Status will break if this correspondence is broken. superBasic is an unresolved problem: there's no analogue in CoinWarmStartBasis::Status. */ enum Status { isFree = 0x00, basic = 0x01, atUpperBound = 0x02, atLowerBound = 0x03, superBasic = 0x04 }; /*! \name Functions to work with variable status Functions to work with the CoinPrePostsolveMatrix::Status enum and related vectors. */ //@{ /// Set row status (i.e., status of artificial for this row) inline void setRowStatus(int sequence, Status status) { unsigned char & st_byte = rowstat_[sequence]; st_byte = static_cast(st_byte & (~7)) ; st_byte = static_cast(st_byte | status) ; } /// Get row status inline Status getRowStatus(int sequence) const {return static_cast (rowstat_[sequence]&7);} /// Check if artificial for this row is basic inline bool rowIsBasic(int sequence) const {return (static_cast (rowstat_[sequence]&7)==basic);} /// Set column status (i.e., status of primal variable) inline void setColumnStatus(int sequence, Status status) { unsigned char & st_byte = colstat_[sequence]; st_byte = static_cast(st_byte & (~7)) ; st_byte = static_cast(st_byte | status) ; # ifdef PRESOLVE_DEBUG switch (status) { case isFree: { if (clo_[sequence] > -PRESOLVE_INF || cup_[sequence] < PRESOLVE_INF) { std::cout << "Bad status: Var " << sequence << " isFree, lb = " << clo_[sequence] << ", ub = " << cup_[sequence] << std::endl ; } break ; } case basic: { break ; } case atUpperBound: { if (cup_[sequence] >= PRESOLVE_INF) { std::cout << "Bad status: Var " << sequence << " atUpperBound, lb = " << clo_[sequence] << ", ub = " << cup_[sequence] << std::endl ; } break ; } case atLowerBound: { if (clo_[sequence] <= -PRESOLVE_INF) { std::cout << "Bad status: Var " << sequence << " atLowerBound, lb = " << clo_[sequence] << ", ub = " << cup_[sequence] << std::endl ; } break ; } case superBasic: { if (clo_[sequence] <= -PRESOLVE_INF && cup_[sequence] >= PRESOLVE_INF) { std::cout << "Bad status: Var " << sequence << " superBasic, lb = " << clo_[sequence] << ", ub = " << cup_[sequence] << std::endl ; } break ; } default: { assert(false) ; break ; } } # endif } /// Get column (structural variable) status inline Status getColumnStatus(int sequence) const {return static_cast (colstat_[sequence]&7);} /// Check if column (structural variable) is basic inline bool columnIsBasic(int sequence) const {return (static_cast (colstat_[sequence]&7)==basic);} /*! \brief Set status of row (artificial variable) to the correct nonbasic status given bounds and current value */ void setRowStatusUsingValue(int iRow); /*! \brief Set status of column (structural variable) to the correct nonbasic status given bounds and current value */ void setColumnStatusUsingValue(int iColumn); /*! \brief Set column (structural variable) status vector */ void setStructuralStatus(const char *strucStatus, int lenParam) ; /*! \brief Set row (artificial variable) status vector */ void setArtificialStatus(const char *artifStatus, int lenParam) ; /*! \brief Set the status of all variables from a basis */ void setStatus(const CoinWarmStartBasis *basis) ; /*! \brief Get status in the form of a CoinWarmStartBasis */ CoinWarmStartBasis *getStatus() ; /*! \brief Return a print string for status of a column (structural variable) */ const char *columnStatusString(int j) const ; /*! \brief Return a print string for status of a row (artificial variable) */ const char *rowStatusString(int i) const ; //@} /*! \name Functions to load problem and solution information These functions can be used to load portions of the problem definition and solution. See also the CoinPresolveMatrix and CoinPostsolveMatrix classes. */ //@{ /// Set the objective function offset for the original system. void setObjOffset(double offset) ; /*! \brief Set the objective sense (max/min) Coded as 1.0 for min, -1.0 for max. */ void setObjSense(double objSense) ; /// Set the primal feasibility tolerance void setPrimalTolerance(double primTol) ; /// Set the dual feasibility tolerance void setDualTolerance(double dualTol) ; /// Set column lower bounds void setColLower(const double *colLower, int lenParam) ; /// Set column upper bounds void setColUpper(const double *colUpper, int lenParam) ; /// Set column solution void setColSolution(const double *colSol, int lenParam) ; /// Set objective coefficients void setCost(const double *cost, int lenParam) ; /// Set reduced costs void setReducedCost(const double *redCost, int lenParam) ; /// Set row lower bounds void setRowLower(const double *rowLower, int lenParam) ; /// Set row upper bounds void setRowUpper(const double *rowUpper, int lenParam) ; /// Set row solution void setRowPrice(const double *rowSol, int lenParam) ; /// Set row activity void setRowActivity(const double *rowAct, int lenParam) ; //@} /*! \name Functions to retrieve problem and solution information */ //@{ /// Get current number of columns inline int getNumCols() { return (ncols_) ; } /// Get current number of rows inline int getNumRows() { return (nrows_) ; } /// Get current number of non-zero coefficients inline int getNumElems() { return (nelems_) ; } /// Get column start vector for column-major packed matrix inline const CoinBigIndex *getColStarts() const { return (mcstrt_) ; } /// Get column length vector for column-major packed matrix inline const int *getColLengths() const { return (hincol_) ; } /// Get vector of row indices for column-major packed matrix inline const int *getRowIndicesByCol() const { return (hrow_) ; } /// Get vector of elements for column-major packed matrix inline const double *getElementsByCol() const { return (colels_) ; } /// Get column lower bounds inline const double *getColLower() const { return (clo_) ; } /// Get column upper bounds inline const double *getColUpper() const { return (cup_) ; } /// Get objective coefficients inline const double *getCost() const { return (cost_) ; } /// Get row lower bounds inline const double *getRowLower() const { return (rlo_) ; } /// Get row upper bounds inline const double *getRowUpper() const { return (rup_) ; } /// Get column solution (primal variable values) inline const double *getColSolution() const { return (sol_) ; } /// Get row activity (constraint lhs values) inline const double *getRowActivity() const { return (acts_) ; } /// Get row solution (dual variables) inline const double *getRowPrice() const { return (rowduals_) ; } /// Get reduced costs inline const double *getReducedCost() const { return (rcosts_) ; } /// Count empty columns inline int countEmptyCols() { int empty = 0 ; for (int i = 0 ; i < ncols_ ; i++) if (hincol_[i] == 0) empty++ ; return (empty) ; } //@} /*! \name Message handling */ //@{ /// Return message handler inline CoinMessageHandler *messageHandler() const { return handler_; } /*! \brief Set message handler The client retains responsibility for the handler --- it will not be destroyed with the \c CoinPrePostsolveMatrix object. */ inline void setMessageHandler(CoinMessageHandler *handler) { if (defaultHandler_ == true) { delete handler_ ; defaultHandler_ = false ; } handler_ = handler ; } /// Return messages inline CoinMessages messages() const { return messages_; } //@} /*! \name Current and Allocated Size During pre- and postsolve, the matrix will change in size. During presolve it will shrink; during postsolve it will grow. Hence there are two sets of size variables, one for the current size and one for the allocated size. (See the general comments for the CoinPrePostsolveMatrix class for more information.) */ //@{ /// current number of columns int ncols_; /// current number of rows int nrows_; /// current number of coefficients CoinBigIndex nelems_; /// Allocated number of columns int ncols0_; /// Allocated number of rows int nrows0_ ; /// Allocated number of coefficients CoinBigIndex nelems0_ ; /*! \brief Allocated size of bulk storage for row indices and coefficients This is the space allocated for hrow_ and colels_. This must be large enough to allow columns to be copied into empty space when they need to be expanded. For efficiency (to minimize the number of times the representation must be compressed) it's recommended that this be at least 2*nelems0_. */ CoinBigIndex bulk0_ ; /// Ratio of bulk0_ to nelems0_; default is 2. double bulkRatio_; //@} /*! \name Problem representation The matrix is the common column-major format: A pair of vectors with positional correspondence to hold coefficients and row indices, and a second pair of vectors giving the starting position and length of each column in the first pair. */ //@{ /// Vector of column start positions in #hrow_, #colels_ CoinBigIndex *mcstrt_; /// Vector of column lengths int *hincol_; /// Row indices (positional correspondence with #colels_) int *hrow_; /// Coefficients (positional correspondence with #hrow_) double *colels_; /// Objective coefficients double *cost_; /// Original objective offset double originalOffset_; /// Column (primal variable) lower bounds double *clo_; /// Column (primal variable) upper bounds double *cup_; /// Row (constraint) lower bounds double *rlo_; /// Row (constraint) upper bounds double *rup_; /// Original column numbers int * originalColumn_; /// Original row numbers int * originalRow_; /// Primal feasibility tolerance double ztolzb_; /// Dual feasibility tolerance double ztoldj_; /// Maximization/minimization double maxmin_; //@} /*! \name Problem solution information The presolve phase will work without any solution information (appropriate for initial optimisation) or with solution information (appropriate for reoptimisation). When solution information is supplied, presolve will maintain it to the best of its ability. #colstat_ is checked to determine the presence/absence of status information. #sol_ is checked for primal solution information, and #rowduals_ for dual solution information. The postsolve phase requires the complete solution information from the presolved problem (status, primal and dual solutions). It will be transformed into a correct solution for the original problem. */ //@{ /*! \brief Vector of primal variable values If #sol_ exists, it is assumed that primal solution information should be updated and that #acts_ also exists. */ double *sol_; /*! \brief Vector of dual variable values If #rowduals_ exists, it is assumed that dual solution information should be updated and that #rcosts_ also exists. */ double *rowduals_; /*! \brief Vector of constraint left-hand-side values (row activity) Produced by evaluating constraints according to #sol_. Updated iff #sol_ exists. */ double *acts_; /*! \brief Vector of reduced costs Produced by evaluating dual constraints according to #rowduals_. Updated iff #rowduals_ exists. */ double *rcosts_; /*! \brief Status of primal variables Coded with CoinPrePostSolveMatrix::Status, one code per char. colstat_ and #rowstat_ MUST be allocated as a single vector. This is to maintain compatibility with ClpPresolve and OsiPresolve, which do it this way. */ unsigned char *colstat_; /*! \brief Status of constraints More accurately, the status of the logical variable associated with the constraint. Coded with CoinPrePostSolveMatrix::Status, one code per char. Note that this must be allocated as a single vector with #colstat_. */ unsigned char *rowstat_; //@} /*! \name Message handling Uses the standard COIN approach: a default handler is installed, and the CoinPrePostsolveMatrix object takes responsibility for it. If the client replaces the handler with one of their own, it becomes their responsibility. */ //@{ /// Message handler CoinMessageHandler *handler_; /// Indicates if the current #handler_ is default (true) or not (false). bool defaultHandler_; /// Standard COIN messages CoinMessage messages_; //@} }; /*! \class presolvehlink \brief Links to aid in packed matrix modification Currently, the matrices held by the CoinPrePostsolveMatrix and CoinPresolveMatrix objects are represented in the same way as a CoinPackedMatrix. In the course of presolve and postsolve transforms, it will happen that a major-dimension vector needs to increase in size. In order to check whether there is enough room to add another coefficient in place, it helps to know the next vector (in memory order) in the bulk storage area. To do that, a linked list of major-dimension vectors is maintained; the "pre" and "suc" fields give the previous and next vector, in memory order (that is, the vector whose mcstrt_ or mrstrt_ entry is next smaller or larger). Consider a column-major matrix with ncols columns. By definition, presolvehlink[ncols].pre points to the column in the last occupied position of the bulk storage arrays. There is no easy way to find the column which occupies the first position (there is no presolvehlink[-1] to consult). If the column that initially occupies the first position is moved for expansion, there is no way to reclaim the space until the bulk storage is compacted. The same holds for the last and first rows of a row-major matrix, of course. */ class presolvehlink { public: int pre, suc; } ; #define NO_LINK -66666666 /*! \relates presolvehlink \brief unlink vector i Remove vector i from the ordering. */ inline void PRESOLVE_REMOVE_LINK(presolvehlink *link, int i) { int ipre = link[i].pre; int isuc = link[i].suc; if (ipre >= 0) { link[ipre].suc = isuc; } if (isuc >= 0) { link[isuc].pre = ipre; } link[i].pre = NO_LINK, link[i].suc = NO_LINK; } /*! \relates presolvehlink \brief insert vector i after vector j Insert vector i between j and j.suc. */ inline void PRESOLVE_INSERT_LINK(presolvehlink *link, int i, int j) { int isuc = link[j].suc; link[j].suc = i; link[i].pre = j; if (isuc >= 0) { link[isuc].pre = i; } link[i].suc = isuc; } /*! \relates presolvehlink \brief relink vector j in place of vector i Replace vector i in the ordering with vector j. This is equivalent to
     int pre = link[i].pre;
     PRESOLVE_REMOVE_LINK(link,i);
     PRESOLVE_INSERT_LINK(link,j,pre);
   
But, this routine will work even if i happens to be first in the order. */ inline void PRESOLVE_MOVE_LINK(presolvehlink *link, int i, int j) { int ipre = link[i].pre; int isuc = link[i].suc; if (ipre >= 0) { link[ipre].suc = j; } if (isuc >= 0) { link[isuc].pre = j; } link[i].pre = NO_LINK, link[i].suc = NO_LINK; } /*! \class CoinPresolveMatrix \brief Augments CoinPrePostsolveMatrix with information about the problem that is only needed during presolve. For problem manipulation, this class adds a row-major matrix representation, linked lists that allow for easy manipulation of the matrix when applying presolve transforms, and vectors to track row and column processing status (changed, needs further processing, change prohibited) For problem representation, this class adds information about variable type (integer or continuous), an objective offset, and a feasibility tolerance. NOTE that the #anyInteger_ and #anyProhibited_ flags are independent of the vectors used to track this information for individual variables (#integerType_ and #rowChanged_ and #colChanged_, respectively). NOTE also that at the end of presolve the column-major and row-major matrix representations are loosely packed (i.e., there may be gaps between columns in the bulk storage arrays). */ class CoinPresolveMatrix : public CoinPrePostsolveMatrix { public: /*! \brief `Native' constructor This constructor creates an empty object which must then be loaded. On the other hand, it doesn't assume that the client is an OsiSolverInterface. */ CoinPresolveMatrix(int ncols_alloc, int nrows_alloc, CoinBigIndex nelems_alloc) ; /*! \brief Clp OSI constructor See Clp code for the definition. */ CoinPresolveMatrix(int ncols0, double maxmin, // end prepost members ClpSimplex * si, // rowrep int nrows, CoinBigIndex nelems, bool doStatus, double nonLinearVariable, double bulkRatio); /*! \brief Update the model held by a Clp OSI */ void update_model(ClpSimplex * si, int nrows0, int ncols0, CoinBigIndex nelems0); /*! \brief Generic OSI constructor See OSI code for the definition. */ CoinPresolveMatrix(int ncols0, double maxmin, // end prepost members OsiSolverInterface * si, // rowrep int nrows, CoinBigIndex nelems, bool doStatus, double nonLinearVariable, const char * prohibited); /*! \brief Update the model held by a generic OSI */ void update_model(OsiSolverInterface * si, int nrows0, int ncols0, CoinBigIndex nelems0); /// Destructor ~CoinPresolveMatrix(); /*! \brief Initialize a CoinPostsolveMatrix object, destroying the CoinPresolveMatrix object. See CoinPostsolveMatrix::assignPresolveToPostsolve. */ friend void assignPresolveToPostsolve (CoinPresolveMatrix *&preObj) ; /*! \name Functions to load the problem representation */ //@{ /*! \brief Load the cofficient matrix. Load the coefficient matrix before loading the other vectors (bounds, objective, variable type) required to define the problem. */ void setMatrix(const CoinPackedMatrix *mtx) ; /// Count number of empty rows inline int countEmptyRows() { int empty = 0 ; for (int i = 0 ; i < nrows_ ; i++) if (hinrow_[i] == 0) empty++ ; return (empty) ; } /*! \brief Set variable type information for a single variable Set \p variableType to 0 for continous, 1 for integer. Does not manipulate the #anyInteger_ flag. */ inline void setVariableType(int i, int variableType) { if (integerType_ == 0) integerType_ = new unsigned char [ncols0_] ; integerType_[i] = static_cast(variableType) ; } /*! \brief Set variable type information for all variables Set \p variableType[i] to 0 for continuous, 1 for integer. Does not manipulate the #anyInteger_ flag. */ void setVariableType(const unsigned char *variableType, int lenParam) ; /*! \brief Set the type of all variables allIntegers should be true to set the type to integer, false to set the type to continuous. */ void setVariableType (bool allIntegers, int lenParam) ; /// Set a flag for presence (true) or absence (false) of integer variables inline void setAnyInteger (bool anyInteger = true) { anyInteger_ = anyInteger ; } //@} /*! \name Functions to retrieve problem information */ //@{ /// Get row start vector for row-major packed matrix inline const CoinBigIndex *getRowStarts() const { return (mrstrt_) ; } /// Get vector of column indices for row-major packed matrix inline const int *getColIndicesByRow() const { return (hcol_) ; } /// Get vector of elements for row-major packed matrix inline const double *getElementsByRow() const { return (rowels_) ; } /*! \brief Check for integrality of the specified variable. Consults the #integerType_ vector if present; fallback is the #anyInteger_ flag. */ inline bool isInteger (int i) const { if (integerType_ == 0) { return (anyInteger_) ; } else if (integerType_[i] == 1) { return (true) ; } else { return (false) ; } } /*! \brief Check if there are any integer variables Consults the #anyInteger_ flag */ inline bool anyInteger () const { return (anyInteger_) ; } /// Picks up any special options inline int presolveOptions() const { return presolveOptions_;} /// Sets any special options inline void setPresolveOptions(int value) { presolveOptions_=value;} //@} /*! \name Matrix storage management links Linked lists, modelled after the linked lists used in OSL factorization. They are used for management of the bulk coefficient and minor index storage areas. */ //@{ /// Linked list for the column-major representation. presolvehlink *clink_; /// Linked list for the row-major representation. presolvehlink *rlink_; //@} /// Objective function offset introduced during presolve double dobias_; /// Adjust objective function constant offset inline void change_bias(double change_amount) { dobias_ += change_amount; #if PRESOLVE_DEBUG assert(fabs(change_amount)<1.0e50); #endif if (change_amount) PRESOLVE_STMT(printf("changing bias by %g to %g\n", change_amount, dobias_)); } /*! \name Row-major representation Common row-major format: A pair of vectors with positional correspondence to hold coefficients and column indices, and a second pair of vectors giving the starting position and length of each row in the first pair. */ //@{ /// Vector of row start positions in #hcol, #rowels_ CoinBigIndex *mrstrt_; /// Vector of row lengths int *hinrow_; /// Coefficients (positional correspondence with #hcol_) double *rowels_; /// Column indices (positional correspondence with #rowels_) int *hcol_; //@} /// Tracks integrality of columns (1 for integer, 0 for continuous) unsigned char *integerType_; /*! \brief Flag to say if any variables are integer Note that this flag is not manipulated by the various \c setVariableType routines. */ bool anyInteger_ ; /// Print statistics for tuning bool tuning_; /// Say we want statistics - also set time void statistics(); /// Start time of presolve double startTime_; /// Bounds can be moved by this to retain feasibility double feasibilityTolerance_; /// Return feasibility tolerance inline double feasibilityTolerance() { return (feasibilityTolerance_) ; } /// Set feasibility tolerance inline void setFeasibilityTolerance (double val) { feasibilityTolerance_ = val ; } /*! \brief Output status: 0 = feasible, 1 = infeasible, 2 = unbounded Actually implemented as single bit flags: 1^0 = infeasible, 1^1 = unbounded. */ int status_; /// Returns problem status (0 = feasible, 1 = infeasible, 2 = unbounded) inline int status() { return (status_) ; } /// Set problem status inline void setStatus(int status) { status_ = (status&0x3) ; } /*! \brief Pass number Used to control the execution of testRedundant (evoked by the implied_free transform). */ int pass_; /// Set pass number inline void setPass (int pass = 0) { pass_ = pass ; } /*! \brief Maximum substitution level Used to control the execution of subst from implied_free */ int maxSubstLevel_; /// Set Maximum substitution level (normally 3) inline void setMaximumSubstitutionLevel (int level) { maxSubstLevel_ = level ; } /*! \name Row and column processing status Information used to determine if rows or columns can be changed and if they require further processing due to changes. There are four major lists: the [row,col]ToDo list, and the [row,col]NextToDo list. In general, a transform processes entries from the ToDo list and adds entries to the NextToDo list. There are two vectors, [row,col]Changed, which track the status of individual rows and columns. */ //@{ /*! \brief Column change status information Coded using the following bits:
  • 0x01: Column has changed
  • 0x02: preprocessing prohibited
  • 0x04: Column has been used
*/ unsigned char * colChanged_; /// Input list of columns to process int * colsToDo_; /// Length of #colsToDo_ int numberColsToDo_; /// Output list of columns to process next int * nextColsToDo_; /// Length of #nextColsToDo_ int numberNextColsToDo_; /*! \brief Row change status information Coded using the following bits:
  • 0x01: Row has changed
  • 0x02: preprocessing prohibited
  • 0x04: Row has been used
*/ unsigned char * rowChanged_; /// Input list of rows to process int * rowsToDo_; /// Length of #rowsToDo_ int numberRowsToDo_; /// Output list of rows to process next int * nextRowsToDo_; /// Length of #nextRowsToDo_ int numberNextRowsToDo_; /** Presolve options 1 set if allow duplicate column tests for integer variables 2 set to allow code to try and fix infeasibilities 4 set to inhibit x+y+z=1 mods */ int presolveOptions_; /*! Flag to say if any rows or columns are marked as prohibited Note that this flag is not manipulated by any of the various \c set*Prohibited routines. */ bool anyProhibited_; //@} /*! \name Functions to manipulate row and column processing status */ //@{ /*! \brief Initialise the column ToDo lists Places all columns in the #colsToDo_ list except for columns marked as prohibited (viz. #colChanged_). */ void initColsToDo () ; /*! \brief Step column ToDo lists Moves columns on the #nextColsToDo_ list to the #colsToDo_ list, emptying #nextColsToDo_. Returns the number of columns transferred. */ int stepColsToDo () ; /// Return the number of columns on the #colsToDo_ list inline int numberColsToDo() { return (numberColsToDo_) ; } /// Has column been changed? inline bool colChanged(int i) const { return (colChanged_[i]&1)!=0; } /// Mark column as not changed inline void unsetColChanged(int i) { colChanged_[i] = static_cast(colChanged_[i] & (~1)) ; } /// Mark column as changed. inline void setColChanged(int i) { colChanged_[i] = static_cast(colChanged_[i] | (1)) ; } /// Mark column as changed and add to list of columns to process next inline void addCol(int i) { if ((colChanged_[i]&1)==0) { colChanged_[i] = static_cast(colChanged_[i] | (1)) ; nextColsToDo_[numberNextColsToDo_++] = i; } } /// Test if column is eligible for preprocessing inline bool colProhibited(int i) const { return (colChanged_[i]&2)!=0; } /*! \brief Test if column is eligible for preprocessing The difference between this method and #colProhibited() is that this method first tests #anyProhibited_ before examining the specific entry for the specified column. */ inline bool colProhibited2(int i) const { if (!anyProhibited_) return false; else return (colChanged_[i]&2)!=0; } /// Mark column as ineligible for preprocessing inline void setColProhibited(int i) { colChanged_[i] = static_cast(colChanged_[i] | (2)) ; } /*! \brief Test if column is marked as used This is for doing faster lookups to see where two columns have entries in common. */ inline bool colUsed(int i) const { return (colChanged_[i]&4)!=0; } /// Mark column as used inline void setColUsed(int i) { colChanged_[i] = static_cast(colChanged_[i] | (4)) ; } /// Mark column as unused inline void unsetColUsed(int i) { colChanged_[i] = static_cast(colChanged_[i] & (~4)) ; } /*! \brief Initialise the row ToDo lists Places all rows in the #rowsToDo_ list except for rows marked as prohibited (viz. #rowChanged_). */ void initRowsToDo () ; /*! \brief Step row ToDo lists Moves rows on the #nextRowsToDo_ list to the #rowsToDo_ list, emptying #nextRowsToDo_. Returns the number of rows transferred. */ int stepRowsToDo () ; /// Return the number of rows on the #rowsToDo_ list inline int numberRowsToDo() { return (numberRowsToDo_) ; } /// Has row been changed? inline bool rowChanged(int i) const { return (rowChanged_[i]&1)!=0; } /// Mark row as not changed inline void unsetRowChanged(int i) { rowChanged_[i] = static_cast(rowChanged_[i] & (~1)) ; } /// Mark row as changed inline void setRowChanged(int i) { rowChanged_[i] = static_cast(rowChanged_[i] | (1)) ; } /// Mark row as changed and add to list of rows to process next inline void addRow(int i) { if ((rowChanged_[i]&1)==0) { rowChanged_[i] = static_cast(rowChanged_[i] | (1)) ; nextRowsToDo_[numberNextRowsToDo_++] = i; } } /// Test if row is eligible for preprocessing inline bool rowProhibited(int i) const { return (rowChanged_[i]&2)!=0; } /*! \brief Test if row is eligible for preprocessing The difference between this method and #rowProhibited() is that this method first tests #anyProhibited_ before examining the specific entry for the specified row. */ inline bool rowProhibited2(int i) const { if (!anyProhibited_) return false; else return (rowChanged_[i]&2)!=0; } /// Mark row as ineligible for preprocessing inline void setRowProhibited(int i) { rowChanged_[i] = static_cast(rowChanged_[i] | (2)) ; } /*! \brief Test if row is marked as used This is for doing faster lookups to see where two rows have entries in common. It can be used anywhere as long as it ends up zeroed out. */ inline bool rowUsed(int i) const { return (rowChanged_[i]&4)!=0; } /// Mark row as used inline void setRowUsed(int i) { rowChanged_[i] = static_cast(rowChanged_[i] | (4)) ; } /// Mark row as unused inline void unsetRowUsed(int i) { rowChanged_[i] = static_cast(rowChanged_[i] & (~4)) ; } /// Check if there are any prohibited rows or columns inline bool anyProhibited() const { return anyProhibited_;} /// Set a flag for presence of prohibited rows or columns inline void setAnyProhibited(bool val = true) { anyProhibited_ = val ; } //@} }; /*! \class CoinPostsolveMatrix \brief Augments CoinPrePostsolveMatrix with information about the problem that is only needed during postsolve. The notable point is that the matrix representation is threaded. The representation is column-major and starts with the standard two pairs of arrays: one pair to hold the row indices and coefficients, the second pair to hold the column starting positions and lengths. But the row indices and coefficients for a column do not necessarily occupy a contiguous block in their respective arrays. Instead, a link array gives the position of the next (row index,coefficient) pair. If the row index and value of a coefficient a occupy position kp in their arrays, then the position of the next coefficient a is found as kq = link[kp]. This threaded representation allows for efficient expansion of columns as rows are reintroduced during postsolve transformations. The basic packed structures are allocated to the expected size of the postsolved matrix, and as new coefficients are added, their location is simply added to the thread for the column. There is no provision to convert the threaded representation to a packed representation. In the context of postsolve, it's not required. (You did keep a copy of the original matrix, eh?) */ class CoinPostsolveMatrix : public CoinPrePostsolveMatrix { public: /*! \brief `Native' constructor This constructor creates an empty object which must then be loaded. On the other hand, it doesn't assume that the client is an OsiSolverInterface. */ CoinPostsolveMatrix(int ncols_alloc, int nrows_alloc, CoinBigIndex nelems_alloc) ; /*! \brief Clp OSI constructor See Clp code for the definition. */ CoinPostsolveMatrix(ClpSimplex * si, int ncols0, int nrows0, CoinBigIndex nelems0, double maxmin_, // end prepost members double *sol, double *acts, unsigned char *colstat, unsigned char *rowstat); /*! \brief Generic OSI constructor See OSI code for the definition. */ CoinPostsolveMatrix(OsiSolverInterface * si, int ncols0, int nrows0, CoinBigIndex nelems0, double maxmin_, // end prepost members double *sol, double *acts, unsigned char *colstat, unsigned char *rowstat); /*! \brief Load an empty CoinPostsolveMatrix from a CoinPresolveMatrix This routine transfers the contents of the CoinPrePostsolveMatrix object from the CoinPresolveMatrix object to the CoinPostsolveMatrix object and completes initialisation of the CoinPostsolveMatrix object. The empty shell of the CoinPresolveMatrix object is destroyed. The routine expects an empty CoinPostsolveMatrix object. If handed a loaded object, a lot of memory will leak. */ void assignPresolveToPostsolve (CoinPresolveMatrix *&preObj) ; /// Destructor ~CoinPostsolveMatrix(); /*! \name Column thread structures As mentioned in the class documentation, the entries for a given column do not necessarily occupy a contiguous block of space. The #link_ array is used to maintain the threading. There is one thread for each column, and a single thread for all free entries in #hrow_ and #colels_. The allocated size of #link_ must be at least as large as the allocated size of #hrow_ and #colels_. */ //@{ /*! \brief First entry in free entries thread */ CoinBigIndex free_list_; /// Allocated size of #link_ int maxlink_; /*! \brief Thread array Within a thread, link_[k] points to the next entry in the thread. */ CoinBigIndex *link_; //@} /*! \name Debugging aids These arrays are allocated only when CoinPresolve is compiled with PRESOLVE_DEBUG defined. They hold codes which track the reason that a column or row is added to the problem during postsolve. */ //@{ char *cdone_; char *rdone_; //@} /// debug void check_nbasic(); }; #define PRESOLVEFINITE(n) (-PRESOLVE_INF < (n) && (n) < PRESOLVE_INF) /*! \defgroup MtxManip Presolve Matrix Manipulation Functions Functions to work with the loosely packed and threaded packed matrix structures used during presolve and postsolve. */ //@{ /*! \relates CoinPrePostsolveMatrix \brief Initialise linked list for major vector order in bulk storage */ void presolve_make_memlists(CoinBigIndex *starts, int *lengths, presolvehlink *link, int n); /*! \relates CoinPrePostsolveMatrix \brief Make sure a major-dimension vector k has room for one more coefficient. You can use this directly, or use the inline wrappers presolve_expand_col and presolve_expand_row */ bool presolve_expand_major(CoinBigIndex *majstrts, double *majels, int *minndxs, int *majlens, presolvehlink *majlinks, int nmaj, int k) ; /*! \relates CoinPrePostsolveMatrix \brief Make sure a column (colx) in a column-major matrix has room for one more coefficient */ inline bool presolve_expand_col(CoinBigIndex *mcstrt, double *colels, int *hrow, int *hincol, presolvehlink *clink, int ncols, int colx) { return presolve_expand_major(mcstrt,colels, hrow,hincol,clink,ncols,colx) ; } /*! \relates CoinPrePostsolveMatrix \brief Make sure a row (rowx) in a row-major matrix has room for one more coefficient */ inline bool presolve_expand_row(CoinBigIndex *mrstrt, double *rowels, int *hcol, int *hinrow, presolvehlink *rlink, int nrows, int rowx) { return presolve_expand_major(mrstrt,rowels, hcol,hinrow,rlink,nrows,rowx) ; } /*! \relates CoinPrePostsolveMatrix \brief Find position of a minor index in a major vector. The routine returns the position \c k in \p minndxs for the specified minor index \p tgt. It will abort if the entry does not exist. Can be used directly or via the inline wrappers presolve_find_row and presolve_find_col. */ CoinBigIndex presolve_find_minor(int tgt, CoinBigIndex ks, CoinBigIndex ke, const int *minndxs); /*! \relates CoinPrePostsolveMatrix \brief Find position of a row in a column in a column-major matrix. The routine returns the position \c k in \p hrow for the specified \p row. It will abort if the entry does not exist. */ inline CoinBigIndex presolve_find_row(int row, CoinBigIndex kcs, CoinBigIndex kce, const int *hrow) { return presolve_find_minor(row,kcs,kce,hrow) ; } /*! \relates CoinPostsolveMatrix \brief Find position of a column in a row in a row-major matrix. The routine returns the position \c k in \p hcol for the specified \p col. It will abort if the entry does not exist. */ inline CoinBigIndex presolve_find_col(int col, CoinBigIndex krs, CoinBigIndex kre, const int *hcol) { return presolve_find_minor(col,krs,kre,hcol) ; } /*! \relates CoinPrePostsolveMatrix \brief Find position of a minor index in a major vector. The routine returns the position \c k in \p minndxs for the specified minor index \p tgt. A return value of \p ke means the entry does not exist. Can be used directly or via the inline wrappers presolve_find_row1 and presolve_find_col1. */ CoinBigIndex presolve_find_minor1(int tgt, CoinBigIndex ks, CoinBigIndex ke, const int *minndxs); /*! \relates CoinPrePostsolveMatrix \brief Find position of a row in a column in a column-major matrix. The routine returns the position \c k in \p hrow for the specified \p row. A return value of \p kce means the entry does not exist. */ inline CoinBigIndex presolve_find_row1(int row, CoinBigIndex kcs, CoinBigIndex kce, const int *hrow) { return presolve_find_minor1(row,kcs,kce,hrow) ; } /*! \relates CoinPrePostsolveMatrix \brief Find position of a column in a row in a row-major matrix. The routine returns the position \c k in \p hcol for the specified \p col. A return value of \p kre means the entry does not exist. */ inline CoinBigIndex presolve_find_col1(int col, CoinBigIndex krs, CoinBigIndex kre, const int *hcol) { return presolve_find_minor1(col,krs,kre,hcol) ; } /*! \relates CoinPostsolveMatrix \brief Find position of a minor index in a major vector in a threaded matrix. The routine returns the position \c k in \p minndxs for the specified minor index \p tgt. It will abort if the entry does not exist. Can be used directly or via the inline wrapper presolve_find_row2. */ CoinBigIndex presolve_find_minor2(int tgt, CoinBigIndex ks, int majlen, const int *minndxs, const CoinBigIndex *majlinks) ; /*! \relates CoinPostsolveMatrix \brief Find position of a row in a column in a column-major threaded matrix. The routine returns the position \c k in \p hrow for the specified \p row. It will abort if the entry does not exist. */ inline CoinBigIndex presolve_find_row2(int row, CoinBigIndex kcs, int collen, const int *hrow, const CoinBigIndex *clinks) { return presolve_find_minor2(row,kcs,collen,hrow,clinks) ; } /*! \relates CoinPostsolveMatrix \brief Find position of a minor index in a major vector in a threaded matrix. The routine returns the position \c k in \p minndxs for the specified minor index \p tgt. It will return -1 if the entry does not exist. Can be used directly or via the inline wrappers presolve_find_row3. */ CoinBigIndex presolve_find_minor3(int tgt, CoinBigIndex ks, int majlen, const int *minndxs, const CoinBigIndex *majlinks) ; /*! \relates CoinPostsolveMatrix \brief Find position of a row in a column in a column-major threaded matrix. The routine returns the position \c k in \p hrow for the specified \p row. It will return -1 if the entry does not exist. */ inline CoinBigIndex presolve_find_row3(int row, CoinBigIndex kcs, int collen, const int *hrow, const CoinBigIndex *clinks) { return presolve_find_minor3(row,kcs,collen,hrow,clinks) ; } /*! \relates CoinPrePostsolveMatrix \brief Delete the entry for a minor index from a major vector. Deletes the entry for \p minndx from the major vector \p majndx. Specifically, the relevant entries are removed from the minor index (\p minndxs) and coefficient (\p els) arrays and the vector length (\p majlens) is decremented. Loose packing is maintained by swapping the last entry in the row into the position occupied by the deleted entry. */ void presolve_delete_from_major(int majndx, int minndx, const CoinBigIndex *majstrts, int *majlens, int *minndxs, double *els) ; // Delete all marked from major (and zero marked) void presolve_delete_many_from_major(int majndx, char * marked, const CoinBigIndex *majstrts, int *majlens, int *minndxs, double *els) ; /*! \relates CoinPrePostsolveMatrix \brief Delete the entry for row \p row from column \p col in a column-major matrix Deletes the entry for \p row from the major vector for \p col. Specifically, the relevant entries are removed from the row index (\p hrow) and coefficient (\p colels) arrays and the vector length (\p hincol) is decremented. Loose packing is maintained by swapping the last entry in the row into the position occupied by the deleted entry. */ inline void presolve_delete_from_col(int row, int col, const CoinBigIndex *mcstrt, int *hincol, int *hrow, double *colels) { presolve_delete_from_major(col,row,mcstrt,hincol,hrow,colels) ; } /*! \relates CoinPrePostsolveMatrix \brief Delete the entry for column \p col from row \p row in a row-major matrix Deletes the entry for \p col from the major vector for \p row. Specifically, the relevant entries are removed from the column index (\p hcol) and coefficient (\p rowels) arrays and the vector length (\p hinrow) is decremented. Loose packing is maintained by swapping the last entry in the column into the position occupied by the deleted entry. */ inline void presolve_delete_from_row(int row, int col, const CoinBigIndex *mrstrt, int *hinrow, int *hcol, double *rowels) { presolve_delete_from_major(row,col,mrstrt,hinrow,hcol,rowels) ; } /*! \relates CoinPostsolveMatrix \brief Delete the entry for a minor index from a major vector in a threaded matrix. Deletes the entry for \p minndx from the major vector \p majndx. Specifically, the relevant entries are removed from the minor index (\p minndxs) and coefficient (\p els) arrays and the vector length (\p majlens) is decremented. The thread for the major vector is relinked around the deleted entry and the space is returned to the free list. */ void presolve_delete_from_major2 (int majndx, int minndx, CoinBigIndex *majstrts, int *majlens, int *minndxs, double *els, int *majlinks, CoinBigIndex *free_listp) ; /*! \relates CoinPostsolveMatrix \brief Delete the entry for row \p row from column \p col in a column-major threaded matrix Deletes the entry for \p row from the major vector for \p col. Specifically, the relevant entries are removed from the row index (\p hrow) and coefficient (\p colels) arrays and the vector length (\p hincol) is decremented. The thread for the major vector is relinked around the deleted entry and the space is returned to the free list. */ inline void presolve_delete_from_col2(int row, int col, CoinBigIndex *mcstrt, int *hincol, int *hrow, double *colels, int *clinks, CoinBigIndex *free_listp) { presolve_delete_from_major2(col,row,mcstrt,hincol,hrow,colels,clinks, free_listp) ; } //@} /*! \defgroup PresolveUtilities Presolve Utility Functions Utilities used by multiple presolve transform objects. */ //@{ /*! \brief Duplicate a major-dimension vector; optionally omit the entry with minor index \p tgt. Designed to copy a major-dimension vector from the paired coefficient (\p elems) and minor index (\p indices) arrays used in the standard packed matrix representation. Copies \p length entries starting at \p offset. If \p tgt is specified, the entry with minor index == \p tgt is omitted from the copy. */ double *presolve_dupmajor(const double *elems, const int *indices, int length, CoinBigIndex offset, int tgt = -1); //@} #endif