ClpSimplex.hpp

// Copyright (C) 2002, International Business Machines
// Corporation and others.  All Rights Reserved.

/* 
   Authors
 
   John Forrest

 */
#ifndef ClpSimplex_H
#define ClpSimplex_H

#include <iostream>
#include <cfloat>
#include "ClpModel.hpp"
#include "ClpMatrixBase.hpp"
#include "ClpSolve.hpp"
class ClpDualRowPivot;
class ClpPrimalColumnPivot;
class ClpFactorization;
class CoinIndexedVector;
class ClpNonLinearCost;
class ClpSimplexProgress;

class ClpSimplex : public ClpModel {
   friend void ClpSimplexUnitTest(const std::string & mpsDir,
                                  const std::string & netlibDir);

public:
  enum Status {
    isFree = 0x00,
    basic = 0x01,
    atUpperBound = 0x02,
    atLowerBound = 0x03,
    superBasic = 0x04,
    isFixed = 0x05
  };
  // For Dual
  enum FakeBound {
    noFake = 0x00,
    bothFake = 0x01,
    upperFake = 0x02,
    lowerFake = 0x03
  };


    ClpSimplex (  );

  ClpSimplex(const ClpSimplex &);
  ClpSimplex(const ClpModel &);
  ClpSimplex (const ClpModel * wholeModel,
              int numberRows, const int * whichRows,
              int numberColumns, const int * whichColumns,
              bool dropNames=true, bool dropIntegers=true);
  ClpSimplex (ClpSimplex * wholeModel,
              int numberColumns, const int * whichColumns);
  void originalModel(ClpSimplex * miniModel);
    ClpSimplex & operator=(const ClpSimplex & rhs);
   ~ClpSimplex (  );
  // Ones below are just ClpModel with some changes
  void loadProblem (  const ClpMatrixBase& matrix,
                     const double* collb, const double* colub,   
                     const double* obj,
                     const double* rowlb, const double* rowub,
                      const double * rowObjective=NULL);
  void loadProblem (  const CoinPackedMatrix& matrix,
                     const double* collb, const double* colub,   
                     const double* obj,
                     const double* rowlb, const double* rowub,
                      const double * rowObjective=NULL);

  void loadProblem (  const int numcols, const int numrows,
                     const CoinBigIndex* start, const int* index,
                     const double* value,
                     const double* collb, const double* colub,   
                     const double* obj,
                      const double* rowlb, const double* rowub,
                      const double * rowObjective=NULL);
  void loadProblem (  const int numcols, const int numrows,
                     const CoinBigIndex* start, const int* index,
                      const double* value,const int * length,
                     const double* collb, const double* colub,   
                     const double* obj,
                      const double* rowlb, const double* rowub,
                      const double * rowObjective=NULL);
  int readMps(const char *filename,
              bool keepNames=false,
              bool ignoreErrors = false);
  void borrowModel(ClpModel & otherModel);

  int initialSolve(ClpSolve & options);
  int initialSolve();
  int initialDualSolve();
  int initialPrimalSolve();
  int dual(int ifValuesPass=0);
  int primal(int ifValuesPass=0);
  int quadraticSLP(int numberPasses,double deltaTolerance);
  int quadraticPrimal(int phase=2);
  void setFactorization( ClpFactorization & factorization);
  void scaling(int mode=1);
  inline int scalingFlag() const {return scalingFlag_;};
  int tightenPrimalBounds(double factor=0.0);
  int crash(double gap,int pivot);
  void setDualRowPivotAlgorithm(ClpDualRowPivot & choice);
  void setPrimalColumnPivotAlgorithm(ClpPrimalColumnPivot & choice);
  int strongBranching(int numberVariables,const int * variables,
                      double * newLower, double * newUpper,
                      double ** outputSolution,
                      int * outputStatus, int * outputIterations,
                      bool stopOnFirstInfeasible=true,
                      bool alwaysFinish=false);

