Changeset 1344


Ignore:
Timestamp:
Mar 16, 2009 6:25:57 AM (11 years ago)
Author:
forrest
Message:

changes to simplex and lots of stuff and start Mumps cholesky

Location:
trunk/Clp/src
Files:
2 added
25 edited

Legend:

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

    r1327 r1344  
    13441344  parameters[numberParameters-1].append("Taucs_dummy");
    13451345#endif
     1346#ifdef MUMPS_BARRIER
     1347  parameters[numberParameters-1].append("Mumps");
     1348#define REAL_BARRIER
     1349#else
     1350  parameters[numberParameters-1].append("Mumps_dummy");   
     1351#endif
    13461352  parameters[numberParameters-1].setLonghelp
    13471353    (
     
    13921398first.  This primitive strategy can be surprsingly effective.  The column order\
    13931399 option is obviously not on costs but easy to code here."
     1400     );
     1401  parameters[numberParameters++]=
     1402    CbcOrClpParam("cplex!Use","Whether to use Cplex!",
     1403                  "off",CPX);
     1404  parameters[numberParameters-1].append("on");
     1405  parameters[numberParameters-1].setLonghelp
     1406    (
     1407     " If the user has Cplex, but wants to use some of Cbc'sheuristics \
     1408then you can!  If this is on, then Cbc will get to the root node and then \
     1409hand over to Cplex.  If heuristics find a solution this can be significantly \
     1410quicker.  You will probably want to switch off Cbc's cuts as Cplex thinks \
     1411they are genuine constraints.  It is also probable that you want to switch \
     1412off preprocessing, although for difficult problems it is worth trying \
     1413both."
    13941414     );
    13951415#endif
  • trunk/Clp/src/CbcOrClpParam.hpp

    r1321 r1344  
    8181    TWOMIRCUTS,PREPROCESS,FPUMP,GREEDY,COMBINE,LOCALTREE,SOS,
    8282    LANDPCUTS,RINS,RESIDCUTS,RENS,DIVINGS,DIVINGC,DIVINGF,DIVINGG,DIVINGL,
    83     DIVINGP,DIVINGV,DINS,PIVOTANDFIX,RANDROUND,NAIVE,ZEROHALFCUTS,
     83    DIVINGP,DIVINGV,DINS,PIVOTANDFIX,RANDROUND,NAIVE,ZEROHALFCUTS,CPX,
    8484   
    8585    DIRECTORY=301,DIRSAMPLE,DIRNETLIB,DIRMIPLIB,IMPORT,EXPORT,RESTORE,SAVE,DUALSIMPLEX,PRIMALSIMPLEX,EITHERSIMPLEX,
  • trunk/Clp/src/ClpFactorization.cpp

    r1321 r1344  
    1010#include "CoinIndexedVector.hpp"
    1111#include "ClpSimplex.hpp"
     12#include "ClpSimplexDual.hpp"
    1213#include "ClpMatrixBase.hpp"
    1314#ifndef SLIM_CLP
     
    12931294                              int solveType, bool valuesPass)
    12941295{
     1296  //if ((model->specialOptions()&16384))
     1297  //printf("factor at %d iterations\n",model->numberIterations());
    12951298  ClpMatrixBase * matrix = model->clpMatrix();
    12961299  int numberRows = model->numberRows();
     
    13011304  factorization_instrument(-1);
    13021305#endif
     1306  bool anyChanged=false;
    13031307  if (coinFactorizationB_) {
    13041308    coinFactorizationB_->setStatus(-99);
     
    14531457      } else {
    14541458        // Change pivotTemp to be correct list
     1459        anyChanged=true;
    14551460        coinFactorizationB_->makeNonSingular(pivotTemp,numberColumns);
    14561461        double * columnLower = model->lowerRegion();
     
    15691574    factorization_instrument(2);
    15701575#endif
     1576    if ( anyChanged&&model->algorithm()<0&&solveType>0) {
     1577      double dummyCost;
     1578      static_cast<ClpSimplexDual *> (model)->changeBounds(3,
     1579                                                          NULL,dummyCost);
     1580    }
    15711581    return coinFactorizationB_->status();
    15721582  }
     
    26162626#endif
    26172627}
    2618 #endif
     2628// Set tolerances to safer of existing and given
     2629void
     2630ClpFactorization::saferTolerances (  double zeroValue,
     2631                                    double pivotValue)
     2632{
     2633  // better to have small tolerance even if slower
     2634  zeroTolerance(CoinMin(zeroTolerance(),zeroValue));
     2635  // better to have large tolerance even if slower
     2636  double newValue;
     2637  if (pivotValue>0.0)
     2638    newValue = pivotValue;
     2639  else
     2640    newValue = -pivotTolerance()*pivotValue;
     2641  pivotTolerance(CoinMin(CoinMax(pivotTolerance(),newValue),0.999));
     2642}
     2643#endif
  • trunk/Clp/src/ClpFactorization.hpp

    r1287 r1344  
    168168    if (coinFactorizationA_) coinFactorizationA_->zeroTolerance(value); else coinFactorizationB_->zeroTolerance(value);
    169169  }
     170  /// Set tolerances to safer of existing and given
     171  void saferTolerances (  double zeroTolerance, double pivotTolerance);
    170172  /**  get sparse threshold */
    171173  inline int sparseThreshold ( ) const
  • trunk/Clp/src/ClpMain.cpp

    r1333 r1344  
    5151#endif
    5252#ifdef TAUCS_BARRIER
     53#define FOREIGN_BARRIER
     54#endif
     55#ifdef MUMPS_BARRIER
    5356#define FOREIGN_BARRIER
    5457#endif
  • trunk/Clp/src/ClpModel.cpp

    r1328 r1344  
    177177    eventHandler_=NULL;
    178178  }
    179   whatsChanged_=0;
     179  whatsChanged_ = 0;
    180180  delete matrix_;
    181181  matrix_=NULL;
     
    291291  delete [] rowObjective_;
    292292  rowObjective_=ClpCopyOfArray(rowObjective,numberRows_);
    293   whatsChanged_=0;
     293  whatsChanged_ = 0;
    294294}
    295295void
     
    796796      rowUpper_ = ClpCopyOfArray ( rhs.rowUpper_, numberRows_ );
    797797      columnLower_ = ClpCopyOfArray ( rhs.columnLower_, numberColumns_ );
    798       int scaleLength = ((specialOptions_&131072)==0) ? 1 : 2;
    799798      columnUpper_ = ClpCopyOfArray ( rhs.columnUpper_, numberColumns_ );
    800       rowScale_ = ClpCopyOfArray(rhs.rowScale_,numberRows_*scaleLength);
    801       columnScale_ = ClpCopyOfArray(rhs.columnScale_,numberColumns_*scaleLength);
     799      rowScale_ = ClpCopyOfArray(rhs.rowScale_,numberRows_*2);
     800      columnScale_ = ClpCopyOfArray(rhs.columnScale_,numberColumns_*2);
    802801      if (rhs.objective_)
    803802        objective_  = rhs.objective_->clone();
     
    855854      ClpDisjointCopyN ( rhs.columnLower_, numberColumns_,columnLower_ );
    856855      assert ((specialOptions_&131072)==0);
     856      abort();
    857857      ClpDisjointCopyN ( rhs.columnUpper_, numberColumns_,columnUpper_ );
    858858      if (rhs.objective_) {
     
    27072707  bool savePrefix =m.messageHandler()->prefix();
    27082708  m.messageHandler()->setPrefix(handler_->prefix());
     2709  m.setSmallElementValue(CoinMax(smallElement_,m.getSmallElementValue()));
    27092710  double time1 = CoinCpuTime(),time2;
    27102711  int status=0;
     
    37773778  }
    37783779  for (i=0;i<numberColumns_;i++) {
    3779     double multiplier = 1.0/columnScale_[i];
     3780    double multiplier = 1.0*inverseColumnScale_[i];
    37803781    columnActivity_[i] *= multiplier;
    37813782    reducedCost_[i] *= columnScale_[i];
     
    38043805    // reverse scaling
    38053806    for (i=0;i<numberRows_;i++)
    3806       rowScale_[i] = 1.0/rowScale_[i];
     3807      rowScale_[i] = 1.0*inverseRowScale_[i];
    38073808    for (i=0;i<numberColumns_;i++)
    3808       columnScale_[i] = 1.0/columnScale_[i];
     3809      columnScale_[i] = 1.0*inverseColumnScale_[i];
    38093810    gutsOfScaling();
    38103811  }
     
    39713972  savedRowScale_ = NULL;
    39723973  savedColumnScale_ = NULL;
     3974}
     3975// Set new row matrix
     3976void
     3977ClpModel::setNewRowCopy(ClpMatrixBase * newCopy)
     3978{
     3979  delete rowCopy_;
     3980  rowCopy_ = newCopy;
    39733981}
    39743982/* Find a network subset.
  • trunk/Clp/src/ClpModel.hpp

    r1329 r1344  
    599599   /// Row Matrix
    600600   inline ClpMatrixBase * rowCopy() const       { return rowCopy_; }
     601   /// Set new row matrix
     602   void setNewRowCopy(ClpMatrixBase * newCopy);
    601603   /// Clp Matrix
    602604   inline ClpMatrixBase * clpMatrix() const     { return matrix_; }
     
    847849      32768 - called from Osi
    848850      65536 - keep arrays around as much as possible (also use maximumR/C)
    849       131072 - scale factor arrays have inverse values at end
     851      131072 - unused
    850852      262144 - extra copy of scaled matrix
    851853      524288 - Clp fast dual
     
    10031005        512 - basis not changed (up to user to set this to 0)
    10041006              top bits may be used internally
     1007        shift by 65336 is 3 all same, 1 all except col bounds     
    10051008  */
    10061009  unsigned int whatsChanged_;
  • trunk/Clp/src/ClpNode.cpp

    r1321 r1344  
    188188  int nFix=0;
    189189  double gap = CoinMax(model->dualObjectiveLimit()-objectiveValue_,1.0e-4);
    190 #define PSEUDO 3
     190#define PSEUDO 1
    191191#if PSEUDO==1||PSEUDO==2
    192192  // Column copy of matrix
    193193  ClpPackedMatrix * matrix = model->clpScaledMatrix();
    194   if (!matrix)
     194  const double *objective = model->costRegion() ;
     195  if (!objective) {
     196    objective=model->objective();
     197    //if (!matrix)
    195198    matrix = dynamic_cast< ClpPackedMatrix*>(model->clpMatrix());
     199  } else if (!matrix) {
     200    matrix = dynamic_cast< ClpPackedMatrix*>(model->clpMatrix());
     201  }
    196202  const double * element = matrix->getElements();
    197203  const int * row = matrix->getIndices();
    198204  const CoinBigIndex * columnStart = matrix->getVectorStarts();
    199205  const int * columnLength = matrix->getVectorLengths();
    200   const double *objective = model->costRegion() ;
    201206  double direction = model->optimizationDirection();
    202207  const double * dual = dualSolution_+numberColumns;
     
    239244  }
    240245#endif
    241 #elif PSEUDO==3
     246#endif
    242247  const double * downPseudo = stuff->downPseudo_;
    243248  const int * numberDown = stuff->numberDown_;
     
    246251  const int * numberUp = stuff->numberUp_;
    247252  const int * numberUpInfeasible = stuff->numberUpInfeasible_;
    248 #endif
    249253  int iInteger=0;
    250254  for (iColumn=0;iColumn<numberColumns;iColumn++) {
     
    260264        double upValue = 0.0;
    261265        double downValue = 0.0;
    262         double value2 = objective ? direction*objective[iColumn] : 0.0;
     266        double value2 = direction*objective[iColumn];
     267        //double dj2=value2;
    263268        if (value2) {
    264269          if (value2>0.0)
     
    274279          if (value2) {
    275280            value2 *= element[j];
     281            //dj2 += value2;
    276282#if PSEUDO==2
    277283            assert (activeWeight[iRow]>0.0||fabs(dual[iRow])<1.0e-6);
     
    284290          }
    285291        }
    286         upValue = CoinMax(upValue,1.0e-8);
    287         downValue = CoinMax(downValue,1.0e-8);
     292        //assert (fabs(dj2)<1.0e-4);
     293        double upValue2 = (upPseudo[iInteger]/(1.0+numberUp[iInteger]));
     294        // Extra for infeasible branches
     295        if (numberUp[iInteger]) {
     296          double ratio = 1.0+static_cast<double>(numberUpInfeasible[iInteger])/
     297            static_cast<double>(numberUp[iInteger]);
     298          upValue2 *= ratio;
     299        }
     300        double downValue2 = (downPseudo[iInteger]/(1.0+numberDown[iInteger]));
     301        if (numberDown[iInteger]) {
     302          double ratio = 1.0+static_cast<double>(numberDownInfeasible[iInteger])/
     303            static_cast<double>(numberDown[iInteger]);
     304          downValue2 *= ratio;
     305        }
     306        //printf("col %d - downPi %g up %g, downPs %g up %g\n",
     307        //     iColumn,upValue,downValue,upValue2,downValue2);
     308        upValue = CoinMax(0.1*upValue,upValue2);
     309        downValue = CoinMax(0.1*downValue,downValue2);
     310        //upValue = CoinMax(upValue,1.0e-8);
     311        //downValue = CoinMax(downValue,1.0e-8);
    288312        upValue *= ceil(value)-value;
    289313        downValue *= value-floor(value);
     
    294318          infeasibility = 0.1*CoinMax(upValue,downValue)+
    295319            0.9*CoinMin(upValue,downValue) + integerTolerance;
     320        estimatedSolution_ += CoinMin(upValue2,downValue2);
    296321#elif PSEUDO==3
    297322        // Extra 100% for infeasible branches
     
    630655  saveOptions_(0),
    631656  solverOptions_(0),
     657  maximumNodes_(0),
    632658  nDepth_(-1),
    633659  nNodes_(0),
     
    659685    saveOptions_(rhs.saveOptions_),
    660686    solverOptions_(rhs.solverOptions_),
     687    maximumNodes_(rhs.maximumNodes_),
    661688    nDepth_(rhs.nDepth_),
    662689    nNodes_(rhs.nNodes_),
     
    689716    saveOptions_ = rhs.saveOptions_;
    690717    solverOptions_ = rhs.solverOptions_;
     718    maximumNodes_ = rhs.maximumNodes_;
    691719    int n = maximumNodes();
    692720    if (n) {
    693721      for (int i=0;i<n;i++)
    694722        delete nodeInfo_[i];
    695       delete [] nodeInfo_;
    696     }
     723    }
     724    delete [] nodeInfo_;
     725    nodeInfo_=NULL;
    697726    nDepth_ = rhs.nDepth_;
    698727    nNodes_ = rhs.nNodes_;
     
    724753    saveOptions_ = 0;
    725754    solverOptions_ = 0;
     755    maximumNodes_ = 0;
    726756    nDepth_ = -1;
    727757    nNodes_ = 0;
     
    747777    for (int i=0;i<n;i++)
    748778      delete nodeInfo_[i];
    749     delete [] nodeInfo_;
    750   }
     779  }
     780  delete [] nodeInfo_;
    751781}
    752782// Return maximum number of nodes
     
    754784ClpNodeStuff::maximumNodes() const
    755785{
    756   int n;
    757   if ((solverOptions_&32)==0)
    758     n = (1<<nDepth_)+1+nDepth_;
    759   else if (nDepth_)
    760     n = 1+1+nDepth_;
    761   else
    762     n = 0;
     786  int n=0;
     787#if 0
     788  if (nDepth_!=-1) {
     789    if ((solverOptions_&32)==0)
     790      n = (1<<nDepth_);
     791    else if (nDepth_)
     792      n = 1;
     793  }
     794  assert (n==maximumNodes_-(1+nDepth_)||n==0);
     795#else
     796  if (nDepth_!=-1) {
     797    n = maximumNodes_-(1+nDepth_);
     798    assert (n>0);
     799  }
     800#endif
    763801  return n;
     802}
     803// Return maximum space for nodes
     804int
     805ClpNodeStuff::maximumSpace() const
     806{
     807  return maximumNodes_;
    764808}
    765809/* Fill with pseudocosts */
  • trunk/Clp/src/ClpNode.hpp

    r1286 r1344  
    192192  /// Return maximum number of nodes
    193193  int maximumNodes() const;
     194  /// Return maximum space for nodes
     195  int maximumSpace() const;
    194196  //@}
    195197 
     
    236238  */
    237239  int solverOptions_;
     240  /// Maximum number of nodes to do
     241  int maximumNodes_;
    238242  /// Number deep
    239243  int nDepth_;
  • trunk/Clp/src/ClpPackedMatrix.cpp

    r1328 r1344  
    3131// dummy
    3232#define coin_prefetch(mem)
     33#endif
     34#ifdef INTEL_MKL
     35#include "mkl_spblas.h"
    3336#endif
    3437
     
    469472                                       double * COIN_RESTRICT spare) const
    470473{
    471   assert (rowScale);
    472474  // get matrix data pointers
    473475  const int * COIN_RESTRICT row = matrix_->getIndices();
    474476  const CoinBigIndex * COIN_RESTRICT columnStart = matrix_->getVectorStarts();
    475477  const double * COIN_RESTRICT elementByColumn = matrix_->getElements();
    476   if (!spare) {
    477     for (int jColumn=0;jColumn<number;jColumn++) {
    478       int iColumn = which[jColumn];
    479       CoinBigIndex j;
    480       CoinBigIndex start=columnStart[iColumn];
    481       CoinBigIndex next=columnStart[iColumn+1];
    482       double value=0.0;
    483       for (j=start;j<next;j++) {
    484         int jRow=row[j];
    485         value += x[jRow]*elementByColumn[j]*rowScale[jRow];
    486       }
    487       y[iColumn] -= value*columnScale[iColumn];
     478  if (!spare||!rowScale) {
     479    if (rowScale) {
     480      for (int jColumn=0;jColumn<number;jColumn++) {
     481        int iColumn = which[jColumn];
     482        CoinBigIndex j;
     483        CoinBigIndex start=columnStart[iColumn];
     484        CoinBigIndex next=columnStart[iColumn+1];
     485        double value=0.0;
     486        for (j=start;j<next;j++) {
     487          int jRow=row[j];
     488          value += x[jRow]*elementByColumn[j]*rowScale[jRow];
     489        }
     490        y[iColumn] -= value*columnScale[iColumn];
     491      }
     492    } else {
     493      for (int jColumn=0;jColumn<number;jColumn++) {
     494        int iColumn = which[jColumn];
     495        CoinBigIndex j;
     496        CoinBigIndex start=columnStart[iColumn];
     497        CoinBigIndex next=columnStart[iColumn+1];
     498        double value=0.0;
     499        for (j=start;j<next;j++) {
     500          int jRow=row[j];
     501          value += x[jRow]*elementByColumn[j];
     502        }
     503        y[iColumn] -= value;
     504      }
    488505    }
    489506  } else {
     
    578595    const double * elementByColumn = matrix_->getElements();
    579596    const double * rowScale = model->rowScale();
     597#if 0
     598    ClpPackedMatrix * scaledMatrix = model->clpScaledMatrix();
     599    if (rowScale&&scaledMatrix) {
     600      rowScale=NULL;
     601      // get matrix data pointers
     602      row = scaledMatrix->getIndices();
     603      columnStart = scaledMatrix->getVectorStarts();
     604      columnLength = scaledMatrix->getVectorLengths();
     605      elementByColumn = scaledMatrix->getElements();
     606    }
     607#endif
    580608    if (packed) {
    581609      // need to expand pi into y
     
    612640        }
    613641      } else {
     642#ifdef CLP_INVESTIGATE
     643        if (model->clpScaledMatrix())
     644          printf("scaledMatrix_ at %d of ClpPackedMatrix\n",__LINE__);
     645#endif
    614646        // scaled
    615647        // modify pi so can collapse to one loop
     
    679711        }
    680712      } else {
     713#ifdef CLP_INVESTIGATE
     714        if (model->clpScaledMatrix())
     715          printf("scaledMatrix_ at %d of ClpPackedMatrix\n",__LINE__);
     716#endif
    681717        // scaled
    682718        if (scalar==-1.0) {
     
    773809  assert (!y->getNumElements());
    774810  assert (numberActiveColumns_>0);
     811  const ClpPackedMatrix * thisMatrix = this;
     812#if 0
     813  ClpPackedMatrix * scaledMatrix = model->clpScaledMatrix();
     814  if (rowScale&&scaledMatrix) {
     815    rowScale=NULL;
     816    // get matrix data pointers
     817    row = scaledMatrix->getIndices();
     818    columnStart = scaledMatrix->getVectorStarts();
     819    elementByColumn = scaledMatrix->getElements();
     820    thisMatrix=scaledMatrix;
     821    //printf("scaledMatrix\n");
     822  } else if (rowScale) {
     823    //printf("no scaledMatrix\n");
     824  } else {
     825    //printf("no rowScale\n");
     826  }
     827#endif
    775828  if (packed) {
    776829    // need to expand pi into y
     
    795848      }
    796849      if (!columnCopy_) {
    797         numberNonZero=gutsOfTransposeTimesUnscaled(pi,columnArray->getIndices(),
    798                                            columnArray->denseVector(),
    799                                            zeroTolerance);
     850        numberNonZero=thisMatrix->gutsOfTransposeTimesUnscaled(pi,
     851                                                               columnArray->getIndices(),
     852                                                               columnArray->denseVector(),
     853                                                               zeroTolerance);
    800854        columnArray->setNumElements(numberNonZero);
    801855        //xA++;
     
    806860      }
    807861    } else {
     862#ifdef CLP_INVESTIGATE
     863      if (model->clpScaledMatrix())
     864        printf("scaledMatrix_ at %d of ClpPackedMatrix\n",__LINE__);
     865#endif
    808866      // scaled
    809867      // modify pi so can collapse to one loop
     
    905963      }
    906964    } else {
     965#ifdef CLP_INVESTIGATE
     966      if (model->clpScaledMatrix())
     967        printf("scaledMatrix_ at %d of ClpPackedMatrix\n",__LINE__);
     968#endif
    907969      // scaled
    908970      if (scalar==-1.0) {
     
    11361198  const CoinBigIndex * COIN_RESTRICT columnStart = matrix_->getVectorStarts();
    11371199  const double * COIN_RESTRICT elementByColumn = matrix_->getElements();
     1200#if 1 //ndef INTEL_MKL
    11381201  double value = 0.0;
    11391202  CoinBigIndex j;
     
    11611224    index[numberNonZero++]=iColumn;
    11621225  }
     1226#else
     1227  char transA='N';
     1228  //int numberRows = matrix_->getNumRows();
     1229  mkl_cspblas_dcsrgemv(&transA,const_cast<int *>(&numberActiveColumns_),
     1230                       const_cast<double *>(elementByColumn),
     1231                       const_cast<int *>(columnStart),
     1232                       const_cast<int *>(row),
     1233                       const_cast<double *>(pi),array);
     1234  int iColumn;
     1235  for (iColumn=0;iColumn<numberActiveColumns_;iColumn++) {
     1236    double value = array[iColumn];
     1237    if (value) {
     1238      array[iColumn]=0.0;
     1239      if (fabs(value)>zeroTolerance) {
     1240        array[numberNonZero]=value;
     1241        index[numberNonZero++]=iColumn;
     1242      }
     1243    }
     1244  }
     1245#endif
    11631246  return numberNonZero;
    11641247}
     
    14011484  assert (!rowArray->packedMode());
    14021485  columnArray->setPacked();
    1403   if (!(flags_&2)&&numberToDo>5) {
    1404     // no gaps and a reasonable number
     1486  ClpPackedMatrix * scaledMatrix = model->clpScaledMatrix();
     1487  int flags=flags_;
     1488  if (rowScale&&scaledMatrix&&!(scaledMatrix->flags()&2)) {
     1489    flags=0;
     1490    rowScale=NULL;
     1491    // get matrix data pointers
     1492    row = scaledMatrix->getIndices();
     1493    columnStart = scaledMatrix->getVectorStarts();
     1494    elementByColumn = scaledMatrix->getElements();
     1495  }
     1496  if (!(flags&2)&&numberToDo) {
     1497    // no gaps
    14051498    if (!rowScale) {
    14061499      int iColumn = which[0];
     
    14251518      array[jColumn]=value;
    14261519    } else {
     1520#ifdef CLP_INVESTIGATE
     1521      if (model->clpScaledMatrix())
     1522        printf("scaledMatrix_ at %d of ClpPackedMatrix\n",__LINE__);
     1523#endif
    14271524      // scaled
    14281525      const double * columnScale = model->columnScale();
     
    14521549      array[jColumn]=value;
    14531550    }
    1454   } else {
     1551  } else if (numberToDo) {
    14551552    // gaps
    14561553    if (!rowScale) {
     
    14671564      }
    14681565    } else {
     1566#ifdef CLP_INVESTIGATE
     1567      if (model->clpScaledMatrix())
     1568        printf("scaledMatrix_ at %d of ClpPackedMatrix - flags %d (%d) n %d\n",
     1569               __LINE__,flags_,model->clpScaledMatrix()->flags(),numberToDo);
     1570#endif
    14691571      // scaled
    14701572      const double * columnScale = model->columnScale();
     
    15571659    const int * whichRow = pi1->getIndices();
    15581660    int i;
     1661    ClpPackedMatrix * scaledMatrix = model->clpScaledMatrix();
     1662    if (rowScale&&scaledMatrix) {
     1663      rowScale=NULL;
     1664      // get matrix data pointers
     1665      row = scaledMatrix->getIndices();
     1666      columnStart = scaledMatrix->getVectorStarts();
     1667      elementByColumn = scaledMatrix->getElements();
     1668    }
    15591669    if (!rowScale) {
    15601670      // modify pi so can collapse to one loop
     
    17401850      }
    17411851    } else {
     1852#ifdef CLP_INVESTIGATE
     1853      if (model->clpScaledMatrix())
     1854        printf("scaledMatrix_ at %d of ClpPackedMatrix\n",__LINE__);
     1855#endif
    17421856      // scaled
    17431857      // can also scale piWeight as not used again
     
    18561970    }
    18571971  } else {
     1972#ifdef CLP_INVESTIGATE
     1973    if (model->clpScaledMatrix())
     1974      printf("scaledMatrix_ at %d of ClpPackedMatrix\n",__LINE__);
     1975#endif
    18581976    // scaled
    18591977    const double * columnScale = model->columnScale();
     
    20162134  }
    20172135}
    2018 //static int scale_stats[5]={0,0,0,0,0};
    2019 // Creates scales for column copy (rowCopy in model may be modified)
     2136#if 0
    20202137int
    2021 ClpPackedMatrix::scale(ClpModel * model, const ClpSimplex * baseModel) const
    2022 {
     2138ClpPackedMatrix::scale2(ClpModel * model) const
     2139{
     2140  ClpSimplex * baseModel=NULL;
    20232141#ifndef NDEBUG
    20242142  //checkFlags();
     
    20342152  }
    20352153  ClpMatrixBase * rowCopyBase=model->rowCopy();
    2036   if (!rowCopyBase) {
    2037     // temporary copy
    2038     rowCopyBase = reverseOrderedCopy();
    2039   }
    2040 #ifndef NDEBUG
    2041   ClpPackedMatrix* rowCopy =
    2042     dynamic_cast< ClpPackedMatrix*>(rowCopyBase);
    2043   // Make sure it is really a ClpPackedMatrix
    2044   assert (rowCopy!=NULL);
    2045 #else
    2046   ClpPackedMatrix* rowCopy =
    2047     static_cast< ClpPackedMatrix*>(rowCopyBase);
    2048 #endif
    2049 
    2050   const int * column = rowCopy->getIndices();
    2051   const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
    2052   const double * element = rowCopy->getElements();
    2053   int scaleLength = ((model->specialOptions()&131072)==0) ? 1 : 2;
    20542154  double * rowScale;
    20552155  double * columnScale;
     
    20592159  double * inverseColumnScale = NULL;
    20602160  if (!model->rowScale()) {
    2061     rowScale = new double [numberRows*scaleLength];
    2062     columnScale = new double [numberColumns*scaleLength];
    2063     if (scaleLength==2) {
    2064       inverseRowScale = rowScale+numberRows;
    2065       inverseColumnScale = columnScale+numberColumns;
    2066     }
     2161    rowScale = new double [numberRows*2];
     2162    columnScale = new double [numberColumns*2];
     2163    inverseRowScale = rowScale+numberRows;
     2164    inverseColumnScale = columnScale+numberColumns;
    20672165    arraysExist=false;
    20682166  } else {
     
    20732171    arraysExist=true;
    20742172  }
     2173  assert (inverseRowScale==rowScale+numberRows);
     2174  assert (inverseColumnScale==columnScale+numberColumns);
    20752175  // we are going to mark bits we are interested in
    20762176  char * usefulRow = new char [numberRows];
     
    20842184  // mark free rows
    20852185  for (iRow=0;iRow<numberRows;iRow++) {
    2086 #ifndef LEAVE_FIXED
     2186#if 0 //ndef LEAVE_FIXED
    20872187    if (rowUpper[iRow]<1.0e20||
    20882188        rowLower[iRow]>-1.0e20)
     
    21012201  double smallest=1.0e50;
    21022202  // get matrix data pointers
    2103   const int * row = matrix_->getIndices();
     2203  int * row = matrix_->getMutableIndices();
    21042204  const CoinBigIndex * columnStart = matrix_->getVectorStarts();
    2105   const int * columnLength = matrix_->getVectorLengths();
    2106   const double * elementByColumn = matrix_->getElements();
     2205  int * columnLength = matrix_->getMutableVectorLengths();
     2206  double * elementByColumn = matrix_->getMutableElements();
     2207  bool deletedElements=false;
    21072208  for (iColumn=0;iColumn<numberColumns;iColumn++) {
    21082209    CoinBigIndex j;
    21092210    char useful=0;
     2211    bool deleteSome = false;
     2212    int start = columnStart[iColumn];
     2213    int end = start + columnLength[iColumn];
    21102214#ifndef LEAVE_FIXED
    21112215    if (columnUpper[iColumn]>
    21122216        columnLower[iColumn]+1.0e-12) {
    21132217#endif
    2114       for (j=columnStart[iColumn];
    2115            j<columnStart[iColumn]+columnLength[iColumn];j++) {
     2218      for (j=start;j<end;j++) {
    21162219        iRow=row[j];
    2117         if(elementByColumn[j]&&usefulRow[iRow]) {
    2118           useful=1;
    2119           largest = CoinMax(largest,fabs(elementByColumn[j]));
    2120           smallest = CoinMin(smallest,fabs(elementByColumn[j]));
     2220        double value = fabs(elementByColumn[j]);
     2221        if (value>1.0e-20) {
     2222          if(usefulRow[iRow]) {
     2223            useful=1;
     2224            largest = CoinMax(largest,fabs(elementByColumn[j]));
     2225            smallest = CoinMin(smallest,fabs(elementByColumn[j]));
     2226          }
     2227        } else {
     2228          // small
     2229          deleteSome=true;
    21212230        }
    21222231      }
    21232232#ifndef LEAVE_FIXED
     2233    } else {
     2234      // just check values
     2235      for (j=start;j<end;j++) {
     2236        double value = fabs(elementByColumn[j]);
     2237        if (value<=1.0e-20) {
     2238          // small
     2239          deleteSome=true;
     2240        }
     2241      }
    21242242    }
    21252243#endif
    21262244    usefulColumn[iColumn]=useful;
     2245    if (deleteSome) {
     2246      deletedElements=true;
     2247      CoinBigIndex put = start;
     2248      for (j=start;j<end;j++) {
     2249        double value = elementByColumn[j];
     2250        if (fabs(value)>1.0e-20) {
     2251          row[put]=row[j];
     2252          elementByColumn[put++]=value;
     2253        }
     2254      }
     2255      columnLength[iColumn]=put-start;
     2256    }
    21272257  }
    21282258  model->messageHandler()->message(CLP_PACKEDSCALE_INITIAL,*model->messagesPointer())
    21292259    <<smallest<<largest
    21302260    <<CoinMessageEol;
    2131   if (smallest>=0.5&&largest<=2.0) {
     2261  if (smallest>=0.5&&largest<=2.0&&!deletedElements) {
    21322262    // don't bother scaling
    21332263    model->messageHandler()->message(CLP_PACKEDSCALE_FORGET,*model->messagesPointer())
     
    21422272    delete [] usefulRow;
    21432273    delete [] usefulColumn;
    2144     if (!model->rowCopy())
    2145       delete rowCopyBase; // get rid of temporary row copy
    21462274    return 1;
    21472275  } else {
    2148       // need to scale
     2276#ifdef CLP_INVESTIGATE
     2277    if (deletedElements)
     2278      printf("DEL_ELS\n");
     2279#endif
     2280    if (!rowCopyBase) {
     2281      // temporary copy
     2282      rowCopyBase = reverseOrderedCopy();
     2283    } else if (deletedElements) {
     2284      rowCopyBase = reverseOrderedCopy();
     2285      model->setNewRowCopy(rowCopyBase);
     2286    }
     2287#ifndef NDEBUG
     2288    ClpPackedMatrix* rowCopy =
     2289      dynamic_cast< ClpPackedMatrix*>(rowCopyBase);
     2290    // Make sure it is really a ClpPackedMatrix
     2291    assert (rowCopy!=NULL);
     2292#else
     2293    ClpPackedMatrix* rowCopy =
     2294      static_cast< ClpPackedMatrix*>(rowCopyBase);
     2295#endif
     2296   
     2297    const int * column = rowCopy->getIndices();
     2298    const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
     2299    const double * element = rowCopy->getElements();
     2300    // need to scale
    21492301    if (largest>1.0e13*smallest) {
    21502302      // safer to have smaller zero tolerance
     
    21812333    bool finished=false;
    21822334    // if scalingMethod 3 then may change
     2335    bool extraDetails = (model->logLevel()>2);
    21832336    while (!finished) {
    21842337      int numberPass=3;
     
    22162369                rowScale[iRow]=1.0/largest;
    22172370#ifdef COIN_DEVELOP
    2218                 if (model->logLevel()>2) {
     2371                if (extraDetails) {
    22192372                  overallLargest = CoinMax(overallLargest,largest);
    22202373                  overallSmallest = CoinMin(overallSmallest,largest);
     
    22362389                  if (usefulColumn[iColumn]) {
    22372390                    double value = fabs(element[j]);
    2238                     // Don't bother with tiny elements
    2239                     if (value>1.0e-20) {
    2240                       value *= columnScale[iColumn];
    2241                       largest = CoinMax(largest,value);
    2242                       smallest = CoinMin(smallest,value);
    2243                     }
     2391                    value *= columnScale[iColumn];
     2392                    largest = CoinMax(largest,value);
     2393                    smallest = CoinMin(smallest,value);
    22442394                  }
    22452395                }
     
    22492399                }
    22502400#ifdef COIN_DEVELOP
    2251                 if (model->logLevel()>2) {
     2401                if (extraDetails) {
    22522402                  overallLargest = CoinMax(largest*rowScale[iRow],overallLargest);
    22532403                  overallSmallest = CoinMin(smallest*rowScale[iRow],overallSmallest);
     
    22782428            rowScale[iRow]=1.0/largest;
    22792429#ifdef COIN_DEVELOP
    2280             if (model->logLevel()>2) {
     2430            if (extraDetails) {
    22812431              overallLargest = CoinMax(overallLargest,largest);
    22822432              overallSmallest = CoinMin(overallSmallest,largest);
     
    23062456                if (usefulColumn[iColumn]) {
    23072457                  double value = fabs(element[j]);
    2308                   // Don't bother with tiny elements
    2309                   if (value>1.0e-20) {
    2310                     value *= columnScale[iColumn];
    2311                     largest = CoinMax(largest,value);
    2312                     smallest = CoinMin(smallest,value);
    2313                   }
     2458                  value *= columnScale[iColumn];
     2459                  largest = CoinMax(largest,value);
     2460                  smallest = CoinMin(smallest,value);
    23142461                }
    23152462              }
     
    23172464              rowScale[iRow]=1.0/sqrt(smallest*largest);
    23182465              //rowScale[iRow]=CoinMax(1.0e-10,CoinMin(1.0e10,rowScale[iRow]));
    2319               if (model->logLevel()>2) {
     2466              if (extraDetails) {
    23202467                overallLargest = CoinMax(largest*rowScale[iRow],overallLargest);
    23212468                overallSmallest = CoinMin(smallest*rowScale[iRow],overallSmallest);
     
    23292476            if (usefulColumn[iColumn]) {
    23302477              double value = fabs(objective[iColumn]);
    2331               // Don't bother with tiny elements
    2332               if (value>1.0e-20) {
    2333                 value *= columnScale[iColumn];
    2334                 largest = CoinMax(largest,value);
    2335                 smallest = CoinMin(smallest,value);
    2336               }
     2478              value *= columnScale[iColumn];
     2479              largest = CoinMax(largest,value);
     2480              smallest = CoinMin(smallest,value);
    23372481            }
    23382482          }
     
    23562500                iRow=row[j];
    23572501                double value = fabs(elementByColumn[j]);
    2358                 // Don't bother with tiny elements
    2359                 if (value>1.0e-20&&usefulRow[iRow]) {
     2502                if (usefulRow[iRow]) {
    23602503                  value *= rowScale[iRow];
    23612504                  largest = CoinMax(largest,value);
     
    24402583        }
    24412584#if 0
    2442         if (model->logLevel()>2) {
     2585        if (extraDetails) {
    24432586          if (finished)
    24442587            printf("equilibrium ratio %g, geometric ratio %g , geo chosen\n",
     
    24942637          // scaledDifference,difference*columnScale[iColumn]);
    24952638        }
    2496         overallSmallest = CoinMin(overallSmallest,smallest*columnScale[iColumn]);
     2639        double value = smallest*columnScale[iColumn];
     2640        if (overallSmallest>value)
     2641          overallSmallest = value;
     2642        //overallSmallest = CoinMin(overallSmallest,smallest*columnScale[iColumn]);
    24972643      }
    24982644    }
     
    25242670#endif
    25252671    if (quadraticObj) {
     2672      if (!rowCopyBase) {
     2673        // temporary copy
     2674        rowCopyBase = reverseOrderedCopy();
     2675      }
     2676#ifndef NDEBUG
     2677      ClpPackedMatrix* rowCopy =
     2678        dynamic_cast< ClpPackedMatrix*>(rowCopyBase);
     2679      // Make sure it is really a ClpPackedMatrix
     2680      assert (rowCopy!=NULL);
     2681#else
     2682      ClpPackedMatrix* rowCopy =
     2683        static_cast< ClpPackedMatrix*>(rowCopyBase);
     2684#endif
     2685      const int * column = rowCopy->getIndices();
     2686      const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
    25262687      CoinPackedMatrix * quadratic = quadraticObj->quadraticObjective();
    25272688      int numberXColumns = quadratic->getNumCols();
     
    25822743    }
    25832744#endif
    2584     if (inverseRowScale) {
    2585       // make copy (could do faster by using previous values)
    2586       // could just do partial
    2587       for (iRow=0;iRow<numberRows;iRow++)
    2588         inverseRowScale[iRow] = 1.0/rowScale[iRow] ;
    2589       for (iColumn=0;iColumn<numberColumns;iColumn++)
    2590         inverseColumnScale[iColumn] = 1.0/columnScale[iColumn] ;
    2591     }
     2745    // make copy (could do faster by using previous values)
     2746    // could just do partial
     2747    for (iRow=0;iRow<numberRows;iRow++)
     2748      inverseRowScale[iRow] = 1.0/rowScale[iRow] ;
     2749    for (iColumn=0;iColumn<numberColumns;iColumn++)
     2750      inverseColumnScale[iColumn] = 1.0/columnScale[iColumn] ;
    25922751    if (!arraysExist) {
    25932752      model->setRowScale(rowScale);
     
    25962755    if (model->rowCopy()) {
    25972756      // need to replace row by row
    2598       double * element = model->rowCopy()->getPackedMatrix()->getMutableElements();
     2757      ClpPackedMatrix* rowCopy =
     2758        static_cast< ClpPackedMatrix*>(model->rowCopy());
     2759      double * element = rowCopy->getMutableElements();
     2760      const int * column = rowCopy->getIndices();
     2761      const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
     2762      // scale row copy
     2763      for (iRow=0;iRow<numberRows;iRow++) {
     2764        CoinBigIndex j;
     2765        double scale = rowScale[iRow];
     2766        double * elementsInThisRow = element + rowStart[iRow];
     2767        const int * columnsInThisRow = column + rowStart[iRow];
     2768        int number = rowStart[iRow+1]-rowStart[iRow];
     2769        assert (number<=numberColumns);
     2770        for (j=0;j<number;j++) {
     2771          int iColumn = columnsInThisRow[j];
     2772          elementsInThisRow[j] *= scale*columnScale[iColumn];
     2773        }
     2774      }
     2775      if ((model->specialOptions()&262144)!=0) {
     2776      //if ((model->specialOptions()&(COIN_CBC_USING_CLP|16384))!=0) {
     2777      //if (model->inCbcBranchAndBound()&&false) {
     2778        // copy without gaps
     2779        CoinPackedMatrix * scaledMatrix = new CoinPackedMatrix(*matrix_,0,0);
     2780        ClpPackedMatrix * scaled = new ClpPackedMatrix(scaledMatrix);
     2781        model->setClpScaledMatrix(scaled);
     2782        // get matrix data pointers
     2783        const int * row = scaledMatrix->getIndices();
     2784        const CoinBigIndex * columnStart = scaledMatrix->getVectorStarts();
     2785#ifndef NDEBUG
     2786        const int * columnLength = scaledMatrix->getVectorLengths();
     2787#endif
     2788        double * elementByColumn = scaledMatrix->getMutableElements();
     2789        for (iColumn=0;iColumn<numberColumns;iColumn++) {
     2790          CoinBigIndex j;
     2791          double scale = columnScale[iColumn];
     2792          assert (columnStart[iColumn+1]==columnStart[iColumn]+columnLength[iColumn]);
     2793          for (j=columnStart[iColumn];
     2794               j<columnStart[iColumn+1];j++) {
     2795            int iRow = row[j];
     2796            elementByColumn[j] *= scale*rowScale[iRow];
     2797          }
     2798        }
     2799      } else {
     2800        //printf("not in b&b\n");
     2801      }
     2802    } else {
     2803      // no row copy
     2804      delete rowCopyBase;
     2805    }
     2806    return 0;
     2807  }
     2808}
     2809#endif
     2810//static int scale_stats[5]={0,0,0,0,0};
     2811// Creates scales for column copy (rowCopy in model may be modified)
     2812int
     2813ClpPackedMatrix::scale(ClpModel * model,const ClpSimplex * baseModel) const
     2814{
     2815  //const ClpSimplex * baseModel=NULL;
     2816  //return scale2(model);
     2817#if 0
     2818  ClpMatrixBase * rowClone=NULL;
     2819  if (model->rowCopy())
     2820    rowClone = model->rowCopy()->clone();
     2821  assert (!model->rowScale());
     2822  assert (!model->columnScale());
     2823  int returnCode=scale2(model);
     2824  if (returnCode)
     2825    return returnCode;
     2826#endif
     2827#ifndef NDEBUG
     2828  //checkFlags();
     2829#endif
     2830  int numberRows = model->numberRows();
     2831  int numberColumns = matrix_->getNumCols();
     2832  model->setClpScaledMatrix(NULL); // get rid of any scaled matrix
     2833  // If empty - return as sanityCheck will trap
     2834  if (!numberRows||!numberColumns) {
     2835    model->setRowScale(NULL);
     2836    model->setColumnScale(NULL);
     2837    return 1;
     2838  }
     2839#if 0
     2840  // start fake
     2841  double * rowScale2 = CoinCopyOfArray(model->rowScale(),numberRows);
     2842  double * columnScale2 = CoinCopyOfArray(model->columnScale(),numberColumns);
     2843  model->setRowScale(NULL);
     2844  model->setColumnScale(NULL);
     2845  model->setNewRowCopy(rowClone);
     2846#endif
     2847  ClpMatrixBase * rowCopyBase=model->rowCopy();
     2848  double * rowScale;
     2849  double * columnScale;
     2850  //assert (!model->rowScale());
     2851  bool arraysExist;
     2852  double * inverseRowScale = NULL;
     2853  double * inverseColumnScale = NULL;
     2854  if (!model->rowScale()) {
     2855    rowScale = new double [numberRows*2];
     2856    columnScale = new double [numberColumns*2];
     2857    inverseRowScale = rowScale+numberRows;
     2858    inverseColumnScale = columnScale+numberColumns;
     2859    arraysExist=false;
     2860  } else {
     2861    rowScale=model->mutableRowScale();
     2862    columnScale=model->mutableColumnScale();
     2863    inverseRowScale = model->mutableInverseRowScale();
     2864    inverseColumnScale = model->mutableInverseColumnScale();
     2865    arraysExist=true;
     2866  }
     2867  assert (inverseRowScale==rowScale+numberRows);
     2868  assert (inverseColumnScale==columnScale+numberColumns);
     2869  // we are going to mark bits we are interested in
     2870  char * usefulColumn = new char [numberColumns];
     2871  double * rowLower = model->rowLower();
     2872  double * rowUpper = model->rowUpper();
     2873  double * columnLower = model->columnLower();
     2874  double * columnUpper = model->columnUpper();
     2875  int iColumn, iRow;
     2876  //#define LEAVE_FIXED
     2877  // mark empty and fixed columns
     2878  // also see if worth scaling
     2879  assert (model->scalingFlag()<=4);
     2880  //  scale_stats[model->scalingFlag()]++;
     2881  double largest=0.0;
     2882  double smallest=1.0e50;
     2883  // get matrix data pointers
     2884  int * row = matrix_->getMutableIndices();
     2885  const CoinBigIndex * columnStart = matrix_->getVectorStarts();
     2886  int * columnLength = matrix_->getMutableVectorLengths();
     2887  double * elementByColumn = matrix_->getMutableElements();
     2888  bool deletedElements=false;
     2889  for (iColumn=0;iColumn<numberColumns;iColumn++) {
     2890    CoinBigIndex j;
     2891    char useful=0;
     2892    bool deleteSome = false;
     2893    int start = columnStart[iColumn];
     2894    int end = start + columnLength[iColumn];
     2895#ifndef LEAVE_FIXED
     2896    if (columnUpper[iColumn]>
     2897        columnLower[iColumn]+1.0e-12) {
     2898#endif
     2899      for (j=start;j<end;j++) {
     2900        iRow=row[j];
     2901        double value = fabs(elementByColumn[j]);
     2902        if (value>1.0e-20) {
     2903          useful=1;
     2904          largest = CoinMax(largest,fabs(elementByColumn[j]));
     2905          smallest = CoinMin(smallest,fabs(elementByColumn[j]));
     2906        } else {
     2907          // small
     2908          deleteSome=true;
     2909        }
     2910      }
     2911#ifndef LEAVE_FIXED
     2912    } else {
     2913      // just check values
     2914      for (j=start;j<end;j++) {
     2915        double value = fabs(elementByColumn[j]);
     2916        if (value<=1.0e-20) {
     2917          // small
     2918          deleteSome=true;
     2919        }
     2920      }
     2921    }
     2922#endif
     2923    usefulColumn[iColumn]=useful;
     2924    if (deleteSome) {
     2925      deletedElements=true;
     2926      CoinBigIndex put = start;
     2927      for (j=start;j<end;j++) {
     2928        double value = elementByColumn[j];
     2929        if (fabs(value)>1.0e-20) {
     2930          row[put]=row[j];
     2931          elementByColumn[put++]=value;
     2932        }
     2933      }
     2934      columnLength[iColumn]=put-start;
     2935    }
     2936  }
     2937  model->messageHandler()->message(CLP_PACKEDSCALE_INITIAL,*model->messagesPointer())
     2938    <<smallest<<largest
     2939    <<CoinMessageEol;
     2940  if (smallest>=0.5&&largest<=2.0&&!deletedElements) {
     2941    // don't bother scaling
     2942    model->messageHandler()->message(CLP_PACKEDSCALE_FORGET,*model->messagesPointer())
     2943      <<CoinMessageEol;
     2944    if (!arraysExist) {
     2945      delete [] rowScale;
     2946      delete [] columnScale;
     2947    } else {
     2948      model->setRowScale(NULL);
     2949      model->setColumnScale(NULL);
     2950    }
     2951    delete [] usefulColumn;
     2952    return 1;
     2953  } else {
     2954#ifdef CLP_INVESTIGATE
     2955    if (deletedElements)
     2956      printf("DEL_ELS\n");
     2957#endif
     2958    if (!rowCopyBase) {
     2959      // temporary copy
     2960      rowCopyBase = reverseOrderedCopy();
     2961    } else if (deletedElements) {
     2962      rowCopyBase = reverseOrderedCopy();
     2963      model->setNewRowCopy(rowCopyBase);
     2964    }
     2965#ifndef NDEBUG
     2966    ClpPackedMatrix* rowCopy =
     2967      dynamic_cast< ClpPackedMatrix*>(rowCopyBase);
     2968    // Make sure it is really a ClpPackedMatrix
     2969    assert (rowCopy!=NULL);
     2970#else
     2971    ClpPackedMatrix* rowCopy =
     2972      static_cast< ClpPackedMatrix*>(rowCopyBase);
     2973#endif
     2974   
     2975    const int * column = rowCopy->getIndices();
     2976    const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
     2977    const double * element = rowCopy->getElements();
     2978    // need to scale
     2979    if (largest>1.0e13*smallest) {
     2980      // safer to have smaller zero tolerance
     2981      double ratio = smallest/largest;
     2982      ClpSimplex * simplex = static_cast<ClpSimplex *> (model);
     2983      double newTolerance = CoinMax(ratio*0.5,1.0e-18);
     2984      if (simplex->zeroTolerance()>newTolerance)
     2985        simplex->setZeroTolerance(newTolerance);
     2986    }
     2987    int scalingMethod = model->scalingFlag();
     2988    if (scalingMethod==4) {
     2989      // As auto
     2990      scalingMethod=3;
     2991    }
     2992    double savedOverallRatio=0.0;
     2993    double tolerance = 5.0*model->primalTolerance();
     2994    double overallLargest=-1.0e-20;
     2995    double overallSmallest=1.0e20;
     2996    bool finished=false;
     2997    // if scalingMethod 3 then may change
     2998    bool extraDetails = (model->logLevel()>2);
     2999    while (!finished) {
     3000      int numberPass=3;
     3001      overallLargest=-1.0e-20;
     3002      overallSmallest=1.0e20;
     3003      ClpFillN ( rowScale, numberRows,1.0);
     3004      ClpFillN ( columnScale, numberColumns,1.0);
     3005      if (scalingMethod==1||scalingMethod==3) {
     3006        // Maximum in each row
     3007        for (iRow=0;iRow<numberRows;iRow++) {
     3008          CoinBigIndex j;
     3009          largest=1.0e-10;
     3010          for (j=rowStart[iRow];j<rowStart[iRow+1];j++) {
     3011            int iColumn = column[j];
     3012            if (usefulColumn[iColumn]) {
     3013              double value = fabs(element[j]);
     3014              largest = CoinMax(largest,value);
     3015              assert (largest<1.0e40);
     3016            }
     3017          }
     3018          rowScale[iRow]=1.0/largest;
     3019#ifdef COIN_DEVELOP
     3020          if (extraDetails) {
     3021            overallLargest = CoinMax(overallLargest,largest);
     3022            overallSmallest = CoinMin(overallSmallest,largest);
     3023          }
     3024#endif
     3025        }
     3026      } else {
     3027#ifdef USE_OBJECTIVE
     3028        // This will be used to help get scale factors
     3029        double * objective = new double[numberColumns];
     3030        CoinMemcpyN(model->costRegion(1),numberColumns,objective);
     3031        double objScale=1.0;
     3032#endif
     3033        while (numberPass) {
     3034          overallLargest=0.0;
     3035          overallSmallest=1.0e50;
     3036          numberPass--;
     3037          // Geometric mean on row scales
     3038          for (iRow=0;iRow<numberRows;iRow++) {
     3039            CoinBigIndex j;
     3040            largest=1.0e-20;
     3041            smallest=1.0e50;
     3042            for (j=rowStart[iRow];j<rowStart[iRow+1];j++) {
     3043              int iColumn = column[j];
     3044              if (usefulColumn[iColumn]) {
     3045                double value = fabs(element[j]);
     3046                value *= columnScale[iColumn];
     3047                largest = CoinMax(largest,value);
     3048                smallest = CoinMin(smallest,value);
     3049              }
     3050            }
     3051           
     3052            rowScale[iRow]=1.0/sqrt(smallest*largest);
     3053            //rowScale[iRow]=CoinMax(1.0e-10,CoinMin(1.0e10,rowScale[iRow]));
     3054            if (extraDetails) {
     3055              overallLargest = CoinMax(largest*rowScale[iRow],overallLargest);
     3056              overallSmallest = CoinMin(smallest*rowScale[iRow],overallSmallest);
     3057            }
     3058          }
     3059#ifdef USE_OBJECTIVE
     3060          largest=1.0e-20;
     3061          smallest=1.0e50;
     3062          for (iColumn=0;iColumn<numberColumns;iColumn++) {
     3063            if (usefulColumn[iColumn]) {
     3064              double value = fabs(objective[iColumn]);
     3065              value *= columnScale[iColumn];
     3066              largest = CoinMax(largest,value);
     3067              smallest = CoinMin(smallest,value);
     3068            }
     3069          }
     3070          objScale=1.0/sqrt(smallest*largest);
     3071#endif
     3072          model->messageHandler()->message(CLP_PACKEDSCALE_WHILE,*model->messagesPointer())
     3073            <<overallSmallest
     3074            <<overallLargest
     3075            <<CoinMessageEol;
     3076          // skip last column round
     3077          if (numberPass==1)
     3078            break;
     3079          // Geometric mean on column scales
     3080          for (iColumn=0;iColumn<numberColumns;iColumn++) {
     3081            if (usefulColumn[iColumn]) {
     3082              CoinBigIndex j;
     3083              largest=1.0e-20;
     3084              smallest=1.0e50;
     3085              for (j=columnStart[iColumn];
     3086                   j<columnStart[iColumn]+columnLength[iColumn];j++) {
     3087                iRow=row[j];
     3088                double value = fabs(elementByColumn[j]);
     3089                value *= rowScale[iRow];
     3090                largest = CoinMax(largest,value);
     3091                smallest = CoinMin(smallest,value);
     3092              }
     3093#ifdef USE_OBJECTIVE
     3094              if (fabs(objective[iColumn])>1.0e-20) {
     3095                double value = fabs(objective[iColumn])*objScale;
     3096                largest = CoinMax(largest,value);
     3097                smallest = CoinMin(smallest,value);
     3098              }
     3099#endif
     3100              columnScale[iColumn]=1.0/sqrt(smallest*largest);
     3101              //columnScale[iColumn]=CoinMax(1.0e-10,CoinMin(1.0e10,columnScale[iColumn]));
     3102            }
     3103          }
     3104        }
     3105#ifdef USE_OBJECTIVE
     3106        delete [] objective;
     3107        printf("obj scale %g - use it if you want\n",objScale);
     3108#endif
     3109      }
     3110      // If ranges will make horrid then scale
     3111      for (iRow=0;iRow<numberRows;iRow++) {
     3112        double difference = rowUpper[iRow]-rowLower[iRow];
     3113        double scaledDifference = difference*rowScale[iRow];
     3114        if (scaledDifference>tolerance&&scaledDifference<1.0e-4) {
     3115          // make gap larger
     3116          rowScale[iRow] *= 1.0e-4/scaledDifference;
     3117          rowScale[iRow]=CoinMax(1.0e-10,CoinMin(1.0e10,rowScale[iRow]));
     3118          //printf("Row %d difference %g scaled diff %g => %g\n",iRow,difference,
     3119          // scaledDifference,difference*rowScale[iRow]);
     3120        }
     3121      }
     3122      // final pass to scale columns so largest is reasonable
     3123      // See what smallest will be if largest is 1.0
     3124      overallSmallest=1.0e50;
     3125      for (iColumn=0;iColumn<numberColumns;iColumn++) {
     3126        if (usefulColumn[iColumn]) {
     3127          CoinBigIndex j;
     3128          largest=1.0e-20;
     3129          smallest=1.0e50;
     3130          for (j=columnStart[iColumn];
     3131               j<columnStart[iColumn]+columnLength[iColumn];j++) {
     3132            iRow=row[j];
     3133            double value = fabs(elementByColumn[j]*rowScale[iRow]);
     3134            largest = CoinMax(largest,value);
     3135            smallest = CoinMin(smallest,value);
     3136          }
     3137          if (overallSmallest*largest>smallest)
     3138            overallSmallest = smallest/largest;
     3139        }
     3140      }
     3141      if (scalingMethod==1||scalingMethod==2) {
     3142        finished=true;
     3143      } else if (savedOverallRatio==0.0&&scalingMethod!=4) {
     3144        savedOverallRatio=overallSmallest;
     3145        scalingMethod=4;
     3146      } else {
     3147        assert (scalingMethod==4);
     3148        if (overallSmallest>2.0*savedOverallRatio) {
     3149          finished=true; // geometric was better
     3150          if (model->scalingFlag()==4) {
     3151            // If in Branch and bound change
     3152            if ((model->specialOptions()&1024)!=0) {
     3153              model->scaling(2);
     3154            }
     3155          }
     3156        } else {
     3157          scalingMethod=1; // redo equilibrium
     3158          if (model->scalingFlag()==4) {
     3159            // If in Branch and bound change
     3160            if ((model->specialOptions()&1024)!=0) {
     3161              model->scaling(1);
     3162            }
     3163          }
     3164        }
     3165#if 0
     3166        if (extraDetails) {
     3167          if (finished)
     3168            printf("equilibrium ratio %g, geometric ratio %g , geo chosen\n",
     3169                   savedOverallRatio,overallSmallest);
     3170          else
     3171            printf("equilibrium ratio %g, geometric ratio %g , equi chosen\n",
     3172                   savedOverallRatio,overallSmallest);
     3173        }
     3174#endif
     3175      }
     3176    }
     3177    //#define RANDOMIZE
     3178#ifdef RANDOMIZE
     3179    // randomize by up to 10%
     3180    for (iRow=0;iRow<numberRows;iRow++) {
     3181      double value = 0.5-randomNumberGenerator_.randomDouble();//between -0.5 to + 0.5
     3182      rowScale[iRow] *= (1.0+0.1*value);
     3183    }
     3184#endif
     3185    overallLargest=1.0;
     3186    if (overallSmallest<1.0e-1)
     3187      overallLargest = 1.0/sqrt(overallSmallest);
     3188    overallLargest = CoinMin(100.0,overallLargest);
     3189    overallSmallest=1.0e50;
     3190    char * usedRow = reinterpret_cast<char *>(inverseRowScale);
     3191    memset(usedRow,0,numberRows);
     3192    //printf("scaling %d\n",model->scalingFlag());
     3193    for (iColumn=0;iColumn<numberColumns;iColumn++) {
     3194      if (columnUpper[iColumn]>
     3195          columnLower[iColumn]+1.0e-12) {
     3196        //if (usefulColumn[iColumn]) {
     3197        CoinBigIndex j;
     3198        largest=1.0e-20;
     3199        smallest=1.0e50;
     3200        for (j=columnStart[iColumn];
     3201             j<columnStart[iColumn]+columnLength[iColumn];j++) {
     3202          iRow=row[j];
     3203          usedRow[iRow]=1;
     3204          double value = fabs(elementByColumn[j]*rowScale[iRow]);
     3205          largest = CoinMax(largest,value);
     3206          smallest = CoinMin(smallest,value);
     3207        }
     3208        columnScale[iColumn]=overallLargest/largest;
     3209        //columnScale[iColumn]=CoinMax(1.0e-10,CoinMin(1.0e10,columnScale[iColumn]));
     3210#ifdef RANDOMIZE
     3211        double value = 0.5-randomNumberGenerator_.randomDouble();//between -0.5 to + 0.5
     3212        columnScale[iColumn] *= (1.0+0.1*value);
     3213#endif
     3214        double difference = columnUpper[iColumn]-columnLower[iColumn];
     3215        if (difference<1.0e-5*columnScale[iColumn]) {
     3216          // make gap larger
     3217          columnScale[iColumn] = difference/1.0e-5;
     3218          //printf("Column %d difference %g scaled diff %g => %g\n",iColumn,difference,
     3219          // scaledDifference,difference*columnScale[iColumn]);
     3220        }
     3221        double value = smallest*columnScale[iColumn];
     3222        if (overallSmallest>value)
     3223          overallSmallest = value;
     3224        //overallSmallest = CoinMin(overallSmallest,smallest*columnScale[iColumn]);
     3225      } else {
     3226        assert(columnScale[iColumn]==1.0);
     3227        //columnScale[iColumn]=1.0;
     3228      }
     3229    }
     3230    for (iRow=0;iRow<numberRows;iRow++) {
     3231      if (!usedRow[iRow]) {
     3232        rowScale[iRow] = 1.0;
     3233      }
     3234    }
     3235    model->messageHandler()->message(CLP_PACKEDSCALE_FINAL,*model->messagesPointer())
     3236      <<overallSmallest
     3237      <<overallLargest
     3238      <<CoinMessageEol;
     3239#if 0
     3240    {
     3241      for (iRow=0;iRow<numberRows;iRow++) {
     3242        assert (rowScale[iRow]==rowScale2[iRow]);
     3243      }
     3244      delete [] rowScale2;
     3245      for (iColumn=0;iColumn<numberColumns;iColumn++) {
     3246        assert (columnScale[iColumn]==columnScale2[iColumn]);
     3247      }
     3248      delete [] columnScale2;
     3249    }
     3250#endif
     3251    if (overallSmallest<1.0e-13) {
     3252      // Change factorization zero tolerance
     3253      double newTolerance = CoinMax(1.0e-15*(overallSmallest/1.0e-13),
     3254                                             1.0e-18);
     3255      ClpSimplex * simplex = static_cast<ClpSimplex *> (model);
     3256      if (simplex->factorization()->zeroTolerance()>newTolerance)
     3257        simplex->factorization()->zeroTolerance(newTolerance);
     3258      newTolerance = CoinMax(overallSmallest*0.5,1.0e-18);
     3259      simplex->setZeroTolerance(newTolerance);
     3260    }
     3261    delete [] usefulColumn;
     3262#ifndef SLIM_CLP
     3263    // If quadratic then make symmetric
     3264    ClpObjective * obj = model->objectiveAsObject();
     3265#ifndef NO_RTTI
     3266    ClpQuadraticObjective * quadraticObj = (dynamic_cast< ClpQuadraticObjective*>(obj));
     3267#else
     3268    ClpQuadraticObjective * quadraticObj = NULL;
     3269    if (obj->type()==2)
     3270      quadraticObj = (static_cast< ClpQuadraticObjective*>(obj));
     3271#endif
     3272    if (quadraticObj) {
     3273      if (!rowCopyBase) {
     3274        // temporary copy
     3275        rowCopyBase = reverseOrderedCopy();
     3276      }
     3277#ifndef NDEBUG
     3278      ClpPackedMatrix* rowCopy =
     3279        dynamic_cast< ClpPackedMatrix*>(rowCopyBase);
     3280      // Make sure it is really a ClpPackedMatrix
     3281      assert (rowCopy!=NULL);
     3282#else
     3283      ClpPackedMatrix* rowCopy =
     3284        static_cast< ClpPackedMatrix*>(rowCopyBase);
     3285#endif
     3286      const int * column = rowCopy->getIndices();
     3287      const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
     3288      CoinPackedMatrix * quadratic = quadraticObj->quadraticObjective();
     3289      int numberXColumns = quadratic->getNumCols();
     3290      if (numberXColumns<numberColumns) {
     3291        // we assume symmetric
     3292        int numberQuadraticColumns=0;
     3293        int i;
     3294        //const int * columnQuadratic = quadratic->getIndices();
     3295        //const int * columnQuadraticStart = quadratic->getVectorStarts();
     3296        const int * columnQuadraticLength = quadratic->getVectorLengths();
     3297        for (i=0;i<numberXColumns;i++) {
     3298          int length=columnQuadraticLength[i];
     3299#ifndef CORRECT_COLUMN_COUNTS
     3300          length=1;
     3301#endif
     3302          if (length)
     3303            numberQuadraticColumns++;
     3304        }
     3305        int numberXRows = numberRows-numberQuadraticColumns;
     3306        numberQuadraticColumns=0;
     3307        for (i=0;i<numberXColumns;i++) {
     3308          int length=columnQuadraticLength[i];
     3309#ifndef CORRECT_COLUMN_COUNTS
     3310          length=1;
     3311#endif
     3312          if (length) {
     3313            rowScale[numberQuadraticColumns+numberXRows] = columnScale[i];
     3314            numberQuadraticColumns++;
     3315          }
     3316        }   
     3317        int numberQuadraticRows=0;
     3318        for (i=0;i<numberXRows;i++) {
     3319          // See if any in row quadratic
     3320          CoinBigIndex j;
     3321          int numberQ=0;
     3322          for (j=rowStart[i];j<rowStart[i+1];j++) {
     3323            int iColumn = column[j];
     3324            if (columnQuadraticLength[iColumn])
     3325              numberQ++;
     3326          }
     3327#ifndef CORRECT_ROW_COUNTS
     3328          numberQ=1;
     3329#endif
     3330          if (numberQ) {
     3331            columnScale[numberQuadraticRows+numberXColumns] = rowScale[i];
     3332            numberQuadraticRows++;
     3333          }
     3334        }
     3335        // and make sure Sj okay
     3336        for (iColumn=numberQuadraticRows+numberXColumns;iColumn<numberColumns;iColumn++) {
     3337          CoinBigIndex j=columnStart[iColumn];
     3338          assert(columnLength[iColumn]==1);
     3339          int iRow=row[j];
     3340          double value = fabs(elementByColumn[j]*rowScale[iRow]);
     3341          columnScale[iColumn]=1.0/value;
     3342        }
     3343      }
     3344    }
     3345#endif
     3346    // make copy (could do faster by using previous values)
     3347    // could just do partial
     3348    for (iRow=0;iRow<numberRows;iRow++)
     3349      inverseRowScale[iRow] = 1.0/rowScale[iRow] ;
     3350    for (iColumn=0;iColumn<numberColumns;iColumn++)
     3351      inverseColumnScale[iColumn] = 1.0/columnScale[iColumn] ;
     3352    if (!arraysExist) {
     3353      model->setRowScale(rowScale);
     3354      model->setColumnScale(columnScale);
     3355    }
     3356    if (model->rowCopy()) {
     3357      // need to replace row by row
     3358      ClpPackedMatrix* rowCopy =
     3359        static_cast< ClpPackedMatrix*>(model->rowCopy());
     3360      double * element = rowCopy->getMutableElements();
     3361      const int * column = rowCopy->getIndices();
     3362      const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
    25993363      // scale row copy
    26003364      for (iRow=0;iRow<numberRows;iRow++) {
  • trunk/Clp/src/ClpPackedMatrix.hpp

    r1321 r1344  
    105105  /** Creates scales for column copy (rowCopy in model may be modified)
    106106      returns non-zero if no scaling done */
    107   virtual int scale(ClpModel * model, const ClpSimplex * baseModel=NULL) const ;
     107  virtual int scale(ClpModel * model,const ClpSimplex * baseModel=NULL) const ;
    108108  /** Scales rowCopy if column copy scaled
    109109      Only called if scales already exist */
  • trunk/Clp/src/ClpPresolve.cpp

    r1321 r1344  
    1515#include "ClpPackedMatrix.hpp"
    1616#include "ClpSimplex.hpp"
     17#include "ClpSimplexOther.hpp"
    1718#ifndef SLIM_CLP
    1819#include "ClpQuadraticObjective.hpp"
     
    272273  originalModel_->times(1.0,originalModel_->primalColumnSolution(),
    273274                        originalModel_->primalRowSolution());
    274   originalModel_->checkSolutionInternal();
     275  originalModel_->checkSolutionInternal(); 
     276  if (originalModel_->sumDualInfeasibilities()>1.0e-1) {
     277    // See if we can fix easily
     278    static_cast<ClpSimplexOther *> (originalModel_)->cleanupAfterPostsolve();
     279  }
    275280  // Messages
    276281  presolvedModel_->messageHandler()->message(COIN_PRESOLVE_POSTSOLVE,
     
    850855    int i;
    851856    for (i=0;i<prob.nrows_;i++) {
    852       if ((prob.rowstat_[i]&7)==1)
     857      if ((prob.rowstat_[i]&7)==1) {
    853858        nr++;
     859      } else if ((prob.rowstat_[i]&7)==2) {
     860        // at ub
     861        assert (prob.acts_[i]>prob.rup_[i]-1.0e-6);
     862      } else if ((prob.rowstat_[i]&7)==3) {
     863        // at lb
     864        assert (prob.acts_[i]<prob.rlo_[i]+1.0e-6);
     865      }
    854866    }
    855867    int nc=0;
  • trunk/Clp/src/ClpSimplex.cpp

    r1342 r1344  
    2929#include "CoinLpIO.hpp"
    3030#include <cfloat>
    31 #ifndef CLP_AUXILIARY_MODEL
    32 #define auxiliaryModel_ false
    33 #endif
    3431
    3532#include <string>
     
    9592  rowActivityWork_(NULL),
    9693  columnActivityWork_(NULL),
    97 #ifdef CLP_AUXILIARY_MODEL
    98   auxiliaryModel_(NULL),
    99 #endif
    10094  numberDualInfeasibilities_(0),
    10195  numberDualInfeasibilitiesWithoutFree_(0),
     
    211205  rowActivityWork_(NULL),
    212206  columnActivityWork_(NULL),
    213 #ifdef CLP_AUXILIARY_MODEL
    214   auxiliaryModel_(NULL),
    215 #endif
    216207  numberDualInfeasibilities_(0),
    217208  numberDualInfeasibilitiesWithoutFree_(0),
     
    364355  rowActivityWork_(NULL),
    365356  columnActivityWork_(NULL),
    366 #ifdef CLP_AUXILIARY_MODEL
    367   auxiliaryModel_(NULL),
    368 #endif
    369357  numberDualInfeasibilities_(0),
    370358  numberDualInfeasibilitiesWithoutFree_(0),
     
    408396  saveStatus_=NULL;
    409397  factorization_ = new ClpFactorization(*rhs->factorization_,-numberRows_);
     398  //factorization_ = new ClpFactorization(*rhs->factorization_,
     399  //                            rhs->factorization_->goDenseThreshold());
    410400  ClpDualRowDantzig * pivot =
    411401    dynamic_cast< ClpDualRowDantzig*>(rhs->dualRowPivot_);
     
    11031093    ClpPackedMatrix* clpMatrix =
    11041094      dynamic_cast< ClpPackedMatrix*>(matrix_);
    1105     if (clpMatrix&&rowScale_&&(clpMatrix->flags()&2)==0) {
     1095    double * saveRowScale = rowScale_;
     1096    //double * saveColumnScale = columnScale_;
     1097    if (scaledMatrix_) {
     1098      rowScale_=NULL;
     1099      clpMatrix = scaledMatrix_;
     1100    }
     1101    if (clpMatrix&&(clpMatrix->flags()&2)==0) {
    11061102      CoinIndexedVector * cVector = columnArray_[0];
    11071103      int * whichColumn = cVector->getIndices();
     
    11311127                                rowScale_,columnScale_,NULL);
    11321128    }
     1129    rowScale_ = saveRowScale;
     1130    //columnScale_ = saveColumnScale;
    11331131    ClpFillN(work,numberRows_,0.0);
    11341132    // Extended duals and check dual infeasibility
     
    18601858  int maximumPivots = factorization_->maximumPivots();
    18611859  int numberDense = factorization_->numberDense();
     1860  bool dontInvert = ((specialOptions_&16384)!=0&&numberIterations_*3>
     1861                     2*maximumIterations());
    18621862  if (numberPivots==maximumPivots||
    18631863      maximumPivots<2) {
     
    18741874    }
    18751875    return 1;
    1876   } else if (factorization_->timeToRefactorize()) {
     1876  } else if (factorization_->timeToRefactorize()&&!dontInvert) {
    18771877    //printf("ret after %d pivots\n",factorization_->pivots());
    18781878    return 1;
     
    19581958  rowActivityWork_(NULL),
    19591959  columnActivityWork_(NULL),
    1960 #ifdef CLP_AUXILIARY_MODEL
    1961   auxiliaryModel_(NULL),
    1962 #endif
    19631960  numberDualInfeasibilities_(0),
    19641961  numberDualInfeasibilitiesWithoutFree_(0),
     
    20672064  rowActivityWork_(NULL),
    20682065  columnActivityWork_(NULL),
    2069 #ifdef CLP_AUXILIARY_MODEL
    2070   auxiliaryModel_(NULL),
    2071 #endif
    20722066  numberDualInfeasibilities_(0),
    20732067  numberDualInfeasibilitiesWithoutFree_(0),
     
    21412135  dontFactorizePivots_ = rhs.dontFactorizePivots_;
    21422136  int numberRows2 = numberRows_+numberExtraRows_;
    2143 #ifdef CLP_AUXILIARY_MODEL
    2144   auxiliaryModel_ = NULL;
    2145 #endif
    21462137  moreSpecialOptions_ = rhs.moreSpecialOptions_;
    21472138  if ((whatsChanged_&1)!=0) {
     
    22152206    saveStatus_ = NULL;
    22162207  }
     2208  delete factorization_;
    22172209  if (rhs.factorization_) {
    2218     delete factorization_;
    22192210    factorization_ = new ClpFactorization(*rhs.factorization_,numberRows_);
    22202211  } else {
     
    23542345  if (!type) {
    23552346    // delete everything
    2356 #ifdef CLP_AUXILIARY_MODEL
    2357     delete auxiliaryModel_;
    2358     auxiliaryModel_ = NULL;
    2359 #endif
    23602347    setEmptyFactorization();
    23612348    delete [] pivotVariable_;
     
    23732360    // delete any size information in methods
    23742361    if (type>1) {
    2375       factorization_->clearArrays();
     2362      //assert (factorization_);
     2363      if (factorization_)
     2364        factorization_->clearArrays();
    23762365      delete [] pivotVariable_;
    23772366      pivotVariable_=NULL;
     
    28442833  }
    28452834  bool oldMatrix = ((startFinishOptions&4)!=0&&(whatsChanged_&1)!=0);
    2846 #ifdef CLP_AUXILIARY_MODEL
    2847   if (auxiliaryModel_) {
    2848     if (auxiliaryModel_->numberRows_==numberRows_&&
    2849         auxiliaryModel_->numberColumns_==numberColumns_&&
    2850         (whatsChanged_&511)==511) {
    2851       oldMatrix=true;
    2852     } else {
    2853       deleteAuxiliaryModel();
    2854       oldMatrix=false;
    2855     }
    2856   }
    2857 #endif
    28582835  if (what==63) {
    28592836    pivotRow_=-1;
     
    29032880  int numberRows2 = numberRows_+numberExtraRows_;
    29042881  int numberTotal = numberRows2+numberColumns_;
    2905   if ((specialOptions_&65536)!=0&&!auxiliaryModel_) {
     2882  if ((specialOptions_&65536)!=0) {
    29062883    assert (!numberExtraRows_);
    29072884    if (!cost_||numberRows2>maximumInternalRows_||
     
    29632940  int i;
    29642941  bool doSanityCheck=true;
    2965   if (what==63&&!auxiliaryModel_) {
     2942  if (what==63) {
    29662943    // We may want to switch stuff off for speed
    29672944    if ((specialOptions_&256)!=0)
     
    29902967    }
    29912968#if 0
    2992     if (what==63&&(specialOptions_&131072)!=0) {
     2969    if (what==63) {
    29932970      int k=rowScale_ ? 1 : 0;
    29942971      if (oldMatrix)
     
    30162993    inverseColumnScale_ = NULL;
    30172994    if (scalingFlag_>0&&!rowScale_) {
    3018       if ((specialOptions_&131072)!=0) {
    3019         inverseRowScale_ = rowScale_+numberRows2;
    3020         inverseColumnScale_ = columnScale_+numberColumns_;
    3021       }
    3022       if ((specialOptions_&65536)!=0&&!auxiliaryModel_) {
     2995      if ((specialOptions_&65536)!=0) {
    30232996        assert (!rowScale_);
    30242997        rowScale_ = savedRowScale_;
     
    30403013      if (matrix_->scale(this))
    30413014        scalingFlag_=-scalingFlag_; // not scaled after all
    3042       if ((specialOptions_&131072)!=0&&rowScale_&&!savedRowScale_) {
    3043         inverseRowScale_ = rowScale_+numberRows2;
    3044         inverseColumnScale_ = columnScale_+numberColumns_;
    3045       }
    30463015      if (rowScale_&&automaticScale_) {
    30473016        // try automatic scaling
     
    30583027          }
    30593028          if (columnLower_[i]>0.0||columnUpper_[i]<0.0) {
    3060             double scale = 1.0/columnScale_[i];
     3029            double scale = 1.0*inverseColumnScale_[i];
    30613030            //printf("%d %g %g %g %g\n",i,scale,lower_[i],upper_[i],largestRhs);
    30623031            if (columnLower_[i]>0)
     
    31613130      matrix_->scaleRowCopy(this);
    31623131    }
     3132    if (rowScale_&&!savedRowScale_) {
     3133      inverseRowScale_ = rowScale_+numberRows2;
     3134      inverseColumnScale_ = columnScale_+numberColumns_;
     3135    }
    31633136    // See if we can try for faster row copy
    31643137    if (makeRowCopy&&!oldMatrix) {
     
    32013174      delete [] solution_;
    32023175      solution_ = new double[numberTotal];
    3203 #ifdef CLP_AUXILIARY_MODEL
    3204     } else if (auxiliaryModel_) {
    3205       assert (!cost_);
    3206       cost_=auxiliaryModel_->cost_;
    3207       assert (!lower_);
    3208       lower_=auxiliaryModel_->lower_;
    3209       assert (!upper_);
    3210       upper_=auxiliaryModel_->upper_;
    3211       assert (!dj_);
    3212       dj_=auxiliaryModel_->dj_;
    3213       assert (!solution_);
    3214       solution_=auxiliaryModel_->solution_;
    3215       assert (!rowScale_);
    3216       assert (!columnScale_);
    3217 #endif
    32183176    }
    32193177    reducedCostWork_ = dj_;
     
    32293187  }
    32303188  if ((what&4)!=0) {
    3231 #ifdef CLP_AUXILIARY_MODEL
    3232     if (!auxiliaryModel_||(what==63&&(auxiliaryModel_->specialOptions_&4)==0))
    3233 #endif
    3234       {
    3235       double direction = optimizationDirection_*objectiveScale_;
    3236       const double * obj = objective();
    3237 #ifdef CLP_AUXILIARY_MODEL
    3238       const double * rowScale = auxiliaryModel_ ? auxiliaryModel_->rowScale_ : rowScale_;
    3239       const double * columnScale = auxiliaryModel_ ? auxiliaryModel_->columnScale_ : columnScale_;
    3240 #else
    3241       const double * rowScale = rowScale_;
    3242       const double * columnScale = columnScale_;
    3243 #endif
    3244       // and also scale by scale factors
    3245       if (rowScale) {
    3246         if (rowObjective_) {
    3247           for (i=0;i<numberRows_;i++)
    3248             rowObjectiveWork_[i] = rowObjective_[i]*direction/rowScale[i];
    3249         } else {
    3250           memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
    3251         }
    3252         // If scaled then do all columns later in one loop
    3253         if (what!=63||auxiliaryModel_) {
    3254           for (i=0;i<numberColumns_;i++) {
    3255             CoinAssert(fabs(obj[i])<1.0e25);
    3256             objectiveWork_[i] = obj[i]*direction*columnScale[i];
    3257           }
    3258         }
    3259       } else {
    3260         if (rowObjective_) {
    3261           for (i=0;i<numberRows_;i++)
    3262             rowObjectiveWork_[i] = rowObjective_[i]*direction;
    3263         } else {
    3264           memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
    3265         }
    3266         for (i=0;i<numberColumns_;i++) {
    3267           CoinAssert(fabs(obj[i])<1.0e25);
    3268           objectiveWork_[i] = obj[i]*direction;
    3269         }
    3270       }
    3271 #ifdef CLP_AUXILIARY_MODEL
    3272       if (auxiliaryModel_) {
    3273         // save costs
    3274         CoinMemcpyN(cost_,numberTotal,auxiliaryModel_->cost_+numberTotal);
    3275       }
    3276     } else {
    3277       // just copy
    3278       CoinMemcpyN(auxiliaryModel_->cost_+numberTotal,numberTotal,cost_);
    3279 #endif
    3280     }
    3281   }
    3282   if ((what&1)!=0) {
    3283 #ifdef CLP_AUXILIARY_MODEL
    3284     const double * rowScale = auxiliaryModel_ ? auxiliaryModel_->rowScale_ : rowScale_;
    3285     const double * columnScale = auxiliaryModel_ ? auxiliaryModel_->columnScale_ : columnScale_;
    3286 #else
     3189    double direction = optimizationDirection_*objectiveScale_;
     3190    const double * obj = objective();
    32873191    const double * rowScale = rowScale_;
    32883192    const double * columnScale = columnScale_;
    3289 #endif
     3193    // and also scale by scale factors
     3194    if (rowScale) {
     3195      if (rowObjective_) {
     3196        for (i=0;i<numberRows_;i++)
     3197          rowObjectiveWork_[i] = rowObjective_[i]*direction/rowScale[i];
     3198      } else {
     3199        memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
     3200      }
     3201      // If scaled then do all columns later in one loop
     3202      if (what!=63) {
     3203        for (i=0;i<numberColumns_;i++) {
     3204          CoinAssert(fabs(obj[i])<1.0e25);
     3205          objectiveWork_[i] = obj[i]*direction*columnScale[i];
     3206        }
     3207      }
     3208    } else {
     3209      if (rowObjective_) {
     3210        for (i=0;i<numberRows_;i++)
     3211          rowObjectiveWork_[i] = rowObjective_[i]*direction;
     3212      } else {
     3213        memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
     3214      }
     3215      for (i=0;i<numberColumns_;i++) {
     3216        CoinAssert(fabs(obj[i])<1.0e25);
     3217        objectiveWork_[i] = obj[i]*direction;
     3218      }
     3219    }
     3220  }
     3221  if ((what&1)!=0) {
     3222    const double * rowScale = rowScale_;
    32903223    // clean up any mismatches on infinity
    32913224    // and fix any variables with tiny gaps
    32923225    double primalTolerance=dblParam_[ClpPrimalTolerance];
    3293 #ifdef CLP_AUXILIARY_MODEL
    3294     if (!auxiliaryModel_) {
    3295 #endif
    3296       if(rowScale) {
    3297         // If scaled then do all columns later in one loop
    3298         if (what!=63) {
    3299           if (!inverseColumnScale_) {
    3300             for (i=0;i<numberColumns_;i++) {
    3301               double multiplier = rhsScale_/columnScale[i];
    3302               double lowerValue = columnLower_[i];
    3303               double upperValue = columnUpper_[i];
    3304               if (lowerValue>-1.0e20) {
    3305                 columnLowerWork_[i]=lowerValue*multiplier;
    3306                 if (upperValue>=1.0e20) {
    3307                   columnUpperWork_[i]=COIN_DBL_MAX;
     3226    if(rowScale) {
     3227      // If scaled then do all columns later in one loop
     3228      if (what!=63) {
     3229        const double * inverseScale = inverseColumnScale_;
     3230        for (i=0;i<numberColumns_;i++) {
     3231          double multiplier = rhsScale_*inverseScale[i];
     3232          double lowerValue = columnLower_[i];
     3233          double upperValue = columnUpper_[i];
     3234          if (lowerValue>-1.0e20) {
     3235            columnLowerWork_[i]=lowerValue*multiplier;
     3236            if (upperValue>=1.0e20) {
     3237              columnUpperWork_[i]=COIN_DBL_MAX;
     3238            } else {
     3239              columnUpperWork_[i]=upperValue*multiplier;
     3240              if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3241                if (columnLowerWork_[i]>=0.0) {
     3242                  columnUpperWork_[i] = columnLowerWork_[i];
     3243                } else if (columnUpperWork_[i]<=0.0) {
     3244                  columnLowerWork_[i] = columnUpperWork_[i];
    33083245                } else {
    3309                   columnUpperWork_[i]=upperValue*multiplier;
    3310                   if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    3311                     if (columnLowerWork_[i]>=0.0) {
    3312                       columnUpperWork_[i] = columnLowerWork_[i];
    3313                     } else if (columnUpperWork_[i]<=0.0) {
    3314                       columnLowerWork_[i] = columnUpperWork_[i];
    3315                     } else {
    3316                       columnUpperWork_[i] = 0.0;
    3317                       columnLowerWork_[i] = 0.0;
    3318                     }
    3319                   }
     3246                  columnUpperWork_[i] = 0.0;
     3247                  columnLowerWork_[i] = 0.0;
    33203248                }
    3321               } else if (upperValue<1.0e20) {
    3322                 columnLowerWork_[i]=-COIN_DBL_MAX;
    3323                 columnUpperWork_[i]=upperValue*multiplier;
    3324               } else {
    3325                 // free
    3326                 columnLowerWork_[i]=-COIN_DBL_MAX;
    3327                 columnUpperWork_[i]=COIN_DBL_MAX;
    33283249              }
    33293250            }
     3251          } else if (upperValue<1.0e20) {
     3252            columnLowerWork_[i]=-COIN_DBL_MAX;
     3253            columnUpperWork_[i]=upperValue*multiplier;
    33303254          } else {
    3331             assert (inverseColumnScale_);
    3332             const double * inverseScale = inverseColumnScale_;
    3333             for (i=0;i<numberColumns_;i++) {
    3334               double multiplier = rhsScale_*inverseScale[i];
    3335               double lowerValue = columnLower_[i];
    3336               double upperValue = columnUpper_[i];
    3337               if (lowerValue>-1.0e20) {
    3338                 columnLowerWork_[i]=lowerValue*multiplier;
    3339                 if (upperValue>=1.0e20) {
    3340                   columnUpperWork_[i]=COIN_DBL_MAX;
    3341                 } else {
    3342                   columnUpperWork_[i]=upperValue*multiplier;
    3343                   if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    3344                     if (columnLowerWork_[i]>=0.0) {
    3345                       columnUpperWork_[i] = columnLowerWork_[i];
    3346                     } else if (columnUpperWork_[i]<=0.0) {
    3347                       columnLowerWork_[i] = columnUpperWork_[i];
    3348                     } else {
    3349                       columnUpperWork_[i] = 0.0;
    3350                       columnLowerWork_[i] = 0.0;
    3351                     }
    3352                   }
    3353                 }
    3354               } else if (upperValue<1.0e20) {
    3355                 columnLowerWork_[i]=-COIN_DBL_MAX;
    3356                 columnUpperWork_[i]=upperValue*multiplier;
     3255            // free
     3256            columnLowerWork_[i]=-COIN_DBL_MAX;
     3257            columnUpperWork_[i]=COIN_DBL_MAX;
     3258          }
     3259        }
     3260      }
     3261      for (i=0;i<numberRows_;i++) {
     3262        double multiplier = rhsScale_*rowScale[i];
     3263        double lowerValue = rowLower_[i];
     3264        double upperValue = rowUpper_[i];
     3265        if (lowerValue>-1.0e20) {
     3266          rowLowerWork_[i]=lowerValue*multiplier;
     3267          if (upperValue>=1.0e20) {
     3268            rowUpperWork_[i]=COIN_DBL_MAX;
     3269          } else {
     3270            rowUpperWork_[i]=upperValue*multiplier;
     3271            if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
     3272              if (rowLowerWork_[i]>=0.0) {
     3273                rowUpperWork_[i] = rowLowerWork_[i];
     3274              } else if (rowUpperWork_[i]<=0.0) {
     3275                rowLowerWork_[i] = rowUpperWork_[i];
    33573276              } else {
    3358                 // free
    3359                 columnLowerWork_[i]=-COIN_DBL_MAX;
    3360                 columnUpperWork_[i]=COIN_DBL_MAX;
     3277                rowUpperWork_[i] = 0.0;
     3278                rowLowerWork_[i] = 0.0;
    33613279              }
    33623280            }
    33633281          }
    3364         }
    3365         for (i=0;i<numberRows_;i++) {
    3366           double multiplier = rhsScale_*rowScale[i];
    3367           double lowerValue = rowLower_[i];
    3368           double upperValue = rowUpper_[i];
    3369           if (lowerValue>-1.0e20) {
    3370             rowLowerWork_[i]=lowerValue*multiplier;
    3371             if (upperValue>=1.0e20) {
    3372               rowUpperWork_[i]=COIN_DBL_MAX;
    3373             } else {
    3374               rowUpperWork_[i]=upperValue*multiplier;
    3375               if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
    3376                 if (rowLowerWork_[i]>=0.0) {
    3377                   rowUpperWork_[i] = rowLowerWork_[i];
    3378                 } else if (rowUpperWork_[i]<=0.0) {
    3379                   rowLowerWork_[i] = rowUpperWork_[i];
    3380                 } else {
    3381                   rowUpperWork_[i] = 0.0;
    3382                   rowLowerWork_[i] = 0.0;
    3383                 }
    3384               }
    3385             }
    3386           } else if (upperValue<1.0e20) {
    3387             rowLowerWork_[i]=-COIN_DBL_MAX;
    3388             rowUpperWork_[i]=upperValue*multiplier;
    3389           } else {
    3390             // free
    3391             rowLowerWork_[i]=-COIN_DBL_MAX;
    3392             rowUpperWork_[i]=COIN_DBL_MAX;
    3393           }
    3394         }
    3395       } else if (rhsScale_!=1.0) {
    3396         for (i=0;i<numberColumns_;i++) {
    3397           double lowerValue = columnLower_[i];
    3398           double upperValue = columnUpper_[i];
    3399           if (lowerValue>-1.0e20) {
    3400             columnLowerWork_[i]=lowerValue*rhsScale_;
    3401             if (upperValue>=1.0e20) {
    3402               columnUpperWork_[i]=COIN_DBL_MAX;
    3403             } else {
    3404               columnUpperWork_[i]=upperValue*rhsScale_;
    3405               if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    3406                 if (columnLowerWork_[i]>=0.0) {
    3407                   columnUpperWork_[i] = columnLowerWork_[i];
    3408                 } else if (columnUpperWork_[i]<=0.0) {
    3409                   columnLowerWork_[i] = columnUpperWork_[i];
    3410                 } else {
    3411                   columnUpperWork_[i] = 0.0;
    3412                   columnLowerWork_[i] = 0.0;
    3413                 }
    3414               }
    3415             }
    3416           } else if (upperValue<1.0e20) {
    3417             columnLowerWork_[i]=-COIN_DBL_MAX;
    3418             columnUpperWork_[i]=upperValue*rhsScale_;
    3419           } else {
    3420             // free
    3421             columnLowerWork_[i]=-COIN_DBL_MAX;
    3422             columnUpperWork_[i]=COIN_DBL_MAX;
    3423           }
    3424         }
    3425         for (i=0;i<numberRows_;i++) {
    3426           double lowerValue = rowLower_[i];
    3427           double upperValue = rowUpper_[i];
    3428           if (lowerValue>-1.0e20) {
    3429             rowLowerWork_[i]=lowerValue*rhsScale_;
    3430             if (upperValue>=1.0e20) {
    3431               rowUpperWork_[i]=COIN_DBL_MAX;
    3432             } else {
    3433               rowUpperWork_[i]=upperValue*rhsScale_;
    3434               if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
    3435                 if (rowLowerWork_[i]>=0.0) {
    3436                   rowUpperWork_[i] = rowLowerWork_[i];
    3437                 } else if (rowUpperWork_[i]<=0.0) {
    3438                   rowLowerWork_[i] = rowUpperWork_[i];
    3439                 } else {
    3440                   rowUpperWork_[i] = 0.0;
    3441                   rowLowerWork_[i] = 0.0;
    3442                 }
    3443               }
    3444             }
    3445           } else if (upperValue<1.0e20) {
    3446             rowLowerWork_[i]=-COIN_DBL_MAX;
    3447             rowUpperWork_[i]=upperValue*rhsScale_;
    3448           } else {
    3449             // free
    3450             rowLowerWork_[i]=-COIN_DBL_MAX;
    3451             rowUpperWork_[i]=COIN_DBL_MAX;
    3452           }
    3453         }
    3454       } else {
    3455         for (i=0;i<numberColumns_;i++) {
    3456           double lowerValue = columnLower_[i];
    3457           double upperValue = columnUpper_[i];
    3458           if (lowerValue>-1.0e20) {
    3459             columnLowerWork_[i]=lowerValue;
    3460             if (upperValue>=1.0e20) {
    3461               columnUpperWork_[i]=COIN_DBL_MAX;
    3462             } else {
    3463               columnUpperWork_[i]=upperValue;
    3464               if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    3465                 if (columnLowerWork_[i]>=0.0) {
    3466                   columnUpperWork_[i] = columnLowerWork_[i];
    3467                 } else if (columnUpperWork_[i]<=0.0) {
    3468                   columnLowerWork_[i] = columnUpperWork_[i];
    3469                 } else {
    3470                   columnUpperWork_[i] = 0.0;
    3471                   columnLowerWork_[i] = 0.0;
    3472                 }
    3473               }
    3474             }
    3475           } else if (upperValue<1.0e20) {
    3476             columnLowerWork_[i]=-COIN_DBL_MAX;
    3477             columnUpperWork_[i]=upperValue;
    3478           } else {
    3479             // free
    3480             columnLowerWork_[i]=-COIN_DBL_MAX;
    3481             columnUpperWork_[i]=COIN_DBL_MAX;
    3482           }
    3483         }
    3484         for (i=0;i<numberRows_;i++) {
    3485           double lowerValue = rowLower_[i];
    3486           double upperValue = rowUpper_[i];
    3487           if (lowerValue>-1.0e20) {
    3488             rowLowerWork_[i]=lowerValue;
    3489             if (upperValue>=1.0e20) {
    3490               rowUpperWork_[i]=COIN_DBL_MAX;
    3491             } else {
    3492               rowUpperWork_[i]=upperValue;
    3493               if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
    3494                 if (rowLowerWork_[i]>=0.0) {
    3495                   rowUpperWork_[i] = rowLowerWork_[i];
    3496                 } else if (rowUpperWork_[i]<=0.0) {
    3497                   rowLowerWork_[i] = rowUpperWork_[i];
    3498                 } else {
    3499                   rowUpperWork_[i] = 0.0;
    3500                   rowLowerWork_[i] = 0.0;
    3501                 }
    3502               }
    3503             }
    3504           } else if (upperValue<1.0e20) {
    3505             rowLowerWork_[i]=-COIN_DBL_MAX;
    3506             rowUpperWork_[i]=upperValue;
    3507           } else {
    3508             // free
    3509             rowLowerWork_[i]=-COIN_DBL_MAX;
    3510             rowUpperWork_[i]=COIN_DBL_MAX;
    3511           }
    3512         }
    3513       }
    3514 #ifdef CLP_AUXILIARY_MODEL
     3282        } else if (upperValue<1.0e20) {
     3283          rowLowerWork_[i]=-COIN_DBL_MAX;
     3284          rowUpperWork_[i]=upperValue*multiplier;
     3285        } else {
     3286          // free
     3287          rowLowerWork_[i]=-COIN_DBL_MAX;
     3288          rowUpperWork_[i]=COIN_DBL_MAX;
     3289        }
     3290      }
     3291    } else if (rhsScale_!=1.0) {
     3292      for (i=0;i<numberColumns_;i++) {
     3293        double lowerValue = columnLower_[i];
     3294        double upperValue = columnUpper_[i];
     3295        if (lowerValue>-1.0e20) {
     3296          columnLowerWork_[i]=lowerValue*rhsScale_;
     3297          if (upperValue>=1.0e20) {
     3298            columnUpperWork_[i]=COIN_DBL_MAX;
     3299          } else {
     3300            columnUpperWork_[i]=upperValue*rhsScale_;
     3301            if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3302              if (columnLowerWork_[i]>=0.0) {
     3303                columnUpperWork_[i] = columnLowerWork_[i];
     3304              } else if (columnUpperWork_[i]<=0.0) {
     3305                columnLowerWork_[i] = columnUpperWork_[i];
     3306              } else {
     3307                columnUpperWork_[i] = 0.0;
     3308                columnLowerWork_[i] = 0.0;
     3309              }
     3310            }
     3311          }
     3312        } else if (upperValue<1.0e20) {
     3313          columnLowerWork_[i]=-COIN_DBL_MAX;
     3314          columnUpperWork_[i]=upperValue*rhsScale_;
     3315        } else {
     3316          // free
     3317          columnLowerWork_[i]=-COIN_DBL_MAX;
     3318          columnUpperWork_[i]=COIN_DBL_MAX;
     3319        }
     3320      }
     3321      for (i=0;i<numberRows_;i++) {
     3322        double lowerValue = rowLower_[i];
     3323        double upperValue = rowUpper_[i];
     3324        if (lowerValue>-1.0e20) {
     3325          rowLowerWork_[i]=lowerValue*rhsScale_;
     3326          if (upperValue>=1.0e20) {
     3327            rowUpperWork_[i]=COIN_DBL_MAX;
     3328          } else {
     3329            rowUpperWork_[i]=upperValue*rhsScale_;
     3330            if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
     3331              if (rowLowerWork_[i]>=0.0) {
     3332                rowUpperWork_[i] = rowLowerWork_[i];
     3333              } else if (rowUpperWork_[i]<=0.0) {
     3334                rowLowerWork_[i] = rowUpperWork_[i];
     3335              } else {
     3336                rowUpperWork_[i] = 0.0;
     3337                rowLowerWork_[i] = 0.0;
     3338              }
     3339            }
     3340          }
     3341        } else if (upperValue<1.0e20) {
     3342          rowLowerWork_[i]=-COIN_DBL_MAX;
     3343          rowUpperWork_[i]=upperValue*rhsScale_;
     3344        } else {
     3345          // free
     3346          rowLowerWork_[i]=-COIN_DBL_MAX;
     3347          rowUpperWork_[i]=COIN_DBL_MAX;
     3348        }
     3349      }
    35153350    } else {
    3516       // auxiliary model
    3517       if (what!=63) {
    3518         // just copy
    3519         CoinMemcpyN(auxiliaryModel_->lower_+numberTotal,numberTotal,lower_);
    3520         CoinMemcpyN(auxiliaryModel_->upper_+numberTotal,numberTotal,upper_);
    3521       } else {
    3522         assert (rhsScale_==1.0);
    3523         assert (objectiveScale_==1.0);
    3524         if ((auxiliaryModel_->specialOptions_&2)==0) {
    3525           // need to do column bounds
    3526           if (columnScale) {
    3527             const double * inverseScale = columnScale+numberColumns_;
    3528             // scaled
    3529             for (i=0;i<numberColumns_;i++) {
    3530               double value;
    3531               value = columnLower_[i];
    3532               if (value>-1.0e20)
    3533                 value *= inverseScale[i];
    3534               lower_[i] = value;
    3535               value = columnUpper_[i];
    3536               if (value<1.0e20)
    3537                 value *= inverseScale[i];
    3538               upper_[i] = value;
    3539             }
    3540           } else {
    3541             for (i=0;i<numberColumns_;i++) {
    3542               lower_[i] = columnLower_[i];
    3543               upper_[i] = columnUpper_[i];
    3544             }
    3545           }
    3546           // save
    3547           CoinMemcpyN(lower_,numberColumns_,auxiliaryModel_->lower_+numberTotal);
    3548           CoinMemcpyN(upper_,numberColumns_,auxiliaryModel_->upper_+numberTotal);
    3549         } else {
    3550           // just copy
    3551           CoinMemcpyN(auxiliaryModel_->lower_+numberTotal,numberColumns_,lower_);
    3552           CoinMemcpyN(auxiliaryModel_->upper_+numberTotal,numberColumns_,upper_);
    3553         }
    3554         if ((auxiliaryModel_->specialOptions_&1)==0) {
    3555           // need to do row bounds
    3556           if (rowScale) {
    3557             // scaled
    3558             for (i=0;i<numberRows_;i++) {
    3559               double value;
    3560               value = rowLower_[i];
    3561               if (value>-1.0e20)
    3562                 value *= rowScale[i];
    3563               lower_[i+numberColumns_] = value;
    3564               value = rowUpper_[i];
    3565               if (value<1.0e20)
    3566                 value *= rowScale[i];
    3567               upper_[i+numberColumns_] = value;
    3568             }
    3569           } else {
    3570             for (i=0;i<numberRows_;i++) {
    3571               lower_[i+numberColumns_] = rowLower_[i];
    3572               upper_[i+numberColumns_] = rowUpper_[i];
    3573             }
    3574           }
    3575           // save
    3576           CoinMemcpyN(                 lower_+numberColumns_,numberRows_,auxiliaryModel_->lower_+numberTotal+numberColumns_);
    3577           CoinMemcpyN(                 upper_+numberColumns_,numberRows_,auxiliaryModel_->upper_+numberTotal+numberColumns_);
    3578         } else {
    3579           // just copy
    3580           CoinMemcpyN(                 auxiliaryModel_->lower_+numberTotal+numberColumns_,
    3581                        numberRows_,lower_+numberColumns_);
    3582           CoinMemcpyN(                 auxiliaryModel_->upper_+numberTotal+numberColumns_,
    3583                        numberRows_,upper_+numberColumns_);
    3584         }
    3585       }
    3586       if (what==63) {
    3587         // do basis
    3588         assert ((auxiliaryModel_->specialOptions_&8)!=0);
    3589         // clean solution trusting basis
    3590         for ( i=0;i<numberTotal;i++) {
    3591           Status status =getStatus(i);
    3592           double value=solution_[i]; // may as well leave if basic (values pass)
    3593           if (status!=basic) {
    3594             if (upper_[i]==lower_[i]) {
    3595               setStatus(i,isFixed);
    3596               value=lower_[i];
    3597             } else if (status==atLowerBound) {
    3598               if (lower_[i]>-1.0e20) {
    3599                 value=lower_[i];
    3600               } else {
    3601                 if (upper_[i]<1.0e20) {
    3602                   value=upper_[i];
    3603                   setStatus(i,atUpperBound);
    3604                 } else {
    3605                   setStatus(i,isFree);
    3606                 }
    3607               }
    3608             } else if (status==atUpperBound) {
    3609               if (upper_[i]<1.0e20) {
    3610               value=upper_[i];
    3611               } else {
    3612                 if (lower_[i]>-1.0e20) {
    3613                   value=lower_[i];
    3614                   setStatus(i,atLowerBound);
    3615                 } else {
    3616                   setStatus(i,isFree);
    3617                 }
    3618               }
    3619             }
    3620           }
    3621           solution_[i]=value;
    3622         }
    3623       }
    3624     }
    3625 #endif
     3351      for (i=0;i<numberColumns_;i++) {
     3352        double lowerValue = columnLower_[i];
     3353        double upperValue = columnUpper_[i];
     3354        if (lowerValue>-1.0e20) {
     3355          columnLowerWork_[i]=lowerValue;
     3356          if (upperValue>=1.0e20) {
     3357            columnUpperWork_[i]=COIN_DBL_MAX;
     3358          } else {
     3359            columnUpperWork_[i]=upperValue;
     3360            if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3361              if (columnLowerWork_[i]>=0.0) {
     3362                columnUpperWork_[i] = columnLowerWork_[i];
     3363              } else if (columnUpperWork_[i]<=0.0) {
     3364                columnLowerWork_[i] = columnUpperWork_[i];
     3365              } else {
     3366                columnUpperWork_[i] = 0.0;
     3367                columnLowerWork_[i] = 0.0;
     3368              }
     3369            }
     3370          }
     3371        } else if (upperValue<1.0e20) {
     3372          columnLowerWork_[i]=-COIN_DBL_MAX;
     3373          columnUpperWork_[i]=upperValue;
     3374        } else {
     3375          // free
     3376          columnLowerWork_[i]=-COIN_DBL_MAX;
     3377          columnUpperWork_[i]=COIN_DBL_MAX;
     3378        }
     3379      }
     3380      for (i=0;i<numberRows_;i++) {
     3381        double lowerValue = rowLower_[i];
     3382        double upperValue = rowUpper_[i];
     3383        if (lowerValue>-1.0e20) {
     3384          rowLowerWork_[i]=lowerValue;
     3385          if (upperValue>=1.0e20) {
     3386            rowUpperWork_[i]=COIN_DBL_MAX;
     3387          } else {
     3388            rowUpperWork_[i]=upperValue;
     3389            if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
     3390              if (rowLowerWork_[i]>=0.0) {
     3391                rowUpperWork_[i] = rowLowerWork_[i];
     3392              } else if (rowUpperWork_[i]<=0.0) {
     3393                rowLowerWork_[i] = rowUpperWork_[i];
     3394              } else {
     3395                rowUpperWork_[i] = 0.0;
     3396                rowLowerWork_[i] = 0.0;
     3397              }
     3398            }
     3399          }
     3400        } else if (upperValue<1.0e20) {
     3401          rowLowerWork_[i]=-COIN_DBL_MAX;
     3402          rowUpperWork_[i]=upperValue;
     3403        } else {
     3404          // free
     3405          rowLowerWork_[i]=-COIN_DBL_MAX;
     3406          rowUpperWork_[i]=COIN_DBL_MAX;
     3407        }
     3408      }
     3409    }
    36263410  }
    36273411  if (what==63) {
     
    36313415    if (direction)
    36323416      direction = 1.0/direction;
    3633 #ifdef CLP_AUXILIARY_MODEL
    3634     if (direction!=1.0&&(!auxiliaryModel_||(auxiliaryModel_->specialOptions_&8)==0))
    3635 #else
    3636     if (direction!=1.0)
    3637 #endif
    3638  {
     3417    if (direction!=1.0) {
    36393418      // reverse all dual signs
    36403419      for (i=0;i<numberColumns_;i++)
     
    36463425      setFakeBound(i,noFake);
    36473426    }
    3648 #ifdef CLP_AUXILIARY_MODEL
    3649     if (!auxiliaryModel_) {
    3650 #endif
    3651       if (rowScale_) {
    3652         const double * obj = objective();
    3653         double direction = optimizationDirection_*objectiveScale_;
    3654         // clean up any mismatches on infinity
    3655         // and fix any variables with tiny gaps
    3656         double primalTolerance=dblParam_[ClpPrimalTolerance];
    3657         // on entry
    3658         if (!inverseColumnScale_) {
    3659           for (i=0;i<numberColumns_;i++) {
    3660             CoinAssert(fabs(obj[i])<1.0e25);
    3661             double scaleFactor = columnScale_[i];
    3662             double multiplier = rhsScale_/scaleFactor;
    3663             scaleFactor *= direction;
    3664             objectiveWork_[i] = obj[i]*scaleFactor;
    3665             reducedCostWork_[i] = reducedCost_[i]*scaleFactor;
    3666             double lowerValue = columnLower_[i];
    3667             double upperValue = columnUpper_[i];
    3668             if (lowerValue>-1.0e20) {
    3669               columnLowerWork_[i]=lowerValue*multiplier;
    3670               if (upperValue>=1.0e20) {
    3671                 columnUpperWork_[i]=COIN_DBL_MAX;
     3427    if (rowScale_) {
     3428      const double * obj = objective();
     3429      double direction = optimizationDirection_*objectiveScale_;
     3430      // clean up any mismatches on infinity
     3431      // and fix any variables with tiny gaps
     3432      double primalTolerance=dblParam_[ClpPrimalTolerance];
     3433      // on entry
     3434      const double * inverseScale = inverseColumnScale_;
     3435      for (i=0;i<numberColumns_;i++) {
     3436        CoinAssert(fabs(obj[i])<1.0e25);
     3437        double scaleFactor = columnScale_[i];
     3438        double multiplier = rhsScale_*inverseScale[i];
     3439        scaleFactor *= direction;
     3440        objectiveWork_[i] = obj[i]*scaleFactor;
     3441        reducedCostWork_[i] = reducedCost_[i]*scaleFactor;
     3442        double lowerValue = columnLower_[i];
     3443        double upperValue = columnUpper_[i];
     3444        if (lowerValue>-1.0e20) {
     3445          columnLowerWork_[i]=lowerValue*multiplier;
     3446          if (upperValue>=1.0e20) {
     3447            columnUpperWork_[i]=COIN_DBL_MAX;
     3448          } else {
     3449            columnUpperWork_[i]=upperValue*multiplier;
     3450            if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3451              if (columnLowerWork_[i]>=0.0) {
     3452                columnUpperWork_[i] = columnLowerWork_[i];
     3453              } else if (columnUpperWork_[i]<=0.0) {
     3454                columnLowerWork_[i] = columnUpperWork_[i];
    36723455              } else {
    3673                 columnUpperWork_[i]=upperValue*multiplier;
    3674                 if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    3675                   if (columnLowerWork_[i]>=0.0) {
    3676                     columnUpperWork_[i] = columnLowerWork_[i];
    3677                   } else if (columnUpperWork_[i]<=0.0) {
    3678                     columnLowerWork_[i] = columnUpperWork_[i];
    3679                   } else {
    3680                     columnUpperWork_[i] = 0.0;
    3681                     columnLowerWork_[i] = 0.0;
    3682                   }
    3683                 }
    3684               }
    3685             } else if (upperValue<1.0e20) {
    3686               columnLowerWork_[i]=-COIN_DBL_MAX;
    3687               columnUpperWork_[i]=upperValue*multiplier;
    3688             } else {
    3689               // free
    3690               columnLowerWork_[i]=-COIN_DBL_MAX;
    3691               columnUpperWork_[i]=COIN_DBL_MAX;
    3692             }
    3693             double value = columnActivity_[i] * multiplier;
    3694             if (fabs(value)>1.0e20) {
    3695               //printf("bad value of %g for column %d\n",value,i);
    3696               setColumnStatus(i,superBasic);
    3697               if (columnUpperWork_[i]<0.0) {
    3698                 value=columnUpperWork_[i];
    3699               } else if (columnLowerWork_[i]>0.0) {
    3700                 value=columnLowerWork_[i];
    3701               } else {
    3702                 value=0.0;
     3456                columnUpperWork_[i] = 0.0;
     3457                columnLowerWork_[i] = 0.0;
    37033458              }
    37043459            }
    3705             columnActivityWork_[i] = value;
    37063460          }
    3707           for (i=0;i<numberRows_;i++) {
    3708             dual_[i] /= rowScale_[i];
    3709             dual_[i] *= objectiveScale_;
    3710             rowReducedCost_[i] = dual_[i];
    3711             double multiplier = rhsScale_*rowScale_[i];
    3712             double value = rowActivity_[i] * multiplier;
    3713             if (fabs(value)>1.0e20) {
    3714               //printf("bad value of %g for row %d\n",value,i);
    3715               setRowStatus(i,superBasic);
    3716               if (rowUpperWork_[i]<0.0) {
    3717                 value=rowUpperWork_[i];
    3718               } else if (rowLowerWork_[i]>0.0) {
    3719                 value=rowLowerWork_[i];
    3720               } else {
    3721                 value=0.0;
    3722               }
    3723             }
    3724             rowActivityWork_[i] = value;
     3461        } else if (upperValue<1.0e20) {
     3462          columnLowerWork_[i]=-COIN_DBL_MAX;
     3463          columnUpperWork_[i]=upperValue*multiplier;
     3464        } else {
     3465          // free
     3466          columnLowerWork_[i]=-COIN_DBL_MAX;
     3467          columnUpperWork_[i]=COIN_DBL_MAX;
     3468        }
     3469        double value = columnActivity_[i] * multiplier;
     3470        if (fabs(value)>1.0e20) {
     3471          //printf("bad value of %g for column %d\n",value,i);
     3472          setColumnStatus(i,superBasic);
     3473          if (columnUpperWork_[i]<0.0) {
     3474            value=columnUpperWork_[i];
     3475          } else if (columnLowerWork_[i]>0.0) {
     3476            value=columnLowerWork_[i];
     3477          } else {
     3478            value=0.0;
    37253479          }
    3726         } else {
    3727           assert (inverseColumnScale_);
    3728           const double * inverseScale = inverseColumnScale_;
    3729           for (i=0;i<numberColumns_;i++) {
    3730             CoinAssert(fabs(obj[i])<1.0e25);
    3731             double scaleFactor = columnScale_[i];
    3732             double multiplier = rhsScale_*inverseScale[i];
    3733             scaleFactor *= direction;
    3734             objectiveWork_[i] = obj[i]*scaleFactor;
    3735             reducedCostWork_[i] = reducedCost_[i]*scaleFactor;
    3736             double lowerValue = columnLower_[i];
    3737             double upperValue = columnUpper_[i];
    3738             if (lowerValue>-1.0e20) {
    3739               columnLowerWork_[i]=lowerValue*multiplier;
    3740               if (upperValue>=1.0e20) {
    3741                 columnUpperWork_[i]=COIN_DBL_MAX;
    3742               } else {
    3743                 columnUpperWork_[i]=upperValue*multiplier;
    3744                 if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    3745                   if (columnLowerWork_[i]>=0.0) {
    3746                     columnUpperWork_[i] = columnLowerWork_[i];
    3747                   } else if (columnUpperWork_[i]<=0.0) {
    3748                     columnLowerWork_[i] = columnUpperWork_[i];
    3749                   } else {
    3750                     columnUpperWork_[i] = 0.0;
    3751                     columnLowerWork_[i] = 0.0;
    3752                   }
    3753                 }
    3754               }
    3755             } else if (upperValue<1.0e20) {
    3756               columnLowerWork_[i]=-COIN_DBL_MAX;
    3757               columnUpperWork_[i]=upperValue*multiplier;
    3758             } else {
    3759               // free
    3760               columnLowerWork_[i]=-COIN_DBL_MAX;
    3761               columnUpperWork_[i]=COIN_DBL_MAX;
    3762             }
    3763             double value = columnActivity_[i] * multiplier;
    3764             if (fabs(value)>1.0e20) {
    3765               //printf("bad value of %g for column %d\n",value,i);
    3766               setColumnStatus(i,superBasic);
    3767               if (columnUpperWork_[i]<0.0) {
    3768                 value=columnUpperWork_[i];
    3769               } else if (columnLowerWork_[i]>0.0) {
    3770                 value=columnLowerWork_[i];
    3771               } else {
    3772                 value=0.0;
    3773               }
    3774             }
    3775             columnActivityWork_[i] = value;
     3480        }
     3481        columnActivityWork_[i] = value;
     3482      }
     3483      inverseScale = inverseRowScale_;
     3484      for (i=0;i<numberRows_;i++) {
     3485        dual_[i] *= inverseScale[i];
     3486        dual_[i] *= objectiveScale_;
     3487        rowReducedCost_[i] = dual_[i];
     3488        double multiplier = rhsScale_*rowScale_[i];
     3489        double value = rowActivity_[i] * multiplier;
     3490        if (fabs(value)>1.0e20) {
     3491          //printf("bad value of %g for row %d\n",value,i);
     3492          setRowStatus(i,superBasic);
     3493          if (rowUpperWork_[i]<0.0) {
     3494            value=rowUpperWork_[i];
     3495          } else if (rowLowerWork_[i]>0.0) {
     3496            value=rowLowerWork_[i];
     3497          } else {
     3498            value=0.0;
    37763499          }
    3777           inverseScale = inverseRowScale_;
    3778           for (i=0;i<numberRows_;i++) {
    3779             dual_[i] *= inverseScale[i];
    3780             dual_[i] *= objectiveScale_;
    3781             rowReducedCost_[i] = dual_[i];
    3782             double multiplier = rhsScale_*rowScale_[i];
    3783             double value = rowActivity_[i] * multiplier;
    3784             if (fabs(value)>1.0e20) {
    3785               //printf("bad value of %g for row %d\n",value,i);
    3786               setRowStatus(i,superBasic);
    3787               if (rowUpperWork_[i]<0.0) {
    3788                 value=rowUpperWork_[i];
    3789               } else if (rowLowerWork_[i]>0.0) {
    3790                 value=rowLowerWork_[i];
    3791               } else {
    3792                 value=0.0;
    3793               }
    3794             }
    3795             rowActivityWork_[i] = value;
     3500        }
     3501        rowActivityWork_[i] = value;
     3502      }
     3503    } else if (objectiveScale_!=1.0||rhsScale_!=1.0) {
     3504      // on entry
     3505      for (i=0;i<numberColumns_;i++) {
     3506        double value = columnActivity_[i];
     3507        value *= rhsScale_;
     3508        if (fabs(value)>1.0e20) {
     3509          //printf("bad value of %g for column %d\n",value,i);
     3510          setColumnStatus(i,superBasic);
     3511          if (columnUpperWork_[i]<0.0) {
     3512            value=columnUpperWork_[i];
     3513          } else if (columnLowerWork_[i]>0.0) {
     3514            value=columnLowerWork_[i];
     3515          } else {
     3516            value=0.0;
    37963517          }
    37973518        }
    3798       } else if (objectiveScale_!=1.0||rhsScale_!=1.0) {
    3799         // on entry
    3800         for (i=0;i<numberColumns_;i++) {
    3801           double value = columnActivity_[i];
    3802           value *= rhsScale_;
    3803           if (fabs(value)>1.0e20) {
    3804             //printf("bad value of %g for column %d\n",value,i);
    3805             setColumnStatus(i,superBasic);
    3806             if (columnUpperWork_[i]<0.0) {
    3807               value=columnUpperWork_[i];
    3808             } else if (columnLowerWork_[i]>0.0) {
    3809               value=columnLowerWork_[i];
    3810             } else {
    3811               value=0.0;
    3812             }
    3813           }
    3814           columnActivityWork_[i] = value;
    3815           reducedCostWork_[i] = reducedCost_[i]*objectiveScale_;
    3816         }
    3817         for (i=0;i<numberRows_;i++) {
    3818           double value = rowActivity_[i];
    3819           value *= rhsScale_;
    3820           if (fabs(value)>1.0e20) {
    3821             //printf("bad value of %g for row %d\n",value,i);
    3822             setRowStatus(i,superBasic);
    3823             if (rowUpperWork_[i]<0.0) {
    3824               value=rowUpperWork_[i];
    3825             } else if (rowLowerWork_[i]>0.0) {
    3826               value=rowLowerWork_[i];
    3827             } else {
    3828               value=0.0;
    3829             }
    3830           }
    3831           rowActivityWork_[i] = value;
    3832           dual_[i] *= objectiveScale_;
    3833           rowReducedCost_[i] = dual_[i];
    3834         }
    3835       } else {
    3836         // on entry
    3837         for (i=0;i<numberColumns_;i++) {
    3838           double value = columnActivity_[i];
    3839           if (fabs(value)>1.0e20) {
    3840             //printf("bad value of %g for column %d\n",value,i);
    3841             setColumnStatus(i,superBasic);
    3842             if (columnUpperWork_[i]<0.0) {
    3843               value=columnUpperWork_[i];
    3844             } else if (columnLowerWork_[i]>0.0) {
    3845               value=columnLowerWork_[i];
    3846             } else {
    3847               value=0.0;
    3848             }
    3849           }
    3850           columnActivityWork_[i] = value;
    3851           reducedCostWork_[i] = reducedCost_[i];
    3852         }
    3853         for (i=0;i<numberRows_;i++) {
    3854           double value = rowActivity_[i];
    3855           if (fabs(value)>1.0e20) {
    3856             //printf("bad value of %g for row %d\n",value,i);
    3857             setRowStatus(i,superBasic);
    3858             if (rowUpperWork_[i]<0.0) {
    3859               value=rowUpperWork_[i];
    3860             } else if (rowLowerWork_[i]>0.0) {
    3861               value=rowLowerWork_[i];
    3862             } else {
    3863               value=0.0;
    3864             }
    3865           }
    3866           rowActivityWork_[i] = value;
    3867           rowReducedCost_[i] = dual_[i];
    3868         }
    3869       }
    3870 #ifdef CLP_AUXILIARY_MODEL
    3871     }
    3872 #endif
     3519        columnActivityWork_[i] = value;
     3520        reducedCostWork_[i] = reducedCost_[i]*objectiveScale_;
     3521      }
     3522      for (i=0;i<numberRows_;i++) {
     3523        double value = rowActivity_[i];
     3524        value *= rhsScale_;
     3525        if (fabs(value)>1.0e20) {
     3526          //printf("bad value of %g for row %d\n",value,i);
     3527          setRowStatus(i,superBasic);
     3528          if (rowUpperWork_[i]<0.0) {
     3529            value=rowUpperWork_[i];
     3530          } else if (rowLowerWork_[i]>0.0) {
     3531            value=rowLowerWork_[i];
     3532          } else {
     3533            value=0.0;
     3534          }
     3535        }
     3536        rowActivityWork_[i] = value;
     3537        dual_[i] *= objectiveScale_;
     3538        rowReducedCost_[i] = dual_[i];
     3539      }
     3540    } else {
     3541      // on entry
     3542      for (i=0;i<numberColumns_;i++) {
     3543        double value = columnActivity_[i];
     3544        if (fabs(value)>1.0e20) {
     3545          //printf("bad value of %g for column %d\n",value,i);
     3546          setColumnStatus(i,superBasic);
     3547          if (columnUpperWork_[i]<0.0) {
     3548            value=columnUpperWork_[i];
     3549          } else if (columnLowerWork_[i]>0.0) {
     3550            value=columnLowerWork_[i];
     3551          } else {
     3552            value=0.0;
     3553          }
     3554        }
     3555        columnActivityWork_[i] = value;
     3556        reducedCostWork_[i] = reducedCost_[i];
     3557      }
     3558      for (i=0;i<numberRows_;i++) {
     3559        double value = rowActivity_[i];
     3560        if (fabs(value)>1.0e20) {
     3561          //printf("bad value of %g for row %d\n",value,i);
     3562          setRowStatus(i,superBasic);
     3563          if (rowUpperWork_[i]<0.0) {
     3564            value=rowUpperWork_[i];
     3565          } else if (rowLowerWork_[i]>0.0) {
     3566            value=rowLowerWork_[i];
     3567          } else {
     3568            value=0.0;
     3569          }
     3570        }
     3571        rowActivityWork_[i] = value;
     3572        rowReducedCost_[i] = dual_[i];
     3573      }
     3574    }
    38733575  }
    38743576 
    3875   if (what==63&&doSanityCheck&&!auxiliaryModel_) {
     3577  if (what==63&&doSanityCheck) {
    38763578    // check rim of problem okay
    38773579    if (!sanityCheck())
     
    39813683    } else {
    39823684      int iRow,iColumn;
    3983 #ifdef CLP_AUXILIARY_MODEL
    3984       if (auxiliaryModel_) {
    3985         for (iRow=0;iRow<4;iRow++) {
    3986           assert (!rowArray_[iRow]);
    3987           rowArray_[iRow]=auxiliaryModel_->rowArray_[iRow];
    3988           // For safety
    3989           memset(rowArray_[iRow]->denseVector(),0,rowArray_[iRow]->capacity()*sizeof(double));
    3990         }
    3991         for (iColumn=0;iColumn<2;iColumn++) {
    3992           assert (!columnArray_[iColumn]);
    3993           columnArray_[iColumn]=auxiliaryModel_->columnArray_[iColumn];
    3994           // For safety
    3995           memset(columnArray_[iColumn]->denseVector(),0,columnArray_[iColumn]->capacity()*sizeof(double));
    3996         }
    3997         // do matrices
    3998         rowCopy_=auxiliaryModel_->rowCopy_;
    3999         ClpMatrixBase * temp = matrix_;
    4000         matrix_=auxiliaryModel_->matrix_;
    4001         auxiliaryModel_->matrix_=temp;
    4002       }
    4003 #endif
    40043685      for (iRow=0;iRow<4;iRow++) {
    40053686        rowArray_[iRow]->clear();
     
    40313712  if ((what&5)!=0)
    40323713    matrix_->generalExpanded(this,9,what); // update costs and bounds if necessary
    4033   if (goodMatrix&&(specialOptions_&65536)!=0&&!auxiliaryModel_) {
     3714  if (goodMatrix&&(specialOptions_&65536)!=0) {
    40343715    int save = maximumColumns_+maximumRows_;
    40353716    CoinMemcpyN(cost_,numberTotal,cost_+save);
     
    40713752  int numberRows2 = numberRows_+numberExtraRows_;
    40723753  int numberTotal = numberRows2+numberColumns_;
    4073   //assert (!auxiliaryModel_);
    4074 #ifdef CLP_AUXILIARY_MODEL
    4075   if (!auxiliaryModel_) {
    4076 #endif
    4077     if ((specialOptions_&65536)!=0&&true) {
    4078       assert (!initial);
    4079       int save = maximumColumns_+maximumRows_;
    4080       CoinMemcpyN(lower_+save,numberTotal,lower_);
    4081       CoinMemcpyN(upper_+save,numberTotal,upper_);
    4082       return;
    4083     }
    4084     const double * rowScale = rowScale_;
    4085     const double * columnScale = columnScale_;
    4086     // clean up any mismatches on infinity
    4087     // and fix any variables with tiny gaps
    4088     double primalTolerance=dblParam_[ClpPrimalTolerance];
    4089     if(rowScale) {
    4090       // If scaled then do all columns later in one loop
    4091       if (!initial) {
    4092         if (!inverseColumnScale_) {
    4093           for (i=0;i<numberColumns_;i++) {
    4094             double multiplier = rhsScale_/columnScale[i];
    4095             double lowerValue = columnLower_[i];
    4096             double upperValue = columnUpper_[i];
    4097             if (lowerValue>-1.0e20) {
    4098               columnLowerWork_[i]=lowerValue*multiplier;
    4099               if (upperValue>=1.0e20) {
    4100                 columnUpperWork_[i]=COIN_DBL_MAX;
    4101               } else {
    4102                 columnUpperWork_[i]=upperValue*multiplier;
    4103                 if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    4104                   if (columnLowerWork_[i]>=0.0) {
    4105                     columnUpperWork_[i] = columnLowerWork_[i];
    4106                   } else if (columnUpperWork_[i]<=0.0) {
    4107                     columnLowerWork_[i] = columnUpperWork_[i];
    4108                   } else {
    4109                     columnUpperWork_[i] = 0.0;
    4110                     columnLowerWork_[i] = 0.0;
    4111                   }
    4112                 }
    4113               }
    4114             } else if (upperValue<1.0e20) {
    4115               columnLowerWork_[i]=-COIN_DBL_MAX;
    4116               columnUpperWork_[i]=upperValue*multiplier;
    4117             } else {
    4118               // free
    4119               columnLowerWork_[i]=-COIN_DBL_MAX;
    4120               columnUpperWork_[i]=COIN_DBL_MAX;
    4121             }
    4122           }
    4123         } else {
    4124           assert (inverseColumnScale_);
    4125           const double * inverseScale = inverseColumnScale_;
    4126           for (i=0;i<numberColumns_;i++) {
    4127             double multiplier = rhsScale_*inverseScale[i];
    4128             double lowerValue = columnLower_[i];
    4129             double upperValue = columnUpper_[i];
    4130             if (lowerValue>-1.0e20) {
    4131               columnLowerWork_[i]=lowerValue*multiplier;
    4132               if (upperValue>=1.0e20) {
    4133                 columnUpperWork_[i]=COIN_DBL_MAX;
    4134               } else {
    4135                 columnUpperWork_[i]=upperValue*multiplier;
    4136                 if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    4137                   if (columnLowerWork_[i]>=0.0) {
    4138                     columnUpperWork_[i] = columnLowerWork_[i];
    4139                   } else if (columnUpperWork_[i]<=0.0) {
    4140                     columnLowerWork_[i] = columnUpperWork_[i];
    4141                   } else {
    4142                     columnUpperWork_[i] = 0.0;
    4143                     columnLowerWork_[i] = 0.0;
    4144                   }
    4145                 }
    4146               }
    4147             } else if (upperValue<1.0e20) {
    4148               columnLowerWork_[i]=-COIN_DBL_MAX;
    4149               columnUpperWork_[i]=upperValue*multiplier;
    4150             } else {
    4151               // free
    4152               columnLowerWork_[i]=-COIN_DBL_MAX;
    4153               columnUpperWork_[i]=COIN_DBL_MAX;
    4154             }
    4155           }
    4156         }
    4157       }
    4158       for (i=0;i<numberRows_;i++) {
    4159         double multiplier = rhsScale_*rowScale[i];
    4160         double lowerValue = rowLower_[i];
    4161         double upperValue = rowUpper_[i];
    4162         if (lowerValue>-1.0e20) {
    4163           rowLowerWork_[i]=lowerValue*multiplier;
    4164           if (upperValue>=1.0e20) {
    4165             rowUpperWork_[i]=COIN_DBL_MAX;
    4166           } else {
    4167             rowUpperWork_[i]=upperValue*multiplier;
    4168             if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
    4169               if (rowLowerWork_[i]>=0.0) {
    4170                 rowUpperWork_[i] = rowLowerWork_[i];
    4171               } else if (rowUpperWork_[i]<=0.0) {
    4172                 rowLowerWork_[i] = rowUpperWork_[i];
    4173               } else {
    4174                 rowUpperWork_[i] = 0.0;
    4175                 rowLowerWork_[i] = 0.0;
    4176               }
    4177             }
    4178           }
    4179         } else if (upperValue<1.0e20) {
    4180           rowLowerWork_[i]=-COIN_DBL_MAX;
    4181           rowUpperWork_[i]=upperValue*multiplier;
    4182         } else {
    4183           // free
    4184           rowLowerWork_[i]=-COIN_DBL_MAX;
    4185           rowUpperWork_[i]=COIN_DBL_MAX;
    4186         }
    4187       }
    4188     } else if (rhsScale_!=1.0) {
     3754  if ((specialOptions_&65536)!=0&&true) {
     3755    assert (!initial);
     3756    int save = maximumColumns_+maximumRows_;
     3757    CoinMemcpyN(lower_+save,numberTotal,lower_);
     3758    CoinMemcpyN(upper_+save,numberTotal,upper_);
     3759    return;
     3760  }
     3761  const double * rowScale = rowScale_;
     3762  // clean up any mismatches on infinity
     3763  // and fix any variables with tiny gaps
     3764  double primalTolerance=dblParam_[ClpPrimalTolerance];
     3765  if(rowScale) {
     3766    // If scaled then do all columns later in one loop
     3767    if (!initial) {
     3768      const double * inverseScale = inverseColumnScale_;
    41893769      for (i=0;i<numberColumns_;i++) {
     3770        double multiplier = rhsScale_*inverseScale[i];
    41903771        double lowerValue = columnLower_[i];
    41913772        double upperValue = columnUpper_[i];
    41923773        if (lowerValue>-1.0e20) {
    4193           columnLowerWork_[i]=lowerValue*rhsScale_;
     3774          columnLowerWork_[i]=lowerValue*multiplier;
    41943775          if (upperValue>=1.0e20) {
    41953776            columnUpperWork_[i]=COIN_DBL_MAX;
    41963777          } else {
    4197             columnUpperWork_[i]=upperValue*rhsScale_;
     3778            columnUpperWork_[i]=upperValue*multiplier;
    41983779            if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    41993780              if (columnLowerWork_[i]>=0.0) {
     
    42093790        } else if (upperValue<1.0e20) {
    42103791          columnLowerWork_[i]=-COIN_DBL_MAX;
    4211           columnUpperWork_[i]=upperValue*rhsScale_;
     3792          columnUpperWork_[i]=upperValue*multiplier;
    42123793        } else {
    42133794          // free
     
    42163797        }
    42173798      }
    4218       for (i=0;i<numberRows_;i++) {
    4219         double lowerValue = rowLower_[i];
    4220         double upperValue = rowUpper_[i];
    4221         if (lowerValue>-1.0e20) {
    4222           rowLowerWork_[i]=lowerValue*rhsScale_;
    4223           if (upperValue>=1.0e20) {
    4224             rowUpperWork_[i]=COIN_DBL_MAX;
    4225           } else {
    4226             rowUpperWork_[i]=upperValue*rhsScale_;
    4227             if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
    4228               if (rowLowerWork_[i]>=0.0) {
    4229                 rowUpperWork_[i] = rowLowerWork_[i];
    4230               } else if (rowUpperWork_[i]<=0.0) {
    4231                 rowLowerWork_[i] = rowUpperWork_[i];
    4232               } else {
    4233                 rowUpperWork_[i] = 0.0;
    4234                 rowLowerWork_[i] = 0.0;
    4235               }
     3799    }
     3800    for (i=0;i<numberRows_;i++) {
     3801      double multiplier = rhsScale_*rowScale[i];
     3802      double lowerValue = rowLower_[i];
     3803      double upperValue = rowUpper_[i];
     3804      if (lowerValue>-1.0e20) {
     3805        rowLowerWork_[i]=lowerValue*multiplier;
     3806        if (upperValue>=1.0e20) {
     3807          rowUpperWork_[i]=COIN_DBL_MAX;
     3808        } else {
     3809          rowUpperWork_[i]=upperValue*multiplier;
     3810          if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
     3811            if (rowLowerWork_[i]>=0.0) {
     3812              rowUpperWork_[i] = rowLowerWork_[i];
     3813            } else if (rowUpperWork_[i]<=0.0) {
     3814              rowLowerWork_[i] = rowUpperWork_[i];
     3815            } else {
     3816              rowUpperWork_[i] = 0.0;
     3817              rowLowerWork_[i] = 0.0;
    42363818            }
    42373819          }
    4238         } else if (upperValue<1.0e20) {
    4239           rowLowerWork_[i]=-COIN_DBL_MAX;
    4240           rowUpperWork_[i]=upperValue*rhsScale_;
     3820        }
     3821      } else if (upperValue<1.0e20) {
     3822        rowLowerWork_[i]=-COIN_DBL_MAX;
     3823        rowUpperWork_[i]=upperValue*multiplier;
     3824      } else {
     3825        // free
     3826        rowLowerWork_[i]=-COIN_DBL_MAX;
     3827        rowUpperWork_[i]=COIN_DBL_MAX;
     3828      }
     3829    }
     3830  } else if (rhsScale_!=1.0) {
     3831    for (i=0;i<numberColumns_;i++) {
     3832      double lowerValue = columnLower_[i];
     3833      double upperValue = columnUpper_[i];
     3834      if (lowerValue>-1.0e20) {
     3835        columnLowerWork_[i]=lowerValue*rhsScale_;
     3836        if (upperValue>=1.0e20) {
     3837          columnUpperWork_[i]=COIN_DBL_MAX;
    42413838        } else {
    4242           // free
    4243           rowLowerWork_[i]=-COIN_DBL_MAX;
    4244           rowUpperWork_[i]=COIN_DBL_MAX;
    4245         }
    4246       }
    4247     } else {
    4248       for (i=0;i<numberColumns_;i++) {
    4249         double lowerValue = columnLower_[i];
    4250         double upperValue = columnUpper_[i];
    4251         if (lowerValue>-1.0e20) {
    4252           columnLowerWork_[i]=lowerValue;
    4253           if (upperValue>=1.0e20) {
    4254             columnUpperWork_[i]=COIN_DBL_MAX;
    4255           } else {
    4256             columnUpperWork_[i]=upperValue;
    4257             if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    4258               if (columnLowerWork_[i]>=0.0) {
    4259                 columnUpperWork_[i] = columnLowerWork_[i];
    4260               } else if (columnUpperWork_[i]<=0.0) {
    4261                 columnLowerWork_[i] = columnUpperWork_[i];
    4262               } else {
    4263                 columnUpperWork_[i] = 0.0;
    4264                 columnLowerWork_[i] = 0.0;
    4265               }
     3839          columnUpperWork_[i]=upperValue*rhsScale_;
     3840          if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3841            if (columnLowerWork_[i]>=0.0) {
     3842              columnUpperWork_[i] = columnLowerWork_[i];
     3843            } else if (columnUpperWork_[i]<=0.0) {
     3844              columnLowerWork_[i] = columnUpperWork_[i];
     3845            } else {
     3846              columnUpperWork_[i] = 0.0;
     3847              columnLowerWork_[i] = 0.0;
    42663848            }
    42673849          }
    4268         } else if (upperValue<1.0e20) {
    4269           columnLowerWork_[i]=-COIN_DBL_MAX;
    4270           columnUpperWork_[i]=upperValue;
     3850        }
     3851      } else if (upperValue<1.0e20) {
     3852        columnLowerWork_[i]=-COIN_DBL_MAX;
     3853        columnUpperWork_[i]=upperValue*rhsScale_;
     3854      } else {
     3855        // free
     3856        columnLowerWork_[i]=-COIN_DBL_MAX;
     3857        columnUpperWork_[i]=COIN_DBL_MAX;
     3858      }
     3859    }
     3860    for (i=0;i<numberRows_;i++) {
     3861      double lowerValue = rowLower_[i];
     3862      double upperValue = rowUpper_[i];
     3863      if (lowerValue>-1.0e20) {
     3864        rowLowerWork_[i]=lowerValue*rhsScale_;
     3865        if (upperValue>=1.0e20) {
     3866          rowUpperWork_[i]=COIN_DBL_MAX;
    42713867        } else {
    4272           // free
    4273           columnLowerWork_[i]=-COIN_DBL_MAX;
    4274           columnUpperWork_[i]=COIN_DBL_MAX;
    4275         }
    4276       }
    4277       for (i=0;i<numberRows_;i++) {
    4278         double lowerValue = rowLower_[i];
    4279         double upperValue = rowUpper_[i];
    4280         if (lowerValue>-1.0e20) {
    4281           rowLowerWork_[i]=lowerValue;
    4282           if (upperValue>=1.0e20) {
    4283             rowUpperWork_[i]=COIN_DBL_MAX;
    4284           } else {
    4285             rowUpperWork_[i]=upperValue;
    4286             if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
    4287               if (rowLowerWork_[i]>=0.0) {
    4288                 rowUpperWork_[i] = rowLowerWork_[i];
    4289               } else if (rowUpperWork_[i]<=0.0) {
    4290                 rowLowerWork_[i] = rowUpperWork_[i];
    4291               } else {
    4292                 rowUpperWork_[i] = 0.0;
    4293                 rowLowerWork_[i] = 0.0;
    4294               }
     3868          rowUpperWork_[i]=upperValue*rhsScale_;
     3869          if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
     3870            if (rowLowerWork_[i]>=0.0) {
     3871              rowUpperWork_[i] = rowLowerWork_[i];
     3872            } else if (rowUpperWork_[i]<=0.0) {
     3873              rowLowerWork_[i] = rowUpperWork_[i];
     3874            } else {
     3875              rowUpperWork_[i] = 0.0;
     3876              rowLowerWork_[i] = 0.0;
    42953877            }
    42963878          }
    4297         } else if (upperValue<1.0e20) {
    4298           rowLowerWork_[i]=-COIN_DBL_MAX;
    4299           rowUpperWork_[i]=upperValue;
     3879        }
     3880      } else if (upperValue<1.0e20) {
     3881        rowLowerWork_[i]=-COIN_DBL_MAX;
     3882        rowUpperWork_[i]=upperValue*rhsScale_;
     3883      } else {
     3884        // free
     3885        rowLowerWork_[i]=-COIN_DBL_MAX;
     3886        rowUpperWork_[i]=COIN_DBL_MAX;
     3887      }
     3888    }
     3889  } else {
     3890    for (i=0;i<numberColumns_;i++) {
     3891      double lowerValue = columnLower_[i];
     3892      double upperValue = columnUpper_[i];
     3893      if (lowerValue>-1.0e20) {
     3894        columnLowerWork_[i]=lowerValue;
     3895        if (upperValue>=1.0e20) {
     3896          columnUpperWork_[i]=COIN_DBL_MAX;
    43003897        } else {
    4301           // free
    4302           rowLowerWork_[i]=-COIN_DBL_MAX;
    4303           rowUpperWork_[i]=COIN_DBL_MAX;
    4304         }
    4305       }
    4306     }
    4307 #ifdef CLP_AUXILIARY_MODEL
    4308   } else {
    4309     // auxiliary model
    4310     if (!initial) {
    4311       // just copy
    4312       CoinMemcpyN(auxiliaryModel_->lower_+numberTotal,numberTotal,lower_);
    4313       CoinMemcpyN(auxiliaryModel_->upper_+numberTotal,numberTotal,upper_);
    4314     } else {
    4315       // clean up any mismatches on infinity
    4316       // and fix any variables with tiny gaps
    4317       const double * rowScale = auxiliaryModel_->rowScale_ ;
    4318       const double * columnScale = auxiliaryModel_->columnScale_ ;
    4319       assert (rhsScale_==1.0);
    4320       assert (objectiveScale_==1.0);
    4321       if ((auxiliaryModel_->specialOptions_&2)==0) {
    4322         // need to do column bounds
    4323         if (columnScale) {
    4324           const double * inverseScale = columnScale+numberColumns_;
    4325           // scaled
    4326           for (i=0;i<numberColumns_;i++) {
    4327             double value;
    4328             value = columnLower_[i];
    4329             if (value>-1.0e20)
    4330               value *= inverseScale[i];
    4331             lower_[i] = value;
    4332             value = columnUpper_[i];
    4333             if (value<1.0e20)
    4334               value *= inverseScale[i];
    4335             upper_[i] = value;
    4336           }
    4337         } else {
    4338           for (i=0;i<numberColumns_;i++) {
    4339             lower_[i] = columnLower_[i];
    4340             upper_[i] = columnUpper_[i];
    4341           }
    4342         }
    4343         // save
    4344  CoinMemcpyN(lower_,numberColumns_,auxiliaryModel_->lower_+numberTotal);
    4345  CoinMemcpyN(upper_,numberColumns_,auxiliaryModel_->upper_+numberTotal);
    4346       } else {
    4347         // just copy
    4348  CoinMemcpyN(auxiliaryModel_->lower_+numberTotal,numberColumns_,lower_);
    4349  CoinMemcpyN(auxiliaryModel_->upper_+numberTotal,numberColumns_,upper_);
    4350       }
    4351       if ((auxiliaryModel_->specialOptions_&1)==0) {
    4352         // need to do row bounds
    4353         if (rowScale) {
    4354           // scaled
    4355           for (i=0;i<numberRows_;i++) {
    4356             double value;
    4357             value = rowLower_[i];
    4358             if (value>-1.0e20)
    4359               value *= rowScale[i];
    4360             lower_[i+numberColumns_] = value;
    4361             value = rowUpper_[i];
    4362             if (value<1.0e20)
    4363               value *= rowScale[i];
    4364             upper_[i+numberColumns_] = value;
    4365           }
    4366         } else {
    4367           for (i=0;i<numberRows_;i++) {
    4368             lower_[i+numberColumns_] = rowLower_[i];
    4369             upper_[i+numberColumns_] = rowUpper_[i];
    4370           }
    4371         }
    4372         // save
    4373  CoinMemcpyN(          lower_+numberColumns_,numberRows_,auxiliaryModel_->lower_+numberTotal+numberColumns_);
    4374  CoinMemcpyN(          upper_+numberColumns_,numberRows_,auxiliaryModel_->upper_+numberTotal+numberColumns_);
    4375       } else {
    4376         // just copy
    4377  CoinMemcpyN(          auxiliaryModel_->lower_+numberTotal+numberColumns_,numberRows_,
    4378               lower_+numberColumns_);
    4379  CoinMemcpyN(          auxiliaryModel_->upper_+numberTotal+numberColumns_,numberRows_,
    4380               upper_+numberColumns_);
    4381       }
    4382     }
    4383     if (initial) {
    4384       // do basis
    4385       assert ((auxiliaryModel_->specialOptions_&8)!=0);
    4386       // clean solution trusting basis
    4387       for ( i=0;i<numberTotal;i++) {
    4388         Status status =getStatus(i);
    4389         double value=solution_[i]; // may as well leave if basic (values pass)
    4390         if (status!=basic) {
    4391           if (upper_[i]==lower_[i]) {
    4392             setStatus(i,isFixed);
    4393             value=lower_[i];
    4394           } else if (status==atLowerBound) {
    4395             if (lower_[i]>-1.0e20) {
    4396               value=lower_[i];
     3898          columnUpperWork_[i]=upperValue;
     3899          if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3900            if (columnLowerWork_[i]>=0.0) {
     3901              columnUpperWork_[i] = columnLowerWork_[i];
     3902            } else if (columnUpperWork_[i]<=0.0) {
     3903              columnLowerWork_[i] = columnUpperWork_[i];
    43973904            } else {
    4398               if (upper_[i]<1.0e20) {
    4399                 value=upper_[i];
    4400                 setStatus(i,atUpperBound);
    4401               } else {
    4402                 setStatus(i,isFree);
    4403               }
    4404             }
    4405           } else if (status==atUpperBound) {
    4406             if (upper_[i]<1.0e20) {
    4407               value=upper_[i];
    4408             } else {
    4409               if (lower_[i]>-1.0e20) {
    4410                 value=lower_[i];
    4411                 setStatus(i,atLowerBound);
    4412               } else {
    4413                 setStatus(i,isFree);
    4414               }
     3905              columnUpperWork_[i] = 0.0;
     3906              columnLowerWork_[i] = 0.0;
    44153907            }
    44163908          }
    44173909        }
    4418         solution_[i]=value;
    4419       }
    4420     }
    4421   }
    4422 #endif
     3910      } else if (upperValue<1.0e20) {
     3911        columnLowerWork_[i]=-COIN_DBL_MAX;
     3912        columnUpperWork_[i]=upperValue;
     3913      } else {
     3914        // free
     3915        columnLowerWork_[i]=-COIN_DBL_MAX;
     3916        columnUpperWork_[i]=COIN_DBL_MAX;
     3917      }
     3918    }
     3919    for (i=0;i<numberRows_;i++) {
     3920      double lowerValue = rowLower_[i];
     3921      double upperValue = rowUpper_[i];
     3922      if (lowerValue>-1.0e20) {
     3923        rowLowerWork_[i]=lowerValue;
     3924        if (upperValue>=1.0e20) {
     3925          rowUpperWork_[i]=COIN_DBL_MAX;
     3926        } else {
     3927          rowUpperWork_[i]=upperValue;
     3928          if (fabs(rowUpperWork_[i]-rowLowerWork_[i])<=primalTolerance) {
     3929            if (rowLowerWork_[i]>=0.0) {
     3930              rowUpperWork_[i] = rowLowerWork_[i];
     3931            } else if (rowUpperWork_[i]<=0.0) {
     3932              rowLowerWork_[i] = rowUpperWork_[i];
     3933            } else {
     3934              rowUpperWork_[i] = 0.0;
     3935              rowLowerWork_[i] = 0.0;
     3936            }
     3937          }
     3938        }
     3939      } else if (upperValue<1.0e20) {
     3940        rowLowerWork_[i]=-COIN_DBL_MAX;
     3941        rowUpperWork_[i]=upperValue;
     3942      } else {
     3943        // free
     3944        rowLowerWork_[i]=-COIN_DBL_MAX;
     3945        rowUpperWork_[i]=COIN_DBL_MAX;
     3946      }
     3947    }
     3948  }
    44233949#ifndef NDEBUG
    44243950  if ((specialOptions_&65536)!=0&&false) {
     
    44393965  int numberRows2 = numberRows_+numberExtraRows_;
    44403966  int numberTotal = numberRows2+numberColumns_;
    4441   //assert (!auxiliaryModel_);
    4442 #ifdef CLP_AUXILIARY_MODEL
    4443   if (!auxiliaryModel_||(initial&&(auxiliaryModel_->specialOptions_&4)==0)) {
    4444 #endif
    4445     if ((specialOptions_&65536)!=0&&true) {
    4446       assert (!initial);
    4447       int save = maximumColumns_+maximumRows_;
    4448       CoinMemcpyN(cost_+save,numberTotal,cost_);
    4449       return;
    4450     }
    4451     double direction = optimizationDirection_*objectiveScale_;
    4452     const double * obj = objective();
    4453 #ifdef CLP_AUXILIARY_MODEL
    4454     const double * rowScale = auxiliaryModel_ ? auxiliaryModel_->rowScale_ : rowScale_;
    4455     const double * columnScale = auxiliaryModel_ ? auxiliaryModel_->columnScale_ : columnScale_;
    4456 #else
    4457     const double * rowScale = rowScale_;
    4458     const double * columnScale = columnScale_;
    4459 #endif
    4460     // and also scale by scale factors
    4461     if (rowScale) {
    4462       if (rowObjective_) {
    4463         for (i=0;i<numberRows_;i++)
    4464           rowObjectiveWork_[i] = rowObjective_[i]*direction/rowScale[i];
    4465       } else {
    4466         memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
    4467       }
    4468       // If scaled then do all columns later in one loop
    4469       if (!initial||auxiliaryModel_) {
    4470         for (i=0;i<numberColumns_;i++) {
    4471           CoinAssert(fabs(obj[i])<1.0e25);
    4472           objectiveWork_[i] = obj[i]*direction*columnScale[i];
    4473         }
    4474       }
     3967  if ((specialOptions_&65536)!=0&&true) {
     3968    assert (!initial);
     3969    int save = maximumColumns_+maximumRows_;
     3970    CoinMemcpyN(cost_+save,numberTotal,cost_);
     3971    return;
     3972  }
     3973  double direction = optimizationDirection_*objectiveScale_;
     3974  const double * obj = objective();
     3975  const double * rowScale = rowScale_;
     3976  const double * columnScale = columnScale_;
     3977  // and also scale by scale factors
     3978  if (rowScale) {
     3979    if (rowObjective_) {
     3980      for (i=0;i<numberRows_;i++)
     3981        rowObjectiveWork_[i] = rowObjective_[i]*direction/rowScale[i];
    44753982    } else {
    4476       if (rowObjective_) {
    4477         for (i=0;i<numberRows_;i++)
    4478           rowObjectiveWork_[i] = rowObjective_[i]*direction;
    4479       } else {
    4480         memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
    4481       }
     3983      memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
     3984    }
     3985    // If scaled then do all columns later in one loop
     3986    if (!initial) {
    44823987      for (i=0;i<numberColumns_;i++) {
    44833988        CoinAssert(fabs(obj[i])<1.0e25);
    4484         objectiveWork_[i] = obj[i]*direction;
    4485       }
    4486     }
    4487 #ifdef CLP_AUXILIARY_MODEL
    4488     if (auxiliaryModel_) {
    4489       // save costs
    4490       CoinMemcpyN(cost_,numberTotal,auxiliaryModel_->cost_+numberTotal);
     3989        objectiveWork_[i] = obj[i]*direction*columnScale[i];
     3990      }
    44913991    }
    44923992  } else {
    4493     // just copy
    4494     CoinMemcpyN(auxiliaryModel_->cost_+numberTotal,numberTotal,cost_);
    4495   }
    4496 #endif
     3993    if (rowObjective_) {
     3994      for (i=0;i<numberRows_;i++)
     3995        rowObjectiveWork_[i] = rowObjective_[i]*direction;
     3996    } else {
     3997      memset(rowObjectiveWork_,0,numberRows_*sizeof(double));
     3998    }
     3999    for (i=0;i<numberColumns_;i++) {
     4000      CoinAssert(fabs(obj[i])<1.0e25);
     4001      objectiveWork_[i] = obj[i]*direction;
     4002    }
     4003  }
    44974004}
    44984005// Does rows and columns and objective
     
    45144021      numberColumns=0;
    45154022  }
    4516 #ifdef CLP_AUXILIARY_MODEL
    4517   if (!auxiliaryModel_) {
    4518 #endif
    4519     int i;
    4520     if (problemStatus_!=1&&problemStatus_!=2) {
    4521       delete [] ray_;
    4522       ray_=NULL;
    4523     }
    4524     // set upperOut_ to furthest away from bound so can use in dual for dualBound_
    4525     upperOut_=1.0;
     4023  int i;
     4024  if (problemStatus_!=1&&problemStatus_!=2) {
     4025    delete [] ray_;
     4026    ray_=NULL;
     4027  }
     4028  // set upperOut_ to furthest away from bound so can use in dual for dualBound_
     4029  upperOut_=1.0;
    45264030#if 0
    4527     {
    4528       int nBad=0;
     4031  {
     4032    int nBad=0;
     4033    for (i=0;i<numberColumns;i++) {
     4034      if (lower_[i]==upper_[i]&&getColumnStatus(i)==basic)
     4035        nBad++;
     4036    }
     4037    if (nBad)
     4038      printf("yy %d basic fixed\n",nBad);
     4039  }
     4040#endif
     4041  // ray may be null if in branch and bound
     4042  if (rowScale_) {
     4043    // Collect infeasibilities
     4044    int numberPrimalScaled=0;
     4045    int numberPrimalUnscaled=0;
     4046    int numberDualScaled=0;
     4047    int numberDualUnscaled=0;
     4048    double scaleC = 1.0/objectiveScale_;
     4049    double scaleR = 1.0/rhsScale_;
     4050    const double * inverseScale = inverseColumnScale_;
     4051    for (i=0;i<numberColumns;i++) {
     4052      double scaleFactor = columnScale_[i];
     4053      double valueScaled = columnActivityWork_[i];
     4054      double lowerScaled = columnLowerWork_[i];
     4055      double upperScaled = columnUpperWork_[i];
     4056      if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     4057        if (valueScaled<lowerScaled-primalTolerance_||
     4058            valueScaled>upperScaled+primalTolerance_)
     4059          numberPrimalScaled++;
     4060        else
     4061          upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     4062      }
     4063      columnActivity_[i] = valueScaled*scaleFactor*scaleR;
     4064      double value = columnActivity_[i];
     4065      if (value<columnLower_[i]-primalTolerance_)
     4066        numberPrimalUnscaled++;
     4067      else if (value>columnUpper_[i]+primalTolerance_)
     4068        numberPrimalUnscaled++;
     4069      double valueScaledDual = reducedCostWork_[i];
     4070      if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     4071        numberDualScaled++;
     4072      if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     4073        numberDualScaled++;
     4074      reducedCost_[i] = (valueScaledDual*scaleC)*inverseScale[i];
     4075      double valueDual = reducedCost_[i];
     4076      if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     4077        numberDualUnscaled++;
     4078      if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     4079        numberDualUnscaled++;
     4080    }
     4081    inverseScale = inverseRowScale_;
     4082    for (i=0;i<numberRows;i++) {
     4083      double scaleFactor = rowScale_[i];
     4084      double valueScaled = rowActivityWork_[i];
     4085      double lowerScaled = rowLowerWork_[i];
     4086      double upperScaled = rowUpperWork_[i];
     4087      if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     4088        if (valueScaled<lowerScaled-primalTolerance_||
     4089            valueScaled>upperScaled+primalTolerance_)
     4090          numberPrimalScaled++;
     4091        else
     4092          upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     4093      }
     4094      rowActivity_[i] = (valueScaled*scaleR)*inverseScale[i];
     4095      double value = rowActivity_[i];
     4096      if (value<rowLower_[i]-primalTolerance_)
     4097        numberPrimalUnscaled++;
     4098      else if (value>rowUpper_[i]+primalTolerance_)
     4099        numberPrimalUnscaled++;
     4100      double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
     4101      if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     4102        numberDualScaled++;
     4103      if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     4104        numberDualScaled++;
     4105      dual_[i] *= scaleFactor*scaleC;
     4106      double valueDual = dual_[i];
     4107      if (rowObjective_)
     4108        valueDual += rowObjective_[i];
     4109      if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     4110        numberDualUnscaled++;
     4111      if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     4112        numberDualUnscaled++;
     4113    }
     4114    if (!problemStatus_&&!secondaryStatus_) {
     4115      // See if we need to set secondary status
     4116      if (numberPrimalUnscaled) {
     4117        if (numberDualUnscaled)
     4118          secondaryStatus_=4;
     4119        else
     4120          secondaryStatus_=2;
     4121      } else {
     4122        if (numberDualUnscaled)
     4123          secondaryStatus_=3;
     4124      }
     4125    }
     4126    if (problemStatus_==2) {
    45294127      for (i=0;i<numberColumns;i++) {
    4530         if (lower_[i]==upper_[i]&&getColumnStatus(i)==basic)
    4531           nBad++;
    4532       }
    4533       if (nBad)
    4534         printf("yy %d basic fixed\n",nBad);
    4535     }
    4536 #endif
    4537     // ray may be null if in branch and bound
    4538     if (rowScale_) {
    4539       // Collect infeasibilities
    4540       int numberPrimalScaled=0;
    4541       int numberPrimalUnscaled=0;
    4542       int numberDualScaled=0;
    4543       int numberDualUnscaled=0;
    4544       double scaleC = 1.0/objectiveScale_;
    4545       double scaleR = 1.0/rhsScale_;
    4546       if (!inverseColumnScale_) {
    4547         for (i=0;i<numberColumns;i++) {
    4548           double scaleFactor = columnScale_[i];
    4549           double valueScaled = columnActivityWork_[i];
    4550           double lowerScaled = columnLowerWork_[i];
    4551           double upperScaled = columnUpperWork_[i];
    4552           if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    4553             if (valueScaled<lowerScaled-primalTolerance_||
    4554                 valueScaled>upperScaled+primalTolerance_)
    4555               numberPrimalScaled++;
    4556             else
    4557               upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    4558           }
    4559           columnActivity_[i] = valueScaled*scaleFactor*scaleR;
    4560           double value = columnActivity_[i];
    4561           if (value<columnLower_[i]-primalTolerance_)
    4562             numberPrimalUnscaled++;
    4563           else if (value>columnUpper_[i]+primalTolerance_)
    4564             numberPrimalUnscaled++;
    4565           double valueScaledDual = reducedCostWork_[i];
    4566           if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    4567             numberDualScaled++;
    4568           if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    4569             numberDualScaled++;
    4570           reducedCost_[i] = (valueScaledDual*scaleC)/scaleFactor;
    4571           double valueDual = reducedCost_[i];
    4572           if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    4573             numberDualUnscaled++;
    4574           if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    4575             numberDualUnscaled++;
    4576         }
    4577         for (i=0;i<numberRows;i++) {
    4578           double scaleFactor = rowScale_[i];
    4579           double valueScaled = rowActivityWork_[i];
    4580           double lowerScaled = rowLowerWork_[i];
    4581           double upperScaled = rowUpperWork_[i];
    4582           if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    4583             if (valueScaled<lowerScaled-primalTolerance_||
    4584                 valueScaled>upperScaled+primalTolerance_)
    4585               numberPrimalScaled++;
    4586             else
    4587               upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    4588           }
    4589           rowActivity_[i] = (valueScaled*scaleR)/scaleFactor;
    4590           double value = rowActivity_[i];
    4591           if (value<rowLower_[i]-primalTolerance_)
    4592             numberPrimalUnscaled++;
    4593           else if (value>rowUpper_[i]+primalTolerance_)
    4594             numberPrimalUnscaled++;
    4595           double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
    4596           if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    4597             numberDualScaled++;
    4598           if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    4599             numberDualScaled++;
    4600           dual_[i] *= scaleFactor*scaleC;
    4601           double valueDual = dual_[i];
    4602           if (rowObjective_)
    4603             valueDual += rowObjective_[i];
    4604           if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    4605             numberDualUnscaled++;
    4606           if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    4607             numberDualUnscaled++;
    4608         }
     4128        ray_[i] *= columnScale_[i];
     4129      }
     4130    } else if (problemStatus_==1&&ray_) {
     4131      for (i=0;i<numberRows;i++) {
     4132        ray_[i] *= rowScale_[i];
     4133      }
     4134    }
     4135  } else if (rhsScale_!=1.0||objectiveScale_!=1.0) {
     4136    // Collect infeasibilities
     4137    int numberPrimalScaled=0;
     4138    int numberPrimalUnscaled=0;
     4139    int numberDualScaled=0;
     4140    int numberDualUnscaled=0;
     4141    double scaleC = 1.0/objectiveScale_;
     4142    double scaleR = 1.0/rhsScale_;
     4143    for (i=0;i<numberColumns;i++) {
     4144      double valueScaled = columnActivityWork_[i];
     4145      double lowerScaled = columnLowerWork_[i];
     4146      double upperScaled = columnUpperWork_[i];
     4147      if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     4148        if (valueScaled<lowerScaled-primalTolerance_||
     4149            valueScaled>upperScaled+primalTolerance_)
     4150          numberPrimalScaled++;
     4151        else
     4152          upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     4153      }
     4154      columnActivity_[i] = valueScaled*scaleR;
     4155      double value = columnActivity_[i];
     4156      if (value<columnLower_[i]-primalTolerance_)
     4157        numberPrimalUnscaled++;
     4158      else if (value>columnUpper_[i]+primalTolerance_)
     4159        numberPrimalUnscaled++;
     4160      double valueScaledDual = reducedCostWork_[i];
     4161      if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     4162        numberDualScaled++;
     4163      if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     4164        numberDualScaled++;
     4165      reducedCost_[i] = valueScaledDual*scaleC;
     4166      double valueDual = reducedCost_[i];
     4167      if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     4168        numberDualUnscaled++;
     4169      if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     4170        numberDualUnscaled++;
     4171    }
     4172    for (i=0;i<numberRows;i++) {
     4173      double valueScaled = rowActivityWork_[i];
     4174      double lowerScaled = rowLowerWork_[i];
     4175      double upperScaled = rowUpperWork_[i];
     4176      if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     4177        if (valueScaled<lowerScaled-primalTolerance_||
     4178            valueScaled>upperScaled+primalTolerance_)
     4179          numberPrimalScaled++;
     4180        else
     4181          upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     4182      }
     4183      rowActivity_[i] = valueScaled*scaleR;
     4184      double value = rowActivity_[i];
     4185      if (value<rowLower_[i]-primalTolerance_)
     4186        numberPrimalUnscaled++;
     4187      else if (value>rowUpper_[i]+primalTolerance_)
     4188        numberPrimalUnscaled++;
     4189      double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
     4190      if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     4191        numberDualScaled++;
     4192      if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     4193        numberDualScaled++;
     4194      dual_[i] *= scaleC;
     4195      double valueDual = dual_[i];
     4196      if (rowObjective_)
     4197        valueDual += rowObjective_[i];
     4198      if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     4199        numberDualUnscaled++;
     4200      if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     4201        numberDualUnscaled++;
     4202    }
     4203    if (!problemStatus_&&!secondaryStatus_) {
     4204      // See if we need to set secondary status
     4205      if (numberPrimalUnscaled) {
     4206        if (numberDualUnscaled)
     4207          secondaryStatus_=4;
     4208        else
     4209          secondaryStatus_=2;
    46094210      } else {
    4610         assert (inverseColumnScale_);
    4611         const double * inverseScale = inverseColumnScale_;
    4612         for (i=0;i<numberColumns;i++) {
    4613           double scaleFactor = columnScale_[i];
    4614           double valueScaled = columnActivityWork_[i];
    4615           double lowerScaled = columnLowerWork_[i];
    4616           double upperScaled = columnUpperWork_[i];
    4617           if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    4618             if (valueScaled<lowerScaled-primalTolerance_||
    4619                 valueScaled>upperScaled+primalTolerance_)
    4620               numberPrimalScaled++;
    4621             else
    4622               upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    4623           }
    4624           columnActivity_[i] = valueScaled*scaleFactor*scaleR;
    4625           double value = columnActivity_[i];
    4626           if (value<columnLower_[i]-primalTolerance_)
    4627             numberPrimalUnscaled++;
    4628           else if (value>columnUpper_[i]+primalTolerance_)
    4629             numberPrimalUnscaled++;
    4630           double valueScaledDual = reducedCostWork_[i];
    4631           if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    4632             numberDualScaled++;
    4633           if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    4634             numberDualScaled++;
    4635           reducedCost_[i] = (valueScaledDual*scaleC)*inverseScale[i];
    4636           double valueDual = reducedCost_[i];
    4637           if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    4638             numberDualUnscaled++;
    4639           if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    4640             numberDualUnscaled++;
    4641         }
    4642         inverseScale = inverseRowScale_;
    4643         for (i=0;i<numberRows;i++) {
    4644           double scaleFactor = rowScale_[i];
    4645           double valueScaled = rowActivityWork_[i];
    4646           double lowerScaled = rowLowerWork_[i];
    4647           double upperScaled = rowUpperWork_[i];
    4648           if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    4649             if (valueScaled<lowerScaled-primalTolerance_||
    4650                 valueScaled>upperScaled+primalTolerance_)
    4651               numberPrimalScaled++;
    4652             else
    4653               upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    4654           }
    4655           rowActivity_[i] = (valueScaled*scaleR)*inverseScale[i];
    4656           double value = rowActivity_[i];
    4657           if (value<rowLower_[i]-primalTolerance_)
    4658             numberPrimalUnscaled++;
    4659           else if (value>rowUpper_[i]+primalTolerance_)
    4660             numberPrimalUnscaled++;
    4661           double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
    4662           if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    4663             numberDualScaled++;
    4664           if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    4665             numberDualScaled++;
    4666           dual_[i] *= scaleFactor*scaleC;
    4667           double valueDual = dual_[i];
    4668           if (rowObjective_)
    4669             valueDual += rowObjective_[i];
    4670           if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    4671             numberDualUnscaled++;
    4672           if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    4673             numberDualUnscaled++;
    4674         }
    4675       }
    4676       if (!problemStatus_&&!secondaryStatus_) {
    4677         // See if we need to set secondary status
    4678         if (numberPrimalUnscaled) {
    4679           if (numberDualUnscaled)
    4680             secondaryStatus_=4;
    4681           else
    4682             secondaryStatus_=2;
    4683         } else {
    4684           if (numberDualUnscaled)
    4685             secondaryStatus_=3;
    4686         }
    4687       }
    4688       if (problemStatus_==2) {
    4689         for (i=0;i<numberColumns;i++) {
    4690           ray_[i] *= columnScale_[i];
    4691         }
    4692       } else if (problemStatus_==1&&ray_) {
    4693         for (i=0;i<numberRows;i++) {
    4694           ray_[i] *= rowScale_[i];
    4695         }
    4696       }
    4697     } else if (rhsScale_!=1.0||objectiveScale_!=1.0) {
    4698       // Collect infeasibilities
    4699       int numberPrimalScaled=0;
    4700       int numberPrimalUnscaled=0;
    4701       int numberDualScaled=0;
    4702       int numberDualUnscaled=0;
    4703       double scaleC = 1.0/objectiveScale_;
    4704       double scaleR = 1.0/rhsScale_;
     4211        if (numberDualUnscaled)
     4212          secondaryStatus_=3;
     4213      }
     4214    }
     4215  } else {
     4216    if (columnActivityWork_) {
    47054217      for (i=0;i<numberColumns;i++) {
    4706         double valueScaled = columnActivityWork_[i];
    4707         double lowerScaled = columnLowerWork_[i];
    4708         double upperScaled = columnUpperWork_[i];
    4709         if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    4710           if (valueScaled<lowerScaled-primalTolerance_||
    4711               valueScaled>upperScaled+primalTolerance_)
    4712             numberPrimalScaled++;
    4713           else
    4714             upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    4715         }
    4716         columnActivity_[i] = valueScaled*scaleR;
    4717         double value = columnActivity_[i];
    4718         if (value<columnLower_[i]-primalTolerance_)
    4719           numberPrimalUnscaled++;
    4720         else if (value>columnUpper_[i]+primalTolerance_)
    4721           numberPrimalUnscaled++;
    4722         double valueScaledDual = reducedCostWork_[i];
    4723         if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    4724           numberDualScaled++;
    4725         if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    4726           numberDualScaled++;
    4727         reducedCost_[i] = valueScaledDual*scaleC;
    4728         double valueDual = reducedCost_[i];
    4729         if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    4730           numberDualUnscaled++;
    4731         if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    4732           numberDualUnscaled++;
     4218        double value = columnActivityWork_[i];
     4219        double lower = columnLowerWork_[i];
     4220        double upper = columnUpperWork_[i];
     4221        if (lower>-1.0e20||upper<1.0e20) {
     4222          if (value>lower&&value<upper)
     4223            upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
     4224        }
     4225        columnActivity_[i] = columnActivityWork_[i];
     4226        reducedCost_[i] = reducedCostWork_[i];
    47334227      }
    47344228      for (i=0;i<numberRows;i++) {
    4735         double valueScaled = rowActivityWork_[i];
    4736         double lowerScaled = rowLowerWork_[i];
    4737         double upperScaled = rowUpperWork_[i];
    4738         if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    4739           if (valueScaled<lowerScaled-primalTolerance_||
    4740               valueScaled>upperScaled+primalTolerance_)
    4741             numberPrimalScaled++;
    4742           else
    4743             upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    4744         }
    4745         rowActivity_[i] = valueScaled*scaleR;
    4746         double value = rowActivity_[i];
    4747         if (value<rowLower_[i]-primalTolerance_)
    4748           numberPrimalUnscaled++;
    4749         else if (value>rowUpper_[i]+primalTolerance_)
    4750           numberPrimalUnscaled++;
    4751         double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
    4752         if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    4753           numberDualScaled++;
    4754         if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    4755           numberDualScaled++;
    4756         dual_[i] *= scaleC;
    4757         double valueDual = dual_[i];
    4758         if (rowObjective_)
    4759           valueDual += rowObjective_[i];
    4760         if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    4761           numberDualUnscaled++;
    4762         if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    4763           numberDualUnscaled++;
    4764       }
    4765       if (!problemStatus_&&!secondaryStatus_) {
    4766         // See if we need to set secondary status
    4767         if (numberPrimalUnscaled) {
    4768           if (numberDualUnscaled)
    4769             secondaryStatus_=4;
    4770           else
    4771             secondaryStatus_=2;
    4772         } else {
    4773           if (numberDualUnscaled)
    4774             secondaryStatus_=3;
    4775         }
    4776       }
    4777     } else {
    4778       if (columnActivityWork_) {
    4779         for (i=0;i<numberColumns;i++) {
    4780           double value = columnActivityWork_[i];
    4781           double lower = columnLowerWork_[i];
    4782           double upper = columnUpperWork_[i];
    4783           if (lower>-1.0e20||upper<1.0e20) {
    4784             if (value>lower&&value<upper)
    4785               upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
    4786           }
    4787           columnActivity_[i] = columnActivityWork_[i];
    4788           reducedCost_[i] = reducedCostWork_[i];
    4789         }
    4790         for (i=0;i<numberRows;i++) {
    4791           double value = rowActivityWork_[i];
    4792           double lower = rowLowerWork_[i];
    4793           double upper = rowUpperWork_[i];
    4794           if (lower>-1.0e20||upper<1.0e20) {
    4795             if (value>lower&&value<upper)
    4796               upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
    4797           }
    4798           rowActivity_[i] = rowActivityWork_[i];
    4799         }
    4800       }
    4801     }
    4802     // switch off scalefactor if auto
    4803     if (automaticScale_) {
    4804       rhsScale_=1.0;
    4805       objectiveScale_=1.0;
    4806     }
    4807     if (optimizationDirection_!=1.0) {
    4808       // and modify all dual signs
    4809       for (i=0;i<numberColumns;i++)
    4810         reducedCost_[i] *= optimizationDirection_;
     4229        double value = rowActivityWork_[i];
     4230        double lower = rowLowerWork_[i];
     4231        double upper = rowUpperWork_[i];
     4232        if (lower>-1.0e20||upper<1.0e20) {
     4233          if (value>lower&&value<upper)
     4234            upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
     4235        }
     4236        rowActivity_[i] = rowActivityWork_[i];
     4237      }
     4238    }
     4239  }
     4240  // switch off scalefactor if auto
     4241  if (automaticScale_) {
     4242    rhsScale_=1.0;
     4243    objectiveScale_=1.0;
     4244  }
     4245  if (optimizationDirection_!=1.0) {
     4246    // and modify all dual signs
     4247    for (i=0;i<numberColumns;i++)
     4248      reducedCost_[i] *= optimizationDirection_;
    48114249      for (i=0;i<numberRows;i++)
    48124250        dual_[i] *= optimizationDirection_;
    4813     }
    4814 #ifdef CLP_AUXILIARY_MODEL
    4815   } else {
    4816     // auxiliary model
    4817     cost_=NULL;
    4818     lower_=NULL;
    4819     upper_=NULL;
    4820     dj_=NULL;
    4821     solution_=NULL;
    4822     int iRow,iColumn;
    4823     for (iRow=0;iRow<4;iRow++) {
    4824       if (rowArray_[iRow]) {
    4825         rowArray_[iRow]->clear();
    4826         //rowArray_[iRow]->checkClear();
    4827       }
    4828       rowArray_[iRow]=NULL;
    4829     }
    4830     for (iColumn=0;iColumn<2;iColumn++) {
    4831       if (columnArray_[iColumn]) {
    4832         columnArray_[iColumn]->clear();
    4833         //columnArray_[iColumn]->checkClear();
    4834       }
    4835       columnArray_[iColumn]=NULL;
    4836     }
    4837     rowCopy_=NULL;
    4838     ClpMatrixBase * temp = matrix_;
    4839     matrix_=auxiliaryModel_->matrix_;
    4840     auxiliaryModel_->matrix_=temp;
    4841     assert ((auxiliaryModel_->specialOptions_&(16+32))==16+32);
    4842     if (problemStatus_!=1&&problemStatus_!=10) {
    4843       int i;
    4844       if (auxiliaryModel_->rowScale_) {
    4845         const double * scale = auxiliaryModel_->columnScale_;
    4846         const double * inverseScale = scale + numberColumns;
    4847         for (i=0;i<numberColumns;i++) {
    4848           columnActivity_[i] = auxiliaryModel_->columnActivityWork_[i]*scale[i];
    4849           reducedCost_[i] = auxiliaryModel_->reducedCostWork_[i]*inverseScale[i];
    4850         }
    4851         scale = auxiliaryModel_->rowScale_;
    4852         inverseScale = scale + numberRows;
    4853         for (i=0;i<numberRows;i++) {
    4854           rowActivity_[i] = auxiliaryModel_->rowActivityWork_[i]*inverseScale[i];
    4855         }
    4856       } else {
    4857         for (i=0;i<numberColumns;i++) {
    4858           columnActivity_[i] = auxiliaryModel_->columnActivityWork_[i];
    4859           reducedCost_[i] = auxiliaryModel_->reducedCostWork_[i];
    4860         }
    4861         for (i=0;i<numberRows;i++) {
    4862           rowActivity_[i] = auxiliaryModel_->rowActivityWork_[i];
    4863         }
    4864       }
    4865       if (optimizationDirection_!=1.0) {
    4866         // and modify reduced costs
    4867         for (i=0;i<numberColumns;i++)
    4868           reducedCost_[i] *= optimizationDirection_;
    4869       }
    4870     } else if (problemStatus_==10) {
    4871       int i;
    4872       if (auxiliaryModel_->rowScale_) {
    4873         const double * scale = auxiliaryModel_->columnScale_;
    4874         const double * inverseScale = scale + numberColumns;
    4875         for (i=0;i<numberColumns;i++) {
    4876           double lower = auxiliaryModel_->columnLowerWork_[i];
    4877           double upper = auxiliaryModel_->columnUpperWork_[i];
    4878           double value = auxiliaryModel_->columnActivityWork_[i];
    4879           if (value>lower&&value<upper) {
    4880             upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
    4881           }
    4882           columnActivity_[i] = value*scale[i];
    4883         }
    4884         scale = auxiliaryModel_->rowScale_;
    4885         inverseScale = scale + numberRows;
    4886         for (i=0;i<numberRows;i++) {
    4887           double lower = auxiliaryModel_->rowLowerWork_[i];
    4888           double upper = auxiliaryModel_->rowUpperWork_[i];
    4889           double value = auxiliaryModel_->rowActivityWork_[i];
    4890           if (value>lower&&value<upper) {
    4891             upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
    4892           }
    4893           rowActivity_[i] = auxiliaryModel_->rowActivityWork_[i]*inverseScale[i];
    4894         }
    4895       } else {
    4896         for (i=0;i<numberColumns;i++) {
    4897           double lower = auxiliaryModel_->columnLowerWork_[i];
    4898           double upper = auxiliaryModel_->columnUpperWork_[i];
    4899           double value = auxiliaryModel_->columnActivityWork_[i];
    4900           if (value>lower&&value<upper) {
    4901             upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
    4902           }
    4903           columnActivity_[i] = value;
    4904         }
    4905         for (i=0;i<numberRows;i++) {
    4906           double lower = auxiliaryModel_->rowLowerWork_[i];
    4907           double upper = auxiliaryModel_->rowUpperWork_[i];
    4908           double value = auxiliaryModel_->rowActivityWork_[i];
    4909           if (value>lower&&value<upper) {
    4910             upperOut_ = CoinMax(upperOut_,CoinMin(value-lower,upper-value));
    4911           }
    4912           rowActivity_[i] = value;
    4913         }
    4914       }
    4915     }
    4916   }
    4917 #endif
     4251  }
    49184252  // scaling may have been turned off
    49194253  scalingFlag_ = abs(scalingFlag_);
     
    56064940  */
    56074941#ifdef COIN_DEVELOP
    5608 #define EXPENSIVE
     4942  //#define EXPENSIVE
    56094943#endif
    56104944#ifdef EXPENSIVE
     
    58135147    return reducedGradient();
    58145148#endif 
    5815   CoinAssert (ifValuesPass>=0&&ifValuesPass<3);
     5149  CoinAssert ((ifValuesPass>=0&&ifValuesPass<3)||
     5150              (ifValuesPass>=12&&ifValuesPass<100)||
     5151              (ifValuesPass>=112&&ifValuesPass<200));
     5152  if (ifValuesPass>=12) {
     5153    int numberProblems = (ifValuesPass-10)%100;
     5154    ifValuesPass = (ifValuesPass<100) ? 1 : 2;
     5155    // Go parallel to do solve
     5156    // Only if all slack basis
     5157    int i;
     5158    for ( i=0;i<numberColumns_;i++) {
     5159      if (getColumnStatus(i)==basic)
     5160        break;
     5161    }
     5162    if (i==numberColumns_) {
     5163      // try if vaguely feasible
     5164      CoinZeroN(rowActivity_,numberRows_);
     5165      const int * row = matrix_->getIndices();
     5166      const CoinBigIndex * columnStart = matrix_->getVectorStarts();
     5167      const int * columnLength = matrix_->getVectorLengths();
     5168      const double * element = matrix_->getElements();
     5169      for (int iColumn=0;iColumn<numberColumns_;iColumn++) {
     5170        CoinBigIndex j;
     5171        double value = columnActivity_[iColumn];
     5172        if (value) {
     5173          CoinBigIndex start = columnStart[iColumn];
     5174          CoinBigIndex end = start + columnLength[iColumn];
     5175          for (j=start; j<end; j++) {
     5176            int iRow=row[j];
     5177            rowActivity_[iRow] += value*element[j];
     5178          }
     5179        }
     5180      }
     5181      checkSolutionInternal();
     5182      if (sumPrimalInfeasibilities_*sqrt(numberRows_)<1.0) {
     5183        // Could do better if can decompose
     5184        // correction to get feasible
     5185        double scaleFactor = 1.0/numberProblems;
     5186        double * correction = new double [numberRows_];
     5187        for (int iRow=0;iRow<numberRows_;iRow++) {
     5188          double value=rowActivity_[iRow];
     5189          if (value>rowUpper_[iRow])
     5190            value = rowUpper_[iRow]-value;
     5191          else if (value<rowLower_[iRow])
     5192            value = rowLower_[iRow]-value;
     5193          else
     5194            value=0.0;
     5195          correction[iRow]=value*scaleFactor;
     5196        }
     5197        int numberColumns = (numberColumns_+numberProblems-1)/numberProblems;
     5198        int * whichRows = new int [numberRows_];
     5199        for (int i=0;i<numberRows_;i++)
     5200          whichRows[i]=i;
     5201        int * whichColumns = new int [numberColumns_];
     5202        ClpSimplex ** model = new ClpSimplex * [numberProblems];
     5203        int startColumn=0;
     5204        double * saveLower = CoinCopyOfArray(rowLower_,numberRows_);
     5205        double * saveUpper = CoinCopyOfArray(rowUpper_,numberRows_);
     5206        for (int i=0;i<numberProblems;i++) {
     5207          int endColumn = CoinMin(startColumn+numberColumns,numberColumns_);
     5208          CoinZeroN(rowActivity_,numberRows_);
     5209          for (int iColumn=startColumn;iColumn<endColumn;iColumn++) {
     5210            whichColumns[iColumn-startColumn]=iColumn;
     5211            CoinBigIndex j;
     5212            double value = columnActivity_[iColumn];
     5213            if (value) {
     5214              CoinBigIndex start = columnStart[iColumn];
     5215              CoinBigIndex end = start + columnLength[iColumn];
     5216              for (j=start; j<end; j++) {
     5217                int iRow=row[j];
     5218                rowActivity_[iRow] += value*element[j];
     5219              }
     5220            }
     5221          }
     5222          // adjust rhs
     5223          for (int iRow=0;iRow<numberRows_;iRow++) {
     5224            double value=rowActivity_[iRow]+correction[iRow];
     5225            if (saveUpper[iRow]<1.0e30)
     5226              rowUpper_[iRow]=value;
     5227            if (saveLower[iRow]>-1.0e30)
     5228              rowLower_[iRow]=value;
     5229          }
     5230          model[i] = new ClpSimplex(this,numberRows_,whichRows,
     5231                                    endColumn-startColumn,whichColumns);
     5232          //#define FEB_TRY
     5233#ifdef FEB_TRY
     5234          model[i]->setPerturbation(perturbation_);
     5235#endif
     5236          startColumn=endColumn;
     5237        }
     5238        memcpy(rowLower_,saveLower,numberRows_*sizeof(double));
     5239        memcpy(rowUpper_,saveUpper,numberRows_*sizeof(double));
     5240        delete [] saveLower;
     5241        delete [] saveUpper;
     5242        delete [] correction;
     5243        // solve (in parallel)
     5244        for (int i=0;i<numberProblems;i++) {
     5245          model[i]->primal(1/*ifValuesPass*/);
     5246        }
     5247        startColumn=0;
     5248        int numberBasic=0;
     5249        // use whichRows as counter
     5250        for (int iRow=0;iRow<numberRows_;iRow++) {
     5251          int startValue=0;
     5252          if (rowUpper_[iRow]>rowLower_[iRow])
     5253            startValue++;
     5254          if (rowUpper_[iRow]>1.0e30)
     5255            startValue++;
     5256          if (rowLower_[iRow]<-1.0e30)
     5257            startValue++;
     5258          whichRows[iRow]=1000*startValue;
     5259        }
     5260        for (int i=0;i<numberProblems;i++) {
     5261          int endColumn = CoinMin(startColumn+numberColumns,numberColumns_);
     5262          ClpSimplex * simplex = model[i];
     5263          const double * solution = simplex->columnActivity_;
     5264          for (int iColumn=startColumn;iColumn<endColumn;iColumn++) {
     5265            columnActivity_[iColumn] = solution[iColumn-startColumn];
     5266            Status status = simplex->getColumnStatus(iColumn-startColumn);
     5267            setColumnStatus(iColumn,status);
     5268            if (status==basic)
     5269              numberBasic++;
     5270          }
     5271          for (int iRow=0;iRow<numberRows_;iRow++) {
     5272            if (simplex->getRowStatus(iRow)==basic)
     5273              whichRows[iRow]++;
     5274          }
     5275          delete model[i];
     5276          startColumn=endColumn;
     5277        }
     5278        delete [] model;
     5279        for (int iRow=0;iRow<numberRows_;iRow++)
     5280          setRowStatus(iRow,superBasic);
     5281        CoinZeroN(rowActivity_,numberRows_);
     5282        for (int iColumn=0;iColumn<numberColumns_;iColumn++) {
     5283          CoinBigIndex j;
     5284          double value = columnActivity_[iColumn];
     5285          if (value) {
     5286            CoinBigIndex start = columnStart[iColumn];
     5287            CoinBigIndex end = start + columnLength[iColumn];
     5288            for (j=start; j<end; j++) {
     5289              int iRow=row[j];
     5290              rowActivity_[iRow] += value*element[j];
     5291            }
     5292          }
     5293        }
     5294        checkSolutionInternal();
     5295        if (numberBasic<numberRows_) {
     5296          int * order = new int [numberRows_];
     5297          for (int iRow=0;iRow<numberRows_;iRow++) {
     5298            setRowStatus(iRow,superBasic);
     5299            int nTimes = whichRows[iRow]%1000;
     5300            if (nTimes)
     5301              nTimes += whichRows[iRow]/500;
     5302            whichRows[iRow]=-nTimes;
     5303            order[iRow]=iRow;
     5304          }
     5305          CoinSort_2(whichRows,whichRows+numberRows_,order);
     5306          int nPut = numberRows_-numberBasic;
     5307          for (int i=0;i<nPut;i++) {
     5308            int iRow = order[i];
     5309            setRowStatus(iRow,basic);
     5310          }
     5311          delete [] order;
     5312        } else if (numberBasic>numberRows_) {
     5313          double * away = new double [numberBasic];
     5314          numberBasic=0;
     5315          for (int iColumn=0;iColumn<numberColumns_;iColumn++) {
     5316            if (getColumnStatus(iColumn)==basic) {
     5317              double value = columnActivity_[iColumn];
     5318              value = CoinMin(value-columnLower_[iColumn],
     5319                              columnUpper_[iColumn]-value);
     5320              away[numberBasic]=value;
     5321              whichColumns[numberBasic++]=iColumn;
     5322            }
     5323          }
     5324          CoinSort_2(away,away+numberBasic,whichColumns);
     5325          int nPut = numberBasic-numberRows_;
     5326          for (int i=0;i<nPut;i++) {
     5327            int iColumn = whichColumns[i];
     5328            double value = columnActivity_[iColumn];
     5329            if (value-columnLower_[iColumn]<
     5330                columnUpper_[iColumn]-value)
     5331              setColumnStatus(iColumn,atLowerBound);
     5332            else
     5333              setColumnStatus(iColumn,atUpperBound);
     5334          }
     5335          delete [] away;
     5336        }
     5337        delete [] whichColumns;
     5338        delete [] whichRows;
     5339      }
     5340    }
     5341  }
    58165342  /*  Note use of "down casting".  The only class the user sees is ClpSimplex.
    58175343      Classes ClpSimplexDual, ClpSimplexPrimal, (ClpSimplexNonlinear)
     
    60345560#include "ClpPredictorCorrector.hpp"
    60355561#include "ClpCholeskyBase.hpp"
    6036 // Preference is WSSMP, TAUCS, UFL (just ordering) then base
     5562// Preference is WSSMP, UFL (just ordering), MUMPS, TAUCS then base
    60375563#if WSSMP_BARRIER
    60385564#include "ClpCholeskyWssmp.hpp"
     
    60405566#elif UFL_BARRIER
    60415567#include "ClpCholeskyUfl.hpp"
     5568#elif MUMPS_BARRIER
     5569#include "ClpCholeskyMumps.hpp"
    60425570#elif TAUCS_BARRIER
    60435571#include "ClpCholeskyTaucs.hpp"
     
    60595587  // If Quadratic we need KKT
    60605588  bool doKKT = (quadraticObj!=NULL);
    6061   // Preference is WSSMP, UFL, TAUCS then base
     5589  // Preference is WSSMP, UFL, MUMPS, TAUCS then base
    60625590#ifdef WSSMP_BARRIER
    60635591 if (!doKKT) {
     
    60835611 ClpCholeskyTaucs * cholesky = new ClpCholeskyTaucs();
    60845612 barrier.setCholesky(cholesky);
     5613#elif MUMPS_BARRIER
     5614 if (!doKKT) {
     5615   ClpCholeskyMumps * cholesky = new ClpCholeskyMumps();
     5616   barrier.setCholesky(cholesky);
     5617 } else {
     5618   printf("***** Unable to do Mumps with KKT\n");
     5619   ClpCholeskyBase * cholesky = new ClpCholeskyBase();
     5620   cholesky->setKKT(true);
     5621   barrier.setCholesky(cholesky);
     5622 }
    60855623#else
    60865624 if (!doKKT) {
     
    70816619        }
    70826620      } else {
     6621        double * inverseColumnScale=columnScale_+numberColumns_;
    70836622        for (i=0;i<numberColumns_;i++) {
    7084           solution_[i] = columnActivity_[i]*(rhsScale_/columnScale_[i]);
     6623          solution_[i] = columnActivity_[i]*(rhsScale_*inverseColumnScale[i]);
    70856624        }
    70866625      }
     
    82297768  if ((startFinishOptions&2)!=0&&(whatsChanged_&(2+512))==2+512)
    82307769    useFactorization=true; // Keep factorization if possible
     7770#if 0
     7771  // seems to be needed if rows deleted later in CbcModel!
     7772  if (!solution_&&scaledMatrix_) {
     7773    // get rid of scaled matrix
     7774    if (scaledMatrix_->getNumRows()!=numberRows_) {
     7775      delete scaledMatrix_;
     7776      scaledMatrix_=NULL;
     7777    }
     7778  }
     7779#endif
    82317780  // sanity check
    82327781  // bad if empty (trap here to avoid using bad matrix_)
     
    83167865      if (perturbation_<100) {
    83177866        if (algorithm_>0&&(objective_->type()<2||!objective_->activated())) {
     7867#ifndef FEB_TRY
    83187868          static_cast<ClpSimplexPrimal *> (this)->perturb(0);
     7869#endif
    83197870        } else if (algorithm_<0) {
    83207871          static_cast<ClpSimplexDual *> (this)->perturb();
     
    83317882    int numberThrownOut = -1;
    83327883    int totalNumberThrownOut=0;
     7884    problemStatus_=-1;
    83337885    // see if we are re-using factorization
    83347886    if (!useFactorization) {
     
    84167968    factorization_->setDenseThreshold(saveThreshold);
    84177969   
    8418     problemStatus_ = -1;
     7970    if (!numberPrimalInfeasibilities_&&!numberDualInfeasibilities_
     7971        &&!ifValuesPass)
     7972      problemStatus_=0;
     7973    else
     7974      assert(problemStatus_ == -1);
    84197975   
    84207976    // number of times we have declared optimality
     
    84407996    getRidOfData=0; // Keep stuff
    84417997    // mark all as current
    8442     whatsChanged_ = 0xffff;
     7998    whatsChanged_ = 0x3ffffff;
    84437999  } else {
    8444     whatsChanged_=0;
     8000    whatsChanged_ &= ~0xffff;
    84458001  }
    84468002  deleteRim(getRidOfData);
     
    86948250    wholeModel->rowCopy_=NULL;
    86958251  }
    8696   whatsChanged_=0;
     8252  whatsChanged_ &= ~0xffff;
    86978253  CoinAssert (wholeModel->matrix_);
    86988254  wholeModel->matrix_ = wholeModel->matrix_->subsetClone(numberRows_,whichRow,
     
    90808636      value = columnScale_[pivot];
    90818637    else
    9082       value = -1.0/rowScale_[pivot-numberColumns_];
     8638      value = -1.0*inverseRowScale_[pivot-numberColumns_];
    90838639  }
    90848640  rowArray1->insert(row,value);
     
    90928648    double * array = columnArray0->denseVector();
    90938649    for (int i=0;i<numberColumns_;i++)
    9094       z[i] = array[i]/columnScale_[i];
     8650      z[i] = array[i]*inverseColumnScale_[i];
    90958651  }
    90968652  if (slack) {
     
    91828738    if (col<numberColumns_) {
    91838739      unpack(rowArray1,col);
    9184       double multiplier = 1.0/columnScale_[col];
     8740      double multiplier = 1.0*inverseColumnScale_[col];
    91858741      int number = rowArray1->getNumElements();
    91868742      int * index = rowArray1->getIndices();
     
    92968852        objectiveWork_[elementIndex] = direction*elementValue;
    92978853      } else {
    9298         assert (!auxiliaryModel_);
    92998854        objectiveWork_[elementIndex] = direction*elementValue
    93008855          *columnScale_[elementIndex];
     
    93258880        rowLowerWork_[elementIndex] = elementValue * rhsScale_;
    93268881      } else {
    9327         assert (!auxiliaryModel_);
    93288882        rowLowerWork_[elementIndex] = elementValue * rhsScale_
    93298883          * rowScale_[elementIndex];
     
    93558909        rowUpperWork_[elementIndex] = elementValue * rhsScale_;
    93568910      } else {
    9357         assert (!auxiliaryModel_);
    93588911        rowUpperWork_[elementIndex] = elementValue * rhsScale_
    93598912          * rowScale_[elementIndex];
     
    93888941        rowLowerWork_[elementIndex] = lowerValue * rhsScale_;
    93898942      } else {
    9390         assert (!auxiliaryModel_);
    93918943        rowLowerWork_[elementIndex] = lowerValue * rhsScale_
    93928944          * rowScale_[elementIndex];
     
    94048956        rowUpperWork_[elementIndex] = upperValue * rhsScale_;
    94058957      } else {
    9406         assert (!auxiliaryModel_);
    94078958        rowUpperWork_[elementIndex] = upperValue * rhsScale_
    94088959          * rowScale_[elementIndex];
     
    94549005        rowLowerWork_[iRow] = rowLower_[iRow] * rhsScale_;
    94559006      } else {
    9456         assert (!auxiliaryModel_);
    94579007        rowLowerWork_[iRow] = rowLower_[iRow] * rhsScale_
    94589008          * rowScale_[iRow];
     
    94639013        rowUpperWork_[iRow] = rowUpper_[iRow] * rhsScale_;
    94649014      } else {
    9465         assert (!auxiliaryModel_);
    94669015        rowUpperWork_[iRow] = rowUpper_[iRow] * rhsScale_
    94679016          * rowScale_[iRow];
     
    94959044        value = elementValue * rhsScale_;
    94969045      } else {
    9497         assert (!auxiliaryModel_);
    94989046        value = elementValue * rhsScale_
    94999047          / columnScale_[elementIndex];
     
    95309078        value = elementValue * rhsScale_;
    95319079      } else {
    9532         assert (!auxiliaryModel_);
    95339080        value = elementValue * rhsScale_
    95349081          / columnScale_[elementIndex];
     
    95659112        lower_[elementIndex] = lowerValue * rhsScale_;
    95669113      } else {
    9567         assert (!auxiliaryModel_);
    95689114        lower_[elementIndex] = lowerValue * rhsScale_
    95699115          / columnScale_[elementIndex];
     
    95849130        upper_[elementIndex] = upperValue * rhsScale_;
    95859131      } else {
    9586         assert (!auxiliaryModel_);
    95879132        upper_[elementIndex] = upperValue * rhsScale_
    95889133          / columnScale_[elementIndex];
     
    96349179        lower_[iColumn] = columnLower_[iColumn] * rhsScale_;
    96359180      } else {
    9636         assert (!auxiliaryModel_);
    96379181        lower_[iColumn] = columnLower_[iColumn] * rhsScale_
    96389182          / columnScale_[iColumn];
     
    96439187        upper_[iColumn] = columnUpper_[iColumn] * rhsScale_;
    96449188      } else {
    9645         assert (!auxiliaryModel_);
    96469189        upper_[iColumn] = columnUpper_[iColumn] * rhsScale_
    96479190          / columnScale_[iColumn];
     
    96779220    double lower = rowLower_[iRow];
    96789221    double upper = rowUpper_[iRow];
     9222    ClpSimplex::Status status = getRowStatus(iRow);
     9223    if (status!=basic) {
     9224      if (lower==upper) {
     9225        status = ClpSimplex::isFixed;
     9226      } else if (primalValue>upper-primalTolerance) {
     9227        status = ClpSimplex::atUpperBound;
     9228      } else if (primalValue<lower+primalTolerance) {
     9229        status = ClpSimplex::atLowerBound;
     9230      }
     9231      setRowStatus(iRow,status);
     9232    }
    96799233    if (primalValue>upper+primalTolerance) {
    96809234      sumPrimalInfeasibilities_ += primalValue-upper-primalTolerance;
     
    96849238      numberPrimalInfeasibilities_ ++;
    96859239    } else {
    9686       switch(getRowStatus(iRow)) {
     9240      switch(status) {
    96879241       
    96889242      case basic:
     
    97299283    double lower = columnLower_[iColumn];
    97309284    double upper = columnUpper_[iColumn];
     9285    ClpSimplex::Status status = getColumnStatus(iColumn);
     9286    if (status!=basic&&lower==upper) {
     9287      status = ClpSimplex::isFixed;
     9288      setColumnStatus(iColumn,ClpSimplex::isFixed);
     9289    }
    97319290    if (primalValue>upper+primalTolerance) {
    97329291      sumPrimalInfeasibilities_ += primalValue-upper-primalTolerance;
     
    97369295      numberPrimalInfeasibilities_ ++;
    97379296    } else {
    9738       switch(getColumnStatus(iColumn)) {
     9297      switch(status) {
    97399298       
    97409299      case basic:
     
    97919350    problemStatus_=-1;
    97929351}
    9793 #ifdef CLP_AUXILIARY_MODEL
    9794 /*
    9795   If you are re-using the same matrix again and again then the setup time
    9796   to do scaling may be significant.  Also you may not want to initialize all values
    9797   or return all values (especially if infeasible).  While an auxiliary model exists
    9798   it will be faster.  If options -1 then model is switched off.  Otherwise switched on
    9799   with following options.
    9800   1 - rhs is constant
    9801   2 - bounds are constant
    9802   4 - objective is constant
    9803   8 - solution in by basis and no djs etc in
    9804   16 - no duals out (but reduced costs)
    9805   32 - no output if infeasible
    9806 */
    9807 void
    9808 ClpSimplex::auxiliaryModel(int options)
    9809 {
    9810   delete auxiliaryModel_;
    9811   auxiliaryModel_=NULL;
    9812   if (options>=0) {
    9813     createRim(63,true,0);
    9814     auxiliaryModel_ = new ClpSimplex(true);
    9815     auxiliaryModel_->specialOptions_=options;
    9816     int i;
    9817     int numberRows2 = numberRows_+numberExtraRows_;
    9818     int numberTotal = numberRows2+numberColumns_;
    9819     auxiliaryModel_->numberRows_=numberRows_;
    9820     auxiliaryModel_->numberColumns_=numberColumns_;
    9821     if (rowScale_) {
    9822       auxiliaryModel_->rowScale_= new double [2*numberRows_];
    9823       for (i=0;i<numberRows_;i++) {
    9824         auxiliaryModel_->rowScale_[i]=rowScale_[i];
    9825         auxiliaryModel_->rowScale_[numberRows_+i]=1.0/rowScale_[i];
    9826       }
    9827       auxiliaryModel_->columnScale_= new double [2*numberColumns_];
    9828       for (i=0;i<numberColumns_;i++) {
    9829         auxiliaryModel_->columnScale_[i]=columnScale_[i];
    9830         auxiliaryModel_->columnScale_[numberColumns_+i]=1.0/columnScale_[i];
    9831       }
    9832     }
    9833     // copy anyway
    9834     auxiliaryModel_->lower_ = new double [2*numberTotal];
    9835     CoinMemcpyN(lower_,numberTotal,auxiliaryModel_->lower_);
    9836     CoinMemcpyN(lower_,numberTotal,auxiliaryModel_->lower_+numberTotal);
    9837     auxiliaryModel_->upper_ = new double [2*numberTotal];
    9838     CoinMemcpyN(upper_,numberTotal,auxiliaryModel_->upper_);
    9839     CoinMemcpyN(upper_,numberTotal,auxiliaryModel_->upper_+numberTotal);
    9840     auxiliaryModel_->cost_ = new double [2*numberTotal];
    9841     CoinMemcpyN(cost_,numberTotal,auxiliaryModel_->cost_);
    9842     CoinMemcpyN(cost_,numberTotal,auxiliaryModel_->cost_+numberTotal);
    9843     auxiliaryModel_->dj_ = new double [2*numberTotal];
    9844     memset(auxiliaryModel_->dj_,0,2*numberTotal*sizeof(double));
    9845     auxiliaryModel_->solution_ = new double [2*numberTotal];
    9846     memset(auxiliaryModel_->solution_,0,2*numberTotal*sizeof(double));
    9847     auxiliaryModel_->reducedCostWork_ = auxiliaryModel_->dj_;
    9848     auxiliaryModel_->rowReducedCost_ = auxiliaryModel_->dj_+numberColumns_;
    9849     auxiliaryModel_->columnActivityWork_ = auxiliaryModel_->solution_;
    9850     auxiliaryModel_->rowActivityWork_ = auxiliaryModel_->solution_+numberColumns_;
    9851     auxiliaryModel_->objectiveWork_ = auxiliaryModel_->cost_;
    9852     auxiliaryModel_->rowObjectiveWork_ = auxiliaryModel_->cost_+numberColumns_;
    9853     auxiliaryModel_->rowLowerWork_ = auxiliaryModel_->lower_+numberColumns_;
    9854     auxiliaryModel_->columnLowerWork_ = auxiliaryModel_->lower_;
    9855     auxiliaryModel_->rowUpperWork_ = auxiliaryModel_->upper_+numberColumns_;
    9856     auxiliaryModel_->columnUpperWork_ = auxiliaryModel_->upper_;
    9857     // delete in model
    9858     delete [] lower_;
    9859     lower_=NULL;
    9860     delete [] upper_;
    9861     upper_=NULL;
    9862     delete [] cost_;
    9863     cost_=NULL;
    9864     delete [] dj_;
    9865     dj_=NULL;
    9866     delete [] solution_;
    9867     solution_=NULL;
    9868     auxiliaryModel_->rowCopy_=rowCopy_;
    9869     ClpMatrixBase * columnCopy = rowCopy_->reverseOrderedCopy();
    9870     rowCopy_=NULL;
    9871     // Point to correct as will be restored
    9872     auxiliaryModel_->matrix_ = matrix_;
    9873     int iRow,iColumn;
    9874     for (iRow=0;iRow<4;iRow++) {
    9875       auxiliaryModel_->rowArray_[iRow]=rowArray_[iRow];
    9876       rowArray_[iRow]=NULL;
    9877     }
    9878     for (iColumn=0;iColumn<2;iColumn++) {
    9879       auxiliaryModel_->columnArray_[iColumn]=columnArray_[iColumn];
    9880       columnArray_[iColumn]=NULL;
    9881     }
    9882     // Set status so won't overwrite solution
    9883     int saveStat=problemStatus_;
    9884     problemStatus_=1;
    9885     deleteRim(1);
    9886     problemStatus_=saveStat;
    9887     auxiliaryModel_->matrix_ = columnCopy;
    9888     // Get rid of scaling
    9889     delete [] rowScale_;
    9890     delete [] columnScale_;
    9891     rowScale_=NULL;
    9892     columnScale_=NULL;
    9893     whatsChanged_ &= ~1;
    9894   }
    9895 }
    9896 #endif
    98979352/*
    98989353  When scaling is on it is possible that the scaled problem
     
    99809435  }
    99819436  return basis;
    9982 }
    9983 #endif
    9984 #ifdef CLP_AUXILIARY_MODEL
    9985 // Switch off e.g. if people using presolve
    9986 void
    9987 ClpSimplex::deleteAuxiliaryModel()
    9988 {
    9989   delete auxiliaryModel_;
    9990   auxiliaryModel_=NULL;
    99919437}
    99929438#endif
     
    103979843      }
    103989844      if (returnCode) {
     9845        bool fixBounds = (info->nNodes_>=0) ? true : false;
    103999846        //check this does everything
    104009847        static_cast<ClpSimplexOther *> (this)->afterCrunch(*small,
     
    104069853            assert (fabs(value-value2)<1.0e-4);
    104079854            columnActivity_[i]=value2;
    10408             columnLower_[i]=value2;
    10409             columnUpper_[i]=value2;
     9855            if (fixBounds) {
     9856              columnLower_[i]=value2;
     9857              columnUpper_[i]=value2;
     9858            }
    104109859          }
    104119860        }
     
    1081610265  ClpNodeStuff * info = reinterpret_cast<ClpNodeStuff *> (stuff);
    1081710266  int nNodes = info->maximumNodes();
     10267  int putNode = info->maximumSpace();
    1081810268  int goodNodes=0;
    1081910269  info->nNodes_=0;
     
    1082410274  double limit = 0.0;
    1082510275  getDblParam(ClpDualObjectiveLimit, limit);
    10826   for (int j=0;j<nNodes;j++) {
     10276  for (int j=0;j<putNode;j++) {
    1082710277    if (nodeInfo[j]) {
    1082810278      nodeInfo[j]->setObjectiveValue(limit);
     
    1097510425    return whichSolution;
    1097610426  }
     10427#ifndef DEBUG
    1097710428  {
    1097810429    int nBasic=0;
     
    1098410435    assert (nBasic==numberRows_);
    1098510436  }
     10437#endif
    1098610438  int returnCode = startFastDual2(info);
    1098710439  if (returnCode) {
     
    1106110513#endif
    1106210514  /* Use nodeInfo for storage
    11063      depth 0 will be nNodes-1, 1 nNodes-2 etc */
    11064   int useDepth=nNodes-1;
     10515     depth 0 will be putNode-1, 1 putNode-2 etc */
     10516  int useDepth=putNode-1;
    1106510517  bool lastFeasible=true;
    1106610518  bool justDive = (info->solverOptions_&32)!=0;
     10519  //printf("putNode %d nDepth %d\n");
    1106710520  while (depth>=0) {
    1106810521    bool stopAtOnce=false;
     
    1109210545      CoinMemcpyN(saveUpper,numberColumns_,columnUpper_);
    1109310546      for (int i=0;i<depth;i++) {
    11094         nodeInfo[nNodes-1-i]->applyNode(this,0);
     10547        nodeInfo[putNode-1-i]->applyNode(this,0);
    1109510548      }
    1109610549      nodeInfo[useDepth]->applyNode(this,1);
     
    1130710760          whichSolution=goodNodes;
    1130810761          goodNodes++;
     10762          if (goodNodes>=nNodes)
     10763            justDive=true; // clean up phase
    1130910764          assert (node->sequence()<0);
    1131010765          bestObjective = objectiveValue-increment;
     
    1136510820          }
    1136610821          goodNodes++;
     10822          if (goodNodes>=nNodes)
     10823            justDive=true; // clean up phase
    1136710824          backtrack=true;
    1136810825        } else {
     
    1137410831    }
    1137510832  }
     10833  //printf("nNodes %d nDepth %d, useDepth %d goodNodes %d\n",
     10834  // nNodes,info->nDepth_,useDepth,goodNodes);
    1137610835#ifdef CHECK_PATH
    1137710836  if (startOptimal) {
     
    1145010909  assert (clpMatrix&&(clpMatrix->flags()&1)==0);
    1145110910#endif
    11452   if (!inverseColumnScale_&&columnScale_) {
    11453     if ((info->solverOptions_&1)!=0) {
    11454       double * temp = new double [2*numberColumns_];
    11455       for (int i=0;i<numberColumns_;i++) {
    11456         double value = columnScale_[i];
    11457         temp[i]=value;
    11458         temp[i+numberColumns_]=1.0/value;
    11459       }
    11460       delete [] columnScale_;
    11461       if (savedColumnScale_&&savedColumnScale_==columnScale_)
    11462         savedColumnScale_=NULL;
    11463       if (savedRowScale_&&savedRowScale_==rowScale_)
    11464         savedRowScale_=NULL;
    11465       columnScale_ = temp;
    11466     }
    11467     if ((info->solverOptions_&4)!=0) {
    11468       double * temp = new double [2*numberRows_];
    11469       for (int i=0;i<numberRows_;i++) {
    11470         double value = rowScale_[i];
    11471         temp[i]=value;
    11472         temp[i+numberRows_]=1.0/value;
    11473       }
    11474       if (savedColumnScale_&&savedColumnScale_==columnScale_)
    11475         savedColumnScale_=NULL;
    11476       if (savedRowScale_&&savedRowScale_==rowScale_)
    11477         savedRowScale_=NULL;
    11478       delete [] rowScale_;
    11479       rowScale_ = temp;
    11480     }
    11481   }
    1148210911  // mark all as current
    11483   whatsChanged_ = 0xffff;
     10912  whatsChanged_ = 0x3ffffff;
    1148410913
    1148510914  // change newLower and newUpper if scaled
     
    1169411123    factorization_->setPersistenceFlag(0);
    1169511124  deleteRim(1);
    11696   whatsChanged_=0;
     11125  whatsChanged_ &= ~0xffff;
    1169711126  assert ((info->solverOptions_&65536)!=0);
    1169811127  info->solverOptions_ &= ~65536;
  • trunk/Clp/src/ClpSimplex.hpp

    r1343 r1344  
    122122  */
    123123  void setPersistenceFlag(int value);
    124 #ifdef CLP_AUXILIARY_MODEL
    125   /**
    126      If you are re-using the same matrix again and again then the setup time
    127      to do scaling may be significant.  Also you may not want to initialize all values
    128      or return all values (especially if infeasible).  While an auxiliary model exists
    129      it will be faster.  If options -1 then model is switched off.  Otherwise switched on
    130      with following options.
    131      1 - rhs is constant
    132      2 - bounds are constant
    133      4 - objective is constant
    134      8 - solution in by basis and no djs etc in
    135      16 - no duals out (but reduced costs)
    136      32 - no output if infeasible
    137   */
    138   void auxiliaryModel(int options);
    139   /// Switch off e.g. if people using presolve
    140   void deleteAuxiliaryModel();
    141   /// See if we have auxiliary model
    142   inline bool usingAuxiliaryModel() const
    143   { return auxiliaryModel_!=NULL;}
    144 #endif
    145124  /// Save a copy of model with certain state - normally without cuts
    146125  void makeBaseModel();
     
    12981277  /// Column activities - working copy
    12991278  double * columnActivityWork_;
    1300 #ifdef CLP_AUXILIARY_MODEL
    1301   /// Auxiliary model
    1302   ClpSimplex * auxiliaryModel_;
    1303 #endif
    13041279  /// Number of dual infeasibilities
    13051280  int numberDualInfeasibilities_;
  • trunk/Clp/src/ClpSimplexDual.cpp

    r1343 r1344  
    127127extern int * nextDescendent;
    128128#endif
    129 #ifndef CLP_AUXILIARY_MODEL
    130 #define auxiliaryModel_ false
    131 #endif
    132129#ifdef NDEBUG
    133130#define NDEBUG_CLP
     
    236233  // Do initial factorization
    237234  // If user asked for perturbation - do it
     235  numberFake_ =0; // Number of variables at fake bounds
     236  numberChanged_ =0; // Number of variables with changed costs
    238237  if (!startup(0,startFinishOptions)) {
    239238    int usePrimal=0;
     
    247246      if (scalingFlag_>0) {
    248247        for (i=0;i<numberRows_;i++) {
    249           dual_[i] = saveDuals[i]/rowScale_[i];
     248          dual_[i] = saveDuals[i]*inverseRowScale_[i];
    250249        }
    251250      } else {
     
    303302
    304303    double objectiveChange;
    305     numberFake_ =0; // Number of variables at fake bounds
    306     numberChanged_ =0; // Number of variables with changed costs
    307     changeBounds(1,NULL,objectiveChange);
     304    assert (!numberFake_);
     305    assert (numberChanged_ ==0);
     306    if (!numberFake_) // if nonzero then adjust
     307      changeBounds(1,NULL,objectiveChange);
    308308   
    309309    if (!ifValuesPass) {
     
    320320#ifdef CLP_INVESTIGATE
    321321      if (numberDualInfeasibilities_)
    322         printf("ZZZ %d primal %d dual - cost %g\n",
     322        printf("ZZZ %d primal %d dual - sumdinf %g\n",
    323323               numberPrimalInfeasibilities_,
    324                numberDualInfeasibilities_,cost_[0]);
     324               numberDualInfeasibilities_,sumDualInfeasibilities_);
    325325#endif
    326326      if (handler_->logLevel()>2) {
     
    360360    if (usePrimal) {
    361361      problemStatus_=10;
    362 #if COIN_DEVELOP>1
     362#if COIN_DEVELOP>2
    363363      printf("returning to use primal (no obj) at %d\n",__LINE__);
    364364#endif
     
    570570  // Save so can see if doing after primal
    571571  int initialStatus=problemStatus_;
     572  if (!returnCode&&!numberDualInfeasibilities_&&
     573      !numberPrimalInfeasibilities_&&perturbation_<101) {
     574    returnCode=1; // to skip gutsOfDual
     575    problemStatus_=0;
     576  }
    572577 
    573578  if (!returnCode)
     
    617622      if (scalingFlag_>0) {
    618623        for (i=0;i<numberRows_;i++) {
    619           dual_[i] = saveDuals[i]/rowScale_[i];
     624          dual_[i] = saveDuals[i]*inverseRowScale_[i];
    620625        }
    621626      } else {
     
    832837  z_thinks=-1;
    833838#endif
    834 #if 0
    835   if (!numberIterations_&&auxiliaryModel_) {
    836     for (int i=0;i<numberColumns_;i++) {
    837       if (!columnLower_[i]==auxiliaryModel_->lowerRegion()[i+numberRows_+numberColumns_])
    838         {printf("%d a\n",i);break;}
    839       if (!columnUpper_[i]==auxiliaryModel_->upperRegion()[i+numberRows_+numberColumns_])
    840         {printf("%d b %g\n",i,auxiliaryModel_->upperRegion()[i+numberRows_+numberColumns_]);break;}
    841     }
    842     for (int i=0 ;i<numberRows_;i++) {
    843       if (!rowLower_[i]==auxiliaryModel_->lowerRegion()[i+numberRows_+2*numberColumns_])
    844         {printf("%d c\n",i);break;}
    845       if (!rowUpper_[i]==auxiliaryModel_->upperRegion()[i+numberRows_+2*numberColumns_])
    846         {printf("%d d\n",i);break;}
    847     }
    848   }
    849 #endif
    850839#ifdef CLP_DEBUG
    851840  int debugIteration=-1;
     
    13231312#ifdef COIN_DEVELOP
    13241313              printf("flag a %g %g\n",btranAlpha,alpha_);
     1314#endif
     1315              //#define FEB_TRY
     1316#ifdef FEB_TRY
     1317              // Make safer?
     1318              factorization_->saferTolerances (1.0e-15,-1.03);
    13251319#endif
    13261320              setFlagged(sequenceOut_);
     
    28222816          if (iSequence<numberColumns_) {
    28232817            if (columnScale_) {
    2824               double multiplier = rhsScale_/columnScale_[iSequence];
     2818              double multiplier = rhsScale_*inverseColumnScale_[iSequence];
    28252819              // lower
    28262820              double value = columnLower_[iSequence];
     
    28762870            double value = columnLower_[iSequence];
    28772871            if (value>-1.0e30) {
    2878               double multiplier = rhsScale_/columnScale_[iSequence];
     2872              double multiplier = rhsScale_*inverseColumnScale_[iSequence];
    28792873              value *= multiplier;
    28802874            }
     
    28852879            double value = columnUpper_[iSequence];
    28862880            if (value<1.0e30) {
    2887               double multiplier = rhsScale_/columnScale_[iSequence];
     2881              double multiplier = rhsScale_*inverseColumnScale_[iSequence];
    28882882              value *= multiplier;
    28892883            }
     
    50755069    numberFake_--;;
    50765070    setFakeBound(iSequence,noFake);
    5077 #ifdef CLP_AUXILIARY_MODEL
    5078     if (auxiliaryModel_) {
    5079       // just copy back
    5080       lower_[iSequence]=auxiliaryModel_->lowerRegion()[iSequence+numberRows_+numberColumns_];
    5081       upper_[iSequence]=auxiliaryModel_->upperRegion()[iSequence+numberRows_+numberColumns_];
    5082       return;
    5083     }
    5084 #endif
    50855071    if (iSequence>=numberColumns_) {
    50865072      // rows
     
    51045090      columnUpperWork_[iSequence]=columnUpper_[iSequence];
    51055091      if (rowScale_) {
    5106         double multiplier = 1.0/columnScale_[iSequence];
     5092        double multiplier = 1.0*inverseColumnScale_[iSequence];
    51075093        if (columnLowerWork_[iSequence]>-1.0e50)
    51085094          columnLowerWork_[iSequence] *= multiplier*rhsScale_;
     
    56295615    // external view - in case really getting optimal
    56305616    columnUpper_[iColumn] = newUpper[i];
     5617    assert (inverseColumnScale_);
    56315618    if (scalingFlag_<=0)
    56325619      upper_[iColumn] = newUpper[i]*rhsScale_;
    56335620    else
    5634       upper_[iColumn] = (newUpper[i]/columnScale_[iColumn])*rhsScale_; // scale
     5621      upper_[iColumn] = (newUpper[i]*inverseColumnScale_[iColumn])*rhsScale_; // scale
    56355622    // Get fake bounds correctly
    56365623    double changeCost;
     
    56965683    // external view - in case really getting optimal
    56975684    columnLower_[iColumn] = newLower[i];
     5685    assert (inverseColumnScale_);
    56985686    if (scalingFlag_<=0)
    56995687      lower_[iColumn] = newLower[i]*rhsScale_;
    57005688    else
    5701       lower_[iColumn] = (newLower[i]/columnScale_[iColumn])*rhsScale_; // scale
     5689      lower_[iColumn] = (newLower[i]*inverseColumnScale_[iColumn])*rhsScale_; // scale
    57025690    // Get fake bounds correctly
    57035691    resetFakeBounds(1);
     
    57915779  if ((startFinishOptions&1)==0) {
    57925780    deleteRim(1);
    5793     whatsChanged_=0;
     5781    whatsChanged_ &= ~0xffff;
    57945782  } else {
    57955783    // Original factorization will have been put back by last loop
     
    57985786    deleteRim(0);
    57995787    // mark all as current
    5800     whatsChanged_ = 0xffff;
     5788    whatsChanged_ = 0x3ffffff;
    58015789  }
    58025790  objectiveValue_ = saveObjectiveValue;
     
    60226010// This does first part of StrongBranching
    60236011ClpFactorization *
    6024 ClpSimplexDual::setupForStrongBranching(char * arrays, int numberRows, int numberColumns)
     6012ClpSimplexDual::setupForStrongBranching(char * arrays, int numberRows,
     6013                                        int numberColumns,bool solveLp)
    60256014{
    6026   algorithm_ = -1;
    6027   // put in standard form (and make row copy)
    6028   // create modifiable copies of model rim and do optional scaling
    6029   int startFinishOptions;
    6030   /*  COIN_CLP_VETTED
    6031       Looks safe for Cbc
    6032   */
    6033   bool safeOption = ((specialOptions_&COIN_CBC_USING_CLP)!=0);
    6034   if((specialOptions_&4096)==0||(auxiliaryModel_&&!safeOption)) {
    6035     startFinishOptions=0;
    6036   } else {
    6037     startFinishOptions=1+2+4;
    6038   }
    6039   createRim(7+8+16+32,true,startFinishOptions);
    6040   // Do initial factorization
    6041   // and set certain stuff
    6042   // We can either set increasing rows so ...IsBasic gives pivot row
    6043   // or we can just increment iBasic one by one
    6044   // for now let ...iBasic give pivot row
    6045   bool useFactorization=false;
    6046   if ((startFinishOptions&2)!=0&&(whatsChanged_&(2+512))==2+512) {
    6047     useFactorization=true; // Keep factorization if possible
    6048     // switch off factorization if bad
    6049     if (pivotVariable_[0]<0||factorization_->numberRows()!=numberRows_)
    6050       useFactorization=false;
    6051   }
    6052   if (!useFactorization) {
    6053     factorization_->setDefaultValues();
    6054 
    6055     int factorizationStatus = internalFactorize(0);
    6056     if (factorizationStatus<0) {
    6057       // some error
    6058       // we should either debug or ignore
     6015  if (solveLp) {
     6016    // solve
     6017    dual(0,7);
     6018    if (problemStatus_==10) {
     6019      ClpSimplex::dual(0,0);
     6020      assert (problemStatus_!=10);
     6021      if (problemStatus_==0) {
     6022        dual(0,7);
     6023        //assert (problemStatus_!=10);
     6024      }
     6025    }
     6026    if (problemStatus_==1)
     6027      return NULL; // say infeasible
     6028    // May be empty
     6029    solveLp = (solution_!=NULL&&problemStatus_==0);
     6030  }
     6031  problemStatus_=0;
     6032  if (!solveLp) {
     6033    algorithm_ = -1;
     6034    // put in standard form (and make row copy)
     6035    // create modifiable copies of model rim and do optional scaling
     6036    int startFinishOptions;
     6037    if((specialOptions_&4096)==0) {
     6038      startFinishOptions=0;
     6039    } else {
     6040      startFinishOptions=1+2+4;
     6041    }
     6042    createRim(7+8+16+32,true,startFinishOptions);
     6043    // Do initial factorization
     6044    // and set certain stuff
     6045    // We can either set increasing rows so ...IsBasic gives pivot row
     6046    // or we can just increment iBasic one by one
     6047    // for now let ...iBasic give pivot row
     6048    bool useFactorization=false;
     6049    if ((startFinishOptions&2)!=0&&(whatsChanged_&(2+512))==2+512) {
     6050      useFactorization=true; // Keep factorization if possible
     6051      // switch off factorization if bad
     6052      if (pivotVariable_[0]<0||factorization_->numberRows()!=numberRows_)
     6053        useFactorization=false;
     6054    }
     6055    if (!useFactorization) {
     6056      factorization_->setDefaultValues();
     6057     
     6058      int factorizationStatus = internalFactorize(0);
     6059      if (factorizationStatus<0) {
     6060        // some error
     6061        // we should either debug or ignore
    60596062#ifndef NDEBUG
    6060       printf("***** ClpDual strong branching factorization error - debug\n");
    6061 #endif
    6062     } else if (factorizationStatus&&factorizationStatus<=numberRows_) {
    6063       handler_->message(CLP_SINGULARITIES,messages_)
    6064         <<factorizationStatus
    6065         <<CoinMessageEol;
     6063        printf("***** ClpDual strong branching factorization error - debug\n");
     6064#endif
     6065      } else if (factorizationStatus&&factorizationStatus<=numberRows_) {
     6066        handler_->message(CLP_SINGULARITIES,messages_)
     6067          <<factorizationStatus
     6068          <<CoinMessageEol;
     6069      }
    60666070    }
    60676071  }
    60686072  // Get fake bounds correctly
    60696073  double dummyChangeCost;
    6070   changeBounds(1,NULL,dummyChangeCost);
     6074  changeBounds(3,NULL,dummyChangeCost);
    60716075  double * arrayD = reinterpret_cast<double *> (arrays);
    60726076  arrayD[0]=objectiveValue()*optimizationDirection_;
     
    60956099          numberRows_+numberColumns_,saveObjective);
    60966100  CoinMemcpyN(pivotVariable_, numberRows_,savePivot);
    6097   return new ClpFactorization(*factorization_,numberRows_);
     6101  ClpFactorization * factorization = factorization_;
     6102  factorization_ = NULL;
     6103  return factorization;
    60986104}
    60996105// This cleans up after strong branching
    61006106void
    6101 ClpSimplexDual::cleanupAfterStrongBranching()
     6107ClpSimplexDual::cleanupAfterStrongBranching(ClpFactorization * factorization)
    61026108{
    61036109  int startFinishOptions;
     
    61056111      Looks safe for Cbc
    61066112  */
    6107   bool safeOption = ((specialOptions_&COIN_CBC_USING_CLP)!=0);
    6108   if((specialOptions_&4096)==0||(auxiliaryModel_&&!safeOption)) {
     6113  if((specialOptions_&4096)==0) {
    61096114    startFinishOptions=0;
    61106115  } else {
     
    61136118  if ((startFinishOptions&1)==0) {
    61146119    deleteRim(1);
    6115     whatsChanged_=0;
    61166120  } else {
    61176121    // Original factorization will have been put back by last loop
    6118     //delete factorization_;
    6119     //factorization_ = new ClpFactorization(saveFactorization);
    6120     deleteRim(0);
     6122    delete factorization_;
     6123    factorization_ = factorization;
     6124    //deleteRim(0);
    61216125    // mark all as current
    6122     whatsChanged_ = 0xffff;
    6123   }
     6126  }
     6127  whatsChanged_ &= ~0xffff;
    61246128}
    61256129/* Checks number of variables at fake bounds.  This is used by fastDual
     
    66356639      const int * thisIndices = column+rowStart[iRow];
    66366640      if (rowScale_) {
    6637         assert (!auxiliaryModel_);
    66386641        // scale row
    66396642        double scale = rowScale_[iRow];
     
    67676770            double value = columnLower_[iSequence];
    67686771            if (value>-1.0e30) {
    6769               double multiplier = rhsScale_/columnScale_[iSequence];
     6772              double multiplier = rhsScale_*inverseColumnScale_[iSequence];
    67706773              value *= multiplier;
    67716774            }
     
    67766779            double value = columnUpper_[iSequence];
    67776780            if (value<1.0e30) {
    6778               double multiplier = rhsScale_/columnScale_[iSequence];
     6781              double multiplier = rhsScale_*inverseColumnScale_[iSequence];
    67796782              value *= multiplier;
    67806783            }
     
    68406843        double upperValue=tempUpper[iSequence];
    68416844        double value = solution_[iSequence];
    6842 #ifndef NDEBUG
    68436845        CoinRelFltEq equal;
    6844 #endif
    68456846        if (fakeStatus==ClpSimplexDual::upperFake) {
    68466847          assert(equal(upper_[iSequence],(lowerValue+dualBound_)));
     
    68726873    if (columnScale_) {
    68736874      for (iSequence=0;iSequence<numberColumns_;iSequence++) {
    6874         double multiplier = rhsScale_/columnScale_[iSequence];
     6875        double multiplier = rhsScale_*inverseColumnScale_[iSequence];
    68756876        // lower
    68766877        double value = columnLower_[iSequence];
  • trunk/Clp/src/ClpSimplexDual.hpp

    r1321 r1344  
    132132                      int startFinishOptions=0);
    133133  /// This does first part of StrongBranching
    134   ClpFactorization * setupForStrongBranching(char * arrays, int numberRows, int numberColumns);
     134  ClpFactorization * setupForStrongBranching(char * arrays, int numberRows,
     135                                             int numberColumns,bool solveLp=false);
    135136  /// This cleans up after strong branching
    136   void cleanupAfterStrongBranching();
     137  void cleanupAfterStrongBranching(ClpFactorization * factorization);
    137138  //@}
    138139
  • trunk/Clp/src/ClpSimplexOther.cpp

    r1331 r1344  
    111111                // can improve
    112112                double movement = (columnScale_==NULL) ? 1.0 :
    113                   rhsScale_/columnScale_[sequenceIncrease];
     113                  rhsScale_*inverseColumnScale_[sequenceIncrease];
    114114                costIncrease = CoinMax(fabs(djValue*movement),costIncrease);
    115115              }
     
    126126                // can improve
    127127                double movement = (columnScale_==NULL) ? 1.0 :
    128                   rhsScale_/columnScale_[sequenceDecrease];
     128                  rhsScale_*inverseColumnScale_[sequenceDecrease];
    129129                costDecrease = CoinMax(fabs(djValue*movement),costDecrease);
    130130              }
     
    159159    }
    160160    double scaleFactor;
    161 #ifdef CLP_AUXILIARY_MODEL
    162     if (!auxiliaryModel_) {
    163 #endif
    164       if (rowScale_) {
    165         if (iSequence<numberColumns_)
    166           scaleFactor = 1.0/(objectiveScale_*columnScale_[iSequence]);
    167         else
    168           scaleFactor = rowScale_[iSequence-numberColumns_]/objectiveScale_;
    169       } else {
    170         scaleFactor = 1.0/objectiveScale_;
    171       }
    172 #ifdef CLP_AUXILIARY_MODEL
     161    if (rowScale_) {
     162      if (iSequence<numberColumns_)
     163        scaleFactor = 1.0/(objectiveScale_*columnScale_[iSequence]);
     164      else
     165        scaleFactor = rowScale_[iSequence-numberColumns_]/objectiveScale_;
    173166    } else {
    174       if (auxiliaryModel_->rowScale()) {
    175         if (iSequence<numberColumns_)
    176           scaleFactor = 1.0/(objectiveScale_*auxiliaryModel_->columnScale()[iSequence]);
    177         else
    178           scaleFactor = auxiliaryModel_->rowScale()[iSequence-numberColumns_]/objectiveScale_;
    179       } else {
    180         scaleFactor = 1.0/objectiveScale_;
    181       }
    182     }
    183 #endif
     167      scaleFactor = 1.0/objectiveScale_;
     168    }
    184169    if (costIncrease<1.0e30)
    185170      costIncrease *= scaleFactor;
     
    11691154  matrix_->transposeTimes(-1.0,dual_,reducedCost_);
    11701155  checkSolutionInternal();
     1156  if (sumDualInfeasibilities_>1.0e-5||sumPrimalInfeasibilities_>1.0e-5) {
     1157    returnCode=1;
     1158#ifdef CLP_INVESTIGATE
     1159    printf("There are %d dual infeasibilities summing to %g ",
     1160           numberDualInfeasibilities_,sumDualInfeasibilities_);
     1161    printf("and %d primal infeasibilities summing to %g\n",
     1162           numberPrimalInfeasibilities_,sumPrimalInfeasibilities_);
     1163#endif
     1164  }
    11711165  // Below will go to ..DEBUG later
    1172 #if 0 //ndef NDEBUG
     1166#if 1 //ndef NDEBUG
    11731167  // Check if correct
    11741168  double * columnActivity = CoinCopyOfArray(columnActivity_,numberColumns_);
     
    11761170  double * reducedCost = CoinCopyOfArray(reducedCost_,numberColumns_);
    11771171  double * dual = CoinCopyOfArray(dual_,numberRows_);
    1178   this->primal();
     1172  this->dual(); //primal();
    11791173  CoinRelFltEq eq(1.0e-5);
    11801174  for (iRow=0;iRow<numberRows_;iRow++) {
     
    15571551    }
    15581552  }
     1553#if 0
     1554  if (small) {
     1555    static int which=0;
     1556    which++;
     1557    char xxxx[20];
     1558    sprintf(xxxx,"bad%d.mps",which);
     1559    small->writeMps(xxxx,0,1);
     1560    sprintf(xxxx,"largebad%d.mps",which);
     1561    writeMps(xxxx,0,1);
     1562    printf("bad%d %x old size %d %d new %d %d\n",which,small,
     1563           numberRows_,numberColumns_,small->numberRows(),small->numberColumns());
     1564#if 0
     1565    for (int i=0;i<numberColumns_;i++)
     1566      printf("Bound %d %g %g\n",i,columnLower_[i],columnUpper_[i]);
     1567    for (int i=0;i<numberRows_;i++)
     1568      printf("Row bound %d %g %g\n",i,rowLower_[i],rowUpper_[i]);
     1569#endif
     1570  }
     1571#endif
    15591572  return small;
    15601573}
     
    15671580                             const int * whichColumn, int nBound)
    15681581{
     1582#ifndef NDEBUG
     1583  for (int i=0;i<small.numberRows();i++)
     1584    assert (whichRow[i]>=0&&whichRow[i]<numberRows_);
     1585  for (int i=0;i<small.numberColumns();i++)
     1586    assert (whichColumn[i]>=0&&whichColumn[i]<numberColumns_);
     1587#endif
    15691588  getbackSolution(small,whichRow,whichColumn);
    15701589  // and deal with status for bounds
     
    32403259  return nelCreate;
    32413260}
     3261// Quick try at cleaning up duals if postsolve gets wrong
     3262void
     3263ClpSimplexOther::cleanupAfterPostsolve()
     3264{
     3265  // First mark singleton equality rows
     3266  char * mark = new char [ numberRows_];
     3267  memset(mark,0,numberRows_);
     3268  const int * row = matrix_->getIndices();
     3269  const CoinBigIndex * columnStart = matrix_->getVectorStarts();
     3270  const int * columnLength = matrix_->getVectorLengths();
     3271  const double * element = matrix_->getElements();
     3272  for (int iColumn=0;iColumn<numberColumns_;iColumn++) {
     3273    for (CoinBigIndex j =columnStart[iColumn];
     3274         j<columnStart[iColumn]+columnLength[iColumn];j++) {
     3275      int iRow = row[j];
     3276      if (mark[iRow])
     3277        mark[iRow]=2;
     3278      else
     3279        mark[iRow]=1;
     3280    }
     3281  }
     3282  // for now just == rows
     3283  for (int iRow=0;iRow<numberRows_;iRow++) {
     3284    if (rowUpper_[iRow]>rowLower_[iRow])
     3285      mark[iRow]=3;
     3286  }
     3287  double dualTolerance=dblParam_[ClpDualTolerance];
     3288  double primalTolerance=dblParam_[ClpPrimalTolerance];
     3289  int numberCleaned=0;
     3290  double maxmin = optimizationDirection_;
     3291  for (int iColumn=0;iColumn<numberColumns_;iColumn++) {
     3292    double dualValue = reducedCost_[iColumn]*maxmin;
     3293    double primalValue = columnActivity_[iColumn];
     3294    double lower = columnLower_[iColumn];
     3295    double upper = columnUpper_[iColumn];
     3296    int way=0;
     3297    switch(getColumnStatus(iColumn)) {
     3298       
     3299    case basic:
     3300      // dual should be zero
     3301      if (dualValue>dualTolerance) {
     3302        way=-1;
     3303      } else if (dualValue<-dualTolerance) {
     3304        way=1;
     3305      }
     3306      break;
     3307    case ClpSimplex::isFixed:
     3308      break;
     3309    case atUpperBound:
     3310      // dual should not be positive
     3311      if (dualValue>dualTolerance) {
     3312        way=-1;
     3313      }
     3314      break;
     3315    case atLowerBound:
     3316      // dual should not be negative
     3317      if (dualValue<-dualTolerance) {
     3318        way=1;
     3319      }
     3320      break;
     3321    case superBasic:
     3322    case isFree:
     3323      if (primalValue<upper-primalTolerance) {
     3324        // dual should not be negative
     3325        if (dualValue<-dualTolerance) {
     3326          way=1;
     3327        }
     3328      }
     3329      if (primalValue>lower+primalTolerance) {
     3330        // dual should not be positive
     3331        if (dualValue>dualTolerance) {
     3332          way=-1;
     3333        }
     3334      }
     3335      break;
     3336    }
     3337    if (way) {
     3338      // see if can find singleton row
     3339      for (CoinBigIndex j =columnStart[iColumn];
     3340           j<columnStart[iColumn]+columnLength[iColumn];j++) {
     3341        int iRow = row[j];
     3342        if (mark[iRow]==1) {
     3343          double value=element[j];
     3344          // dj - addDual*value == 0.0
     3345          double addDual = dualValue/value;
     3346          dual_[iRow] += addDual;
     3347          reducedCost_[iColumn]=0.0;
     3348          numberCleaned++;
     3349          break;
     3350        }
     3351      }
     3352    }
     3353  }
     3354  delete [] mark;
     3355#ifdef CLP_INVESTIGATE
     3356  printf("cleanupAfterPostsolve cleaned up %d columns\n",numberCleaned);
     3357#endif
     3358  // Redo
     3359  memcpy(reducedCost_,this->objective(),numberColumns_*sizeof(double));
     3360  matrix_->transposeTimes(-1.0,dual_,reducedCost_);
     3361  checkSolutionInternal();
     3362}
  • trunk/Clp/src/ClpSimplexOther.hpp

    r1128 r1344  
    184184                   const int * whichRows, const int * whichColumns,
    185185                   int nBound);
     186  /// Quick try at cleaning up duals if postsolve gets wrong
     187  void cleanupAfterPostsolve();
    186188  /** Tightens integer bounds - returns number tightened or -1 if infeasible
    187189  */
  • trunk/Clp/src/ClpSimplexPrimal.cpp

    r1321 r1344  
    232232    // This says whether to restore things etc
    233233    int factorType=0;
    234     if (problemStatus_<0&&perturbation_<100) {
     234    if (problemStatus_<0&&perturbation_<100&&!ifValuesPass) {
    235235      perturb(0);
    236236      // Can't get here if values pass
     
    448448            break;
    449449          }
     450          //#define FEB_TRY
     451#ifdef FEB_TRY
     452          if (perturbation_<100)
     453            perturb(0);
     454#endif
    450455        }
    451456      }
     
    19091914                               CoinIndexedVector * spareColumn2)
    19101915{
     1916 
     1917  ClpMatrixBase * saveMatrix = matrix_;
     1918  double * saveRowScale = rowScale_;
     1919  if (scaledMatrix_) {
     1920    rowScale_=NULL;
     1921    matrix_ = scaledMatrix_;
     1922  }
    19111923  sequenceIn_ = primalColumnPivot_->pivotColumn(updates,spareRow1,
    19121924                                               spareRow2,spareColumn1,
    19131925                                               spareColumn2);
     1926  if (scaledMatrix_) {
     1927    matrix_ = saveMatrix;
     1928    rowScale_ = saveRowScale;
     1929  }
    19141930  if (sequenceIn_>=0) {
    19151931    valueIn_=solution_[sequenceIn_];
     
    22242240    else
    22252241      perturbation = 1.0e-1;
     2242    if (perturbation_>50&&perturbation_<60) {
     2243      // reduce
     2244      while (perturbation_>50) {
     2245        perturbation_--;
     2246        perturbation *= 0.25;
     2247      }
     2248    }
    22262249  } else if (perturbation_<100) {
    22272250    perturbation = pow(10.0,perturbation_);
     
    26612684    if (largestDualError_>1.0e-5)
    26622685      checkValue=1.0e-1;
    2663     if (!ifValuesPass&&solveType_==1&&(saveDj*dualIn_<1.0e-20||
     2686    double test2 = dualTolerance_;
     2687    double test1 = 1.0e-20;
     2688#if 0 //def FEB_TRY
     2689    if (factorization_->pivots()<1) {
     2690      test1 = -1.0e-4;
     2691      if ((saveDj<0.0&&dualIn_<-1.0e-5*dualTolerance_)||
     2692          (saveDj>0.0&&dualIn_>1.0e-5*dualTolerance_))
     2693        test2=0.0; // allow through
     2694    }
     2695#endif
     2696    if (!ifValuesPass&&solveType_==1&&(saveDj*dualIn_<test1||
    26642697        fabs(saveDj-dualIn_)>checkValue*(1.0+fabs(saveDj))||
    2665                         fabs(dualIn_)<dualTolerance_)) {
     2698                        fabs(dualIn_)<test2)) {
    26662699      char x = isColumn(sequenceIn_) ? 'C' :'R';
    26672700      handler_->message(CLP_PRIMAL_DJ,messages_)
     
    26822715      } else {
    26832716        // take on more relaxed criterion
    2684         if (saveDj*dualIn_<1.0e-20||
     2717        if (saveDj*dualIn_<test1||
    26852718            fabs(saveDj-dualIn_)>2.0e-1*(1.0+fabs(dualIn_))||
    2686             fabs(dualIn_)<dualTolerance_) {
     2719            fabs(dualIn_)<test2) {
    26872720          // need to reject something
    26882721          char x = isColumn(sequenceIn_) ? 'C' :'R';
     
    26912724            <<CoinMessageEol;
    26922725          setFlagged(sequenceIn_);
     2726#ifdef FEB_TRY
     2727          // Make safer?
     2728          factorization_->saferTolerances (1.0e-15,-1.03);
     2729#endif
    26932730          progress_.clearBadTimes();
    26942731          lastBadIteration_ = numberIterations_; // say be more cautious
  • trunk/Clp/src/ClpSolve.cpp

    r1343 r1344  
    4444#ifdef TAUCS_BARRIER
    4545#include "ClpCholeskyTaucs.hpp"
     46#define FAST_BARRIER
     47#endif
     48#ifdef MUMPS_BARRIER
     49#include "ClpCholeskyMumps.hpp"
    4650#define FAST_BARRIER
    4751#endif
     
    19551959      {
    19561960        ClpCholeskyTaucs * cholesky = new ClpCholeskyTaucs();
     1961        barrier.setCholesky(cholesky);
     1962        assert (!doKKT);
     1963      }
     1964      break;
     1965#endif
     1966#ifdef MUMPS_BARRIER
     1967    case 6:
     1968      {
     1969        ClpCholeskyMumps * cholesky = new ClpCholeskyMumps();
    19571970        barrier.setCholesky(cholesky);
    19581971        assert (!doKKT);
  • trunk/Clp/src/IdiSolve.cpp

    r1321 r1344  
    7070Idiot::objval(int nrows, int ncols, double * rowsol , double * colsol,
    7171              double * pi, double * djs, const double * cost ,
    72