Changeset 2322 for trunk


Ignore:
Timestamp:
Mar 31, 2018 12:42:43 PM (19 months ago)
Author:
forrest
Message:

changes to improve ray handling

Location:
trunk/Clp/src
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • trunk/Clp/src/ClpModel.hpp

    r2271 r2322  
    10391039         524288 - Clp fast dual
    10401040         1048576 - don't need to finish dual (can return 3)
    1041          2097152 - zero costs!
     1041         2097152 - ray even if >2 pivots AND if problem is "crunched"
    10421042         4194304 - don't scale integer variables
    10431043         8388608 - Idiot when not really sure about it
     1044         16777216 - zero costs!
    10441045         NOTE - many applications can call Clp but there may be some short cuts
    10451046                which are taken which are not guaranteed safe from all applications.
  • trunk/Clp/src/ClpSimplex.cpp

    r2315 r2322  
    718718     }
    719719#if CAN_HAVE_ZERO_OBJ>1
    720      if ((specialOptions_&2097152)==0) {
     720     if ((specialOptions_&16777216)==0) {
    721721#endif
    722722     computeDuals(givenDuals);
     
    31243124          if (distanceUp < -primalTolerance) {
    31253125               double infeasibility = -distanceUp;
     3126               if (getStatus(iSequence) != basic)
     3127                 moreSpecialOptions_ &= ~8;  // Say superbasic variables exist
    31263128               sumPrimalInfeasibilities_ += infeasibility - primalTolerance_;
    31273129               if (infeasibility > relaxedToleranceP)
     
    31353137          } else if (distanceDown < -primalTolerance) {
    31363138               double infeasibility = -distanceDown;
     3139               if (getStatus(iSequence) != basic)
     3140                 moreSpecialOptions_ &= ~8;  // Say superbasic variables exist
    31373141               sumPrimalInfeasibilities_ += infeasibility - primalTolerance_;
    31383142               if (infeasibility > relaxedToleranceP)
     
    56475651          problemStatus_ = 0; // ignore
    56485652     if (problemStatus_==1&&((specialOptions_&(1024 | 4096)) == 0 || (specialOptions_ & 32) != 0)
    5649          &&numberFake_) {
     5653         &&(static_cast<ClpSimplexDual *>(this))->checkFakeBounds()) {
    56505654       problemStatus_ = 10; // clean up in primal as fake bounds
    56515655     }
     
    56975701          baseIteration_ = 0;
    56985702          bool inCbcOrOther = (specialOptions_ & 0x03000000) != 0;
    5699           if (inCbcOrOther) {
     5703          if (inCbcOrOther && (specialOptions_&32)==0) {
    57005704            delete [] ray_;
    57015705            ray_ = NULL;
     
    57905794     //if (problemStatus_==1&&lastAlgorithm==1)
    57915795     //returnCode=10; // so will do primal after postsolve
     5796#ifdef CHECK_RAY
     5797     if (problemStatus_==1&&ray_) {
     5798       double * lower = rowLower_;
     5799       double * upper = rowUpper_;
     5800       double * solution = rowActivity_;
     5801       double * dj = dual_;
     5802       assert(ray_);
     5803       double largestBad=0.0;
     5804       double largestBadDj=0.0;
     5805       for (int iRow = 0; iRow < numberRows_; iRow++) {
     5806         if (upper[iRow]==lower[iRow])
     5807           continue;
     5808         if (solution[iRow]<lower[iRow]+primalTolerance_) {
     5809           largestBadDj=CoinMax(largestBadDj,-dj[iRow]);
     5810           largestBad=CoinMax(largestBad,ray_[iRow]);
     5811         } else if (solution[iRow]>upper[iRow]-primalTolerance_) {
     5812           largestBadDj=CoinMax(largestBadDj,dj[iRow]);
     5813           largestBad=CoinMax(largestBad,-ray_[iRow]);
     5814         }
     5815       }
     5816       double * result = new double[numberColumns_];
     5817       CoinFillN ( result, numberColumns_, 0.0);
     5818       this->matrix()->transposeTimes(ray_, result);
     5819       lower = columnLower_;
     5820       upper = columnUpper_;
     5821       solution = columnActivity_;
     5822       dj = reducedCost_;
     5823       for (int iColumn = 0; iColumn < numberColumns_; iColumn++) {
     5824         // should have been ..Times -1.0
     5825         result[iColumn]=-result[iColumn];
     5826         if (upper[iColumn]==lower[iColumn])
     5827           continue;
     5828         if (solution[iColumn]<lower[iColumn]+primalTolerance_) {
     5829           largestBadDj=CoinMax(largestBadDj,-dj[iColumn]);
     5830           largestBad=CoinMax(largestBad,result[iColumn]);
     5831         } else if (solution[iColumn]>upper[iColumn]-primalTolerance_) {
     5832           largestBadDj=CoinMax(largestBadDj,dj[iColumn]);
     5833           largestBad=CoinMax(largestBad,-result[iColumn]);
     5834         }
     5835       }
     5836       if (largestBad>1.0e-5||largestBadDj>1.0e-5) {
     5837         if (numberPrimalInfeasibilities_==1)
     5838           printf("BAD ");
     5839         printf("bad ray %g bad dj %g\n",largestBad,largestBadDj);
     5840       }
     5841       delete [] result;
     5842     }
     5843#endif
    57925844     return returnCode;
    57935845}
     
    90489100          } else {
    90499101               // using previous factorization - we assume fine
    9050                if ((moreSpecialOptions_ & 8) == 0) {
     9102               if ((moreSpecialOptions_ & 16777216) == 0) {
    90519103                    // but we need to say not optimal
    90529104                    numberPrimalInfeasibilities_ = 1;
     
    1089610948     info->nNodes_ = 0;
    1089710949     // say can declare optimal
    10898      moreSpecialOptions_ |= 8;
     10950     moreSpecialOptions_ |= 16777216;
    1089910951     int saveMaxIterations = maximumIterations();
    1090010952     setMaximumIterations((((moreSpecialOptions_&2048)==0) ? 200 : 2000)
     
    1160611658     assert (nodeInfo);
    1160711659     // say can declare optimal
    11608      moreSpecialOptions_ |= 8;
     11660     moreSpecialOptions_ |= 16777216;
    1160911661     double limit = 0.0;
    1161011662     getDblParam(ClpDualObjectiveLimit, limit);
  • trunk/Clp/src/ClpSimplex.hpp

    r2271 r2322  
    12731273         2 bit - if presolved problem infeasible return
    12741274         4 bit - keep arrays like upper_ around
    1275          8 bit - if factorization kept can still declare optimal at once
     1275         8 bit - no free or superBasic variables
    12761276         16 bit - if checking replaceColumn accuracy before updating
    12771277         32 bit - say optimal if primal feasible!
     
    12901290         524288 bit - stop when primal feasible
    12911291         1048576 bit - stop when primal feasible after n-1000000 iterations
     1292         2097152 bit - no primal in fastDual2 if feasible
     1293         4194304 bit - tolerances have been changed by code
     1294         8388608 bit - tolerances are dynamic (at first)
     1295         16777216 bit - if factorization kept can still declare optimal at once
    12921296     */
    12931297     inline int moreSpecialOptions() const {
     
    13221326         4194304 bit - tolerances have been changed by code
    13231327         8388608 bit - tolerances are dynamic (at first)
     1328         16777216 bit - if factorization kept can still declare optimal at once
    13241329     */
    13251330     inline void setMoreSpecialOptions(int value) {
  • trunk/Clp/src/ClpSimplexDual.cpp

    r2293 r2322  
    510510            //         i,cost_[i],solution_[i],lower_[i],upper_[i]);
    511511            //}
    512           }
    513           if ((specialOptions_&2097152)!=0&&problemStatus_==1&&!ray_&&
    514               !numberRayTries && numberIterations_) {
     512          } else if ((specialOptions_&(32|2097152))!=0&&
     513                     problemStatus_==1&&!ray_&&
     514                     !numberRayTries && numberIterations_) {
    515515            numberRayTries=1;
    516516            problemStatus_=-1;
     
    610610     algorithm_ = -1;
    611611     moreSpecialOptions_ &= ~16; // clear check replaceColumn accuracy
     612     delete [] ray_;
     613     ray_ = NULL;
    612614     // save data
    613615     ClpDataSave data = saveData();
     
    12811283                      acceptablePivot_ = - fabs(acceptablePivot_); // stop early exit
    12821284#if CAN_HAVE_ZERO_OBJ>1
    1283                     if ((specialOptions_&2097152)!=0)
     1285                    if ((specialOptions_&16777216)!=0)
    12841286                      theta_=0.0;
    12851287#endif
     
    14611463                    if (candidate == -1) {
    14621464#if CLP_CAN_HAVE_ZERO_OBJ>1
    1463                     if ((specialOptions_&2097152)==0) {
     1465                    if ((specialOptions_&16777216)==0) {
    14641466#endif
    14651467                         // make sure incoming doesn't count
     
    18321834                         printf("** no column pivot\n");
    18331835#endif
     1836                    delete [] ray_;
     1837                    ray_ = NULL;
    18341838                    if ((factorization_->pivots() < 2
    18351839                         ||((specialOptions_&2097152)!=0&&factorization_->pivots()<50))
     
    18371841                         //&&goodAccuracy()) {
    18381842                         // If not in branch and bound etc save ray
    1839                          delete [] ray_;
    18401843                         if ((specialOptions_&(1024 | 4096)) == 0 || (specialOptions_ & (32|2097152)) != 0) {
    18411844                              // create ray anyway
     
    53925395     if (!numberPrimalInfeasibilities_ && ((!numberDualInfeasibilitiesWithoutFree_ &&
    53935396                                            numberDualInfeasibilities_)||
    5394                                            (moreSpecialOptions_&2097152)!=0))
     5397                                           (moreSpecialOptions_&16777216)!=0))
    53955398          problemStatus_ = 10;
    53965399     // dual bound coming in
     
    55115514                                   (numberChanged_ || (specialOptions_ & 4096) == 0)) {
    55125515#if CLP_CAN_HAVE_ZERO_OBJ
    5513                            if ((specialOptions_&2097152)==0) {
     5516                           if ((specialOptions_&16777216)==0) {
    55145517#endif
    55155518                              doOriginalTolerance = 2;
     
    57935796                         dualTolerance_ = 1.0e50;
    57945797#if CLP_CAN_HAVE_ZERO_OBJ>1
    5795                     if ((specialOptions_&2097152)==0) {
     5798                    if ((specialOptions_&16777216)==0) {
    57965799#endif
    57975800                    updateDualsInDual(rowArray_[0], columnArray_[0], rowArray_[1],
     
    58355838                         dualTolerance_ = 1.0e50;
    58365839#if CLP_CAN_HAVE_ZERO_OBJ>1
    5837                     if ((specialOptions_&2097152)==0) {
     5840                    if ((specialOptions_&16777216)==0) {
    58385841#endif
    58395842                    updateDualsInDual(rowArray_[0], columnArray_[0], rowArray_[1],
  • trunk/Clp/src/ClpSimplexOther.cpp

    r2271 r2322  
    39123912                    //columnArray_[0]->cleanAndPackSafe(1.0e-60);
    39133913#if CLP_CAN_HAVE_ZERO_OBJ
    3914                     if ((specialOptions_&2097152)==0) {
     3914                    if ((specialOptions_&16777216)==0) {
    39153915#endif
    39163916                      nswapped = reinterpret_cast<ClpSimplexDual *> ( this)->updateDualsInDual(rowArray_[0], columnArray_[0],
     
    40954095                    // update primal solution
    40964096#if CLP_CAN_HAVE_ZERO_OBJ
    4097                     if ((specialOptions_&2097152)!=0)
     4097                    if ((specialOptions_&16777216)!=0)
    40984098                      theta_=0.0;
    40994099#endif
     
    41644164                    objectiveChange = 0.0;
    41654165#if CLP_CAN_HAVE_ZERO_OBJ
    4166                     if ((specialOptions_&2097152)==0) {
     4166                    if ((specialOptions_&16777216)==0) {
    41674167#endif
    41684168                      for (int i=0;i<numberTotal;i++)
     
    41824182#if CLP_CAN_HAVE_ZERO_OBJ>1
    41834183#ifdef COIN_REUSE_RANDOM
    4184                          if ((specialOptions_&2097152)!=0) {
     4184                         if ((specialOptions_&16777216)!=0) {
    41854185                           dj_[sequenceOut_] = 1.0e-9*(1.0+CoinDrand48());;
    41864186                         }
     
    41924192#if CLP_CAN_HAVE_ZERO_OBJ>1
    41934193#ifdef COIN_REUSE_RANDOM
    4194                          if ((specialOptions_&2097152)!=0) {
     4194                         if ((specialOptions_&16777216)!=0) {
    41954195                           dj_[sequenceOut_] = -1.0e-9*(1.0+CoinDrand48());;
    41964196                         }
  • trunk/Clp/src/ClpSimplexPrimal.cpp

    r2293 r2322  
    14051405                         if (infeasibilityCost_ >= MAX_INFEASIBILITY_COST ||
    14061406                                   numberDualInfeasibilities_ == 0) {
    1407                               // we are infeasible - use as ray
    1408                               delete [] ray_;
    1409                               ray_ = new double [numberRows_];
    1410                               // swap sign
    1411                               for (int i=0;i<numberRows_;i++)
    1412                                 ray_[i] = -dual_[i];
     1407                              // we are infeasible - use as ray
     1408                              delete [] ray_;
     1409                              //if ((specialOptions_&(32|0x01000000))!=0x01000000) {
     1410                              if ((specialOptions_&0x01000000)==0) {
     1411                                ray_ = new double [numberRows_];
     1412                                double * saveObjective = CoinCopyOfArray(objective(),
     1413                                                                         numberColumns_);
     1414                                memset(objective(),0,numberColumns_*sizeof(double));
     1415                                infeasibilityCost_ = 1.0;
     1416                                createRim(4);
     1417                                nonLinearCost_->checkInfeasibilities(primalTolerance_);
     1418                                gutsOfSolution(NULL, NULL, false);
     1419                                memcpy(objective(),saveObjective,numberColumns_*sizeof(double));
     1420                                delete [] saveObjective;
     1421                                // swap sign
     1422                                for (int i=0;i<numberRows_;i++)
     1423                                  ray_[i] = -dual_[i];
     1424                              } else {
     1425                                ray_ = NULL;
     1426                              }
    14131427#ifdef PRINT_RAY_METHOD
    14141428                              printf("Primal creating infeasibility ray\n");
  • trunk/Clp/src/OsiClp/OsiClpSolverInterface.cpp

    r2315 r2322  
    74867486  }
    74877487}
    7488 
    74897488// Crunch down model
    74907489void
     
    78567855            //printf("%d basic (%d non-unit) %d non-basic\n",n,nn,k);
    78577856            modelPtr_->ray_=ray;
     7857#ifdef CHECK_RAY
     7858            {
     7859              double * lower = modelPtr_->rowLower_;
     7860              double * upper = modelPtr_->rowUpper_;
     7861              double * solution = modelPtr_->rowActivity_;
     7862              double * dj = modelPtr_->dual_;
     7863              double largestBad=0.0;
     7864              double largestBadDj=0.0;
     7865              for (int iRow = 0; iRow < modelPtr_->numberRows_; iRow++) {
     7866                if (upper[iRow]==lower[iRow])
     7867                  continue;
     7868                if (solution[iRow]<lower[iRow]+modelPtr_->primalTolerance_) {
     7869                  largestBadDj=CoinMax(largestBadDj,-dj[iRow]);
     7870                  largestBad=CoinMax(largestBad,ray[iRow]);
     7871                } else if (solution[iRow]>upper[iRow]-modelPtr_->primalTolerance_) {
     7872                  largestBadDj=CoinMax(largestBadDj,dj[iRow]);
     7873                  largestBad=CoinMax(largestBad,-ray[iRow]);
     7874                }
     7875              }
     7876              double * result = new double[modelPtr_->numberColumns_];
     7877              CoinFillN ( result, modelPtr_->numberColumns_, 0.0);
     7878              modelPtr_->matrix()->transposeTimes(ray, result);
     7879              lower = modelPtr_->columnLower_;
     7880              upper = modelPtr_->columnUpper_;
     7881              solution = modelPtr_->columnActivity_;
     7882              dj = modelPtr_->reducedCost_;
     7883              for (int iColumn = 0; iColumn < modelPtr_->numberColumns_; iColumn++) {
     7884                // should have been ..Times -1.0
     7885                result[iColumn]=-result[iColumn];
     7886                if (upper[iColumn]==lower[iColumn])
     7887                  continue;
     7888                if (solution[iColumn]<lower[iColumn]+modelPtr_->primalTolerance_) {
     7889                  largestBadDj=CoinMax(largestBadDj,-dj[iColumn]);
     7890                  largestBad=CoinMax(largestBad,result[iColumn]);
     7891                } else if (solution[iColumn]>upper[iColumn]-modelPtr_->primalTolerance_) {
     7892                  largestBadDj=CoinMax(largestBadDj,dj[iColumn]);
     7893                  largestBad=CoinMax(largestBad,-result[iColumn]);
     7894                }
     7895              }
     7896              if (largestBad>1.0e-5||largestBadDj>1.0e-5) {
     7897                if (modelPtr_->numberPrimalInfeasibilities_==1)
     7898                  printf("BAD ");
     7899                printf("bad ray %g bad dj %g\n",largestBad,largestBadDj);
     7900              }
     7901              delete [] result;
     7902            }
     7903#endif
    78587904            modelPtr_->directionOut_=small->directionOut_;
    78597905            if (small->sequenceOut_<small->numberColumns_)
  • trunk/Clp/src/unitTest.cpp

    r2271 r2322  
    351351     else
    352352          dirNetlib = dirsep == '/' ? "../../Data/Netlib/" : "..\\..\\Data\\Netlib\\";
    353 #if 0 //FACTORIZATION_STATISTICS==0
     353#if FACTORIZATION_STATISTICS==0
    354354     if (!empty.numberRows()) {
    355355          testingMessage( "Testing ClpSimplex\n" );
     
    15501550               for (int i = 0; i < numberRows; i++)
    15511551                    solution.setRowStatus(i, ClpSimplex::basic);
     1552               solution.objective()[1]=-2.0;
    15521553               solution.primal(1);
    15531554               assert (solution.secondaryStatus() == 102); // Came out at end of pass
     
    20762077                    // test infeasible and ray
    20772078                    solution.columnUpper()[0] = 0.0;
    2078 #ifdef DUAL
     2079#if 1 //def DUAL
    20792080                    solution.setDualBound(100.0);
    20802081                    solution.dual();
Note: See TracChangeset for help on using the changeset viewer.