Ignore:
Timestamp:
Jan 6, 2019 2:43:06 PM (4 months ago)
Author:
unxusr
Message:

formatting

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Clp/src/AbcSimplexPrimal.cpp

    r2166 r2385  
    108108   Variable will change theta if currentValue - currentTheta*alpha < 0.0
    109109*/
    110 bool userChoiceValid1(const AbcSimplex * model,
    111                       int sequenceOut,
    112                       double currentValue,
    113                       double currentTheta,
    114                       double alpha,
    115                       double realAlpha);
     110bool userChoiceValid1(const AbcSimplex *model,
     111  int sequenceOut,
     112  double currentValue,
     113  double currentTheta,
     114  double alpha,
     115  double realAlpha);
    116116/* This returns true if chosen in/out pair valid.
    117117   The main thing to check would be variable flipping bounds may be
     
    119119   If you return false sequenceIn_ will be flagged as ineligible.
    120120*/
    121 bool userChoiceValid2(const AbcSimplex * model);
     121bool userChoiceValid2(const AbcSimplex *model);
    122122/* If a good pivot then you may wish to unflag some variables.
    123123 */
    124 void userChoiceWasGood(AbcSimplex * model);
     124void userChoiceWasGood(AbcSimplex *model);
    125125#endif
    126126#ifdef TRY_SPLIT_VALUES_PASS
    127 static double valuesChunk=-10.0;
    128 static double valuesRatio=0.1;
    129 static int valuesStop=-1;
    130 static int keepValuesPass=-1;
    131 static int doOrdinaryVariables=-1;
     127static double valuesChunk = -10.0;
     128static double valuesRatio = 0.1;
     129static int valuesStop = -1;
     130static int keepValuesPass = -1;
     131static int doOrdinaryVariables = -1;
    132132#endif
    133133// primal
    134 int AbcSimplexPrimal::primal (int ifValuesPass , int /*startFinishOptions*/)
     134int AbcSimplexPrimal::primal(int ifValuesPass, int /*startFinishOptions*/)
    135135{
    136  
     136
    137137  /*
    138138    Method
     
    214214   
    215215  */
    216  
     216
    217217  algorithm_ = +1;
    218218  moreSpecialOptions_ &= ~16; // clear check replaceColumn accuracy
    219 #if ABC_PARALLEL>0
     219#if ABC_PARALLEL > 0
    220220  if (parallelMode()) {
    221221    // extra regions
    222222    // for moment allow for ordered factorization
    223     int length=2*numberRows_+abcFactorization_->maximumPivots();
     223    int length = 2 * numberRows_ + abcFactorization_->maximumPivots();
    224224    for (int i = 0; i < 6; i++) {
    225225      delete rowArray_[i];
     
    236236  dualTolerance_ = dblParam_[ClpDualTolerance];
    237237  primalTolerance_ = dblParam_[ClpPrimalTolerance];
    238   currentDualTolerance_=dualTolerance_;
     238  currentDualTolerance_ = dualTolerance_;
    239239  //nextCleanNonBasicIteration_=baseIteration_+numberRows_;
    240240  // Save so can see if doing after dual
     
    243243  int initialNegDjs = -1;
    244244  // copy bounds to perturbation
    245   CoinAbcMemcpy(abcPerturbation_,abcLower_,numberTotal_);
    246   CoinAbcMemcpy(abcPerturbation_+numberTotal_,abcUpper_,numberTotal_);
     245  CoinAbcMemcpy(abcPerturbation_, abcLower_, numberTotal_);
     246  CoinAbcMemcpy(abcPerturbation_ + numberTotal_, abcUpper_, numberTotal_);
    247247#if 0
    248248  if (numberRows_>80000&&numberRows_<90000) {
     
    266266  statusOfProblemInPrimal(0);
    267267  /*if (!startup(ifValuesPass))*/ {
    268    
     268
    269269    // Set average theta
    270270    abcNonLinearCost_->setAverageTheta(1.0e3);
    271    
     271
    272272    // Say no pivot has occurred (for steepest edge and updates)
    273273    pivotRow_ = -2;
    274    
     274
    275275    // This says whether to restore things etc
    276276    int factorType = 0;
    277277    if (problemStatus_ < 0 && perturbation_ < 100 && !ifValuesPass) {
    278278      perturb(0);
    279       if (perturbation_!=100) {
    280         // Can't get here if values pass
    281         assert (!ifValuesPass);
    282         gutsOfPrimalSolution(3);
    283         if (handler_->logLevel() > 2) {
    284           handler_->message(CLP_SIMPLEX_STATUS, messages_)
    285             << numberIterations_ << objectiveValue();
    286           handler_->printing(sumPrimalInfeasibilities_ > 0.0)
    287             << sumPrimalInfeasibilities_ << numberPrimalInfeasibilities_;
    288           handler_->printing(sumDualInfeasibilities_ > 0.0)
    289             << sumDualInfeasibilities_ << numberDualInfeasibilities_;
    290           handler_->printing(numberDualInfeasibilitiesWithoutFree_
    291                              < numberDualInfeasibilities_)
    292             << numberDualInfeasibilitiesWithoutFree_;
    293           handler_->message() << CoinMessageEol;
    294         }
     279      if (perturbation_ != 100) {
     280        // Can't get here if values pass
     281        assert(!ifValuesPass);
     282        gutsOfPrimalSolution(3);
     283        if (handler_->logLevel() > 2) {
     284          handler_->message(CLP_SIMPLEX_STATUS, messages_)
     285            << numberIterations_ << objectiveValue();
     286          handler_->printing(sumPrimalInfeasibilities_ > 0.0)
     287            << sumPrimalInfeasibilities_ << numberPrimalInfeasibilities_;
     288          handler_->printing(sumDualInfeasibilities_ > 0.0)
     289            << sumDualInfeasibilities_ << numberDualInfeasibilities_;
     290          handler_->printing(numberDualInfeasibilitiesWithoutFree_
     291            < numberDualInfeasibilities_)
     292            << numberDualInfeasibilitiesWithoutFree_;
     293          handler_->message() << CoinMessageEol;
     294        }
    295295      }
    296296    }
     
    298298    abcProgress_.fillFromModel(this);
    299299    abcProgress_.startCheck();
    300     bool needNewWeights=false;
    301     double pivotTolerance=abcFactorization_->pivotTolerance();
     300    bool needNewWeights = false;
     301    double pivotTolerance = abcFactorization_->pivotTolerance();
    302302    /*
    303303      Status of problem:
     
    316316      // If getting nowhere - why not give it a kick
    317317#if 1
    318       if (perturbation_ < 101 && numberIterations_ ==6078/*> 2 * (numberRows_ + numberColumns_)*/ && (specialOptions_ & 4) == 0
    319           && initialStatus != 10) {
    320         perturb(1);
     318      if (perturbation_ < 101 && numberIterations_ == 6078 /*> 2 * (numberRows_ + numberColumns_)*/ && (specialOptions_ & 4) == 0
     319        && initialStatus != 10) {
     320        perturb(1);
    321321      }
    322322#endif
    323323      // If we have done no iterations - special
    324324      if (lastGoodIteration_ == numberIterations_ && factorType)
    325         factorType = 3;
    326       if(pivotTolerance<abcFactorization_->pivotTolerance()) {
    327         // force factorization
    328         pivotTolerance = abcFactorization_->pivotTolerance();
    329         factorType=5;
    330       }
    331      
     325        factorType = 3;
     326      if (pivotTolerance < abcFactorization_->pivotTolerance()) {
     327        // force factorization
     328        pivotTolerance = abcFactorization_->pivotTolerance();
     329        factorType = 5;
     330      }
     331
    332332      // may factorize, checks if problem finished
    333333      if (factorType)
    334         statusOfProblemInPrimal(factorType + (needNewWeights ? 10 : 0));
    335       needNewWeights=false;
    336       if (problemStatus_==10 && (moreSpecialOptions_ & 2048) != 0) {
    337         problemStatus_=0;
     334        statusOfProblemInPrimal(factorType + (needNewWeights ? 10 : 0));
     335      needNewWeights = false;
     336      if (problemStatus_ == 10 && (moreSpecialOptions_ & 2048) != 0) {
     337        problemStatus_ = 0;
    338338      }
    339339      if (initialStatus == 10) {
    340         initialStatus=-1;
    341         // cleanup phase
    342         if(initialIterations != numberIterations_) {
    343           if (numberDualInfeasibilities_ > 10000 && numberDualInfeasibilities_ > 10 * initialNegDjs) {
    344             // getting worse - try perturbing
    345             if (perturbation_ < 101 && (specialOptions_ & 4) == 0) {
    346               perturb(1);
    347               statusOfProblemInPrimal(factorType);
    348             }
    349           }
    350         } else {
    351           // save number of negative djs
    352           if (!numberPrimalInfeasibilities_)
    353             initialNegDjs = numberDualInfeasibilities_;
    354           // make sure weight won't be changed
    355           if (infeasibilityCost_ == 1.0e10)
    356             infeasibilityCost_ = 1.000001e10;
    357         }
    358       }
    359      
     340        initialStatus = -1;
     341        // cleanup phase
     342        if (initialIterations != numberIterations_) {
     343          if (numberDualInfeasibilities_ > 10000 && numberDualInfeasibilities_ > 10 * initialNegDjs) {
     344            // getting worse - try perturbing
     345            if (perturbation_ < 101 && (specialOptions_ & 4) == 0) {
     346              perturb(1);
     347              statusOfProblemInPrimal(factorType);
     348            }
     349          }
     350        } else {
     351          // save number of negative djs
     352          if (!numberPrimalInfeasibilities_)
     353            initialNegDjs = numberDualInfeasibilities_;
     354          // make sure weight won't be changed
     355          if (infeasibilityCost_ == 1.0e10)
     356            infeasibilityCost_ = 1.000001e10;
     357        }
     358      }
     359
    360360      // Say good factorization
    361361      factorType = 1;
    362      
     362
    363363      // Say no pivot has occurred (for steepest edge and updates)
    364364      pivotRow_ = -2;
    365      
     365
    366366      // exit if victory declared
    367367      if (problemStatus_ >= 0) {
    368368#ifdef CLP_USER_DRIVEN
    369         int status =
    370           eventHandler_->event(ClpEventHandler::endInPrimal);
    371         if (status>=0&&status<10) {
    372           // carry on
    373           problemStatus_=-1;
    374           if (status==0)
    375             continue; // re-factorize
    376         } else if (status>=10) {
    377           problemStatus_=status-10;
    378           break;
    379         } else {
    380           break;
    381         }
     369        int status = eventHandler_->event(ClpEventHandler::endInPrimal);
     370        if (status >= 0 && status < 10) {
     371          // carry on
     372          problemStatus_ = -1;
     373          if (status == 0)
     374            continue; // re-factorize
     375        } else if (status >= 10) {
     376          problemStatus_ = status - 10;
     377          break;
     378        } else {
     379          break;
     380        }
    382381#else
    383         break;
    384 #endif
    385       } 
    386      
     382        break;
     383#endif
     384      }
     385
    387386      // test for maximum iterations
    388387      if (hitMaximumIterations() || (ifValuesPass == 2 && firstFree_ < 0)) {
    389         problemStatus_ = 3;
    390         break;
    391       }
    392      
     388        problemStatus_ = 3;
     389        break;
     390      }
     391
    393392      if (firstFree_ < 0) {
    394         if (ifValuesPass) {
    395           // end of values pass
    396           ifValuesPass = 0;
    397           needNewWeights=true;
    398           int status = eventHandler_->event(ClpEventHandler::endOfValuesPass);
    399           if (status >= 0) {
    400             problemStatus_ = 5;
    401             secondaryStatus_ = ClpEventHandler::endOfValuesPass;
    402             break;
    403           }
    404           //#define FEB_TRY
     393        if (ifValuesPass) {
     394          // end of values pass
     395          ifValuesPass = 0;
     396          needNewWeights = true;
     397          int status = eventHandler_->event(ClpEventHandler::endOfValuesPass);
     398          if (status >= 0) {
     399            problemStatus_ = 5;
     400            secondaryStatus_ = ClpEventHandler::endOfValuesPass;
     401            break;
     402          }
     403          //#define FEB_TRY
    405404#if 1 //def FEB_TRY
    406           if (perturbation_ < 100
     405          if (perturbation_ < 100
    407406#ifdef TRY_SPLIT_VALUES_PASS
    408               &&valuesStop<0
    409 #endif
    410               )
    411             perturb(0);
    412 #endif
    413         }
    414       }
    415       if (abcNonLinearCost_->numberInfeasibilities() > 0 && !abcProgress_.initialWeight_ && !ifValuesPass && infeasibilityCost_ == 1.0e10 
     407            && valuesStop < 0
     408#endif
     409          )
     410            perturb(0);
     411#endif
     412        }
     413      }
     414      if (abcNonLinearCost_->numberInfeasibilities() > 0 && !abcProgress_.initialWeight_ && !ifValuesPass && infeasibilityCost_ == 1.0e10
    416415#ifdef TRY_SPLIT_VALUES_PASS
    417           && valuesStop<0
    418 #endif
    419           ) {
    420         // first time infeasible - start up weight computation
    421         // compute with original costs
    422         CoinAbcMemcpy(djSaved_,abcDj_,numberTotal_);
    423         CoinAbcMemcpy(abcCost_,costSaved_,numberTotal_);
    424         CoinAbcGatherFrom(abcCost_,costBasic_,abcPivotVariable_,numberRows_);
    425         gutsOfPrimalSolution(1);
    426         int numberSame = 0;
    427         int numberDifferent = 0;
    428         int numberZero = 0;
    429         int numberFreeSame = 0;
    430         int numberFreeDifferent = 0;
    431         int numberFreeZero = 0;
    432         int n = 0;
    433         for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
    434           if (getInternalStatus(iSequence) != basic && !flagged(iSequence)) {
    435             // not basic
    436             double distanceUp = abcUpper_[iSequence] - abcSolution_[iSequence];
    437             double distanceDown = abcSolution_[iSequence] - abcLower_[iSequence];
    438             double feasibleDj = abcDj_[iSequence];
    439             double infeasibleDj = djSaved_[iSequence] - feasibleDj;
    440             double value = feasibleDj * infeasibleDj;
    441             if (distanceUp > primalTolerance_) {
    442               // Check if "free"
    443               if (distanceDown > primalTolerance_) {
    444                 // free
    445                 if (value > dualTolerance_) {
    446                   numberFreeSame++;
    447                 } else if(value < -dualTolerance_) {
    448                   numberFreeDifferent++;
    449                   abcDj_[n++] = feasibleDj / infeasibleDj;
    450                 } else {
    451                   numberFreeZero++;
    452                 }
    453               } else {
    454                 // should not be negative
    455                 if (value > dualTolerance_) {
    456                   numberSame++;
    457                 } else if(value < -dualTolerance_) {
    458                   numberDifferent++;
    459                   abcDj_[n++] = feasibleDj / infeasibleDj;
    460                 } else {
    461                   numberZero++;
    462                 }
    463               }
    464             } else if (distanceDown > primalTolerance_) {
    465               // should not be positive
    466               if (value > dualTolerance_) {
    467                 numberSame++;
    468               } else if(value < -dualTolerance_) {
    469                 numberDifferent++;
    470                 abcDj_[n++] = feasibleDj / infeasibleDj;
    471               } else {
    472                 numberZero++;
    473               }
    474             }
    475       }
    476           abcProgress_.initialWeight_ = -1.0;
    477     }
     416        && valuesStop < 0
     417#endif
     418      ) {
     419        // first time infeasible - start up weight computation
     420        // compute with original costs
     421        CoinAbcMemcpy(djSaved_, abcDj_, numberTotal_);
     422        CoinAbcMemcpy(abcCost_, costSaved_, numberTotal_);
     423        CoinAbcGatherFrom(abcCost_, costBasic_, abcPivotVariable_, numberRows_);
     424        gutsOfPrimalSolution(1);
     425        int numberSame = 0;
     426        int numberDifferent = 0;
     427        int numberZero = 0;
     428        int numberFreeSame = 0;
     429        int numberFreeDifferent = 0;
     430        int numberFreeZero = 0;
     431        int n = 0;
     432        for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
     433          if (getInternalStatus(iSequence) != basic && !flagged(iSequence)) {
     434            // not basic
     435            double distanceUp = abcUpper_[iSequence] - abcSolution_[iSequence];
     436            double distanceDown = abcSolution_[iSequence] - abcLower_[iSequence];
     437            double feasibleDj = abcDj_[iSequence];
     438            double infeasibleDj = djSaved_[iSequence] - feasibleDj;
     439            double value = feasibleDj * infeasibleDj;
     440            if (distanceUp > primalTolerance_) {
     441              // Check if "free"
     442              if (distanceDown > primalTolerance_) {
     443                // free
     444                if (value > dualTolerance_) {
     445                  numberFreeSame++;
     446                } else if (value < -dualTolerance_) {
     447                  numberFreeDifferent++;
     448                  abcDj_[n++] = feasibleDj / infeasibleDj;
     449                } else {
     450                  numberFreeZero++;
     451                }
     452              } else {
     453                // should not be negative
     454                if (value > dualTolerance_) {
     455                  numberSame++;
     456                } else if (value < -dualTolerance_) {
     457                  numberDifferent++;
     458                  abcDj_[n++] = feasibleDj / infeasibleDj;
     459                } else {
     460                  numberZero++;
     461                }
     462              }
     463            } else if (distanceDown > primalTolerance_) {
     464              // should not be positive
     465              if (value > dualTolerance_) {
     466                numberSame++;
     467              } else if (value < -dualTolerance_) {
     468                numberDifferent++;
     469                abcDj_[n++] = feasibleDj / infeasibleDj;
     470              } else {
     471                numberZero++;
     472              }
     473            }
     474          }
     475          abcProgress_.initialWeight_ = -1.0;
     476        }
    478477#ifdef CLP_USEFUL_PRINTOUT
    479         printf("XXXX %d same, %d different, %d zero, -- free %d %d %d\n",
    480                numberSame,numberDifferent,numberZero,
    481                numberFreeSame,numberFreeDifferent,numberFreeZero);
    482 #endif
    483         // we want most to be same
    484         if (n) {
    485           std::sort(abcDj_, abcDj_ + n);
     478        printf("XXXX %d same, %d different, %d zero, -- free %d %d %d\n",
     479          numberSame, numberDifferent, numberZero,
     480          numberFreeSame, numberFreeDifferent, numberFreeZero);
     481#endif
     482        // we want most to be same
     483        if (n) {
     484          std::sort(abcDj_, abcDj_ + n);
    486485#ifdef CLP_USEFUL_PRINTOUT
    487           double most = 0.95;
    488           int which = static_cast<int> ((1.0 - most) * static_cast<double> (n));
    489           double take2 = -abcDj_[which] * infeasibilityCost_;
    490           printf("XXXX inf cost %g take %g (range %g %g)\n",infeasibilityCost_,take2,-abcDj_[0]*infeasibilityCost_,-abcDj_[n-1]*infeasibilityCost_);
    491 #endif
    492           double take = -abcDj_[0] * infeasibilityCost_;
    493           // was infeasibilityCost_ = CoinMin(CoinMax(1000.0 * take, 1.0e8), 1.0000001e10);
    494           infeasibilityCost_ = CoinMin(CoinMax(1000.0 * take, 1.0e3), 1.0000001e10);
     486          double most = 0.95;
     487          int which = static_cast< int >((1.0 - most) * static_cast< double >(n));
     488          double take2 = -abcDj_[which] * infeasibilityCost_;
     489          printf("XXXX inf cost %g take %g (range %g %g)\n", infeasibilityCost_, take2, -abcDj_[0] * infeasibilityCost_, -abcDj_[n - 1] * infeasibilityCost_);
     490#endif
     491          double take = -abcDj_[0] * infeasibilityCost_;
     492          // was infeasibilityCost_ = CoinMin(CoinMax(1000.0 * take, 1.0e8), 1.0000001e10);
     493          infeasibilityCost_ = CoinMin(CoinMax(1000.0 * take, 1.0e3), 1.0000001e10);
    495494#ifdef CLP_USEFUL_PRINTOUT
    496           printf("XXXX changing weight to %g\n",infeasibilityCost_);
    497 #endif
    498           // may need to force redo of weights
    499           needNewWeights=true;
    500         }
    501         abcNonLinearCost_->checkInfeasibilities(0.0);
    502         gutsOfPrimalSolution(3);
     495          printf("XXXX changing weight to %g\n", infeasibilityCost_);
     496#endif
     497          // may need to force redo of weights
     498          needNewWeights = true;
     499        }
     500        abcNonLinearCost_->checkInfeasibilities(0.0);
     501        gutsOfPrimalSolution(3);
    503502      }
    504503      // Check event
    505504      {
    506         int status = eventHandler_->event(ClpEventHandler::endOfFactorization);
    507         if (status >= 0) {
    508           problemStatus_ = 5;
    509           secondaryStatus_ = ClpEventHandler::endOfFactorization;
    510           break;
    511         }
     505        int status = eventHandler_->event(ClpEventHandler::endOfFactorization);
     506        if (status >= 0) {
     507          problemStatus_ = 5;
     508          secondaryStatus_ = ClpEventHandler::endOfFactorization;
     509          break;
     510        }
    512511      }
    513512      // Iterate
     
    532531  unflag();
    533532  abcProgress_.clearTimesFlagged();
    534 #if ABC_PARALLEL>0
     533#if ABC_PARALLEL > 0
    535534  if (parallelMode()) {
    536535    for (int i = 0; i < 6; i++) {
     
    544543  //finish(startFinishOptions);
    545544  restoreData(data);
    546   setObjectiveValue(abcNonLinearCost_->feasibleReportCost()+objectiveOffset_);
     545  setObjectiveValue(abcNonLinearCost_->feasibleReportCost() + objectiveOffset_);
    547546  // clean up
    548   if (problemStatus_==10)
     547  if (problemStatus_ == 10)
    549548    abcNonLinearCost_->checkInfeasibilities(0.0);
    550549  delete abcNonLinearCost_;
    551   abcNonLinearCost_=NULL;
     550  abcNonLinearCost_ = NULL;
    552551#if 0
    553552  if (numberRows_>80000&&numberRows_<90000) {
     
    575574  +3 max iterations
    576575*/
    577 int
    578 AbcSimplexPrimal::whileIterating(int valuesOption)
     576int AbcSimplexPrimal::whileIterating(int valuesOption)
    579577{
    580578#if 1
    581579#define DELAYED_UPDATE
    582   arrayForBtran_=0;
     580  arrayForBtran_ = 0;
    583581  //setUsedArray(arrayForBtran_);
    584   arrayForFtran_=1;
     582  arrayForFtran_ = 1;
    585583  setUsedArray(arrayForFtran_);
    586   arrayForFlipBounds_=2;
     584  arrayForFlipBounds_ = 2;
    587585  setUsedArray(arrayForFlipBounds_);
    588   arrayForTableauRow_=3;
     586  arrayForTableauRow_ = 3;
    589587  setUsedArray(arrayForTableauRow_);
    590588  //arrayForDualColumn_=4;
    591589  //setUsedArray(arrayForDualColumn_);
    592 #if ABC_PARALLEL<2
    593   arrayForReplaceColumn_=4; //4;
     590#if ABC_PARALLEL < 2
     591  arrayForReplaceColumn_ = 4; //4;
    594592#else
    595   arrayForReplaceColumn_=6; //4;
     593  arrayForReplaceColumn_ = 6; //4;
    596594  setUsedArray(arrayForReplaceColumn_);
    597595#endif
     
    605603  if (valuesOption > 1)
    606604    superBasicType = 3;
    607   int numberStartingInfeasibilities=abcNonLinearCost_->numberInfeasibilities();
     605  int numberStartingInfeasibilities = abcNonLinearCost_->numberInfeasibilities();
    608606  // status stays at -1 while iterating, >=0 finished, -2 to invert
    609607  // status -3 to go to top without an invert
    610   int numberFlaggedStart =abcProgress_.timesFlagged();
     608  int numberFlaggedStart = abcProgress_.timesFlagged();
    611609  while (problemStatus_ == -1) {
    612610    if (!ifValuesPass) {
     
    616614      // NOTE rowArray_[0] is used by computeDuals which is a
    617615      // slow way of getting duals but might be used
    618       int saveSequence=sequenceIn_;
     616      int saveSequence = sequenceIn_;
    619617      primalColumn(&usefulArray_[arrayForFtran_], &usefulArray_[arrayForTableauRow_],
    620                    &usefulArray_[arrayForFlipBounds_]);
    621       if (saveSequence>=0&&saveSequence!=sequenceOut_) {
    622         if (getInternalStatus(saveSequence)==basic)
    623           abcDj_[saveSequence]=0.0;
    624       }
    625 #if ABC_NORMAL_DEBUG>0
    626       if (handler_->logLevel()==63) {
    627         for (int i=0;i<numberTotal_;i++) {
    628           if (fabs(abcDj_[i])>1.0e2)
    629             assert (getInternalStatus(i)!=basic);
    630         }
    631         double largestCost=0.0;
    632         double largestDj=0.0;
    633         double largestGoodDj=0.0;
    634         int iLargest=-1;
    635         int numberInfeasibilities=0;
    636         double sum=0.0;
    637         for (int i=0;i<numberTotal_;i++) {
    638           if (getInternalStatus(i)==isFixed)
    639             continue;
    640           largestCost=CoinMax(largestCost,fabs(abcCost_[i]));
    641           double dj=abcDj_[i];
    642           if (getInternalStatus(i)==atLowerBound)
    643             dj=-dj;
    644           largestDj=CoinMax(largestDj,fabs(dj));
    645           if(largestGoodDj<dj) {
    646             largestGoodDj=dj;
    647             iLargest=i;
    648           }
    649           if (getInternalStatus(i)==basic) {
    650             if (abcSolution_[i]>abcUpper_[i]+primalTolerance_) {
    651               numberInfeasibilities++;
    652               sum += abcSolution_[i]-abcUpper_[i]-primalTolerance_;
    653             } else if (abcSolution_[i]<abcLower_[i]-primalTolerance_) {
    654               numberInfeasibilities++;
    655               sum -= abcSolution_[i]-abcLower_[i]+primalTolerance_;
    656             }
    657           }
    658         }
    659         char xx;
    660         int kLargest;
    661         if (iLargest>=numberRows_) {
    662           xx='C';
    663           kLargest=iLargest-numberRows_;
    664         } else {
    665           xx='R';
    666           kLargest=iLargest;
    667         }
    668         printf("largest cost %g, largest dj %g best %g (%d==%c%d) - %d infeasible (sum %g) nonlininf %d\n",
    669                largestCost,largestDj,largestGoodDj,iLargest,xx,kLargest,numberInfeasibilities,sum,
    670                abcNonLinearCost_->numberInfeasibilities());
    671         assert (getInternalStatus(iLargest)!=basic);
     618        &usefulArray_[arrayForFlipBounds_]);
     619      if (saveSequence >= 0 && saveSequence != sequenceOut_) {
     620        if (getInternalStatus(saveSequence) == basic)
     621          abcDj_[saveSequence] = 0.0;
     622      }
     623#if ABC_NORMAL_DEBUG > 0
     624      if (handler_->logLevel() == 63) {
     625        for (int i = 0; i < numberTotal_; i++) {
     626          if (fabs(abcDj_[i]) > 1.0e2)
     627            assert(getInternalStatus(i) != basic);
     628        }
     629        double largestCost = 0.0;
     630        double largestDj = 0.0;
     631        double largestGoodDj = 0.0;
     632        int iLargest = -1;
     633        int numberInfeasibilities = 0;
     634        double sum = 0.0;
     635        for (int i = 0; i < numberTotal_; i++) {
     636          if (getInternalStatus(i) == isFixed)
     637            continue;
     638          largestCost = CoinMax(largestCost, fabs(abcCost_[i]));
     639          double dj = abcDj_[i];
     640          if (getInternalStatus(i) == atLowerBound)
     641            dj = -dj;
     642          largestDj = CoinMax(largestDj, fabs(dj));
     643          if (largestGoodDj < dj) {
     644            largestGoodDj = dj;
     645            iLargest = i;
     646          }
     647          if (getInternalStatus(i) == basic) {
     648            if (abcSolution_[i] > abcUpper_[i] + primalTolerance_) {
     649              numberInfeasibilities++;
     650              sum += abcSolution_[i] - abcUpper_[i] - primalTolerance_;
     651            } else if (abcSolution_[i] < abcLower_[i] - primalTolerance_) {
     652              numberInfeasibilities++;
     653              sum -= abcSolution_[i] - abcLower_[i] + primalTolerance_;
     654            }
     655          }
     656        }
     657        char xx;
     658        int kLargest;
     659        if (iLargest >= numberRows_) {
     660          xx = 'C';
     661          kLargest = iLargest - numberRows_;
     662        } else {
     663          xx = 'R';
     664          kLargest = iLargest;
     665        }
     666        printf("largest cost %g, largest dj %g best %g (%d==%c%d) - %d infeasible (sum %g) nonlininf %d\n",
     667          largestCost, largestDj, largestGoodDj, iLargest, xx, kLargest, numberInfeasibilities, sum,
     668          abcNonLinearCost_->numberInfeasibilities());
     669        assert(getInternalStatus(iLargest) != basic);
    672670      }
    673671#endif
    674672    } else {
    675673      // in values pass
    676       for (int i=0;i<4;i++)
    677         multipleSequenceIn_[i]=-1;
     674      for (int i = 0; i < 4; i++)
     675        multipleSequenceIn_[i] = -1;
    678676      int sequenceIn = nextSuperBasic(superBasicType, &usefulArray_[arrayForFlipBounds_]);
    679677      if (valuesOption > 1)
    680         superBasicType = 2;
     678        superBasicType = 2;
    681679      if (sequenceIn < 0) {
    682         // end of values pass - initialize weights etc
    683         sequenceIn_=-1;
    684         handler_->message(CLP_END_VALUES_PASS, messages_)
    685           << numberIterations_;
    686         stateOfProblem_ &= ~VALUES_PASS;
     680        // end of values pass - initialize weights etc
     681        sequenceIn_ = -1;
     682        handler_->message(CLP_END_VALUES_PASS, messages_)
     683          << numberIterations_;
     684        stateOfProblem_ &= ~VALUES_PASS;
    687685#ifdef TRY_SPLIT_VALUES_PASS
    688         valuesStop=numberIterations_+doOrdinaryVariables;
    689 #endif
    690         abcPrimalColumnPivot_->saveWeights(this, 5);
    691         problemStatus_ = -2; // factorize now
    692         pivotRow_ = -1; // say no weights update
    693         returnCode = -4;
    694         // Clean up
    695         for (int i = 0; i < numberTotal_; i++) {
    696           if (getInternalStatus(i) == atLowerBound || getInternalStatus(i) == isFixed)
    697             abcSolution_[i] = abcLower_[i];
    698           else if (getInternalStatus(i) == atUpperBound)
    699             abcSolution_[i] = abcUpper_[i];
    700         }
    701         break;
     686        valuesStop = numberIterations_ + doOrdinaryVariables;
     687#endif
     688        abcPrimalColumnPivot_->saveWeights(this, 5);
     689        problemStatus_ = -2; // factorize now
     690        pivotRow_ = -1; // say no weights update
     691        returnCode = -4;
     692        // Clean up
     693        for (int i = 0; i < numberTotal_; i++) {
     694          if (getInternalStatus(i) == atLowerBound || getInternalStatus(i) == isFixed)
     695            abcSolution_[i] = abcLower_[i];
     696          else if (getInternalStatus(i) == atUpperBound)
     697            abcSolution_[i] = abcUpper_[i];
     698        }
     699        break;
    702700      } else {
    703         // normal
    704         sequenceIn_ = sequenceIn;
    705         valueIn_ = abcSolution_[sequenceIn_];
    706         lowerIn_ = abcLower_[sequenceIn_];
    707         upperIn_ = abcUpper_[sequenceIn_];
    708         dualIn_ = abcDj_[sequenceIn_];
    709         // see if any more
    710         if (maximumIterations()==100000) {
    711           multipleSequenceIn_[0]=sequenceIn;
    712           for (int i=1;i<4;i++) {
    713             int sequenceIn = nextSuperBasic(superBasicType, &usefulArray_[arrayForFlipBounds_]);
    714             if (sequenceIn>=0)
    715               multipleSequenceIn_[i]=sequenceIn;
    716             else
    717               break;
    718           }
    719         }
     701        // normal
     702        sequenceIn_ = sequenceIn;
     703        valueIn_ = abcSolution_[sequenceIn_];
     704        lowerIn_ = abcLower_[sequenceIn_];
     705        upperIn_ = abcUpper_[sequenceIn_];
     706        dualIn_ = abcDj_[sequenceIn_];
     707        // see if any more
     708        if (maximumIterations() == 100000) {
     709          multipleSequenceIn_[0] = sequenceIn;
     710          for (int i = 1; i < 4; i++) {
     711            int sequenceIn = nextSuperBasic(superBasicType, &usefulArray_[arrayForFlipBounds_]);
     712            if (sequenceIn >= 0)
     713              multipleSequenceIn_[i] = sequenceIn;
     714            else
     715              break;
     716          }
     717        }
    720718      }
    721719    }
     
    725723    if (sequenceIn_ >= 0) {
    726724      // we found a pivot column
    727       assert (!flagged(sequenceIn_));
     725      assert(!flagged(sequenceIn_));
    728726      //#define MULTIPLE_PRICE
    729727      // do second half of iteration
    730       if (multipleSequenceIn_[1]==-1||maximumIterations()!=100000) {
    731         returnCode = pivotResult(ifValuesPass);
     728      if (multipleSequenceIn_[1] == -1 || maximumIterations() != 100000) {
     729        returnCode = pivotResult(ifValuesPass);
    732730      } else {
    733         if (multipleSequenceIn_[0]<0)
    734           multipleSequenceIn_[0]=sequenceIn_;
    735         returnCode = pivotResult4(ifValuesPass);
     731        if (multipleSequenceIn_[0] < 0)
     732          multipleSequenceIn_[0] = sequenceIn_;
     733        returnCode = pivotResult4(ifValuesPass);
    736734#ifdef MULTIPLE_PRICE
    737         if (sequenceIn_>=0)
    738           returnCode = pivotResult(ifValuesPass);
    739 #endif
    740       }
    741       if(numberStartingInfeasibilities&&!abcNonLinearCost_->numberInfeasibilities()) {
    742         //if (abcFactorization_->pivots()>200&&numberIterations_>2*(numberRows_+numberColumns_))
    743         if (abcFactorization_->pivots()>2&&numberIterations_>(numberRows_+numberColumns_)&&(stateOfProblem_&PESSIMISTIC)!=0)
    744           returnCode=-2; // refactorize - maybe just after n iterations
     735        if (sequenceIn_ >= 0)
     736          returnCode = pivotResult(ifValuesPass);
     737#endif
     738      }
     739      if (numberStartingInfeasibilities && !abcNonLinearCost_->numberInfeasibilities()) {
     740        //if (abcFactorization_->pivots()>200&&numberIterations_>2*(numberRows_+numberColumns_))
     741        if (abcFactorization_->pivots() > 2 && numberIterations_ > (numberRows_ + numberColumns_) && (stateOfProblem_ & PESSIMISTIC) != 0)
     742          returnCode = -2; // refactorize - maybe just after n iterations
    745743      }
    746744      if (returnCode < -1 && returnCode > -5) {
    747         problemStatus_ = -2; //
     745        problemStatus_ = -2; //
    748746      } else if (returnCode == -5) {
    749         if (abcProgress_.timesFlagged()>10+numberFlaggedStart)
    750           problemStatus_ =-2;
    751         if ((moreSpecialOptions_ & 16) == 0 && abcFactorization_->pivots()) {
    752           moreSpecialOptions_ |= 16;
    753           problemStatus_ = -2;
    754         }
    755         // otherwise something flagged - continue;
     747        if (abcProgress_.timesFlagged() > 10 + numberFlaggedStart)
     748          problemStatus_ = -2;
     749        if ((moreSpecialOptions_ & 16) == 0 && abcFactorization_->pivots()) {
     750          moreSpecialOptions_ |= 16;
     751          problemStatus_ = -2;
     752        }
     753        // otherwise something flagged - continue;
    756754      } else if (returnCode == 2) {
    757         problemStatus_ = -5; // looks unbounded
     755        problemStatus_ = -5; // looks unbounded
    758756      } else if (returnCode == 4) {
    759         problemStatus_ = -2; // looks unbounded but has iterated
     757        problemStatus_ = -2; // looks unbounded but has iterated
    760758      } else if (returnCode != -1) {
    761         assert(returnCode == 3);
    762         if (problemStatus_ != 5)
    763           problemStatus_ = 3;
     759        assert(returnCode == 3);
     760        if (problemStatus_ != 5)
     761          problemStatus_ = 3;
    764762      }
    765763    } else {
    766764      // no pivot column
    767 #if ABC_NORMAL_DEBUG>3
     765#if ABC_NORMAL_DEBUG > 3
    768766      if (handler_->logLevel() & 32)
    769         printf("** no column pivot\n");
     767        printf("** no column pivot\n");
    770768#endif
    771769      if (abcNonLinearCost_->numberInfeasibilities())
    772         problemStatus_ = -4; // might be infeasible
     770        problemStatus_ = -4; // might be infeasible
    773771      // Force to re-factorize early next time
    774772      int numberPivots = abcFactorization_->pivots();
     
    776774#ifdef CLP_USER_DRIVEN
    777775      // If large number of pivots trap later?
    778       if (problemStatus_==-1 && numberPivots<13000) {
    779         int status = eventHandler_->event(ClpEventHandler::noCandidateInPrimal);
    780         if (status>=0&&status<10) {
    781           // carry on
    782           problemStatus_=-1;
    783           if (status==0)
    784             break;
    785         } else if (status>=10) {
    786           problemStatus_=status-10;
    787           break;
    788         } else {
    789           forceFactorization_ = CoinMin(forceFactorization_, (numberPivots + 1) >> 1);
    790           break;
    791         }
     776      if (problemStatus_ == -1 && numberPivots < 13000) {
     777        int status = eventHandler_->event(ClpEventHandler::noCandidateInPrimal);
     778        if (status >= 0 && status < 10) {
     779          // carry on
     780          problemStatus_ = -1;
     781          if (status == 0)
     782            break;
     783        } else if (status >= 10) {
     784          problemStatus_ = status - 10;
     785          break;
     786        } else {
     787          forceFactorization_ = CoinMin(forceFactorization_, (numberPivots + 1) >> 1);
     788          break;
     789        }
    792790      }
    793791#else
     
    802800}
    803801/* Checks if finished.  Updates status */
    804 void
    805 AbcSimplexPrimal::statusOfProblemInPrimal(int type)
     802void AbcSimplexPrimal::statusOfProblemInPrimal(int type)
    806803{
    807804  int saveFirstFree = firstFree_;
     
    813810#endif
    814811  //int saveType=type;
    815   if (type>3) {
     812  if (type > 3) {
    816813    // force factorization
    817     numberPivots=9999999;
    818     if (type>10)
    819       type -=10;
    820     else 
     814    numberPivots = 9999999;
     815    if (type > 10)
     816      type -= 10;
     817    else
    821818      type &= 3;
    822819  }
    823820#ifndef TRY_ABC_GUS
    824   bool doFactorization =  (type != 3&&(numberPivots>dontFactorizePivots_||numberIterations_==baseIteration_||problemStatus_==10));
     821  bool doFactorization = (type != 3 && (numberPivots > dontFactorizePivots_ || numberIterations_ == baseIteration_ || problemStatus_ == 10));
    825822#else
    826   bool doFactorization =  (type != 3&&(numberPivots>dontFactorizePivots_||numberIterations_==baseIteration_));
    827 #endif
    828   bool doWeights = doFactorization||problemStatus_==10;
     823  bool doFactorization = (type != 3 && (numberPivots > dontFactorizePivots_ || numberIterations_ == baseIteration_));
     824#endif
     825  bool doWeights = doFactorization || problemStatus_ == 10;
    829826  if (type == 2) {
    830827    // trouble - restore solution
     
    841838#endif
    842839#if CLP_CAUTION
    843   double lastAverageInfeasibility = sumDualInfeasibilities_ /
    844     static_cast<double>(numberDualInfeasibilities_ + 1);
    845 #endif
    846   if (numberIterations_&&type) {
     840  double lastAverageInfeasibility = sumDualInfeasibilities_ / static_cast< double >(numberDualInfeasibilities_ + 1);
     841#endif
     842  if (numberIterations_ && type) {
    847843    lastSumInfeasibility = abcNonLinearCost_->sumInfeasibilities();
    848844    lastNumberInfeasibility = abcNonLinearCost_->numberInfeasibilities();
    849845  } else {
    850     lastAverageInfeasibility=1.0e10;
    851   }
    852   bool ifValuesPass=(stateOfProblem_&VALUES_PASS)!=0;
    853   bool takenAction=false;
    854   double sumInfeasibility=0.0;
     846    lastAverageInfeasibility = 1.0e10;
     847  }
     848  bool ifValuesPass = (stateOfProblem_ & VALUES_PASS) != 0;
     849  bool takenAction = false;
     850  double sumInfeasibility = 0.0;
    855851  if (problemStatus_ > -3 || problemStatus_ == -4) {
    856852    // factorize
     
    861857    if (doFactorization)
    862858      abcPrimalColumnPivot_->saveWeights(this, 1);
    863    
     859
    864860    if (!type) {
    865861      // be optimistic
     
    871867    if (doFactorization) {
    872868      // is factorization okay?
    873       int solveType = ifValuesPass ? 11 :1;
     869      int solveType = ifValuesPass ? 11 : 1;
    874870      if (!type)
    875         solveType++;
     871        solveType++;
    876872      int factorStatus = internalFactorize(solveType);
    877873      if (factorStatus) {
    878         if (type != 1 || largestPrimalError_ > 1.0e3
    879             || largestDualError_ > 1.0e3) {
    880 #if ABC_NORMAL_DEBUG>0
    881           printf("Bad initial basis\n");
    882 #endif
    883           internalFactorize(2);
    884         } else {
    885           // no - restore previous basis
    886           // Keep any flagged variables
    887           for (int i = 0; i < numberTotal_; i++) {
    888             if (flagged(i))
    889               internalStatusSaved_[i] |= 64; //say flagged
    890           }
    891           restoreGoodStatus(1);
    892           if (numberPivots <= 1) {
    893             // throw out something
    894             if (sequenceIn_ >= 0 && getInternalStatus(sequenceIn_) != basic) {
    895               setFlagged(sequenceIn_);
    896             } else if (sequenceOut_ >= 0 && getInternalStatus(sequenceOut_) != basic) {
    897               setFlagged(sequenceOut_);
    898             }
    899             abcProgress_.incrementTimesFlagged();
    900             double newTolerance = CoinMax(0.5 + 0.499 * randomNumberGenerator_.randomDouble(), abcFactorization_->pivotTolerance());
    901             abcFactorization_->pivotTolerance(newTolerance);
    902           } else {
    903             // Go to safe
    904             abcFactorization_->pivotTolerance(0.99);
    905           }
    906           forceFactorization_ = 1; // a bit drastic but ..
    907           type = 2;
    908           abcNonLinearCost_->checkInfeasibilities();
    909           if (internalFactorize(2) != 0) {
    910             largestPrimalError_ = 1.0e4; // force other type
    911           }
    912         }
    913         changeMade_++; // say change made
     874        if (type != 1 || largestPrimalError_ > 1.0e3
     875          || largestDualError_ > 1.0e3) {
     876#if ABC_NORMAL_DEBUG > 0
     877          printf("Bad initial basis\n");
     878#endif
     879          internalFactorize(2);
     880        } else {
     881          // no - restore previous basis
     882          // Keep any flagged variables
     883          for (int i = 0; i < numberTotal_; i++) {
     884            if (flagged(i))
     885              internalStatusSaved_[i] |= 64; //say flagged
     886          }
     887          restoreGoodStatus(1);
     888          if (numberPivots <= 1) {
     889            // throw out something
     890            if (sequenceIn_ >= 0 && getInternalStatus(sequenceIn_) != basic) {
     891              setFlagged(sequenceIn_);
     892            } else if (sequenceOut_ >= 0 && getInternalStatus(sequenceOut_) != basic) {
     893              setFlagged(sequenceOut_);
     894            }
     895            abcProgress_.incrementTimesFlagged();
     896            double newTolerance = CoinMax(0.5 + 0.499 * randomNumberGenerator_.randomDouble(), abcFactorization_->pivotTolerance());
     897            abcFactorization_->pivotTolerance(newTolerance);
     898          } else {
     899            // Go to safe
     900            abcFactorization_->pivotTolerance(0.99);
     901          }
     902          forceFactorization_ = 1; // a bit drastic but ..
     903          type = 2;
     904          abcNonLinearCost_->checkInfeasibilities();
     905          if (internalFactorize(2) != 0) {
     906            largestPrimalError_ = 1.0e4; // force other type
     907          }
     908        }
     909        changeMade_++; // say change made
    914910      }
    915911    }
    916912    if (problemStatus_ != -4)
    917913      problemStatus_ = -3;
    918     if(abcProgress_.realInfeasibility_[0]<1.0e-1 &&
    919        primalTolerance_==1.0e-7&&abcProgress_.iterationNumber_[0]>0&&
    920        abcProgress_.iterationNumber_[CLP_PROGRESS-1]-abcProgress_.iterationNumber_[0]>25) {
     914    if (abcProgress_.realInfeasibility_[0] < 1.0e-1 && primalTolerance_ == 1.0e-7 && abcProgress_.iterationNumber_[0] > 0 && abcProgress_.iterationNumber_[CLP_PROGRESS - 1] - abcProgress_.iterationNumber_[0] > 25) {
    921915      // default - so user did not set
    922916      int iP;
    923       double minAverage=COIN_DBL_MAX;
    924       double maxAverage=0.0;
    925       for (iP=0;iP<CLP_PROGRESS;iP++) {
    926         int n=abcProgress_.numberInfeasibilities_[iP];
    927         if (!n) {
    928           break;
    929         } else {
    930           double average=abcProgress_.realInfeasibility_[iP];
    931           if (average>0.1)
    932             break;
    933           average /= static_cast<double>(n);
    934           minAverage=CoinMin(minAverage,average);
    935           maxAverage=CoinMax(maxAverage,average);
    936         }
    937       }
    938       if (iP==CLP_PROGRESS&&minAverage<1.0e-5&&maxAverage<1.0e-3) {
    939         // change tolerance
    940 #if CBC_USEFUL_PRINTING>0
    941         printf("CCchanging tolerance\n");
    942 #endif
    943         primalTolerance_=1.0e-6;
    944         dblParam_[ClpPrimalTolerance]=1.0e-6;
    945         moreSpecialOptions_ |= 4194304;
     917      double minAverage = COIN_DBL_MAX;
     918      double maxAverage = 0.0;
     919      for (iP = 0; iP < CLP_PROGRESS; iP++) {
     920        int n = abcProgress_.numberInfeasibilities_[iP];
     921        if (!n) {
     922          break;
     923        } else {
     924          double average = abcProgress_.realInfeasibility_[iP];
     925          if (average > 0.1)
     926            break;
     927          average /= static_cast< double >(n);
     928          minAverage = CoinMin(minAverage, average);
     929          maxAverage = CoinMax(maxAverage, average);
     930        }
     931      }
     932      if (iP == CLP_PROGRESS && minAverage < 1.0e-5 && maxAverage < 1.0e-3) {
     933        // change tolerance
     934#if CBC_USEFUL_PRINTING > 0
     935        printf("CCchanging tolerance\n");
     936#endif
     937        primalTolerance_ = 1.0e-6;
     938        dblParam_[ClpPrimalTolerance] = 1.0e-6;
     939        moreSpecialOptions_ |= 4194304;
    946940      }
    947941    }
     
    950944    // put back original costs and then check
    951945    // createRim(4); // costs do not change
    952     if (ifValuesPass&&numberIterations_==baseIteration_) {
     946    if (ifValuesPass && numberIterations_ == baseIteration_) {
    953947      abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
    954948      lastSumInfeasibility = abcNonLinearCost_->largestInfeasibility();
     
    960954    }
    961955#endif
    962     if (ifValuesPass&&numberIterations_==baseIteration_) {
    963       double * save = new double[numberRows_];
     956    if (ifValuesPass && numberIterations_ == baseIteration_) {
     957      double *save = new double[numberRows_];
    964958      int numberOut = 1;
    965959      while (numberOut) {
    966         for (int iRow = 0; iRow < numberRows_; iRow++) {
    967           int iPivot = abcPivotVariable_[iRow];
    968           save[iRow] = abcSolution_[iPivot];
    969         }
    970         gutsOfPrimalSolution(2);
    971         double badInfeasibility = abcNonLinearCost_->largestInfeasibility();
    972         numberOut = 0;
    973         // But may be very large rhs etc
    974         double useError = CoinMin(largestPrimalError_,
    975                                   1.0e5 / CoinAbcMaximumAbsElement(abcSolution_, numberTotal_));
    976         if ((lastSumInfeasibility < incomingInfeasibility_ || badInfeasibility >
    977              (CoinMax(10.0, 100.0 * lastSumInfeasibility)))
    978             && (badInfeasibility > 10.0 ||useError > 1.0e-3)) {
    979           //printf("Original largest infeas %g, now %g, primalError %g\n",
    980           //     lastSumInfeasibility,abcNonLinearCost_->largestInfeasibility(),
    981           //     largestPrimalError_);
    982           // throw out up to 1000 structurals
    983           int * sort = new int[numberRows_];
    984           // first put back solution and store difference
    985           for (int iRow = 0; iRow < numberRows_; iRow++) {
    986             int iPivot = abcPivotVariable_[iRow];
    987             double difference = fabs(abcSolution_[iPivot] - save[iRow]);
    988             abcSolution_[iPivot] = save[iRow];
    989             save[iRow] = difference;
    990           }
    991           abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
    992           if (handler_->logLevel()>1)
    993             printf("Largest infeasibility %g\n",abcNonLinearCost_->largestInfeasibility());
    994           int numberBasic = 0;
    995           for (int iRow = 0; iRow < numberRows_; iRow++) {
    996             int iPivot = abcPivotVariable_[iRow];
    997             if (iPivot >= numberRows_) {
    998               // column
    999               double difference = save[iRow];
    1000               if (difference > 1.0e-4) {
    1001                 sort[numberOut] = iRow;
    1002                 save[numberOut++] = -difference;
    1003                 if (getInternalStatus(iPivot) == basic)
    1004                   numberBasic++;
    1005               }
    1006             }
    1007           }
    1008           if (!numberBasic) {
    1009             //printf("no errors on basic - going to all slack - numberOut %d\n",numberOut);
    1010             // allow
    1011             numberOut = 0;
    1012           }
    1013           CoinSort_2(save, save + numberOut, sort);
    1014           numberOut = CoinMin(1000, numberOut);
    1015           // for now bring in any slack
    1016           int jRow=0;
    1017           for (int i = 0; i < numberOut; i++) {
    1018             int iRow = sort[i];
    1019             int iColumn = abcPivotVariable_[iRow];
    1020             assert (getInternalStatus(iColumn)==basic);
    1021             setInternalStatus(iColumn, superBasic);
    1022             while (getInternalStatus(jRow)==basic)
    1023               jRow++;
    1024             setInternalStatus(jRow, basic);
    1025             abcPivotVariable_[iRow] = jRow;
    1026             if (fabs(abcSolution_[iColumn]) > 1.0e10) {
    1027               if (abcUpper_[iColumn] < 0.0) {
    1028                 abcSolution_[iColumn] = abcUpper_[iColumn];
    1029               } else if (abcLower_[iColumn] > 0.0) {
    1030                 abcSolution_[iColumn] = abcLower_[iColumn];
    1031               } else {
    1032                 abcSolution_[iColumn] = 0.0;
    1033               }
    1034             }
    1035           }
    1036           delete [] sort;
    1037         }
    1038         if (numberOut) {
    1039           double savePivotTolerance = abcFactorization_->pivotTolerance();
    1040           abcFactorization_->pivotTolerance(0.99);
    1041           int factorStatus = internalFactorize(12);
    1042           abcFactorization_->pivotTolerance(savePivotTolerance);
    1043           assert (!factorStatus);
    1044         }
    1045       }
    1046       delete [] save;
     960        for (int iRow = 0; iRow < numberRows_; iRow++) {
     961          int iPivot = abcPivotVariable_[iRow];
     962          save[iRow] = abcSolution_[iPivot];
     963        }
     964        gutsOfPrimalSolution(2);
     965        double badInfeasibility = abcNonLinearCost_->largestInfeasibility();
     966        numberOut = 0;
     967        // But may be very large rhs etc
     968        double useError = CoinMin(largestPrimalError_,
     969          1.0e5 / CoinAbcMaximumAbsElement(abcSolution_, numberTotal_));
     970        if ((lastSumInfeasibility < incomingInfeasibility_ || badInfeasibility > (CoinMax(10.0, 100.0 * lastSumInfeasibility)))
     971          && (badInfeasibility > 10.0 || useError > 1.0e-3)) {
     972          //printf("Original largest infeas %g, now %g, primalError %g\n",
     973          //     lastSumInfeasibility,abcNonLinearCost_->largestInfeasibility(),
     974          //     largestPrimalError_);
     975          // throw out up to 1000 structurals
     976          int *sort = new int[numberRows_];
     977          // first put back solution and store difference
     978          for (int iRow = 0; iRow < numberRows_; iRow++) {
     979            int iPivot = abcPivotVariable_[iRow];
     980            double difference = fabs(abcSolution_[iPivot] - save[iRow]);
     981            abcSolution_[iPivot] = save[iRow];
     982            save[iRow] = difference;
     983          }
     984          abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
     985          if (handler_->logLevel() > 1)
     986            printf("Largest infeasibility %g\n", abcNonLinearCost_->largestInfeasibility());
     987          int numberBasic = 0;
     988          for (int iRow = 0; iRow < numberRows_; iRow++) {
     989            int iPivot = abcPivotVariable_[iRow];
     990            if (iPivot >= numberRows_) {
     991              // column
     992              double difference = save[iRow];
     993              if (difference > 1.0e-4) {
     994                sort[numberOut] = iRow;
     995                save[numberOut++] = -difference;
     996                if (getInternalStatus(iPivot) == basic)
     997                  numberBasic++;
     998              }
     999            }
     1000          }
     1001          if (!numberBasic) {
     1002            //printf("no errors on basic - going to all slack - numberOut %d\n",numberOut);
     1003            // allow
     1004            numberOut = 0;
     1005          }
     1006          CoinSort_2(save, save + numberOut, sort);
     1007          numberOut = CoinMin(1000, numberOut);
     1008          // for now bring in any slack
     1009          int jRow = 0;
     1010          for (int i = 0; i < numberOut; i++) {
     1011            int iRow = sort[i];
     1012            int iColumn = abcPivotVariable_[iRow];
     1013            assert(getInternalStatus(iColumn) == basic);
     1014            setInternalStatus(iColumn, superBasic);
     1015            while (getInternalStatus(jRow) == basic)
     1016              jRow++;
     1017            setInternalStatus(jRow, basic);
     1018            abcPivotVariable_[iRow] = jRow;
     1019            if (fabs(abcSolution_[iColumn]) > 1.0e10) {
     1020              if (abcUpper_[iColumn] < 0.0) {
     1021                abcSolution_[iColumn] = abcUpper_[iColumn];
     1022              } else if (abcLower_[iColumn] > 0.0) {
     1023                abcSolution_[iColumn] = abcLower_[iColumn];
     1024              } else {
     1025                abcSolution_[iColumn] = 0.0;
     1026              }
     1027            }
     1028          }
     1029          delete[] sort;
     1030        }
     1031        if (numberOut) {
     1032          double savePivotTolerance = abcFactorization_->pivotTolerance();
     1033          abcFactorization_->pivotTolerance(0.99);
     1034          int factorStatus = internalFactorize(12);
     1035          abcFactorization_->pivotTolerance(savePivotTolerance);
     1036          assert(!factorStatus);
     1037        }
     1038      }
     1039      delete[] save;
    10471040    }
    10481041    gutsOfPrimalSolution(3);
    1049     sumInfeasibility =  abcNonLinearCost_->sumInfeasibilities();
     1042    sumInfeasibility = abcNonLinearCost_->sumInfeasibilities();
    10501043    int reason2 = 0;
    10511044#if CLP_CAUTION
    1052 #if CLP_CAUTION==2
     1045#if CLP_CAUTION == 2
    10531046    double test2 = 1.0e5;
    10541047#else
    10551048    double test2 = 1.0e-1;
    10561049#endif
    1057     if (!lastSumInfeasibility && sumInfeasibility>1.0e3*primalTolerance_ &&
    1058         lastAverageInfeasibility < test2 && numberPivots > 10&&!ifValuesPass)
     1050    if (!lastSumInfeasibility && sumInfeasibility > 1.0e3 * primalTolerance_ && lastAverageInfeasibility < test2 && numberPivots > 10 && !ifValuesPass)
    10591051      reason2 = 3;
    1060     if (lastSumInfeasibility < 1.0e-6 && sumInfeasibility > 1.0e-3 &&
    1061         numberPivots > 10&&!ifValuesPass)
     1052    if (lastSumInfeasibility < 1.0e-6 && sumInfeasibility > 1.0e-3 && numberPivots > 10 && !ifValuesPass)
    10621053      reason2 = 4;
    10631054#endif
    10641055    if ((sumInfeasibility > 1.0e7 && sumInfeasibility > 100.0 * lastSumInfeasibility
    1065          && abcFactorization_->pivotTolerance() < 0.11) ||
    1066         (largestPrimalError_ > 1.0e10 && largestDualError_ > 1.0e10))
     1056          && abcFactorization_->pivotTolerance() < 0.11)
     1057      || (largestPrimalError_ > 1.0e10 && largestDualError_ > 1.0e10))
    10671058      reason2 = 2;
    10681059    if (reason2) {
    1069       takenAction=true;
     1060      takenAction = true;
    10701061      problemStatus_ = tentativeStatus;
    10711062      doFactorization = true;
    10721063      if (numberPivots) {
    1073         // go back
    1074         // trouble - restore solution
    1075         restoreGoodStatus(1);
    1076         sequenceIn_=-1;
    1077         sequenceOut_=-1;
    1078         abcNonLinearCost_->checkInfeasibilities();
    1079         if (reason2 < 3) {
    1080           // Go to safe
    1081           abcFactorization_->pivotTolerance(CoinMin(0.99, 1.01 * abcFactorization_->pivotTolerance()));
    1082           forceFactorization_ = 1; // a bit drastic but ..
    1083         } else if (forceFactorization_ < 0) {
    1084           forceFactorization_ = CoinMin(numberPivots / 4, 100);
    1085         } else {
    1086           forceFactorization_ = CoinMin(forceFactorization_,
    1087                                         CoinMax(3, numberPivots / 4));
    1088         }
    1089         pivotRow_ = -1; // say no weights update
    1090         changeMade_++; // say change made
    1091         if (numberPivots == 1) {
    1092           // throw out something
    1093           if (sequenceIn_ >= 0 && getInternalStatus(sequenceIn_) != basic) {
    1094             setFlagged(sequenceIn_);
    1095           } else if (sequenceOut_ >= 0 && getInternalStatus(sequenceOut_) != basic) {
    1096             setFlagged(sequenceOut_);
    1097           }
    1098           abcProgress_.incrementTimesFlagged();
    1099         }
    1100         type = 2; // so will restore weights
    1101         if (internalFactorize(2) != 0) {
    1102           largestPrimalError_ = 1.0e4; // force other type
    1103         }
    1104         numberPivots = 0;
    1105         gutsOfPrimalSolution(3);
    1106         sumInfeasibility = abcNonLinearCost_->sumInfeasibilities();
     1064        // go back
     1065        // trouble - restore solution
     1066        restoreGoodStatus(1);
     1067        sequenceIn_ = -1;
     1068        sequenceOut_ = -1;
     1069        abcNonLinearCost_->checkInfeasibilities();
     1070        if (reason2 < 3) {
     1071          // Go to safe
     1072          abcFactorization_->pivotTolerance(CoinMin(0.99, 1.01 * abcFactorization_->pivotTolerance()));
     1073          forceFactorization_ = 1; // a bit drastic but ..
     1074        } else if (forceFactorization_ < 0) {
     1075          forceFactorization_ = CoinMin(numberPivots / 4, 100);
     1076        } else {
     1077          forceFactorization_ = CoinMin(forceFactorization_,
     1078            CoinMax(3, numberPivots / 4));
     1079        }
     1080        pivotRow_ = -1; // say no weights update
     1081        changeMade_++; // say change made
     1082        if (numberPivots == 1) {
     1083          // throw out something
     1084          if (sequenceIn_ >= 0 && getInternalStatus(sequenceIn_) != basic) {
     1085            setFlagged(sequenceIn_);
     1086          } else if (sequenceOut_ >= 0 && getInternalStatus(sequenceOut_) != basic) {
     1087            setFlagged(sequenceOut_);
     1088          }
     1089          abcProgress_.incrementTimesFlagged();
     1090        }
     1091        type = 2; // so will restore weights
     1092        if (internalFactorize(2) != 0) {
     1093          largestPrimalError_ = 1.0e4; // force other type
     1094        }
     1095        numberPivots = 0;
     1096        gutsOfPrimalSolution(3);
     1097        sumInfeasibility = abcNonLinearCost_->sumInfeasibilities();
    11071098      }
    11081099    }
     
    11171108  // If in primal and small dj give up
    11181109  if ((specialOptions_ & 1024) != 0 && !numberPrimalInfeasibilities_ && numberDualInfeasibilities_) {
    1119     double average = sumDualInfeasibilities_ / (static_cast<double> (numberDualInfeasibilities_));
     1110    double average = sumDualInfeasibilities_ / (static_cast< double >(numberDualInfeasibilities_));
    11201111    if (numberIterations_ > 300 && average < 1.0e-4) {
    11211112      numberDualInfeasibilities_ = 0;
     
    11251116  // Check if looping
    11261117  int loop;
    1127   ifValuesPass=firstFree_>=0;
     1118  ifValuesPass = firstFree_ >= 0;
    11281119  if (type != 2 && !ifValuesPass)
    11291120    loop = abcProgress_.looping();
     
    11311122    loop = -1;
    11321123  if ((moreSpecialOptions_ & 2048) != 0 && !numberPrimalInfeasibilities_ && numberDualInfeasibilities_) {
    1133     double average = sumDualInfeasibilities_ / (static_cast<double> (numberDualInfeasibilities_));
    1134     if (abcProgress_.lastIterationNumber(2)==numberIterations_&&
    1135         ((abcProgress_.timesFlagged()>2&&average < 1.0e-1)||
    1136          abcProgress_.timesFlagged()>20)) {
     1124    double average = sumDualInfeasibilities_ / (static_cast< double >(numberDualInfeasibilities_));
     1125    if (abcProgress_.lastIterationNumber(2) == numberIterations_ && ((abcProgress_.timesFlagged() > 2 && average < 1.0e-1) || abcProgress_.timesFlagged() > 20)) {
    11371126      numberDualInfeasibilities_ = 0;
    11381127      sumDualInfeasibilities_ = 0.0;
    1139       problemStatus_=3;
    1140       loop=0;
     1128      problemStatus_ = 3;
     1129      loop = 0;
    11411130    }
    11421131  }
     
    11541143    }
    11551144    problemStatus_ = 10; // instead - try other algorithm
    1156     return ;
     1145    return;
    11571146  } else if (loop < -1) {
    11581147    // Is it time for drastic measures
    1159     if (abcNonLinearCost_->numberInfeasibilities() && abcProgress_.badTimes() > 5 &&
    1160         abcProgress_.oddState() < 10 && abcProgress_.oddState() >= 0) {
     1148    if (abcNonLinearCost_->numberInfeasibilities() && abcProgress_.badTimes() > 5 && abcProgress_.oddState() < 10 && abcProgress_.oddState() >= 0) {
    11611149      abcProgress_.newOddState();
    11621150      abcNonLinearCost_->zapCosts();
     
    11741162  }
    11751163  abcProgress_.modifyObjective(abcNonLinearCost_->feasibleReportCost());
    1176   if (!lastNumberInfeasibility && sumInfeasibility && numberPivots > 1&&!ifValuesPass&&
    1177       !takenAction&&
    1178       abcProgress_.lastObjective(0)>=abcProgress_.lastObjective(CLP_PROGRESS-1)) {
     1164  if (!lastNumberInfeasibility && sumInfeasibility && numberPivots > 1 && !ifValuesPass && !takenAction && abcProgress_.lastObjective(0) >= abcProgress_.lastObjective(CLP_PROGRESS - 1)) {
    11791165    // first increase minimumThetaMovement_;
    11801166    // be more careful
    11811167    //abcFactorization_->saferTolerances (-0.99, -1.002);
    1182 #if ABC_NORMAL_DEBUG>0
     1168#if ABC_NORMAL_DEBUG > 0
    11831169    if (handler_->logLevel() == 63)
    1184       printf("thought feasible but sum is %g force %d minimum theta %g\n", sumInfeasibility, 
    1185              forceFactorization_,minimumThetaMovement_);
     1170      printf("thought feasible but sum is %g force %d minimum theta %g\n", sumInfeasibility,
     1171        forceFactorization_, minimumThetaMovement_);
    11861172#endif
    11871173    stateOfProblem_ |= PESSIMISTIC;
    1188     if (minimumThetaMovement_<1.0e-10) {
     1174    if (minimumThetaMovement_ < 1.0e-10) {
    11891175      minimumThetaMovement_ *= 1.1;
    11901176    } else if (0) {
    11911177      int maxFactor = abcFactorization_->maximumPivots();
    11921178      if (maxFactor > 10) {
    1193         if (forceFactorization_ < 0)
    1194           forceFactorization_ = maxFactor;
    1195         forceFactorization_ = CoinMax(1, (forceFactorization_ >> 2));
    1196         if (handler_->logLevel() == 63)
    1197           printf("Reducing factorization frequency to %d\n", forceFactorization_);
     1179        if (forceFactorization_ < 0)
     1180          forceFactorization_ = maxFactor;
     1181        forceFactorization_ = CoinMax(1, (forceFactorization_ >> 2));
     1182        if (handler_->logLevel() == 63)
     1183          printf("Reducing factorization frequency to %d\n", forceFactorization_);
    11981184      }
    11991185    }
     
    12051191  //problemStatus_=-1;;
    12061192  progressFlag_ = 0; //reset progress flag
    1207  
     1193
    12081194  handler_->message(CLP_SIMPLEX_STATUS, messages_)
    1209     << numberIterations_ << abcNonLinearCost_->feasibleReportCost()+objectiveOffset_;
     1195    << numberIterations_ << abcNonLinearCost_->feasibleReportCost() + objectiveOffset_;
    12101196  handler_->printing(abcNonLinearCost_->numberInfeasibilities() > 0)
    12111197    << abcNonLinearCost_->sumInfeasibilities() << abcNonLinearCost_->numberInfeasibilities();
     
    12131199    << sumDualInfeasibilities_ << numberDualInfeasibilities_;
    12141200  handler_->printing(numberDualInfeasibilitiesWithoutFree_
    1215                      < numberDualInfeasibilities_)
     1201    < numberDualInfeasibilities_)
    12161202    << numberDualInfeasibilitiesWithoutFree_;
    12171203  handler_->message() << CoinMessageEol;
     
    12441230    lastObj3 += infeasibilityCost_ * 2.0 * lastInf3;
    12451231    if (lastObj < thisObj - 1.0e-5 * CoinMax(fabs(thisObj), fabs(lastObj)) - 1.0e-7
    1246         && firstFree_ < 0 && thisInf >= lastInf) {
    1247 #if ABC_NORMAL_DEBUG>0
     1232      && firstFree_ < 0 && thisInf >= lastInf) {
     1233#if ABC_NORMAL_DEBUG > 0
    12481234      if (handler_->logLevel() == 63)
    1249         printf("lastobj %g this %g force %d\n", lastObj, thisObj, forceFactorization_);
     1235        printf("lastobj %g this %g force %d\n", lastObj, thisObj, forceFactorization_);
    12501236#endif
    12511237      int maxFactor = abcFactorization_->maximumPivots();
    12521238      if (maxFactor > 10) {
    1253         if (forceFactorization_ < 0)
    1254           forceFactorization_ = maxFactor;
    1255         forceFactorization_ = CoinMax(1, (forceFactorization_ >> 2));
    1256 #if ABC_NORMAL_DEBUG>0
    1257         if (handler_->logLevel() == 63)
    1258           printf("Reducing factorization frequency to %d\n", forceFactorization_);
     1239        if (forceFactorization_ < 0)
     1240          forceFactorization_ = maxFactor;
     1241        forceFactorization_ = CoinMax(1, (forceFactorization_ >> 2));
     1242#if ABC_NORMAL_DEBUG > 0
     1243        if (handler_->logLevel() == 63)
     1244          printf("Reducing factorization frequency to %d\n", forceFactorization_);
    12591245#endif
    12601246      }
    12611247    } else if (lastObj3 < thisObj - 1.0e-5 * CoinMax(fabs(thisObj), fabs(lastObj3)) - 1.0e-7
    1262                && firstFree_ < 0 && thisInf >= lastInf) {
    1263 #if ABC_NORMAL_DEBUG>0
     1248      && firstFree_ < 0 && thisInf >= lastInf) {
     1249#if ABC_NORMAL_DEBUG > 0
    12641250      if (handler_->logLevel() == 63)
    1265         printf("lastobj3 %g this3 %g force %d\n", lastObj3, thisObj, forceFactorization_);
     1251        printf("lastobj3 %g this3 %g force %d\n", lastObj3, thisObj, forceFactorization_);
    12661252#endif
    12671253      int maxFactor = abcFactorization_->maximumPivots();
    12681254      if (maxFactor > 10) {
    1269         if (forceFactorization_ < 0)
    1270           forceFactorization_ = maxFactor;
    1271         forceFactorization_ = CoinMax(1, (forceFactorization_ * 2) / 3);
    1272 #if ABC_NORMAL_DEBUG>0
    1273         if (handler_->logLevel() == 63)
    1274           printf("Reducing factorization frequency to %d\n", forceFactorization_);
    1275 #endif
    1276       }
    1277     } else if(lastInf < testValue || trueInfeasibility == 1.123456e10) {
     1255        if (forceFactorization_ < 0)
     1256          forceFactorization_ = maxFactor;
     1257        forceFactorization_ = CoinMax(1, (forceFactorization_ * 2) / 3);
     1258#if ABC_NORMAL_DEBUG > 0
     1259        if (handler_->logLevel() == 63)
     1260          printf("Reducing factorization frequency to %d\n", forceFactorization_);
     1261#endif
     1262      }
     1263    } else if (lastInf < testValue || trueInfeasibility == 1.123456e10) {
    12781264      if (infeasibilityCost_ < 1.0e14) {
    1279         infeasibilityCost_ *= 1.5;
    1280         // reset looping criterion
    1281         abcProgress_.reset();
    1282 #if ABC_NORMAL_DEBUG>0
    1283         if (handler_->logLevel() == 63)
    1284           printf("increasing weight to %g\n", infeasibilityCost_);
    1285 #endif
    1286         gutsOfPrimalSolution(3);
     1265        infeasibilityCost_ *= 1.5;
     1266        // reset looping criterion
     1267        abcProgress_.reset();
     1268#if ABC_NORMAL_DEBUG > 0
     1269        if (handler_->logLevel() == 63)
     1270          printf("increasing weight to %g\n", infeasibilityCost_);
     1271#endif
     1272        gutsOfPrimalSolution(3);
    12871273      }
    12881274    }
     
    12921278#if CLP_CAUTION
    12931279  // If twice nearly there ...
    1294   if (lastAverageInfeasibility<2.0*dualTolerance_) {
    1295     double averageInfeasibility = sumDualInfeasibilities_ /
    1296       static_cast<double>(numberDualInfeasibilities_ + 1);
    1297     printf("last av %g now %g\n",lastAverageInfeasibility,
    1298            averageInfeasibility);
    1299     if (averageInfeasibility<2.0*dualTolerance_)
     1280  if (lastAverageInfeasibility < 2.0 * dualTolerance_) {
     1281    double averageInfeasibility = sumDualInfeasibilities_ / static_cast< double >(numberDualInfeasibilities_ + 1);
     1282    printf("last av %g now %g\n", lastAverageInfeasibility,
     1283      averageInfeasibility);
     1284    if (averageInfeasibility < 2.0 * dualTolerance_)
    13001285      sumOfRelaxedDualInfeasibilities_ = 0.0;
    13011286  }
    13021287#endif
    13031288  // give code benefit of doubt
    1304   if (sumOfRelaxedDualInfeasibilities_ == 0.0 &&
    1305       sumOfRelaxedPrimalInfeasibilities_ == 0.0) {
     1289  if (sumOfRelaxedDualInfeasibilities_ == 0.0 && sumOfRelaxedPrimalInfeasibilities_ == 0.0) {
    13061290    // say been feasible
    13071291    abcState_ |= CLP_ABC_BEEN_FEASIBLE;
     
    13181302      double error = CoinMin(1.0e-2, largestPrimalError_);
    13191303      // allow tolerance at least slightly bigger than standard
    1320       relaxedToleranceP = relaxedToleranceP +  error;
     1304      relaxedToleranceP = relaxedToleranceP + error;
    13211305      int ninfeas = abcNonLinearCost_->numberInfeasibilities();
    13221306      double sum = abcNonLinearCost_->sumInfeasibilities();
    1323       double average = sum / static_cast<double> (ninfeas);
    1324 #if ABC_NORMAL_DEBUG>3
     1307      double average = sum / static_cast< double >(ninfeas);
     1308#if ABC_NORMAL_DEBUG > 3
    13251309      if (handler_->logLevel() > 0)
    1326         printf("nonLinearCost says infeasible %d summing to %g\n",
    1327                ninfeas, sum);
     1310        printf("nonLinearCost says infeasible %d summing to %g\n",
     1311          ninfeas, sum);
    13281312#endif
    13291313      if (average > relaxedToleranceP) {
    1330         sumOfRelaxedPrimalInfeasibilities_ = sum;
    1331         numberPrimalInfeasibilities_ = ninfeas;
    1332         sumPrimalInfeasibilities_ = sum;
    1333 #if ABC_NORMAL_DEBUG>3
    1334         bool unflagged =
    1335 #endif
    1336           unflag();
    1337         abcProgress_.clearTimesFlagged();
    1338 #if ABC_NORMAL_DEBUG>3
    1339         if (unflagged && handler_->logLevel() > 0)
    1340           printf(" - but flagged variables\n");
     1314        sumOfRelaxedPrimalInfeasibilities_ = sum;
     1315        numberPrimalInfeasibilities_ = ninfeas;
     1316        sumPrimalInfeasibilities_ = sum;
     1317#if ABC_NORMAL_DEBUG > 3
     1318        bool unflagged =
     1319#endif
     1320          unflag();
     1321        abcProgress_.clearTimesFlagged();
     1322#if ABC_NORMAL_DEBUG > 3
     1323        if (unflagged && handler_->logLevel() > 0)
     1324          printf(" - but flagged variables\n");
    13411325#endif
    13421326      }
     
    13471331  if ((dualFeasible() || problemStatus_ == -4) && !ifValuesPass) {
    13481332    // see if extra helps
    1349     if (abcNonLinearCost_->numberInfeasibilities() &&
    1350         (abcNonLinearCost_->sumInfeasibilities() > 1.0e-3 || sumOfRelaxedPrimalInfeasibilities_)
    1351         && !alwaysOptimal) {
     1333    if (abcNonLinearCost_->numberInfeasibilities() && (abcNonLinearCost_->sumInfeasibilities() > 1.0e-3 || sumOfRelaxedPrimalInfeasibilities_)
     1334      && !alwaysOptimal) {
    13521335      //may need infeasiblity cost changed
    13531336      // we can see if we can construct a ray
     
    13551338      // put non-basics to bounds in case tolerance moved
    13561339      // put back original costs
    1357       CoinAbcMemcpy(abcCost_,costSaved_,numberTotal_);;
     1340      CoinAbcMemcpy(abcCost_, costSaved_, numberTotal_);
     1341      ;
    13581342      abcNonLinearCost_->checkInfeasibilities(0.0);
    13591343      gutsOfPrimalSolution(3);
    13601344      // put back original costs
    1361       CoinAbcMemcpy(abcCost_,costSaved_,numberTotal_);;
     1345      CoinAbcMemcpy(abcCost_, costSaved_, numberTotal_);
     1346      ;
    13621347      abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
    13631348      gutsOfPrimalSolution(3);
    13641349      // may have fixed infeasibilities - double check
    13651350      if (abcNonLinearCost_->numberInfeasibilities() == 0) {
    1366         // carry on
    1367         problemStatus_ = -1;
    1368         abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
     1351        // carry on
     1352        problemStatus_ = -1;
     1353        abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
    13691354      } else {
    1370         if ((infeasibilityCost_ >= 1.0e18 ||
    1371              numberDualInfeasibilities_ == 0) && perturbation_ == 101
    1372             && (specialOptions_&8192)==0) {
    1373           goToDual = unPerturb(); // stop any further perturbation
     1355        if ((infeasibilityCost_ >= 1.0e18 || numberDualInfeasibilities_ == 0) && perturbation_ == 101
     1356          && (specialOptions_ & 8192) == 0) {
     1357          goToDual = unPerturb(); // stop any further perturbation
    13741358#ifndef TRY_ABC_GUS
    1375           if (abcNonLinearCost_->sumInfeasibilities() > 1.0e-1)
    1376             goToDual = false;
    1377 #endif
    1378           numberDualInfeasibilities_ = 1; // carry on
    1379           problemStatus_ = -1;
    1380         } else if (numberDualInfeasibilities_ == 0 && largestDualError_ > 1.0e-2
     1359          if (abcNonLinearCost_->sumInfeasibilities() > 1.0e-1)
     1360            goToDual = false;
     1361#endif
     1362          numberDualInfeasibilities_ = 1; // carry on
     1363          problemStatus_ = -1;
     1364        } else if (numberDualInfeasibilities_ == 0 && largestDualError_ > 1.0e-2
    13811365#ifndef TRY_ABC_GUS
    1382                    &&((moreSpecialOptions_ & 256) == 0&&(specialOptions_ & 8192) == 0)
     1366          && ((moreSpecialOptions_ & 256) == 0 && (specialOptions_ & 8192) == 0)
    13831367#else
    1384                    &&(specialOptions_ & 8192) == 0
    1385 #endif
    1386                    ) {
    1387           goToDual = true;
    1388           abcFactorization_->pivotTolerance(CoinMax(0.5, abcFactorization_->pivotTolerance()));
    1389         }
    1390         if (!goToDual) {
    1391           if (infeasibilityCost_ >= 1.0e20 ||
    1392               numberDualInfeasibilities_ == 0) {
    1393             // we are infeasible - use as ray
    1394             delete [] ray_;
    1395             ray_ = new double [numberRows_];
    1396             CoinMemcpyN(dual_, numberRows_, ray_);
    1397             // and get feasible duals
    1398             infeasibilityCost_ = 0.0;
    1399             CoinAbcMemcpy(abcCost_,costSaved_,numberTotal_);;
    1400             abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
    1401             gutsOfPrimalSolution(3);
    1402             // so will exit
    1403             infeasibilityCost_ = 1.0e30;
    1404             // reset infeasibilities
    1405             sumPrimalInfeasibilities_ = abcNonLinearCost_->sumInfeasibilities();;
    1406             numberPrimalInfeasibilities_ =
    1407               abcNonLinearCost_->numberInfeasibilities();
    1408           }
    1409           if (infeasibilityCost_ < 1.0e20) {
    1410             infeasibilityCost_ *= 5.0;
    1411             // reset looping criterion
    1412             abcProgress_.reset();
    1413             changeMade_++; // say change made
    1414             handler_->message(CLP_PRIMAL_WEIGHT, messages_)
    1415               << infeasibilityCost_
    1416               << CoinMessageEol;
    1417             // put back original costs and then check
    1418             CoinAbcMemcpy(abcCost_,costSaved_,numberTotal_);;
    1419             abcNonLinearCost_->checkInfeasibilities(0.0);
    1420             gutsOfPrimalSolution(3);
    1421             problemStatus_ = -1; //continue
     1368          && (specialOptions_ & 8192) == 0
     1369#endif
     1370        ) {
     1371          goToDual = true;
     1372          abcFactorization_->pivotTolerance(CoinMax(0.5, abcFactorization_->pivotTolerance()));
     1373        }
     1374        if (!goToDual) {
     1375          if (infeasibilityCost_ >= 1.0e20 || numberDualInfeasibilities_ == 0) {
     1376            // we are infeasible - use as ray
     1377            delete[] ray_;
     1378            ray_ = new double[numberRows_];
     1379            CoinMemcpyN(dual_, numberRows_, ray_);
     1380            // and get feasible duals
     1381            infeasibilityCost_ = 0.0;
     1382            CoinAbcMemcpy(abcCost_, costSaved_, numberTotal_);
     1383            ;
     1384            abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
     1385            gutsOfPrimalSolution(3);
     1386            // so will exit
     1387            infeasibilityCost_ = 1.0e30;
     1388            // reset infeasibilities
     1389            sumPrimalInfeasibilities_ = abcNonLinearCost_->sumInfeasibilities();
     1390            ;
     1391            numberPrimalInfeasibilities_ = abcNonLinearCost_->numberInfeasibilities();
     1392          }
     1393          if (infeasibilityCost_ < 1.0e20) {
     1394            infeasibilityCost_ *= 5.0;
     1395            // reset looping criterion
     1396            abcProgress_.reset();
     1397            changeMade_++; // say change made
     1398            handler_->message(CLP_PRIMAL_WEIGHT, messages_)
     1399              << infeasibilityCost_
     1400              << CoinMessageEol;
     1401            // put back original costs and then check
     1402            CoinAbcMemcpy(abcCost_, costSaved_, numberTotal_);
     1403            ;
     1404            abcNonLinearCost_->checkInfeasibilities(0.0);
     1405            gutsOfPrimalSolution(3);
     1406            problemStatus_ = -1; //continue
    14221407#ifndef TRY_ABC_GUS
    1423             goToDual = false;
     1408            goToDual = false;
    14241409#else
    1425             if((specialOptions_&8192)==0&&!sumOfRelaxedDualInfeasibilities_)
    1426               goToDual=true;
    1427 #endif
    1428           } else {
    1429             // say infeasible
    1430             problemStatus_ = 1;
    1431           }
    1432         }
     1410            if ((specialOptions_ & 8192) == 0 && !sumOfRelaxedDualInfeasibilities_)
     1411              goToDual = true;
     1412#endif
     1413          } else {
     1414            // say infeasible
     1415            problemStatus_ = 1;
     1416          }
     1417        }
    14331418      }
    14341419    } else {
    14351420      // may be optimal
    14361421      if (perturbation_ == 101) {
    1437         goToDual = unPerturb(); // stop any further perturbation
     1422        goToDual = unPerturb(); // stop any further perturbation
    14381423#ifndef TRY_ABC_GUS
    1439         if ((numberRows_ > 20000 || numberDualInfeasibilities_) && !numberTimesOptimal_)
     1424        if ((numberRows_ > 20000 || numberDualInfeasibilities_) && !numberTimesOptimal_)
    14401425#else
    1441           if ((specialOptions_&8192)!=0)
    1442 #endif
    1443           goToDual = false; // Better to carry on a bit longer
    1444         lastCleaned_ = -1; // carry on
     1426        if ((specialOptions_ & 8192) != 0)
     1427#endif
     1428          goToDual = false; // Better to carry on a bit longer
     1429        lastCleaned_ = -1; // carry on
    14451430      }
    14461431      bool unflagged = (unflag() != 0);
    14471432      abcProgress_.clearTimesFlagged();
    1448       if ( lastCleaned_ != numberIterations_ || unflagged) {
    1449         handler_->message(CLP_PRIMAL_OPTIMAL, messages_)
    1450           << primalTolerance_
    1451           << CoinMessageEol;
    1452         if (numberTimesOptimal_ < 4) {
    1453           numberTimesOptimal_++;
    1454           changeMade_++; // say change made
    1455           if (numberTimesOptimal_ == 1) {
    1456             // better to have small tolerance even if slower
    1457             abcFactorization_->zeroTolerance(CoinMin(abcFactorization_->zeroTolerance(), 1.0e-15));
    1458           }
    1459           lastCleaned_ = numberIterations_;
    1460           if (primalTolerance_ != dblParam_[ClpPrimalTolerance])
    1461             handler_->message(CLP_PRIMAL_ORIGINAL, messages_)
    1462               << CoinMessageEol;
    1463           double oldTolerance = primalTolerance_;
    1464           primalTolerance_ = dblParam_[ClpPrimalTolerance];
    1465           // put back original costs and then check
     1433      if (lastCleaned_ != numberIterations_ || unflagged) {
     1434        handler_->message(CLP_PRIMAL_OPTIMAL, messages_)
     1435          << primalTolerance_
     1436          << CoinMessageEol;
     1437        if (numberTimesOptimal_ < 4) {
     1438          numberTimesOptimal_++;
     1439          changeMade_++; // say change made
     1440          if (numberTimesOptimal_ == 1) {
     1441            // better to have small tolerance even if slower
     1442            abcFactorization_->zeroTolerance(CoinMin(abcFactorization_->zeroTolerance(), 1.0e-15));
     1443          }
     1444          lastCleaned_ = numberIterations_;
     1445          if (primalTolerance_ != dblParam_[ClpPrimalTolerance])
     1446            handler_->message(CLP_PRIMAL_ORIGINAL, messages_)
     1447              << CoinMessageEol;
     1448          double oldTolerance = primalTolerance_;
     1449          primalTolerance_ = dblParam_[ClpPrimalTolerance];
     1450          // put back original costs and then check
    14661451#if 0 //ndef NDEBUG
    14671452          double largestDifference=0.0;
     
    14711456            printf("largest change in cost %g\n",largestDifference);
    14721457#endif
    1473           CoinAbcMemcpy(abcCost_,costSaved_,numberTotal_);;
    1474           abcNonLinearCost_->checkInfeasibilities(oldTolerance);
    1475           gutsOfPrimalSolution(3);
    1476           if (sumOfRelaxedDualInfeasibilities_ == 0.0 &&
    1477               sumOfRelaxedPrimalInfeasibilities_ == 0.0) {
    1478             // say optimal (with these bounds etc)
    1479             numberDualInfeasibilities_ = 0;
    1480             sumDualInfeasibilities_ = 0.0;
    1481             numberPrimalInfeasibilities_ = 0;
    1482             sumPrimalInfeasibilities_ = 0.0;
    1483             goToDual=false;
    1484           }
    1485           if (dualFeasible() && !abcNonLinearCost_->numberInfeasibilities() && lastCleaned_ >= 0)
    1486             problemStatus_ = 0;
    1487           else
    1488             problemStatus_ = -1;
    1489         } else {
    1490           problemStatus_ = 0; // optimal
    1491           if (lastCleaned_ < numberIterations_) {
    1492             handler_->message(CLP_SIMPLEX_GIVINGUP, messages_)
    1493               << CoinMessageEol;
    1494           }
    1495         }
     1458          CoinAbcMemcpy(abcCost_, costSaved_, numberTotal_);
     1459          ;
     1460          abcNonLinearCost_->checkInfeasibilities(oldTolerance);
     1461          gutsOfPrimalSolution(3);
     1462          if (sumOfRelaxedDualInfeasibilities_ == 0.0 && sumOfRelaxedPrimalInfeasibilities_ == 0.0) {
     1463            // say optimal (with these bounds etc)
     1464            numberDualInfeasibilities_ = 0;
     1465            sumDualInfeasibilities_ = 0.0;
     1466            numberPrimalInfeasibilities_ = 0;
     1467            sumPrimalInfeasibilities_ = 0.0;
     1468            goToDual = false;
     1469          }
     1470          if (dualFeasible() && !abcNonLinearCost_->numberInfeasibilities() && lastCleaned_ >= 0)
     1471            problemStatus_ = 0;
     1472          else
     1473            problemStatus_ = -1;
     1474        } else {
     1475          problemStatus_ = 0; // optimal
     1476          if (lastCleaned_ < numberIterations_) {
     1477            handler_->message(CLP_SIMPLEX_GIVINGUP, messages_)
     1478              << CoinMessageEol;
     1479          }
     1480        }
    14961481      } else {
    1497         if (!alwaysOptimal || !sumOfRelaxedPrimalInfeasibilities_)
    1498           problemStatus_ = 0; // optimal
    1499         else
    1500           problemStatus_ = 1; // infeasible
     1482        if (!alwaysOptimal || !sumOfRelaxedPrimalInfeasibilities_)
     1483          problemStatus_ = 0; // optimal
     1484        else
     1485          problemStatus_ = 1; // infeasible
    15011486      }
    15021487    }
     
    15051490    if (problemStatus_ == -5) {
    15061491      if (abcNonLinearCost_->numberInfeasibilities()) {
    1507         if (infeasibilityCost_ > 1.0e18 && perturbation_ == 101) {
    1508           // back off weight
    1509           infeasibilityCost_ = 1.0e13;
    1510           // reset looping criterion
    1511           abcProgress_.reset();
    1512           unPerturb(); // stop any further perturbation
    1513         }
    1514         //we need infeasiblity cost changed
    1515         if (infeasibilityCost_ < 1.0e20) {
    1516           infeasibilityCost_ *= 5.0;
    1517           // reset looping criterion
    1518           abcProgress_.reset();
    1519           changeMade_++; // say change made
    1520           handler_->message(CLP_PRIMAL_WEIGHT, messages_)
    1521             << infeasibilityCost_
    1522             << CoinMessageEol;
    1523           // put back original costs and then check
    1524           CoinAbcMemcpy(abcCost_,costSaved_,numberTotal_);;
    1525           gutsOfPrimalSolution(3);
    1526           problemStatus_ = -1; //continue
    1527         } else {
    1528           // say infeasible
    1529           problemStatus_ = 1;
    1530           // we are infeasible - use as ray
    1531           delete [] ray_;
    1532           ray_ = new double [numberRows_];
    1533           CoinMemcpyN(dual_, numberRows_, ray_);
    1534         }
     1492        if (infeasibilityCost_ > 1.0e18 && perturbation_ == 101) {
     1493          // back off weight
     1494          infeasibilityCost_ = 1.0e13;
     1495          // reset looping criterion
     1496          abcProgress_.reset();
     1497          unPerturb(); // stop any further perturbation
     1498        }
     1499        //we need infeasiblity cost changed
     1500        if (infeasibilityCost_ < 1.0e20) {
     1501          infeasibilityCost_ *= 5.0;
     1502          // reset looping criterion
     1503          abcProgress_.reset();
     1504          changeMade_++; // say change made
     1505          handler_->message(CLP_PRIMAL_WEIGHT, messages_)
     1506            << infeasibilityCost_
     1507            << CoinMessageEol;
     1508          // put back original costs and then check
     1509          CoinAbcMemcpy(abcCost_, costSaved_, numberTotal_);
     1510          ;
     1511          gutsOfPrimalSolution(3);
     1512          problemStatus_ = -1; //continue
     1513        } else {
     1514          // say infeasible
     1515          problemStatus_ = 1;
     1516          // we are infeasible - use as ray
     1517          delete[] ray_;
     1518          ray_ = new double[numberRows_];
     1519          CoinMemcpyN(dual_, numberRows_, ray_);
     1520        }
    15351521      } else {
    1536         // say unbounded
    1537         problemStatus_ = 2;
     1522        // say unbounded
     1523        problemStatus_ = 2;
    15381524      }
    15391525    } else {
    15401526      // carry on
    15411527      problemStatus_ = -1;
    1542       if(type == 3 && !ifValuesPass) {
    1543         //bool unflagged =
    1544         unflag();
    1545         abcProgress_.clearTimesFlagged();
    1546         if (sumDualInfeasibilities_ < 1.0e-3 ||
    1547             (sumDualInfeasibilities_ / static_cast<double> (numberDualInfeasibilities_)) < 1.0e-5 ||
    1548             abcProgress_.lastIterationNumber(0) == numberIterations_) {
    1549           if (!numberPrimalInfeasibilities_) {
    1550             if (numberTimesOptimal_ < 4) {
    1551               numberTimesOptimal_++;
    1552               changeMade_++; // say change made
    1553             } else {
    1554               problemStatus_ = 0;
    1555               secondaryStatus_ = 5;
    1556             }
    1557           }
    1558         }
     1528      if (type == 3 && !ifValuesPass) {
     1529        //bool unflagged =
     1530        unflag();
     1531        abcProgress_.clearTimesFlagged();
     1532        if (sumDualInfeasibilities_ < 1.0e-3 || (sumDualInfeasibilities_ / static_cast< double >(numberDualInfeasibilities_)) < 1.0e-5 || abcProgress_.lastIterationNumber(0) == numberIterations_) {
     1533          if (!numberPrimalInfeasibilities_) {
     1534            if (numberTimesOptimal_ < 4) {
     1535              numberTimesOptimal_++;
     1536              changeMade_++; // say change made
     1537            } else {
     1538              problemStatus_ = 0;
     1539              secondaryStatus_ = 5;
     1540            }
     1541          }
     1542        }
    15591543      }
    15601544    }
     
    15621546  if (problemStatus_ == 0) {
    15631547    double objVal = (abcNonLinearCost_->feasibleCost()
    1564                      + objective_->nonlinearOffset());
     1548      + objective_->nonlinearOffset());
    15651549    objVal /= (objectiveScale_ * rhsScale_);
    15661550    double tol = 1.0e-10 * CoinMax(fabs(objVal), fabs(objectiveValue_)) + 1.0e-8;
    15671551    if (fabs(objVal - objectiveValue_) > tol) {
    1568 #if ABC_NORMAL_DEBUG>3
     1552#if ABC_NORMAL_DEBUG > 3
    15691553      if (handler_->logLevel() > 0)
    1570         printf("nonLinearCost has feasible obj of %g, objectiveValue_ is %g\n",
    1571                objVal, objectiveValue_);
     1554        printf("nonLinearCost has feasible obj of %g, objectiveValue_ is %g\n",
     1555          objVal, objectiveValue_);
    15721556#endif
    15731557      objectiveValue_ = objVal;
     
    15931577    moreSpecialOptions_ &= ~16; // clear check accuracy flag
    15941578#ifndef TRY_ABC_GUS
    1595   if ((moreSpecialOptions_ & 256) != 0||saveSum||(specialOptions_ & 8192) != 0)
    1596     goToDual=false;
     1579  if ((moreSpecialOptions_ & 256) != 0 || saveSum || (specialOptions_ & 8192) != 0)
     1580    goToDual = false;
    15971581#else
    15981582  if ((specialOptions_ & 8192) != 0)
    1599     goToDual=false;
    1600 #endif
    1601   if (goToDual || (numberIterations_ > 1000+baseIteration_ && largestPrimalError_ > 1.0e6
    1602                    && largestDualError_ > 1.0e6)) {
     1583    goToDual = false;
     1584#endif
     1585  if (goToDual || (numberIterations_ > 1000 + baseIteration_ && largestPrimalError_ > 1.0e6 && largestDualError_ > 1.0e6)) {
    16031586    problemStatus_ = 10; // try dual
    16041587    // See if second call
    1605     if ((moreSpecialOptions_ & 256) != 0||abcNonLinearCost_->sumInfeasibilities()>1.0e20) {
     1588    if ((moreSpecialOptions_ & 256) != 0 || abcNonLinearCost_->sumInfeasibilities() > 1.0e20) {
    16061589      numberPrimalInfeasibilities_ = abcNonLinearCost_->numberInfeasibilities();
    16071590      sumPrimalInfeasibilities_ = abcNonLinearCost_->sumInfeasibilities();
    16081591      // say infeasible
    1609       if (numberPrimalInfeasibilities_&&(abcState_&CLP_ABC_BEEN_FEASIBLE)==0)
    1610         problemStatus_ = 1;
     1592      if (numberPrimalInfeasibilities_ && (abcState_ & CLP_ABC_BEEN_FEASIBLE) == 0)
     1593        problemStatus_ = 1;
    16111594    }
    16121595  }
     
    16171600  }
    16181601#ifdef TRY_SPLIT_VALUES_PASS
    1619   if (valuesChunk>0) {
    1620     bool doFixing=firstFree_<0;
    1621     if (numberIterations_==baseIteration_&&
    1622         numberDualInfeasibilitiesWithoutFree_+1000<numberDualInfeasibilities_) {
    1623       doFixing=true;
    1624       int numberSuperBasic=0;
     1602  if (valuesChunk > 0) {
     1603    bool doFixing = firstFree_ < 0;
     1604    if (numberIterations_ == baseIteration_ && numberDualInfeasibilitiesWithoutFree_ + 1000 < numberDualInfeasibilities_) {
     1605      doFixing = true;
     1606      int numberSuperBasic = 0;
    16251607      for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
    1626         if (getInternalStatus(iSequence) == superBasic)
    1627           numberSuperBasic++;
    1628       }
    1629       keepValuesPass=static_cast<int>((1.01*numberSuperBasic)/valuesChunk);
    1630       doOrdinaryVariables=static_cast<int>(valuesRatio*keepValuesPass);
    1631     } else if (valuesStop>0) {
    1632       if (numberIterations_>=valuesStop||problemStatus_>=0) {
    1633         gutsOfSolution(3);
    1634         abcNonLinearCost_->refreshFromPerturbed(primalTolerance_);
    1635         gutsOfSolution(3);
    1636         int numberSuperBasic=0;
    1637         for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
    1638           if (getInternalStatus(iSequence) == isFixed) {
    1639             if (abcUpper_[iSequence]>abcLower_[iSequence]) {
    1640               setInternalStatus(iSequence,superBasic);
    1641               numberSuperBasic++;
    1642             }
    1643           }
    1644         }
    1645         if (numberSuperBasic) {
    1646           stateOfProblem_ |= VALUES_PASS;
    1647           problemStatus_=-1;
    1648           gutsOfSolution(3);
    1649         } else {
    1650           doFixing=false;
    1651         }
    1652         valuesStop=-1;
     1608        if (getInternalStatus(iSequence) == superBasic)
     1609          numberSuperBasic++;
     1610      }
     1611      keepValuesPass = static_cast< int >((1.01 * numberSuperBasic) / valuesChunk);
     1612      doOrdinaryVariables = static_cast< int >(valuesRatio * keepValuesPass);
     1613    } else if (valuesStop > 0) {
     1614      if (numberIterations_ >= valuesStop || problemStatus_ >= 0) {
     1615        gutsOfSolution(3);
     1616        abcNonLinearCost_->refreshFromPerturbed(primalTolerance_);
     1617        gutsOfSolution(3);
     1618        int numberSuperBasic = 0;
     1619        for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
     1620          if (getInternalStatus(iSequence) == isFixed) {
     1621            if (abcUpper_[iSequence] > abcLower_[iSequence]) {
     1622              setInternalStatus(iSequence, superBasic);
     1623              numberSuperBasic++;
     1624            }
     1625          }
     1626        }
     1627        if (numberSuperBasic) {
     1628          stateOfProblem_ |= VALUES_PASS;
     1629          problemStatus_ = -1;
     1630          gutsOfSolution(3);
     1631        } else {
     1632          doFixing = false;
     1633        }
     1634        valuesStop = -1;
    16531635      } else {
    1654         doFixing=false;
    1655       }
    1656     } 
     1636        doFixing = false;
     1637      }
     1638    }
    16571639    if (doFixing) {
    16581640      abcNonLinearCost_->refreshFromPerturbed(primalTolerance_);
    1659         int numberSuperBasic=0;
     1641      int numberSuperBasic = 0;
    16601642      for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
    1661         if (getInternalStatus(iSequence) == isFixed) {
    1662           if (abcUpper_[iSequence]>abcLower_[iSequence])
    1663             setInternalStatus(iSequence,superBasic);
    1664         }
    1665         if (getInternalStatus(iSequence) == superBasic)
    1666           numberSuperBasic++;
    1667       }
    1668       int numberFixed=0;
    1669       firstFree_=-1;
     1643        if (getInternalStatus(iSequence) == isFixed) {
     1644          if (abcUpper_[iSequence] > abcLower_[iSequence])
     1645            setInternalStatus(iSequence, superBasic);
     1646        }
     1647        if (getInternalStatus(iSequence) == superBasic)
     1648          numberSuperBasic++;
     1649      }
     1650      int numberFixed = 0;
     1651      firstFree_ = -1;
    16701652      if (numberSuperBasic) {
    1671         double threshold = static_cast<double>(keepValuesPass)/numberSuperBasic;
    1672         if (1.1*keepValuesPass>numberSuperBasic)
    1673           threshold=1.1;
    1674         for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
    1675           if (getInternalStatus(iSequence) == superBasic) {
    1676             if (randomNumberGenerator_.randomDouble()<threshold) {
    1677               if (firstFree_<0)
    1678                 firstFree_=iSequence;
    1679             } else {
    1680               numberFixed++;
    1681               abcLower_[iSequence]=abcSolution_[iSequence];
    1682               abcUpper_[iSequence]=abcSolution_[iSequence];
    1683               setInternalStatus(iSequence,isFixed);
    1684             }
    1685           }
    1686         }
     1653        double threshold = static_cast< double >(keepValuesPass) / numberSuperBasic;
     1654        if (1.1 * keepValuesPass > numberSuperBasic)
     1655          threshold = 1.1;
     1656        for (int iSequence = 0; iSequence < numberTotal_; iSequence++) {
     1657          if (getInternalStatus(iSequence) == superBasic) {
     1658            if (randomNumberGenerator_.randomDouble() < threshold) {
     1659              if (firstFree_ < 0)
     1660                firstFree_ = iSequence;
     1661            } else {
     1662              numberFixed++;
     1663              abcLower_[iSequence] = abcSolution_[iSequence];
     1664              abcUpper_[iSequence] = abcSolution_[iSequence];
     1665              setInternalStatus(iSequence, isFixed);
     1666            }
     1667          }
     1668        }
    16871669      }
    16881670      if (!numberFixed) {
    1689         // end
    1690         valuesChunk=-1;
    1691       }
    1692     }
    1693   }
    1694 #endif
    1695   if (doFactorization||doWeights) {
     1671        // end
     1672        valuesChunk = -1;
     1673      }
     1674    }
     1675  }
     1676#endif
     1677  if (doFactorization || doWeights) {
    16961678    // restore weights (if saved) - also recompute infeasibility list
    16971679    if (tentativeStatus > -3)
     
    17021684}
    17031685// Computes solutions - 1 do duals, 2 do primals, 3 both
    1704 int
    1705 AbcSimplex::gutsOfPrimalSolution(int type)
     1686int AbcSimplex::gutsOfPrimalSolution(int type)
    17061687{
    17071688  //work space
    17081689  int whichArray[2];
    1709   for (int i=0;i<2;i++)
    1710     whichArray[i]=getAvailableArray();
    1711   CoinIndexedVector * array1 = &usefulArray_[whichArray[0]];
    1712   CoinIndexedVector * array2 = &usefulArray_[whichArray[1]];
     1690  for (int i = 0; i < 2; i++)
     1691    whichArray[i] = getAvailableArray();
     1692  CoinIndexedVector *array1 = &usefulArray_[whichArray[0]];
     1693  CoinIndexedVector *array2 = &usefulArray_[whichArray[1]];
    17131694  // do work and check
    1714   int numberRefinements=0;
    1715   if ((type&2)!=0) {
    1716     numberRefinements=computePrimals(array1,array2);
     1695  int numberRefinements = 0;
     1696  if ((type & 2) != 0) {
     1697    numberRefinements = computePrimals(array1, array2);
    17171698    if (algorithm_ > 0 && abcNonLinearCost_ != NULL) {
    17181699      // primal algorithm
     
    17201701      abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
    17211702      if (abcNonLinearCost_->numberInfeasibilities())
    1722         if (handler_->detail(CLP_SIMPLEX_NONLINEAR, messages_) < 100) {
    1723           handler_->message(CLP_SIMPLEX_NONLINEAR, messages_)
    1724             << abcNonLinearCost_->changeInCost()
    1725             << abcNonLinearCost_->numberInfeasibilities()
    1726             << CoinMessageEol;
    1727         }
     1703        if (handler_->detail(CLP_SIMPLEX_NONLINEAR, messages_) < 100) {
     1704          handler_->message(CLP_SIMPLEX_NONLINEAR, messages_)
     1705            << abcNonLinearCost_->changeInCost()
     1706            << abcNonLinearCost_->numberInfeasibilities()
     1707            << CoinMessageEol;
     1708        }
    17281709    }
    17291710    checkPrimalSolution(true);
    17301711  }
    1731   if ((type&1)!=0
    1732 #if CAN_HAVE_ZERO_OBJ>1
    1733       &&(specialOptions_&2097152)==0
    1734 #endif
    1735       ) {
    1736     numberRefinements+=computeDuals(NULL,array1,array2);
     1712  if ((type & 1) != 0
     1713#if CAN_HAVE_ZERO_OBJ > 1
     1714    && (specialOptions_ & 2097152) == 0
     1715#endif
     1716  ) {
     1717    numberRefinements += computeDuals(NULL, array1, array2);
    17371718    checkDualSolution();
    17381719  }
    1739   for (int i=0;i<2;i++)
     1720  for (int i = 0; i < 2; i++)
    17401721    setAvailableArray(whichArray[i]);
    1741   rawObjectiveValue_ +=sumNonBasicCosts_;
     1722  rawObjectiveValue_ += sumNonBasicCosts_;
    17421723  //computeObjective(); // ? done in checkDualSolution??
    17431724  setClpSimplexObjectiveValue();
    1744   if (handler_->logLevel() > 3 || (largestPrimalError_ > 1.0e-2 ||
    1745                                    largestDualError_ > 1.0e-2))
     1725  if (handler_->logLevel() > 3 || (largestPrimalError_ > 1.0e-2 || largestDualError_ > 1.0e-2))
    17461726    handler_->message(CLP_SIMPLEX_ACCURACY, messages_)
    17471727      << largestPrimalError_
    17481728      << largestDualError_
    17491729      << CoinMessageEol;
    1750   if ((largestPrimalError_ > 1.0e1||largestDualError_>1.0e1)
    1751       && numberRows_ > 100 && numberIterations_) {
    1752 #if ABC_NORMAL_DEBUG>0
     1730  if ((largestPrimalError_ > 1.0e1 || largestDualError_ > 1.0e1)
     1731    && numberRows_ > 100 && numberIterations_) {
     1732#if ABC_NORMAL_DEBUG > 0
    17531733    printf("Large errors - primal %g dual %g\n",
    1754            largestPrimalError_,largestDualError_);
     1734      largestPrimalError_, largestDualError_);
    17551735#endif
    17561736    // Change factorization tolerance
     
    17671747  On exit rhsArray will have changes in costs of basic variables
    17681748*/
    1769 void
    1770 AbcSimplexPrimal::primalRow(CoinIndexedVector * rowArray,
    1771                             CoinIndexedVector * rhsArray,
    1772                             CoinIndexedVector * spareArray,
    1773                             int valuesPass)
     1749void AbcSimplexPrimal::primalRow(CoinIndexedVector *rowArray,
     1750  CoinIndexedVector *rhsArray,
     1751  CoinIndexedVector *spareArray,
     1752  int valuesPass)
    17741753{
    17751754#if 1
    1776   for (int iRow=0;iRow<numberRows_;iRow++) {
    1777     int iPivot=abcPivotVariable_[iRow];
    1778     assert (costBasic_[iRow]==abcCost_[iPivot]);
    1779     assert (lowerBasic_[iRow]==abcLower_[iPivot]);
    1780     assert (upperBasic_[iRow]==abcUpper_[iPivot]);
    1781     assert (solutionBasic_[iRow]==abcSolution_[iPivot]);
     1755  for (int iRow = 0; iRow < numberRows_; iRow++) {
     1756    int iPivot = abcPivotVariable_[iRow];
     1757    assert(costBasic_[iRow] == abcCost_[iPivot]);
     1758    assert(lowerBasic_[iRow] == abcLower_[iPivot]);
     1759    assert(upperBasic_[iRow] == abcUpper_[iPivot]);
     1760    assert(solutionBasic_[iRow] == abcSolution_[iPivot]);
    17821761  }
    17831762#endif
     
    17851764  if (valuesPass) {
    17861765    dualIn_ = abcCost_[sequenceIn_];
    1787    
    1788     double * work = rowArray->denseVector();
     1766
     1767    double *work = rowArray->denseVector();
    17891768    int number = rowArray->getNumElements();
    1790     int * which = rowArray->getIndices();
    1791    
     1769    int *which = rowArray->getIndices();
     1770
    17921771    int iIndex;
    17931772    for (iIndex = 0; iIndex < number; iIndex++) {
    1794      
     1773
    17951774      int iRow = which[iIndex];
    17961775      double alpha = work[iRow];
     
    18051784      // towards nearest bound
    18061785      if (valueIn_ - lowerIn_ < upperIn_ - valueIn_) {
    1807         directionIn_ = -1;
    1808         dualIn_ = dualTolerance_;
     1786        directionIn_ = -1;
     1787        dualIn_ = dualTolerance_;
    18091788      } else {
    1810         directionIn_ = 1;
    1811         dualIn_ = -dualTolerance_;
    1812       }
    1813     }
    1814   }
    1815  
     1789        directionIn_ = 1;
     1790        dualIn_ = -dualTolerance_;
     1791      }
     1792    }
     1793  }
     1794
    18161795  // sequence stays as row number until end
    18171796  pivotRow_ = -1;
    18181797  int numberRemaining = 0;
    1819  
     1798
    18201799  double totalThru = 0.0; // for when variables flip
    18211800  // Allow first few iterations to take tiny
     
    18331812  double lastPivot = 0.0;
    18341813  double lastTheta = 1.0e50;
    1835  
     1814
    18361815  // use spareArrays to put ones looked at in
    18371816  // First one is list of candidates
     
    18391818  // Second array has current set of pivot candidates
    18401819  // with a backup list saved in double * part of indexed vector
    1841  
     1820
    18421821  // pivot elements
    1843   double * spare;
     1822  double *spare;
    18441823  // indices
    1845   int * index;
     1824  int *index;
    18461825  spareArray->clear();
    18471826  spare = spareArray->denseVector();
    18481827  index = spareArray->getIndices();
    1849  
     1828
    18501829  // we also need somewhere for effective rhs
    1851   double * rhs = rhsArray->denseVector();
     1830  double *rhs = rhsArray->denseVector();
    18521831  // and we can use indices to point to alpha
    18531832  // that way we can store fabs(alpha)
    1854   int * indexPoint = rhsArray->getIndices();
     1833  int *indexPoint = rhsArray->getIndices();
    18551834  //int numberFlip=0; // Those which may change if flips
    1856  
     1835
    18571836  /*
    18581837    First we get a list of possible pivots.  We can also see if the
     
    18641843   
    18651844  */
    1866  
     1845
    18671846  // do first pass to get possibles
    18681847  // We can also see if unbounded
    1869  
    1870   double * work = rowArray->denseVector();
     1848
     1849  double *work = rowArray->denseVector();
    18711850  int number = rowArray->getNumElements();
    1872   int * which = rowArray->getIndices();
    1873  
     1851  int *which = rowArray->getIndices();
     1852
    18741853  // we need to swap sign if coming in from ub
    18751854  double way = directionIn_;
     
    18791858  else
    18801859    maximumMovement = CoinMin(1.0e30, valueIn_ - lowerIn_);
    1881  
     1860
    18821861  double averageTheta = abcNonLinearCost_->averageTheta();
    18831862  double tentativeTheta = CoinMin(10.0 * averageTheta, maximumMovement);
     
    18921871  // but make a bit more pessimistic
    18931872  dualCheck = CoinMax(dualCheck - 100.0 * dualTolerance_, 0.99 * dualCheck);
    1894  
     1873
    18951874  int iIndex;
    18961875  //#define CLP_DEBUG
    1897 #if ABC_NORMAL_DEBUG>3
     1876#if ABC_NORMAL_DEBUG > 3
    18981877  if (numberIterations_ == 1369 || numberIterations_ == -3840) {
    18991878    double dj = abcCost_[sequenceIn_];
    19001879    printf("cost in on %d is %g, dual in %g\n", sequenceIn_, dj, dualIn_);
    19011880    for (iIndex = 0; iIndex < number; iIndex++) {
    1902      
     1881
    19031882      int iRow = which[iIndex];
    19041883      double alpha = work[iRow];
     
    19061885      dj -= alpha * costBasic_[iRow];
    19071886      printf("row %d var %d current %g %g %g, alpha %g so sol => %g (cost %g, dj %g)\n",
    1908              iRow, iPivot, lowerBasic_[iRow], solutionBasic_[iRow], upperBasic_[iRow],
    1909              alpha, solutionBasic_[iRow] - 1.0e9 * alpha, costBasic_[iRow], dj);
     1887        iRow, iPivot, lowerBasic_[iRow], solutionBasic_[iRow], upperBasic_[iRow],
     1888        alpha, solutionBasic_[iRow] - 1.0e9 * alpha, costBasic_[iRow], dj);
    19101889    }
    19111890  }
     
    19201899#endif
    19211900    for (iIndex = 0; iIndex < number; iIndex++) {
    1922      
     1901
    19231902      int iRow = which[iIndex];
    19241903      double alpha = work[iRow];
     
    19311910      // do computation same way as later on in primal
    19321911      if (alpha > 0.0) {
    1933         // basic variable going towards lower bound
    1934         double bound = lowerBasic_[iRow];
    1935         // must be exactly same as when used
    1936         double change = tentativeTheta * alpha;
    1937         possible = (oldValue - change) <= bound + primalTolerance_;
    1938         oldValue -= bound;
     1912        // basic variable going towards lower bound
     1913        double bound = lowerBasic_[iRow];
     1914        // must be exactly same as when used
     1915        double change = tentativeTheta * alpha;
     1916        possible = (oldValue - change) <= bound + primalTolerance_;
     1917        oldValue -= bound;
    19391918      } else {
    1940         // basic variable going towards upper bound
    1941         double bound = upperBasic_[iRow];
    1942         // must be exactly same as when used
    1943         double change = tentativeTheta * alpha;
    1944         possible = (oldValue - change) >= bound - primalTolerance_;
    1945         oldValue = bound - oldValue;
    1946         alpha = - alpha;
     1919        // basic variable going towards upper bound
     1920        double bound = upperBasic_[iRow];
     1921        // must be exactly same as when used
     1922        double change = tentativeTheta * alpha;
     1923        possible = (oldValue - change) >= bound - primalTolerance_;
     1924        oldValue = bound - oldValue;
     1925        alpha = -alpha;
    19471926      }
    19481927      double value;
    19491928      //assert (oldValue >= -10.0*tolerance);
    19501929      if (possible) {
    1951         value = oldValue - upperTheta * alpha;
     1930        value = oldValue - upperTheta * alpha;
    19521931#ifdef CLP_USER_DRIVEN1
    1953         if(!userChoiceValid1(this,iPivot,oldValue,
    1954                              upperTheta,alpha,work[iRow]*way))
    1955           value =0.0; // say can't use
    1956 #endif
    1957         if (value < -primalTolerance_ && alpha >= acceptablePivot) {
    1958           upperTheta = (oldValue + primalTolerance_) / alpha;
    1959         }
    1960         // add to list
    1961         spare[numberRemaining] = alpha;
    1962         rhs[numberRemaining] = oldValue;
    1963         indexPoint[numberRemaining] = iIndex;
    1964         index[numberRemaining++] = iRow;
    1965         totalThru += alpha;
    1966         setActive(iRow);
    1967         //} else if (value<primalTolerance_*1.002) {
    1968         // May change if is a flip
    1969         //indexRhs[numberFlip++]=iRow;
    1970       }
    1971     }
    1972     if (upperTheta < maximumMovement && totalThru*infeasibilityCost_ >= 1.0001 * dualCheck) {
     1932        if (!userChoiceValid1(this, iPivot, oldValue,
     1933              upperTheta, alpha, work[iRow] * way))
     1934          value = 0.0; // say can't use
     1935#endif
     1936        if (value < -primalTolerance_ && alpha >= acceptablePivot) {
     1937          upperTheta = (oldValue + primalTolerance_) / alpha;
     1938        }
     1939        // add to list
     1940        spare[numberRemaining] = alpha;
     1941        rhs[numberRemaining] = oldValue;
     1942        indexPoint[numberRemaining] = iIndex;
     1943        index[numberRemaining++] = iRow;
     1944        totalThru += alpha;
     1945        setActive(iRow);
     1946        //} else if (value<primalTolerance_*1.002) {
     1947        // May change if is a flip
     1948        //indexRhs[numberFlip++]=iRow;
     1949      }
     1950    }
     1951    if (upperTheta < maximumMovement && totalThru * infeasibilityCost_ >= 1.0001 * dualCheck) {
    19731952      // Can pivot here
    19741953      break;
     
    19821961  }
    19831962  totalThru = 0.0;
    1984  
     1963
    19851964  theta_ = maximumMovement;
    1986  
     1965
    19871966  bool goBackOne = false;
    19881967  if (objective_->type() > 1)
    19891968    dualIn_ = saveDj;
    1990  
     1969
    19911970  //printf("%d remain out of %d\n",numberRemaining,number);
    19921971  int iTry = 0;
     
    19941973  if (numberRemaining && upperTheta < maximumMovement) {
    19951974    // First check if previously chosen one will work
    1996    
     1975
    19971976    // first get ratio with tolerance
    1998     for ( ; iTry < MAXTRY; iTry++) {
    1999      
     1977    for (; iTry < MAXTRY; iTry++) {
     1978
    20001979      upperTheta = maximumMovement;
    20011980      int iBest = -1;
    20021981      for (iIndex = 0; iIndex < numberRemaining; iIndex++) {
    2003        
    2004         double alpha = spare[iIndex];
    2005         double oldValue = rhs[iIndex];
    2006         double value = oldValue - upperTheta * alpha;
    2007        
     1982
     1983        double alpha = spare[iIndex];
     1984        double oldValue = rhs[iIndex];
     1985        double value = oldValue - upperTheta * alpha;
     1986
    20081987#ifdef CLP_USER_DRIVEN1
    2009         int sequenceOut=abcPivotVariable_[index[iIndex]];
    2010         if(!userChoiceValid1(this,sequenceOut,oldValue,
    2011                              upperTheta,alpha, 0.0))
    2012           value =0.0; // say can't use
    2013 #endif
    2014         if (value < -primalTolerance_ && alpha >= acceptablePivot) {
    2015           upperTheta = (oldValue + primalTolerance_) / alpha;
    2016           iBest = iIndex; // just in case weird numbers
    2017         }
    2018       }
    2019      
     1988        int sequenceOut = abcPivotVariable_[index[iIndex]];
     1989        if (!userChoiceValid1(this, sequenceOut, oldValue,
     1990              upperTheta, alpha, 0.0))
     1991          value = 0.0; // say can't use
     1992#endif
     1993        if (value < -primalTolerance_ && alpha >= acceptablePivot) {
     1994          upperTheta = (oldValue + primalTolerance_) / alpha;
     1995          iBest = iIndex; // just in case weird numbers
     1996        }
     1997      }
     1998
    20201999      // now look at best in this lot
    20212000      // But also see how infeasible small pivots will make
     
    20242003      pivotRow_ = -1;
    20252004      for (iIndex = 0; iIndex < numberRemaining; iIndex++) {
    2026        
    2027         int iRow = index[iIndex];
    2028         double alpha = spare[iIndex];
    2029         double oldValue = rhs[iIndex];
    2030         double value = oldValue - upperTheta * alpha;
    2031        
    2032         if (value <= 0 || iBest == iIndex) {
    2033           // how much would it cost to go thru and modify bound
    2034           double trueAlpha = way * work[iRow];
    2035           totalThru += abcNonLinearCost_->changeInCost(iRow, trueAlpha, rhs[iIndex]);
    2036           setActive(iRow);
    2037           if (alpha > bestPivot) {
    2038             bestPivot = alpha;
    2039             theta_ = oldValue / bestPivot;
    2040             pivotRow_ = iRow;
    2041           } else if (alpha < acceptablePivot
     2005
     2006        int iRow = index[iIndex];
     2007        double alpha = spare[iIndex];
     2008        double oldValue = rhs[iIndex];
     2009        double value = oldValue - upperTheta * alpha;
     2010
     2011        if (value <= 0 || iBest == iIndex) {
     2012          // how much would it cost to go thru and modify bound
     2013          double trueAlpha = way * work[iRow];
     2014          totalThru += abcNonLinearCost_->changeInCost(iRow, trueAlpha, rhs[iIndex]);
     2015          setActive(iRow);
     2016          if (alpha > bestPivot) {
     2017            bestPivot = alpha;
     2018            theta_ = oldValue / bestPivot;
     2019            pivotRow_ = iRow;
     2020          } else if (alpha < acceptablePivot
    20422021#ifdef CLP_USER_DRIVEN1
    2043                      ||!userChoiceValid1(this,abcPivotVariable_[index[iIndex]],
    2044                                          oldValue,upperTheta,alpha,0.0)
    2045 #endif
    2046                      ) {
    2047             if (value < -primalTolerance_)
    2048               sumInfeasibilities += -value - primalTolerance_;
    2049           }
    2050         }
    2051       }
    2052       if (bestPivot < 0.1 * bestEverPivot &&
    2053           bestEverPivot > 1.0e-6 && bestPivot < 1.0e-3) {
    2054         // back to previous one
    2055         goBackOne = true;
    2056         break;
     2022            || !userChoiceValid1(this, abcPivotVariable_[index[iIndex]],
     2023                 oldValue, upperTheta, alpha, 0.0)
     2024#endif
     2025          ) {
     2026            if (value < -primalTolerance_)
     2027              sumInfeasibilities += -value - primalTolerance_;
     2028          }
     2029        }
     2030      }
     2031      if (bestPivot < 0.1 * bestEverPivot && bestEverPivot > 1.0e-6 && bestPivot < 1.0e-3) {
     2032        // back to previous one
     2033        goBackOne = true;
     2034        break;
    20572035      } else if (pivotRow_ == -1 && upperTheta > largeValue_) {
    2058         if (lastPivot > acceptablePivot) {
    2059           // back to previous one
    2060           goBackOne = true;
    2061         } else {
    2062           // can only get here if all pivots so far too small
    2063         }
    2064         break;
     2036        if (lastPivot > acceptablePivot) {
     2037          // back to previous one
     2038          goBackOne = true;
     2039        } else {
     2040          // can only get here if all pivots so far too small
     2041        }
     2042        break;
    20652043      } else if (totalThru >= dualCheck) {
    2066         if (sumInfeasibilities > primalTolerance_ && !abcNonLinearCost_->numberInfeasibilities()) {
    2067           // Looks a bad choice
    2068           if (lastPivot > acceptablePivot) {
    2069             goBackOne = true;
    2070           } else {
    2071             // say no good
    2072             dualIn_ = 0.0;
    2073           }
    2074         }
    2075         break; // no point trying another loop
     2044        if (sumInfeasibilities > primalTolerance_ && !abcNonLinearCost_->numberInfeasibilities()) {
     2045          // Looks a bad choice
     2046          if (lastPivot > acceptablePivot) {
     2047            goBackOne = true;
     2048          } else {
     2049            // say no good
     2050            dualIn_ = 0.0;
     2051          }
     2052        }
     2053        break; // no point trying another loop
    20762054      } else {
    2077         lastPivotRow = pivotRow_;
    2078         lastTheta = theta_;
    2079         if (bestPivot > bestEverPivot)
    2080           bestEverPivot = bestPivot;
     2055        lastPivotRow = pivotRow_;
     2056        lastTheta = theta_;
     2057        if (bestPivot > bestEverPivot)
     2058          bestEverPivot = bestPivot;
    20812059      }
    20822060    }
     
    20952073      tentativeTheta = 1.0e50;
    20962074      for (iIndex = 0; iIndex < number; iIndex++) {
    2097         int iRow = which[iIndex];
    2098         double alpha = work[iRow];
    2099         alpha *= way;
    2100         double oldValue = solutionBasic_[iRow];
    2101         // get where in bound sequence
    2102         // note that after this alpha is actually fabs(alpha)
    2103         if (alpha > 0.0) {
    2104           // basic variable going towards lower bound
    2105           double bound = lowerBasic_[iRow];
    2106           oldValue -= bound;
    2107         } else {
    2108           // basic variable going towards upper bound
    2109           double bound = upperBasic_[iRow];
    2110           oldValue = bound - oldValue;
    2111           alpha = - alpha;
    2112         }
    2113         if (oldValue - tentativeTheta * alpha < 0.0) {
    2114           tentativeTheta = oldValue / alpha;
    2115         }
     2075        int iRow = which[iIndex];
     2076        double alpha = work[iRow];
     2077        alpha *= way;
     2078        double oldValue = solutionBasic_[iRow];
     2079        // get where in bound sequence
     2080        // note that after this alpha is actually fabs(alpha)
     2081        if (alpha > 0.0) {
     2082          // basic variable going towards lower bound
     2083          double bound = lowerBasic_[iRow];
     2084          oldValue -= bound;
     2085        } else {
     2086          // basic variable going towards upper bound
     2087          double bound = upperBasic_[iRow];
     2088          oldValue = bound - oldValue;
     2089          alpha = -alpha;
     2090        }
     2091        if (oldValue - tentativeTheta * alpha < 0.0) {
     2092          tentativeTheta = oldValue / alpha;
     2093        }
    21162094      }
    21172095      // If free in then see if we can get to 0.0
    21182096      if (lowerIn_ < -1.0e20 && upperIn_ > 1.0e20) {
    2119         if (dualIn_ * valueIn_ > 0.0) {
    2120           if (fabs(valueIn_) < 1.0e-2 && (tentativeTheta < fabs(valueIn_) || tentativeTheta > 1.0e20)) {
    2121             tentativeTheta = fabs(valueIn_);
    2122           }
    2123         }
     2097        if (dualIn_ * valueIn_ > 0.0) {
     2098          if (fabs(valueIn_) < 1.0e-2 && (tentativeTheta < fabs(valueIn_) || tentativeTheta > 1.0e20)) {
     2099            tentativeTheta = fabs(valueIn_);
     2100          }
     2101        }
    21242102      }
    21252103      if (tentativeTheta < 1.0e10)
    2126         valueOut_ = valueIn_ + way * tentativeTheta;
     2104        valueOut_ = valueIn_ + way * tentativeTheta;
    21272105    }
    21282106  }
     
    21582136      theta_ = minimumTheta;
    21592137      for (iIndex = 0; iIndex < numberRemaining - numberRemaining; iIndex++) {
    2160         largestInfeasibility = CoinMax(largestInfeasibility,
    2161                                        -(rhs[iIndex] - spare[iIndex] * theta_));
     2138        largestInfeasibility = CoinMax(largestInfeasibility,
     2139          -(rhs[iIndex] - spare[iIndex] * theta_));
    21622140      }
    21632141      //#define CLP_DEBUG
    2164 #if ABC_NORMAL_DEBUG>3
     2142#if ABC_NORMAL_DEBUG > 3
    21652143      if (largestInfeasibility > primalTolerance_ && (handler_->logLevel() & 32) > -1)
    2166         printf("Primal tolerance increased from %g to %g\n",
    2167                primalTolerance_, largestInfeasibility);
     2144        printf("Primal tolerance increased from %g to %g\n",
     2145          primalTolerance_, largestInfeasibility);
    21682146#endif
    21692147      //#undef CLP_DEBUG
     
    21732151    if (theta_ > tentativeTheta) {
    21742152      for (iIndex = 0; iIndex < number; iIndex++)
    2175         setActive(which[iIndex]);
     2153        setActive(which[iIndex]);
    21762154    }
    21772155    if (way < 0.0)
    2178       theta_ = - theta_;
     2156      theta_ = -theta_;
    21792157    double newValue = valueOut_ - theta_ * alpha_;
    21802158    // If  4 bit set - Force outgoing variables to exact bound (primal)
    21812159    if (alpha_ * way < 0.0) {
    2182       directionOut_ = -1;    // to upper bound
     2160      directionOut_ = -1; // to upper bound
    21832161      if (fabs(theta_) > 1.0e-6 || (specialOptions_ & 4) != 0) {
    2184         upperOut_ = abcNonLinearCost_->nearest(pivotRow_, newValue);
     2162        upperOut_ = abcNonLinearCost_->nearest(pivotRow_, newValue);
    21852163      } else {
    2186         upperOut_ = newValue;
     2164        upperOut_ = newValue;
    21872165      }
    21882166    } else {
    2189       directionOut_ = 1;    // to lower bound
     2167      directionOut_ = 1; // to lower bound
    21902168      if (fabs(theta_) > 1.0e-6 || (specialOptions_ & 4) != 0) {
    2191         lowerOut_ = abcNonLinearCost_->nearest(pivotRow_, newValue);
     2169        lowerOut_ = abcNonLinearCost_->nearest(pivotRow_, newValue);
    21922170      } else {
    2193         lowerOut_ = newValue;
     2171        lowerOut_ = newValue;
    21942172      }
    21952173    }
     
    22052183    alpha_ = 0.0;
    22062184    if (way < 0.0) {
    2207       directionOut_ = 1;    // to lower bound
     2185      directionOut_ = 1; // to lower bound
    22082186      theta_ = lowerOut_ - valueOut_;
    22092187    } else {
    2210       directionOut_ = -1;    // to upper bound
     2188      directionOut_ = -1; // to upper bound
    22112189      theta_ = upperOut_ - valueOut_;
    22122190    }
    22132191  }
    2214  
     2192
    22152193  double theta1 = CoinMax(theta_, 1.0e-12);
    22162194  double theta2 = numberIterations_ * abcNonLinearCost_->averageTheta();
    22172195  // Set average theta
    2218   abcNonLinearCost_->setAverageTheta((theta1 + theta2) / (static_cast<double> (numberIterations_ + 1)));
     2196  abcNonLinearCost_->setAverageTheta((theta1 + theta2) / (static_cast< double >(numberIterations_ + 1)));
    22192197  //if (numberIterations_%1000==0)
    22202198  //printf("average theta is %g\n",abcNonLinearCost_->averageTheta());
    2221  
     2199
    22222200  // clear arrays
    22232201  CoinZeroN(spare, numberRemaining);
    22242202  CoinZeroN(rhs, numberRemaining);
    2225  
     2203
    22262204  // put back original bounds etc
    22272205  CoinMemcpyN(index, numberRemaining, rhsArray->getIndices());
     
    22302208  abcNonLinearCost_->goBackAll(rhsArray);
    22312209  rhsArray->clear();
    2232  
    22332210}
    22342211/*
     
    22392216  On exit rhsArray will have changes in costs of basic variables
    22402217*/
    2241 void
    2242 AbcSimplexPrimal::primalRow(CoinIndexedVector * rowArray,
    2243                             CoinIndexedVector * rhsArray,
    2244                             CoinIndexedVector * spareArray,
    2245                             pivotStruct & stuff)
     2218void AbcSimplexPrimal::primalRow(CoinIndexedVector *rowArray,
     2219  CoinIndexedVector *rhsArray,
     2220  CoinIndexedVector *spareArray,
     2221  pivotStruct &stuff)
    22462222{
    2247   int valuesPass=stuff.valuesPass_;
    2248   double dualIn=stuff.dualIn_;
    2249   double lowerIn=stuff.lowerIn_;
    2250   double upperIn=stuff.upperIn_;
    2251   double valueIn=stuff.valueIn_;
    2252   int sequenceIn=stuff.sequenceIn_;
    2253   int directionIn=stuff.directionIn_;
    2254   int pivotRow=-1;
    2255   int sequenceOut=-1;
    2256   double dualOut=0.0;
    2257   double lowerOut=0.0;
    2258   double upperOut=0.0;
    2259   double valueOut=0.0;
     2223  int valuesPass = stuff.valuesPass_;
     2224  double dualIn = stuff.dualIn_;
     2225  double lowerIn = stuff.lowerIn_;
     2226  double upperIn = stuff.upperIn_;
     2227  double valueIn = stuff.valueIn_;
     2228  int sequenceIn = stuff.sequenceIn_;
     2229  int directionIn = stuff.directionIn_;
     2230  int pivotRow = -1;
     2231  int sequenceOut = -1;
     2232  double dualOut = 0.0;
     2233  double lowerOut = 0.0;
     2234  double upperOut = 0.0;
     2235  double valueOut = 0.0;
    22602236  double theta;
    22612237  double alpha;
    2262   int directionOut=0;
     2238  int directionOut = 0;
    22632239  if (valuesPass) {
    22642240    dualIn = abcCost_[sequenceIn];
    2265    
    2266     double * work = rowArray->denseVector();
     2241
     2242    double *work = rowArray->denseVector();
    22672243    int number = rowArray->getNumElements();
    2268     int * which = rowArray->getIndices();
    2269    
     2244    int *which = rowArray->getIndices();
     2245
    22702246    int iIndex;
    22712247    for (iIndex = 0; iIndex < number; iIndex++) {
    2272      
     2248
    22732249      int iRow = which[iIndex];
    22742250      double alpha = work[iRow];
     
    22832259      // towards nearest bound
    22842260      if (valueIn - lowerIn < upperIn - valueIn) {
    2285         directionIn = -1;
    2286         dualIn = dualTolerance_;
     2261        directionIn = -1;
     2262        dualIn = dualTolerance_;
    22872263      } else {
    2288         directionIn = 1;
    2289         dualIn = -dualTolerance_;
    2290       }
    2291     }
    2292   }
    2293  
     2264        directionIn = 1;
     2265        dualIn = -dualTolerance_;
     2266      }
     2267    }
     2268  }
     2269
    22942270  // sequence stays as row number until end
    22952271  pivotRow = -1;
    22962272  int numberRemaining = 0;
    2297  
     2273
    22982274  double totalThru = 0.0; // for when variables flip
    22992275  // Allow first few iterations to take tiny
     
    23112287  double lastPivot = 0.0;
    23122288  double lastTheta = 1.0e50;
    2313  
     2289
    23142290  // use spareArrays to put ones looked at in
    23152291  // First one is list of candidates
     
    23172293  // Second array has current set of pivot candidates
    23182294  // with a backup list saved in double * part of indexed vector
    2319  
     2295
    23202296  // pivot elements
    2321   double * spare;
     2297  double *spare;
    23222298  // indices
    2323   int * index;
     2299  int *index;
    23242300  spareArray->clear();
    23252301  spare = spareArray->denseVector();
    23262302  index = spareArray->getIndices();
    2327  
     2303
    23282304  // we also need somewhere for effective rhs
    2329   double * rhs = rhsArray->denseVector();
     2305  double *rhs = rhsArray->denseVector();
    23302306  // and we can use indices to point to alpha
    23312307  // that way we can store fabs(alpha)
    2332   int * indexPoint = rhsArray->getIndices();
     2308  int *indexPoint = rhsArray->getIndices();
    23332309  //int numberFlip=0; // Those which may change if flips
    2334  
     2310
    23352311  /*
    23362312    First we get a list of possible pivots.  We can also see if the
     
    23422318   
    23432319  */
    2344  
     2320
    23452321  // do first pass to get possibles
    23462322  // We can also see if unbounded
    2347  
    2348   double * work = rowArray->denseVector();
     2323
     2324  double *work = rowArray->denseVector();
    23492325  int number = rowArray->getNumElements();
    2350   int * which = rowArray->getIndices();
    2351  
     2326  int *which = rowArray->getIndices();
     2327
    23522328  // we need to swap sign if coming in from ub
    23532329  double way = directionIn;
     
    23572333  else
    23582334    maximumMovement = CoinMin(1.0e30, valueIn - lowerIn);
    2359  
     2335
    23602336  double averageTheta = abcNonLinearCost_->averageTheta();
    23612337  double tentativeTheta = CoinMin(10.0 * averageTheta, maximumMovement);
     
    23702346  // but make a bit more pessimistic
    23712347  dualCheck = CoinMax(dualCheck - 100.0 * dualTolerance_, 0.99 * dualCheck);
    2372  
     2348
    23732349  int iIndex;
    23742350  while (true) {
     
    23882364      // do computation same way as later on in primal
    23892365      if (alpha > 0.0) {
    2390         // basic variable going towards lower bound
    2391         double bound = lowerBasic_[iRow];
    2392         // must be exactly same as when used
    2393         double change = tentativeTheta * alpha;
    2394         possible = (oldValue - change) <= bound + primalTolerance_;
    2395         oldValue -= bound;
     2366        // basic variable going towards lower bound
     2367        double bound = lowerBasic_[iRow];
     2368        // must be exactly same as when used
     2369        double change = tentativeTheta * alpha;
     2370        possible = (oldValue - change) <= bound + primalTolerance_;
     2371        oldValue -= bound;
    23962372      } else {
    2397         // basic variable going towards upper bound
    2398         double bound = upperBasic_[iRow];
    2399         // must be exactly same as when used
    2400         double change = tentativeTheta * alpha;
    2401         possible = (oldValue - change) >= bound - primalTolerance_;
    2402         oldValue = bound - oldValue;
    2403         alpha = - alpha;
     2373        // basic variable going towards upper bound
     2374        double bound = upperBasic_[iRow];
     2375        // must be exactly same as when used
     2376        double change = tentativeTheta * alpha;
     2377        possible = (oldValue - change) >= bound - primalTolerance_;
     2378        oldValue = bound - oldValue;
     2379        alpha = -alpha;
    24042380      }
    24052381      double value;
    24062382      if (possible) {
    2407         value = oldValue - upperTheta * alpha;
     2383        value = oldValue - upperTheta * alpha;
    24082384#ifdef CLP_USER_DRIVEN1
    2409         if(!userChoiceValid1(this,iPivot,oldValue,
    2410                              upperTheta,alpha,work[iRow]*way))
    2411           value =0.0; // say can't use
    2412 #endif
    2413         if (value < -primalTolerance_ && alpha >= acceptablePivot) {
    2414           upperTheta = (oldValue + primalTolerance_) / alpha;
    2415         }
    2416         // add to list
    2417         spare[numberRemaining] = alpha;
    2418         rhs[numberRemaining] = oldValue;
    2419         indexPoint[numberRemaining] = iIndex;
    2420         index[numberRemaining++] = iRow;
    2421         totalThru += alpha;
    2422         setActive(iRow);
    2423       }
    2424     }
    2425     if (upperTheta < maximumMovement && totalThru*infeasibilityCost_ >= 1.0001 * dualCheck) {
     2385        if (!userChoiceValid1(this, iPivot, oldValue,
     2386              upperTheta, alpha, work[iRow] * way))
     2387          value = 0.0; // say can't use
     2388#endif
     2389        if (value < -primalTolerance_ && alpha >= acceptablePivot) {
     2390          upperTheta = (oldValue + primalTolerance_) / alpha;
     2391        }
     2392        // add to list
     2393        spare[numberRemaining] = alpha;
     2394        rhs[numberRemaining] = oldValue;
     2395        indexPoint[numberRemaining] = iIndex;
     2396        index[numberRemaining++] = iRow;
     2397        totalThru += alpha;
     2398        setActive(iRow);
     2399      }
     2400    }
     2401    if (upperTheta < maximumMovement && totalThru * infeasibilityCost_ >= 1.0001 * dualCheck) {
    24262402      // Can pivot here
    24272403      break;
     
    24352411  }
    24362412  totalThru = 0.0;
    2437  
     2413
    24382414  theta = maximumMovement;
    2439  
     2415
    24402416  bool goBackOne = false;
    2441  
     2417
    24422418  //printf("%d remain out of %d\n",numberRemaining,number);
    24432419  int iTry = 0;
    24442420  if (numberRemaining && upperTheta < maximumMovement) {
    24452421    // First check if previously chosen one will work
    2446    
     2422
    24472423    // first get ratio with tolerance
    2448     for ( ; iTry < MAXTRY; iTry++) {
    2449      
     2424    for (; iTry < MAXTRY; iTry++) {
     2425
    24502426      upperTheta = maximumMovement;
    24512427      int iBest = -1;
    24522428      for (iIndex = 0; iIndex < numberRemaining; iIndex++) {
    2453        
    2454         double alpha = spare[iIndex];
    2455         double oldValue = rhs[iIndex];
    2456         double value = oldValue - upperTheta * alpha;
    2457        
     2429
     2430        double alpha = spare[iIndex];
     2431        double oldValue = rhs[iIndex];
     2432        double value = oldValue - upperTheta * alpha;
     2433
    24582434#ifdef CLP_USER_DRIVEN1
    2459         int sequenceOut=abcPivotVariable_[index[iIndex]];
    2460         if(!userChoiceValid1(this,sequenceOut,oldValue,
    2461                              upperTheta,alpha, 0.0))
    2462           value =0.0; // say can't use
    2463 #endif
    2464         if (value < -primalTolerance_ && alpha >= acceptablePivot) {
    2465           upperTheta = (oldValue + primalTolerance_) / alpha;
    2466           iBest = iIndex; // just in case weird numbers
    2467         }
    2468       }
    2469      
     2435        int sequenceOut = abcPivotVariable_[index[iIndex]];
     2436        if (!userChoiceValid1(this, sequenceOut, oldValue,
     2437              upperTheta, alpha, 0.0))
     2438          value = 0.0; // say can't use
     2439#endif
     2440        if (value < -primalTolerance_ && alpha >= acceptablePivot) {
     2441          upperTheta = (oldValue + primalTolerance_) / alpha;
     2442          iBest = iIndex; // just in case weird numbers
     2443        }
     2444      }
     2445
    24702446      // now look at best in this lot
    24712447      // But also see how infeasible small pivots will make
     
    24742450      pivotRow = -1;
    24752451      for (iIndex = 0; iIndex < numberRemaining; iIndex++) {
    2476        
    2477         int iRow = index[iIndex];
    2478         double alpha = spare[iIndex];
    2479         double oldValue = rhs[iIndex];
    2480         double value = oldValue - upperTheta * alpha;
    2481        
    2482         if (value <= 0 || iBest == iIndex) {
    2483           // how much would it cost to go thru and modify bound
    2484           double trueAlpha = way * work[iRow];
    2485           totalThru += abcNonLinearCost_->changeInCost(iRow, trueAlpha, rhs[iIndex]);
    2486           setActive(iRow);
    2487           if (alpha > bestPivot) {
    2488             bestPivot = alpha;
    2489             theta = oldValue / bestPivot;
    2490             pivotRow = iRow;
    2491           } else if (alpha < acceptablePivot
     2452
     2453        int iRow = index[iIndex];
     2454        double alpha = spare[iIndex];
     2455        double oldValue = rhs[iIndex];
     2456        double value = oldValue - upperTheta * alpha;
     2457
     2458        if (value <= 0 || iBest == iIndex) {
     2459          // how much would it cost to go thru and modify bound
     2460          double trueAlpha = way * work[iRow];
     2461          totalThru += abcNonLinearCost_->changeInCost(iRow, trueAlpha, rhs[iIndex]);
     2462          setActive(iRow);
     2463          if (alpha > bestPivot) {
     2464            bestPivot = alpha;
     2465            theta = oldValue / bestPivot;
     2466            pivotRow = iRow;
     2467          } else if (alpha < acceptablePivot
    24922468#ifdef CLP_USER_DRIVEN1
    2493                      ||!userChoiceValid1(this,abcPivotVariable_[index[iIndex]],
    2494                                          oldValue,upperTheta,alpha,0.0)
    2495 #endif
    2496                      ) {
    2497             if (value < -primalTolerance_)
    2498               sumInfeasibilities += -value - primalTolerance_;
    2499           }
    2500         }
    2501       }
    2502       if (bestPivot < 0.1 * bestEverPivot &&
    2503           bestEverPivot > 1.0e-6 && bestPivot < 1.0e-3) {
    2504         // back to previous one
    2505         goBackOne = true;
    2506         break;
     2469            || !userChoiceValid1(this, abcPivotVariable_[index[iIndex]],
     2470                 oldValue, upperTheta, alpha, 0.0)
     2471#endif
     2472          ) {
     2473            if (value < -primalTolerance_)
     2474              sumInfeasibilities += -value - primalTolerance_;
     2475          }
     2476        }
     2477      }
     2478      if (bestPivot < 0.1 * bestEverPivot && bestEverPivot > 1.0e-6 && bestPivot < 1.0e-3) {
     2479        // back to previous one
     2480        goBackOne = true;
     2481        break;
    25072482      } else if (pivotRow == -1 && upperTheta > largeValue_) {
    2508         if (lastPivot > acceptablePivot) {
    2509           // back to previous one
    2510           goBackOne = true;
    2511         } else {
    2512           // can only get here if all pivots so far too small
    2513         }
    2514         break;
     2483        if (lastPivot > acceptablePivot) {
     2484          // back to previous one
     2485          goBackOne = true;
     2486        } else {
     2487          // can only get here if all pivots so far too small
     2488        }
     2489        break;
    25152490      } else if (totalThru >= dualCheck) {
    2516         if (sumInfeasibilities > primalTolerance_ && !abcNonLinearCost_->numberInfeasibilities()) {
    2517           // Looks a bad choice
    2518           if (lastPivot > acceptablePivot) {
    2519             goBackOne = true;
    2520           } else {
    2521             // say no good
    2522             dualIn = 0.0;
    2523           }
    2524         }
    2525         break; // no point trying another loop
     2491        if (sumInfeasibilities > primalTolerance_ && !abcNonLinearCost_->numberInfeasibilities()) {
     2492          // Looks a bad choice
     2493          if (lastPivot > acceptablePivot) {
     2494            goBackOne = true;
     2495          } else {
     2496            // say no good
     2497            dualIn = 0.0;
     2498          }
     2499        }
     2500        break; // no point trying another loop
    25262501      } else {
    2527         lastPivotRow = pivotRow;
    2528         lastTheta = theta;
    2529         if (bestPivot > bestEverPivot)
    2530           bestEverPivot = bestPivot;
     2502        lastPivotRow = pivotRow;
     2503        lastTheta = theta;
     2504        if (bestPivot > bestEverPivot)
     2505          bestEverPivot = bestPivot;
    25312506      }
    25322507    }
     
    25452520      tentativeTheta = 1.0e50;
    25462521      for (iIndex = 0; iIndex < number; iIndex++) {
    2547         int iRow = which[iIndex];
    2548         double alpha = work[iRow];
    2549         alpha *= way;
    2550         double oldValue = solutionBasic_[iRow];
    2551         // get where in bound sequence
    2552         // note that after this alpha is actually fabs(alpha)
    2553         if (alpha > 0.0) {
    2554           // basic variable going towards lower bound
    2555           double bound = lowerBasic_[iRow];
    2556           oldValue -= bound;
    2557         } else {
    2558           // basic variable going towards upper bound
    2559           double bound = upperBasic_[iRow];
    2560           oldValue = bound - oldValue;
    2561           alpha = - alpha;
    2562         }
    2563         if (oldValue - tentativeTheta * alpha < 0.0) {
    2564           tentativeTheta = oldValue / alpha;
    2565         }
     2522        int iRow = which[iIndex];
     2523        double alpha = work[iRow];
     2524        alpha *= way;
     2525        double oldValue = solutionBasic_[iRow];
     2526        // get where in bound sequence
     2527        // note that after this alpha is actually fabs(alpha)
     2528        if (alpha > 0.0) {
     2529          // basic variable going towards lower bound
     2530          double bound = lowerBasic_[iRow];
     2531          oldValue -= bound;
     2532        } else {
     2533          // basic variable going towards upper bound
     2534          double bound = upperBasic_[iRow];
     2535          oldValue = bound - oldValue;
     2536          alpha = -alpha;
     2537        }
     2538        if (oldValue - tentativeTheta * alpha < 0.0) {
     2539          tentativeTheta = oldValue / alpha;
     2540        }
    25662541      }
    25672542      // If free in then see if we can get to 0.0
    25682543      if (lowerIn < -1.0e20 && upperIn > 1.0e20) {
    2569         if (dualIn * valueIn > 0.0) {
    2570           if (fabs(valueIn) < 1.0e-2 && (tentativeTheta < fabs(valueIn) || tentativeTheta > 1.0e20)) {
    2571             tentativeTheta = fabs(valueIn);
    2572           }
    2573         }
     2544        if (dualIn * valueIn > 0.0) {
     2545          if (fabs(valueIn) < 1.0e-2 && (tentativeTheta < fabs(valueIn) || tentativeTheta > 1.0e20)) {
     2546            tentativeTheta = fabs(valueIn);
     2547          }
     2548        }
    25742549      }
    25752550      if (tentativeTheta < 1.0e10)
    2576         valueOut = valueIn + way * tentativeTheta;
     2551        valueOut = valueIn + way * tentativeTheta;
    25772552    }
    25782553  }
     
    26122587    if (theta > tentativeTheta) {
    26132588      for (iIndex = 0; iIndex < number; iIndex++)
    2614         setActive(which[iIndex]);
     2589        setActive(which[iIndex]);
    26152590    }
    26162591    if (way < 0.0)
    2617       theta = - theta;
     2592      theta = -theta;
    26182593    double newValue = valueOut - theta * alpha;
    26192594    // If  4 bit set - Force outgoing variables to exact bound (primal)
    26202595    if (alpha * way < 0.0) {
    2621       directionOut = -1;    // to upper bound
     2596      directionOut = -1; // to upper bound
    26222597      if (fabs(theta) > 1.0e-6 || (specialOptions_ & 4) != 0) {
    2623         upperOut = abcNonLinearCost_->nearest(pivotRow, newValue);
     2598        upperOut = abcNonLinearCost_->nearest(pivotRow, newValue);
    26242599      } else {
    2625         upperOut = newValue;
     2600        upperOut = newValue;
    26262601      }
    26272602    } else {
    2628       directionOut = 1;    // to lower bound
     2603      directionOut = 1; // to lower bound
    26292604      if (fabs(theta) > 1.0e-6 || (specialOptions_ & 4) != 0) {
    2630         lowerOut = abcNonLinearCost_->nearest(pivotRow, newValue);
     2605        lowerOut = abcNonLinearCost_->nearest(pivotRow, newValue);
    26312606      } else {
    2632         lowerOut = newValue;
     2607        lowerOut = newValue;
    26332608      }
    26342609    }
     
    26442619    alpha = 0.0;
    26452620    if (way < 0.0) {
    2646       directionOut = 1;    // to lower bound
     2621      directionOut = 1; // to lower bound
    26472622      theta = lowerOut - valueOut;
    26482623    } else {
    2649       directionOut = -1;    // to upper bound
     2624      directionOut = -1; // to upper bound
    26502625      theta = upperOut - valueOut;
    26512626    }
    26522627  }
    2653  
    2654  
     2628
    26552629  // clear arrays
    26562630  CoinZeroN(spare, numberRemaining);
    26572631  CoinZeroN(rhs, numberRemaining);
    2658  
     2632
    26592633  // put back original bounds etc
    26602634  CoinMemcpyN(index, numberRemaining, rhsArray->getIndices());
     
    26622636  abcNonLinearCost_->goBackAll(rhsArray);
    26632637  rhsArray->clear();
    2664   stuff.theta_=theta;
    2665   stuff.alpha_=alpha;
    2666   stuff.dualIn_=dualIn;
    2667   stuff.dualOut_=dualOut;
    2668   stuff.lowerOut_=lowerOut;
    2669   stuff.upperOut_=upperOut;
    2670   stuff.valueOut_=valueOut;
    2671   stuff.sequenceOut_=sequenceOut;
    2672   stuff.directionOut_=directionOut;
    2673   stuff.pivotRow_=pivotRow;
     2638  stuff.theta_ = theta;
     2639  stuff.alpha_ = alpha;
     2640  stuff.dualIn_ = dualIn;
     2641  stuff.dualOut_ = dualOut;
     2642  stuff.lowerOut_ = lowerOut;
     2643  stuff.upperOut_ = upperOut;
     2644  stuff.valueOut_ = valueOut;
     2645  stuff.sequenceOut_ = sequenceOut;
     2646  stuff.directionOut_ = directionOut;
     2647  stuff.pivotRow_ = pivotRow;
    26742648}
    26752649/*
     
    26802654  For easy problems we can just choose one of the first columns we look at
    26812655*/
    2682 void
    2683 AbcSimplexPrimal::primalColumn(CoinPartitionedVector * updates,
    2684                                CoinPartitionedVector * spareRow2,
    2685                                CoinPartitionedVector * spareColumn1)
     2656void AbcSimplexPrimal::primalColumn(CoinPartitionedVector *updates,
     2657  CoinPartitionedVector *spareRow2,
     2658  CoinPartitionedVector *spareColumn1)
    26862659{
    2687   for (int i=0;i<4;i++)
    2688     multipleSequenceIn_[i]=-1;
     2660  for (int i = 0; i < 4; i++)
     2661    multipleSequenceIn_[i] = -1;
    26892662  sequenceIn_ = abcPrimalColumnPivot_->pivotColumn(updates,
    2690                                                    spareRow2, spareColumn1);
     2663    spareRow2, spareColumn1);
    26912664  if (sequenceIn_ >= 0) {
    26922665    valueIn_ = abcSolution_[sequenceIn_];
     
    27062679   After rowArray will have list of cost changes
    27072680*/
    2708 int
    2709 AbcSimplexPrimal::updatePrimalsInPrimal(CoinIndexedVector * rowArray,
    2710                                         double theta,
    2711                                         double & objectiveChange,
    2712                                         int valuesPass)
     2681int AbcSimplexPrimal::updatePrimalsInPrimal(CoinIndexedVector *rowArray,
     2682  double theta,
     2683  double &objectiveChange,
     2684  int valuesPass)
    27132685{
    27142686  // Cost on pivot row may change - may need to change dualIn
     
    27162688  if (pivotRow_ >= 0)
    27172689    oldCost = abcCost_[sequenceOut_];
    2718   double * work = rowArray->denseVector();
     2690  double *work = rowArray->denseVector();
    27192691  int number = rowArray->getNumElements();
    2720   int * which = rowArray->getIndices();
    2721  
     2692  int *which = rowArray->getIndices();
     2693
    27222694  int newNumber = 0;
    27232695  abcNonLinearCost_->setChangeInCost(0.0);
     
    27282700  if (!valuesPass) {
    27292701    for (iIndex = 0; iIndex < number; iIndex++) {
    2730      
     2702
    27312703      int iRow = which[iIndex];
    27322704      double alpha = work[iRow];
     
    27412713      // check if not active then okay
    27422714      if (!active(iRow) && (specialOptions_ & 4) == 0 && pivotRow_ != -1) {
    2743         // But make sure one going out is feasible
    2744         if (change > 0.0) {
    2745           // going down
    2746           if (value <= abcLower_[iPivot] + primalTolerance_) {
    2747             if (iPivot == sequenceOut_ && value > abcLower_[iPivot] - relaxedTolerance)
    2748               value = abcLower_[iPivot];
    2749             //double difference = abcNonLinearCost_->setOneBasic(iRow, value);
    2750             //assert (!difference || fabs(change) > 1.0e9);
    2751           }
    2752         } else {
    2753           // going up
    2754           if (value >= abcUpper_[iPivot] - primalTolerance_) {
    2755             if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
    2756               value = abcUpper_[iPivot];
    2757             //double difference = abcNonLinearCost_->setOneBasic(iRow, value);
    2758             //assert (!difference || fabs(change) > 1.0e9);
    2759           }
    2760         }
     2715        // But make sure one going out is feasible
     2716        if (change > 0.0) {
     2717          // going down
     2718          if (value <= abcLower_[iPivot] + primalTolerance_) {
     2719            if (iPivot == sequenceOut_ && value > abcLower_[iPivot] - relaxedTolerance)
     2720              value = abcLower_[iPivot];
     2721            //double difference = abcNonLinearCost_->setOneBasic(iRow, value);
     2722            //assert (!difference || fabs(change) > 1.0e9);
     2723          }
     2724        } else {
     2725          // going up
     2726          if (value >= abcUpper_[iPivot] - primalTolerance_) {
     2727            if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
     2728              value = abcUpper_[iPivot];
     2729            //double difference = abcNonLinearCost_->setOneBasic(iRow, value);
     2730            //assert (!difference || fabs(change) > 1.0e9);
     2731          }
     2732        }
    27612733      }
    27622734#endif
    27632735      if (active(iRow) || theta_ < 0.0) {
    2764         clearActive(iRow);
    2765         // But make sure one going out is feasible
    2766         if (change > 0.0) {
    2767           // going down
    2768           if (value <= abcLower_[iPivot] + primalTolerance_) {
    2769             if (iPivot == sequenceOut_ && value >= abcLower_[iPivot] - relaxedTolerance)
    2770               value = abcLower_[iPivot];
    2771             double difference = abcNonLinearCost_->setOneBasic(iRow, value);
    2772             if (difference) {
    2773               work[iRow] = difference;
    2774               //change reduced cost on this
    2775               //abcDj_[iPivot] = -difference;
    2776               which[newNumber++] = iRow;
    2777             }
    2778           }
    2779         } else {
    2780           // going up
    2781           if (value >= abcUpper_[iPivot] - primalTolerance_) {
    2782             if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
    2783               value = abcUpper_[iPivot];
    2784             double difference = abcNonLinearCost_->setOneBasic(iRow, value);
    2785             if (difference) {
    2786               work[iRow] = difference;
    2787               //change reduced cost on this
    2788               //abcDj_[iPivot] = -difference;
    2789               which[newNumber++] = iRow;
    2790             }
    2791           }
    2792         }
     2736        clearActive(iRow);
     2737        // But make sure one going out is feasible
     2738        if (change > 0.0) {
     2739          // going down
     2740          if (value <= abcLower_[iPivot] + primalTolerance_) {
     2741            if (iPivot == sequenceOut_ && value >= abcLower_[iPivot] - relaxedTolerance)
     2742              value = abcLower_[iPivot];
     2743            double difference = abcNonLinearCost_->setOneBasic(iRow, value);
     2744            if (difference) {
     2745              work[iRow] = difference;
     2746              //change reduced cost on this
     2747              //abcDj_[iPivot] = -difference;
     2748              which[newNumber++] = iRow;
     2749            }
     2750          }
     2751        } else {
     2752          // going up
     2753          if (value >= abcUpper_[iPivot] - primalTolerance_) {
     2754            if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
     2755              value = abcUpper_[iPivot];
     2756            double difference = abcNonLinearCost_->setOneBasic(iRow, value);
     2757            if (difference) {
     2758              work[iRow] = difference;
     2759              //change reduced cost on this
     2760              //abcDj_[iPivot] = -difference;
     2761              which[newNumber++] = iRow;
     2762            }
     2763          }
     2764        }
    27932765      }
    27942766    }
     
    27962768    // values pass so look at all
    27972769    for (iIndex = 0; iIndex < number; iIndex++) {
    2798      
     2770
    27992771      int iRow = which[iIndex];
    28002772      double alpha = work[iRow];
     
    28092781      // But make sure one going out is feasible
    28102782      if (change > 0.0) {
    2811         // going down
    2812         if (value <= abcLower_[iPivot] + primalTolerance_) {
    2813           if (iPivot == sequenceOut_ && value > abcLower_[iPivot] - relaxedTolerance)
    2814             value = abcLower_[iPivot];
    2815           double difference = abcNonLinearCost_->setOneBasic(iRow, value);
    2816           if (difference) {
    2817             work[iRow] = difference;
    2818             //change reduced cost on this
    2819             abcDj_[iPivot] = -difference;
    2820             which[newNumber++] = iRow;
    2821           }
    2822         }
     2783        // going down
     2784        if (value <= abcLower_[iPivot] + primalTolerance_) {
     2785          if (iPivot == sequenceOut_ && value > abcLower_[iPivot] - relaxedTolerance)
     2786            value = abcLower_[iPivot];
     2787          double difference = abcNonLinearCost_->setOneBasic(iRow, value);
     2788          if (difference) {
     2789            work[iRow] = difference;
     2790            //change reduced cost on this
     2791            abcDj_[iPivot] = -difference;
     2792            which[newNumber++] = iRow;
     2793          }
     2794        }
    28232795      } else {
    2824         // going up
    2825         if (value >= abcUpper_[iPivot] - primalTolerance_) {
    2826           if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
    2827             value = abcUpper_[iPivot];
    2828           double difference = abcNonLinearCost_->setOneBasic(iRow, value);
    2829           if (difference) {
    2830             work[iRow] = difference;
    2831             //change reduced cost on this
    2832             abcDj_[iPivot] = -difference;
    2833             which[newNumber++] = iRow;
    2834           }
    2835         }
     2796        // going up
     2797        if (value >= abcUpper_[iPivot] - primalTolerance_) {
     2798          if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
     2799            value = abcUpper_[iPivot];
     2800          double difference = abcNonLinearCost_->setOneBasic(iRow, value);
     2801          if (difference) {
     2802            work[iRow] = difference;
     2803            //change reduced cost on this
     2804            abcDj_[iPivot] = -difference;
     2805            which[newNumber++] = iRow;
     2806          }
     2807        }
    28362808      }
    28372809    }
     
    28522824}
    28532825// Perturbs problem
    2854 void
    2855 AbcSimplexPrimal::perturb(int /*type*/)
     2826void AbcSimplexPrimal::perturb(int /*type*/)
    28562827{
    28572828  if (perturbation_ > 100)
     
    28682839  double largestPositive;
    28692840  matrix_->rangeOfElements(smallestNegative, largestNegative,
    2870                            smallestPositive, largestPositive);
     2841    smallestPositive, largestPositive);
    28712842  smallestPositive = CoinMin(fabs(smallestNegative), smallestPositive);
    28722843  largestPositive = CoinMax(fabs(largestNegative), largestPositive);
    28732844  double elementRatio = largestPositive / smallestPositive;
    2874   if ((!numberIterations_ ||initialSumInfeasibilities_==1.23456789e-5)&& perturbation_ == 50) {
     2845  if ((!numberIterations_ || initialSumInfeasibilities_ == 1.23456789e-5) && perturbation_ == 50) {
    28752846    // See if we need to perturb
    28762847    int numberTotal = CoinMax(numberRows_, numberColumns_);
    2877     double * sort = new double[numberTotal];
     2848    double *sort = new double[numberTotal];
    28782849    int nFixed = 0;
    28792850    for (i = 0; i < numberRows_; i++) {
     
    28822853      double value = 0.0;
    28832854      if (lo && lo < 1.0e20) {
    2884         if (up && up < 1.0e20) {
    2885           value = 0.5 * (lo + up);
    2886           if (lo == up)
    2887             nFixed++;
    2888         } else {
    2889           value = lo;
    2890         }
     2855        if (up && up < 1.0e20) {
     2856          value = 0.5 * (lo + up);
     2857          if (lo == up)
     2858            nFixed++;
     2859        } else {
     2860          value = lo;
     2861        }
    28912862      } else {
    2892         if (up && up < 1.0e20)
    2893           value = up;
     2863        if (up && up < 1.0e20)
     2864          value = up;
    28942865      }
    28952866      sort[i] = value;
     
    29002871    for (i = 1; i < numberRows_; i++) {
    29012872      if (last != sort[i])
    2902         number++;
     2873        number++;
    29032874      last = sort[i];
    29042875    }
     
    29092880    if (number * 4 > numberRows_ || elementRatio > 1.0e12) {
    29102881      perturbation_ = 100;
    2911       delete [] sort;
     2882      delete[] sort;
    29122883      // Make sure feasible bounds
    29132884      if (abcNonLinearCost_) {
    2914         abcNonLinearCost_->refresh();
    2915         abcNonLinearCost_->checkInfeasibilities();
    2916         //abcNonLinearCost_->feasibleBounds();
     2885        abcNonLinearCost_->refresh();
     2886        abcNonLinearCost_->checkInfeasibilities();
     2887        //abcNonLinearCost_->feasibleBounds();
    29172888      }
    29182889      moveToBasic();
     
    29252896      nFixed = 0;
    29262897      for (i = 0; i < numberColumns_; i++) {
    2927         double lo = fabs(columnLower_[i]);
    2928         double up = fabs(columnUpper_[i]);
    2929         double value = 0.0;
    2930         if (lo && lo < 1.0e20) {
    2931           if (up && up < 1.0e20) {
    2932             value = 0.5 * (lo + up);
    2933             if (lo == up)
    2934               nFixed++;
    2935           } else {
    2936             value = lo;
    2937           }
    2938         } else {
    2939           if (up && up < 1.0e20)
    2940             value = up;
    2941         }
    2942         sort[i] = value;
     2898        double lo = fabs(columnLower_[i]);
     2899        double up = fabs(columnUpper_[i]);
     2900        double value = 0.0;
     2901        if (lo && lo < 1.0e20) {
     2902          if (up && up < 1.0e20) {
     2903            value = 0.5 * (lo + up);
     2904            if (lo == up)
     2905              nFixed++;
     2906          } else {
     2907            value = lo;
     2908          }
     2909        } else {
     2910          if (up && up < 1.0e20)
     2911            value = up;
     2912        }
     2913        sort[i] = value;
    29432914      }
    29442915      std::sort(sort, sort + numberColumns_);
     
    29462917      last = sort[0];
    29472918      for (i = 1; i < numberColumns_; i++) {
    2948         if (last != sort[i])
    2949           number++;
    2950         last = sort[i];
     2919        if (last != sort[i])
     2920          number++;
     2921        last = sort[i];
    29512922      }
    29522923      //printf("cratio number diff bounds %g (%d %d fixed)\n",((double)number)/((double) numberColumns_),
     
    29542925    }
    29552926#endif
    2956     delete [] sort;
     2927    delete[] sort;
    29572928    if (number * 4 > numberColumns_) {
    29582929      perturbation_ = 100;
    29592930      // Make sure feasible bounds
    29602931      if (abcNonLinearCost_) {
    2961         abcNonLinearCost_->refresh();
    2962         abcNonLinearCost_->checkInfeasibilities();
    2963         //abcNonLinearCost_->feasibleBounds();
     2932        abcNonLinearCost_->refresh();
     2933        abcNonLinearCost_->checkInfeasibilities();
     2934        //abcNonLinearCost_->feasibleBounds();
    29642935      }
    29652936      moveToBasic();
     
    29732944  // maximum fraction of rhs/bounds to perturb
    29742945  double maximumFraction = 1.0e-5;
    2975   double overallMultiplier= (perturbation_==50||perturbation_>54) ? 2.0 : 0.2;
     2946  double overallMultiplier = (perturbation_ == 50 || perturbation_ > 54) ? 2.0 : 0.2;
    29762947#ifdef HEAVY_PERTURBATION
    2977     if (perturbation_==50)
    2978       perturbation_=HEAVY_PERTURBATION;
     2948  if (perturbation_ == 50)
     2949    perturbation_ = HEAVY_PERTURBATION;
    29792950#endif
    29802951  if (perturbation_ >= 50) {
     
    29822953    for (i = 0; i < numberColumns_ + numberRows_; i++) {
    29832954      if (abcUpper_[i] > abcLower_[i] + primalTolerance_) {
    2984         double lowerValue, upperValue;
    2985         if (abcLower_[i] > -1.0e20)
    2986           lowerValue = fabs(abcLower_[i]);
    2987         else
    2988           lowerValue = 0.0;
    2989         if (abcUpper_[i] < 1.0e20)
    2990           upperValue = fabs(abcUpper_[i]);
    2991         else
    2992           upperValue = 0.0;
    2993         double value = CoinMax(fabs(lowerValue), fabs(upperValue));
    2994         value = CoinMin(value, abcUpper_[i] - abcLower_[i]);
     2955        double lowerValue, upperValue;
     2956        if (abcLower_[i] > -1.0e20)
     2957          lowerValue = fabs(abcLower_[i]);
     2958        else
     2959          lowerValue = 0.0;
     2960        if (abcUpper_[i] < 1.0e20)
     2961          upperValue = fabs(abcUpper_[i]);
     2962        else
     2963          upperValue = 0.0;
     2964        double value = CoinMax(fabs(lowerValue), fabs(upperValue));
     2965        value = CoinMin(value, abcUpper_[i] - abcLower_[i]);
    29952966#if 1
    2996         if (value) {
    2997           perturbation += value;
    2998           numberNonZero++;
    2999         }
     2967        if (value) {
     2968          perturbation += value;
     2969          numberNonZero++;
     2970        }
    30002971#else
    3001         perturbation = CoinMax(perturbation, value);
     2972        perturbation = CoinMax(perturbation, value);
    30022973#endif
    30032974      }
    30042975    }
    30052976    if (numberNonZero)
    3006       perturbation /= static_cast<double> (numberNonZero);
     2977      perturbation /= static_cast< double >(numberNonZero);
    30072978    else
    30082979      perturbation = 1.0e-1;
     
    30102981      // reduce
    30112982      while (perturbation_ < 55) {
    3012         perturbation_++;
    3013         perturbation *= 0.25;
    3014         bias *= 0.25;
     2983        perturbation_++;
     2984        perturbation *= 0.25;
     2985        bias *= 0.25;
    30152986      }
    30162987      perturbation_ = 50;
     
    30182989      // increase
    30192990      while (perturbation_ > 55) {
    3020         overallMultiplier *= 1.2;
    3021         perturbation_--;
    3022         perturbation *= 4.0;
     2991        overallMultiplier *= 1.2;
     2992        perturbation_--;
     2993        perturbation *= 4.0;
    30232994      }
    30242995      perturbation_ = 50;
     
    30453016    maximumFraction *= 0.1;
    30463017  }
    3047   if (savePerturbation==51) {
    3048     perturbation = CoinMin(0.1,perturbation);
    3049     maximumFraction *=0.1;
     3018  if (savePerturbation == 51) {
     3019    perturbation = CoinMin(0.1, perturbation);
     3020    maximumFraction *= 0.1;
    30503021  }
    30513022  //if (number != numberRows_)
     
    30583029  //const double * COIN_RESTRICT perturbationArray = perturbationSaved_;
    30593030  // Make sure feasible bounds
    3060   if (abcNonLinearCost_&&true) {
     3031  if (abcNonLinearCost_ && true) {
    30613032    abcNonLinearCost_->refresh();
    30623033    abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
     
    30643035  }
    30653036  double tolerance = 100.0 * primalTolerance_;
    3066   int numberChanged=0;
     3037  int numberChanged = 0;
    30673038  // Set bit if fixed
    3068   for (int i=0;i<numberRows_;i++) {
    3069     if (rowLower_[i]!=rowUpper_[i])
     3039  for (int i = 0; i < numberRows_; i++) {
     3040    if (rowLower_[i] != rowUpper_[i])
    30703041      internalStatus_[i] &= ~128;
    30713042    else
    30723043      internalStatus_[i] |= 128;
    30733044  }
    3074   for (int i=0;i<numberColumns_;i++) {
    3075     if (columnLower_[i]!=columnUpper_[i])
    3076       internalStatus_[i+numberRows_] &= ~128;
     3045  for (int i = 0; i < numberColumns_; i++) {
     3046    if (columnLower_[i] != columnUpper_[i])
     3047      internalStatus_[i + numberRows_] &= ~128;
    30773048    else
    3078       internalStatus_[i+numberRows_] |= 128;
     3049      internalStatus_[i + numberRows_] |= 128;
    30793050  }
    30803051  //double multiplier = perturbation*maximumFraction;
    30813052  for (iSequence = 0; iSequence < numberRows_ + numberColumns_; iSequence++) {
    30823053    if (getInternalStatus(iSequence) == basic) {
    3083       if ((internalStatus_[i] &128)!=0)
    3084         continue;
     3054      if ((internalStatus_[i] & 128) != 0)
     3055        continue;
    30853056      double lowerValue = abcLower_[iSequence];
    30863057      double upperValue = abcUpper_[iSequence];
    30873058      if (upperValue > lowerValue + tolerance) {
    3088         double solutionValue = abcSolution_[iSequence];
    3089         double difference = upperValue - lowerValue;
    3090         difference = CoinMin(difference, perturbation);
    3091         difference = CoinMin(difference, fabs(solutionValue) + 1.0);
    3092         double value = maximumFraction * (difference + bias);
    3093         value = CoinMin(value, 0.1);
    3094         value = CoinMax(value,primalTolerance_);
    3095         double perturbationValue=overallMultiplier*randomNumberGenerator_.randomDouble();
    3096         value *= perturbationValue;
    3097         if (solutionValue - lowerValue <= primalTolerance_) {
    3098           abcLower_[iSequence] -= value;
    3099           // change correct saved value
    3100           if (abcNonLinearCost_->getCurrentStatus(iSequence)==CLP_FEASIBLE)
    3101             abcPerturbation_[iSequence]=abcLower_[iSequence];
    3102           else
    3103             abcPerturbation_[iSequence+numberTotal_]=abcLower_[iSequence];
    3104         } else if (upperValue - solutionValue <= primalTolerance_) {
    3105           abcUpper_[iSequence] += value;
    3106           // change correct saved value
    3107           if (abcNonLinearCost_->getCurrentStatus(iSequence)==CLP_FEASIBLE)
    3108             abcPerturbation_[iSequence+numberTotal_]=abcUpper_[iSequence];
    3109           else
    3110             abcPerturbation_[iSequence]=abcUpper_[iSequence];
    3111         } else {
     3059        double solutionValue = abcSolution_[iSequence];
     3060        double difference = upperValue - lowerValue;
     3061        difference = CoinMin(difference, perturbation);
     3062        difference = CoinMin(difference, fabs(solutionValue) + 1.0);
     3063        double value = maximumFraction * (difference + bias);
     3064        value = CoinMin(value, 0.1);
     3065        value = CoinMax(value, primalTolerance_);
     3066        double perturbationValue = overallMultiplier * randomNumberGenerator_.randomDouble();
     3067        value *= perturbationValue;
     3068        if (solutionValue - lowerValue <= primalTolerance_) {
     3069          abcLower_[iSequence] -= value;
     3070          // change correct saved value
     3071          if (abcNonLinearCost_->getCurrentStatus(iSequence) == CLP_FEASIBLE)
     3072            abcPerturbation_[iSequence] = abcLower_[iSequence];
     3073          else
     3074            abcPerturbation_[iSequence + numberTotal_] = abcLower_[iSequence];
     3075        } else if (upperValue - solutionValue <= primalTolerance_) {
     3076          abcUpper_[iSequence] += value;
     3077          // change correct saved value
     3078          if (abcNonLinearCost_->getCurrentStatus(iSequence) == CLP_FEASIBLE)
     3079            abcPerturbation_[iSequence + numberTotal_] = abcUpper_[iSequence];
     3080          else
     3081            abcPerturbation_[iSequence] = abcUpper_[iSequence];
     3082        } else {
    31123083#if 0
    31133084          if (iSequence >= numberColumns_) {
     
    31213092          }
    31223093#else
    3123           value = 0.0;
    3124 #endif
    3125         }
    3126         if (value) {
    3127           numberChanged++;
    3128           if (printOut)
    3129             printf("col %d lower from %g to %g, upper from %g to %g\n",
    3130                    iSequence, abcLower_[iSequence], lowerValue, abcUpper_[iSequence], upperValue);
    3131           if (solutionValue) {
    3132             largest = CoinMax(largest, value);
    3133             if (value > (fabs(solutionValue) + 1.0)*largestPerCent)
    3134               largestPerCent = value / (fabs(solutionValue) + 1.0);
    3135           } else {
    3136             largestZero = CoinMax(largestZero, value);
    3137           }
    3138         }
     3094          value = 0.0;
     3095#endif
     3096        }
     3097        if (value) {
     3098          numberChanged++;
     3099          if (printOut)
     3100            printf("col %d lower from %g to %g, upper from %g to %g\n",
     3101              iSequence, abcLower_[iSequence], lowerValue, abcUpper_[iSequence], upperValue);
     3102          if (solutionValue) {
     3103            largest = CoinMax(largest, value);
     3104            if (value > (fabs(solutionValue) + 1.0) * largestPerCent)
     3105              largestPerCent = value / (fabs(solutionValue) + 1.0);
     3106          } else {
     3107            largestZero = CoinMax(largestZero, value);
     3108          }
     3109        }
    31393110      }
    31403111    }
     
    31443115    for (iSequence = 0; iSequence < maximumAbcNumberRows_ + numberColumns_; iSequence++) {
    31453116      if (getInternalStatus(iSequence) != basic) {
    3146         if ((internalStatus_[i] &128)!=0)
    3147           continue;
    3148         double lowerValue = abcLower_[iSequence];
    3149         double upperValue = abcUpper_[iSequence];
    3150         if (upperValue > lowerValue + tolerance) {
    3151           double solutionValue = abcSolution_[iSequence];
    3152           double difference = upperValue - lowerValue;
    3153           difference = CoinMin(difference, perturbation);
    3154           difference = CoinMin(difference, fabs(solutionValue) + 1.0);
    3155           double value = maximumFraction * (difference + bias);
    3156           value = CoinMin(value, 0.1);
    3157           value = CoinMax(value,primalTolerance_);
    3158           double perturbationValue=overallMultiplier*randomNumberGenerator_.randomDouble();
    3159           value *= perturbationValue;
    3160           if (solutionValue - lowerValue <= primalTolerance_) {
    3161             abcLower_[iSequence] -= value;
    3162             // change correct saved value
    3163             if (abcNonLinearCost_->getCurrentStatus(iSequence)==CLP_FEASIBLE)
    3164               abcPerturbation_[iSequence]=abcLower_[iSequence];
    3165             else
    3166               abcPerturbation_[iSequence+numberTotal_]=abcLower_[iSequence];
    3167           } else if (upperValue - solutionValue <= primalTolerance_) {
    3168             abcUpper_[iSequence] += value;
    3169             // change correct saved value
    3170             if (abcNonLinearCost_->getCurrentStatus(iSequence)==CLP_FEASIBLE)
    3171               abcPerturbation_[iSequence+numberTotal_]=abcUpper_[iSequence];
    3172             else
    3173               abcPerturbation_[iSequence]=abcUpper_[iSequence];
    3174           } else {
    3175             value = 0.0;
    3176           }
    3177           if (value) {
    3178             if (printOut)
    3179               printf("col %d lower from %g to %g, upper from %g to %g\n",
    3180                      iSequence, abcLower_[iSequence], lowerValue, abcUpper_[iSequence], upperValue);
    3181             if (solutionValue) {
    3182               largest = CoinMax(largest, value);
    3183               if (value > (fabs(solutionValue) + 1.0)*largestPerCent)
    3184                 largestPerCent = value / (fabs(solutionValue) + 1.0);
    3185             } else {
    3186               largestZero = CoinMax(largestZero, value);
    3187             }
    3188           }
    3189         }
     3117        if ((internalStatus_[i] & 128) != 0)
     3118          continue;
     3119        double lowerValue = abcLower_[iSequence];
     3120        double upperValue = abcUpper_[iSequence];
     3121        if (upperValue > lowerValue + tolerance) {
     3122          double solutionValue = abcSolution_[iSequence];
     3123          double difference = upperValue - lowerValue;
     3124          difference = CoinMin(difference, perturbation);
     3125          difference = CoinMin(difference, fabs(solutionValue) + 1.0);
     3126          double value = maximumFraction * (difference + bias);
     3127          value = CoinMin(value, 0.1);
     3128          value = CoinMax(value, primalTolerance_);
     3129          double perturbationValue = overallMultiplier * randomNumberGenerator_.randomDouble();
     3130          value *= perturbationValue;
     3131          if (solutionValue - lowerValue <= primalTolerance_) {
     3132            abcLower_[iSequence] -= value;
     3133            // change correct saved value
     3134            if (abcNonLinearCost_->getCurrentStatus(iSequence) == CLP_FEASIBLE)
     3135              abcPerturbation_[iSequence] = abcLower_[iSequence];
     3136            else
     3137              abcPerturbation_[iSequence + numberTotal_] = abcLower_[iSequence];
     3138          } else if (upperValue - solutionValue <= primalTolerance_) {
     3139            abcUpper_[iSequence] += value;
     3140            // change correct saved value
     3141            if (abcNonLinearCost_->getCurrentStatus(iSequence) == CLP_FEASIBLE)
     3142              abcPerturbation_[iSequence + numberTotal_] = abcUpper_[iSequence];
     3143            else
     3144              abcPerturbation_[iSequence] = abcUpper_[iSequence];
     3145          } else {
     3146            value = 0.0;
     3147          }
     3148          if (value) {
     3149            if (printOut)
     3150              printf("col %d lower from %g to %g, upper from %g to %g\n",
     3151                iSequence, abcLower_[iSequence], lowerValue, abcUpper_[iSequence], upperValue);
     3152            if (solutionValue) {
     3153              largest = CoinMax(largest, value);
     3154              if (value > (fabs(solutionValue) + 1.0) * largestPerCent)
     3155                largestPerCent = value / (fabs(solutionValue) + 1.0);
     3156            } else {
     3157              largestZero = CoinMax(largestZero, value);
     3158            }
     3159          }
     3160        }
    31903161      }
    31913162    }
     
    31943165  for (i = 0; i < numberColumns_ + numberRows_; i++) {
    31953166    internalStatus_[i] &= ~128;
    3196     switch(getInternalStatus(i)) {
    3197      
     3167    switch (getInternalStatus(i)) {
     3168
    31983169    case basic:
    31993170      break;
     
    32323203}
    32333204// un perturb
    3234 bool
    3235 AbcSimplexPrimal::unPerturb()
     3205bool AbcSimplexPrimal::unPerturb()
    32363206{
    32373207  if (perturbation_ != 101)
     
    32403210  copyFromSaved();
    32413211  // copy bounds to perturbation
    3242   CoinAbcMemcpy(abcPerturbation_,abcLower_,numberTotal_);
    3243   CoinAbcMemcpy(abcPerturbation_+numberTotal_,abcUpper_,numberTotal_);
     3212  CoinAbcMemcpy(abcPerturbation_, abcLower_, numberTotal_);
     3213  CoinAbcMemcpy(abcPerturbation_ + numberTotal_, abcUpper_, numberTotal_);
    32443214  //sanityCheck();
    32453215  // unflag
     
    32553225  abcNonLinearCost_->checkInfeasibilities(primalTolerance_);
    32563226  // Try using dual
    3257   return abcNonLinearCost_->sumInfeasibilities()!=0.0;
    3258  
     3227  return abcNonLinearCost_->sumInfeasibilities() != 0.0;
    32593228}
    32603229// Unflag all variables and return number unflagged
    3261 int
    3262 AbcSimplexPrimal::unflag()
     3230int AbcSimplexPrimal::unflag()
    32633231{
    32643232  int i;
     
    32733241      // only say if reasonable dj
    32743242      if (fabs(abcDj_[i]) > relaxedToleranceD)
    3275         numberFlagged++;
    3276     }
    3277   }
    3278 #if ABC_NORMAL_DEBUG>0
     3243        numberFlagged++;
     3244    }
     3245  }
     3246#if ABC_NORMAL_DEBUG > 0
    32793247  if (handler_->logLevel() > 2 && numberFlagged && objective_->type() > 1)
    32803248    printf("%d unflagged\n", numberFlagged);
     
    32833251}
    32843252// Do not change infeasibility cost and always say optimal
    3285 void
    3286 AbcSimplexPrimal::alwaysOptimal(bool onOff)
     3253void AbcSimplexPrimal::alwaysOptimal(bool onOff)
    32873254{
    32883255  if (onOff)
     
    32913258    specialOptions_ &= ~1;
    32923259}
    3293 bool
    3294 AbcSimplexPrimal::alwaysOptimal() const
     3260bool AbcSimplexPrimal::alwaysOptimal() const
    32953261{
    32963262  return (specialOptions_ & 1) != 0;
    32973263}
    32983264// Flatten outgoing variables i.e. - always to exact bound
    3299 void
    3300 AbcSimplexPrimal::exactOutgoing(bool onOff)
     3265void AbcSimplexPrimal::exactOutgoing(bool onOff)
    33013266{
    33023267  if (onOff)
     
    33053270    specialOptions_ &= ~4;
    33063271}
    3307 bool
    3308 AbcSimplexPrimal::exactOutgoing() const
     3272bool AbcSimplexPrimal::exactOutgoing() const
    33093273{
    33103274  return (specialOptions_ & 4) != 0;
     
    33203284  +3 max iterations (iteration done)
    33213285*/
    3322 int
    3323 AbcSimplexPrimal::pivotResult(int ifValuesPass)
     3286int AbcSimplexPrimal::pivotResult(int ifValuesPass)
    33243287{
    3325  
     3288
    33263289  bool roundAgain = true;
    33273290  int returnCode = -1;
    3328  
     3291
    33293292  // loop round if user setting and doing refactorization
    33303293  while (roundAgain) {
     
    33343297    sequenceOut_ = -1;
    33353298    usefulArray_[arrayForFtran_].clear();
    3336    
     3299
    33373300    // we found a pivot column
    33383301    // update the incoming column
     
    33453308    if ((moreSpecialOptions_ & 512) == 0) {
    33463309#endif
    3347       primalRow(&usefulArray_[arrayForFtran_], &usefulArray_[arrayForFlipBounds_], 
    3348                 &usefulArray_[arrayForTableauRow_],
    3349                 ifValuesPass);
     3310      primalRow(&usefulArray_[arrayForFtran_], &usefulArray_[arrayForFlipBounds_],
     3311        &usefulArray_[arrayForTableauRow_],
     3312        ifValuesPass);
    33503313#ifdef CLP_USER_DRIVEN
    33513314      // user can tell which use it is
    33523315      int status = eventHandler_->event(ClpEventHandler::pivotRow);
    33533316      if (status >= 0) {
    3354         problemStatus_ = 5;
    3355         secondaryStatus_ = ClpEventHandler::pivotRow;
    3356         break;
     3317        problemStatus_ = 5;
     3318        secondaryStatus_ = ClpEventHandler::pivotRow;
     3319        break;
    33573320      }
    33583321    } else {
    33593322      int status = eventHandler_->event(ClpEventHandler::pivotRow);
    33603323      if (status >= 0) {
    3361         problemStatus_ = 5;
    3362         secondaryStatus_ = ClpEventHandler::pivotRow;
    3363         break;
     3324        problemStatus_ = 5;
     3325        secondaryStatus_ = ClpEventHandler::pivotRow;
     3326        break;
    33643327      }
    33653328    }
     
    33693332      //assert (fabs(alpha_)>=1.0e-5||(objective_->type()<2||!objective_->activated())||pivotRow_==-2);
    33703333      if (pivotRow_ == -1 || (pivotRow_ >= 0 && fabs(alpha_) < 1.0e-5)) {
    3371         if(fabs(dualIn_) < 1.0e2 * dualTolerance_ && objective_->type() < 2) {
    3372           // try other way
    3373           directionIn_ = -directionIn_;
    3374           primalRow(&usefulArray_[arrayForFtran_], &usefulArray_[arrayForFlipBounds_],
    3375                     &usefulArray_[arrayForTableauRow_],
    3376                     0);
    3377         }
    3378         if (pivotRow_ == -1 || (pivotRow_ >= 0 && fabs(alpha_) < 1.0e-5)) {
    3379           // reject it
    3380           char x = isColumn(sequenceIn_) ? 'C' : 'R';
    3381           handler_->message(CLP_SIMPLEX_FLAG, messages_)
    3382             << x << sequenceWithin(sequenceIn_)
    3383             << CoinMessageEol;
    3384           setFlagged(sequenceIn_);
    3385           abcProgress_.incrementTimesFlagged();
    3386           abcProgress_.clearBadTimes();
    3387           lastBadIteration_ = numberIterations_; // say be more cautious
    3388           clearAll();
    3389           pivotRow_ = -1;
    3390           returnCode = -5;
    3391           break;
    3392         }
     3334        if (fabs(dualIn_) < 1.0e2 * dualTolerance_ && objective_->type() < 2) {
     3335          // try other way
     3336          directionIn_ = -directionIn_;
     3337          primalRow(&usefulArray_[arrayForFtran_], &usefulArray_[arrayForFlipBounds_],
     3338            &usefulArray_[arrayForTableauRow_],
     3339            0);
     3340        }
     3341        if (pivotRow_ == -1 || (pivotRow_ >= 0 && fabs(alpha_) < 1.0e-5)) {
     3342          // reject it
     3343          char x = isColumn(sequenceIn_) ? 'C' : 'R';
     3344          handler_->message(CLP_SIMPLEX_FLAG, messages_)
     3345            << x << sequenceWithin(sequenceIn_)
     3346            << CoinMessageEol;
     3347          setFlagged(sequenceIn_);
     3348          abcProgress_.incrementTimesFlagged();
     3349          abcProgress_.clearBadTimes();
     3350          lastBadIteration_ = numberIterations_; // say be more cautious
     3351          clearAll();
     3352          pivotRow_ = -1;
     3353          returnCode = -5;
     3354          break;
     3355        }
    33933356      }
    33943357    }
     
    34063369    }
    34073370#endif
    3408     if (!ifValuesPass && (saveDj * dualIn_ < test1 ||
    3409                                              fabs(saveDj - dualIn_) > checkValue*(1.0 + fabs(saveDj)) ||
    3410                           fabs(dualIn_) < test2)) {
    3411       if (!(saveDj * dualIn_ > 0.0 && CoinMin(fabs(saveDj), fabs(dualIn_)) >
    3412             1.0e5)) {
    3413         char x = isColumn(sequenceIn_) ? 'C' : 'R';
    3414         handler_->message(CLP_PRIMAL_DJ, messages_)
    3415           << x << sequenceWithin(sequenceIn_)
    3416           << saveDj << dualIn_
    3417           << CoinMessageEol;
    3418         if(lastGoodIteration_ != numberIterations_) {
    3419           clearAll();
    3420           pivotRow_ = -1; // say no weights update
    3421           returnCode = -4;
    3422           if(lastGoodIteration_ + 1 == numberIterations_) {
    3423             // not looking wonderful - try cleaning bounds
    3424             // put non-basics to bounds in case tolerance moved
    3425             abcNonLinearCost_->checkInfeasibilities(0.0);
    3426           }
    3427           sequenceOut_ = -1;
    3428           sequenceIn_=-1;
    3429           if (abcFactorization_->pivots()<10&&abcFactorization_->pivotTolerance()<0.25)
    3430             abcFactorization_->saferTolerances(1.0,-1.03);
    3431           break;
    3432         } else {
    3433           // take on more relaxed criterion
    3434           if (saveDj * dualIn_ < test1 ||
    3435                                  fabs(saveDj - dualIn_) > 2.0e-1 * (1.0 + fabs(dualIn_)) ||
    3436               fabs(dualIn_) < test2) {
    3437             // need to reject something
    3438             char x = isColumn(sequenceIn_) ? 'C' : 'R';
    3439             handler_->message(CLP_SIMPLEX_FLAG, messages_)
    3440               << x << sequenceWithin(sequenceIn_)
    3441               << CoinMessageEol;
    3442             setFlagged(sequenceIn_);
    3443             abcProgress_.incrementTimesFlagged();
    3444 #if 1 //def FEB_TRY
    3445             // Make safer?
    3446             double tolerance=abcFactorization_->pivotTolerance();
    3447             abcFactorization_->saferTolerances (1.0, -1.03);
    3448 #endif
    3449             abcProgress_.clearBadTimes();
    3450             lastBadIteration_ = numberIterations_; // say be more cautious
    3451             clearAll();
    3452             pivotRow_ = -1;
    3453             returnCode = -5;
    3454             if(tolerance<abcFactorization_->pivotTolerance())
    3455               returnCode=-4;
    3456             sequenceOut_ = -1;
    3457             sequenceIn_=-1;
    3458             break;
    3459           }
    3460         }
     3371    if (!ifValuesPass && (saveDj * dualIn_ < test1 || fabs(saveDj - dualIn_) > checkValue * (1.0 + fabs(saveDj)) || fabs(dualIn_) < test2)) {
     3372      if (!(saveDj * dualIn_ > 0.0 && CoinMin(fabs(saveDj), fabs(dualIn_)) > 1.0e5)) {
     3373        char x = isColumn(sequenceIn_) ? 'C' : 'R';
     3374        handler_->message(CLP_PRIMAL_DJ, messages_)
     3375          << x << sequenceWithin(sequenceIn_)
     3376          << saveDj << dualIn_
     3377          << CoinMessageEol;
     3378        if (lastGoodIteration_ != numberIterations_) {
     3379          clearAll();
     3380          pivotRow_ = -1; // say no weights update
     3381          returnCode = -4;
     3382          if (lastGoodIteration_ + 1 == numberIterations_) {
     3383            // not looking wonderful - try cleaning bounds
     3384            // put non-basics to bounds in case tolerance moved
     3385            abcNonLinearCost_->checkInfeasibilities(0.0);
     3386          }
     3387          sequenceOut_ = -1;
     3388          sequenceIn_ = -1;
     3389          if (abcFactorization_->pivots() < 10 && abcFactorization_->pivotTolerance() < 0.25)
     3390            abcFactorization_->saferTolerances(1.0, -1.03);
     3391          break;
     3392        } else {
     3393          // take on more relaxed criterion
     3394          if (saveDj * dualIn_ < test1 || fabs(saveDj - dualIn_) > 2.0e-1 * (1.0 + fabs(dualIn_)) || fabs(dualIn_) < test2) {
     3395            // need to reject something
     3396            char x = isColumn(sequenceIn_) ? 'C' : 'R';
     3397            handler_->message(CLP_SIMPLEX_FLAG, messages_)
     3398              << x << sequenceWithin(sequenceIn_)
     3399              << CoinMessageEol;
     3400            setFlagged(sequenceIn_);
     3401            abcProgress_.incrementTimesFlagged();
     3402#if 1 //def FEB_TRY \
     3403  // Make safer?
     3404            double tolerance = abcFactorization_->pivotTolerance();
     3405            abcFactorization_->saferTolerances(1.0, -1.03);
     3406#endif
     3407            abcProgress_.clearBadTimes();
     3408            lastBadIteration_ = numberIterations_; // say be more cautious
     3409            clearAll();
     3410            pivotRow_ = -1;
     3411            returnCode = -5;
     3412            if (tolerance < abcFactorization_->pivotTolerance())
     3413              returnCode = -4;
     3414            sequenceOut_ = -1;
     3415            sequenceIn_ = -1;
     3416            break;
     3417          }
     3418        }
    34613419      } else {
    3462         //printf("%d %g %g\n",numberIterations_,saveDj,dualIn_);
     3420        //printf("%d %g %g\n",numberIterations_,saveDj,dualIn_);
    34633421      }
    34643422    }
     
    34723430      //abcFactorization_->checkReplacePart1a(&usefulArray_[arrayForReplaceColumn_],pivotRow_);
    34733431      //usefulArray_[arrayForReplaceColumn_].print();
    3474       ftAlpha_=abcFactorization_->checkReplacePart1(&usefulArray_[arrayForReplaceColumn_],pivotRow_);
    3475       int updateStatus=abcFactorization_->checkReplacePart2(pivotRow_,btranAlpha_,alpha_,ftAlpha_);
     3432      ftAlpha_ = abcFactorization_->checkReplacePart1(&usefulArray_[arrayForReplaceColumn_], pivotRow_);
     3433      int updateStatus = abcFactorization_->checkReplacePart2(pivotRow_, btranAlpha_, alpha_, ftAlpha_);
    34763434      abcFactorization_->replaceColumnPart3(this,
    3477                                             &usefulArray_[arrayForReplaceColumn_],
    3478                                             &usefulArray_[arrayForFtran_],
    3479                                             pivotRow_,
    3480                                             ftAlpha_?ftAlpha_:alpha_);
    3481      
     3435        &usefulArray_[arrayForReplaceColumn_],
     3436        &usefulArray_[arrayForFtran_],
     3437        pivotRow_,
     3438        ftAlpha_ ? ftAlpha_ : alpha_);
     3439
    34823440      // if no pivots, bad update but reasonable alpha - take and invert
    3483       if (updateStatus == 2 &&
    3484           lastGoodIteration_ == numberIterations_ && fabs(alpha_) > 1.0e-5)
    3485         updateStatus = 4;
     3441      if (updateStatus == 2 && lastGoodIteration_ == numberIterations_ && fabs(alpha_) > 1.0e-5)
     3442        updateStatus = 4;
    34863443      if (updateStatus == 1 || updateStatus == 4) {
    3487         // slight error
    3488         if (abcFactorization_->pivots() > 5 || updateStatus == 4) {
    3489           returnCode = -3;
    3490         }
     3444        // slight error
     3445        if (abcFactorization_->pivots() > 5 || updateStatus == 4) {
     3446          returnCode = -3;
     3447        }
    34913448      } else if (updateStatus == 2) {
    3492         // major error
    3493         // better to have small tolerance even if slower
    3494         abcFactorization_->zeroTolerance(CoinMin(abcFactorization_->zeroTolerance(), 1.0e-15));
    3495         int maxFactor = abcFactorization_->maximumPivots();
    3496         if (maxFactor > 10) {
    3497           if (forceFactorization_ < 0)
    3498             forceFactorization_ = maxFactor;
    3499           forceFactorization_ = CoinMax(1, (forceFactorization_ >> 1));
    3500         }
    3501         // later we may need to unwind more e.g. fake bounds
    3502         if(lastGoodIteration_ != numberIterations_) {
    3503           clearAll();
    3504           pivotRow_ = -1;
    3505           sequenceIn_ = -1;
    3506           sequenceOut_ = -1;
    3507           returnCode = -4;
    3508           break;
    3509         } else {
    3510           // need to reject something
    3511           char x = isColumn(sequenceIn_) ? 'C' : 'R';
    3512           handler_->message(CLP_SIMPLEX_FLAG, messages_)
    3513             << x << sequenceWithin(sequenceIn_)
    3514             << CoinMessageEol;
    3515           setFlagged(sequenceIn_);
    3516           abcProgress_.incrementTimesFlagged();
    3517           abcProgress_.clearBadTimes();
    3518           lastBadIteration_ = numberIterations_; // say be more cautious
    3519           clearAll();
    3520           pivotRow_ = -1;
    3521           sequenceIn_ = -1;
    3522           sequenceOut_ = -1;
    3523           returnCode = -5;
    3524           break;
    3525          
    3526         }
     3449        // major error
     3450        // better to have small tolerance even if slower
     3451        abcFactorization_->zeroTolerance(CoinMin(abcFactorization_->zeroTolerance(), 1.0e-15));
     3452        int maxFactor = abcFactorization_->maximumPivots();
     3453        if (maxFactor > 10) {
     3454          if (forceFactorization_ < 0)
     3455            forceFactorization_ = maxFactor;
     3456          forceFactorization_ = CoinMax(1, (forceFactorization_ >> 1));
     3457        }
     3458        // later we may need to unwind more e.g. fake bounds
     3459        if (lastGoodIteration_ != numberIterations_) {
     3460          clearAll();
     3461          pivotRow_ = -1;
     3462          sequenceIn_ = -1;
     3463          sequenceOut_ = -1;
     3464          returnCode = -4;
     3465          break;
     3466        } else {
     3467          // need to reject something
     3468          char x = isColumn(sequenceIn_) ? 'C' : 'R';
     3469          handler_->message(CLP_SIMPLEX_FLAG, messages_)
     3470            << x << sequenceWithin(sequenceIn_)
     3471            << CoinMessageEol;
     3472          setFlagged(sequenceIn_);
     3473          abcProgress_.incrementTimesFlagged();
     3474          abcProgress_.clearBadTimes();
     3475          lastBadIteration_ = numberIterations_; // say be more cautious
     3476          clearAll();
     3477          pivotRow_ = -1;
     3478          sequenceIn_ = -1;
     3479          sequenceOut_ = -1;
     3480          returnCode = -5;
     3481          break;
     3482        }
    35273483      } else if (updateStatus == 3) {
    3528         // out of memory
    3529         // increase space if not many iterations
    3530         if (abcFactorization_->pivots() <
    3531             0.5 * abcFactorization_->maximumPivots() &&
    3532             abcFactorization_->pivots() < 200)
    3533           abcFactorization_->areaFactor(
    3534                                         abcFactorization_->areaFactor() * 1.1);
    3535         returnCode = -2; // factorize now
     3484        // out of memory
     3485        // increase space if not many iterations
     3486        if (abcFactorization_->pivots() < 0.5 * abcFactorization_->maximumPivots() && abcFactorization_->pivots() < 200)
     3487          abcFactorization_->areaFactor(
     3488            abcFactorization_->areaFactor() * 1.1);
     3489        returnCode = -2; // factorize now
    35363490      } else if (updateStatus == 5) {
    3537         problemStatus_ = -2; // factorize now
     3491        problemStatus_ = -2; // factorize now
    35383492      }
    35393493      // here do part of steepest - ready for next iteration
    35403494      if (!ifValuesPass)
    3541         abcPrimalColumnPivot_->updateWeights(&usefulArray_[arrayForFtran_]);
     3495        abcPrimalColumnPivot_->updateWeights(&usefulArray_[arrayForFtran_]);
    35423496    } else {
    35433497      if (pivotRow_ == -1) {
    3544         // no outgoing row is valid
    3545         if (valueOut_ != COIN_DBL_MAX) {
    3546           double objectiveChange = 0.0;
    3547           theta_ = valueOut_ - valueIn_;
    3548           updatePrimalsInPrimal(&usefulArray_[arrayForFtran_], theta_, objectiveChange, ifValuesPass);
    3549           abcSolution_[sequenceIn_] += theta_;
    3550         }
     3498        // no outgoing row is valid
     3499        if (valueOut_ != COIN_DBL_MAX) {
     3500          double objectiveChange = 0.0;
     3501          theta_ = valueOut_ - valueIn_;
     3502          updatePrimalsInPrimal(&usefulArray_[arrayForFtran_], theta_, objectiveChange, ifValuesPass);
     3503          abcSolution_[sequenceIn_] += theta_;
     3504        }
    35513505#ifdef CLP_USER_DRIVEN1
    3552         /* Note if valueOut_ < COIN_DBL_MAX and
     3506        /* Note if valueOut_ < COIN_DBL_MAX and
    35533507           theta_ reasonable then this may be a valid sub flip */
    3554         if(!userChoiceValid2(this)) {
    3555           if (abcFactorization_->pivots()<5) {
    3556             // flag variable
    3557             char x = isColumn(sequenceIn_) ? 'C' : 'R';
    3558             handler_->message(CLP_SIMPLEX_FLAG, messages_)
    3559               << x << sequenceWithin(sequenceIn_)
    3560               << CoinMessageEol;
    3561             setFlagged(sequenceIn_);
    3562             abcProgress_.incrementTimesFlagged();
    3563             abcProgress_.clearBadTimes();
    3564             roundAgain = true;
    3565             continue;
    3566           } else {
    3567             // try refactorizing first
    3568             returnCode = 4; //say looks odd but has iterated
    3569             break;
    3570           }
    3571         }
    3572 #endif
    3573         if (!abcFactorization_->pivots() && acceptablePivot_ <= 1.0e-8) {
    3574           returnCode = 2; //say looks unbounded
    3575           // do ray
    3576           primalRay(&usefulArray_[arrayForFtran_]);
    3577         } else {
    3578           acceptablePivot_ = 1.0e-8;
    3579           returnCode = 4; //say looks unbounded but has iterated
    3580         }
    3581         break;
     3508        if (!userChoiceValid2(this)) {
     3509          if (abcFactorization_->pivots() < 5) {
     3510            // flag variable
     3511            char x = isColumn(sequenceIn_) ? 'C' : 'R';
     3512            handler_->message(CLP_SIMPLEX_FLAG, messages_)
     3513              << x << sequenceWithin(sequenceIn_)
     3514              << CoinMessageEol;
     3515            setFlagged(sequenceIn_);
     3516            abcProgress_.incrementTimesFlagged();
     3517            abcProgress_.clearBadTimes();
     3518            roundAgain = true;
     3519            continue;
     3520          } else {
     3521            // try refactorizing first
     3522            returnCode = 4; //say looks odd but has iterated
     3523            break;
     3524          }
     3525        }
     3526#endif
     3527        if (!abcFactorization_->pivots() && acceptablePivot_ <= 1.0e-8) {
     3528          returnCode = 2; //say looks unbounded
     3529          // do ray
     3530          primalRay(&usefulArray_[arrayForFtran_]);
     3531        } else {
     3532          acceptablePivot_ = 1.0e-8;
     3533          returnCode = 4; //say looks unbounded but has iterated
     3534        }
     3535        break;
    35823536      } else {
    3583         // flipping from bound to bound
    3584       }
    3585     }
    3586    
     3537        // flipping from bound to bound
     3538      }
     3539    }
     3540
    35873541    double oldCost = 0.0;
    35883542    if (sequenceOut_ >= 0)
    35893543      oldCost = abcCost_[sequenceOut_];
    35903544    // update primal solution
    3591    
     3545
    35923546    double objectiveChange = 0.0;
    35933547    // after this usefulArray_[arrayForFtran_] is not empty - used to update djs
    35943548#ifdef CLP_USER_DRIVEN
    3595     if (theta_<0.0) {
    3596       if (theta_>=-1.0e-12)
    3597         theta_=0.0;
     3549    if (theta_ < 0.0) {
     3550      if (theta_ >= -1.0e-12)
     3551        theta_ = 0.0;
    35983552      //else
    35993553      //printf("negative theta %g\n",theta_);
     
    36013555#endif
    36023556    updatePrimalsInPrimal(&usefulArray_[arrayForFtran_], theta_, objectiveChange, ifValuesPass);
    3603    
     3557
    36043558    double oldValue = valueIn_;
    36053559    if (directionIn_ == -1) {
    36063560      // as if from upper bound
    36073561      if (sequenceIn_ != sequenceOut_) {
    3608         // variable becoming basic
    3609         valueIn_ -= fabs(theta_);
     3562        // variable becoming basic
     3563        valueIn_ -= fabs(theta_);
    36103564      } else {
    3611         valueIn_ = lowerIn_;
     3565        valueIn_ = lowerIn_;
    36123566      }
    36133567    } else {
    36143568      // as if from lower bound
    36153569      if (sequenceIn_ != sequenceOut_) {
    3616         // variable becoming basic
    3617         valueIn_ += fabs(theta_);
     3570        // variable becoming basic
     3571        valueIn_ += fabs(theta_);
    36183572      } else {
    3619         valueIn_ = upperIn_;
     3573        valueIn_ = upperIn_;
    36203574      }
    36213575    }
     
    36243578    if (sequenceIn_ != sequenceOut_) {
    36253579      if (directionOut_ > 0) {
    3626         valueOut_ = lowerOut_;
     3580        valueOut_ = lowerOut_;
    36273581      } else {
    3628         valueOut_ = upperOut_;
    3629       }
    3630       if(valueOut_ < abcLower_[sequenceOut_] - primalTolerance_)
    3631         valueOut_ = abcLower_[sequenceOut_] - 0.9 * primalTolerance_;
     3582        valueOut_ = upperOut_;
     3583      }
     3584      if (valueOut_ < abcLower_[sequenceOut_] - primalTolerance_)
     3585        valueOut_ = abcLower_[sequenceOut_] - 0.9 * primalTolerance_;
    36323586      else if (valueOut_ > abcUpper_[sequenceOut_] + primalTolerance_)
    3633         valueOut_ = abcUpper_[sequenceOut_] + 0.9 * primalTolerance_;
     3587        valueOut_ = abcUpper_[sequenceOut_] + 0.9 * primalTolerance_;
    36343588      // may not be exactly at bound and bounds may have changed
    36353589      // Make sure outgoing looks feasible
     
    36513605      // maximum iterations or equivalent
    36523606      returnCode = 3;
    3653     } else if(numberIterations_ == lastGoodIteration_
    3654               + 2 * abcFactorization_->maximumPivots()) {
     3607    } else if (numberIterations_ == lastGoodIteration_ + 2 * abcFactorization_->maximumPivots()) {
    36553608      // done a lot of flips - be safe
    36563609      returnCode = -2; // refactorize
     
    36603613      int status = eventHandler_->event(ClpEventHandler::endOfIteration);
    36613614      if (status >= 0) {
    3662         problemStatus_ = 5;
    3663         secondaryStatus_ = ClpEventHandler::endOfIteration;
    3664         returnCode = 3;
     3615        problemStatus_ = 5;
     3616        secondaryStatus_ = ClpEventHandler::endOfIteration;
     3617        returnCode = 3;
    36653618      }
    36663619    }
     
    36783631  +3 max iterations (iteration done)
    36793632*/
    3680 int
    3681 AbcSimplexPrimal::pivotResult4(int ifValuesPass)
     3633int AbcSimplexPrimal::pivotResult4(int ifValuesPass)
    36823634{
    3683  
     3635
    36843636  int returnCode = -1;
    3685   int numberMinor=0;
    3686   int numberDone=0;
    3687   CoinIndexedVector * vector[16];
     3637  int numberMinor = 0;
     3638  int numberDone = 0;
     3639  CoinIndexedVector *vector[16];
    36883640  pivotStruct stuff[4];
    3689   vector[0]=&usefulArray_[arrayForFtran_];
    3690   vector[1]=&usefulArray_[arrayForFlipBounds_];
    3691   vector[2]=&usefulArray_[arrayForTableauRow_];
    3692   vector[3]=&usefulArray_[arrayForBtran_];
     3641  vector[0] = &usefulArray_[arrayForFtran_];
     3642  vector[1] = &usefulArray_[arrayForFlipBounds_];
     3643  vector[2] = &usefulArray_[arrayForTableauRow_];
     3644  vector[3] = &usefulArray_[arrayForBtran_];
    36933645  /*
    36943646    For pricing we need to get a vector with difference in costs
     
    36993651  //double * saveCosts=CoinCopyOfArray(costBasic_,numberRows_);
    37003652  //double saveCostsIn[4];
    3701   for (int iMinor=0;iMinor<4;iMinor++) {
    3702     int sequenceIn=multipleSequenceIn_[iMinor];
    3703     if (sequenceIn<0)
     3653  for (int iMinor = 0; iMinor < 4; iMinor++) {
     3654    int sequenceIn = multipleSequenceIn_[iMinor];
     3655    if (sequenceIn < 0)
    37043656      break;
    3705     stuff[iMinor].valuesPass_=ifValuesPass;
    3706     stuff[iMinor].lowerIn_=abcLower_[sequenceIn];
    3707     stuff[iMinor].upperIn_=abcUpper_[sequenceIn];
    3708     stuff[iMinor].valueIn_=abcSolution_[sequenceIn];
    3709     stuff[iMinor].sequenceIn_=sequenceIn;
     3657    stuff[iMinor].valuesPass_ = ifValuesPass;
     3658    stuff[iMinor].lowerIn_ = abcLower_[sequenceIn];
     3659    stuff[iMinor].upperIn_ = abcUpper_[sequenceIn];
     3660    stuff[iMinor].valueIn_ = abcSolution_[sequenceIn];
     3661    stuff[iMinor].sequenceIn_ = sequenceIn;
    37103662    //saveCostsIn[iMinor]=abcCost_[sequenceIn];
    37113663    numberMinor++;
    37123664    if (iMinor) {
    3713       vector[4*iMinor]=rowArray_[2*iMinor-2];
    3714       vector[4*iMinor+1]=rowArray_[2*iMinor-1];
    3715       vector[4*iMinor+2]=columnArray_[2*iMinor-2];
    3716       vector[4*iMinor+3]=columnArray_[2*iMinor-1];
    3717     }
    3718     for (int i=0;i<4;i++)
    3719       vector[4*iMinor+i]->checkClear();
    3720     unpack(*vector[4*iMinor],sequenceIn);
    3721   }
    3722   int numberLeft=numberMinor;
     3665      vector[4 * iMinor] = rowArray_[2 * iMinor - 2];
     3666      vector[4 * iMinor + 1] = rowArray_[2 * iMinor - 1];
     3667      vector[4 * iMinor + 2] = columnArray_[2 * iMinor - 2];
     3668      vector[4 * iMinor + 3] = columnArray_[2 * iMinor - 1];
     3669    }
     3670    for (int i = 0; i < 4; i++)
     3671      vector[4 * iMinor + i]->checkClear();
     3672    unpack(*vector[4 * iMinor], sequenceIn);
     3673  }
     3674  int numberLeft = numberMinor;
    37233675  // parallel (with cpu)
    3724   for (int iMinor=1;iMinor<numberMinor;iMinor++) {
     3676  for (int iMinor = 1; iMinor < numberMinor; iMinor++) {
    37253677    // update the incoming columns
    3726     cilk_spawn abcFactorization_->updateColumnFT(*vector[4*iMinor],*vector[4*iMinor+3],iMinor);
    3727   }
    3728   abcFactorization_->updateColumnFT(*vector[0],*vector[+3],0);
     3678    cilk_spawn abcFactorization_->updateColumnFT(*vector[4 * iMinor], *vector[4 * iMinor + 3], iMinor);
     3679  }
     3680  abcFactorization_->updateColumnFT(*vector[0], *vector[+3], 0);
    37293681  cilk_sync;
    3730   for (int iMinor=0;iMinor<numberMinor;iMinor++) {
     3682  for (int iMinor = 0; iMinor < numberMinor; iMinor++) {
    37313683    // find best (or first if values pass)
    3732     int numberDo=1;
    3733     int jMinor=0;
    3734     while (jMinor<numberDo) {
    3735       int sequenceIn=stuff[jMinor].sequenceIn_;
    3736       double dj=abcDj_[sequenceIn];
    3737       bool bad=false;
     3684    int numberDo = 1;
     3685    int jMinor = 0;
     3686    while (jMinor < numberDo) {
     3687      int sequenceIn = stuff[jMinor].sequenceIn_;
     3688      double dj = abcDj_[sequenceIn];
     3689      bool bad = false;
    37383690      if (!bad) {
    3739         stuff[jMinor].dualIn_=dj;
    3740         stuff[jMinor].saveDualIn_=dj;
    3741         if (dj<0.0)
    3742           stuff[jMinor].directionIn_=1;
    3743         else
    3744           stuff[jMinor].directionIn_=-1;
    3745         jMinor++;
     3691        stuff[jMinor].dualIn_ = dj;
     3692        stuff[jMinor].saveDualIn_ = dj;
     3693        if (dj < 0.0)
     3694          stuff[jMinor].directionIn_ = 1;
     3695        else
     3696          stuff[jMinor].directionIn_ = -1;
     3697        jMinor++;
    37463698      } else {
    3747         numberDo--;
    3748         numberLeft--;
    3749         // throw away
    3750         for (int i=0;i<4;i++) {
    3751           vector[4*jMinor+i]->clear();
    3752           vector[4*jMinor+i]=vector[4*numberDo+i];
    3753           vector[4*numberDo+i]=NULL;
    3754         }
    3755         stuff[jMinor]=stuff[numberDo];
    3756       }
    3757     }
    3758     for (int jMinor=0;jMinor<numberDo;jMinor++) {
     3699        numberDo--;
     3700        numberLeft--;
     3701        // throw away
     3702        for (int i = 0; i < 4; i++) {
     3703          vector[4 * jMinor + i]->clear();
     3704          vector[4 * jMinor + i] = vector[4 * numberDo + i];
     3705          vector[4 * numberDo + i] = NULL;
     3706        }
     3707        stuff[jMinor] = stuff[numberDo];
     3708      }
     3709    }
     3710    for (int jMinor = 0; jMinor < numberDo; jMinor++) {
    37593711      // do ratio test and re-compute dj
    3760       primalRow(vector[4*jMinor],vector[4*jMinor+1],vector[4*jMinor+2],stuff[jMinor]);
     3712      primalRow(vector[4 * jMinor], vector[4 * jMinor + 1], vector[4 * jMinor + 2], stuff[jMinor]);
    37613713    }
    37623714    // choose best
    3763     int iBest=-1;
    3764     double bestMovement=-COIN_DBL_MAX;
    3765     for (int jMinor=0;jMinor<numberDo;jMinor++) {
    3766       double movement=stuff[jMinor].theta_*fabs(stuff[jMinor].dualIn_);
    3767       if (movement>bestMovement) {
    3768         bestMovement=movement;
    3769         iBest=jMinor;
     3715    int iBest = -1;
     3716    double bestMovement = -COIN_DBL_MAX;
     3717    for (int jMinor = 0; jMinor < numberDo; jMinor++) {
     3718      double movement = stuff[jMinor].theta_ * fabs(stuff[jMinor].dualIn_);
     3719      if (movement > bestMovement) {
     3720        bestMovement = movement;
     3721        iBest = jMinor;
    37703722      }
    37713723    }
     
    37743726      iBest=0;
    37753727#endif
    3776     if (iBest>=0) {
    3777       dualIn_=stuff[iBest].dualIn_;
    3778       dualOut_=stuff[iBest].dualOut_;
    3779       lowerOut_=stuff[iBest].lowerOut_;
    3780       upperOut_=stuff[iBest].upperOut_;
    3781       valueOut_=stuff[iBest].valueOut_;
    3782       sequenceOut_=stuff[iBest].sequenceOut_;
    3783       sequenceIn_=stuff[iBest].sequenceIn_;
    3784       lowerIn_=stuff[iBest].lowerIn_;
    3785       upperIn_=stuff[iBest].upperIn_;
    3786       valueIn_=stuff[iBest].valueIn_;
    3787       directionIn_=stuff[iBest].directionIn_;
    3788       directionOut_=stuff[iBest].directionOut_;
    3789       pivotRow_=stuff[iBest].pivotRow_;
    3790       theta_=stuff[iBest].theta_;
    3791       alpha_=stuff[iBest].alpha_;
     3728    if (iBest >= 0) {
     3729      dualIn_ = stuff[iBest].dualIn_;
     3730      dualOut_ = stuff[iBest].dualOut_;
     3731      lowerOut_ = stuff[iBest].lowerOut_;
     3732      upperOut_ = stuff[iBest].upperOut_;
     3733      valueOut_ = stuff[iBest].valueOut_;
     3734      sequenceOut_ = stuff[iBest].sequenceOut_;
     3735      sequenceIn_ = stuff[iBest].sequenceIn_;
     3736      lowerIn_ = stuff[iBest].lowerIn_;
     3737      upperIn_ = stuff[iBest].upperIn_;
     3738      valueIn_ = stuff[iBest].valueIn_;
     3739      directionIn_ = stuff[iBest].directionIn_;
     3740      directionOut_ = stuff[iBest].directionOut_;
     3741      pivotRow_ = stuff[iBest].pivotRow_;
     3742      theta_ = stuff[iBest].theta_;
     3743      alpha_ = stuff[iBest].alpha_;
    37923744#ifdef MULTIPLE_PRICE
    3793       for (int i=0;i<4*numberLeft;i++)
    3794         vector[i]->clear();
     3745      for (int i = 0; i < 4 * numberLeft; i++)
     3746        vector[i]->clear();
    37953747      return 0;
    37963748#endif
     
    37993751      double theta2 = numberIterations_ * abcNonLinearCost_->averageTheta();
    38003752      // Set average theta
    3801       abcNonLinearCost_->setAverageTheta((theta1 + theta2) / (static_cast<double> (numberIterations_ + 1)));
     3753      abcNonLinearCost_->setAverageTheta((theta1 + theta2) / (static_cast< double >(numberIterations_ + 1)));
    38023754      if (pivotRow_ == -1 || (pivotRow_ >= 0 && fabs(alpha_) < 1.0e-5)) {
    3803         if(fabs(dualIn_) < 1.0e2 * dualTolerance_) {
    3804           // try other way
    3805           stuff[iBest].directionIn_ = -directionIn_;
    3806           stuff[iBest].valuesPass_=0;
    3807           primalRow(vector[4*iBest],vector[4*iBest+1],vector[4*iBest+2],stuff[iBest]);
    3808           dualIn_=stuff[iBest].dualIn_;
    3809           dualOut_=stuff[iBest].dualOut_;
    3810           lowerOut_=stuff[iBest].lowerOut_;
    3811           upperOut_=stuff[iBest].upperOut_;
    3812           valueOut_=stuff[iBest].valueOut_;
    3813           sequenceOut_=stuff[iBest].sequenceOut_;
    3814           sequenceIn_=stuff[iBest].sequenceIn_;
    3815           lowerIn_=stuff[iBest].lowerIn_;
    3816           upperIn_=stuff[iBest].upperIn_;
    3817           valueIn_=stuff[iBest].valueIn_;
    3818           directionIn_=stuff[iBest].directionIn_;
    3819           directionOut_=stuff[iBest].directionOut_;
    3820           pivotRow_=stuff[iBest].pivotRow_;
    3821           theta_=stuff[iBest].theta_;
    3822           alpha_=stuff[iBest].alpha_;
    3823         }
    3824         if (pivotRow_ == -1 || (pivotRow_ >= 0 && fabs(alpha_) < 1.0e-5)) {
    3825           // reject it
    3826           char x = isColumn(sequenceIn_) ? 'C' : 'R';
    3827           handler_->message(CLP_SIMPLEX_FLAG, messages_)
    3828             << x << sequenceWithin(sequenceIn_)
    3829             << CoinMessageEol;
    3830           setFlagged(sequenceIn_);
    3831           abcProgress_.incrementTimesFlagged();
    3832           abcProgress_.clearBadTimes();
    3833           lastBadIteration_ = numberIterations_; // say be more cautious
    3834           clearAll();
    3835           pivotRow_ = -1;
    3836           returnCode = -5;
    3837           break;
    3838         }
    3839       }
    3840       CoinIndexedVector * bestUpdate=NULL;
     3755        if (fabs(dualIn_) < 1.0e2 * dualTolerance_) {
     3756          // try other way
     3757          stuff[iBest].directionIn_ = -directionIn_;
     3758          stuff[iBest].valuesPass_ = 0;
     3759          primalRow(vector[4 * iBest], vector[4 * iBest + 1], vector[4 * iBest + 2], stuff[iBest]);
     3760          dualIn_ = stuff[iBest].dualIn_;
     3761          dualOut_ = stuff[iBest].dualOut_;
     3762          lowerOut_ = stuff[iBest].lowerOut_;
     3763          upperOut_ = stuff[iBest].upperOut_;
     3764          valueOut_ = stuff[iBest].valueOut_;
     3765          sequenceOut_ = stuff[iBest].sequenceOut_;
     3766          sequenceIn_ = stuff[iBest].sequenceIn_;
     3767          lowerIn_ = stuff[iBest].lowerIn_;
     3768          upperIn_ = stuff[iBest].upperIn_;
     3769          valueIn_ = stuff[iBest].valueIn_;
     3770          directionIn_ = stuff[iBest].directionIn_;
     3771          directionOut_ = stuff[iBest].directionOut_;
     3772          pivotRow_ = stuff[iBest].pivotRow_;
     3773          theta_ = stuff[iBest].theta_;
     3774          alpha_ = stuff[iBest].alpha_;
     3775        }
     3776        if (pivotRow_ == -1 || (pivotRow_ >= 0 && fabs(alpha_) < 1.0e-5)) {
     3777          // reject it
     3778          char x = isColumn(sequenceIn_) ? 'C' : 'R';
     3779          handler_->message(CLP_SIMPLEX_FLAG, messages_)
     3780            << x << sequenceWithin(sequenceIn_)
     3781            << CoinMessageEol;
     3782          setFlagged(sequenceIn_);
     3783          abcProgress_.incrementTimesFlagged();
     3784          abcProgress_.clearBadTimes();
     3785          lastBadIteration_ = numberIterations_; // say be more cautious
     3786          clearAll();
     3787          pivotRow_ = -1;
     3788          returnCode = -5;
     3789          break;
     3790        }
     3791      }
     3792      CoinIndexedVector *bestUpdate = NULL;
    38413793      if (pivotRow_ >= 0) {
    3842         // Normal pivot
    3843         // if stable replace in basis
    3844         // check update
    3845         CoinIndexedVector * tempVector[4];
    3846         for (int i=0;i<4;i++)
    3847           tempVector[i] = vector[4*iBest+i];
    3848         returnCode = cilk_spawn doFTUpdate(tempVector);
    3849         bestUpdate=tempVector[0];
    3850         // after this bestUpdate is not empty - used to update djs ??
    3851         cilk_spawn updatePrimalsInPrimal(*bestUpdate, theta_, ifValuesPass);
    3852         numberDone++;
    3853         numberLeft--;
    3854         // throw away
    3855         for (int i=0;i<4;i++) {
    3856           vector[4*iBest+i]=vector[4*numberLeft+i];
    3857           vector[4*numberLeft+i]=NULL;
    3858         }
    3859         stuff[iBest]=stuff[numberLeft];
    3860         // update pi and other vectors and FT stuff
    3861         // parallel (can go 8 way?)
    3862         for (int jMinor=0;jMinor<numberLeft;jMinor++) {
    3863           cilk_spawn updatePartialUpdate(*vector[4*jMinor+3]);
    3864         }
    3865         // parallel
    3866         for (int jMinor=0;jMinor<numberLeft;jMinor++) {
    3867           stuff[jMinor].dualIn_=cilk_spawn updateMinorCandidate(*bestUpdate,*vector[4*jMinor],stuff[jMinor].sequenceIn_);
    3868         }
    3869         cilk_sync;
    3870         // throw away
    3871         for (int i=1;i<4;i++)
    3872           tempVector[i]->clear();
    3873         if (returnCode<-3) {
    3874           clearAll();
    3875           break;
    3876         }
    3877         // end Normal pivot
     3794        // Normal pivot
     3795        // if stable replace in basis
     3796        // check update
     3797        CoinIndexedVector *tempVector[4];
     3798        for (int i = 0; i < 4; i++)
     3799          tempVector[i] = vector[4 * iBest + i];
     3800        returnCode = cilk_spawn doFTUpdate(tempVector);
     3801        bestUpdate = tempVector[0];
     3802        // after this bestUpdate is not empty - used to update djs ??
     3803        cilk_spawn updatePrimalsInPrimal(*bestUpdate, theta_, ifValuesPass);
     3804        numberDone++;
     3805        numberLeft--;
     3806        // throw away
     3807        for (int i = 0; i < 4; i++) {
     3808          vector[4 * iBest + i] = vector[4 * numberLeft + i];
     3809          vector[4 * numberLeft + i] = NULL;
     3810        }
     3811        stuff[iBest] = stuff[numberLeft];
     3812        // update pi and other vectors and FT stuff
     3813        // parallel (can go 8 way?)
     3814        for (int jMinor = 0; jMinor < numberLeft; jMinor++) {
     3815          cilk_spawn updatePartialUpdate(*vector[4 * jMinor + 3]);
     3816        }
     3817        // parallel
     3818        for (int jMinor = 0; jMinor < numberLeft; jMinor++) {
     3819          stuff[jMinor].dualIn_ = cilk_spawn updateMinorCandidate(*bestUpdate, *vector[4 * jMinor], stuff[jMinor].sequenceIn_);
     3820        }
     3821        cilk_sync;
     3822        // throw away
     3823        for (int i = 1; i < 4; i++)
     3824          tempVector[i]->clear();
     3825        if (returnCode < -3) {
     3826          clearAll();
     3827          break;
     3828        }
     3829        // end Normal pivot
    38783830      } else {
    3879         // start Flip
    3880         vector[4*iBest+3]->clear();
    3881         if (pivotRow_ == -1) {
    3882           // no outgoing row is valid
    3883           if (valueOut_ != COIN_DBL_MAX) {
    3884             theta_ = valueOut_ - valueIn_;
    3885             updatePrimalsInPrimal(*vector[4*iBest], theta_, ifValuesPass);
    3886             abcSolution_[sequenceIn_] += theta_;
    3887           }
    3888           if (!abcFactorization_->pivots() && acceptablePivot_ <= 1.0e-8) {
    3889             returnCode = 2; //say looks unbounded
    3890             // do ray
    3891             primalRay(vector[4*iBest]);
    3892           } else {
    3893             acceptablePivot_ = 1.0e-8;
    3894             returnCode = 4; //say looks unbounded but has iterated
    3895           }
    3896           break;
    3897         } else {
    3898           // flipping from bound to bound
    3899           bestUpdate=vector[4*iBest];
    3900           // after this bestUpdate is not empty - used to update djs ??
    3901           updatePrimalsInPrimal(*bestUpdate, theta_, ifValuesPass);
    3902           // throw away best BUT remember we are updating costs somehow
    3903           numberDone++;
    3904           numberLeft--;
    3905           // throw away
    3906           for (int i=0;i<4;i++) {
    3907             if (i)
    3908               vector[4*iBest+i]->clear();
    3909             vector[4*iBest+i]=vector[4*numberLeft+i];
    3910             vector[4*numberLeft+i]=NULL;
    3911           }
    3912           stuff[iBest]=stuff[numberLeft];
    3913           // parallel
    3914           for (int jMinor=0;jMinor<numberLeft;jMinor++) {
    3915             stuff[jMinor].dualIn_=cilk_spawn updateMinorCandidate(*bestUpdate,*vector[4*jMinor],stuff[jMinor].sequenceIn_);
    3916           }
    3917           cilk_sync;
    3918         }
    3919         // end Flip
     3831        // start Flip
     3832        vector[4 * iBest + 3]->clear();
     3833        if (pivotRow_ == -1) {
     3834          // no outgoing row is valid
     3835          if (valueOut_ != COIN_DBL_MAX) {
     3836            theta_ = valueOut_ - valueIn_;
     3837            updatePrimalsInPrimal(*vector[4 * iBest], theta_, ifValuesPass);
     3838            abcSolution_[sequenceIn_] += theta_;
     3839          }
     3840          if (!abcFactorization_->pivots() && acceptablePivot_ <= 1.0e-8) {
     3841            returnCode = 2; //say looks unbounded
     3842            // do ray
     3843            primalRay(vector[4 * iBest]);
     3844          } else {
     3845            acceptablePivot_ = 1.0e-8;
     3846            returnCode = 4; //say looks unbounded but has iterated
     3847          }
     3848          break;
     3849        } else {
     3850          // flipping from bound to bound
     3851          bestUpdate = vector[4 * iBest];
     3852          // after this bestUpdate is not empty - used to update djs ??
     3853          updatePrimalsInPrimal(*bestUpdate, theta_, ifValuesPass);
     3854          // throw away best BUT remember we are updating costs somehow
     3855          numberDone++;
     3856          numberLeft--;
     3857          // throw away
     3858          for (int i = 0; i < 4; i++) {
     3859            if (i)
     3860              vector[4 * iBest + i]->clear();
     3861            vector[4 * iBest + i] = vector[4 * numberLeft + i];
     3862            vector[4 * numberLeft + i] = NULL;
     3863          }
     3864          stuff[iBest] = stuff[numberLeft];
     3865          // parallel
     3866          for (int jMinor = 0; jMinor < numberLeft; jMinor++) {
     3867            stuff[jMinor].dualIn_ = cilk_spawn updateMinorCandidate(*bestUpdate, *vector[4 * jMinor], stuff[jMinor].sequenceIn_);
     3868          }
     3869          cilk_sync;
     3870        }
     3871        // end Flip
    39203872      }
    39213873      bestUpdate->clear(); // as only works in values pass
    39223874
    39233875      if (directionIn_ == -1) {
    3924         // as if from upper bound
    3925         if (sequenceIn_ != sequenceOut_) {
    3926           // variable becoming basic
    3927           valueIn_ -= fabs(theta_);
    3928         } else {
    3929           valueIn_ = lowerIn_;
    3930         }
     3876        // as if from upper bound
     3877        if (sequenceIn_ != sequenceOut_) {
     3878          // variable becoming basic
     3879          valueIn_ -= fabs(theta_);
     3880        } else {
     3881          valueIn_ = lowerIn_;
     3882        }
    39313883      } else {
    3932         // as if from lower bound
    3933         if (sequenceIn_ != sequenceOut_) {
    3934           // variable becoming basic
    3935           valueIn_ += fabs(theta_);
    3936         } else {
    3937           valueIn_ = upperIn_;
    3938         }
     3884        // as if from lower bound
     3885        if (sequenceIn_ != sequenceOut_) {
     3886          // variable becoming basic
     3887          valueIn_ += fabs(theta_);
     3888        } else {
     3889          valueIn_ = upperIn_;
     3890        }
    39393891      }
    39403892      // outgoing
    39413893      if (sequenceIn_ != sequenceOut_) {
    3942         if (directionOut_ > 0) {
    3943           valueOut_ = lowerOut_;
    3944         } else {
    3945           valueOut_ = upperOut_;
    3946         }
    3947         if(valueOut_ < abcLower_[sequenceOut_] - primalTolerance_)
    3948           valueOut_ = abcLower_[sequenceOut_] - 0.9 * primalTolerance_;
    3949         else if (valueOut_ > abcUpper_[sequenceOut_] + primalTolerance_)
    3950           valueOut_ = abcUpper_[sequenceOut_] + 0.9 * primalTolerance_;
    3951         // may not be exactly at bound and bounds may have changed
    3952         // Make sure outgoing looks feasible
    3953         directionOut_ = abcNonLinearCost_->setOneOutgoing(pivotRow_, valueOut_);
    3954         abcSolution_[sequenceOut_] = valueOut_;
     3894        if (directionOut_ > 0) {
     3895          valueOut_ = lowerOut_;
     3896        } else {
     3897          valueOut_ = upperOut_;
     3898        }
     3899        if (valueOut_ < abcLower_[sequenceOut_] - primalTolerance_)
     3900          valueOut_ = abcLower_[sequenceOut_] - 0.9 * primalTolerance_;
     3901        else if (valueOut_ > abcUpper_[sequenceOut_] + primalTolerance_)
     3902          valueOut_ = abcUpper_[sequenceOut_] + 0.9 * primalTolerance_;
     3903        // may not be exactly at bound and bounds may have changed
     3904        // Make sure outgoing looks feasible
     3905        directionOut_ = abcNonLinearCost_->setOneOutgoing(pivotRow_, valueOut_);
     3906        abcSolution_[sequenceOut_] = valueOut_;
    39553907      }
    39563908      // change cost and bounds on incoming if primal
     
    39593911      swapPrimalStuff();
    39603912      if (whatNext == 1) {
    3961         returnCode = -2; // refactorize
     3913        returnCode = -2; // refactorize
    39623914      } else if (whatNext == 2) {
    3963         // maximum iterations or equivalent
    3964         returnCode = 3;
    3965       } else if(numberIterations_ == lastGoodIteration_
    3966                 + 2 * abcFactorization_->maximumPivots()) {
    3967         // done a lot of flips - be safe
    3968         returnCode = -2; // refactorize
     3915        // maximum iterations or equivalent
     3916        returnCode = 3;
     3917      } else if (numberIterations_ == lastGoodIteration_ + 2 * abcFactorization_->maximumPivots()) {
     3918        // done a lot of flips - be safe
     3919        returnCode = -2; // refactorize
    39693920      }
    39703921      // Check event
    39713922      {
    3972         int status = eventHandler_->event(ClpEventHandler::endOfIteration);
    3973         if (status >= 0) {
    3974           problemStatus_ = 5;
    3975           secondaryStatus_ = ClpEventHandler::endOfIteration;
    3976           returnCode = 3;
    3977         }
    3978       }
    3979     }
    3980     if (returnCode!=-1)
     3923        int status = eventHandler_->event(ClpEventHandler::endOfIteration);
     3924        if (status >= 0) {
     3925          problemStatus_ = 5;
     3926          secondaryStatus_ = ClpEventHandler::endOfIteration;
     3927          returnCode = 3;
     3928        }
     3929      }
     3930    }
     3931    if (returnCode != -1)
    39813932      break;
    39823933  }
    3983   for (int i=0;i<16;i++) {
     3934  for (int i = 0; i < 16; i++) {
    39843935    if (vector[i])
    39853936      vector[i]->clear();
     
    39913942}
    39923943// Do FT update as separate function for minor iterations (nonzero return code on problems)
    3993 int
    3994 AbcSimplexPrimal::doFTUpdate(CoinIndexedVector * vector[4])
     3944int AbcSimplexPrimal::doFTUpdate(CoinIndexedVector *vector[4])
    39953945{
    3996   ftAlpha_=abcFactorization_->checkReplacePart1(&usefulArray_[arrayForReplaceColumn_],
    3997                                                 vector[3],pivotRow_);
    3998   int updateStatus=abcFactorization_->checkReplacePart2(pivotRow_,btranAlpha_,alpha_,ftAlpha_);
     3946  ftAlpha_ = abcFactorization_->checkReplacePart1(&usefulArray_[arrayForReplaceColumn_],
     3947    vector[3], pivotRow_);
     3948  int updateStatus = abcFactorization_->checkReplacePart2(pivotRow_, btranAlpha_, alpha_, ftAlpha_);
    39993949  if (!updateStatus)
    40003950    abcFactorization_->replaceColumnPart3(this,
    4001                                           &usefulArray_[arrayForReplaceColumn_],
    4002                                           vector[0],
    4003                                           vector[3],
    4004                                           pivotRow_,
    4005                                           ftAlpha_?ftAlpha_:alpha_);
     3951      &usefulArray_[arrayForReplaceColumn_],
     3952      vector[0],
     3953      vector[3],
     3954      pivotRow_,
     3955      ftAlpha_ ? ftAlpha_ : alpha_);
    40063956  else
    40073957    vector[3]->clear();
    4008   int returnCode=0;
     3958  int returnCode = 0;
    40093959  // if no pivots, bad update but reasonable alpha - take and invert
    4010   if (updateStatus == 2 &&
    4011       lastGoodIteration_ == numberIterations_ && fabs(alpha_) > 1.0e-5)
     3960  if (updateStatus == 2 && lastGoodIteration_ == numberIterations_ && fabs(alpha_) > 1.0e-5)
    40123961    updateStatus = 4;
    40133962  if (updateStatus == 1 || updateStatus == 4) {
     
    40233972    if (maxFactor > 10) {
    40243973      if (forceFactorization_ < 0)
    4025         forceFactorization_ = maxFactor;
     3974        forceFactorization_ = maxFactor;
    40263975      forceFactorization_ = CoinMax(1, (forceFactorization_ >> 1));
    40273976    }
    40283977    // later we may need to unwind more e.g. fake bounds
    4029     if(lastGoodIteration_ != numberIterations_) {
     3978    if (lastGoodIteration_ != numberIterations_) {
    40303979      pivotRow_ = -1;
    40313980      sequenceIn_ = -1;
     
    40363985      char x = isColumn(sequenceIn_) ? 'C' : 'R';
    40373986      handler_->message(CLP_SIMPLEX_FLAG, messages_)
    4038         << x << sequenceWithin(sequenceIn_)
    4039         << CoinMessageEol;
     3987        << x << sequenceWithin(sequenceIn_)
     3988        << CoinMessageEol;
    40403989      setFlagged(sequenceIn_);
    40413990      abcProgress_.incrementTimesFlagged();
     
    40463995      sequenceOut_ = -1;
    40473996      returnCode = -5;
    4048      
    40493997    }
    40503998  } else if (updateStatus == 3) {
    40513999    // out of memory
    40524000    // increase space if not many iterations
    4053     if (abcFactorization_->pivots() <
    4054         0.5 * abcFactorization_->maximumPivots() &&
    4055         abcFactorization_->pivots() < 200)
     4001    if (abcFactorization_->pivots() < 0.5 * abcFactorization_->maximumPivots() && abcFactorization_->pivots() < 200)
    40564002      abcFactorization_->areaFactor(
    4057                                     abcFactorization_->areaFactor() * 1.1);
     4003        abcFactorization_->areaFactor() * 1.1);
    40584004    returnCode = -2; // factorize now
    40594005  } else if (updateStatus == 5) {
    40604006    problemStatus_ = -2; // factorize now
    4061     returnCode=-2;
     4007    returnCode = -2;
    40624008  }
    40634009  return returnCode;
     
    40664012   costs are changed
    40674013*/
    4068 void
    4069 AbcSimplexPrimal::updatePrimalsInPrimal(CoinIndexedVector & rowArray,
    4070                                         double theta,bool valuesPass)
     4014void AbcSimplexPrimal::updatePrimalsInPrimal(CoinIndexedVector &rowArray,
     4015  double theta, bool valuesPass)
    40714016{
    40724017  // Cost on pivot row may change - may need to change dualIn
     
    40744019  //if (pivotRow_ >= 0)
    40754020  //oldCost = abcCost_[sequenceOut_];
    4076   double * work = rowArray.denseVector();
     4021  double *work = rowArray.denseVector();
    40774022  int number = rowArray.getNumElements();
    4078   int * which = rowArray.getIndices();
     4023  int *which = rowArray.getIndices();
    40794024  // allow for case where bound+tolerance == bound
    40804025  //double tolerance = 0.999999*primalTolerance_;
     
    40834028  if (!valuesPass) {
    40844029    for (iIndex = 0; iIndex < number; iIndex++) {
    4085      
     4030
    40864031      int iRow = which[iIndex];
    40874032      double alpha = work[iRow];
     
    40944039      solutionBasic_[iRow] = value;
    40954040      if (active(iRow) || theta_ < 0.0) {
    4096         clearActive(iRow);
    4097         // But make sure one going out is feasible
    4098         if (change > 0.0) {
    4099           // going down
    4100           if (value <= abcLower_[iPivot] + primalTolerance_) {
    4101             if (iPivot == sequenceOut_ && value >= abcLower_[iPivot] - relaxedTolerance)
    4102               value = abcLower_[iPivot];
    4103             abcNonLinearCost_->setOneBasic(iRow, value);
    4104           }
    4105         } else {
    4106           // going up
    4107           if (value >= abcUpper_[iPivot] - primalTolerance_) {
    4108             if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
    4109               value = abcUpper_[iPivot];
    4110             abcNonLinearCost_->setOneBasic(iRow, value);
    4111           }
    4112         }
     4041        clearActive(iRow);
     4042        // But make sure one going out is feasible
     4043        if (change > 0.0) {
     4044          // going down
     4045          if (value <= abcLower_[iPivot] + primalTolerance_) {
     4046            if (iPivot == sequenceOut_ && value >= abcLower_[iPivot] - relaxedTolerance)
     4047              value = abcLower_[iPivot];
     4048            abcNonLinearCost_->setOneBasic(iRow, value);
     4049          }
     4050        } else {
     4051          // going up
     4052          if (value >= abcUpper_[iPivot] - primalTolerance_) {
     4053            if (iPivot == sequenceOut_ && value < abcUpper_[iPivot] + relaxedTolerance)
     4054              value = abcUpper_[iPivot];
     4055            abcNonLinearCost_->setOneBasic(iRow, value);
     4056          }
     4057        }
    41134058      }
    41144059    }
     
    41164061    // values pass so look at all
    41174062    for (iIndex = 0; iIndex < number; iIndex++) {
    4118      
     4063
    41194064      int iRow = which[iIndex];
    41204065      double alpha = work[iRow];