  int pivot();

  int primalPivotResult();
  
  int dualPivotResult();

  int startup(int ifValuesPass);
  void finish();
  
  bool statusOfProblem();


  inline bool primalFeasible() const
         { return (numberPrimalInfeasibilities_==0);};
  inline bool dualFeasible() const
         { return (numberDualInfeasibilities_==0);};
  inline ClpFactorization * factorization() const 
          { return factorization_;};
  bool sparseFactorization() const;
  void setSparseFactorization(bool value);
  int factorizationFrequency() const;
  void setFactorizationFrequency(int value);
  inline double dualBound() const
          { return dualBound_;};
  void setDualBound(double value);
  inline double infeasibilityCost() const
          { return infeasibilityCost_;};
  void setInfeasibilityCost(double value);
  inline int perturbation() const
    { return perturbation_;};
  void setPerturbation(int value);
  inline int algorithm() const 
  {return algorithm_; } ;
  inline void setAlgorithm(int value)
  {algorithm_=value; } ;
  inline double sumDualInfeasibilities() const 
          { return sumDualInfeasibilities_;} ;
  inline int numberDualInfeasibilities() const 
          { return numberDualInfeasibilities_;} ;
  inline double sumPrimalInfeasibilities() const 
          { return sumPrimalInfeasibilities_;} ;
  inline void setSumPrimalInfeasibilities(double value)
          { sumPrimalInfeasibilities_=value;} ;
  inline int numberPrimalInfeasibilities() const 
          { return numberPrimalInfeasibilities_;} ;
  int saveModel(const char * fileName);
  int restoreModel(const char * fileName);
  
  void checkSolution();
  inline CoinIndexedVector * rowArray(int index) const
  { return rowArray_[index];};
  inline CoinIndexedVector * columnArray(int index) const
  { return columnArray_[index];};

  /******************** End of most useful part **************/
  int getSolution (  const double * rowActivities,
                     const double * columnActivities);
  int getSolution ();
  int createPiecewiseLinearCosts(const int * starts,
                   const double * lower, const double * gradient);
  void returnModel(ClpSimplex & otherModel);

    ClpDataSave saveData() ;
    void restoreData(ClpDataSave saved);
  int internalFactorize(int solveType);
  void cleanStatus();
  int factorize();
  void computeDuals(double * givenDjs);
  void computePrimals (  const double * rowActivities,
                     const double * columnActivities);
  void unpack(CoinIndexedVector * rowArray) const ;
  void unpack(CoinIndexedVector * rowArray,int sequence) const;
  void unpackPacked(CoinIndexedVector * rowArray) ;
  void unpackPacked(CoinIndexedVector * rowArray,int sequence);
  
  int housekeeping(double objectiveChange);
  void checkPrimalSolution(const double * rowActivities=NULL,
                           const double * columnActivies=NULL);
  void checkDualSolution();
  void setValuesPassAction(float incomingInfeasibility,
                           float allowedInfeasibility);

   void times(double scalar,
                       const double * x, double * y) const;
    void transposeTimes(double scalar,
                                const double * x, double * y) const ;


