Ignore:
Timestamp:
Aug 7, 2002 11:55:35 AM (18 years ago)
Author:
forrest
Message:

Breaking out whileIterating

File:
1 edited

Legend:

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

    r12 r14  
    325325    }
    326326
    327     // status stays at -1 while iterating, >=0 finished, -2 to invert
    328     // status -3 to go to top without an invert
    329     while (problemStatus_==-1) {
     327    // Do iterations
     328    whileIterating();
     329  }
     330
     331  // at present we are leaving factorization around
     332  // maybe we should empty it
     333  deleteRim();
     334  handler_->message(CLP_SIMPLEX_FINISHED+problemStatus_,messages_)
     335    <<objectiveValue()
     336    <<OsiMessageEol;
     337  // Restore any saved stuff
     338  perturbation_ = savePerturbation;
     339  factorization_->sparseThreshold(saveSparse);
     340  dualBound_ = saveDualBound_;
     341  return problemStatus_;
     342}
     343void
     344ClpSimplexDual::whileIterating()
     345{
     346  // status stays at -1 while iterating, >=0 finished, -2 to invert
     347  // status -3 to go to top without an invert
     348  while (problemStatus_==-1) {
    330349#ifdef CLP_DEBUG
    331       {
    332         int i;
    333         for (i=0;i<4;i++) {
    334           rowArray_[i]->checkClear();
    335         }   
    336         for (i=0;i<2;i++) {
    337           columnArray_[i]->checkClear();
    338         }   
    339       }     
     350    {
     351      int i;
     352      for (i=0;i<4;i++) {
     353        rowArray_[i]->checkClear();
     354      }   
     355      for (i=0;i<2;i++) {
     356        columnArray_[i]->checkClear();
     357      }   
     358    }     
    340359#endif
    341360#if CLP_DEBUG>2
    342       // very expensive
    343       if (numberIterations_>0&&numberIterations_<-801) {
    344         handler_->setLogLevel(63);
    345         double saveValue = objectiveValue_;
    346         double * saveRow1 = new double[numberRows_];
    347         double * saveRow2 = new double[numberRows_];
    348         memcpy(saveRow1,rowReducedCost_,numberRows_*sizeof(double));
    349         memcpy(saveRow2,rowActivityWork_,numberRows_*sizeof(double));
    350         double * saveColumn1 = new double[numberColumns_];
    351         double * saveColumn2 = new double[numberColumns_];
    352         memcpy(saveColumn1,reducedCostWork_,numberColumns_*sizeof(double));
    353         memcpy(saveColumn2,columnActivityWork_,numberColumns_*sizeof(double));
    354         gutsOfSolution(rowActivityWork_,columnActivityWork_);
    355         printf("xxx %d old obj %g, recomputed %g, sum dual inf %g\n",
    356                numberIterations_,
    357                saveValue,objectiveValue_,sumDualInfeasibilities_);
    358         if (saveValue>objectiveValue_+1.0e-2)
    359           printf("**bad**\n");
    360         memcpy(rowReducedCost_,saveRow1,numberRows_*sizeof(double));
    361         memcpy(rowActivityWork_,saveRow2,numberRows_*sizeof(double));
    362         memcpy(reducedCostWork_,saveColumn1,numberColumns_*sizeof(double));
    363         memcpy(columnActivityWork_,saveColumn2,numberColumns_*sizeof(double));
    364         delete [] saveRow1;
    365         delete [] saveRow2;
    366         delete [] saveColumn1;
    367         delete [] saveColumn2;
    368         objectiveValue_=saveValue;
     361    // very expensive
     362    if (numberIterations_>0&&numberIterations_<-801) {
     363      handler_->setLogLevel(63);
     364      double saveValue = objectiveValue_;
     365      double * saveRow1 = new double[numberRows_];
     366      double * saveRow2 = new double[numberRows_];
     367      memcpy(saveRow1,rowReducedCost_,numberRows_*sizeof(double));
     368      memcpy(saveRow2,rowActivityWork_,numberRows_*sizeof(double));
     369      double * saveColumn1 = new double[numberColumns_];
     370      double * saveColumn2 = new double[numberColumns_];
     371      memcpy(saveColumn1,reducedCostWork_,numberColumns_*sizeof(double));
     372      memcpy(saveColumn2,columnActivityWork_,numberColumns_*sizeof(double));
     373      gutsOfSolution(rowActivityWork_,columnActivityWork_);
     374      printf("xxx %d old obj %g, recomputed %g, sum dual inf %g\n",
     375             numberIterations_,
     376             saveValue,objectiveValue_,sumDualInfeasibilities_);
     377      if (saveValue>objectiveValue_+1.0e-2)
     378        printf("**bad**\n");
     379      memcpy(rowReducedCost_,saveRow1,numberRows_*sizeof(double));
     380      memcpy(rowActivityWork_,saveRow2,numberRows_*sizeof(double));
     381      memcpy(reducedCostWork_,saveColumn1,numberColumns_*sizeof(double));
     382      memcpy(columnActivityWork_,saveColumn2,numberColumns_*sizeof(double));
     383      delete [] saveRow1;
     384      delete [] saveRow2;
     385      delete [] saveColumn1;
     386      delete [] saveColumn2;
     387      objectiveValue_=saveValue;
     388    }
     389#endif
     390#ifdef CLP_DEBUG
     391    {
     392      int iSequence, number=numberRows_+numberColumns_;
     393      for (iSequence=0;iSequence<number;iSequence++) {
     394        double lowerValue=lower_[iSequence];
     395        double upperValue=upper_[iSequence];
     396        double value=solution_[iSequence];
     397        if(getStatus(iSequence)!=ClpSimplex::basic) {
     398          assert(lowerValue>-1.0e20);
     399          assert(upperValue<1.0e20);
     400        }
     401        switch(getStatus(iSequence)) {
     402         
     403        case ClpSimplex::basic:
     404          break;
     405        case ClpSimplex::isFree:
     406        case ClpSimplex::superBasic:
     407          break;
     408        case ClpSimplex::atUpperBound:
     409          assert (fabs(value-upperValue)<=primalTolerance_) ;
     410          break;
     411        case ClpSimplex::atLowerBound:
     412          assert (fabs(value-lowerValue)<=primalTolerance_) ;
     413          break;
     414        }
    369415      }
    370 #endif
     416    }
     417    if(numberIterations_==debugIteration) {
     418      printf("dodgy iteration coming up\n");
     419    }
     420#endif
     421    // choose row to go out
     422    dualRow();
     423    if (pivotRow_>=0) {
     424      // we found a pivot row
     425      handler_->message(CLP_SIMPLEX_PIVOTROW,messages_)
     426        <<pivotRow_
     427        <<OsiMessageEol;
     428      // check accuracy of weights
     429      dualRowPivot_->checkAccuracy();
     430      // get sign for finding row of tableau
     431      rowArray_[0]->insert(pivotRow_,directionOut_);
     432      factorization_->updateColumnTranspose(rowArray_[1],rowArray_[0]);
     433      // put row of tableau in rowArray[0] and columnArray[0]
     434      matrix_->transposeTimes(this,-1.0,
     435                              rowArray_[0],columnArray_[1],columnArray_[0]);
     436      // rowArray has pi equivalent
     437      // do ratio test
     438      dualColumn(rowArray_[0],columnArray_[0],columnArray_[1],
     439                 rowArray_[3]);
     440      if (sequenceIn_>=0) {
     441        // normal iteration
     442        // update the incoming column
     443        unpack(rowArray_[1]);
     444        factorization_->updateColumn(rowArray_[2],rowArray_[1],true);
     445        // and update dual weights (can do in parallel - with extra array)
     446        dualRowPivot_->updateWeights(rowArray_[0],rowArray_[2],
     447                                     rowArray_[1]);
     448        // see if update stable
     449        double btranAlpha = -alpha_*directionOut_; // for check
     450        alpha_=(*rowArray_[1])[pivotRow_];
    371451#ifdef CLP_DEBUG
    372       {
    373         int iSequence, number=numberRows_+numberColumns_;
    374         for (iSequence=0;iSequence<number;iSequence++) {
    375           double lowerValue=lower_[iSequence];
    376           double upperValue=upper_[iSequence];
    377           double value=solution_[iSequence];
    378           if(getStatus(iSequence)!=ClpSimplex::basic) {
    379             assert(lowerValue>-1.0e20);
    380             assert(upperValue<1.0e20);
    381           }
    382           switch(getStatus(iSequence)) {
    383            
    384           case ClpSimplex::basic:
     452        if ((handler_->logLevel()&32))
     453          printf("btran alpha %g, ftran alpha %g\n",btranAlpha,alpha_);
     454#endif
     455        if (fabs(btranAlpha)<1.0e-12||fabs(alpha_)<1.0e-12||
     456            fabs(btranAlpha-alpha_)>1.0e-7*(1.0+fabs(alpha_))) {
     457          handler_->message(CLP_DUAL_CHECK,messages_)
     458            <<btranAlpha
     459            <<alpha_
     460            <<OsiMessageEol;
     461          dualRowPivot_->unrollWeights();
     462          if (factorization_->pivots()) {
     463            problemStatus_=-2; // factorize now
     464            rowArray_[0]->clear();
     465            rowArray_[1]->clear();
     466            columnArray_[0]->clear();
    385467            break;
    386           case ClpSimplex::isFree:
    387           case ClpSimplex::superBasic:
    388             break;
    389           case ClpSimplex::atUpperBound:
    390             assert (fabs(value-upperValue)<=primalTolerance_) ;
    391             break;
    392           case ClpSimplex::atLowerBound:
    393             assert (fabs(value-lowerValue)<=primalTolerance_) ;
    394               break;
    395           }
    396         }
    397       }
    398       if(numberIterations_==debugIteration) {
    399         printf("dodgy iteration coming up\n");
    400       }
    401 #endif
    402       // choose row to go out
    403       dualRow();
    404       if (pivotRow_>=0) {
    405         // we found a pivot row
    406         handler_->message(CLP_SIMPLEX_PIVOTROW,messages_)
    407           <<pivotRow_
    408           <<OsiMessageEol;
    409         // check accuracy of weights
    410         dualRowPivot_->checkAccuracy();
    411         // get sign for finding row of tableau
    412         rowArray_[0]->insert(pivotRow_,directionOut_);
    413         factorization_->updateColumnTranspose(rowArray_[1],rowArray_[0]);
    414         // put row of tableau in rowArray[0] and columnArray[0]
    415         matrix_->transposeTimes(this,-1.0,
    416                                 rowArray_[0],columnArray_[1],columnArray_[0]);
    417         // rowArray has pi equivalent
    418         // do ratio test
    419         dualColumn(rowArray_[0],columnArray_[0],columnArray_[1],
    420                    rowArray_[3]);
    421         if (sequenceIn_>=0) {
    422           // normal iteration
    423           // update the incoming column
    424           unpack(rowArray_[1]);
    425           factorization_->updateColumn(rowArray_[2],rowArray_[1],true);
    426           // and update dual weights (can do in parallel - with extra array)
    427           dualRowPivot_->updateWeights(rowArray_[0],rowArray_[2],
    428                                        rowArray_[1]);
    429           // see if update stable
    430           double btranAlpha = -alpha_*directionOut_; // for check
    431           alpha_=(*rowArray_[1])[pivotRow_];
    432 #ifdef CLP_DEBUG
    433           if ((handler_->logLevel()&32))
    434             printf("btran alpha %g, ftran alpha %g\n",btranAlpha,alpha_);
    435 #endif
    436           if (fabs(btranAlpha)<1.0e-12||fabs(alpha_)<1.0e-12||
    437               fabs(btranAlpha-alpha_)>1.0e-7*(1.0+fabs(alpha_))) {
    438             handler_->message(CLP_DUAL_CHECK,messages_)
    439               <<btranAlpha
    440               <<alpha_
    441               <<OsiMessageEol;
    442             dualRowPivot_->unrollWeights();
    443             if (factorization_->pivots()) {
    444               problemStatus_=-2; // factorize now
    445               rowArray_[0]->clear();
    446               rowArray_[1]->clear();
    447               columnArray_[0]->clear();
    448               break;
    449             } else {
    450               // take on more relaxed criterion
    451               if (fabs(btranAlpha)<1.0e-12||fabs(alpha_)<1.0e-12||
    452                   fabs(btranAlpha-alpha_)>1.0e-4*(1.0+fabs(alpha_))) {
    453                 // need to reject something
    454                 char x = isColumn(sequenceOut_) ? 'C' :'R';
    455                 handler_->message(CLP_SIMPLEX_FLAG,messages_)
    456                   <<x<<sequenceWithin(sequenceOut_)
    457                   <<OsiMessageEol;
    458                 setFlagged(sequenceOut_);
    459                 lastBadIteration_ = numberIterations_; // say be more cautious
    460                 rowArray_[0]->clear();
    461                 rowArray_[1]->clear();
    462                 columnArray_[0]->clear();
    463                 continue;
    464               }
    465             }
    466           }
    467           // update duals BEFORE replaceColumn so can do updateColumn
    468           objectiveChange=0.0;
    469           // do duals first as variables may flip bounds
    470           // rowArray_[0] and columnArray_[0] may have flips
    471           // so use rowArray_[3] for work array from here on
    472           int nswapped =
    473             updateDualsInDual(rowArray_[0],columnArray_[0],rowArray_[2],theta_,
    474                         objectiveChange);
    475           // which will change basic solution
    476           if (nswapped) {
    477 #ifdef CLP_DEBUG
    478           if ((handler_->logLevel()&16))
    479             printf("old dualOut_ %g, v %g, l %g, u %g - new ",
    480                    dualOut_,valueOut_,lowerOut_,upperOut_);
    481             double oldOut=dualOut_;
    482 #endif
    483             factorization_->updateColumn(rowArray_[3],rowArray_[2],false);
    484             dualRowPivot_->updatePrimalSolution(rowArray_[2],
    485                                                 1.0,objectiveChange);
    486 
    487             // recompute dualOut_
    488             valueOut_ = solution_[sequenceOut_];
    489             if (directionOut_<0) {
    490               dualOut_ = valueOut_ - upperOut_;
    491             } else {
    492               dualOut_ = lowerOut_ - valueOut_;
    493             }
    494 #ifdef CLP_DEBUG
    495             if ((handler_->logLevel()&16))
    496               printf("%g\n",dualOut_);
    497             assert(dualOut_<=oldOut);
    498 #endif
    499             if(dualOut_<0.0&&factorization_->pivots()) {
    500               // going backwards - factorize
    501               dualRowPivot_->unrollWeights();
    502               problemStatus_=-2; // factorize now
    503               break;
    504             }
    505           }
    506           // amount primal will move
    507           double movement = -dualOut_*directionOut_/alpha_;
    508           // so objective should increase by fabs(dj)*movement
    509           // but we already have objective change - so check will be good
    510           if (objectiveChange+fabs(movement*dualIn_)<-1.0e-5) {
    511 #ifdef CLP_DEBUG
    512             if (handler_->logLevel()&32)
    513               printf("movement %g, swap change %g, rest %g  * %g\n",
    514                      objectiveChange+fabs(movement*dualIn_),
    515                      objectiveChange,movement,dualIn_);
    516 #endif
    517             if(factorization_->pivots()>5) {
    518               // going backwards - factorize
    519               dualRowPivot_->unrollWeights();
    520               problemStatus_=-2; // factorize now
    521               break;
    522             }
    523           }
    524           // if stable replace in basis
    525           int updateStatus = factorization_->replaceColumn(rowArray_[2],
    526                                                            pivotRow_,
    527                                                            alpha_);
    528           if (updateStatus==1) {
    529             // slight error
    530             if (factorization_->pivots()>5)
    531               problemStatus_=-2; // factorize now
    532           } else if (updateStatus==2) {
    533             // major error
    534             dualRowPivot_->unrollWeights();
    535             // later we may need to unwind more e.g. fake bounds
    536             if (factorization_->pivots()) {
    537               problemStatus_=-2; // factorize now
    538               break;
    539             } else {
     468          } else {
     469            // take on more relaxed criterion
     470            if (fabs(btranAlpha)<1.0e-12||fabs(alpha_)<1.0e-12||
     471                fabs(btranAlpha-alpha_)>1.0e-4*(1.0+fabs(alpha_))) {
    540472              // need to reject something
    541473              char x = isColumn(sequenceOut_) ? 'C' :'R';
     
    550482              continue;
    551483            }
    552           } else if (updateStatus==3) {
    553             // out of memory
    554             // increase space if not many iterations
    555             if (factorization_->pivots()<
    556                 0.5*factorization_->maximumPivots()&&
    557                 factorization_->pivots()<200)
    558               factorization_->areaFactor(
    559                                          factorization_->areaFactor() * 1.1);
     484          }
     485        }
     486        // update duals BEFORE replaceColumn so can do updateColumn
     487        double objectiveChange=0.0;
     488        // do duals first as variables may flip bounds
     489        // rowArray_[0] and columnArray_[0] may have flips
     490        // so use rowArray_[3] for work array from here on
     491        int nswapped =
     492          updateDualsInDual(rowArray_[0],columnArray_[0],rowArray_[2],theta_,
     493                            objectiveChange);
     494        // which will change basic solution
     495        if (nswapped) {
     496#ifdef CLP_DEBUG
     497          if ((handler_->logLevel()&16))
     498            printf("old dualOut_ %g, v %g, l %g, u %g - new ",
     499                   dualOut_,valueOut_,lowerOut_,upperOut_);
     500          double oldOut=dualOut_;
     501#endif
     502          factorization_->updateColumn(rowArray_[3],rowArray_[2],false);
     503          dualRowPivot_->updatePrimalSolution(rowArray_[2],
     504                                              1.0,objectiveChange);
     505         
     506          // recompute dualOut_
     507          valueOut_ = solution_[sequenceOut_];
     508          if (directionOut_<0) {
     509            dualOut_ = valueOut_ - upperOut_;
     510          } else {
     511            dualOut_ = lowerOut_ - valueOut_;
     512          }
     513#ifdef CLP_DEBUG
     514          if ((handler_->logLevel()&16))
     515            printf("%g\n",dualOut_);
     516          assert(dualOut_<=oldOut);
     517#endif
     518          if(dualOut_<0.0&&factorization_->pivots()) {
     519            // going backwards - factorize
     520            dualRowPivot_->unrollWeights();
    560521            problemStatus_=-2; // factorize now
    561           }
    562           // update primal solution
    563           if (theta_<0.0) {
    564 #ifdef CLP_DEBUG
    565             if (handler_->logLevel()&32)
    566               printf("negative theta %g\n",theta_);
    567 #endif
    568             theta_=0.0;
    569           }
    570           // do actual flips
    571           flipBounds(rowArray_[0],columnArray_[0],theta_);
    572           dualRowPivot_->updatePrimalSolution(rowArray_[1],
    573                                               movement,
    574                                               objectiveChange);
    575 #ifdef CLP_DEBUG
    576           double oldobj=objectiveValue_;
    577 #endif
    578           int whatNext=housekeeping(objectiveChange);
    579           // and set bounds correctly
    580           originalBound(sequenceIn_);
    581           changeBound(sequenceOut_);
    582 #ifdef CLP_DEBUG
    583           if (objectiveValue_<oldobj-1.0e-5&&(handler_->logLevel()&16))
    584             printf("obj backwards %g %g\n",objectiveValue_,oldobj);
    585 #endif
    586           if (whatNext==1) {
    587             problemStatus_ =-2; // refactorize
    588           } else if (whatNext==2) {
    589             // maximum iterations or equivalent
    590             problemStatus_= 3;
    591522            break;
    592523          }
    593         } else {
    594           // no incoming column is valid
     524        }
     525        // amount primal will move
     526        double movement = -dualOut_*directionOut_/alpha_;
     527        // so objective should increase by fabs(dj)*movement
     528        // but we already have objective change - so check will be good
     529        if (objectiveChange+fabs(movement*dualIn_)<-1.0e-5) {
    595530#ifdef CLP_DEBUG
    596531          if (handler_->logLevel()&32)
    597             printf("** no column pivot\n");
    598 #endif
    599           if (factorization_->pivots()<5) {
    600             problemStatus_=-4; //say looks infeasible
    601             // create ray anyway
    602             delete [] ray_;
    603             ray_ = new double [ numberRows_];
    604             CoinDisjointCopyN(rowArray_[0]->denseVector(),numberRows_,ray_);
    605           }
    606           rowArray_[0]->clear();
    607           columnArray_[0]->clear();
     532            printf("movement %g, swap change %g, rest %g  * %g\n",
     533                   objectiveChange+fabs(movement*dualIn_),
     534                   objectiveChange,movement,dualIn_);
     535#endif
     536          if(factorization_->pivots()>5) {
     537            // going backwards - factorize
     538            dualRowPivot_->unrollWeights();
     539            problemStatus_=-2; // factorize now
     540            break;
     541          }
     542        }
     543        // if stable replace in basis
     544        int updateStatus = factorization_->replaceColumn(rowArray_[2],
     545                                                         pivotRow_,
     546                                                         alpha_);
     547        if (updateStatus==1) {
     548          // slight error
     549          if (factorization_->pivots()>5)
     550            problemStatus_=-2; // factorize now
     551        } else if (updateStatus==2) {
     552          // major error
     553          dualRowPivot_->unrollWeights();
     554          // later we may need to unwind more e.g. fake bounds
     555          if (factorization_->pivots()) {
     556            problemStatus_=-2; // factorize now
     557            break;
     558          } else {
     559            // need to reject something
     560            char x = isColumn(sequenceOut_) ? 'C' :'R';
     561            handler_->message(CLP_SIMPLEX_FLAG,messages_)
     562              <<x<<sequenceWithin(sequenceOut_)
     563              <<OsiMessageEol;
     564            setFlagged(sequenceOut_);
     565            lastBadIteration_ = numberIterations_; // say be more cautious
     566            rowArray_[0]->clear();
     567            rowArray_[1]->clear();
     568            columnArray_[0]->clear();
     569            continue;
     570          }
     571        } else if (updateStatus==3) {
     572          // out of memory
     573          // increase space if not many iterations
     574          if (factorization_->pivots()<
     575              0.5*factorization_->maximumPivots()&&
     576              factorization_->pivots()<200)
     577            factorization_->areaFactor(
     578                                       factorization_->areaFactor() * 1.1);
     579          problemStatus_=-2; // factorize now
     580        }
     581        // update primal solution
     582        if (theta_<0.0) {
     583#ifdef CLP_DEBUG
     584          if (handler_->logLevel()&32)
     585            printf("negative theta %g\n",theta_);
     586#endif
     587          theta_=0.0;
     588        }
     589        // do actual flips
     590        flipBounds(rowArray_[0],columnArray_[0],theta_);
     591        dualRowPivot_->updatePrimalSolution(rowArray_[1],
     592                                            movement,
     593                                            objectiveChange);
     594#ifdef CLP_DEBUG
     595        double oldobj=objectiveValue_;
     596#endif
     597        int whatNext=housekeeping(objectiveChange);
     598        // and set bounds correctly
     599        originalBound(sequenceIn_);
     600        changeBound(sequenceOut_);
     601#ifdef CLP_DEBUG
     602        if (objectiveValue_<oldobj-1.0e-5&&(handler_->logLevel()&16))
     603          printf("obj backwards %g %g\n",objectiveValue_,oldobj);
     604#endif
     605        if (whatNext==1) {
     606          problemStatus_ =-2; // refactorize
     607        } else if (whatNext==2) {
     608          // maximum iterations or equivalent
     609          problemStatus_= 3;
    608610          break;
    609611        }
    610612      } else {
    611         // no pivot row
     613        // no incoming column is valid
    612614#ifdef CLP_DEBUG
    613615        if (handler_->logLevel()&32)
    614           printf("** no row pivot\n");
    615 #endif
    616         if (!factorization_->pivots()) {
    617           // may have crept through - so may be optimal
    618           //problemStatus_=-5; //say looks unbounded
    619           problemStatus_=0;
    620           // check any flagged variables
    621           int iRow;
    622           for (iRow=0;iRow<numberRows_;iRow++) {
    623             int iPivot=pivotVariable_[iRow];
    624             if (flagged(iPivot))
    625               break;
    626           }
    627           if (iRow<numberRows_) {
    628 #ifdef CLP_DEBUG
    629             std::cerr<<"Flagged variables at end - infeasible?"<<std::endl;
    630 #endif
    631             problemStatus_=-4; //say looks infeasible
    632             // create ray anyway
    633             delete [] ray_;
    634             ray_ = new double [ numberRows_];
    635             CoinDisjointCopyN(rowArray_[0]->denseVector(),numberRows_,ray_);
    636           }
    637         }
     616          printf("** no column pivot\n");
     617#endif
     618        if (factorization_->pivots()<5) {
     619          problemStatus_=-4; //say looks infeasible
     620          // create ray anyway
     621          delete [] ray_;
     622          ray_ = new double [ numberRows_];
     623          CoinDisjointCopyN(rowArray_[0]->denseVector(),numberRows_,ray_);
     624        }
     625        rowArray_[0]->clear();
     626        columnArray_[0]->clear();
    638627        break;
    639628      }
    640     }
    641   }
    642 
    643   // at present we are leaving factorization around
    644   // maybe we should empty it
    645   deleteRim();
    646   handler_->message(CLP_SIMPLEX_FINISHED+problemStatus_,messages_)
    647     <<objectiveValue()
    648     <<OsiMessageEol;
    649   // Restore any saved stuff
    650   perturbation_ = savePerturbation;
    651   factorization_->sparseThreshold(saveSparse);
    652   dualBound_ = saveDualBound_;
    653   return problemStatus_;
     629    } else {
     630      // no pivot row
     631#ifdef CLP_DEBUG
     632      if (handler_->logLevel()&32)
     633        printf("** no row pivot\n");
     634#endif
     635      if (!factorization_->pivots()) {
     636        // may have crept through - so may be optimal
     637        //problemStatus_=-5; //say looks unbounded
     638        problemStatus_=0;
     639        // check any flagged variables
     640        int iRow;
     641        for (iRow=0;iRow<numberRows_;iRow++) {
     642          int iPivot=pivotVariable_[iRow];
     643          if (flagged(iPivot))
     644            break;
     645        }
     646        if (iRow<numberRows_) {
     647#ifdef CLP_DEBUG
     648          std::cerr<<"Flagged variables at end - infeasible?"<<std::endl;
     649#endif
     650          problemStatus_=-4; //say looks infeasible
     651          // create ray anyway
     652          delete [] ray_;
     653          ray_ = new double [ numberRows_];
     654          CoinDisjointCopyN(rowArray_[0]->denseVector(),numberRows_,ray_);
     655        }
     656      }
     657      break;
     658    }
     659  }
    654660}
    655661/* The duals are updated by the given arrays.
    656     Returns number of infeasibilities.
    657     rowArray and columnarray will have flipped
    658     The output vector has movement (row length array) */
     662   Returns number of infeasibilities.
     663   rowArray and columnarray will have flipped
     664   The output vector has movement (row length array) */
    659665int
    660666ClpSimplexDual::updateDualsInDual(OsiIndexedVector * rowArray,
    661                                 OsiIndexedVector * columnArray,
    662                                 OsiIndexedVector * outputArray,
    663                                 double theta,
    664                                 double & objectiveChange)
     667                                  OsiIndexedVector * columnArray,
     668                                  OsiIndexedVector * outputArray,
     669                                  double theta,
     670                                  double & objectiveChange)
    665671{
    666 
     672 
    667673  outputArray->clear();
    668674 
     
    678684                       columnArray->getNumElements()==numberColumns_);
    679685  int numberAtFake=0;
    680 
     686 
    681687  // use a tighter tolerance except for all being okay
    682688  double tolerance = dualTolerance_;
    683 
     689 
    684690  double changeObj=0.0;
    685 
     691 
    686692  int iSection;
    687693 
Note: See TracChangeset for help on using the changeset viewer.