Changeset 8 for branches


Ignore:
Timestamp:
Aug 5, 2002 11:57:49 AM (17 years ago)
Author:
forrest
Message:

Improving reliability using IBM Burlington problems (Primal)

Location:
branches/devel-1
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • branches/devel-1/ClpMessage.cpp

    r2 r8  
    2323  {CLP_SIMPLEX_STATUS,5,1,"%d  Objective %g%? Primal infeas %g (%d)%? Dual infeas %g (%d)%? without free dual infeas (%d)"},
    2424  {CLP_DUAL_BOUNDS,5,3,"Looking optimal checking bounds with %g"},
    25   {CLP_SIMPLEX_ACCURACY,6,1,"Primal error %g, dual error %g"},
     25  {CLP_SIMPLEX_ACCURACY,6,3,"Primal error %g, dual error %g"},
    2626  {CLP_SIMPLEX_BADFACTOR,7,1,"Singular factorization of basis - status %d"},
    2727  {CLP_SIMPLEX_BOUNDTIGHTEN,8,1,"Bounds were tightened %d times"},
    2828  {CLP_SIMPLEX_INFEASIBILITIES,9,1,"%d infeasibilities"},
    29   {CLP_SIMPLEX_FLAG,10,1,"Flagging variable %c%d"},
     29  {CLP_SIMPLEX_FLAG,10,3,"Flagging variable %c%d"},
    3030  {CLP_SIMPLEX_GIVINGUP,11,2,"Stopping as close enough"},
    3131  {CLP_DUAL_CHECKB,12,3,"Looking infeasible checking bounds with %g"},
  • branches/devel-1/ClpNonLinearCost.cpp

    r2 r8  
    3131  changeCost_(0.0),
    3232  largestInfeasibility_(0.0),
     33  sumInfeasibilities_(0.0),
    3334  convex_(true)
    3435{
     
    5455  changeCost_=0.0;
    5556  double infeasibilityCost = model_->infeasibilityCost();
     57  sumInfeasibilities_=0.0;
    5658  largestInfeasibility_=0.0;
    5759
     
    115117  double infeasibilityCost = model_->infeasibilityCost();
    116118  largestInfeasibility_=0.0;
     119  sumInfeasibilities_=0.0;
    117120
    118121  int iSequence;
     
    196199  changeCost_(0.0),
    197200  largestInfeasibility_(0.0),
     201  sumInfeasibilities_(0.0),
    198202  convex_(true)
    199203
     
    215219    changeCost_ = rhs.changeCost_;
    216220    largestInfeasibility_ = rhs.largestInfeasibility_;
     221    sumInfeasibilities_ = rhs.sumInfeasibilities_;
    217222    convex_ = rhs.convex_;
    218223  }
     
    268273    changeCost_ = rhs.changeCost_;
    269274    largestInfeasibility_ = rhs.largestInfeasibility_;
     275    sumInfeasibilities_ = rhs.sumInfeasibilities_;
    270276    convex_ = rhs.convex_;
    271277  }
     
    284290  changeCost_=0.0;
    285291  largestInfeasibility_ = 0.0;
     292  sumInfeasibilities_ = 0.0;
    286293  double primalTolerance = model_->currentPrimalTolerance();
    287294 
     
    328335      // slot in here
    329336      if (value<lower[iSequence]-primalTolerance) {
    330         largestInfeasibility_ = max(largestInfeasibility_,
    331                                     lower[iSequence]-value);
     337        value = lower[iSequence]-value;
     338        sumInfeasibilities_ += value;
     339        largestInfeasibility_ = max(largestInfeasibility_,value);
    332340        changeCost_ -= lower[iSequence]*
    333341          (cost_[iRange]-cost[iSequence]);
    334342        numberInfeasibilities_++;
    335343      } else if (value>upper[iSequence]+primalTolerance) {
    336         largestInfeasibility_ = max(largestInfeasibility_,
    337                                     value-upper[iSequence]);
     344        value = value-upper[iSequence];
     345        sumInfeasibilities_ += value;
     346        largestInfeasibility_ = max(largestInfeasibility_,value);
    338347        changeCost_ -= upper[iSequence]*
    339348          (cost_[iRange]-cost[iSequence]);
  • branches/devel-1/ClpPackedMatrix.cpp

    r2 r8  
    9898  const double * elementByColumn = matrix_->getElements();
    9999  int numberColumns = matrix_->getNumCols();
     100  //memset(y,0,matrix_->getNumRows()*sizeof(double));
    100101  for (iColumn=0;iColumn<numberColumns;iColumn++) {
    101102    int j;
     
    121122  const double * elementByColumn = matrix_->getElements();
    122123  int numberColumns = matrix_->getNumCols();
     124  //memset(y,0,numberColumns*sizeof(double));
    123125  for (iColumn=0;iColumn<numberColumns;iColumn++) {
    124126    int j;
     
    145147  const double * elementByColumn = matrix_->getElements();
    146148  int numberColumns = matrix_->getNumCols();
     149  //memset(y,0,matrix_->getNumRows()*sizeof(double));
    147150  for (iColumn=0;iColumn<numberColumns;iColumn++) {
    148151    int j;
     
    176179  const double * elementByColumn = matrix_->getElements();
    177180  int numberColumns = matrix_->getNumCols();
     181  //memset(y,0,numberColumns*sizeof(double));
    178182  for (iColumn=0;iColumn<numberColumns;iColumn++) {
    179183    int j;
  • branches/devel-1/ClpPrimalColumnSteepest.cpp

    r7 r8  
    536536  // we can't really trust infeasibilities if there is dual error
    537537  double checkTolerance = 1.0e-8;
    538   if (!model_->factorization()->pivots())
     538  if (!model_->factorization()->pivots()&&
     539      model_->numberIterations()<model_->lastBadIteration()+200)
    539540    checkTolerance = 1.0e-6;
    540541  if (model_->largestDualError()>checkTolerance)
     
    664665      }
    665666      alternateWeights_->setNumElements(number);
     667#ifdef CLP_DEBUG
     668      // Can happen but I should clean up
    666669      assert(found>=0);
     670#endif
    667671      pivotSequence_ = found;
    668672      delete [] temp;
  • branches/devel-1/ClpSimplex.cpp

    r2 r8  
    125125  infeasibilityCost_(1.0e7),
    126126  nonLinearCost_(NULL),
    127   specialOptions_(0)
     127  specialOptions_(0),
     128  lastBadIteration_(-999999)
    128129
    129130{
     
    179180    for (iRow=0;iRow<numberRows_;iRow++) {
    180181      int iPivot=pivotVariable_[iRow];
    181       if (iPivot>=numberColumns_) {
    182         // slack
    183         save[iRow] = rowActivityWork_[iPivot-numberColumns_];
    184       } else {
    185         // column
    186         save[iRow] = columnActivityWork_[iPivot];
    187       }
     182      save[iRow] = solution_[iPivot];
    188183    }
    189184  }
     
    203198  }
    204199  if (valuesPass) {
     200#ifdef CLP_DEBUG
    205201    std::cout<<"Largest given infeasibility "<<oldValue
    206202             <<" now "<<nonLinearCost_->largestInfeasibility()<<std::endl;
     203#endif
    207204    int numberOut=0;
    208205    if (oldValue<1.0
     
    212209      int iRow;
    213210      int * sort = new int[numberRows_];
     211      // first put back solution and store difference
    214212      for (iRow=0;iRow<numberRows_;iRow++) {
    215213        int iPivot=pivotVariable_[iRow];
    216 
     214        double difference = fabs(solution_[iPivot]=save[iRow]);
     215        solution_[iPivot]=save[iRow];
     216        save[iRow]=difference;
     217      }
     218      for (iRow=0;iRow<numberRows_;iRow++) {
     219        int iPivot=pivotVariable_[iRow];
     220       
    217221        if (iPivot<numberColumns_) {
    218222          // column
    219           double difference= fabs(save[iRow] - columnActivityWork_[iPivot]);
     223          double difference= save[iRow];
    220224          if (difference>1.0e-4) {
    221225            sort[numberOut]=iPivot;
     
    245249  objectiveValue_ += objectiveModification;
    246250  checkDualSolution();
    247   if (handler_->logLevel()>3||(largestPrimalError_>1.0e-4||
    248                                largestDualError_>1.0e-4))
     251  if (handler_->logLevel()>3||(largestPrimalError_>1.0e-2||
     252                               largestDualError_>1.0e-2))
    249253    handler_->message(CLP_SIMPLEX_ACCURACY,messages_)
    250254      <<largestPrimalError_
     
    278282  for (iRow=0;iRow<numberRows_;iRow++) {
    279283    int iPivot=pivotVariable_[iRow];
    280     if (iPivot>=numberColumns_) {
    281       // slack
    282       rowActivityWork_[iPivot-numberColumns_] = 0.0;
    283     } else {
    284       // column
    285       columnActivityWork_[iPivot] = 0.0;
    286     }
     284    solution_[iPivot] = 0.0;
    287285  }
    288286  times(-1.0,columnActivityWork_,array);
     
    305303    for (iRow=0;iRow<numberRows_;iRow++) {
    306304      int iPivot=pivotVariable_[iRow];
    307       if (iPivot>=numberColumns_) {
    308         // slack
    309         rowActivityWork_[iPivot-numberColumns_] = array[iRow];
    310       } else {
    311         // column
    312         columnActivityWork_[iPivot] = array[iRow];
    313       }
     305      solution_[iPivot] = array[iRow];
    314306    }
    315307    // check Ax == b  (for all)
     
    359351  for (iRow=0;iRow<numberRows_;iRow++) {
    360352    int iPivot=pivotVariable_[iRow];
    361     if (iPivot>=numberColumns_) {
    362       // slack
    363       rowActivityWork_[iPivot-numberColumns_] = array[iRow];
    364     } else {
    365       // column
    366       columnActivityWork_[iPivot] = array[iRow];
    367      }
     353    solution_[iPivot] = array[iRow];
    368354  }
    369355  delete [] previous;
     
    530516    if (!valuesPass) {
    531517      // not values pass so set to bounds
     518      bool allSlack=true;
    532519      if (status_) {
     520        for (iRow=0;iRow<numberRows_;iRow++) {
     521          if (getRowStatus(iRow)!=ClpSimplex::basic) {
     522            allSlack=false;
     523            break;
     524          }
     525        }
     526      }
     527      if (!allSlack) {
    533528        // set values from warm start (if sensible)
    534529        int numberBasic=0;
     
    677672      // values pass has less coding
    678673      // make row activities correct
    679       times(-1.0,columnActivityWork_,rowActivityWork_);
     674      memset(rowActivityWork_,0,numberRows_*sizeof(double));
     675      times(1.0,columnActivityWork_,rowActivityWork_);
    680676      if (status_) {
    681677        // set values from warm start (if sensible)
     
    790786              ClpSimplex::basic) {
    791787            // take out
    792             double lower=columnLowerWork_[iColumn];
    793             double upper=columnUpperWork_[iColumn];
    794             double value=columnActivityWork_[iColumn];
    795             if (lower>-largeValue_||upper<largeValue_) {
    796               if (fabs(value-lower)<fabs(value-upper)) {
    797                 setColumnStatus(iColumn,ClpSimplex::atLowerBound);
    798                 columnActivityWork_[iColumn]=lower;
     788            if (!valuesPass) {
     789              double lower=columnLowerWork_[iColumn];
     790              double upper=columnUpperWork_[iColumn];
     791              double value=columnActivityWork_[iColumn];
     792              if (lower>-largeValue_||upper<largeValue_) {
     793                if (fabs(value-lower)<fabs(value-upper)) {
     794                  setColumnStatus(iColumn,ClpSimplex::atLowerBound);
     795                  columnActivityWork_[iColumn]=lower;
     796                } else {
     797                  setColumnStatus(iColumn,ClpSimplex::atUpperBound);
     798                  columnActivityWork_[iColumn]=upper;
     799                }
    799800              } else {
    800                 setColumnStatus(iColumn,ClpSimplex::atUpperBound);
    801                 columnActivityWork_[iColumn]=upper;
     801                setColumnStatus(iColumn,ClpSimplex::isFree);
    802802              }
    803803            } else {
    804               setColumnStatus(iColumn,ClpSimplex::isFree);
     804              setColumnStatus(iColumn,ClpSimplex::superBasic);
    805805            }
    806806          }
     
    10591059  infeasibilityCost_(1.0e7),
    10601060  nonLinearCost_(NULL),
    1061   specialOptions_(0)
     1061  specialOptions_(0),
     1062  lastBadIteration_(-999999)
    10621063{
    10631064  int i;
     
    11431144  infeasibilityCost_(1.0e7),
    11441145  nonLinearCost_(NULL),
    1145   specialOptions_(0)
     1146  specialOptions_(0),
     1147  lastBadIteration_(-999999)
    11461148{
    11471149  int i;
     
    12701272  infeasibilityCost_ = rhs.infeasibilityCost_;
    12711273  specialOptions_ = rhs.specialOptions_;
     1274  lastBadIteration_ = rhs.lastBadIteration_;
    12721275  if (rhs.nonLinearCost_!=NULL)
    12731276    nonLinearCost_ = new ClpNonLinearCost(*rhs.nonLinearCost_);
     
    13921395    if (infeasibility>largestSolutionError_)
    13931396      largestSolutionError_=infeasibility;
    1394     if (columnLowerWork_[iColumn]!=columnUpperWork_[iColumn])
     1397    if (columnUpperWork_[iColumn]-columnLowerWork_[iColumn]>primalTolerance_)
    13951398      clearFixed(iColumn);
    13961399    else
     
    18031806ClpSimplex::setBasis ( const OsiWarmStartBasis & basis)
    18041807{
    1805   if (basis.getNumStructural()) {
    1806     // transform basis to status arrays
    1807     int iRow,iColumn;
    1808     if (!status_) {
    1809       /*
    1810         get status arrays
    1811         OsiWarmStartBasis would seem to have overheads and we will need
    1812         extra bits anyway.
    1813       */
    1814       status_ = new unsigned char [numberColumns_+numberRows_];
    1815       memset(status_,0,(numberColumns_+numberRows_)*sizeof(char));
    1816     }
    1817     OsiWarmStartBasis basis2 = basis;
    1818     // resize if necessary
    1819     basis2.resize(numberRows_,numberColumns_);
    1820     // move status
    1821     for (iRow=0;iRow<numberRows_;iRow++) {
    1822       setRowStatus(iRow,
    1823                    (Status) basis2.getArtifStatus(iRow));
    1824     }
    1825     for (iColumn=0;iColumn<numberColumns_;iColumn++) {
    1826       setColumnStatus(iColumn,
    1827                       (Status) basis2.getStructStatus(iColumn));
    1828     }
    1829   } else {
    1830     // null basis so just delete any status arrays
    1831     delete [] status_;
    1832     status_=NULL;
     1808  // transform basis to status arrays
     1809  int iRow,iColumn;
     1810  if (!status_) {
     1811    /*
     1812      get status arrays
     1813      OsiWarmStartBasis would seem to have overheads and we will need
     1814      extra bits anyway.
     1815    */
     1816    status_ = new unsigned char [numberColumns_+numberRows_];
     1817    memset(status_,0,(numberColumns_+numberRows_)*sizeof(char));
     1818  }
     1819  OsiWarmStartBasis basis2 = basis;
     1820  // resize if necessary
     1821  basis2.resize(numberRows_,numberColumns_);
     1822  // move status
     1823  for (iRow=0;iRow<numberRows_;iRow++) {
     1824    setRowStatus(iRow,
     1825                 (Status) basis2.getArtifStatus(iRow));
     1826  }
     1827  for (iColumn=0;iColumn<numberColumns_;iColumn++) {
     1828    setColumnStatus(iColumn,
     1829                    (Status) basis2.getStructStatus(iColumn));
    18331830  }
    18341831}
  • branches/devel-1/ClpSimplexDual.cpp

    r7 r8  
    449449                  <<OsiMessageEol;
    450450                setFlagged(sequenceOut_);
     451                lastBadIteration_ = numberIterations_; // say be more cautious
    451452                rowArray_[0]->clear();
    452453                rowArray_[1]->clear();
     
    535536                <<OsiMessageEol;
    536537              setFlagged(sequenceOut_);
     538              lastBadIteration_ = numberIterations_; // say be more cautious
    537539              rowArray_[0]->clear();
    538540              rowArray_[1]->clear();
  • branches/devel-1/ClpSimplexPrimal.cpp

    r7 r8  
    263263    if (internalFactorize(0+10*ifValuesPass))
    264264      return 1; // some error
    265 
    266265    // for this we need clean basis so it is after factorize
    267266    numberThrownOut=gutsOfSolution(rowActivityWork_,columnActivityWork_,
     
    289288  // This says whether to restore things etc
    290289  int factorType=0;
     290  // Save iteration number
     291  int saveNumber = -1;
    291292  /*
    292293    Status of problem:
     
    319320      perturb();
    320321#endif
     322    // If we have done no iterations - special
     323    if (saveNumber==numberIterations_)
     324      factorType=3;
    321325    // may factorize, checks if problem finished
    322326    statusOfProblemInPrimal(lastCleaned,factorType);
     
    329333
    330334    // Save iteration number
    331     int saveNumber = numberIterations_;
     335    saveNumber = numberIterations_;
    332336
    333337    // status stays at -1 while iterating, >=0 finished, -2 to invert
     
    446450                <<OsiMessageEol;
    447451              setFlagged(sequenceIn_);
     452              lastBadIteration_ = numberIterations_; // say be more cautious
    448453              rowArray_[1]->clear();
    449454              pivotRow_=-1;
     
    458463#endif
    459464        if (saveDj*dualIn_<1.0e-20||
    460             fabs(saveDj-dualIn_)>1.0e-7*(1.0+fabs(dualIn_))) {
     465            fabs(saveDj-dualIn_)>1.0e-5*(1.0+fabs(saveDj))) {
    461466          handler_->message(CLP_PRIMAL_DJ,messages_)
    462467            <<saveDj<<dualIn_
     
    477482                <<OsiMessageEol;
    478483              setFlagged(sequenceIn_);
     484              lastBadIteration_ = numberIterations_; // say be more cautious
    479485              rowArray_[1]->clear();
    480486              pivotRow_=-1;
     
    505511                <<OsiMessageEol;
    506512              setFlagged(sequenceIn_);
     513              lastBadIteration_ = numberIterations_; // say be more cautious
    507514              rowArray_[1]->clear();
    508515              pivotRow_=-1;
     
    619626  }
    620627
     628  // if infeasible get real values
     629  if (problemStatus_) {
     630    infeasibilityCost_=0.0;
     631    createRim(7);
     632    nonLinearCost_->checkInfeasibilities(true);
     633    sumPrimalInfeasibilities_=nonLinearCost_->sumInfeasibilities();
     634    numberPrimalInfeasibilities_= nonLinearCost_->numberInfeasibilities();
     635  }
    621636  // at present we are leaving factorization around
    622637  // maybe we should empty it
     
    690705  // we may wish to say it is optimal even if infeasible
    691706  bool alwaysOptimal = (specialOptions_&1)!=0;
    692   if (dualFeasible()||problemStatus_==-4) {
     707  if (dualFeasible()||problemStatus_==-4||type==3) {
    693708    if (nonLinearCost_->numberInfeasibilities()&&!alwaysOptimal) {
    694709      //may need infeasiblity cost changed
  • branches/devel-1/Makefile

    r7 r8  
    88# I seem to need this at present
    99OPTFLAG := -O2
    10 OptLevel := -g
    11 OPTFLAG := -g
     10#OptLevel := -g
     11#OPTFLAG := -g
    1212ifeq ($(OptLevel),-g)
    1313    CXXFLAGS += -DCLP_DEBUG
  • branches/devel-1/include/ClpModel.hpp

    r2 r8  
    8787  void deleteRows(int number, const int * which);
    8888  /// Deletes columns
     89 
    8990  void deleteColumns(int number, const int * which);
    9091  /** Borrow model.  This is so we dont have to copy large amounts
  • branches/devel-1/include/ClpNonLinearCost.hpp

    r2 r8  
    134134  inline double changeInCost() const
    135135  {return changeCost_;};
     136  /// Sum of infeasibilities
     137  inline double sumInfeasibilities() const
     138  {return sumInfeasibilities_;};
    136139  /// Largest infeasibility
    137140  inline double largestInfeasibility() const
     
    168171  /// Largest infeasibility
    169172  double largestInfeasibility_;
     173  /// Sum of infeasibilities
     174  double sumInfeasibilities_;
    170175  /// If all non-linear costs convex
    171176  bool convex_;
  • branches/devel-1/include/ClpSimplex.hpp

    r2 r8  
    442442  inline bool flagged(int sequence) const
    443443  {return (((status_[sequence]>>6)&1)!=0);};
     444  /// So we know when to be cautious
     445  inline int lastBadIteration() const
     446  {return lastBadIteration_;};
    444447  //@}
    445448
     
    611614  /// For advanced options
    612615  unsigned int specialOptions_;
     616  /// So we know when to be cautious
     617  int lastBadIteration_;
    613618  //@}
    614619};
Note: See TracChangeset for help on using the changeset viewer.