  inline double columnPrimalInfeasibility() const 
          { return columnPrimalInfeasibility_;} ;
  inline int columnPrimalSequence() const 
          { return columnPrimalSequence_;} ;
  inline double rowPrimalInfeasibility() const 
          { return rowPrimalInfeasibility_;} ;
  inline int rowPrimalSequence() const 
          { return rowPrimalSequence_;} ;
  inline double columnDualInfeasibility() const 
          { return columnDualInfeasibility_;} ;
  inline int columnDualSequence() const 
          { return columnDualSequence_;} ;
  inline double rowDualInfeasibility() const 
          { return rowDualInfeasibility_;} ;
  inline int rowDualSequence() const 
          { return rowDualSequence_;} ;
  inline double primalToleranceToGetOptimal() const 
          { return primalToleranceToGetOptimal_;} ;
  inline double remainingDualInfeasibility() const 
          { return remainingDualInfeasibility_;} ;
  inline double largeValue() const 
          { return largeValue_;} ;
  void setLargeValue( double value) ;
  inline double largestPrimalError() const
          { return largestPrimalError_;} ;
  inline double largestDualError() const
          { return largestDualError_;} ;
  inline double largestSolutionError() const
          { return largestSolutionError_;} ;
  inline int * pivotVariable() const
          { return pivotVariable_;};
  inline double objectiveScale() const 
          { return objectiveScale_;} ;
  inline void setObjectiveScale(double value)
          { objectiveScale_ = value;} ;
  inline double currentDualTolerance() const 
          { return dualTolerance_;} ;
  inline void setCurrentDualTolerance(double value)
          { dualTolerance_ = value;} ;
  inline double currentPrimalTolerance() const 
          { return primalTolerance_;} ;
  inline void setCurrentPrimalTolerance(double value)
          { primalTolerance_ = value;} ;
  inline int numberRefinements() const 
          { return numberRefinements_;} ;
  void setNumberRefinements( int value) ;
  inline double alpha() const { return alpha_;};
  inline double dualIn() const { return dualIn_;};
  inline int pivotRow() const{ return pivotRow_;};
  double valueIncomingDual() const;

  protected:
  int gutsOfSolution ( double * givenDuals,
                       const double * givenPrimals,
                       bool valuesPass=false);
  void gutsOfDelete(int type);
  void gutsOfCopy(const ClpSimplex & rhs);
  bool createRim(int what,bool makeRowCopy=false);
  void deleteRim(int getRidOfFactorizationData=2);
  bool sanityCheck();
  public:
  inline double * solutionRegion(int section)
  { if (!section) return rowActivityWork_; else return columnActivityWork_;};
  inline double * djRegion(int section)
  { if (!section) return rowReducedCost_; else return reducedCostWork_;};
  inline double * lowerRegion(int section)
  { if (!section) return rowLowerWork_; else return columnLowerWork_;};
  inline double * upperRegion(int section)
  { if (!section) return rowUpperWork_; else return columnUpperWork_;};
  inline double * costRegion(int section)
  { if (!section) return rowObjectiveWork_; else return objectiveWork_;};
  inline double * solutionRegion()
  { return solution_;};
  inline double * djRegion()
  { return dj_;};
  inline double * lowerRegion()
  { return lower_;};
  inline double * upperRegion()
  { return upper_;};
  inline double * costRegion()
  { return cost_;};
  inline Status getStatus(int sequence) const
  {return static_cast<Status> (status_[sequence]&7);};
  inline void setStatus(int sequence, Status status)
  {
    unsigned char & st_byte = status_[sequence];
    st_byte &= ~7;
    st_byte |= status;
  };
  void setInitialDenseFactorization(bool onOff);
  bool  initialDenseFactorization() const;
  inline int sequenceIn() const
  {return sequenceIn_;};
  inline int sequenceOut() const
  {return sequenceOut_;};
  inline void  setSequenceIn(int sequence)
  { sequenceIn_=sequence;};
  inline void  setSequenceOut(int sequence)
  { sequenceOut_=sequence;};
  inline int directionIn() const
  {return directionIn_;};
  inline int directionOut() const
  {return directionOut_;};
  inline void  setDirectionIn(int direction)
  { directionIn_=direction;};
  inline void  setDirectionOut(int direction)
  { directionOut_=direction;};
  inline int isColumn(int sequence) const
  { return sequence<numberColumns_ ? 1 : 0;};
  inline int sequenceWithin(int sequence) const
  { return sequence<numberColumns_ ? sequence : sequence-numberColumns_;};
  inline double solution(int sequence)
  { return solution_[sequence];};
  inline double & solutionAddress(int sequence)
  { return solution_[sequence];};
  inline double reducedCost(int sequence)
   { return dj_[sequence];};
  inline double & reducedCostAddress(int sequence)
   { return dj_[sequence];};
  inline double lower(int sequence)
  { return lower_[sequence];};
  inline double & lowerAddress(int sequence)
  { return lower_[sequence];};
  inline double upper(int sequence)
  { return upper_[sequence];};
  inline double & upperAddress(int sequence)
  { return upper_[sequence];};
  inline double cost(int sequence)
  { return cost_[sequence];};
  inline double & costAddress(int sequence)
  { return cost_[sequence];};
  inline double originalLower(int iSequence) const
  { if (iSequence<numberColumns_) return columnLower_[iSequence]; else
    return rowLower_[iSequence-numberColumns_];};
  inline double originalUpper(int iSequence) const
  { if (iSequence<numberColumns_) return columnUpper_[iSequence]; else
    return rowUpper_[iSequence-numberColumns_];};
  inline double theta() const
  { return theta_;};
  const double * rowScale() const {return rowScale_;};
  const double * columnScale() const {return columnScale_;};
  void setRowScale(double * scale) { rowScale_ = scale;};
  void setColumnScale(double * scale) { columnScale_ = scale;};
  inline ClpNonLinearCost * nonLinearCost() const
  { return nonLinearCost_;};

  inline void setFakeBound(int sequence, FakeBound fakeBound)
  {
    unsigned char & st_byte = status_[sequence];
    st_byte &= ~24;
    st_byte |= fakeBound<<3;
  };
  inline FakeBound getFakeBound(int sequence) const
  {return static_cast<FakeBound> ((status_[sequence]>>3)&3);};
  inline void setRowStatus(int sequence, Status status)
  {
    unsigned char & st_byte = status_[sequence+numberColumns_];
    st_byte &= ~7;
    st_byte |= status;
  };
  inline Status getRowStatus(int sequence) const
  {return static_cast<Status> (status_[sequence+numberColumns_]&7);};
  inline void setColumnStatus(int sequence, Status status)
  {
    unsigned char & st_byte = status_[sequence];
    st_byte &= ~7;
    st_byte |= status;
  };
  inline Status getColumnStatus(int sequence) const
  {return static_cast<Status> (status_[sequence]&7);};
  inline void setPivoted( int sequence)
  { status_[sequence] |= 32;};
  inline void clearPivoted( int sequence)
  { status_[sequence] &= ~32; };
  inline bool pivoted(int sequence) const
  {return (((status_[sequence]>>5)&1)!=0);};
  inline void setFlagged( int sequence)
  {
    status_[sequence] |= 64;
  };
  inline void clearFlagged( int sequence)
  {
    status_[sequence] &= ~64;
  };
  inline bool flagged(int sequence) const
  {return ((status_[sequence]&64)!=0);};
  inline void setActive( int iRow)
  {
    status_[iRow] |= 128;
  };
  inline void clearActive( int iRow)
  {
    status_[iRow] &= ~128;
  };
  inline bool active(int iRow) const
  {return ((status_[iRow]&128)!=0);};
  void createStatus() ;
  inline void allSlackBasis()
  { createStatus();};
    
  inline int lastBadIteration() const
  {return lastBadIteration_;};
  inline int progressFlag() const
  {return progressFlag_;};
  inline void forceFactorization(int value)
  { forceFactorization_ = value;};
  inline double rawObjectiveValue() const
  { return objectiveValue_;};

protected:


  double columnPrimalInfeasibility_;
  double rowPrimalInfeasibility_;
  int columnPrimalSequence_;
  int rowPrimalSequence_;
  double columnDualInfeasibility_;
  double rowDualInfeasibility_;
  int columnDualSequence_;
  int rowDualSequence_;
  double primalToleranceToGetOptimal_;
  double remainingDualInfeasibility_;
  double largeValue_;
  double objectiveScale_;
  double largestPrimalError_;
  double largestDualError_;
  double largestSolutionError_;
  double dualBound_;
  double alpha_;
  double theta_;
  double lowerIn_;
  double valueIn_;
  double upperIn_;
  double dualIn_;
  double lowerOut_;
  double valueOut_;
  double upperOut_;
  double dualOut_;
  double dualTolerance_;
  double primalTolerance_;
  double sumDualInfeasibilities_;
  double sumPrimalInfeasibilities_;
  double infeasibilityCost_;
  double sumOfRelaxedDualInfeasibilities_;
  double sumOfRelaxedPrimalInfeasibilities_;
  double * lower_;
  double * rowLowerWork_;
  double * columnLowerWork_;
  double * upper_;
  double * rowUpperWork_;
  double * columnUpperWork_;
  double * cost_;
  double * rowObjectiveWork_;
  double * objectiveWork_;
  CoinIndexedVector * rowArray_[6];
  CoinIndexedVector * columnArray_[6];
  int sequenceIn_;
  int directionIn_;
  int sequenceOut_;
  int directionOut_;
  int pivotRow_;
  int lastGoodIteration_;
  double * dj_;
  double * rowReducedCost_;
  double * reducedCostWork_;
  double * solution_;
  double * rowActivityWork_;
  double * columnActivityWork_;
  int numberDualInfeasibilities_;
  int numberDualInfeasibilitiesWithoutFree_;
  int numberPrimalInfeasibilities_;
  int numberRefinements_;
  ClpDualRowPivot * dualRowPivot_;
  ClpPrimalColumnPivot * primalColumnPivot_;
  int * pivotVariable_;
  ClpFactorization * factorization_;
  // ****** get working simply then make coding more efficient
  // on full matrix operations
  double * rowScale_;
  double * savedSolution_;
  double * columnScale_;
  int scalingFlag_;
  int numberTimesOptimal_;
  int changeMade_;
  int algorithm_;
  int forceFactorization_;
  int perturbation_;
  unsigned char * saveStatus_;
  ClpNonLinearCost * nonLinearCost_;
  unsigned int specialOptions_;
  int lastBadIteration_;
  int numberFake_;
  int progressFlag_;
  int firstFree_;
  float incomingInfeasibility_;
  float allowedInfeasibility_;
  ClpSimplexProgress * progress_;
};
//#############################################################################
void
ClpSimplexUnitTest(const std::string & mpsDir,
                   const std::string & netlibDir);


class ClpSimplexProgress {

public:



    ClpSimplexProgress (  );

    ClpSimplexProgress ( ClpSimplex * model );

  ClpSimplexProgress(const ClpSimplexProgress &);

    ClpSimplexProgress & operator=(const ClpSimplexProgress & rhs);
   ~ClpSimplexProgress (  );

    int looping (  );
  void startCheck();
  int cycle(int in, int out,int wayIn,int wayOut); 

  double lastObjective(int back=1) const;
  void modifyObjective(double value);
  int lastIterationNumber(int back=1) const;
  inline void newOddState()
  { oddState_= - oddState_-1;};
  inline void endOddState()
  { oddState_=abs(oddState_);};
  inline void clearOddState() 
  { oddState_=0;};
  inline int oddState() const
  { return oddState_;};
  inline int badTimes() const
  { return numberBadTimes_;};
  inline void clearBadTimes()
  { numberBadTimes_=0;};


#define CLP_PROGRESS 5


  double objective_[CLP_PROGRESS];
  double infeasibility_[CLP_PROGRESS];
#define CLP_CYCLE 12

  double obj_[CLP_CYCLE];
  int in_[CLP_CYCLE];
  int out_[CLP_CYCLE];
  char way_[CLP_CYCLE];
  ClpSimplex * model_;
  int numberInfeasibilities_[CLP_PROGRESS];
  int iterationNumber_[CLP_PROGRESS];
  int numberTimes_;
  int numberBadTimes_;
  int oddState_;
};
#endif

---