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

formatting

File:
1 edited

Legend:

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

    <
    r2274 r2385  
    3535double factorization_instrument(int type)
    3636{
    37      static int times[10] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
    38      static double startTime = 0.0;
    39      static double totalTimes [10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
    40      if (type < 0) {
    41           assert (!startTime);
    42           startTime = CoinCpuTime();
    43           return 0.0;
    44      } else if (type > 0) {
    45           times[type]++;
    46           double difference = CoinCpuTime() - startTime;
    47           totalTimes[type] += difference;
    48           startTime = 0.0;
    49           return difference;
    50      } else {
    51           // report
    52           const char * types[10] = {
    53                "", "fac=rhs_etc", "factorize", "replace", "update_FT",
    54                "update", "update_transpose", "gosparse", "getWeights!", "update2_FT"
    55           };
    56           double total = 0.0;
    57           for (int i = 1; i < 10; i++) {
    58                if (times[i]) {
    59                     printf("%s was called %d times taking %g seconds\n",
    60                            types[i], times[i], totalTimes[i]);
    61                     total += totalTimes[i];
    62                     times[i] = 0;
    63                     totalTimes[i] = 0.0;
    64                }
    65           }
    66           return total;
    67      }
     37  static int times[10] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 };
     38  static double startTime = 0.0;
     39  static double totalTimes[10] = { 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 };
     40  if (type < 0) {
     41    assert(!startTime);
     42    startTime = CoinCpuTime();
     43    return 0.0;
     44  } else if (type > 0) {
     45    times[type]++;
     46    double difference = CoinCpuTime() - startTime;
     47    totalTimes[type] += difference;
     48    startTime = 0.0;
     49    return difference;
     50  } else {
     51    // report
     52    const char *types[10] = {
     53      "", "fac=rhs_etc", "factorize", "replace", "update_FT",
     54      "update", "update_transpose", "gosparse", "getWeights!", "update2_FT"
     55    };
     56    double total = 0.0;
     57    for (int i = 1; i < 10; i++) {
     58      if (times[i]) {
     59        printf("%s was called %d times taking %g seconds\n",
     60          types[i], times[i], totalTimes[i]);
     61        total += totalTimes[i];
     62        times[i] = 0;
     63        totalTimes[i] = 0.0;
     64      }
     65    }
     66    return total;
     67  }
    6868}
    6969#endif
     
    7676// Default Constructor
    7777//-------------------------------------------------------------------
    78 ClpFactorization::ClpFactorization () :
    79     CoinFactorization()
    80 {
    81 #ifndef SLIM_CLP
    82      networkBasis_ = NULL;
     78ClpFactorization::ClpFactorization()
     79  : CoinFactorization()
     80{
     81#ifndef SLIM_CLP
     82  networkBasis_ = NULL;
    8383#endif
    8484}
     
    8787// Copy constructor
    8888//-------------------------------------------------------------------
    89 ClpFactorization::ClpFactorization (const ClpFactorization & rhs,
    90                                     int dummyDenseIfSmaller) :
    91     CoinFactorization(rhs)
    92 {
    93 #ifndef SLIM_CLP
    94      if (rhs.networkBasis_)
    95           networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
    96      else
    97           networkBasis_ = NULL;
     89ClpFactorization::ClpFactorization(const ClpFactorization &rhs,
     90  int dummyDenseIfSmaller)
     91  : CoinFactorization(rhs)
     92{
     93#ifndef SLIM_CLP
     94  if (rhs.networkBasis_)
     95    networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
     96  else
     97    networkBasis_ = NULL;
    9898#endif
    9999}
    100100
    101 ClpFactorization::ClpFactorization (const CoinFactorization & rhs) :
    102     CoinFactorization(rhs)
    103 {
    104 #ifndef SLIM_CLP
    105      networkBasis_ = NULL;
     101ClpFactorization::ClpFactorization(const CoinFactorization &rhs)
     102  : CoinFactorization(rhs)
     103{
     104#ifndef SLIM_CLP
     105  networkBasis_ = NULL;
    106106#endif
    107107}
     
    110110// Destructor
    111111//-------------------------------------------------------------------
    112 ClpFactorization::~ClpFactorization ()
    113 {
    114 #ifndef SLIM_CLP
    115      delete networkBasis_;
     112ClpFactorization::~ClpFactorization()
     113{
     114#ifndef SLIM_CLP
     115  delete networkBasis_;
    116116#endif
    117117}
     
    121121//-------------------------------------------------------------------
    122122ClpFactorization &
    123 ClpFactorization::operator=(const ClpFactorization& rhs)
    124 {
    125 #ifdef CLP_FACTORIZATION_INSTRUMENT
    126      factorization_instrument(-1);
    127 #endif
    128      if (this != &rhs) {
    129           CoinFactorization::operator=(rhs);
    130 #ifndef SLIM_CLP
    131           delete networkBasis_;
    132           if (rhs.networkBasis_)
    133                networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
    134           else
    135                networkBasis_ = NULL;
    136 #endif
    137      }
    138 #ifdef CLP_FACTORIZATION_INSTRUMENT
    139      factorization_instrument(1);
    140 #endif
    141      return *this;
     123ClpFactorization::operator=(const ClpFactorization &rhs)
     124{
     125#ifdef CLP_FACTORIZATION_INSTRUMENT
     126  factorization_instrument(-1);
     127#endif
     128  if (this != &rhs) {
     129    CoinFactorization::operator=(rhs);
     130#ifndef SLIM_CLP
     131    delete networkBasis_;
     132    if (rhs.networkBasis_)
     133      networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
     134    else
     135      networkBasis_ = NULL;
     136#endif
     137  }
     138#ifdef CLP_FACTORIZATION_INSTRUMENT
     139  factorization_instrument(1);
     140#endif
     141  return *this;
    142142}
    143143#if 0
     
    155155}
    156156#endif
    157 int
    158 ClpFactorization::factorize ( ClpSimplex * model,
    159                               int solveType, bool valuesPass)
    160 {
    161 #ifdef CLP_FACTORIZATION_INSTRUMENT
    162      factorization_instrument(-1);
    163 #endif
    164      ClpMatrixBase * matrix = model->clpMatrix();
    165      int numberRows = model->numberRows();
    166      int numberColumns = model->numberColumns();
    167      if (!numberRows)
    168           return 0;
    169      // If too many compressions increase area
    170      if (numberPivots_ > 1 && numberCompressions_ * 10 > numberPivots_ + 10) {
    171           areaFactor_ *= 1.1;
    172      }
    173      //int numberPivots=numberPivots_;
     157int ClpFactorization::factorize(ClpSimplex *model,
     158  int solveType, bool valuesPass)
     159{
     160#ifdef CLP_FACTORIZATION_INSTRUMENT
     161  factorization_instrument(-1);
     162#endif
     163  ClpMatrixBase *matrix = model->clpMatrix();
     164  int numberRows = model->numberRows();
     165  int numberColumns = model->numberColumns();
     166  if (!numberRows)
     167    return 0;
     168  // If too many compressions increase area
     169  if (numberPivots_ > 1 && numberCompressions_ * 10 > numberPivots_ + 10) {
     170    areaFactor_ *= 1.1;
     171  }
     172  //int numberPivots=numberPivots_;
    174173#if 0
    175174     if (model->algorithm() > 0)
     
    177176#endif
    178177#ifndef SLIM_CLP
    179      if (!networkBasis_ || doCheck) {
    180 #endif
    181           status_ = -99;
    182           int * pivotVariable = model->pivotVariable();
    183           //returns 0 -okay, -1 singular, -2 too many in basis, -99 memory */
    184           while (status_ < -98) {
     178  if (!networkBasis_ || doCheck) {
     179#endif
     180    status_ = -99;
     181    int *pivotVariable = model->pivotVariable();
     182    //returns 0 -okay, -1 singular, -2 too many in basis, -99 memory */
     183    while (status_ < -98) {
    185184
    186                int i;
    187                int numberBasic = 0;
    188                int numberRowBasic;
    189                // Move pivot variables across if they look good
    190                int * pivotTemp = model->rowArray(0)->getIndices();
    191                assert (!model->rowArray(0)->getNumElements());
    192                if (!matrix->rhsOffset(model)) {
     185      int i;
     186      int numberBasic = 0;
     187      int numberRowBasic;
     188      // Move pivot variables across if they look good
     189      int *pivotTemp = model->rowArray(0)->getIndices();
     190      assert(!model->rowArray(0)->getNumElements());
     191      if (!matrix->rhsOffset(model)) {
    193192#if 0
    194193                    if (numberSave > 0) {
     
    259258                    }
    260259#endif
    261                     // Seems to prefer things in order so quickest
    262                     // way is to go though like this
    263                     for (i = 0; i < numberRows; i++) {
    264                          if (model->getRowStatus(i) == ClpSimplex::basic)
    265                               pivotTemp[numberBasic++] = i;
    266                     }
    267                     numberRowBasic = numberBasic;
    268                     /* Put column basic variables into pivotVariable
     260        // Seems to prefer things in order so quickest
     261        // way is to go though like this
     262        for (i = 0; i < numberRows; i++) {
     263          if (model->getRowStatus(i) == ClpSimplex::basic)
     264            pivotTemp[numberBasic++] = i;
     265        }
     266        numberRowBasic = numberBasic;
     267        /* Put column basic variables into pivotVariable
    269268                       This is done by ClpMatrixBase to allow override for gub
    270269                    */
    271                     matrix->generalExpanded(model, 0, numberBasic);
    272                } else {
    273                     // Long matrix - do a different way
    274                     bool fullSearch = false;
    275                     for (i = 0; i < numberRows; i++) {
    276                          int iPivot = pivotVariable[i];
    277                          if (iPivot >= numberColumns) {
    278                               pivotTemp[numberBasic++] = iPivot - numberColumns;
    279                          }
    280                     }
    281                     numberRowBasic = numberBasic;
    282                     for (i = 0; i < numberRows; i++) {
    283                          int iPivot = pivotVariable[i];
    284                          if (iPivot < numberColumns) {
    285                               if (iPivot >= 0) {
    286                                    pivotTemp[numberBasic++] = iPivot;
    287                               } else {
    288                                    // not full basis
    289                                    fullSearch = true;
    290                                    break;
    291                               }
    292                          }
    293                     }
    294                     if (fullSearch) {
    295                          // do slow way
    296                          numberBasic = 0;
    297                          for (i = 0; i < numberRows; i++) {
    298                               if (model->getRowStatus(i) == ClpSimplex::basic)
    299                                    pivotTemp[numberBasic++] = i;
    300                          }
    301                          numberRowBasic = numberBasic;
    302                          /* Put column basic variables into pivotVariable
     270        matrix->generalExpanded(model, 0, numberBasic);
     271      } else {
     272        // Long matrix - do a different way
     273        bool fullSearch = false;
     274        for (i = 0; i < numberRows; i++) {
     275          int iPivot = pivotVariable[i];
     276          if (iPivot >= numberColumns) {
     277            pivotTemp[numberBasic++] = iPivot - numberColumns;
     278          }
     279        }
     280        numberRowBasic = numberBasic;
     281        for (i = 0; i < numberRows; i++) {
     282          int iPivot = pivotVariable[i];
     283          if (iPivot < numberColumns) {
     284            if (iPivot >= 0) {
     285              pivotTemp[numberBasic++] = iPivot;
     286            } else {
     287              // not full basis
     288              fullSearch = true;
     289              break;
     290            }
     291          }
     292        }
     293        if (fullSearch) {
     294          // do slow way
     295          numberBasic = 0;
     296          for (i = 0; i < numberRows; i++) {
     297            if (model->getRowStatus(i) == ClpSimplex::basic)
     298              pivotTemp[numberBasic++] = i;
     299          }
     300          numberRowBasic = numberBasic;
     301          /* Put column basic variables into pivotVariable
    303302                            This is done by ClpMatrixBase to allow override for gub
    304303                         */
    305                          matrix->generalExpanded(model, 0, numberBasic);
    306                     }
    307                }
    308                if (numberBasic > model->maximumBasic()) {
     304          matrix->generalExpanded(model, 0, numberBasic);
     305        }
     306      }
     307      if (numberBasic > model->maximumBasic()) {
    309308#if 0 // ndef NDEBUG
    310309                    printf("%d basic - should only be %d\n",
    311310                           numberBasic, numberRows);
    312311#endif
    313                     // Take out some
    314                     numberBasic = numberRowBasic;
    315                     for (int i = 0; i < numberColumns; i++) {
    316                          if (model->getColumnStatus(i) == ClpSimplex::basic) {
    317                               if (numberBasic < numberRows)
    318                                    numberBasic++;
    319                               else
    320                                    model->setColumnStatus(i, ClpSimplex::superBasic);
    321                          }
    322                     }
    323                     numberBasic = numberRowBasic;
    324                     matrix->generalExpanded(model, 0, numberBasic);
    325                }
    326 #ifndef SLIM_CLP
    327                // see if matrix a network
     312        // Take out some
     313        numberBasic = numberRowBasic;
     314        for (int i = 0; i < numberColumns; i++) {
     315          if (model->getColumnStatus(i) == ClpSimplex::basic) {
     316            if (numberBasic < numberRows)
     317              numberBasic++;
     318            else
     319              model->setColumnStatus(i, ClpSimplex::superBasic);
     320          }
     321        }
     322        numberBasic = numberRowBasic;
     323        matrix->generalExpanded(model, 0, numberBasic);
     324      }
     325#ifndef SLIM_CLP
     326      // see if matrix a network
    328327#ifndef NO_RTTI
    329                ClpNetworkMatrix* networkMatrix =
    330                     dynamic_cast< ClpNetworkMatrix*>(model->clpMatrix());
     328      ClpNetworkMatrix *networkMatrix = dynamic_cast< ClpNetworkMatrix * >(model->clpMatrix());
    331329#else
    332                ClpNetworkMatrix* networkMatrix = NULL;
    333                if (model->clpMatrix()->type() == 11)
    334                     networkMatrix =
    335                          static_cast< ClpNetworkMatrix*>(model->clpMatrix());
    336 #endif
    337                // If network - still allow ordinary factorization first time for laziness
    338                if (networkMatrix)
    339                     biasLU_ = 0; // All to U if network
    340                //int saveMaximumPivots = maximumPivots();
    341                delete networkBasis_;
    342                networkBasis_ = NULL;
    343                if (networkMatrix && !doCheck)
    344                     maximumPivots(1);
    345 #endif
    346                //printf("L, U, R %d %d %d\n",numberElementsL(),numberElementsU(),numberElementsR());
    347                while (status_ == -99) {
    348                     // maybe for speed will be better to leave as many regions as possible
    349                     gutsOfDestructor();
    350                     gutsOfInitialize(2);
    351                     CoinBigIndex numberElements = numberRowBasic;
     330      ClpNetworkMatrix *networkMatrix = NULL;
     331      if (model->clpMatrix()->type() == 11)
     332        networkMatrix = static_cast< ClpNetworkMatrix * >(model->clpMatrix());
     333#endif
     334      // If network - still allow ordinary factorization first time for laziness
     335      if (networkMatrix)
     336        biasLU_ = 0; // All to U if network
     337      //int saveMaximumPivots = maximumPivots();
     338      delete networkBasis_;
     339      networkBasis_ = NULL;
     340      if (networkMatrix && !doCheck)
     341        maximumPivots(1);
     342#endif
     343      //printf("L, U, R %d %d %d\n",numberElementsL(),numberElementsU(),numberElementsR());
     344      while (status_ == -99) {
     345        // maybe for speed will be better to leave as many regions as possible
     346        gutsOfDestructor();
     347        gutsOfInitialize(2);
     348        CoinBigIndex numberElements = numberRowBasic;
    352349
    353                     // compute how much in basis
     350        // compute how much in basis
    354351
    355                     int i;
    356                     // can change for gub
    357                     int numberColumnBasic = numberBasic - numberRowBasic;
     352        int i;
     353        // can change for gub
     354        int numberColumnBasic = numberBasic - numberRowBasic;
    358355
    359                     numberElements += matrix->countBasis(model,
    360                                                          pivotTemp + numberRowBasic,
    361                                                          numberRowBasic,
    362                                                          numberColumnBasic);
    363                     // and recompute as network side say different
    364                     if (model->numberIterations())
    365                          numberRowBasic = numberBasic - numberColumnBasic;
    366                     numberElements = 3 * numberBasic + 3 * numberElements + 20000;
     356        numberElements += matrix->countBasis(model,
     357          pivotTemp + numberRowBasic,
     358          numberRowBasic,
     359          numberColumnBasic);
     360        // and recompute as network side say different
     361        if (model->numberIterations())
     362          numberRowBasic = numberBasic - numberColumnBasic;
     363        numberElements = 3 * numberBasic + 3 * numberElements + 20000;
    367364#if 0
    368365                    // If iteration not zero then may be compressed
     
    371368                               2 * numberElements );
    372369#else
    373                     getAreas ( numberRows,
    374                                numberRowBasic + numberColumnBasic, numberElements,
    375                                2 * numberElements );
    376 #endif
    377                     //fill
    378                     // Fill in counts so we can skip part of preProcess
    379                     int * numberInRow = numberInRow_.array();
    380                     int * numberInColumn = numberInColumn_.array();
    381                     CoinZeroN ( numberInRow, numberRows_ + 1 );
    382                     CoinZeroN ( numberInColumn, maximumColumnsExtra_ + 1 );
    383                     double * elementU = elementU_.array();
    384                     int * indexRowU = indexRowU_.array();
    385                     CoinBigIndex * startColumnU = startColumnU_.array();
    386                     for (i = 0; i < numberRowBasic; i++) {
    387                          int iRow = pivotTemp[i];
    388                          indexRowU[i] = iRow;
    389                          startColumnU[i] = i;
    390                          elementU[i] = slackValue_;
    391                          numberInRow[iRow] = 1;
    392                          numberInColumn[i] = 1;
    393                     }
    394                     startColumnU[numberRowBasic] = numberRowBasic;
    395                     // can change for gub so redo
    396                     numberColumnBasic = numberBasic - numberRowBasic;
    397                     matrix->fillBasis(model,
    398                                       pivotTemp + numberRowBasic,
    399                                       numberColumnBasic,
    400                                       indexRowU,
    401                                       startColumnU + numberRowBasic,
    402                                       numberInRow,
    403                                       numberInColumn + numberRowBasic,
    404                                       elementU);
     370      getAreas(numberRows,
     371        numberRowBasic + numberColumnBasic, numberElements,
     372        2 * numberElements);
     373#endif
     374        //fill
     375        // Fill in counts so we can skip part of preProcess
     376        int *numberInRow = numberInRow_.array();
     377        int *numberInColumn = numberInColumn_.array();
     378        CoinZeroN(numberInRow, numberRows_ + 1);
     379        CoinZeroN(numberInColumn, maximumColumnsExtra_ + 1);
     380        double *elementU = elementU_.array();
     381        int *indexRowU = indexRowU_.array();
     382        CoinBigIndex *startColumnU = startColumnU_.array();
     383        for (i = 0; i < numberRowBasic; i++) {
     384          int iRow = pivotTemp[i];
     385          indexRowU[i] = iRow;
     386          startColumnU[i] = i;
     387          elementU[i] = slackValue_;
     388          numberInRow[iRow] = 1;
     389          numberInColumn[i] = 1;
     390        }
     391        startColumnU[numberRowBasic] = numberRowBasic;
     392        // can change for gub so redo
     393        numberColumnBasic = numberBasic - numberRowBasic;
     394        matrix->fillBasis(model,
     395          pivotTemp + numberRowBasic,
     396          numberColumnBasic,
     397          indexRowU,
     398          startColumnU + numberRowBasic,
     399          numberInRow,
     400          numberInColumn + numberRowBasic,
     401          elementU);
    405402#if 0
    406403                    {
     
    411408                    }
    412409#endif
    413                     // recompute number basic
    414                     numberBasic = numberRowBasic + numberColumnBasic;
    415                     if (numberBasic)
    416                          numberElements = startColumnU[numberBasic-1]
    417                                           + numberInColumn[numberBasic-1];
    418                     else
    419                          numberElements = 0;
    420                     lengthU_ = numberElements;
    421                     //saveFactorization("dump.d");
    422                     if (biasLU_ >= 3 || numberRows_ != numberColumns_)
    423                          preProcess ( 2 );
    424                     else
    425                          preProcess ( 3 ); // no row copy
    426                     factor (  );
    427                     if (status_ == -99) {
    428                          // get more memory
    429                          areaFactor(2.0 * areaFactor());
    430                     } else if (status_ == -1 && model->numberIterations() == 0 &&
    431                                denseThreshold_) {
    432                          // Round again without dense
    433                          denseThreshold_ = 0;
    434                          status_ = -99;
    435                     }
    436                }
    437                // If we get here status is 0 or -1
    438                if (status_ == 0) {
    439                     // We may need to tamper with order and redo - e.g. network with side
    440                     int useNumberRows = numberRows;
    441                     // **** we will also need to add test in dual steepest to do
    442                     // as we do for network
    443                     matrix->generalExpanded(model, 12, useNumberRows);
    444                     const int * permuteBack = permuteBack_.array();
    445                     const int * back = pivotColumnBack_.array();
    446                     //int * pivotTemp = pivotColumn_.array();
    447                     //ClpDisjointCopyN ( pivotVariable, numberRows , pivotTemp  );
    448                     // Redo pivot order
    449                     for (i = 0; i < numberRowBasic; i++) {
    450                          int k = pivotTemp[i];
    451                          // so rowIsBasic[k] would be permuteBack[back[i]]
    452                          pivotVariable[permuteBack[back[i]]] = k + numberColumns;
    453                     }
    454                     for (; i < useNumberRows; i++) {
    455                          int k = pivotTemp[i];
    456                          // so rowIsBasic[k] would be permuteBack[back[i]]
    457                          pivotVariable[permuteBack[back[i]]] = k;
    458                     }
     410        // recompute number basic
     411        numberBasic = numberRowBasic + numberColumnBasic;
     412        if (numberBasic)
     413          numberElements = startColumnU[numberBasic - 1]
     414            + numberInColumn[numberBasic - 1];
     415        else
     416          numberElements = 0;
     417        lengthU_ = numberElements;
     418        //saveFactorization("dump.d");
     419        if (biasLU_ >= 3 || numberRows_ != numberColumns_)
     420          preProcess(2);
     421        else
     422          preProcess(3); // no row copy
     423        factor();
     424        if (status_ == -99) {
     425          // get more memory
     426          areaFactor(2.0 * areaFactor());
     427        } else if (status_ == -1 && model->numberIterations() == 0 && denseThreshold_) {
     428          // Round again without dense
     429          denseThreshold_ = 0;
     430          status_ = -99;
     431        }
     432      }
     433      // If we get here status is 0 or -1
     434      if (status_ == 0) {
     435        // We may need to tamper with order and redo - e.g. network with side
     436        int useNumberRows = numberRows;
     437        // **** we will also need to add test in dual steepest to do
     438        // as we do for network
     439        matrix->generalExpanded(model, 12, useNumberRows);
     440        const int *permuteBack = permuteBack_.array();
     441        const int *back = pivotColumnBack_.array();
     442        //int * pivotTemp = pivotColumn_.array();
     443        //ClpDisjointCopyN ( pivotVariable, numberRows , pivotTemp  );
     444        // Redo pivot order
     445        for (i = 0; i < numberRowBasic; i++) {
     446          int k = pivotTemp[i];
     447          // so rowIsBasic[k] would be permuteBack[back[i]]
     448          pivotVariable[permuteBack[back[i]]] = k + numberColumns;
     449        }
     450        for (; i < useNumberRows; i++) {
     451          int k = pivotTemp[i];
     452          // so rowIsBasic[k] would be permuteBack[back[i]]
     453          pivotVariable[permuteBack[back[i]]] = k;
     454        }
    459455#if 0
    460456                    if (numberSave >= 0) {
     
    467463                    }
    468464#endif
    469                     // Set up permutation vector
    470                     // these arrays start off as copies of permute
    471                     // (and we could use permute_ instead of pivotColumn (not back though))
    472                     ClpDisjointCopyN ( permute_.array(), useNumberRows , pivotColumn_.array()  );
    473                     ClpDisjointCopyN ( permuteBack_.array(), useNumberRows , pivotColumnBack_.array()  );
    474 #ifndef SLIM_CLP
    475                     if (networkMatrix) {
    476                          maximumPivots(CoinMax(2000, maximumPivots()));
    477                          // redo arrays
    478                          for (int iRow = 0; iRow < 4; iRow++) {
    479                               int length = model->numberRows() + maximumPivots();
    480                               if (iRow == 3 || model->objectiveAsObject()->type() > 1)
    481                                    length += model->numberColumns();
    482                               model->rowArray(iRow)->reserve(length);
    483                          }
    484                          // create network factorization
    485                          if (doCheck)
    486                               delete networkBasis_; // temp
    487                          networkBasis_ = new ClpNetworkBasis(model, numberRows_,
    488                                                              pivotRegion_.array(),
    489                                                              permuteBack_.array(),
    490                                                              startColumnU_.array(),
    491                                                              numberInColumn_.array(),
    492                                                              indexRowU_.array(),
    493                                                              elementU_.array());
    494                          // kill off arrays in ordinary factorization
    495                          if (!doCheck) {
    496                               gutsOfDestructor();
    497                               // but make sure numberRows_ set
    498                               numberRows_ = model->numberRows();
    499                               status_ = 0;
     465        // Set up permutation vector
     466        // these arrays start off as copies of permute
     467        // (and we could use permute_ instead of pivotColumn (not back though))
     468        ClpDisjointCopyN(permute_.array(), useNumberRows, pivotColumn_.array());
     469        ClpDisjointCopyN(permuteBack_.array(), useNumberRows, pivotColumnBack_.array());
     470#ifndef SLIM_CLP
     471        if (networkMatrix) {
     472          maximumPivots(CoinMax(2000, maximumPivots()));
     473          // redo arrays
     474          for (int iRow = 0; iRow < 4; iRow++) {
     475            int length = model->numberRows() + maximumPivots();
     476            if (iRow == 3 || model->objectiveAsObject()->type() > 1)
     477              length += model->numberColumns();
     478            model->rowArray(iRow)->reserve(length);
     479          }
     480          // create network factorization
     481          if (doCheck)
     482            delete networkBasis_; // temp
     483          networkBasis_ = new ClpNetworkBasis(model, numberRows_,
     484            pivotRegion_.array(),
     485            permuteBack_.array(),
     486            startColumnU_.array(),
     487            numberInColumn_.array(),
     488            indexRowU_.array(),
     489            elementU_.array());
     490          // kill off arrays in ordinary factorization
     491          if (!doCheck) {
     492            gutsOfDestructor();
     493            // but make sure numberRows_ set
     494            numberRows_ = model->numberRows();
     495            status_ = 0;
    500496#if 0
    501497                              // but put back permute arrays so odd things will work
     
    509505                              }
    510506#endif
    511                          }
    512                     } else {
    513 #endif
    514                          // See if worth going sparse and when
    515                          if (numberFtranCounts_ > 100) {
    516                               ftranCountInput_ = CoinMax(ftranCountInput_, 1.0);
    517                               ftranAverageAfterL_ = CoinMax(ftranCountAfterL_ / ftranCountInput_, 1.0);
    518                               ftranAverageAfterR_ = CoinMax(ftranCountAfterR_ / ftranCountAfterL_, 1.0);
    519                               ftranAverageAfterU_ = CoinMax(ftranCountAfterU_ / ftranCountAfterR_, 1.0);
    520                               if (btranCountInput_ && btranCountAfterU_ && btranCountAfterR_) {
    521                                    btranAverageAfterU_ = CoinMax(btranCountAfterU_ / btranCountInput_, 1.0);
    522                                    btranAverageAfterR_ = CoinMax(btranCountAfterR_ / btranCountAfterU_, 1.0);
    523                                    btranAverageAfterL_ = CoinMax(btranCountAfterL_ / btranCountAfterR_, 1.0);
    524                               } else {
    525                                    // we have not done any useful btrans (values pass?)
    526                                    btranAverageAfterU_ = 1.0;
    527                                    btranAverageAfterR_ = 1.0;
    528                                    btranAverageAfterL_ = 1.0;
    529                               }
    530                          }
    531                          // scale back
     507          }
     508        } else {
     509#endif
     510          // See if worth going sparse and when
     511          if (numberFtranCounts_ > 100) {
     512            ftranCountInput_ = CoinMax(ftranCountInput_, 1.0);
     513            ftranAverageAfterL_ = CoinMax(ftranCountAfterL_ / ftranCountInput_, 1.0);
     514            ftranAverageAfterR_ = CoinMax(ftranCountAfterR_ / ftranCountAfterL_, 1.0);
     515            ftranAverageAfterU_ = CoinMax(ftranCountAfterU_ / ftranCountAfterR_, 1.0);
     516            if (btranCountInput_ && btranCountAfterU_ && btranCountAfterR_) {
     517              btranAverageAfterU_ = CoinMax(btranCountAfterU_ / btranCountInput_, 1.0);
     518              btranAverageAfterR_ = CoinMax(btranCountAfterR_ / btranCountAfterU_, 1.0);
     519              btranAverageAfterL_ = CoinMax(btranCountAfterL_ / btranCountAfterR_, 1.0);
     520            } else {
     521              // we have not done any useful btrans (values pass?)
     522              btranAverageAfterU_ = 1.0;
     523              btranAverageAfterR_ = 1.0;
     524              btranAverageAfterL_ = 1.0;
     525            }
     526          }
     527          // scale back
    532528
    533                          ftranCountInput_ *= 0.8;
    534                          ftranCountAfterL_ *= 0.8;
    535                          ftranCountAfterR_ *= 0.8;
    536                          ftranCountAfterU_ *= 0.8;
    537                          btranCountInput_ *= 0.8;
    538                          btranCountAfterU_ *= 0.8;
    539                          btranCountAfterR_ *= 0.8;
    540                          btranCountAfterL_ *= 0.8;
    541 #ifndef SLIM_CLP
    542                     }
    543 #endif
    544                } else if (status_ == -1 && (solveType == 0 || solveType == 2)) {
    545                     // This needs redoing as it was merged coding - does not need array
    546                     int numberTotal = numberRows + numberColumns;
    547                     int * isBasic = new int [numberTotal];
    548                     int * rowIsBasic = isBasic + numberColumns;
    549                     int * columnIsBasic = isBasic;
    550                     for (i = 0; i < numberTotal; i++)
    551                          isBasic[i] = -1;
    552                     for (i = 0; i < numberRowBasic; i++) {
    553                          int iRow = pivotTemp[i];
    554                          rowIsBasic[iRow] = 1;
    555                     }
    556                     for (; i < numberBasic; i++) {
    557                          int iColumn = pivotTemp[i];
    558                          columnIsBasic[iColumn] = 1;
    559                     }
    560                     numberBasic = 0;
    561                     for (i = 0; i < numberRows; i++)
    562                          pivotVariable[i] = -1;
    563                     // mark as basic or non basic
    564                     const int * pivotColumn = pivotColumn_.array();
    565                     for (i = 0; i < numberRows; i++) {
    566                          if (rowIsBasic[i] >= 0) {
    567                               if (pivotColumn[numberBasic] >= 0) {
    568                                    rowIsBasic[i] = pivotColumn[numberBasic];
    569                               } else {
    570                                    rowIsBasic[i] = -1;
    571                                    model->setRowStatus(i, ClpSimplex::superBasic);
    572                               }
    573                               numberBasic++;
    574                          }
    575                     }
    576                     for (i = 0; i < numberColumns; i++) {
    577                          if (columnIsBasic[i] >= 0) {
    578                               if (pivotColumn[numberBasic] >= 0)
    579                                    columnIsBasic[i] = pivotColumn[numberBasic];
    580                               else
    581                                    columnIsBasic[i] = -1;
    582                               numberBasic++;
    583                          }
    584                     }
    585                     // leave pivotVariable in useful form for cleaning basis
    586                     int * pivotVariable = model->pivotVariable();
    587                     for (i = 0; i < numberRows; i++) {
    588                          pivotVariable[i] = -1;
    589                     }
     529          ftranCountInput_ *= 0.8;
     530          ftranCountAfterL_ *= 0.8;
     531          ftranCountAfterR_ *= 0.8;
     532          ftranCountAfterU_ *= 0.8;
     533          btranCountInput_ *= 0.8;
     534          btranCountAfterU_ *= 0.8;
     535          btranCountAfterR_ *= 0.8;
     536          btranCountAfterL_ *= 0.8;
     537#ifndef SLIM_CLP
     538        }
     539#endif
     540      } else if (status_ == -1 && (solveType == 0 || solveType == 2)) {
     541        // This needs redoing as it was merged coding - does not need array
     542        int numberTotal = numberRows + numberColumns;
     543        int *isBasic = new int[numberTotal];
     544        int *rowIsBasic = isBasic + numberColumns;
     545        int *columnIsBasic = isBasic;
     546        for (i = 0; i < numberTotal; i++)
     547          isBasic[i] = -1;
     548        for (i = 0; i < numberRowBasic; i++) {
     549          int iRow = pivotTemp[i];
     550          rowIsBasic[iRow] = 1;
     551        }
     552        for (; i < numberBasic; i++) {
     553          int iColumn = pivotTemp[i];
     554          columnIsBasic[iColumn] = 1;
     555        }
     556        numberBasic = 0;
     557        for (i = 0; i < numberRows; i++)
     558          pivotVariable[i] = -1;
     559        // mark as basic or non basic
     560        const int *pivotColumn = pivotColumn_.array();
     561        for (i = 0; i < numberRows; i++) {
     562          if (rowIsBasic[i] >= 0) {
     563            if (pivotColumn[numberBasic] >= 0) {
     564              rowIsBasic[i] = pivotColumn[numberBasic];
     565            } else {
     566              rowIsBasic[i] = -1;
     567              model->setRowStatus(i, ClpSimplex::superBasic);
     568            }
     569            numberBasic++;
     570          }
     571        }
     572        for (i = 0; i < numberColumns; i++) {
     573          if (columnIsBasic[i] >= 0) {
     574            if (pivotColumn[numberBasic] >= 0)
     575              columnIsBasic[i] = pivotColumn[numberBasic];
     576            else
     577              columnIsBasic[i] = -1;
     578            numberBasic++;
     579          }
     580        }
     581        // leave pivotVariable in useful form for cleaning basis
     582        int *pivotVariable = model->pivotVariable();
     583        for (i = 0; i < numberRows; i++) {
     584          pivotVariable[i] = -1;
     585        }
    590586
    591                     for (i = 0; i < numberRows; i++) {
    592                          if (model->getRowStatus(i) == ClpSimplex::basic) {
    593                               int iPivot = rowIsBasic[i];
    594                               if (iPivot >= 0)
    595                                    pivotVariable[iPivot] = i + numberColumns;
    596                          }
    597                     }
    598                     for (i = 0; i < numberColumns; i++) {
    599                          if (model->getColumnStatus(i) == ClpSimplex::basic) {
    600                               int iPivot = columnIsBasic[i];
    601                               if (iPivot >= 0)
    602                                    pivotVariable[iPivot] = i;
    603                          }
    604                     }
    605                     delete [] isBasic;
    606                     double * columnLower = model->lowerRegion();
    607                     double * columnUpper = model->upperRegion();
    608                     double * columnActivity = model->solutionRegion();
    609                     double * rowLower = model->lowerRegion(0);
    610                     double * rowUpper = model->upperRegion(0);
    611                     double * rowActivity = model->solutionRegion(0);
    612                     //redo basis - first take ALL columns out
    613                     int iColumn;
    614                     double largeValue = model->largeValue();
    615                     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    616                          if (model->getColumnStatus(iColumn) == ClpSimplex::basic) {
    617                               // take out
    618                               if (!valuesPass) {
    619                                    double lower = columnLower[iColumn];
    620                                    double upper = columnUpper[iColumn];
    621                                    double value = columnActivity[iColumn];
    622                                    if (lower > -largeValue || upper < largeValue) {
    623                                         if (fabs(value - lower) < fabs(value - upper)) {
    624                                              model->setColumnStatus(iColumn, ClpSimplex::atLowerBound);
    625                                              columnActivity[iColumn] = lower;
    626                                         } else {
    627                                              model->setColumnStatus(iColumn, ClpSimplex::atUpperBound);
    628                                              columnActivity[iColumn] = upper;
    629                                         }
    630                                    } else {
    631                                         model->setColumnStatus(iColumn, ClpSimplex::isFree);
    632                                    }
    633                               } else {
    634                                    model->setColumnStatus(iColumn, ClpSimplex::superBasic);
    635                               }
    636                          }
    637                     }
    638                     int iRow;
    639                     for (iRow = 0; iRow < numberRows; iRow++) {
    640                          int iSequence = pivotVariable[iRow];
    641                          if (iSequence >= 0) {
    642                               // basic
    643                               if (iSequence >= numberColumns) {
    644                                    // slack in - leave
    645                                    //assert (iSequence-numberColumns==iRow);
    646                               } else {
    647                                    assert(model->getRowStatus(iRow) != ClpSimplex::basic);
    648                                    // put back structural
    649                                    model->setColumnStatus(iSequence, ClpSimplex::basic);
    650                               }
    651                          } else {
    652                               // put in slack
    653                               model->setRowStatus(iRow, ClpSimplex::basic);
    654                          }
    655                     }
    656                     // Put back any key variables for gub (status_ not touched)
    657                     matrix->generalExpanded(model, 1, status_);
    658                     // signal repeat
    659                     status_ = -99;
    660                     // set fixed if they are
    661                     for (iRow = 0; iRow < numberRows; iRow++) {
    662                          if (model->getRowStatus(iRow) != ClpSimplex::basic ) {
    663                               if (rowLower[iRow] == rowUpper[iRow]) {
    664                                    rowActivity[iRow] = rowLower[iRow];
    665                                    model->setRowStatus(iRow, ClpSimplex::isFixed);
    666                               }
    667                          }
    668                     }
    669                     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    670                          if (model->getColumnStatus(iColumn) != ClpSimplex::basic ) {
    671                               if (columnLower[iColumn] == columnUpper[iColumn]) {
    672                                    columnActivity[iColumn] = columnLower[iColumn];
    673                                    model->setColumnStatus(iColumn, ClpSimplex::isFixed);
    674                               }
    675                          }
    676                     }
    677                }
    678           }
    679 #ifndef SLIM_CLP
    680      } else {
    681           // network - fake factorization - do nothing
    682           status_ = 0;
    683           numberPivots_ = 0;
    684      }
    685 #endif
    686 #ifndef SLIM_CLP
    687      if (!status_) {
    688           // take out part if quadratic
    689           if (model->algorithm() == 2) {
    690                ClpObjective * obj = model->objectiveAsObject();
     587        for (i = 0; i < numberRows; i++) {
     588          if (model->getRowStatus(i) == ClpSimplex::basic) {
     589            int iPivot = rowIsBasic[i];
     590            if (iPivot >= 0)
     591              pivotVariable[iPivot] = i + numberColumns;
     592          }
     593        }
     594        for (i = 0; i < numberColumns; i++) {
     595          if (model->getColumnStatus(i) == ClpSimplex::basic) {
     596            int iPivot = columnIsBasic[i];
     597            if (iPivot >= 0)
     598              pivotVariable[iPivot] = i;
     599          }
     600        }
     601        delete[] isBasic;
     602        double *columnLower = model->lowerRegion();
     603        double *columnUpper = model->upperRegion();
     604        double *columnActivity = model->solutionRegion();
     605        double *rowLower = model->lowerRegion(0);
     606        double *rowUpper = model->upperRegion(0);
     607        double *rowActivity = model->solutionRegion(0);
     608        //redo basis - first take ALL columns out
     609        int iColumn;
     610        double largeValue = model->largeValue();
     611        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     612          if (model->getColumnStatus(iColumn) == ClpSimplex::basic) {
     613            // take out
     614            if (!valuesPass) {
     615              double lower = columnLower[iColumn];
     616              double upper = columnUpper[iColumn];
     617              double value = columnActivity[iColumn];
     618              if (lower > -largeValue || upper < largeValue) {
     619                if (fabs(value - lower) < fabs(value - upper)) {
     620                  model->setColumnStatus(iColumn, ClpSimplex::atLowerBound);
     621                  columnActivity[iColumn] = lower;
     622                } else {
     623                  model->setColumnStatus(iColumn, ClpSimplex::atUpperBound);
     624                  columnActivity[iColumn] = upper;
     625                }
     626              } else {
     627                model->setColumnStatus(iColumn, ClpSimplex::isFree);
     628              }
     629            } else {
     630              model->setColumnStatus(iColumn, ClpSimplex::superBasic);
     631            }
     632          }
     633        }
     634        int iRow;
     635        for (iRow = 0; iRow < numberRows; iRow++) {
     636          int iSequence = pivotVariable[iRow];
     637          if (iSequence >= 0) {
     638            // basic
     639            if (iSequence >= numberColumns) {
     640              // slack in - leave
     641              //assert (iSequence-numberColumns==iRow);
     642            } else {
     643              assert(model->getRowStatus(iRow) != ClpSimplex::basic);
     644              // put back structural
     645              model->setColumnStatus(iSequence, ClpSimplex::basic);
     646            }
     647          } else {
     648            // put in slack
     649            model->setRowStatus(iRow, ClpSimplex::basic);
     650          }
     651        }
     652        // Put back any key variables for gub (status_ not touched)
     653        matrix->generalExpanded(model, 1, status_);
     654        // signal repeat
     655        status_ = -99;
     656        // set fixed if they are
     657        for (iRow = 0; iRow < numberRows; iRow++) {
     658          if (model->getRowStatus(iRow) != ClpSimplex::basic) {
     659            if (rowLower[iRow] == rowUpper[iRow]) {
     660              rowActivity[iRow] = rowLower[iRow];
     661              model->setRowStatus(iRow, ClpSimplex::isFixed);
     662            }
     663          }
     664        }
     665        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     666          if (model->getColumnStatus(iColumn) != ClpSimplex::basic) {
     667            if (columnLower[iColumn] == columnUpper[iColumn]) {
     668              columnActivity[iColumn] = columnLower[iColumn];
     669              model->setColumnStatus(iColumn, ClpSimplex::isFixed);
     670            }
     671          }
     672        }
     673      }
     674    }
     675#ifndef SLIM_CLP
     676  } else {
     677    // network - fake factorization - do nothing
     678    status_ = 0;
     679    numberPivots_ = 0;
     680  }
     681#endif
     682#ifndef SLIM_CLP
     683  if (!status_) {
     684    // take out part if quadratic
     685    if (model->algorithm() == 2) {
     686      ClpObjective *obj = model->objectiveAsObject();
    691687#ifndef NDEBUG
    692                ClpQuadraticObjective * quadraticObj = (dynamic_cast< ClpQuadraticObjective*>(obj));
    693                assert (quadraticObj);
     688      ClpQuadraticObjective *quadraticObj = (dynamic_cast< ClpQuadraticObjective * >(obj));
     689      assert(quadraticObj);
    694690#else
    695                ClpQuadraticObjective * quadraticObj = (static_cast< ClpQuadraticObjective*>(obj));
    696 #endif
    697                CoinPackedMatrix * quadratic = quadraticObj->quadraticObjective();
    698                int numberXColumns = quadratic->getNumCols();
    699                assert (numberXColumns < numberColumns);
    700                int base = numberColumns - numberXColumns;
    701                int * which = new int [numberXColumns];
    702                int * pivotVariable = model->pivotVariable();
    703                int * permute = pivotColumn();
    704                int i;
    705                int n = 0;
    706                for (i = 0; i < numberRows; i++) {
    707                     int iSj = pivotVariable[i] - base;
    708                     if (iSj >= 0 && iSj < numberXColumns)
    709                          which[n++] = permute[i];
    710                }
    711                if (n)
    712                     emptyRows(n, which);
    713                delete [] which;
    714           }
    715      }
    716 #endif
    717 #ifdef CLP_FACTORIZATION_INSTRUMENT
    718      factorization_instrument(2);
    719 #endif
    720      return status_;
     691      ClpQuadraticObjective *quadraticObj = (static_cast< ClpQuadraticObjective * >(obj));
     692#endif
     693      CoinPackedMatrix *quadratic = quadraticObj->quadraticObjective();
     694      int numberXColumns = quadratic->getNumCols();
     695      assert(numberXColumns < numberColumns);
     696      int base = numberColumns - numberXColumns;
     697      int *which = new int[numberXColumns];
     698      int *pivotVariable = model->pivotVariable();
     699      int *permute = pivotColumn();
     700      int i;
     701      int n = 0;
     702      for (i = 0; i < numberRows; i++) {
     703        int iSj = pivotVariable[i] - base;
     704        if (iSj >= 0 && iSj < numberXColumns)
     705          which[n++] = permute[i];
     706      }
     707      if (n)
     708        emptyRows(n, which);
     709      delete[] which;
     710    }
     711  }
     712#endif
     713#ifdef CLP_FACTORIZATION_INSTRUMENT
     714  factorization_instrument(2);
     715#endif
     716  return status_;
    721717}
    722718/* Replaces one Column to basis,
     
    727723   after factorization and thereafter re-factorize
    728724   partial update already in U */
    729 int
    730 ClpFactorization::replaceColumn ( const ClpSimplex * model,
    731                                   CoinIndexedVector * regionSparse,
    732                                   CoinIndexedVector * tableauColumn,
    733                                   int pivotRow,
    734                                   double pivotCheck ,
    735                                   bool checkBeforeModifying,
    736                                   double acceptablePivot)
    737 {
    738      int returnCode;
    739 #ifndef SLIM_CLP
    740      if (!networkBasis_) {
    741 #endif
    742 #ifdef CLP_FACTORIZATION_INSTRUMENT
    743           factorization_instrument(-1);
    744 #endif
    745           // see if FT
    746           if (doForrestTomlin_) {
    747              returnCode = CoinFactorization::replaceColumn(regionSparse,
    748                             pivotRow,
    749                             pivotCheck,
    750                             checkBeforeModifying,
    751                             acceptablePivot);
    752           } else {
    753                returnCode = CoinFactorization::replaceColumnPFI(tableauColumn,
    754                             pivotRow, pivotCheck); // Note array
    755           }
    756 #ifdef CLP_FACTORIZATION_INSTRUMENT
    757           factorization_instrument(3);
     725int ClpFactorization::replaceColumn(const ClpSimplex *model,
     726  CoinIndexedVector *regionSparse,
     727  CoinIndexedVector *tableauColumn,
     728  int pivotRow,
     729  double pivotCheck,
     730  bool checkBeforeModifying,
     731  double acceptablePivot)
     732{
     733  int returnCode;
     734#ifndef SLIM_CLP
     735  if (!networkBasis_) {
     736#endif
     737#ifdef CLP_FACTORIZATION_INSTRUMENT
     738    factorization_instrument(-1);
     739#endif
     740    // see if FT
     741    if (doForrestTomlin_) {
     742      returnCode = CoinFactorization::replaceColumn(regionSparse,
     743        pivotRow,
     744        pivotCheck,
     745        checkBeforeModifying,
     746        acceptablePivot);
     747    } else {
     748      returnCode = CoinFactorization::replaceColumnPFI(tableauColumn,
     749        pivotRow, pivotCheck); // Note array
     750    }
     751#ifdef CLP_FACTORIZATION_INSTRUMENT
     752    factorization_instrument(3);
    758753#endif
    759754
    760755#ifndef SLIM_CLP
    761      } else {
    762           if (doCheck) {
    763                returnCode = CoinFactorization::replaceColumn(regionSparse,
    764                             pivotRow,
    765                             pivotCheck,
    766                             checkBeforeModifying,
    767                             acceptablePivot);
    768                networkBasis_->replaceColumn(regionSparse,
    769                                             pivotRow);
    770           } else {
    771                // increase number of pivots
    772                numberPivots_++;
    773                returnCode = networkBasis_->replaceColumn(regionSparse,
    774                             pivotRow);
    775           }
    776      }
    777 #endif
    778      return returnCode;
     756  } else {
     757    if (doCheck) {
     758      returnCode = CoinFactorization::replaceColumn(regionSparse,
     759        pivotRow,
     760        pivotCheck,
     761        checkBeforeModifying,
     762        acceptablePivot);
     763      networkBasis_->replaceColumn(regionSparse,
     764        pivotRow);
     765    } else {
     766      // increase number of pivots
     767      numberPivots_++;
     768      returnCode = networkBasis_->replaceColumn(regionSparse,
     769        pivotRow);
     770    }
     771  }
     772#endif
     773  return returnCode;
    779774}
    780775
     
    782777   number returned is negative if no room
    783778   region1 starts as zero and is zero at end */
    784 int
    785 ClpFactorization::updateColumnFT ( CoinIndexedVector * regionSparse,
    786                                    CoinIndexedVector * regionSparse2)
     779int ClpFactorization::updateColumnFT(CoinIndexedVector *regionSparse,
     780  CoinIndexedVector *regionSparse2)
    787781{
    788782#ifdef CLP_DEBUG
    789      regionSparse->checkClear();
    790 #endif
    791      if (!numberRows_)
    792           return 0;
    793 #ifndef SLIM_CLP
    794      if (!networkBasis_) {
    795 #endif
    796 #ifdef CLP_FACTORIZATION_INSTRUMENT
    797           factorization_instrument(-1);
    798 #endif
    799           collectStatistics_ = true;
    800           int returnCode = CoinFactorization::updateColumnFT(regionSparse,
    801                            regionSparse2);
    802           collectStatistics_ = false;
    803 #ifdef CLP_FACTORIZATION_INSTRUMENT
    804           factorization_instrument(4);
    805 #endif
    806           return returnCode;
    807 #ifndef SLIM_CLP
    808      } else {
     783  regionSparse->checkClear();
     784#endif
     785  if (!numberRows_)
     786    return 0;
     787#ifndef SLIM_CLP
     788  if (!networkBasis_) {
     789#endif
     790#ifdef CLP_FACTORIZATION_INSTRUMENT
     791    factorization_instrument(-1);
     792#endif
     793    collectStatistics_ = true;
     794    int returnCode = CoinFactorization::updateColumnFT(regionSparse,
     795      regionSparse2);
     796    collectStatistics_ = false;
     797#ifdef CLP_FACTORIZATION_INSTRUMENT
     798    factorization_instrument(4);
     799#endif
     800    return returnCode;
     801#ifndef SLIM_CLP
     802  } else {
    809803#ifdef CHECK_NETWORK
    810           CoinIndexedVector * save = new CoinIndexedVector(*regionSparse2);
    811           double * check = new double[numberRows_];
    812           int returnCode = CoinFactorization::updateColumnFT(regionSparse,
    813                            regionSparse2);
    814           networkBasis_->updateColumn(regionSparse, save, -1);
    815           int i;
    816           double * array = regionSparse2->denseVector();
    817           int * indices = regionSparse2->getIndices();
    818           int n = regionSparse2->getNumElements();
    819           memset(check, 0, numberRows_ * sizeof(double));
    820           double * array2 = save->denseVector();
    821           int * indices2 = save->getIndices();
    822           int n2 = save->getNumElements();
    823           assert (n == n2);
    824           if (save->packedMode()) {
    825                for (i = 0; i < n; i++) {
    826                     check[indices[i]] = array[i];
    827                }
    828                for (i = 0; i < n; i++) {
    829                     double value2 = array2[i];
    830                     assert (check[indices2[i]] == value2);
    831                }
    832           } else {
    833                for (i = 0; i < numberRows_; i++) {
    834                     double value1 = array[i];
    835                     double value2 = array2[i];
    836                     assert (value1 == value2);
    837                }
    838           }
    839           delete save;
    840           delete [] check;
    841           return returnCode;
     804    CoinIndexedVector *save = new CoinIndexedVector(*regionSparse2);
     805    double *check = new double[numberRows_];
     806    int returnCode = CoinFactorization::updateColumnFT(regionSparse,
     807      regionSparse2);
     808    networkBasis_->updateColumn(regionSparse, save, -1);
     809    int i;
     810    double *array = regionSparse2->denseVector();
     811    int *indices = regionSparse2->getIndices();
     812    int n = regionSparse2->getNumElements();
     813    memset(check, 0, numberRows_ * sizeof(double));
     814    double *array2 = save->denseVector();
     815    int *indices2 = save->getIndices();
     816    int n2 = save->getNumElements();
     817    assert(n == n2);
     818    if (save->packedMode()) {
     819      for (i = 0; i < n; i++) {
     820        check[indices[i]] = array[i];
     821      }
     822      for (i = 0; i < n; i++) {
     823        double value2 = array2[i];
     824        assert(check[indices2[i]] == value2);
     825      }
     826    } else {
     827      for (i = 0; i < numberRows_; i++) {
     828        double value1 = array[i];
     829        double value2 = array2[i];
     830        assert(value1 == value2);
     831      }
     832    }
     833    delete save;
     834    delete[] check;
     835    return returnCode;
    842836#else
    843           networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
    844           return 1;
    845 #endif
    846      }
     837    networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
     838    return 1;
     839#endif
     840  }
    847841#endif
    848842}
     
    850844   number returned is negative if no room
    851845   region1 starts as zero and is zero at end */
    852 int
    853 ClpFactorization::updateColumn ( CoinIndexedVector * regionSparse,
    854                                  CoinIndexedVector * regionSparse2,
    855                                  bool noPermute) const
     846int ClpFactorization::updateColumn(CoinIndexedVector *regionSparse,
     847  CoinIndexedVector *regionSparse2,
     848  bool noPermute) const
    856849{
    857850#ifdef CLP_DEBUG
    858      if (!noPermute)
    859           regionSparse->checkClear();
    860 #endif
    861      if (!numberRows_)
    862           return 0;
    863 #ifndef SLIM_CLP
    864      if (!networkBasis_) {
    865 #endif
    866 #ifdef CLP_FACTORIZATION_INSTRUMENT
    867           factorization_instrument(-1);
    868 #endif
    869           collectStatistics_ = true;
    870           int returnCode = CoinFactorization::updateColumn(regionSparse,
    871                            regionSparse2,
    872                            noPermute);
    873           collectStatistics_ = false;
    874 #ifdef CLP_FACTORIZATION_INSTRUMENT
    875           factorization_instrument(5);
    876 #endif
    877           return returnCode;
    878 #ifndef SLIM_CLP
    879      } else {
     851  if (!noPermute)
     852    regionSparse->checkClear();
     853#endif
     854  if (!numberRows_)
     855    return 0;
     856#ifndef SLIM_CLP
     857  if (!networkBasis_) {
     858#endif
     859#ifdef CLP_FACTORIZATION_INSTRUMENT
     860    factorization_instrument(-1);
     861#endif
     862    collectStatistics_ = true;
     863    int returnCode = CoinFactorization::updateColumn(regionSparse,
     864      regionSparse2,
     865      noPermute);
     866    collectStatistics_ = false;
     867#ifdef CLP_FACTORIZATION_INSTRUMENT
     868    factorization_instrument(5);
     869#endif
     870    return returnCode;
     871#ifndef SLIM_CLP
     872  } else {
    880873#ifdef CHECK_NETWORK
    881           CoinIndexedVector * save = new CoinIndexedVector(*regionSparse2);
    882           double * check = new double[numberRows_];
    883           int returnCode = CoinFactorization::updateColumn(regionSparse,
    884                            regionSparse2,
    885                            noPermute);
    886           networkBasis_->updateColumn(regionSparse, save, -1);
    887           int i;
    888           double * array = regionSparse2->denseVector();
    889           int * indices = regionSparse2->getIndices();
    890           int n = regionSparse2->getNumElements();
    891           memset(check, 0, numberRows_ * sizeof(double));
    892           double * array2 = save->denseVector();
    893           int * indices2 = save->getIndices();
    894           int n2 = save->getNumElements();
    895           assert (n == n2);
    896           if (save->packedMode()) {
    897                for (i = 0; i < n; i++) {
    898                     check[indices[i]] = array[i];
    899                }
    900                for (i = 0; i < n; i++) {
    901                     double value2 = array2[i];
    902                     assert (check[indices2[i]] == value2);
    903                }
    904           } else {
    905                for (i = 0; i < numberRows_; i++) {
    906                     double value1 = array[i];
    907                     double value2 = array2[i];
    908                     assert (value1 == value2);
    909                }
    910           }
    911           delete save;
    912           delete [] check;
    913           return returnCode;
     874    CoinIndexedVector *save = new CoinIndexedVector(*regionSparse2);
     875    double *check = new double[numberRows_];
     876    int returnCode = CoinFactorization::updateColumn(regionSparse,
     877      regionSparse2,
     878      noPermute);
     879    networkBasis_->updateColumn(regionSparse, save, -1);
     880    int i;
     881    double *array = regionSparse2->denseVector();
     882    int *indices = regionSparse2->getIndices();
     883    int n = regionSparse2->getNumElements();
     884    memset(check, 0, numberRows_ * sizeof(double));
     885    double *array2 = save->denseVector();
     886    int *indices2 = save->getIndices();
     887    int n2 = save->getNumElements();
     888    assert(n == n2);
     889    if (save->packedMode()) {
     890      for (i = 0; i < n; i++) {
     891        check[indices[i]] = array[i];
     892      }
     893      for (i = 0; i < n; i++) {
     894        double value2 = array2[i];
     895        assert(check[indices2[i]] == value2);
     896      }
     897    } else {
     898      for (i = 0; i < numberRows_; i++) {
     899        double value1 = array[i];
     900        double value2 = array2[i];
     901        assert(value1 == value2);
     902      }
     903    }
     904    delete save;
     905    delete[] check;
     906    return returnCode;
    914907#else
    915           networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
    916           return 1;
    917 #endif
    918      }
     908    networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
     909    return 1;
     910#endif
     911  }
    919912#endif
    920913}
     
    924917   Also updates region3
    925918   region1 starts as zero and is zero at end */
    926 int
    927 ClpFactorization::updateTwoColumnsFT ( CoinIndexedVector * regionSparse1,
    928                                        CoinIndexedVector * regionSparse2,
    929                                        CoinIndexedVector * regionSparse3,
    930                                        bool noPermuteRegion3)
    931 {
    932      int returnCode = updateColumnFT(regionSparse1, regionSparse2);
    933      updateColumn(regionSparse1, regionSparse3, noPermuteRegion3);
    934      return returnCode;
     919int ClpFactorization::updateTwoColumnsFT(CoinIndexedVector *regionSparse1,
     920  CoinIndexedVector *regionSparse2,
     921  CoinIndexedVector *regionSparse3,
     922  bool noPermuteRegion3)
     923{
     924  int returnCode = updateColumnFT(regionSparse1, regionSparse2);
     925  updateColumn(regionSparse1, regionSparse3, noPermuteRegion3);
     926  return returnCode;
    935927}
    936928/* Updates one column (FTRAN) from region2
    937929   number returned is negative if no room
    938930   region1 starts as zero and is zero at end */
    939 int
    940 ClpFactorization::updateColumnForDebug ( CoinIndexedVector * regionSparse,
    941           CoinIndexedVector * regionSparse2,
    942           bool noPermute) const
    943 {
    944      if (!noPermute)
    945           regionSparse->checkClear();
    946      if (!numberRows_)
    947           return 0;
    948      collectStatistics_ = false;
    949      int returnCode = CoinFactorization::updateColumn(regionSparse,
    950                       regionSparse2,
    951                       noPermute);
    952      return returnCode;
     931int ClpFactorization::updateColumnForDebug(CoinIndexedVector *regionSparse,
     932  CoinIndexedVector *regionSparse2,
     933  bool noPermute) const
     934{
     935  if (!noPermute)
     936    regionSparse->checkClear();
     937  if (!numberRows_)
     938    return 0;
     939  collectStatistics_ = false;
     940  int returnCode = CoinFactorization::updateColumn(regionSparse,
     941    regionSparse2,
     942    noPermute);
     943  return returnCode;
    953944}
    954945/* Updates one column (BTRAN) from region2
    955946   region1 starts as zero and is zero at end */
    956 int
    957 ClpFactorization::updateColumnTranspose ( CoinIndexedVector * regionSparse,
    958           CoinIndexedVector * regionSparse2) const
    959 {
    960      if (!numberRows_)
    961           return 0;
    962 #ifndef SLIM_CLP
    963      if (!networkBasis_) {
    964 #endif
    965 #ifdef CLP_FACTORIZATION_INSTRUMENT
    966           factorization_instrument(-1);
    967 #endif
    968           collectStatistics_ = true;
    969           int returnCode = CoinFactorization::updateColumnTranspose(regionSparse,
    970                            regionSparse2);
    971           collectStatistics_ = false;
    972 #ifdef CLP_FACTORIZATION_INSTRUMENT
    973           factorization_instrument(6);
    974 #endif
    975           return returnCode;
    976 #ifndef SLIM_CLP
    977      } else {
     947int ClpFactorization::updateColumnTranspose(CoinIndexedVector *regionSparse,
     948  CoinIndexedVector *regionSparse2) const
     949{
     950  if (!numberRows_)
     951    return 0;
     952#ifndef SLIM_CLP
     953  if (!networkBasis_) {
     954#endif
     955#ifdef CLP_FACTORIZATION_INSTRUMENT
     956    factorization_instrument(-1);
     957#endif
     958    collectStatistics_ = true;
     959    int returnCode = CoinFactorization::updateColumnTranspose(regionSparse,
     960      regionSparse2);
     961    collectStatistics_ = false;
     962#ifdef CLP_FACTORIZATION_INSTRUMENT
     963    factorization_instrument(6);
     964#endif
     965    return returnCode;
     966#ifndef SLIM_CLP
     967  } else {
    978968#ifdef CHECK_NETWORK
    979           CoinIndexedVector * save = new CoinIndexedVector(*regionSparse2);
    980           double * check = new double[numberRows_];
    981           int returnCode = CoinFactorization::updateColumnTranspose(regionSparse,
    982                            regionSparse2);
    983           networkBasis_->updateColumnTranspose(regionSparse, save);
    984           int i;
    985           double * array = regionSparse2->denseVector();
    986           int * indices = regionSparse2->getIndices();
    987           int n = regionSparse2->getNumElements();
    988           memset(check, 0, numberRows_ * sizeof(double));
    989           double * array2 = save->denseVector();
    990           int * indices2 = save->getIndices();
    991           int n2 = save->getNumElements();
    992           assert (n == n2);
    993           if (save->packedMode()) {
    994                for (i = 0; i < n; i++) {
    995                     check[indices[i]] = array[i];
    996                }
    997                for (i = 0; i < n; i++) {
    998                     double value2 = array2[i];
    999                     assert (check[indices2[i]] == value2);
    1000                }
    1001           } else {
    1002                for (i = 0; i < numberRows_; i++) {
    1003                     double value1 = array[i];
    1004                     double value2 = array2[i];
    1005                     assert (value1 == value2);
    1006                }
    1007           }
    1008           delete save;
    1009           delete [] check;
    1010           return returnCode;
     969    CoinIndexedVector *save = new CoinIndexedVector(*regionSparse2);
     970    double *check = new double[numberRows_];
     971    int returnCode = CoinFactorization::updateColumnTranspose(regionSparse,
     972      regionSparse2);
     973    networkBasis_->updateColumnTranspose(regionSparse, save);
     974    int i;
     975    double *array = regionSparse2->denseVector();
     976    int *indices = regionSparse2->getIndices();
     977    int n = regionSparse2->getNumElements();
     978    memset(check, 0, numberRows_ * sizeof(double));
     979    double *array2 = save->denseVector();
     980    int *indices2 = save->getIndices();
     981    int n2 = save->getNumElements();
     982    assert(n == n2);
     983    if (save->packedMode()) {
     984      for (i = 0; i < n; i++) {
     985        check[indices[i]] = array[i];
     986      }
     987      for (i = 0; i < n; i++) {
     988        double value2 = array2[i];
     989        assert(check[indices2[i]] == value2);
     990      }
     991    } else {
     992      for (i = 0; i < numberRows_; i++) {
     993        double value1 = array[i];
     994        double value2 = array2[i];
     995        assert(value1 == value2);
     996      }
     997    }
     998    delete save;
     999    delete[] check;
     1000    return returnCode;
    10111001#else
    1012           return networkBasis_->updateColumnTranspose(regionSparse, regionSparse2);
    1013 #endif
    1014      }
     1002    return networkBasis_->updateColumnTranspose(regionSparse, regionSparse2);
     1003#endif
     1004  }
    10151005#endif
    10161006}
     
    10191009   Note - if regionSparse2 packed on input - will be packed on output - same for 3
    10201010*/
    1021 void
    1022 ClpFactorization::updateTwoColumnsTranspose ( CoinIndexedVector * regionSparse,
    1023                                   CoinIndexedVector * regionSparse2,
    1024                                   CoinIndexedVector * regionSparse3) const
    1025 {
    1026      if (!numberRows_)
    1027           return 0;
    1028 #ifndef SLIM_CLP
    1029      if (!networkBasis_) {
    1030 #endif
    1031 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1032           factorization_instrument(-1);
    1033 #endif
    1034           collectStatistics_ = true;
    1035           CoinFactorization::updateTwoColumnsTranspose(regionSparse,
    1036                                                        regionSparse2,regionSparse3);
    1037           collectStatistics_ = false;
    1038 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1039           factorization_instrument(6);
    1040 #endif
    1041           return returnCode;
    1042 #ifndef SLIM_CLP
    1043      } else {
    1044        factorization__->updateColumnTranspose(regionSparse, regionSparse2);
    1045        factorization__->updateColumnTranspose(regionSparse, regionSparse3);
    1046      }
     1011void ClpFactorization::updateTwoColumnsTranspose(CoinIndexedVector *regionSparse,
     1012  CoinIndexedVector *regionSparse2,
     1013  CoinIndexedVector *regionSparse3) const
     1014{
     1015  if (!numberRows_)
     1016    return 0;
     1017#ifndef SLIM_CLP
     1018  if (!networkBasis_) {
     1019#endif
     1020#ifdef CLP_FACTORIZATION_INSTRUMENT
     1021    factorization_instrument(-1);
     1022#endif
     1023    collectStatistics_ = true;
     1024    CoinFactorization::updateTwoColumnsTranspose(regionSparse,
     1025      regionSparse2, regionSparse3);
     1026    collectStatistics_ = false;
     1027#ifdef CLP_FACTORIZATION_INSTRUMENT
     1028    factorization_instrument(6);
     1029#endif
     1030    return returnCode;
     1031#ifndef SLIM_CLP
     1032  } else {
     1033    factorization__->updateColumnTranspose(regionSparse, regionSparse2);
     1034    factorization__->updateColumnTranspose(regionSparse, regionSparse3);
     1035  }
    10471036#endif
    10481037}
    10491038/* makes a row copy of L for speed and to allow very sparse problems */
    1050 void
    1051 ClpFactorization::goSparse()
    1052 {
    1053 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1054      factorization_instrument(-1);
    1055 #endif
    1056 #ifndef SLIM_CLP
    1057      if (!networkBasis_)
    1058 #endif
    1059           CoinFactorization::goSparse();
    1060 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1061      factorization_instrument(7);
     1039void ClpFactorization::goSparse()
     1040{
     1041#ifdef CLP_FACTORIZATION_INSTRUMENT
     1042  factorization_instrument(-1);
     1043#endif
     1044#ifndef SLIM_CLP
     1045  if (!networkBasis_)
     1046#endif
     1047    CoinFactorization::goSparse();
     1048#ifdef CLP_FACTORIZATION_INSTRUMENT
     1049  factorization_instrument(7);
    10621050#endif
    10631051}
    10641052// Cleans up i.e. gets rid of network basis
    1065 void
    1066 ClpFactorization::cleanUp()
    1067 {
    1068 #ifndef SLIM_CLP
    1069      delete networkBasis_;
    1070      networkBasis_ = NULL;
    1071 #endif
    1072      resetStatistics();
     1053void ClpFactorization::cleanUp()
     1054{
     1055#ifndef SLIM_CLP
     1056  delete networkBasis_;
     1057  networkBasis_ = NULL;
     1058#endif
     1059  resetStatistics();
    10731060}
    10741061/// Says whether to redo pivot order
    1075 bool
    1076 ClpFactorization::needToReorder() const
     1062bool ClpFactorization::needToReorder() const
    10771063{
    10781064#ifdef CHECK_NETWORK
    1079      return true;
    1080 #endif
    1081 #ifndef SLIM_CLP
    1082      if (!networkBasis_)
    1083 #endif
    1084           return true;
    1085 #ifndef SLIM_CLP
    1086      else
    1087           return false;
     1065  return true;
     1066#endif
     1067#ifndef SLIM_CLP
     1068  if (!networkBasis_)
     1069#endif
     1070    return true;
     1071#ifndef SLIM_CLP
     1072  else
     1073    return false;
    10881074#endif
    10891075}
    10901076// Get weighted row list
    1091 void
    1092 ClpFactorization::getWeights(int * weights) const
    1093 {
    1094 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1095      factorization_instrument(-1);
    1096 #endif
    1097 #ifndef SLIM_CLP
    1098      if (networkBasis_) {
    1099           // Network - just unit
    1100           for (int i = 0; i < numberRows_; i++)
    1101                weights[i] = 1;
    1102           return;
    1103      }
    1104 #endif
    1105      int * numberInRow = numberInRow_.array();
    1106      int * numberInColumn = numberInColumn_.array();
    1107      int * permuteBack = pivotColumnBack_.array();
    1108      int * indexRowU = indexRowU_.array();
    1109      const CoinBigIndex * startColumnU = startColumnU_.array();
    1110      const CoinBigIndex * startRowL = startRowL_.array();
    1111      if (!startRowL || !numberInRow_.array()) {
    1112           int * temp = new int[numberRows_];
    1113           memset(temp, 0, numberRows_ * sizeof(int));
    1114           int i;
    1115           for (i = 0; i < numberRows_; i++) {
    1116                // one for pivot
    1117                temp[i]++;
    1118                CoinBigIndex j;
    1119                for (j = startColumnU[i]; j < startColumnU[i] + numberInColumn[i]; j++) {
    1120                     int iRow = indexRowU[j];
    1121                     temp[iRow]++;
    1122                }
    1123           }
    1124           CoinBigIndex * startColumnL = startColumnL_.array();
    1125           int * indexRowL = indexRowL_.array();
    1126           for (i = baseL_; i < baseL_ + numberL_; i++) {
    1127                CoinBigIndex j;
    1128                for (j = startColumnL[i]; j < startColumnL[i+1]; j++) {
    1129                     int iRow = indexRowL[j];
    1130                     temp[iRow]++;
    1131                }
    1132           }
    1133           for (i = 0; i < numberRows_; i++) {
    1134                int number = temp[i];
    1135                int iPermute = permuteBack[i];
    1136                weights[iPermute] = number;
    1137           }
    1138           delete [] temp;
    1139      } else {
    1140           int i;
    1141           for (i = 0; i < numberRows_; i++) {
    1142                int number = startRowL[i+1] - startRowL[i] + numberInRow[i] + 1;
    1143                //number = startRowL[i+1]-startRowL[i]+1;
    1144                //number = numberInRow[i]+1;
    1145                int iPermute = permuteBack[i];
    1146                weights[iPermute] = number;
    1147           }
    1148      }
    1149 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1150      factorization_instrument(8);
     1077void ClpFactorization::getWeights(int *weights) const
     1078{
     1079#ifdef CLP_FACTORIZATION_INSTRUMENT
     1080  factorization_instrument(-1);
     1081#endif
     1082#ifndef SLIM_CLP
     1083  if (networkBasis_) {
     1084    // Network - just unit
     1085    for (int i = 0; i < numberRows_; i++)
     1086      weights[i] = 1;
     1087    return;
     1088  }
     1089#endif
     1090  int *numberInRow = numberInRow_.array();
     1091  int *numberInColumn = numberInColumn_.array();
     1092  int *permuteBack = pivotColumnBack_.array();
     1093  int *indexRowU = indexRowU_.array();
     1094  const CoinBigIndex *startColumnU = startColumnU_.array();
     1095  const CoinBigIndex *startRowL = startRowL_.array();
     1096  if (!startRowL || !numberInRow_.array()) {
     1097    int *temp = new int[numberRows_];
     1098    memset(temp, 0, numberRows_ * sizeof(int));
     1099    int i;
     1100    for (i = 0; i < numberRows_; i++) {
     1101      // one for pivot
     1102      temp[i]++;
     1103      CoinBigIndex j;
     1104      for (j = startColumnU[i]; j < startColumnU[i] + numberInColumn[i]; j++) {
     1105        int iRow = indexRowU[j];
     1106        temp[iRow]++;
     1107      }
     1108    }
     1109    CoinBigIndex *startColumnL = startColumnL_.array();
     1110    int *indexRowL = indexRowL_.array();
     1111    for (i = baseL_; i < baseL_ + numberL_; i++) {
     1112      CoinBigIndex j;
     1113      for (j = startColumnL[i]; j < startColumnL[i + 1]; j++) {
     1114        int iRow = indexRowL[j];
     1115        temp[iRow]++;
     1116      }
     1117    }
     1118    for (i = 0; i < numberRows_; i++) {
     1119      int number = temp[i];
     1120      int iPermute = permuteBack[i];
     1121      weights[iPermute] = number;
     1122    }
     1123    delete[] temp;
     1124  } else {
     1125    int i;
     1126    for (i = 0; i < numberRows_; i++) {
     1127      int number = startRowL[i + 1] - startRowL[i] + numberInRow[i] + 1;
     1128      //number = startRowL[i+1]-startRowL[i]+1;
     1129      //number = numberInRow[i]+1;
     1130      int iPermute = permuteBack[i];
     1131      weights[iPermute] = number;
     1132    }
     1133  }
     1134#ifdef CLP_FACTORIZATION_INSTRUMENT
     1135  factorization_instrument(8);
    11511136#endif
    11521137}
    11531138#else
    1154 #if COIN_BIG_INDEX==0
     1139#if COIN_BIG_INDEX == 0
    11551140// This one allows multiple factorizations
    11561141#if CLP_MULTIPLE_FACTORIZATIONS == 1
     
    11741159// Default Constructor
    11751160//-------------------------------------------------------------------
    1176 ClpFactorization::ClpFactorization ()
    1177 {
    1178 #ifndef SLIM_CLP
    1179      networkBasis_ = NULL;
    1180 #endif
    1181      //coinFactorizationA_ = NULL;
    1182      coinFactorizationA_ = new CoinFactorization() ;
    1183      coinFactorizationB_ = NULL;
    1184      //coinFactorizationB_ = new CoinOtherFactorization();
    1185      forceB_ = 0;
    1186      goOslThreshold_ = -1;
    1187      goDenseThreshold_ = -1;
    1188      goSmallThreshold_ = -1;
    1189      doStatistics_ = true;
    1190      memset(&shortestAverage_,0,3*(sizeof(double)+sizeof(int)));
     1161ClpFactorization::ClpFactorization()
     1162{
     1163#ifndef SLIM_CLP
     1164  networkBasis_ = NULL;
     1165#endif
     1166  //coinFactorizationA_ = NULL;
     1167  coinFactorizationA_ = new CoinFactorization();
     1168  coinFactorizationB_ = NULL;
     1169  //coinFactorizationB_ = new CoinOtherFactorization();
     1170  forceB_ = 0;
     1171  goOslThreshold_ = -1;
     1172  goDenseThreshold_ = -1;
     1173  goSmallThreshold_ = -1;
     1174  doStatistics_ = true;
     1175  memset(&shortestAverage_, 0, 3 * (sizeof(double) + sizeof(int)));
    11911176}
    11921177
     
    11941179// Copy constructor
    11951180//-------------------------------------------------------------------
    1196 ClpFactorization::ClpFactorization (const ClpFactorization & rhs,
    1197                                     int denseIfSmaller)
    1198 {
    1199 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1200      factorization_instrument(-1);
    1201 #endif
    1202 #ifndef SLIM_CLP
    1203      if (rhs.networkBasis_)
    1204           networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
    1205      else
    1206           networkBasis_ = NULL;
    1207 #endif
    1208      forceB_ = rhs.forceB_;
    1209      goOslThreshold_ = rhs.goOslThreshold_;
    1210      goDenseThreshold_ = rhs.goDenseThreshold_;
    1211      goSmallThreshold_ = rhs.goSmallThreshold_;
    1212      doStatistics_ = rhs.doStatistics_;
    1213      int goDense = 0;
     1181ClpFactorization::ClpFactorization(const ClpFactorization &rhs,
     1182  int denseIfSmaller)
     1183{
     1184#ifdef CLP_FACTORIZATION_INSTRUMENT
     1185  factorization_instrument(-1);
     1186#endif
     1187#ifndef SLIM_CLP
     1188  if (rhs.networkBasis_)
     1189    networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
     1190  else
     1191    networkBasis_ = NULL;
     1192#endif
     1193  forceB_ = rhs.forceB_;
     1194  goOslThreshold_ = rhs.goOslThreshold_;
     1195  goDenseThreshold_ = rhs.goDenseThreshold_;
     1196  goSmallThreshold_ = rhs.goSmallThreshold_;
     1197  doStatistics_ = rhs.doStatistics_;
     1198  int goDense = 0;
    12141199#ifdef CLP_REUSE_ETAS
    1215      model_=rhs.model_;
    1216 #endif
    1217 #if COIN_BIG_INDEX==0
    1218      if (denseIfSmaller > 0 && denseIfSmaller <= goDenseThreshold_) {
    1219           CoinDenseFactorization * denseR =
    1220                dynamic_cast<CoinDenseFactorization *>(rhs.coinFactorizationB_);
    1221           if (!denseR)
    1222                goDense = 1;
    1223      }
    1224      if (denseIfSmaller > 0 && !rhs.coinFactorizationB_) {
    1225           if (denseIfSmaller <= goDenseThreshold_)
    1226                goDense = 1;
    1227           else if (denseIfSmaller <= goSmallThreshold_)
    1228                goDense = 2;
    1229           else if (denseIfSmaller <= goOslThreshold_)
    1230                goDense = 3;
    1231      } else if (denseIfSmaller < 0) {
    1232           if (-denseIfSmaller <= goDenseThreshold_)
    1233                goDense = 1;
    1234           else if (-denseIfSmaller <= goSmallThreshold_)
    1235                goDense = 2;
    1236           else if (-denseIfSmaller <= goOslThreshold_)
    1237                goDense = 3;
    1238      }
    1239 #endif
    1240      if (rhs.coinFactorizationA_ && !goDense)
    1241           coinFactorizationA_ = new CoinFactorization(*(rhs.coinFactorizationA_));
    1242      else
    1243           coinFactorizationA_ = NULL;
    1244      if (rhs.coinFactorizationB_ && (denseIfSmaller >= 0 || !goDense))
    1245           coinFactorizationB_ = rhs.coinFactorizationB_->clone();
    1246      else
    1247           coinFactorizationB_ = NULL;
    1248 #if COIN_BIG_INDEX==0
    1249      if (goDense) {
    1250           delete coinFactorizationB_;
    1251           if (goDense == 1)
    1252                coinFactorizationB_ = new CoinDenseFactorization();
    1253           else if (goDense == 2)
    1254                coinFactorizationB_ = new CoinSimpFactorization();
    1255           else
    1256                coinFactorizationB_ = new CoinOslFactorization();
    1257           if (rhs.coinFactorizationA_) {
    1258                coinFactorizationB_->maximumPivots(rhs.coinFactorizationA_->maximumPivots());
    1259                coinFactorizationB_->pivotTolerance(rhs.coinFactorizationA_->pivotTolerance());
    1260                coinFactorizationB_->zeroTolerance(rhs.coinFactorizationA_->zeroTolerance());
    1261           } else {
    1262                assert (coinFactorizationB_);
    1263                coinFactorizationB_->maximumPivots(rhs.coinFactorizationB_->maximumPivots());
    1264                coinFactorizationB_->pivotTolerance(rhs.coinFactorizationB_->pivotTolerance());
    1265                coinFactorizationB_->zeroTolerance(rhs.coinFactorizationB_->zeroTolerance());
    1266           }
    1267      }
    1268 #endif
    1269      assert (!coinFactorizationA_ || !coinFactorizationB_);
    1270 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1271      factorization_instrument(1);
    1272 #endif
    1273      memcpy(&shortestAverage_,&rhs.shortestAverage_,3*(sizeof(double)+sizeof(int)));
     1200  model_ = rhs.model_;
     1201#endif
     1202#if COIN_BIG_INDEX == 0
     1203  if (denseIfSmaller > 0 && denseIfSmaller <= goDenseThreshold_) {
     1204    CoinDenseFactorization *denseR = dynamic_cast< CoinDenseFactorization * >(rhs.coinFactorizationB_);
     1205    if (!denseR)
     1206      goDense = 1;
     1207  }
     1208  if (denseIfSmaller > 0 && !rhs.coinFactorizationB_) {
     1209    if (denseIfSmaller <= goDenseThreshold_)
     1210      goDense = 1;
     1211    else if (denseIfSmaller <= goSmallThreshold_)
     1212      goDense = 2;
     1213    else if (denseIfSmaller <= goOslThreshold_)
     1214      goDense = 3;
     1215  } else if (denseIfSmaller < 0) {
     1216    if (-denseIfSmaller <= goDenseThreshold_)
     1217      goDense = 1;
     1218    else if (-denseIfSmaller <= goSmallThreshold_)
     1219      goDense = 2;
     1220    else if (-denseIfSmaller <= goOslThreshold_)
     1221      goDense = 3;
     1222  }
     1223#endif
     1224  if (rhs.coinFactorizationA_ && !goDense)
     1225    coinFactorizationA_ = new CoinFactorization(*(rhs.coinFactorizationA_));
     1226  else
     1227    coinFactorizationA_ = NULL;
     1228  if (rhs.coinFactorizationB_ && (denseIfSmaller >= 0 || !goDense))
     1229    coinFactorizationB_ = rhs.coinFactorizationB_->clone();
     1230  else
     1231    coinFactorizationB_ = NULL;
     1232#if COIN_BIG_INDEX == 0
     1233  if (goDense) {
     1234    delete coinFactorizationB_;
     1235    if (goDense == 1)
     1236      coinFactorizationB_ = new CoinDenseFactorization();
     1237    else if (goDense == 2)
     1238      coinFactorizationB_ = new CoinSimpFactorization();
     1239    else
     1240      coinFactorizationB_ = new CoinOslFactorization();
     1241    if (rhs.coinFactorizationA_) {
     1242      coinFactorizationB_->maximumPivots(rhs.coinFactorizationA_->maximumPivots());
     1243      coinFactorizationB_->pivotTolerance(rhs.coinFactorizationA_->pivotTolerance());
     1244      coinFactorizationB_->zeroTolerance(rhs.coinFactorizationA_->zeroTolerance());
     1245    } else {
     1246      assert(coinFactorizationB_);
     1247      coinFactorizationB_->maximumPivots(rhs.coinFactorizationB_->maximumPivots());
     1248      coinFactorizationB_->pivotTolerance(rhs.coinFactorizationB_->pivotTolerance());
     1249      coinFactorizationB_->zeroTolerance(rhs.coinFactorizationB_->zeroTolerance());
     1250    }
     1251  }
     1252#endif
     1253  assert(!coinFactorizationA_ || !coinFactorizationB_);
     1254#ifdef CLP_FACTORIZATION_INSTRUMENT
     1255  factorization_instrument(1);
     1256#endif
     1257  memcpy(&shortestAverage_, &rhs.shortestAverage_, 3 * (sizeof(double) + sizeof(int)));
    12741258}
    12751259
    1276 ClpFactorization::ClpFactorization (const CoinFactorization & rhs)
    1277 {
    1278 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1279      factorization_instrument(-1);
    1280 #endif
    1281 #ifndef SLIM_CLP
    1282      networkBasis_ = NULL;
    1283 #endif
    1284      coinFactorizationA_ = new CoinFactorization(rhs);
    1285      coinFactorizationB_ = NULL;
    1286 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1287      factorization_instrument(1);
    1288 #endif
    1289      forceB_ = 0;
    1290      goOslThreshold_ = -1;
    1291      goDenseThreshold_ = -1;
    1292      goSmallThreshold_ = -1;
    1293      doStatistics_ = true;
    1294      assert (!coinFactorizationA_ || !coinFactorizationB_);
    1295      memset(&shortestAverage_,0,3*(sizeof(double)+sizeof(int)));
     1260ClpFactorization::ClpFactorization(const CoinFactorization &rhs)
     1261{
     1262#ifdef CLP_FACTORIZATION_INSTRUMENT
     1263  factorization_instrument(-1);
     1264#endif
     1265#ifndef SLIM_CLP
     1266  networkBasis_ = NULL;
     1267#endif
     1268  coinFactorizationA_ = new CoinFactorization(rhs);
     1269  coinFactorizationB_ = NULL;
     1270#ifdef CLP_FACTORIZATION_INSTRUMENT
     1271  factorization_instrument(1);
     1272#endif
     1273  forceB_ = 0;
     1274  goOslThreshold_ = -1;
     1275  goDenseThreshold_ = -1;
     1276  goSmallThreshold_ = -1;
     1277  doStatistics_ = true;
     1278  assert(!coinFactorizationA_ || !coinFactorizationB_);
     1279  memset(&shortestAverage_, 0, 3 * (sizeof(double) + sizeof(int)));
    12961280}
    12971281
    1298 ClpFactorization::ClpFactorization (const CoinOtherFactorization & rhs)
    1299 {
    1300 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1301      factorization_instrument(-1);
    1302 #endif
    1303 #ifndef SLIM_CLP
    1304      networkBasis_ = NULL;
    1305 #endif
    1306      coinFactorizationA_ = NULL;
    1307      coinFactorizationB_ = rhs.clone();
    1308      //coinFactorizationB_ = new CoinOtherFactorization(rhs);
    1309      forceB_ = 0;
    1310      goOslThreshold_ = -1;
    1311      goDenseThreshold_ = -1;
    1312      goSmallThreshold_ = -1;
    1313      doStatistics_ = true;
    1314 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1315      factorization_instrument(1);
    1316 #endif
    1317      assert (!coinFactorizationA_ || !coinFactorizationB_);
    1318      memset(&shortestAverage_,0,3*(sizeof(double)+sizeof(int)));
     1282ClpFactorization::ClpFactorization(const CoinOtherFactorization &rhs)
     1283{
     1284#ifdef CLP_FACTORIZATION_INSTRUMENT
     1285  factorization_instrument(-1);
     1286#endif
     1287#ifndef SLIM_CLP
     1288  networkBasis_ = NULL;
     1289#endif
     1290  coinFactorizationA_ = NULL;
     1291  coinFactorizationB_ = rhs.clone();
     1292  //coinFactorizationB_ = new CoinOtherFactorization(rhs);
     1293  forceB_ = 0;
     1294  goOslThreshold_ = -1;
     1295  goDenseThreshold_ = -1;
     1296  goSmallThreshold_ = -1;
     1297  doStatistics_ = true;
     1298#ifdef CLP_FACTORIZATION_INSTRUMENT
     1299  factorization_instrument(1);
     1300#endif
     1301  assert(!coinFactorizationA_ || !coinFactorizationB_);
     1302  memset(&shortestAverage_, 0, 3 * (sizeof(double) + sizeof(int)));
    13191303}
    13201304
     
    13221306// Destructor
    13231307//-------------------------------------------------------------------
    1324 ClpFactorization::~ClpFactorization ()
    1325 {
    1326 #ifndef SLIM_CLP
    1327      delete networkBasis_;
    1328 #endif
    1329      delete coinFactorizationA_;
    1330      delete coinFactorizationB_;
     1308ClpFactorization::~ClpFactorization()
     1309{
     1310#ifndef SLIM_CLP
     1311  delete networkBasis_;
     1312#endif
     1313  delete coinFactorizationA_;
     1314  delete coinFactorizationB_;
    13311315}
    13321316
     
    13351319//-------------------------------------------------------------------
    13361320ClpFactorization &
    1337 ClpFactorization::operator=(const ClpFactorization& rhs)
    1338 {
    1339 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1340      factorization_instrument(-1);
    1341 #endif
    1342      if (this != &rhs) {
    1343 #ifndef SLIM_CLP
    1344           delete networkBasis_;
    1345           if (rhs.networkBasis_)
    1346                networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
    1347           else
    1348                networkBasis_ = NULL;
    1349 #endif
    1350           forceB_ = rhs.forceB_;
     1321ClpFactorization::operator=(const ClpFactorization &rhs)
     1322{
     1323#ifdef CLP_FACTORIZATION_INSTRUMENT
     1324  factorization_instrument(-1);
     1325#endif
     1326  if (this != &rhs) {
     1327#ifndef SLIM_CLP
     1328    delete networkBasis_;
     1329    if (rhs.networkBasis_)
     1330      networkBasis_ = new ClpNetworkBasis(*(rhs.networkBasis_));
     1331    else
     1332      networkBasis_ = NULL;
     1333#endif
     1334    forceB_ = rhs.forceB_;
    13511335#ifdef CLP_REUSE_ETAS
    1352           model_=rhs.model_;
    1353 #endif
    1354           goOslThreshold_ = rhs.goOslThreshold_;
    1355           goDenseThreshold_ = rhs.goDenseThreshold_;
    1356           goSmallThreshold_ = rhs.goSmallThreshold_;
    1357           doStatistics_ = rhs.doStatistics_;
    1358           memcpy(&shortestAverage_,&rhs.shortestAverage_,3*(sizeof(double)+sizeof(int)));
    1359           if (rhs.coinFactorizationA_) {
    1360                if (coinFactorizationA_)
    1361                     *coinFactorizationA_ = *(rhs.coinFactorizationA_);
    1362                else
    1363                     coinFactorizationA_ = new CoinFactorization(*rhs.coinFactorizationA_);
    1364           } else {
    1365                delete coinFactorizationA_;
    1366                coinFactorizationA_ = NULL;
    1367           }
    1368 #if COIN_BIG_INDEX==0
    1369           if (rhs.coinFactorizationB_) {
    1370                if (coinFactorizationB_) {
    1371                     CoinDenseFactorization * denseR = dynamic_cast<CoinDenseFactorization *>(rhs.coinFactorizationB_);
    1372                     CoinDenseFactorization * dense = dynamic_cast<CoinDenseFactorization *>(coinFactorizationB_);
    1373                     CoinOslFactorization * oslR = dynamic_cast<CoinOslFactorization *>(rhs.coinFactorizationB_);
    1374                     CoinOslFactorization * osl = dynamic_cast<CoinOslFactorization *>(coinFactorizationB_);
    1375                     CoinSimpFactorization * simpR = dynamic_cast<CoinSimpFactorization *>(rhs.coinFactorizationB_);
    1376                     CoinSimpFactorization * simp = dynamic_cast<CoinSimpFactorization *>(coinFactorizationB_);
    1377                     if (dense && denseR) {
    1378                          *dense = *denseR;
    1379                     } else if (osl && oslR) {
    1380                          *osl = *oslR;
    1381                     } else if (simp && simpR) {
    1382                          *simp = *simpR;
    1383                     } else {
    1384                          delete coinFactorizationB_;
    1385                          coinFactorizationB_ = rhs.coinFactorizationB_->clone();
    1386                     }
    1387                } else {
    1388                     coinFactorizationB_ = rhs.coinFactorizationB_->clone();
    1389                }
    1390           } else {
    1391                delete coinFactorizationB_;
    1392                coinFactorizationB_ = NULL;
    1393           }
    1394 #endif
    1395      }
    1396 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1397      factorization_instrument(1);
    1398 #endif
    1399      assert (!coinFactorizationA_ || !coinFactorizationB_);
    1400      return *this;
     1336    model_ = rhs.model_;
     1337#endif
     1338    goOslThreshold_ = rhs.goOslThreshold_;
     1339    goDenseThreshold_ = rhs.goDenseThreshold_;
     1340    goSmallThreshold_ = rhs.goSmallThreshold_;
     1341    doStatistics_ = rhs.doStatistics_;
     1342    memcpy(&shortestAverage_, &rhs.shortestAverage_, 3 * (sizeof(double) + sizeof(int)));
     1343    if (rhs.coinFactorizationA_) {
     1344      if (coinFactorizationA_)
     1345        *coinFactorizationA_ = *(rhs.coinFactorizationA_);
     1346      else
     1347        coinFactorizationA_ = new CoinFactorization(*rhs.coinFactorizationA_);
     1348    } else {
     1349      delete coinFactorizationA_;
     1350      coinFactorizationA_ = NULL;
     1351    }
     1352#if COIN_BIG_INDEX == 0
     1353    if (rhs.coinFactorizationB_) {
     1354      if (coinFactorizationB_) {
     1355        CoinDenseFactorization *denseR = dynamic_cast< CoinDenseFactorization * >(rhs.coinFactorizationB_);
     1356        CoinDenseFactorization *dense = dynamic_cast< CoinDenseFactorization * >(coinFactorizationB_);
     1357        CoinOslFactorization *oslR = dynamic_cast< CoinOslFactorization * >(rhs.coinFactorizationB_);
     1358        CoinOslFactorization *osl = dynamic_cast< CoinOslFactorization * >(coinFactorizationB_);
     1359        CoinSimpFactorization *simpR = dynamic_cast< CoinSimpFactorization * >(rhs.coinFactorizationB_);
     1360        CoinSimpFactorization *simp = dynamic_cast< CoinSimpFactorization * >(coinFactorizationB_);
     1361        if (dense && denseR) {
     1362          *dense = *denseR;
     1363        } else if (osl && oslR) {
     1364          *osl = *oslR;
     1365        } else if (simp && simpR) {
     1366          *simp = *simpR;
     1367        } else {
     1368          delete coinFactorizationB_;
     1369          coinFactorizationB_ = rhs.coinFactorizationB_->clone();
     1370        }
     1371      } else {
     1372        coinFactorizationB_ = rhs.coinFactorizationB_->clone();
     1373      }
     1374    } else {
     1375      delete coinFactorizationB_;
     1376      coinFactorizationB_ = NULL;
     1377    }
     1378#endif
     1379  }
     1380#ifdef CLP_FACTORIZATION_INSTRUMENT
     1381  factorization_instrument(1);
     1382#endif
     1383  assert(!coinFactorizationA_ || !coinFactorizationB_);
     1384  return *this;
    14011385}
    14021386// Go over to dense code
    1403 void
    1404 ClpFactorization::goDenseOrSmall(int numberRows)
    1405 {
    1406      if (!forceB_) {
    1407           if (numberRows <= goDenseThreshold_) {
    1408                delete coinFactorizationA_;
    1409                delete coinFactorizationB_;
    1410                coinFactorizationA_ = NULL;
    1411                coinFactorizationB_ = new CoinDenseFactorization();
    1412                //printf("going dense\n");
    1413 #if COIN_BIG_INDEX==0
    1414           } else if (numberRows <= goSmallThreshold_) {
    1415                delete coinFactorizationA_;
    1416                delete coinFactorizationB_;
    1417                coinFactorizationA_ = NULL;
    1418                coinFactorizationB_ = new CoinSimpFactorization();
    1419                //printf("going small\n");
    1420           } else if (numberRows <= goOslThreshold_) {
    1421                delete coinFactorizationA_;
    1422                delete coinFactorizationB_;
    1423                coinFactorizationA_ = NULL;
    1424                coinFactorizationB_ = new CoinOslFactorization();
    1425                //printf("going small\n");
    1426 #endif
    1427           }
    1428      }
    1429      assert (!coinFactorizationA_ || !coinFactorizationB_);
     1387void ClpFactorization::goDenseOrSmall(int numberRows)
     1388{
     1389  if (!forceB_) {
     1390    if (numberRows <= goDenseThreshold_) {
     1391      delete coinFactorizationA_;
     1392      delete coinFactorizationB_;
     1393      coinFactorizationA_ = NULL;
     1394      coinFactorizationB_ = new CoinDenseFactorization();
     1395      //printf("going dense\n");
     1396#if COIN_BIG_INDEX == 0
     1397    } else if (numberRows <= goSmallThreshold_) {
     1398      delete coinFactorizationA_;
     1399      delete coinFactorizationB_;
     1400      coinFactorizationA_ = NULL;
     1401      coinFactorizationB_ = new CoinSimpFactorization();
     1402      //printf("going small\n");
     1403    } else if (numberRows <= goOslThreshold_) {
     1404      delete coinFactorizationA_;
     1405      delete coinFactorizationB_;
     1406      coinFactorizationA_ = NULL;
     1407      coinFactorizationB_ = new CoinOslFactorization();
     1408      //printf("going small\n");
     1409#endif
     1410    }
     1411  }
     1412  assert(!coinFactorizationA_ || !coinFactorizationB_);
    14301413}
    14311414// If nonzero force use of 1,dense 2,small 3,osl
    1432 void
    1433 ClpFactorization::forceOtherFactorization(int which)
    1434 {
    1435      delete coinFactorizationB_;
    1436      forceB_ = 0;
    1437      coinFactorizationB_ = NULL;
    1438      if (which > 0 && which < 4) {
    1439           delete coinFactorizationA_;
    1440           coinFactorizationA_ = NULL;
    1441           forceB_ = which;
    1442           switch (which) {
    1443           case 1:
    1444                coinFactorizationB_ = new CoinDenseFactorization();
    1445                goDenseThreshold_ = COIN_INT_MAX;
    1446                break;
    1447 #if COIN_BIG_INDEX==0
    1448           case 2:
    1449                coinFactorizationB_ = new CoinSimpFactorization();
    1450                goSmallThreshold_ = COIN_INT_MAX;
    1451                break;
    1452           case 3:
    1453                coinFactorizationB_ = new CoinOslFactorization();
    1454                goOslThreshold_ = COIN_INT_MAX;
    1455                break;
    1456 #endif
    1457           }
    1458      } else if (!coinFactorizationA_) {
    1459           coinFactorizationA_ = new CoinFactorization();
    1460           goOslThreshold_ = -1;
    1461           goDenseThreshold_ = -1;
    1462           goSmallThreshold_ = -1;
    1463      }
     1415void ClpFactorization::forceOtherFactorization(int which)
     1416{
     1417  delete coinFactorizationB_;
     1418  forceB_ = 0;
     1419  coinFactorizationB_ = NULL;
     1420  if (which > 0 && which < 4) {
     1421    delete coinFactorizationA_;
     1422    coinFactorizationA_ = NULL;
     1423    forceB_ = which;
     1424    switch (which) {
     1425    case 1:
     1426      coinFactorizationB_ = new CoinDenseFactorization();
     1427      goDenseThreshold_ = COIN_INT_MAX;
     1428      break;
     1429#if COIN_BIG_INDEX == 0
     1430    case 2:
     1431      coinFactorizationB_ = new CoinSimpFactorization();
     1432      goSmallThreshold_ = COIN_INT_MAX;
     1433      break;
     1434    case 3:
     1435      coinFactorizationB_ = new CoinOslFactorization();
     1436      goOslThreshold_ = COIN_INT_MAX;
     1437      break;
     1438#endif
     1439    }
     1440  } else if (!coinFactorizationA_) {
     1441    coinFactorizationA_ = new CoinFactorization();
     1442    goOslThreshold_ = -1;
     1443    goDenseThreshold_ = -1;
     1444    goSmallThreshold_ = -1;
     1445  }
    14641446}
    14651447#ifdef CLP_FACTORIZATION_NEW_TIMING
     
    14831465extern int currentLengthU;
    14841466extern int currentTakeoutU;
    1485  extern double averageLengthR;
    1486  extern double averageLengthL;
    1487  extern double averageLengthU;
    1488  extern double scaledLengthDense;
    1489  extern double scaledLengthDenseSquared;
    1490  extern double scaledLengthL;
    1491  extern double scaledLengthR;
    1492  extern double scaledLengthU;
     1467extern double averageLengthR;
     1468extern double averageLengthL;
     1469extern double averageLengthU;
     1470extern double scaledLengthDense;
     1471extern double scaledLengthDenseSquared;
     1472extern double scaledLengthL;
     1473extern double scaledLengthR;
     1474extern double scaledLengthU;
    14931475extern int startLengthU;
    14941476extern int endLengthU;
     
    15051487#ifdef CLP_STATIC
    15061488static double shortestAverage_;
    1507 static double totalInR_=0.0;
    1508 static double totalInIncreasingU_=0.0;
     1489static double totalInR_ = 0.0;
     1490static double totalInIncreasingU_ = 0.0;
    15091491//static int lastR=0;
    15101492//static int lastU=0;
    1511 static int lastNumberPivots_=0;
    1512 static int effectiveStartNumberU_=0;
     1493static int lastNumberPivots_ = 0;
     1494static int effectiveStartNumberU_ = 0;
    15131495#else
    15141496//#define shortestAverage shortestAverage_
     
    15201502#endif
    15211503#ifdef CLP_USEFUL_PRINTOUT
    1522 static bool readTwiddle=false;
    1523 static double weightIncU=1.0;
    1524 static double weightR=2.0;
    1525 static double weightRest=1.0;
    1526 static double weightFactL=30.0;
    1527 static double weightFactDense=0.1;
    1528 static double weightNrows=10.0;
    1529 static double increaseNeeded=1.1;
     1504static bool readTwiddle = false;
     1505static double weightIncU = 1.0;
     1506static double weightR = 2.0;
     1507static double weightRest = 1.0;
     1508static double weightFactL = 30.0;
     1509static double weightFactDense = 0.1;
     1510static double weightNrows = 10.0;
     1511static double increaseNeeded = 1.1;
    15301512static double constWeightIterate = 1.0;
    15311513static double weightNrowsIterate = 3.0;
     
    15351517#define weightRest 1.0
    15361518#define weightFactL 30.0
    1537 #define weightFactDense 0.1 
     1519#define weightFactDense 0.1
    15381520#define weightNrows 10.0
    15391521#define increaseNeeded 1.1
    1540 #define constWeightIterate   1.0
    1541 #define weightNrowsIterate   3.0
    1542 #endif
    1543 bool
    1544 ClpFactorization::timeToRefactorize() const
     1522#define constWeightIterate 1.0
     1523#define weightNrowsIterate 3.0
     1524#endif
     1525bool ClpFactorization::timeToRefactorize() const
    15451526{
    15461527  if (coinFactorizationA_) {
    1547     bool reFactor = (coinFactorizationA_->pivots() * 3 > coinFactorizationA_->maximumPivots() * 2 &&
    1548                      coinFactorizationA_->numberElementsR() * 3 > (coinFactorizationA_->numberElementsL() +
    1549                                                                    coinFactorizationA_->numberElementsU()) * 2 + 1000 &&
    1550                      !coinFactorizationA_->numberDense());
    1551     reFactor=false;
    1552     bool reFactor3=false;
    1553     int numberPivots=coinFactorizationA_->pivots();
     1528    bool reFactor = (coinFactorizationA_->pivots() * 3 > coinFactorizationA_->maximumPivots() * 2 && coinFactorizationA_->numberElementsR() * 3 > (coinFactorizationA_->numberElementsL() + coinFactorizationA_->numberElementsU()) * 2 + 1000 && !coinFactorizationA_->numberDense());
     1529    reFactor = false;
     1530    bool reFactor3 = false;
     1531    int numberPivots = coinFactorizationA_->pivots();
    15541532    //if (coinFactorizationA_->pivots()<2)
    1555     if (numberPivots>lastNumberPivots_) {
     1533    if (numberPivots > lastNumberPivots_) {
    15561534      if (!lastNumberPivots_) {
    15571535        //lastR=0;
    1558         //lastU=endLengthU;
    1559         totalInR_=0.0;
    1560         totalInIncreasingU_=0.0;
    1561         shortestAverage_=COIN_DBL_MAX;
     1536        //lastU=endLengthU;
     1537        totalInR_ = 0.0;
     1538        totalInIncreasingU_ = 0.0;
     1539        shortestAverage_ = COIN_DBL_MAX;
    15621540#ifdef CLP_USEFUL_PRINTOUT
    1563         if (!readTwiddle) {
    1564           readTwiddle=true;
    1565           char * environ = getenv("CLP_TWIDDLE");
    1566           if (environ) {
    1567             sscanf(environ,"%lg %lg %lg %lg %lg %lg %lg %lg %lg",
    1568       &weightIncU,&weightR,&weightRest,&weightFactL,
    1569       &weightFactDense,&weightNrows,&increaseNeeded,
    1570       &constWeightIterate,&weightNrowsIterate);
    1571           }
    1572             printf("weightIncU %g, weightR %g, weightRest %g, weightFactL %g, weightFactDense %g, weightNrows %g increaseNeeded %g constWeightIterate %g weightNrowsIterate %g\n",
    1573       weightIncU,weightR,weightRest,weightFactL,
    1574         weightFactDense,weightNrows,increaseNeeded,
    1575       constWeightIterate,weightNrowsIterate);
    1576         }
    1577 #endif
    1578       }
    1579       lastNumberPivots_=numberPivots;
    1580       int numberDense=coinFactorizationA_->numberDense();
    1581       double nnd=numberDense*numberDense;
    1582       int lengthL=coinFactorizationA_->numberElementsL();
    1583       int lengthR=coinFactorizationA_->numberElementsR();
     1541        if (!readTwiddle) {
     1542          readTwiddle = true;
     1543          char *environ = getenv("CLP_TWIDDLE");
     1544          if (environ) {
     1545            sscanf(environ, "%lg %lg %lg %lg %lg %lg %lg %lg %lg",
     1546              &weightIncU, &weightR, &weightRest, &weightFactL,
     1547              &weightFactDense, &weightNrows, &increaseNeeded,
     1548              &constWeightIterate, &weightNrowsIterate);
     1549          }
     1550          printf("weightIncU %g, weightR %g, weightRest %g, weightFactL %g, weightFactDense %g, weightNrows %g increaseNeeded %g constWeightIterate %g weightNrowsIterate %g\n",
     1551            weightIncU, weightR, weightRest, weightFactL,
     1552            weightFactDense, weightNrows, increaseNeeded,
     1553            constWeightIterate, weightNrowsIterate);
     1554        }
     1555#endif
     1556      }
     1557      lastNumberPivots_ = numberPivots;
     1558      int numberDense = coinFactorizationA_->numberDense();
     1559      double nnd = numberDense * numberDense;
     1560      int lengthL = coinFactorizationA_->numberElementsL();
     1561      int lengthR = coinFactorizationA_->numberElementsR();
    15841562      int numberRows = coinFactorizationA_->numberRows();
    1585       int lengthU=coinFactorizationA_->numberElementsU()-
    1586         (numberRows-numberDense);
     1563      int lengthU = coinFactorizationA_->numberElementsU() - (numberRows - numberDense);
    15871564      totalInR_ += lengthR;
    1588       int effectiveU=lengthU-effectiveStartNumberU_;
     1565      int effectiveU = lengthU - effectiveStartNumberU_;
    15891566      totalInIncreasingU_ += effectiveU;
    15901567      //lastR=lengthR;
    15911568      //lastU=lengthU;
    1592       double rest=lengthL+0.05*nnd;
    1593       double constWeightFactor = weightFactL*lengthL+weightFactDense*nnd
    1594         + weightNrows*numberRows;
    1595       double constWeightIterateX = constWeightIterate*(lengthL+endLengthU_)
    1596         + weightNrowsIterate*numberRows;
    1597       double variableWeight = weightIncU*totalInIncreasingU_+
    1598                                weightR*totalInR_+weightRest*rest;
    1599       double average=constWeightIterateX+
    1600         (constWeightFactor+variableWeight)/static_cast<double>(numberPivots);
     1569      double rest = lengthL + 0.05 * nnd;
     1570      double constWeightFactor = weightFactL * lengthL + weightFactDense * nnd
     1571        + weightNrows * numberRows;
     1572      double constWeightIterateX = constWeightIterate * (lengthL + endLengthU_)
     1573        + weightNrowsIterate * numberRows;
     1574      double variableWeight = weightIncU * totalInIncreasingU_ + weightR * totalInR_ + weightRest * rest;
     1575      double average = constWeightIterateX + (constWeightFactor + variableWeight) / static_cast< double >(numberPivots);
    16011576#if 0
    16021577      if ((numberPivots%20)==0&&!ifPrint3)
     
    16051580      lengthU,lengthL,lengthR,nnd,average);
    16061581#endif
    1607       shortestAverage_=CoinMin(shortestAverage_,average);
    1608       if (average>increaseNeeded*shortestAverage_&&
    1609         coinFactorizationA_->pivots()>30) {
    1610         //printf("PIVX %d nrow %d startU %d now %d L %d R %d dense %g average %g\n",
    1611         //numberPivots,numberRows,effectiveStartNumberU_,
    1612         //lengthU,lengthL,lengthR,nnd,average);
    1613         reFactor3=true;
    1614       }
    1615   }
    1616     if (reFactor|| reFactor3) {
     1582      shortestAverage_ = CoinMin(shortestAverage_, average);
     1583      if (average > increaseNeeded * shortestAverage_ && coinFactorizationA_->pivots() > 30) {
     1584        //printf("PIVX %d nrow %d startU %d now %d L %d R %d dense %g average %g\n",
     1585        //numberPivots,numberRows,effectiveStartNumberU_,
     1586        //lengthU,lengthL,lengthR,nnd,average);
     1587        reFactor3 = true;
     1588      }
     1589    }
     1590    if (reFactor || reFactor3) {
    16171591#if 0
    16181592        printf("%c%cftranCountInput_ %g  ,ftranCountAfterL_ %g  ,ftranCountAfterR_ %g  ,ftranCountAfterU_ %g  ,btranCountInput_ %g  ,btranCountAfterU_ %g  ,btranCountAfterR_ %g  ,btranCountAfterL_ %g  ,numberFtranCounts_ %d  ,numberBtranCounts_ %d  ,ftranAverageAfterL_ %g  ,ftranAverageAfterR_ %g  ,ftranAverageAfterU_ %g  ,btranAverageAfterU_ %g  ,btranAverageAfterR_ %g  ,btranAverageAfterL_ %g\n"
     
    16361610               ,coinFactorizationA_->btranAverageAfterL_);
    16371611#endif
    1638       reFactor=true;
     1612      reFactor = true;
    16391613    }
    16401614    return reFactor;
     
    16431617  }
    16441618}
    1645 #if CLP_FACTORIZATION_NEW_TIMING>1
    1646 void
    1647 ClpFactorization::statsRefactor(char when) const
    1648 {
    1649   int numberPivots=coinFactorizationA_->pivots();
    1650   int numberDense=coinFactorizationA_->numberDense();
    1651   double nnd=numberDense*numberDense;
    1652   int lengthL=coinFactorizationA_->numberElementsL();
    1653   int lengthR=coinFactorizationA_->numberElementsR();
     1619#if CLP_FACTORIZATION_NEW_TIMING > 1
     1620void ClpFactorization::statsRefactor(char when) const
     1621{
     1622  int numberPivots = coinFactorizationA_->pivots();
     1623  int numberDense = coinFactorizationA_->numberDense();
     1624  double nnd = numberDense * numberDense;
     1625  int lengthL = coinFactorizationA_->numberElementsL();
     1626  int lengthR = coinFactorizationA_->numberElementsR();
    16541627  int numberRows = coinFactorizationA_->numberRows();
    1655   int lengthU=coinFactorizationA_->numberElementsU()-
    1656     (numberRows-numberDense);
    1657   double rest=lengthL+0.05*nnd;
    1658   double constWeightFactor = weightFactL*lengthL+weightFactDense*nnd
    1659     + weightNrows*numberRows;
    1660   double constWeightIterateX = constWeightIterate*(lengthL+endLengthU_)
    1661     + weightNrowsIterate*numberRows;
    1662   double variableWeight = weightIncU*totalInIncreasingU_+
    1663     weightR*totalInR_+weightRest*rest;
    1664   double average=constWeightIterateX+
    1665     (constWeightFactor+variableWeight)/static_cast<double>(numberPivots);
     1628  int lengthU = coinFactorizationA_->numberElementsU() - (numberRows - numberDense);
     1629  double rest = lengthL + 0.05 * nnd;
     1630  double constWeightFactor = weightFactL * lengthL + weightFactDense * nnd
     1631    + weightNrows * numberRows;
     1632  double constWeightIterateX = constWeightIterate * (lengthL + endLengthU_)
     1633    + weightNrowsIterate * numberRows;
     1634  double variableWeight = weightIncU * totalInIncreasingU_ + weightR * totalInR_ + weightRest * rest;
     1635  double average = constWeightIterateX + (constWeightFactor + variableWeight) / static_cast< double >(numberPivots);
    16661636  printf("PIV%c %d nrow %d startU %d now %d L %d R %d dense %g average %g - shortest %g\n",
    1667          when,numberPivots,numberRows,effectiveStartNumberU_,
    1668          lengthU,lengthL,lengthR,nnd,average,shortestAverage_);
     1637    when, numberPivots, numberRows, effectiveStartNumberU_,
     1638    lengthU, lengthL, lengthR, nnd, average, shortestAverage_);
    16691639}
    16701640#endif
    16711641#else
    1672 bool
    1673 ClpFactorization::timeToRefactorize() const
    1674 {
    1675     if (coinFactorizationA_) {
    1676     return (coinFactorizationA_->pivots() * 3 > coinFactorizationA_->maximumPivots() * 2 &&
    1677       coinFactorizationA_->numberElementsR() * 3 > (coinFactorizationA_->numberElementsL() +
    1678       coinFactorizationA_->numberElementsU()) * 2 + 1000 &&
    1679       !coinFactorizationA_->numberDense());
     1642bool ClpFactorization::timeToRefactorize() const
     1643{
     1644  if (coinFactorizationA_) {
     1645    return (coinFactorizationA_->pivots() * 3 > coinFactorizationA_->maximumPivots() * 2 && coinFactorizationA_->numberElementsR() * 3 > (coinFactorizationA_->numberElementsL() + coinFactorizationA_->numberElementsU()) * 2 + 1000 && !coinFactorizationA_->numberDense());
    16801646  } else {
    16811647    return coinFactorizationB_->pivots() > coinFactorizationB_->numberRows() / 2.45 + 20;
    16821648  }
    16831649#endif
    1684 int
    1685 ClpFactorization::factorize ( ClpSimplex * model,
    1686                               int solveType, bool valuesPass)
     1650int ClpFactorization::factorize(ClpSimplex *model,
     1651  int solveType, bool valuesPass)
    16871652{
    16881653#ifdef CLP_REUSE_ETAS
    1689      model_= model;
    1690 #endif
    1691      //if ((model->specialOptions()&16384))
    1692      //printf("factor at %d iterations\n",model->numberIterations());
    1693      ClpMatrixBase * matrix = model->clpMatrix();
    1694      int numberRows = model->numberRows();
    1695      int numberColumns = model->numberColumns();
    1696      if (!numberRows)
    1697           return 0;
    1698 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1699      factorization_instrument(-1);
    1700      if (!timeInUpdate) {
    1701       printf("ZZZZ,timeInFactor,nRows,nrows_2,\
     1654  model_ = model;
     1655#endif
     1656  //if ((model->specialOptions()&16384))
     1657  //printf("factor at %d iterations\n",model->numberIterations());
     1658  ClpMatrixBase *matrix = model->clpMatrix();
     1659  int numberRows = model->numberRows();
     1660  int numberColumns = model->numberColumns();
     1661  if (!numberRows)
     1662    return 0;
     1663#ifdef CLP_FACTORIZATION_INSTRUMENT
     1664  factorization_instrument(-1);
     1665  if (!timeInUpdate) {
     1666    printf("ZZZZ,timeInFactor,nRows,nrows_2,\
    17021667startU,endU,lengthL,dense,dense_2,nslacks\n");
    1703       printf("YYYY,time,added,dense,dense_2,averageR,averageL,averageU,\
     1668    printf("YYYY,time,added,dense,dense_2,averageR,averageL,averageU,\
    17041669average2R,average2L,average2U,\
    17051670scaledDense,scaledDense_2,scaledL,scaledR,scaledU\n");
    1706       printf("WWWW,time,size0,ratio1,ratio05,ratio02\n");
    1707     }
    1708 #endif
    1709 #ifdef CLP_FACTORIZATION_INSTRUMENT
    1710      externalTimeStart=CoinCpuTime();
    1711      memset(average,0,sizeof(average));
    1712      shortest=COIN_DBL_MAX;
    1713      ifPrint=0;
    1714      if (!model->numberIterations())
    1715        printf("maxfactor %d\n",coinFactorizationA_->maximumPivots());
    1716 #endif
    1717      bool anyChanged = false;
    1718      if (coinFactorizationB_) {
    1719           coinFactorizationB_->setStatus(-99);
    1720           int * pivotVariable = model->pivotVariable();
    1721           //returns 0 -okay, -1 singular, -2 too many in basis */
    1722           // allow dense
    1723           int solveMode = 2;
    1724           if (model->numberIterations())
    1725                solveMode += 8;
    1726           if (valuesPass)
    1727                solveMode += 4;
    1728           coinFactorizationB_->setSolveMode(solveMode);
    1729           while (status() < -98) {
     1671    printf("WWWW,time,size0,ratio1,ratio05,ratio02\n");
     1672  }
     1673#endif
     1674#ifdef CLP_FACTORIZATION_INSTRUMENT
     1675  externalTimeStart = CoinCpuTime();
     1676  memset(average, 0, sizeof(average));
     1677  shortest = COIN_DBL_MAX;
     1678  ifPrint = 0;
     1679  if (!model->numberIterations())
     1680    printf("maxfactor %d\n", coinFactorizationA_->maximumPivots());
     1681#endif
     1682  bool anyChanged = false;
     1683  if (coinFactorizationB_) {
     1684    coinFactorizationB_->setStatus(-99);
     1685    int *pivotVariable = model->pivotVariable();
     1686    //returns 0 -okay, -1 singular, -2 too many in basis */
     1687    // allow dense
     1688    int solveMode = 2;
     1689    if (model->numberIterations())
     1690      solveMode += 8;
     1691    if (valuesPass)
     1692      solveMode += 4;
     1693    coinFactorizationB_->setSolveMode(solveMode);
     1694    while (status() < -98) {
    17301695
    1731                int i;
    1732                int numberBasic = 0;
    1733                int numberRowBasic;
    1734                // Move pivot variables across if they look good
    1735                int * pivotTemp = model->rowArray(0)->getIndices();
    1736                assert (!model->rowArray(0)->getNumElements());
    1737                if (!matrix->rhsOffset(model)) {
    1738                     // Seems to prefer things in order so quickest
    1739                     // way is to go though like this
    1740                     for (i = 0; i < numberRows; i++) {
    1741                          if (model->getRowStatus(i) == ClpSimplex::basic)
    1742                               pivotTemp[numberBasic++] = i;
    1743                     }
    1744                     numberRowBasic = numberBasic;
    1745                     /* Put column basic variables into pivotVariable
     1696      int i;
     1697      int numberBasic = 0;
     1698      int numberRowBasic;
     1699      // Move pivot variables across if they look good
     1700      int *pivotTemp = model->rowArray(0)->getIndices();
     1701      assert(!model->rowArray(0)->getNumElements());
     1702      if (!matrix->rhsOffset(model)) {
     1703        // Seems to prefer things in order so quickest
     1704        // way is to go though like this
     1705        for (i = 0; i < numberRows; i++) {
     1706          if (model->getRowStatus(i) == ClpSimplex::basic)
     1707            pivotTemp[numberBasic++] = i;
     1708        }
     1709        numberRowBasic = numberBasic;
     1710        /* Put column basic variables into pivotVariable
    17461711                       This is done by ClpMatrixBase to allow override for gub
    17471712                    */
    1748                     matrix->generalExpanded(model, 0, numberBasic);
    1749                } else {
    1750                     // Long matrix - do a different way
    1751                     bool fullSearch = false;
    1752                     for (i = 0; i < numberRows; i++) {
    1753                          int iPivot = pivotVariable[i];
    1754                          if (iPivot >= numberColumns) {
    1755                               pivotTemp[numberBasic++] = iPivot - numberColumns;
    1756                          }
    1757                     }
    1758                     numberRowBasic = numberBasic;
    1759                     for (i = 0; i < numberRows; i++) {
    1760                          int iPivot = pivotVariable[i];
    1761                          if (iPivot < numberColumns) {
    1762                               if (iPivot >= 0) {
    1763                                    pivotTemp[numberBasic++] = iPivot;
    1764                               } else {
    1765                                    // not full basis
    1766                                    fullSearch = true;
    1767                                    break;
    1768                               }
    1769                          }
    1770                     }
    1771                     if (fullSearch) {
    1772                          // do slow way
    1773                          numberBasic = 0;
    1774                          for (i = 0; i < numberRows; i++) {
    1775                               if (model->getRowStatus(i) == ClpSimplex::basic)
    1776                                    pivotTemp[numberBasic++] = i;
    1777                          }
    1778                          numberRowBasic = numberBasic;
    1779                          /* Put column basic variables into pivotVariable
     1713        matrix->generalExpanded(model, 0, numberBasic);
     1714      } else {
     1715        // Long matrix - do a different way
     1716        bool fullSearch = false;
     1717        for (i = 0; i < numberRows; i++) {
     1718          int iPivot = pivotVariable[i];
     1719          if (iPivot >= numberColumns) {
     1720            pivotTemp[numberBasic++] = iPivot - numberColumns;
     1721          }
     1722        }
     1723        numberRowBasic = numberBasic;
     1724        for (i = 0; i < numberRows; i++) {
     1725          int iPivot = pivotVariable[i];
     1726          if (iPivot < numberColumns) {
     1727            if (iPivot >= 0) {
     1728              pivotTemp[numberBasic++] = iPivot;
     1729            } else {
     1730              // not full basis
     1731              fullSearch = true;
     1732              break;
     1733            }
     1734          }
     1735        }
     1736        if (fullSearch) {
     1737          // do slow way
     1738          numberBasic = 0;
     1739          for (i = 0; i < numberRows; i++) {
     1740            if (model->getRowStatus(i) == ClpSimplex::basic)
     1741              pivotTemp[numberBasic++] = i;
     1742          }
     1743          numberRowBasic = numberBasic;
     1744          /* Put column basic variables into pivotVariable
    17801745                            This is done by ClpMatrixBase to allow override for gub
    17811746                         */
    1782                          matrix->generalExpanded(model, 0, numberBasic);
    1783                     }
    1784                }
    1785                if (numberBasic > model->maximumBasic()) {
    1786                     // Take out some
    1787                     numberBasic = numberRowBasic;
    1788                     for (int i = 0; i < numberColumns; i++) {
    1789                          if (model->getColumnStatus(i) == ClpSimplex::basic) {
    1790                               if (numberBasic < numberRows)
    1791                                    numberBasic++;
    1792                               else
    1793                                    model->setColumnStatus(i, ClpSimplex::superBasic);
    1794                          }
    1795                     }
    1796                     numberBasic = numberRowBasic;
    1797                     matrix->generalExpanded(model, 0, numberBasic);
    1798                } else if (numberBasic < numberRows) {
    1799                     // add in some
    1800                     int needed = numberRows - numberBasic;
    1801                     // move up columns
    1802                     for (i = numberBasic - 1; i >= numberRowBasic; i--)
    1803                          pivotTemp[i+needed] = pivotTemp[i];
    1804                     numberRowBasic = 0;
    1805                     numberBasic = numberRows;
    1806                     for (i = 0; i < numberRows; i++) {
    1807                          if (model->getRowStatus(i) == ClpSimplex::basic) {
    1808                               pivotTemp[numberRowBasic++] = i;
    1809                          } else if (needed) {
    1810                               needed--;
    1811                               model->setRowStatus(i, ClpSimplex::basic);
    1812                               pivotTemp[numberRowBasic++] = i;
    1813                          }
    1814                     }
    1815                }
    1816                int numberElements = numberRowBasic;
     1747          matrix->generalExpanded(model, 0, numberBasic);
     1748        }
     1749      }
     1750      if (numberBasic > model->maximumBasic()) {
     1751        // Take out some
     1752        numberBasic = numberRowBasic;
     1753        for (int i = 0; i < numberColumns; i++) {
     1754          if (model->getColumnStatus(i) == ClpSimplex::basic) {
     1755            if (numberBasic < numberRows)
     1756              numberBasic++;
     1757            else
     1758              model->setColumnStatus(i, ClpSimplex::superBasic);
     1759          }
     1760        }
     1761        numberBasic = numberRowBasic;
     1762        matrix->generalExpanded(model, 0, numberBasic);
     1763      } else if (numberBasic < numberRows) {
     1764        // add in some
     1765        int needed = numberRows - numberBasic;
     1766        // move up columns
     1767        for (i = numberBasic - 1; i >= numberRowBasic; i--)
     1768          pivotTemp[i + needed] = pivotTemp[i];
     1769        numberRowBasic = 0;
     1770        numberBasic = numberRows;
     1771        for (i = 0; i < numberRows; i++) {
     1772          if (model->getRowStatus(i) == ClpSimplex::basic) {
     1773            pivotTemp[numberRowBasic++] = i;
     1774          } else if (needed) {
     1775            needed--;
     1776            model->setRowStatus(i, ClpSimplex::basic);
     1777            pivotTemp[numberRowBasic++] = i;
     1778          }
     1779        }
     1780      }
     1781      int numberElements = numberRowBasic;
    18171782
    1818                // compute how much in basis
    1819                // can change for gub
    1820                int numberColumnBasic = numberBasic - numberRowBasic;
     1783      // compute how much in basis
     1784      // can change for gub
     1785      int numberColumnBasic = numberBasic - numberRowBasic;
    18211786
    1822                numberElements += matrix->countBasis(pivotTemp + numberRowBasic,
    1823                                                     numberColumnBasic);
     1787      numberElements += matrix->countBasis(pivotTemp + numberRowBasic,
     1788        numberColumnBasic);
    18241789#ifndef NDEBUG
    18251790//#define CHECK_CLEAN_BASIS
    18261791#ifdef CHECK_CLEAN_BASIS
    1827                int saveNumberElements = numberElements;
    1828 #endif
    1829 #endif
    1830                // Not needed for dense
    1831                numberElements = 3 * numberBasic + 3 * numberElements + 20000;
    1832                int numberIterations = model->numberIterations();
    1833                coinFactorizationB_->setUsefulInformation(&numberIterations, 0);
    1834                coinFactorizationB_->getAreas ( numberRows,
    1835                                                numberRowBasic + numberColumnBasic, numberElements,
    1836                                                2 * numberElements );
    1837                // Fill in counts so we can skip part of preProcess
    1838                // This is NOT needed for dense but would be needed for later versions
    1839                CoinFactorizationDouble * elementU;
    1840                int * indexRowU;
    1841                int * startColumnU;
    1842                int * numberInRow;
    1843                int * numberInColumn;
    1844                elementU = coinFactorizationB_->elements();
    1845                indexRowU = coinFactorizationB_->indices();
    1846                startColumnU = coinFactorizationB_->starts();
     1792      int saveNumberElements = numberElements;
     1793#endif
     1794#endif
     1795      // Not needed for dense
     1796      numberElements = 3 * numberBasic + 3 * numberElements + 20000;
     1797      int numberIterations = model->numberIterations();
     1798      coinFactorizationB_->setUsefulInformation(&numberIterations, 0);
     1799      coinFactorizationB_->getAreas(numberRows,
     1800        numberRowBasic + numberColumnBasic, numberElements,
     1801        2 * numberElements);
     1802      // Fill in counts so we can skip part of preProcess
     1803      // This is NOT needed for dense but would be needed for later versions
     1804      CoinFactorizationDouble *elementU;
     1805      int *indexRowU;
     1806      int *startColumnU;
     1807      int *numberInRow;
     1808      int *numberInColumn;
     1809      elementU = coinFactorizationB_->elements();
     1810      indexRowU = coinFactorizationB_->indices();
     1811      startColumnU = coinFactorizationB_->starts();
    18471812#ifdef CHECK_CLEAN_BASIS
    1848                for (int i=0;i<saveNumberElements;i++) {
    1849                  elementU[i]=0.0;
    1850                  indexRowU[i]=-1;
    1851                }
    1852                for (int i=0;i<numberRows;i++)
    1853                  startColumnU[i]=-1;
     1813      for (int i = 0; i < saveNumberElements; i++) {
     1814        elementU[i] = 0.0;
     1815        indexRowU[i] = -1;
     1816      }
     1817      for (int i = 0; i < numberRows; i++)
     1818        startColumnU[i] = -1;
    18541819#endif
    18551820#ifndef COIN_FAST_CODE
    1856                double slackValue;
    1857                slackValue = coinFactorizationB_->slackValue();
     1821      double slackValue;
     1822      slackValue = coinFactorizationB_->slackValue();
    18581823#else
    18591824#define slackValue -1.0
    18601825#endif
    1861                numberInRow = coinFactorizationB_->numberInRow();
    1862                numberInColumn = coinFactorizationB_->numberInColumn();
    1863                CoinZeroN ( numberInRow, numberRows  );
    1864                CoinZeroN ( numberInColumn, numberRows );
    1865                for (i = 0; i < numberRowBasic; i++) {
    1866                     int iRow = pivotTemp[i];
    1867                     // Change pivotTemp to correct sequence
    1868                     pivotTemp[i] = iRow + numberColumns;
    1869                     indexRowU[i] = iRow;
    1870                     startColumnU[i] = i;
    1871                     elementU[i] = slackValue;
    1872                     numberInRow[iRow] = 1;
    1873                     numberInColumn[i] = 1;
    1874                }
    1875                startColumnU[numberRowBasic] = numberRowBasic;
    1876                // can change for gub so redo
    1877                numberColumnBasic = numberBasic - numberRowBasic;
    1878                matrix->fillBasis(model,
    1879                                  pivotTemp + numberRowBasic,
    1880                                  numberColumnBasic,
    1881                                  indexRowU,
    1882                                  startColumnU + numberRowBasic,
    1883                                  numberInRow,
    1884                                  numberInColumn + numberRowBasic,
    1885                                  elementU);
    1886                // recompute number basic
    1887                numberBasic = numberRowBasic + numberColumnBasic;
    1888                for (i = numberBasic; i < numberRows; i++)
    1889                     pivotTemp[i] = -1; // mark not there
    1890                if (numberBasic)
    1891                     numberElements = startColumnU[numberBasic-1]
    1892                                      + numberInColumn[numberBasic-1];
    1893                else
    1894                     numberElements = 0;
     1826      numberInRow = coinFactorizationB_->numberInRow();
     1827      numberInColumn = coinFactorizationB_->numberInColumn();
     1828      CoinZeroN(numberInRow, numberRows);
     1829      CoinZeroN(numberInColumn, numberRows);
     1830      for (i = 0; i < numberRowBasic; i++) {
     1831        int iRow = pivotTemp[i];
     1832        // Change pivotTemp to correct sequence
     1833        pivotTemp[i] = iRow + numberColumns;
     1834        indexRowU[i] = iRow;
     1835        startColumnU[i] = i;
     1836        elementU[i] = slackValue;
     1837        numberInRow[iRow] = 1;
     1838        numberInColumn[i] = 1;
     1839      }
     1840      startColumnU[numberRowBasic] = numberRowBasic;
     1841      // can change for gub so redo
     1842      numberColumnBasic = numberBasic - numberRowBasic;
     1843      matrix->fillBasis(model,
     1844        pivotTemp + numberRowBasic,
     1845        numberColumnBasic,
     1846        indexRowU,
     1847        startColumnU + numberRowBasic,
     1848        numberInRow,
     1849        numberInColumn + numberRowBasic,
     1850        elementU);
     1851      // recompute number basic
     1852      numberBasic = numberRowBasic + numberColumnBasic;
     1853      for (i = numberBasic; i < numberRows; i++)
     1854        pivotTemp[i] = -1; // mark not there
     1855      if (numberBasic)
     1856        numberElements = startColumnU[numberBasic - 1]
     1857          + numberInColumn[numberBasic - 1];
     1858      else
     1859        numberElements = 0;
    18951860#ifdef CHECK_CLEAN_BASIS
    1896                assert (!startColumnU[0]);
    1897                int lastStart=0;
    1898                for (int i=0;i<numberRows;i++) {
    1899                  assert (startColumnU[i+1]>lastStart);
    1900                  lastStart=startColumnU[i+1];
    1901                }
    1902                assert (lastStart==saveNumberElements);
    1903                for (int i=0;i<saveNumberElements;i++) {
    1904                  assert(elementU[i]);
    1905                  assert(indexRowU[i]>=0&&indexRowU[i]<numberRows);
    1906                }
    1907 #endif
    1908                coinFactorizationB_->preProcess ( );
    1909                coinFactorizationB_->factor (  );
    1910                if (coinFactorizationB_->status() == -1 &&
    1911                          (coinFactorizationB_->solveMode() % 3) != 0) {
    1912                     int solveMode = coinFactorizationB_->solveMode();
    1913                     solveMode -= solveMode % 3; // so bottom will be 0
    1914                     coinFactorizationB_->setSolveMode(solveMode);
    1915                     coinFactorizationB_->setStatus(-99);
    1916                }
    1917                if (coinFactorizationB_->status() == -99)
    1918                     continue;
    1919                // If we get here status is 0 or -1
    1920                if (coinFactorizationB_->status() == 0 && numberBasic == numberRows) {
    1921                     coinFactorizationB_->postProcess(pivotTemp, pivotVariable);
    1922                } else if (solveType == 0 || solveType == 2/*||solveType==1*/) {
    1923                     // Change pivotTemp to be correct list
    1924                     anyChanged = true;
    1925                     coinFactorizationB_->makeNonSingular(pivotTemp, numberColumns);
    1926                     double * columnLower = model->lowerRegion();
    1927                     double * columnUpper = model->upperRegion();
    1928                     double * columnActivity = model->solutionRegion();
    1929                     double * rowLower = model->lowerRegion(0);
    1930                     double * rowUpper = model->upperRegion(0);
    1931                     double * rowActivity = model->solutionRegion(0);
    1932                     //redo basis - first take ALL out
    1933                     int iColumn;
    1934                     double largeValue = model->largeValue();
    1935                     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1936                          if (model->getColumnStatus(iColumn) == ClpSimplex::basic) {
    1937                               // take out
    1938                               if (!valuesPass) {
    1939                                    double lower = columnLower[iColumn];
    1940                                    double upper = columnUpper[iColumn];
    1941                                    double value = columnActivity[iColumn];
    1942                                    if (lower > -largeValue || upper < largeValue) {
    1943                                         if (fabs(value - lower) < fabs(value - upper)) {
    1944                                              model->setColumnStatus(iColumn, ClpSimplex::atLowerBound);
    1945                                              columnActivity[iColumn] = lower;
    1946                                         } else {
    1947                                              model->setColumnStatus(iColumn, ClpSimplex::atUpperBound);
    1948                                              columnActivity[iColumn] = upper;
    1949                                         }
    1950                                    } else {
    1951                                         model->setColumnStatus(iColumn, ClpSimplex::isFree);
    1952                                    }
    1953                               } else {
    1954                                    model->setColumnStatus(iColumn, ClpSimplex::superBasic);
    1955                               }
    1956                          }
    1957                     }
    1958                     int iRow;
    1959                     for (iRow = 0; iRow < numberRows; iRow++) {
    1960                          if (model->getRowStatus(iRow) == ClpSimplex::basic) {
    1961                               // take out
    1962                               if (!valuesPass) {
    1963                                    double lower = columnLower[iRow];
    1964                                    double upper = columnUpper[iRow];
    1965                                    double value = columnActivity[iRow];
    1966                                    if (lower > -largeValue || upper < largeValue) {
    1967                                         if (fabs(value - lower) < fabs(value - upper)) {
    1968                                              model->setRowStatus(iRow, ClpSimplex::atLowerBound);
    1969                                              columnActivity[iRow] = lower;
    1970                                         } else {
    1971                                              model->setRowStatus(iRow, ClpSimplex::atUpperBound);
    1972                                              columnActivity[iRow] = upper;
    1973                                         }
    1974                                    } else {
    1975                                         model->setRowStatus(iRow, ClpSimplex::isFree);
    1976                                    }
    1977                               } else {
    1978                                    model->setRowStatus(iRow, ClpSimplex::superBasic);
    1979                               }
    1980                          }
    1981                     }
    1982                     for (iRow = 0; iRow < numberRows; iRow++) {
    1983                          int iSequence = pivotTemp[iRow];
    1984                          assert (iSequence >= 0);
    1985                          // basic
    1986                          model->setColumnStatus(iSequence, ClpSimplex::basic);
    1987                     }
    1988                     // signal repeat
    1989                     coinFactorizationB_->setStatus(-99);
    1990                     // set fixed if they are
    1991                     for (iRow = 0; iRow < numberRows; iRow++) {
    1992                          if (model->getRowStatus(iRow) != ClpSimplex::basic ) {
    1993                               if (rowLower[iRow] == rowUpper[iRow]) {
    1994                                    rowActivity[iRow] = rowLower[iRow];
    1995                                    model->setRowStatus(iRow, ClpSimplex::isFixed);
    1996                               }
    1997                          }
    1998                     }
    1999                     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    2000                          if (model->getColumnStatus(iColumn) != ClpSimplex::basic ) {
    2001                               if (columnLower[iColumn] == columnUpper[iColumn]) {
    2002                                    columnActivity[iColumn] = columnLower[iColumn];
    2003                                    model->setColumnStatus(iColumn, ClpSimplex::isFixed);
    2004                               }
    2005                          }
    2006                     }
    2007                }
    2008           }
     1861      assert(!startColumnU[0]);
     1862      int lastStart = 0;
     1863      for (int i = 0; i < numberRows; i++) {
     1864        assert(startColumnU[i + 1] > lastStart);
     1865        lastStart = startColumnU[i + 1];
     1866      }
     1867      assert(lastStart == saveNumberElements);
     1868      for (int i = 0; i < saveNumberElements; i++) {
     1869        assert(elementU[i]);
     1870        assert(indexRowU[i] >= 0 && indexRowU[i] < numberRows);
     1871      }
     1872#endif
     1873      coinFactorizationB_->preProcess();
     1874      coinFactorizationB_->factor();
     1875      if (coinFactorizationB_->status() == -1 && (coinFactorizationB_->solveMode() % 3) != 0) {
     1876        int solveMode = coinFactorizationB_->solveMode();
     1877        solveMode -= solveMode % 3; // so bottom will be 0
     1878        coinFactorizationB_->setSolveMode(solveMode);
     1879        coinFactorizationB_->setStatus(-99);
     1880      }
     1881      if (coinFactorizationB_->status() == -99)
     1882        continue;
     1883      // If we get here status is 0 or -1
     1884      if (coinFactorizationB_->status() == 0 && numberBasic == numberRows) {
     1885        coinFactorizationB_->postProcess(pivotTemp, pivotVariable);
     1886      } else if (solveType == 0 || solveType == 2 /*||solveType==1*/) {
     1887        // Change pivotTemp to be correct list
     1888        anyChanged = true;
     1889        coinFactorizationB_->makeNonSingular(pivotTemp, numberColumns);
     1890        double *columnLower = model->lowerRegion();
     1891        double *columnUpper = model->upperRegion();
     1892        double *columnActivity = model->solutionRegion();
     1893        double *rowLower = model->lowerRegion(0);
     1894        double *rowUpper = model->upperRegion(0);
     1895        double *rowActivity = model->solutionRegion(0);
     1896        //redo basis - first take ALL out
     1897        int iColumn;
     1898        double largeValue = model->largeValue();
     1899        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1900          if (model->getColumnStatus(iColumn) == ClpSimplex::basic) {
     1901            // take out
     1902            if (!valuesPass) {
     1903              double lower = columnLower[iColumn];
     1904              double upper = columnUpper[iColumn];
     1905              double value = columnActivity[iColumn];
     1906              if (lower > -largeValue || upper < largeValue) {
     1907                if (fabs(value - lower) < fabs(value - upper)) {
     1908                  model->setColumnStatus(iColumn, ClpSimplex::atLowerBound);
     1909                  columnActivity[iColumn] = lower;
     1910                } else {
     1911                  model->setColumnStatus(iColumn, ClpSimplex::atUpperBound);
     1912                  columnActivity[iColumn] = upper;
     1913                }
     1914              } else {
     1915                model->setColumnStatus(iColumn, ClpSimplex::isFree);
     1916              }
     1917            } else {
     1918              model->setColumnStatus(iColumn, ClpSimplex::superBasic);
     1919            }
     1920          }
     1921        }
     1922        int iRow;
     1923        for (iRow = 0; iRow < numberRows; iRow++) {
     1924          if (model->getRowStatus(iRow) == ClpSimplex::basic) {
     1925            // take out
     1926            if (!valuesPass) {
     1927              double lower = columnLower[iRow];
     1928              double upper = columnUpper[iRow];
     1929              double value = columnActivity[iRow];
     1930              if (lower > -largeValue || upper < largeValue) {
     1931                if (fabs(value - lower) < fabs(value - upper)) {
     1932                  model->setRowStatus(iRow, ClpSimplex::atLowerBound);
     1933                  columnActivity[iRow] = lower;
     1934                } else {
     1935                  model->setRowStatus(iRow, ClpSimplex::atUpperBound);
     1936                  columnActivity[iRow] = upper;
     1937                }
     1938              } else {
     1939                model->setRowStatus(iRow, ClpSimplex::isFree);
     1940              }
     1941            } else {
     1942              model->setRowStatus(iRow, ClpSimplex::superBasic);
     1943            }
     1944          }
     1945        }
     1946        for (iRow = 0; iRow < numberRows; iRow++) {
     1947          int iSequence = pivotTemp[iRow];
     1948          assert(iSequence >= 0);
     1949          // basic
     1950          model->setColumnStatus(iSequence, ClpSimplex::basic);
     1951        }
     1952        // signal repeat
     1953        coinFactorizationB_->setStatus(-99);
     1954        // set fixed if they are
     1955        for (iRow = 0; iRow < numberRows; iRow++) {
     1956          if (model->getRowStatus(iRow) != ClpSimplex::basic) {
     1957            if (rowLower[iRow] == rowUpper[iRow]) {
     1958              rowActivity[iRow] = rowLower[iRow];
     1959              model->setRowStatus(iRow, ClpSimplex::isFixed);
     1960            }
     1961          }
     1962        }
     1963        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1964          if (model->getColumnStatus(iColumn) != ClpSimplex::basic) {
     1965            if (columnLower[iColumn] == columnUpper[iColumn]) {
     1966              columnActivity[iColumn] = columnLower[iColumn];
     1967              model->setColumnStatus(iColumn, ClpSimplex::isFixed);
     1968            }
     1969          }
     1970        }
     1971      }
     1972    }
    20091973#ifdef CLP_DEBUG
    2010           // check basic
    2011           CoinIndexedVector region1(2 * numberRows);
    2012           CoinIndexedVector region2B(2 * numberRows);
    2013           int iPivot;
    2014           double * arrayB = region2B.denseVector();
    2015           int i;
    2016           for (iPivot = 0; iPivot < numberRows; iPivot++) {
    2017                int iSequence = pivotVariable[iPivot];
    2018                model->unpack(&region2B, iSequence);
    2019                coinFactorizationB_->updateColumn(&region1, &region2B);
    2020                if (fabs(arrayB[iPivot] - 1.0) < 1.0e-4) {
    2021                     // OK?
    2022                     arrayB[iPivot] = 0.0;
    2023                } else {
    2024                     assert (fabs(arrayB[iPivot]) < 1.0e-4);
    2025                     for (i = 0; i < numberRows; i++) {
    2026                          if (fabs(arrayB[i] - 1.0) < 1.0e-4)
    2027                               break;
    2028                     }
    2029                     assert (i < numberRows);
    2030                     printf("variable on row %d landed up on row %d\n", iPivot, i);
    2031                     arrayB[i] = 0.0;
    2032                }
    2033                for (i = 0; i < numberRows; i++)
    2034                     assert (fabs(arrayB[i]) < 1.0e-4);
    2035                region2B.clear();
    2036           }
    2037 #endif
    2038 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2039           factorization_instrument(2);
    2040 #endif
    2041           if ( anyChanged && model->algorithm() < 0 && solveType > 0) {
    2042                double dummyCost;
    2043                static_cast<ClpSimplexDual *> (model)->changeBounds(3,
    2044                          NULL, dummyCost);
    2045           }
    2046           return coinFactorizationB_->status();
    2047      }
    2048      // If too many compressions increase area
    2049      if (coinFactorizationA_->pivots() > 1 && coinFactorizationA_->numberCompressions() * 10 > coinFactorizationA_->pivots() + 10) {
    2050           coinFactorizationA_->areaFactor( coinFactorizationA_->areaFactor() * 1.1);
    2051      }
    2052      //int numberPivots=coinFactorizationA_->pivots();
     1974    // check basic
     1975    CoinIndexedVector region1(2 * numberRows);
     1976    CoinIndexedVector region2B(2 * numberRows);
     1977    int iPivot;
     1978    double *arrayB = region2B.denseVector();
     1979    int i;
     1980    for (iPivot = 0; iPivot < numberRows; iPivot++) {
     1981      int iSequence = pivotVariable[iPivot];
     1982      model->unpack(&region2B, iSequence);
     1983      coinFactorizationB_->updateColumn(&region1, &region2B);
     1984      if (fabs(arrayB[iPivot] - 1.0) < 1.0e-4) {
     1985        // OK?
     1986        arrayB[iPivot] = 0.0;
     1987      } else {
     1988        assert(fabs(arrayB[iPivot]) < 1.0e-4);
     1989        for (i = 0; i < numberRows; i++) {
     1990          if (fabs(arrayB[i] - 1.0) < 1.0e-4)
     1991            break;
     1992        }
     1993        assert(i < numberRows);
     1994        printf("variable on row %d landed up on row %d\n", iPivot, i);
     1995        arrayB[i] = 0.0;
     1996      }
     1997      for (i = 0; i < numberRows; i++)
     1998        assert(fabs(arrayB[i]) < 1.0e-4);
     1999      region2B.clear();
     2000    }
     2001#endif
     2002#ifdef CLP_FACTORIZATION_INSTRUMENT
     2003    factorization_instrument(2);
     2004#endif
     2005    if (anyChanged && model->algorithm() < 0 && solveType > 0) {
     2006      double dummyCost;
     2007      static_cast< ClpSimplexDual * >(model)->changeBounds(3,
     2008        NULL, dummyCost);
     2009    }
     2010    return coinFactorizationB_->status();
     2011  }
     2012  // If too many compressions increase area
     2013  if (coinFactorizationA_->pivots() > 1 && coinFactorizationA_->numberCompressions() * 10 > coinFactorizationA_->pivots() + 10) {
     2014    coinFactorizationA_->areaFactor(coinFactorizationA_->areaFactor() * 1.1);
     2015  }
     2016  //int numberPivots=coinFactorizationA_->pivots();
    20532017#if 0
    20542018     if (model->algorithm() > 0)
     
    20562020#endif
    20572021#ifndef SLIM_CLP
    2058      if (!networkBasis_ || doCheck) {
    2059 #endif
    2060           coinFactorizationA_->setStatus(-99);
    2061           int * pivotVariable = model->pivotVariable();
    2062           int nTimesRound = 0;
    2063           //returns 0 -okay, -1 singular, -2 too many in basis, -99 memory */
    2064           while (coinFactorizationA_->status() < -98) {
    2065                nTimesRound++;
     2022  if (!networkBasis_ || doCheck) {
     2023#endif
     2024    coinFactorizationA_->setStatus(-99);
     2025    int *pivotVariable = model->pivotVariable();
     2026    int nTimesRound = 0;
     2027    //returns 0 -okay, -1 singular, -2 too many in basis, -99 memory */
     2028    while (coinFactorizationA_->status() < -98) {
     2029      nTimesRound++;
    20662030
    2067                int i;
    2068                int numberBasic = 0;
    2069                int numberRowBasic;
    2070                // Move pivot variables across if they look good
    2071                int * pivotTemp = model->rowArray(0)->getIndices();
    2072                assert (!model->rowArray(0)->getNumElements());
    2073                if (!matrix->rhsOffset(model)) {
     2031      int i;
     2032      int numberBasic = 0;
     2033      int numberRowBasic;
     2034      // Move pivot variables across if they look good
     2035      int *pivotTemp = model->rowArray(0)->getIndices();
     2036      assert(!model->rowArray(0)->getNumElements());
     2037      if (!matrix->rhsOffset(model)) {
    20742038#if 0
    20752039                    if (numberSave > 0) {
     
    21402104                    }
    21412105#endif
    2142                     // Seems to prefer things in order so quickest
    2143                     // way is to go though like this
    2144                     for (i = 0; i < numberRows; i++) {
    2145                          if (model->getRowStatus(i) == ClpSimplex::basic)
    2146                               pivotTemp[numberBasic++] = i;
    2147                     }
    2148                     numberRowBasic = numberBasic;
    2149                     /* Put column basic variables into pivotVariable
     2106        // Seems to prefer things in order so quickest
     2107        // way is to go though like this
     2108        for (i = 0; i < numberRows; i++) {
     2109          if (model->getRowStatus(i) == ClpSimplex::basic)
     2110            pivotTemp[numberBasic++] = i;
     2111        }
     2112        numberRowBasic = numberBasic;
     2113        /* Put column basic variables into pivotVariable
    21502114                       This is done by ClpMatrixBase to allow override for gub
    21512115                    */
    2152                     matrix->generalExpanded(model, 0, numberBasic);
    2153                } else {
    2154                     // Long matrix - do a different way
    2155                     bool fullSearch = false;
    2156                     for (i = 0; i < numberRows; i++) {
    2157                          int iPivot = pivotVariable[i];
    2158                          if (iPivot >= numberColumns) {
    2159                               pivotTemp[numberBasic++] = iPivot - numberColumns;
    2160                          }
    2161                     }
    2162                     numberRowBasic = numberBasic;
    2163                     for (i = 0; i < numberRows; i++) {
    2164                          int iPivot = pivotVariable[i];
    2165                          if (iPivot < numberColumns) {
    2166                               if (iPivot >= 0) {
    2167                                    pivotTemp[numberBasic++] = iPivot;
    2168                               } else {
    2169                                    // not full basis
    2170                                    fullSearch = true;
    2171                                    break;
    2172                               }
    2173                          }
    2174                     }
    2175                     if (fullSearch) {
    2176                          // do slow way
    2177                          numberBasic = 0;
    2178                          for (i = 0; i < numberRows; i++) {
    2179                               if (model->getRowStatus(i) == ClpSimplex::basic)
    2180                                    pivotTemp[numberBasic++] = i;
    2181                          }
    2182                          numberRowBasic = numberBasic;
    2183                          /* Put column basic variables into pivotVariable
     2116        matrix->generalExpanded(model, 0, numberBasic);
     2117      } else {
     2118        // Long matrix - do a different way
     2119        bool fullSearch = false;
     2120        for (i = 0; i < numberRows; i++) {
     2121          int iPivot = pivotVariable[i];
     2122          if (iPivot >= numberColumns) {
     2123            pivotTemp[numberBasic++] = iPivot - numberColumns;
     2124          }
     2125        }
     2126        numberRowBasic = numberBasic;
     2127        for (i = 0; i < numberRows; i++) {
     2128          int iPivot = pivotVariable[i];
     2129          if (iPivot < numberColumns) {
     2130            if (iPivot >= 0) {
     2131              pivotTemp[numberBasic++] = iPivot;
     2132            } else {
     2133              // not full basis
     2134              fullSearch = true;
     2135              break;
     2136            }
     2137          }
     2138        }
     2139        if (fullSearch) {
     2140          // do slow way
     2141          numberBasic = 0;
     2142          for (i = 0; i < numberRows; i++) {
     2143            if (model->getRowStatus(i) == ClpSimplex::basic)
     2144              pivotTemp[numberBasic++] = i;
     2145          }
     2146          numberRowBasic = numberBasic;
     2147          /* Put column basic variables into pivotVariable
    21842148                            This is done by ClpMatrixBase to allow override for gub
    21852149                         */
    2186                          matrix->generalExpanded(model, 0, numberBasic);
    2187                     }
    2188                }
    2189                if (numberBasic > model->maximumBasic()) {
     2150          matrix->generalExpanded(model, 0, numberBasic);
     2151        }
     2152      }
     2153      if (numberBasic > model->maximumBasic()) {
    21902154#if 0 // ndef NDEBUG
    21912155                    printf("%d basic - should only be %d\n",
    21922156                           numberBasic, numberRows);
    21932157#endif
    2194                     // Take out some
    2195                     numberBasic = numberRowBasic;
    2196                     for (int i = 0; i < numberColumns; i++) {
    2197                          if (model->getColumnStatus(i) == ClpSimplex::basic) {
    2198                               if (numberBasic < numberRows)
    2199                                    numberBasic++;
    2200                               else
    2201                                    model->setColumnStatus(i, ClpSimplex::superBasic);
    2202                          }
    2203                     }
    2204                     numberBasic = numberRowBasic;
    2205                     matrix->generalExpanded(model, 0, numberBasic);
    2206                }
    2207 #ifndef SLIM_CLP
    2208                // see if matrix a network
     2158        // Take out some
     2159        numberBasic = numberRowBasic;
     2160        for (int i = 0; i < numberColumns; i++) {
     2161          if (model->getColumnStatus(i) == ClpSimplex::basic) {
     2162            if (numberBasic < numberRows)
     2163              numberBasic++;
     2164            else
     2165              model->setColumnStatus(i, ClpSimplex::superBasic);
     2166          }
     2167        }
     2168        numberBasic = numberRowBasic;
     2169        matrix->generalExpanded(model, 0, numberBasic);
     2170      }
     2171#ifndef SLIM_CLP
     2172      // see if matrix a network
    22092173#ifndef NO_RTTI
    2210                ClpNetworkMatrix* networkMatrix =
    2211                     dynamic_cast< ClpNetworkMatrix*>(model->clpMatrix());
     2174      ClpNetworkMatrix *networkMatrix = dynamic_cast< ClpNetworkMatrix * >(model->clpMatrix());
    22122175#else
    2213 ClpNetworkMatrix* networkMatrix = NULL;
    2214 if (model->clpMatrix()->type() == 11)
    2215      networkMatrix =
    2216           static_cast< ClpNetworkMatrix*>(model->clpMatrix());
    2217 #endif
    2218                // If network - still allow ordinary factorization first time for laziness
    2219                if (networkMatrix)
    2220                     coinFactorizationA_->setBiasLU(0); // All to U if network
    2221                //int saveMaximumPivots = maximumPivots();
    2222                delete networkBasis_;
    2223                networkBasis_ = NULL;
    2224                if (networkMatrix && !doCheck)
    2225                     maximumPivots(1);
    2226 #endif
    2227                //printf("L, U, R %d %d %d\n",numberElementsL(),numberElementsU(),numberElementsR());
    2228                while (coinFactorizationA_->status() == -99) {
    2229                     // maybe for speed will be better to leave as many regions as possible
    2230                     coinFactorizationA_->gutsOfDestructor();
    2231                     coinFactorizationA_->gutsOfInitialize(2);
    2232                     int numberElements = numberRowBasic;
     2176      ClpNetworkMatrix *networkMatrix = NULL;
     2177      if (model->clpMatrix()->type() == 11)
     2178        networkMatrix = static_cast< ClpNetworkMatrix * >(model->clpMatrix());
     2179#endif
     2180      // If network - still allow ordinary factorization first time for laziness
     2181      if (networkMatrix)
     2182        coinFactorizationA_->setBiasLU(0); // All to U if network
     2183      //int saveMaximumPivots = maximumPivots();
     2184      delete networkBasis_;
     2185      networkBasis_ = NULL;
     2186      if (networkMatrix && !doCheck)
     2187        maximumPivots(1);
     2188#endif
     2189      //printf("L, U, R %d %d %d\n",numberElementsL(),numberElementsU(),numberElementsR());
     2190      while (coinFactorizationA_->status() == -99) {
     2191        // maybe for speed will be better to leave as many regions as possible
     2192        coinFactorizationA_->gutsOfDestructor();
     2193        coinFactorizationA_->gutsOfInitialize(2);
     2194        int numberElements = numberRowBasic;
    22332195
    2234                     // compute how much in basis
     2196        // compute how much in basis
    22352197
    2236                     int i;
    2237                     // can change for gub
    2238                     int numberColumnBasic = numberBasic - numberRowBasic;
     2198        int i;
     2199        // can change for gub
     2200        int numberColumnBasic = numberBasic - numberRowBasic;
    22392201
    2240                     numberElements += matrix->countBasis( pivotTemp + numberRowBasic,
    2241                                                           numberColumnBasic);
    2242                     // and recompute as network side say different
    2243                     if (model->numberIterations())
    2244                          numberRowBasic = numberBasic - numberColumnBasic;
    2245                     numberElements = 3 * numberBasic + 3 * numberElements + 20000;
    2246                     coinFactorizationA_->getAreas ( numberRows,
    2247                                                     numberRowBasic + numberColumnBasic, numberElements,
    2248                                                     2 * numberElements );
    2249                     //fill
    2250                     // Fill in counts so we can skip part of preProcess
    2251                     int * numberInRow = coinFactorizationA_->numberInRow();
    2252                     int * numberInColumn = coinFactorizationA_->numberInColumn();
    2253                     CoinZeroN ( numberInRow, coinFactorizationA_->numberRows() + 1 );
    2254                     CoinZeroN ( numberInColumn, coinFactorizationA_->maximumColumnsExtra() + 1 );
    2255                     CoinFactorizationDouble * elementU = coinFactorizationA_->elementU();
    2256                     int * indexRowU = coinFactorizationA_->indexRowU();
    2257                     int * startColumnU = coinFactorizationA_->startColumnU();
     2202        numberElements += matrix->countBasis(pivotTemp + numberRowBasic,
     2203          numberColumnBasic);
     2204        // and recompute as network side say different
     2205        if (model->numberIterations())
     2206          numberRowBasic = numberBasic - numberColumnBasic;
     2207        numberElements = 3 * numberBasic + 3 * numberElements + 20000;
     2208        coinFactorizationA_->getAreas(numberRows,
     2209          numberRowBasic + numberColumnBasic, numberElements,
     2210          2 * numberElements);
     2211        //fill
     2212        // Fill in counts so we can skip part of preProcess
     2213        int *numberInRow = coinFactorizationA_->numberInRow();
     2214        int *numberInColumn = coinFactorizationA_->numberInColumn();
     2215        CoinZeroN(numberInRow, coinFactorizationA_->numberRows() + 1);
     2216        CoinZeroN(numberInColumn, coinFactorizationA_->maximumColumnsExtra() + 1);
     2217        CoinFactorizationDouble *elementU = coinFactorizationA_->elementU();
     2218        int *indexRowU = coinFactorizationA_->indexRowU();
     2219        int *startColumnU = coinFactorizationA_->startColumnU();
    22582220#ifndef COIN_FAST_CODE
    2259                     double slackValue = coinFactorizationA_->slackValue();
    2260 #endif
    2261                     for (i = 0; i < numberRowBasic; i++) {
    2262                          int iRow = pivotTemp[i];
    2263                          indexRowU[i] = iRow;
    2264                          startColumnU[i] = i;
    2265                          elementU[i] = slackValue;
    2266                          numberInRow[iRow] = 1;
    2267                          numberInColumn[i] = 1;
    2268                     }
    2269                     startColumnU[numberRowBasic] = numberRowBasic;
    2270                     // can change for gub so redo
    2271                     numberColumnBasic = numberBasic - numberRowBasic;
    2272                     matrix->fillBasis(model,
    2273                                       pivotTemp + numberRowBasic,
    2274                                       numberColumnBasic,
    2275                                       indexRowU,
    2276                                       startColumnU + numberRowBasic,
    2277                                       numberInRow,
    2278                                       numberInColumn + numberRowBasic,
    2279                                       elementU);
     2221        double slackValue = coinFactorizationA_->slackValue();
     2222#endif
     2223        for (i = 0; i < numberRowBasic; i++) {
     2224          int iRow = pivotTemp[i];
     2225          indexRowU[i] = iRow;
     2226          startColumnU[i] = i;
     2227          elementU[i] = slackValue;
     2228          numberInRow[iRow] = 1;
     2229          numberInColumn[i] = 1;
     2230        }
     2231        startColumnU[numberRowBasic] = numberRowBasic;
     2232        // can change for gub so redo
     2233        numberColumnBasic = numberBasic - numberRowBasic;
     2234        matrix->fillBasis(model,
     2235          pivotTemp + numberRowBasic,
     2236          numberColumnBasic,
     2237          indexRowU,
     2238          startColumnU + numberRowBasic,
     2239          numberInRow,
     2240          numberInColumn + numberRowBasic,
     2241          elementU);
    22802242#if 0
    22812243                    {
     
    22862248                    }
    22872249#endif
    2288                     // recompute number basic
    2289                     numberBasic = numberRowBasic + numberColumnBasic;
    2290                     if (numberBasic)
    2291                          numberElements = startColumnU[numberBasic-1]
    2292                                           + numberInColumn[numberBasic-1];
    2293                     else
    2294                          numberElements = 0;
     2250        // recompute number basic
     2251        numberBasic = numberRowBasic + numberColumnBasic;
     2252        if (numberBasic)
     2253          numberElements = startColumnU[numberBasic - 1]
     2254            + numberInColumn[numberBasic - 1];
     2255        else
     2256          numberElements = 0;
    22952257#ifdef CLP_FACTORIZATION_NEW_TIMING
    2296                     lastNumberPivots_=0;
    2297                     effectiveStartNumberU_=numberElements-numberRows;
    2298                     //printf("%d slacks,%d in U at beginning\n",
    2299                     //numberRowBasic,numberElements);
    2300 #endif
    2301                     coinFactorizationA_->setNumberElementsU(numberElements);
    2302                     //saveFactorization("dump.d");
    2303                     if (coinFactorizationA_->biasLU() >= 3 || coinFactorizationA_->numberRows() != coinFactorizationA_->numberColumns())
    2304                          coinFactorizationA_->preProcess ( 2 );
    2305                     else
    2306                          coinFactorizationA_->preProcess ( 3 ); // no row copy
    2307                     coinFactorizationA_->factor (  );
     2258        lastNumberPivots_ = 0;
     2259        effectiveStartNumberU_ = numberElements - numberRows;
     2260        //printf("%d slacks,%d in U at beginning\n",
     2261        //numberRowBasic,numberElements);
     2262#endif
     2263        coinFactorizationA_->setNumberElementsU(numberElements);
     2264        //saveFactorization("dump.d");
     2265        if (coinFactorizationA_->biasLU() >= 3 || coinFactorizationA_->numberRows() != coinFactorizationA_->numberColumns())
     2266          coinFactorizationA_->preProcess(2);
     2267        else
     2268          coinFactorizationA_->preProcess(3); // no row copy
     2269        coinFactorizationA_->factor();
    23082270#ifdef CLP_FACTORIZATION_NEW_TIMING
    2309                     endLengthU_ = coinFactorizationA_->numberElements() -
    2310                       coinFactorizationA_->numberDense()*coinFactorizationA_->numberDense()
    2311                       -coinFactorizationA_->numberElementsL();
    2312 #endif
    2313                     if (coinFactorizationA_->status() == -99) {
    2314                          // get more memory
    2315                          coinFactorizationA_->areaFactor(2.0 * coinFactorizationA_->areaFactor());
    2316                     } else if (coinFactorizationA_->status() == -1 &&
    2317                                (model->numberIterations() == 0 || nTimesRound > 2) &&
    2318                                coinFactorizationA_->denseThreshold()) {
    2319                          // Round again without dense
    2320                          coinFactorizationA_->setDenseThreshold(0);
    2321                          coinFactorizationA_->setStatus(-99);
    2322                     }
    2323                }
    2324                // If we get here status is 0 or -1
    2325                if (coinFactorizationA_->status() == 0) {
    2326                     // We may need to tamper with order and redo - e.g. network with side
    2327                     int useNumberRows = numberRows;
    2328                     // **** we will also need to add test in dual steepest to do
    2329                     // as we do for network
    2330                     matrix->generalExpanded(model, 12, useNumberRows);
    2331                     const int * permuteBack = coinFactorizationA_->permuteBack();
    2332                     const int * back = coinFactorizationA_->pivotColumnBack();
    2333                     //int * pivotTemp = pivotColumn_.array();
    2334                     //ClpDisjointCopyN ( pivotVariable, numberRows , pivotTemp  );
     2271        endLengthU_ = coinFactorizationA_->numberElements() - coinFactorizationA_->numberDense() * coinFactorizationA_->numberDense()
     2272          - coinFactorizationA_->numberElementsL();
     2273#endif
     2274        if (coinFactorizationA_->status() == -99) {
     2275          // get more memory
     2276          coinFactorizationA_->areaFactor(2.0 * coinFactorizationA_->areaFactor());
     2277        } else if (coinFactorizationA_->status() == -1 && (model->numberIterations() == 0 || nTimesRound > 2) && coinFactorizationA_->denseThreshold()) {
     2278          // Round again without dense
     2279          coinFactorizationA_->setDenseThreshold(0);
     2280          coinFactorizationA_->setStatus(-99);
     2281        }
     2282      }
     2283      // If we get here status is 0 or -1
     2284      if (coinFactorizationA_->status() == 0) {
     2285        // We may need to tamper with order and redo - e.g. network with side
     2286        int useNumberRows = numberRows;
     2287        // **** we will also need to add test in dual steepest to do
     2288        // as we do for network
     2289        matrix->generalExpanded(model, 12, useNumberRows);
     2290        const int *permuteBack = coinFactorizationA_->permuteBack();
     2291        const int *back = coinFactorizationA_->pivotColumnBack();
     2292        //int * pivotTemp = pivotColumn_.array();
     2293        //ClpDisjointCopyN ( pivotVariable, numberRows , pivotTemp  );
    23352294#ifndef NDEBUG
    2336                     CoinFillN(pivotVariable, numberRows, -1);
    2337 #endif
    2338                     // Redo pivot order
    2339                     for (i = 0; i < numberRowBasic; i++) {
    2340                          int k = pivotTemp[i];
    2341                          // so rowIsBasic[k] would be permuteBack[back[i]]
    2342                          int j = permuteBack[back[i]];
    2343                          assert (pivotVariable[j] == -1);
    2344                          pivotVariable[j] = k + numberColumns;
    2345                     }
    2346                     for (; i < useNumberRows; i++) {
    2347                          int k = pivotTemp[i];
    2348                          // so rowIsBasic[k] would be permuteBack[back[i]]
    2349                          int j = permuteBack[back[i]];
    2350                          assert (pivotVariable[j] == -1);
    2351                          pivotVariable[j] = k;
    2352                     }
     2295        CoinFillN(pivotVariable, numberRows, -1);
     2296#endif
     2297        // Redo pivot order
     2298        for (i = 0; i < numberRowBasic; i++) {
     2299          int k = pivotTemp[i];
     2300          // so rowIsBasic[k] would be permuteBack[back[i]]
     2301          int j = permuteBack[back[i]];
     2302          assert(pivotVariable[j] == -1);
     2303          pivotVariable[j] = k + numberColumns;
     2304        }
     2305        for (; i < useNumberRows; i++) {
     2306          int k = pivotTemp[i];
     2307          // so rowIsBasic[k] would be permuteBack[back[i]]
     2308          int j = permuteBack[back[i]];
     2309          assert(pivotVariable[j] == -1);
     2310          pivotVariable[j] = k;
     2311        }
    23532312#if 0
    23542313                    if (numberSave >= 0) {
     
    23612320                    }
    23622321#endif
    2363                     // Set up permutation vector
    2364                     // these arrays start off as copies of permute
    2365                     // (and we could use permute_ instead of pivotColumn (not back though))
    2366                     ClpDisjointCopyN ( coinFactorizationA_->permute(), useNumberRows , coinFactorizationA_->pivotColumn()  );
    2367                     ClpDisjointCopyN ( coinFactorizationA_->permuteBack(), useNumberRows , coinFactorizationA_->pivotColumnBack()  );
    2368 #ifndef SLIM_CLP
    2369                     if (networkMatrix) {
    2370                          maximumPivots(CoinMax(2000, maximumPivots()));
    2371                          // redo arrays
    2372                          for (int iRow = 0; iRow < 4; iRow++) {
    2373                               int length = model->numberRows() + maximumPivots();
    2374                               if (iRow == 3 || model->objectiveAsObject()->type() > 1)
    2375                                    length += model->numberColumns();
    2376                               model->rowArray(iRow)->reserve(length);
    2377                          }
    2378                          // create network factorization
    2379                          if (doCheck)
    2380                               delete networkBasis_; // temp
    2381                          networkBasis_ = new ClpNetworkBasis(model, coinFactorizationA_->numberRows(),
    2382                                                              coinFactorizationA_->pivotRegion(),
    2383                                                              coinFactorizationA_->permuteBack(),
    2384                                                              coinFactorizationA_->startColumnU(),
    2385                                                              coinFactorizationA_->numberInColumn(),
    2386                                                              coinFactorizationA_->indexRowU(),
    2387                                                              coinFactorizationA_->elementU());
    2388                          // kill off arrays in ordinary factorization
    2389                          if (!doCheck) {
    2390                               coinFactorizationA_->gutsOfDestructor();
    2391                               // but make sure coinFactorizationA_->numberRows() set
    2392                               coinFactorizationA_->setNumberRows(model->numberRows());
    2393                               coinFactorizationA_->setStatus(0);
     2322        // Set up permutation vector
     2323        // these arrays start off as copies of permute
     2324        // (and we could use permute_ instead of pivotColumn (not back though))
     2325        ClpDisjointCopyN(coinFactorizationA_->permute(), useNumberRows, coinFactorizationA_->pivotColumn());
     2326        ClpDisjointCopyN(coinFactorizationA_->permuteBack(), useNumberRows, coinFactorizationA_->pivotColumnBack());
     2327#ifndef SLIM_CLP
     2328        if (networkMatrix) {
     2329          maximumPivots(CoinMax(2000, maximumPivots()));
     2330          // redo arrays
     2331          for (int iRow = 0; iRow < 4; iRow++) {
     2332            int length = model->numberRows() + maximumPivots();
     2333            if (iRow == 3 || model->objectiveAsObject()->type() > 1)
     2334              length += model->numberColumns();
     2335            model->rowArray(iRow)->reserve(length);
     2336          }
     2337          // create network factorization
     2338          if (doCheck)
     2339            delete networkBasis_; // temp
     2340          networkBasis_ = new ClpNetworkBasis(model, coinFactorizationA_->numberRows(),
     2341            coinFactorizationA_->pivotRegion(),
     2342            coinFactorizationA_->permuteBack(),
     2343            coinFactorizationA_->startColumnU(),
     2344            coinFactorizationA_->numberInColumn(),
     2345            coinFactorizationA_->indexRowU(),
     2346            coinFactorizationA_->elementU());
     2347          // kill off arrays in ordinary factorization
     2348          if (!doCheck) {
     2349            coinFactorizationA_->gutsOfDestructor();
     2350            // but make sure coinFactorizationA_->numberRows() set
     2351            coinFactorizationA_->setNumberRows(model->numberRows());
     2352            coinFactorizationA_->setStatus(0);
    23942353#if 0
    23952354                              // but put back permute arrays so odd things will work
     
    24032362                              }
    24042363#endif
    2405                          }
    2406                     } else {
    2407 #endif
    2408                          // See if worth going sparse and when
    2409                          coinFactorizationA_->checkSparse();
    2410 #ifndef SLIM_CLP
    2411                     }
    2412 #endif
    2413                } else if (coinFactorizationA_->status() == -1 && (solveType == 0 || solveType == 2)) {
    2414                     // This needs redoing as it was merged coding - does not need array
     2364          }
     2365        } else {
     2366#endif
     2367          // See if worth going sparse and when
     2368          coinFactorizationA_->checkSparse();
     2369#ifndef SLIM_CLP
     2370        }
     2371#endif
     2372      } else if (coinFactorizationA_->status() == -1 && (solveType == 0 || solveType == 2)) {
     2373        // This needs redoing as it was merged coding - does not need array
    24152374#if 1
    2416                     int numberTotal = numberRows + numberColumns;
    2417                     int * isBasic = new int [numberTotal];
    2418                     int * rowIsBasic = isBasic + numberColumns;
    2419                     int * columnIsBasic = isBasic;
    2420                     for (i = 0; i < numberTotal; i++)
    2421                          isBasic[i] = -1;
    2422                     for (i = 0; i < numberRowBasic; i++) {
    2423                          int iRow = pivotTemp[i];
    2424                          rowIsBasic[iRow] = 1;
    2425                     }
    2426                     for (; i < numberBasic; i++) {
    2427                          int iColumn = pivotTemp[i];
    2428                          columnIsBasic[iColumn] = 1;
    2429                     }
    2430                     numberBasic = 0;
    2431                     for (i = 0; i < numberRows; i++)
    2432                          pivotVariable[i] = -1;
    2433                     // mark as basic or non basic
    2434                     const int * pivotColumn = coinFactorizationA_->pivotColumn();
    2435                     for (i = 0; i < numberRows; i++) {
    2436                          if (rowIsBasic[i] >= 0) {
    2437                               if (pivotColumn[numberBasic] >= 0) {
    2438                                    rowIsBasic[i] = pivotColumn[numberBasic];
    2439                               } else {
    2440                                    rowIsBasic[i] = -1;
    2441                                    model->setRowStatus(i, ClpSimplex::superBasic);
    2442                               }
    2443                               numberBasic++;
    2444                          }
    2445                     }
    2446                     for (i = 0; i < numberColumns; i++) {
    2447                          if (columnIsBasic[i] >= 0) {
    2448                               if (pivotColumn[numberBasic] >= 0)
    2449                                    columnIsBasic[i] = pivotColumn[numberBasic];
    2450                               else
    2451                                    columnIsBasic[i] = -1;
    2452                               numberBasic++;
    2453                          }
    2454                     }
    2455                     // leave pivotVariable in useful form for cleaning basis
    2456                     int * pivotVariable = model->pivotVariable();
    2457                     for (i = 0; i < numberRows; i++) {
    2458                          pivotVariable[i] = -1;
    2459                     }
     2375        int numberTotal = numberRows + numberColumns;
     2376        int *isBasic = new int[numberTotal];
     2377        int *rowIsBasic = isBasic + numberColumns;
     2378        int *columnIsBasic = isBasic;
     2379        for (i = 0; i < numberTotal; i++)
     2380          isBasic[i] = -1;
     2381        for (i = 0; i < numberRowBasic; i++) {
     2382          int iRow = pivotTemp[i];
     2383          rowIsBasic[iRow] = 1;
     2384        }
     2385        for (; i < numberBasic; i++) {
     2386          int iColumn = pivotTemp[i];
     2387          columnIsBasic[iColumn] = 1;
     2388        }
     2389        numberBasic = 0;
     2390        for (i = 0; i < numberRows; i++)
     2391          pivotVariable[i] = -1;
     2392        // mark as basic or non basic
     2393        const int *pivotColumn = coinFactorizationA_->pivotColumn();
     2394        for (i = 0; i < numberRows; i++) {
     2395          if (rowIsBasic[i] >= 0) {
     2396            if (pivotColumn[numberBasic] >= 0) {
     2397              rowIsBasic[i] = pivotColumn[numberBasic];
     2398            } else {
     2399              rowIsBasic[i] = -1;
     2400              model->setRowStatus(i, ClpSimplex::superBasic);
     2401            }
     2402            numberBasic++;
     2403          }
     2404        }
     2405        for (i = 0; i < numberColumns; i++) {
     2406          if (columnIsBasic[i] >= 0) {
     2407            if (pivotColumn[numberBasic] >= 0)
     2408              columnIsBasic[i] = pivotColumn[numberBasic];
     2409            else
     2410              columnIsBasic[i] = -1;
     2411            numberBasic++;
     2412          }
     2413        }
     2414        // leave pivotVariable in useful form for cleaning basis
     2415        int *pivotVariable = model->pivotVariable();
     2416        for (i = 0; i < numberRows; i++) {
     2417          pivotVariable[i] = -1;
     2418        }
    24602419
    2461                     for (i = 0; i < numberRows; i++) {
    2462                          if (model->getRowStatus(i) == ClpSimplex::basic) {
    2463                               int iPivot = rowIsBasic[i];
    2464                               if (iPivot >= 0)
    2465                                    pivotVariable[iPivot] = i + numberColumns;
    2466                          }
    2467                     }
    2468                     for (i = 0; i < numberColumns; i++) {
    2469                          if (model->getColumnStatus(i) == ClpSimplex::basic) {
    2470                               int iPivot = columnIsBasic[i];
    2471                               if (iPivot >= 0)
    2472                                    pivotVariable[iPivot] = i;
    2473                          }
    2474                     }
    2475                     delete [] isBasic;
     2420        for (i = 0; i < numberRows; i++) {
     2421          if (model->getRowStatus(i) == ClpSimplex::basic) {
     2422            int iPivot = rowIsBasic[i];
     2423            if (iPivot >= 0)
     2424              pivotVariable[iPivot] = i + numberColumns;
     2425          }
     2426        }
     2427        for (i = 0; i < numberColumns; i++) {
     2428          if (model->getColumnStatus(i) == ClpSimplex::basic) {
     2429            int iPivot = columnIsBasic[i];
     2430            if (iPivot >= 0)
     2431              pivotVariable[iPivot] = i;
     2432          }
     2433        }
     2434        delete[] isBasic;
    24762435#else
    2477                     {
    2478                       //int * lastColumn = lastColumn_.array(); // -1 or pivot row
    2479                       int * lastRow = coinFactorizationA_->lastRow(); // -1 or pivot sequence (inside sequence)
    2480                       for ( int i=0;i<numberRows;i++) {
    2481                         int iSeq = lastRow[i];
    2482                         if (iSeq >=0) {
    2483                           pivotVariable[i]=pivotTemp[iSeq];
    2484                           model->setRowStatus(i, ClpSimplex::superBasic);
    2485                         } else {
    2486                           pivotVariable[i]=i+numberColumns;
    2487                         }
    2488                       }
    2489                     }
    2490 #endif
    2491                     double * columnLower = model->lowerRegion();
    2492                     double * columnUpper = model->upperRegion();
    2493                     double * columnActivity = model->solutionRegion();
    2494                     double * rowLower = model->lowerRegion(0);
    2495                     double * rowUpper = model->upperRegion(0);
    2496                     double * rowActivity = model->solutionRegion(0);
    2497                     //redo basis - first take ALL columns out
    2498                     int iColumn;
    2499                     double largeValue = model->largeValue();
    2500                     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    2501                          if (model->getColumnStatus(iColumn) == ClpSimplex::basic) {
    2502                               // take out
    2503                               if (!valuesPass) {
    2504                                    double lower = columnLower[iColumn];
    2505                                    double upper = columnUpper[iColumn];
    2506                                    double value = columnActivity[iColumn];
    2507                                    if (lower > -largeValue || upper < largeValue) {
    2508                                         if (fabs(value - lower) < fabs(value - upper)) {
    2509                                              model->setColumnStatus(iColumn, ClpSimplex::atLowerBound);
    2510                                              columnActivity[iColumn] = lower;
    2511                                         } else {
    2512                                              model->setColumnStatus(iColumn, ClpSimplex::atUpperBound);
    2513                                              columnActivity[iColumn] = upper;
    2514                                         }
    2515                                    } else {
    2516                                         model->setColumnStatus(iColumn, ClpSimplex::isFree);
    2517                                    }
    2518                               } else {
    2519                                    model->setColumnStatus(iColumn, ClpSimplex::superBasic);
    2520                               }
    2521                          }
    2522                     }
    2523                     int iRow;
    2524                     for (iRow = 0; iRow < numberRows; iRow++) {
    2525                          int iSequence = pivotVariable[iRow];
    2526                          if (iSequence >= 0) {
    2527                               // basic
    2528                               if (iSequence >= numberColumns) {
    2529                                    // slack in - leave
    2530                                    //assert (iSequence-numberColumns==iRow);
    2531                               } else {
    2532                                    assert(model->getRowStatus(iRow) != ClpSimplex::basic);
    2533                                    // put back structural
    2534                                    model->setColumnStatus(iSequence, ClpSimplex::basic);
    2535                               }
    2536                          } else {
    2537                               // put in slack
    2538                               model->setRowStatus(iRow, ClpSimplex::basic);
    2539                          }
    2540                     }
    2541                     // Put back any key variables for gub
    2542                     int dummy;
    2543                     matrix->generalExpanded(model, 1, dummy);
    2544                     // signal repeat
    2545                     coinFactorizationA_->setStatus(-99);
    2546                     // set fixed if they are
    2547                     for (iRow = 0; iRow < numberRows; iRow++) {
    2548                          if (model->getRowStatus(iRow) != ClpSimplex::basic ) {
    2549                               if (rowLower[iRow] == rowUpper[iRow]) {
    2550                                    rowActivity[iRow] = rowLower[iRow];
    2551                                    model->setRowStatus(iRow, ClpSimplex::isFixed);
    2552                               }
    2553                          }
    2554                     }
    2555                     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    2556                          if (model->getColumnStatus(iColumn) != ClpSimplex::basic ) {
    2557                               if (columnLower[iColumn] == columnUpper[iColumn]) {
    2558                                    columnActivity[iColumn] = columnLower[iColumn];
    2559                                    model->setColumnStatus(iColumn, ClpSimplex::isFixed);
    2560                               }
    2561                          }
    2562                     }
    2563                }
    2564           }
    2565 #ifndef SLIM_CLP
    2566      } else {
    2567           // network - fake factorization - do nothing
    2568           coinFactorizationA_->setStatus(0);
    2569           coinFactorizationA_->setPivots(0);
    2570      }
    2571 #endif
    2572 #ifndef SLIM_CLP
    2573      if (!coinFactorizationA_->status()) {
    2574           // take out part if quadratic
    2575           if (model->algorithm() == 2) {
    2576                ClpObjective * obj = model->objectiveAsObject();
     2436        {
     2437          //int * lastColumn = lastColumn_.array(); // -1 or pivot row
     2438          int *lastRow = coinFactorizationA_->lastRow(); // -1 or pivot sequence (inside sequence)
     2439          for (int i = 0; i < numberRows; i++) {
     2440            int iSeq = lastRow[i];
     2441            if (iSeq >= 0) {
     2442              pivotVariable[i] = pivotTemp[iSeq];
     2443              model->setRowStatus(i, ClpSimplex::superBasic);
     2444            } else {
     2445              pivotVariable[i] = i + numberColumns;
     2446            }
     2447          }
     2448        }
     2449#endif
     2450        double *columnLower = model->lowerRegion();
     2451        double *columnUpper = model->upperRegion();
     2452        double *columnActivity = model->solutionRegion();
     2453        double *rowLower = model->lowerRegion(0);
     2454        double *rowUpper = model->upperRegion(0);
     2455        double *rowActivity = model->solutionRegion(0);
     2456        //redo basis - first take ALL columns out
     2457        int iColumn;
     2458        double largeValue = model->largeValue();
     2459        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     2460          if (model->getColumnStatus(iColumn) == ClpSimplex::basic) {
     2461            // take out
     2462            if (!valuesPass) {
     2463              double lower = columnLower[iColumn];
     2464              double upper = columnUpper[iColumn];
     2465              double value = columnActivity[iColumn];
     2466              if (lower > -largeValue || upper < largeValue) {
     2467                if (fabs(value - lower) < fabs(value - upper)) {
     2468                  model->setColumnStatus(iColumn, ClpSimplex::atLowerBound);
     2469                  columnActivity[iColumn] = lower;
     2470                } else {
     2471                  model->setColumnStatus(iColumn, ClpSimplex::atUpperBound);
     2472                  columnActivity[iColumn] = upper;
     2473                }
     2474              } else {
     2475                model->setColumnStatus(iColumn, ClpSimplex::isFree);
     2476              }
     2477            } else {
     2478              model->setColumnStatus(iColumn, ClpSimplex::superBasic);
     2479            }
     2480          }
     2481        }
     2482        int iRow;
     2483        for (iRow = 0; iRow < numberRows; iRow++) {
     2484          int iSequence = pivotVariable[iRow];
     2485          if (iSequence >= 0) {
     2486            // basic
     2487            if (iSequence >= numberColumns) {
     2488              // slack in - leave
     2489              //assert (iSequence-numberColumns==iRow);
     2490            } else {
     2491              assert(model->getRowStatus(iRow) != ClpSimplex::basic);
     2492              // put back structural
     2493              model->setColumnStatus(iSequence, ClpSimplex::basic);
     2494            }
     2495          } else {
     2496            // put in slack
     2497            model->setRowStatus(iRow, ClpSimplex::basic);
     2498          }
     2499        }
     2500        // Put back any key variables for gub
     2501        int dummy;
     2502        matrix->generalExpanded(model, 1, dummy);
     2503        // signal repeat
     2504        coinFactorizationA_->setStatus(-99);
     2505        // set fixed if they are
     2506        for (iRow = 0; iRow < numberRows; iRow++) {
     2507          if (model->getRowStatus(iRow) != ClpSimplex::basic) {
     2508            if (rowLower[iRow] == rowUpper[iRow]) {
     2509              rowActivity[iRow] = rowLower[iRow];
     2510              model->setRowStatus(iRow, ClpSimplex::isFixed);
     2511            }
     2512          }
     2513        }
     2514        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     2515          if (model->getColumnStatus(iColumn) != ClpSimplex::basic) {
     2516            if (columnLower[iColumn] == columnUpper[iColumn]) {
     2517              columnActivity[iColumn] = columnLower[iColumn];
     2518              model->setColumnStatus(iColumn, ClpSimplex::isFixed);
     2519            }
     2520          }
     2521        }
     2522      }
     2523    }
     2524#ifndef SLIM_CLP
     2525  } else {
     2526    // network - fake factorization - do nothing
     2527    coinFactorizationA_->setStatus(0);
     2528    coinFactorizationA_->setPivots(0);
     2529  }
     2530#endif
     2531#ifndef SLIM_CLP
     2532  if (!coinFactorizationA_->status()) {
     2533    // take out part if quadratic
     2534    if (model->algorithm() == 2) {
     2535      ClpObjective *obj = model->objectiveAsObject();
    25772536#ifndef NDEBUG
    2578                ClpQuadraticObjective * quadraticObj = (dynamic_cast< ClpQuadraticObjective*>(obj));
    2579                assert (quadraticObj);
     2537      ClpQuadraticObjective *quadraticObj = (dynamic_cast< ClpQuadraticObjective * >(obj));
     2538      assert(quadraticObj);
    25802539#else
    2581 ClpQuadraticObjective * quadraticObj = (static_cast< ClpQuadraticObjective*>(obj));
    2582 #endif
    2583                CoinPackedMatrix * quadratic = quadraticObj->quadraticObjective();
    2584                int numberXColumns = quadratic->getNumCols();
    2585                assert (numberXColumns < numberColumns);
    2586                int base = numberColumns - numberXColumns;
    2587                int * which = new int [numberXColumns];
    2588                int * pivotVariable = model->pivotVariable();
    2589                int * permute = pivotColumn();
    2590                int i;
    2591                int n = 0;
    2592                for (i = 0; i < numberRows; i++) {
    2593                     int iSj = pivotVariable[i] - base;
    2594                     if (iSj >= 0 && iSj < numberXColumns)
    2595                          which[n++] = permute[i];
    2596                }
    2597                if (n)
    2598                     coinFactorizationA_->emptyRows(n, which);
    2599                delete [] which;
    2600           }
    2601      }
    2602 #endif
    2603 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2604      factorization_instrument(2);
    2605 #endif
    2606      return coinFactorizationA_->status();
     2540      ClpQuadraticObjective *quadraticObj = (static_cast< ClpQuadraticObjective * >(obj));
     2541#endif
     2542      CoinPackedMatrix *quadratic = quadraticObj->quadraticObjective();
     2543      int numberXColumns = quadratic->getNumCols();
     2544      assert(numberXColumns < numberColumns);
     2545      int base = numberColumns - numberXColumns;
     2546      int *which = new int[numberXColumns];
     2547      int *pivotVariable = model->pivotVariable();
     2548      int *permute = pivotColumn();
     2549      int i;
     2550      int n = 0;
     2551      for (i = 0; i < numberRows; i++) {
     2552        int iSj = pivotVariable[i] - base;
     2553        if (iSj >= 0 && iSj < numberXColumns)
     2554          which[n++] = permute[i];
     2555      }
     2556      if (n)
     2557        coinFactorizationA_->emptyRows(n, which);
     2558      delete[] which;
     2559    }
     2560  }
     2561#endif
     2562#ifdef CLP_FACTORIZATION_INSTRUMENT
     2563  factorization_instrument(2);
     2564#endif
     2565  return coinFactorizationA_->status();
    26072566}
    26082567/* Replaces one Column in basis,
     
    26132572   after factorization and thereafter re-factorize
    26142573   partial update already in U */
    2615 int
    2616 ClpFactorization::replaceColumn ( const ClpSimplex * model,
    2617                                   CoinIndexedVector * regionSparse,
    2618                                   CoinIndexedVector * tableauColumn,
    2619                                   int pivotRow,
    2620                                   double pivotCheck ,
    2621                                   bool checkBeforeModifying,
    2622                                   double acceptablePivot)
    2623 {
    2624 #ifndef SLIM_CLP
    2625      if (!networkBasis_) {
     2574int ClpFactorization::replaceColumn(const ClpSimplex *model,
     2575  CoinIndexedVector *regionSparse,
     2576  CoinIndexedVector *tableauColumn,
     2577  int pivotRow,
     2578  double pivotCheck,
     2579  bool checkBeforeModifying,
     2580  double acceptablePivot)
     2581{
     2582#ifndef SLIM_CLP
     2583  if (!networkBasis_) {
    26262584#endif
    26272585#ifdef CLP_FACTORIZATION_NEW_TIMING
    26282586#ifdef CLP_FACTORIZATION_INSTRUMENT
    2629           factorization_instrument(-1);
    2630 #endif
    2631           int nOld=0;
    2632           int nNew=0;
    2633           int seq;
    2634 #if COIN_BIG_INDEX==0 && CLP_POOL_MATRIX==0       
    2635           const CoinPackedMatrix * matrix=model->matrix();
    2636           const int * columnLength = matrix->getVectorLengths();
     2587    factorization_instrument(-1);
     2588#endif
     2589    int nOld = 0;
     2590    int nNew = 0;
     2591    int seq;
     2592#if COIN_BIG_INDEX == 0 && CLP_POOL_MATRIX == 0
     2593    const CoinPackedMatrix *matrix = model->matrix();
     2594    const int *columnLength = matrix->getVectorLengths();
    26372595#else
    2638           const int * columnLength = model->clpMatrix()->getVectorLengths();
    2639 #endif
    2640           seq=model->sequenceIn();
    2641           if (seq>=0&&seq<model->numberColumns()+model->numberRows()) {
    2642             if (seq<model->numberColumns()) {
    2643               nNew=columnLength[seq];
    2644             } else {
    2645               nNew=1;
    2646             }
    2647           }
    2648           seq=model->sequenceOut();
    2649           if (seq>=0&&seq<model->numberColumns()+model->numberRows()) {
    2650             if (seq<model->numberColumns()) {
    2651               nOld=columnLength[seq];
    2652             } else {
    2653               nOld=1;
    2654             }
    2655           }
    2656           effectiveStartNumberU_ += nNew-nOld;
    2657 #endif
    2658           int returnCode;
    2659           // see if FT
    2660           if (!coinFactorizationA_ || coinFactorizationA_->forrestTomlin()) {
    2661                if (coinFactorizationA_) {
    2662 #if ABC_USE_COIN_FACTORIZATION<2
    2663                  returnCode =
    2664                    coinFactorizationA_->replaceColumn(regionSparse,
    2665                                                       pivotRow,
    2666                                                       pivotCheck,
    2667                                                       checkBeforeModifying,
    2668                                                       acceptablePivot);
     2596    const int *columnLength = model->clpMatrix()->getVectorLengths();
     2597#endif
     2598    seq = model->sequenceIn();
     2599    if (seq >= 0 && seq < model->numberColumns() + model->numberRows()) {
     2600      if (seq < model->numberColumns()) {
     2601        nNew = columnLength[seq];
     2602      } else {
     2603        nNew = 1;
     2604      }
     2605    }
     2606    seq = model->sequenceOut();
     2607    if (seq >= 0 && seq < model->numberColumns() + model->numberRows()) {
     2608      if (seq < model->numberColumns()) {
     2609        nOld = columnLength[seq];
     2610      } else {
     2611        nOld = 1;
     2612      }
     2613    }
     2614    effectiveStartNumberU_ += nNew - nOld;
     2615#endif
     2616    int returnCode;
     2617    // see if FT
     2618    if (!coinFactorizationA_ || coinFactorizationA_->forrestTomlin()) {
     2619      if (coinFactorizationA_) {
     2620#if ABC_USE_COIN_FACTORIZATION < 2
     2621        returnCode = coinFactorizationA_->replaceColumn(regionSparse,
     2622          pivotRow,
     2623          pivotCheck,
     2624          checkBeforeModifying,
     2625          acceptablePivot);
    26692626#else
    2670                  // fake btran alpha until I understand
    2671                  double btranAlpha=model->alpha();
    2672                  double ftAlpha =
    2673                    coinFactorizationA_->checkReplacePart1(regionSparse,
    2674                                                           pivotRow);
    2675                  returnCode =
    2676                    coinFactorizationA_->checkReplacePart2(pivotRow,
    2677                                                           btranAlpha,
    2678                                                           model->alpha(),
    2679                                                           ftAlpha,
    2680                                                           acceptablePivot);
    2681                  if (returnCode<2)
    2682                    coinFactorizationA_->replaceColumnPart3(regionSparse,
    2683                                                            pivotRow,
    2684                                                            model->alpha());
    2685 #endif
    2686                } else {
    2687                     bool tab = coinFactorizationB_->wantsTableauColumn();
     2627        // fake btran alpha until I understand
     2628        double btranAlpha = model->alpha();
     2629        double ftAlpha = coinFactorizationA_->checkReplacePart1(regionSparse,
     2630          pivotRow);
     2631        returnCode = coinFactorizationA_->checkReplacePart2(pivotRow,
     2632          btranAlpha,
     2633          model->alpha(),
     2634          ftAlpha,
     2635          acceptablePivot);
     2636        if (returnCode < 2)
     2637          coinFactorizationA_->replaceColumnPart3(regionSparse,
     2638            pivotRow,
     2639            model->alpha());
     2640#endif
     2641      } else {
     2642        bool tab = coinFactorizationB_->wantsTableauColumn();
    26882643#ifdef CLP_REUSE_ETAS
    2689                     int tempInfo[2];
    2690                     tempInfo[1] = model_->sequenceOut();
     2644        int tempInfo[2];
     2645        tempInfo[1] = model_->sequenceOut();
    26912646#else
    2692                     int tempInfo[1];
    2693 #endif
    2694                     tempInfo[0] = model->numberIterations();
    2695                     coinFactorizationB_->setUsefulInformation(tempInfo, 1);
    2696                     returnCode =
    2697                          coinFactorizationB_->replaceColumn(tab ? tableauColumn : regionSparse,
    2698                                                             pivotRow,
    2699                                                             pivotCheck,
    2700                                                             checkBeforeModifying,
    2701                                                             acceptablePivot);
     2647        int tempInfo[1];
     2648#endif
     2649        tempInfo[0] = model->numberIterations();
     2650        coinFactorizationB_->setUsefulInformation(tempInfo, 1);
     2651        returnCode = coinFactorizationB_->replaceColumn(tab ? tableauColumn : regionSparse,
     2652          pivotRow,
     2653          pivotCheck,
     2654          checkBeforeModifying,
     2655          acceptablePivot);
    27022656#ifdef CLP_DEBUG
    2703                     // check basic
    2704                     int numberRows = coinFactorizationB_->numberRows();
    2705                     CoinIndexedVector region1(2 * numberRows);
    2706                     CoinIndexedVector region2A(2 * numberRows);
    2707                     CoinIndexedVector region2B(2 * numberRows);
    2708                     int iPivot;
    2709                     double * arrayB = region2B.denseVector();
    2710                     int * pivotVariable = model->pivotVariable();
    2711                     int i;
    2712                     for (iPivot = 0; iPivot < numberRows; iPivot++) {
    2713                          int iSequence = pivotVariable[iPivot];
    2714                          if (iPivot == pivotRow)
    2715                               iSequence = model->sequenceIn();
    2716                          model->unpack(&region2B, iSequence);
    2717                          coinFactorizationB_->updateColumn(&region1, &region2B);
    2718                          assert (fabs(arrayB[iPivot] - 1.0) < 1.0e-4);
    2719                          arrayB[iPivot] = 0.0;
    2720                          for (i = 0; i < numberRows; i++)
    2721                               assert (fabs(arrayB[i]) < 1.0e-4);
    2722                          region2B.clear();
    2723                     }
    2724 #endif
    2725                }
    2726           } else {
    2727                returnCode = coinFactorizationA_->replaceColumnPFI(tableauColumn,
    2728                             pivotRow, pivotCheck); // Note array
    2729           }
    2730 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2731           factorization_instrument(3);
    2732 #endif
    2733           return returnCode;
     2657        // check basic
     2658        int numberRows = coinFactorizationB_->numberRows();
     2659        CoinIndexedVector region1(2 * numberRows);
     2660        CoinIndexedVector region2A(2 * numberRows);
     2661        CoinIndexedVector region2B(2 * numberRows);
     2662        int iPivot;
     2663        double *arrayB = region2B.denseVector();
     2664        int *pivotVariable = model->pivotVariable();
     2665        int i;
     2666        for (iPivot = 0; iPivot < numberRows; iPivot++) {
     2667          int iSequence = pivotVariable[iPivot];
     2668          if (iPivot == pivotRow)
     2669            iSequence = model->sequenceIn();
     2670          model->unpack(&region2B, iSequence);
     2671          coinFactorizationB_->updateColumn(&region1, &region2B);
     2672          assert(fabs(arrayB[iPivot] - 1.0) < 1.0e-4);
     2673          arrayB[iPivot] = 0.0;
     2674          for (i = 0; i < numberRows; i++)
     2675            assert(fabs(arrayB[i]) < 1.0e-4);
     2676          region2B.clear();
     2677        }
     2678#endif
     2679      }
     2680    } else {
     2681      returnCode = coinFactorizationA_->replaceColumnPFI(tableauColumn,
     2682        pivotRow, pivotCheck); // Note array
     2683    }
     2684#ifdef CLP_FACTORIZATION_INSTRUMENT
     2685    factorization_instrument(3);
     2686#endif
     2687    return returnCode;
    27342688
    27352689#ifndef SLIM_CLP
    2736      } else {
    2737           if (doCheck) {
    2738                int returnCode = coinFactorizationA_->replaceColumn(regionSparse,
    2739                                 pivotRow,
    2740                                 pivotCheck,
    2741                                 checkBeforeModifying,
    2742                                 acceptablePivot);
    2743                networkBasis_->replaceColumn(regionSparse,
    2744                                             pivotRow);
    2745                return returnCode;
    2746           } else {
    2747                // increase number of pivots
    2748                coinFactorizationA_->setPivots(coinFactorizationA_->pivots() + 1);
    2749                return networkBasis_->replaceColumn(regionSparse,
    2750                                                    pivotRow);
    2751           }
    2752      }
     2690  } else {
     2691    if (doCheck) {
     2692      int returnCode = coinFactorizationA_->replaceColumn(regionSparse,
     2693        pivotRow,
     2694        pivotCheck,
     2695        checkBeforeModifying,
     2696        acceptablePivot);
     2697      networkBasis_->replaceColumn(regionSparse,
     2698        pivotRow);
     2699      return returnCode;
     2700    } else {
     2701      // increase number of pivots
     2702      coinFactorizationA_->setPivots(coinFactorizationA_->pivots() + 1);
     2703      return networkBasis_->replaceColumn(regionSparse,
     2704        pivotRow);
     2705    }
     2706  }
    27532707#endif
    27542708}
     
    27572711   number returned is negative if no room
    27582712   region1 starts as zero and is zero at end */
    2759 int
    2760 ClpFactorization::updateColumnFT ( CoinIndexedVector * regionSparse,
    2761                                    CoinIndexedVector * regionSparse2)
     2713int ClpFactorization::updateColumnFT(CoinIndexedVector *regionSparse,
     2714  CoinIndexedVector *regionSparse2)
    27622715{
    27632716#ifdef CLP_DEBUG
    2764      regionSparse->checkClear();
    2765 #endif
    2766      if (!numberRows())
    2767           return 0;
    2768 #ifndef SLIM_CLP
    2769      if (!networkBasis_) {
    2770 #endif
    2771 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2772           factorization_instrument(-1);
    2773 #endif
    2774           int returnCode;
    2775           if (coinFactorizationA_) {
    2776                coinFactorizationA_->setCollectStatistics(true);
    2777                returnCode = coinFactorizationA_->updateColumnFT(regionSparse,
    2778                             regionSparse2);
    2779                coinFactorizationA_->setCollectStatistics(false);
    2780           } else {
     2717  regionSparse->checkClear();
     2718#endif
     2719  if (!numberRows())
     2720    return 0;
     2721#ifndef SLIM_CLP
     2722  if (!networkBasis_) {
     2723#endif
     2724#ifdef CLP_FACTORIZATION_INSTRUMENT
     2725    factorization_instrument(-1);
     2726#endif
     2727    int returnCode;
     2728    if (coinFactorizationA_) {
     2729      coinFactorizationA_->setCollectStatistics(true);
     2730      returnCode = coinFactorizationA_->updateColumnFT(regionSparse,
     2731        regionSparse2);
     2732      coinFactorizationA_->setCollectStatistics(false);
     2733    } else {
    27812734#ifdef CLP_REUSE_ETAS
    2782               int tempInfo[2];
    2783               tempInfo[0] = model_->numberIterations();
    2784               tempInfo[1] = model_->sequenceIn();
    2785               coinFactorizationB_->setUsefulInformation(tempInfo, 2);
    2786 #endif
    2787               returnCode = coinFactorizationB_->updateColumnFT(regionSparse,
    2788                             regionSparse2);
    2789           }
    2790 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2791           factorization_instrument(4);
    2792 #endif
    2793           return returnCode;
    2794 #ifndef SLIM_CLP
    2795      } else {
     2735      int tempInfo[2];
     2736      tempInfo[0] = model_->numberIterations();
     2737      tempInfo[1] = model_->sequenceIn();
     2738      coinFactorizationB_->setUsefulInformation(tempInfo, 2);
     2739#endif
     2740      returnCode = coinFactorizationB_->updateColumnFT(regionSparse,
     2741        regionSparse2);
     2742    }
     2743#ifdef CLP_FACTORIZATION_INSTRUMENT
     2744    factorization_instrument(4);
     2745#endif
     2746    return returnCode;
     2747#ifndef SLIM_CLP
     2748  } else {
    27962749#ifdef CHECK_NETWORK
    2797           CoinIndexedVector * save = new CoinIndexedVector(*regionSparse2);
    2798           double * check = new double[coinFactorizationA_->numberRows()];
    2799           int returnCode = coinFactorizationA_->updateColumnFT(regionSparse,
    2800                            regionSparse2);
    2801           networkBasis_->updateColumn(regionSparse, save, -1);
    2802           int i;
    2803           double * array = regionSparse2->denseVector();
    2804           int * indices = regionSparse2->getIndices();
    2805           int n = regionSparse2->getNumElements();
    2806           memset(check, 0, coinFactorizationA_->numberRows()*sizeof(double));
    2807           double * array2 = save->denseVector();
    2808           int * indices2 = save->getIndices();
    2809           int n2 = save->getNumElements();
    2810           assert (n == n2);
    2811           if (save->packedMode()) {
    2812                for (i = 0; i < n; i++) {
    2813                     check[indices[i]] = array[i];
    2814                }
    2815                for (i = 0; i < n; i++) {
    2816                     double value2 = array2[i];
    2817                     assert (check[indices2[i]] == value2);
    2818                }
    2819           } else {
    2820                int numberRows = coinFactorizationA_->numberRows();
    2821                for (i = 0; i < numberRows; i++) {
    2822                     double value1 = array[i];
    2823                     double value2 = array2[i];
    2824                     assert (value1 == value2);
    2825                }
    2826           }
    2827           delete save;
    2828           delete [] check;
    2829           return returnCode;
     2750    CoinIndexedVector *save = new CoinIndexedVector(*regionSparse2);
     2751    double *check = new double[coinFactorizationA_->numberRows()];
     2752    int returnCode = coinFactorizationA_->updateColumnFT(regionSparse,
     2753      regionSparse2);
     2754    networkBasis_->updateColumn(regionSparse, save, -1);
     2755    int i;
     2756    double *array = regionSparse2->denseVector();
     2757    int *indices = regionSparse2->getIndices();
     2758    int n = regionSparse2->getNumElements();
     2759    memset(check, 0, coinFactorizationA_->numberRows() * sizeof(double));
     2760    double *array2 = save->denseVector();
     2761    int *indices2 = save->getIndices();
     2762    int n2 = save->getNumElements();
     2763    assert(n == n2);
     2764    if (save->packedMode()) {
     2765      for (i = 0; i < n; i++) {
     2766        check[indices[i]] = array[i];
     2767      }
     2768      for (i = 0; i < n; i++) {
     2769        double value2 = array2[i];
     2770        assert(check[indices2[i]] == value2);
     2771      }
     2772    } else {
     2773      int numberRows = coinFactorizationA_->numberRows();
     2774      for (i = 0; i < numberRows; i++) {
     2775        double value1 = array[i];
     2776        double value2 = array2[i];
     2777        assert(value1 == value2);
     2778      }
     2779    }
     2780    delete save;
     2781    delete[] check;
     2782    return returnCode;
    28302783#else
    2831 networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
    2832 return 1;
    2833 #endif
    2834      }
     2784    networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
     2785    return 1;
     2786#endif
     2787  }
    28352788#endif
    28362789}
     
    28382791   number returned is negative if no room
    28392792   region1 starts as zero and is zero at end */
    2840 int
    2841 ClpFactorization::updateColumn ( CoinIndexedVector * regionSparse,
    2842                                  CoinIndexedVector * regionSparse2,
    2843                                  bool noPermute) const
     2793int ClpFactorization::updateColumn(CoinIndexedVector *regionSparse,
     2794  CoinIndexedVector *regionSparse2,
     2795  bool noPermute) const
    28442796{
    28452797#ifdef CLP_DEBUG
    2846      if (!noPermute)
    2847           regionSparse->checkClear();
    2848 #endif
    2849      if (!numberRows())
    2850           return 0;
    2851 #ifndef SLIM_CLP
    2852      if (!networkBasis_) {
    2853 #endif
    2854 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2855           factorization_instrument(-1);
    2856 #endif
    2857           int returnCode;
    2858           if (coinFactorizationA_) {
    2859                coinFactorizationA_->setCollectStatistics(doStatistics_);
    2860                returnCode = coinFactorizationA_->updateColumn(regionSparse,
    2861                             regionSparse2,
    2862                             noPermute);
    2863                coinFactorizationA_->setCollectStatistics(false);
    2864           } else {
    2865                returnCode = coinFactorizationB_->updateColumn(regionSparse,
    2866                             regionSparse2,
    2867                             noPermute);
    2868           }
    2869 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2870           factorization_instrument(5);
    2871 #endif
    2872           //#define PRINT_VECTOR
     2798  if (!noPermute)
     2799    regionSparse->checkClear();
     2800#endif
     2801  if (!numberRows())
     2802    return 0;
     2803#ifndef SLIM_CLP
     2804  if (!networkBasis_) {
     2805#endif
     2806#ifdef CLP_FACTORIZATION_INSTRUMENT
     2807    factorization_instrument(-1);
     2808#endif
     2809    int returnCode;
     2810    if (coinFactorizationA_) {
     2811      coinFactorizationA_->setCollectStatistics(doStatistics_);
     2812      returnCode = coinFactorizationA_->updateColumn(regionSparse,
     2813        regionSparse2,
     2814        noPermute);
     2815      coinFactorizationA_->setCollectStatistics(false);
     2816    } else {
     2817      returnCode = coinFactorizationB_->updateColumn(regionSparse,
     2818        regionSparse2,
     2819        noPermute);
     2820    }
     2821#ifdef CLP_FACTORIZATION_INSTRUMENT
     2822    factorization_instrument(5);
     2823#endif
     2824    //#define PRINT_VECTOR
    28732825#ifdef PRINT_VECTOR
    2874           printf("Update\n");
    2875           regionSparse2->print();
    2876 #endif
    2877           return returnCode;
    2878 #ifndef SLIM_CLP
    2879      } else {
     2826    printf("Update\n");
     2827    regionSparse2->print();
     2828#endif
     2829    return returnCode;
     2830#ifndef SLIM_CLP
     2831  } else {
    28802832#ifdef CHECK_NETWORK
    2881           CoinIndexedVector * save = new CoinIndexedVector(*regionSparse2);
    2882           double * check = new double[coinFactorizationA_->numberRows()];
    2883           int returnCode = coinFactorizationA_->updateColumn(regionSparse,
    2884                            regionSparse2,
    2885                            noPermute);
    2886           networkBasis_->updateColumn(regionSparse, save, -1);
    2887           int i;
    2888           double * array = regionSparse2->denseVector();
    2889           int * indices = regionSparse2->getIndices();
    2890           int n = regionSparse2->getNumElements();
    2891           memset(check, 0, coinFactorizationA_->numberRows()*sizeof(double));
    2892           double * array2 = save->denseVector();
    2893           int * indices2 = save->getIndices();
    2894           int n2 = save->getNumElements();
    2895           assert (n == n2);
    2896           if (save->packedMode()) {
    2897                for (i = 0; i < n; i++) {
    2898                     check[indices[i]] = array[i];
    2899                }
    2900                for (i = 0; i < n; i++) {
    2901                     double value2 = array2[i];
    2902                     assert (check[indices2[i]] == value2);
    2903                }
    2904           } else {
    2905                int numberRows = coinFactorizationA_->numberRows();
    2906                for (i = 0; i < numberRows; i++) {
    2907                     double value1 = array[i];
    2908                     double value2 = array2[i];
    2909                     assert (value1 == value2);
    2910                }
    2911           }
    2912           delete save;
    2913           delete [] check;
    2914           return returnCode;
     2833    CoinIndexedVector *save = new CoinIndexedVector(*regionSparse2);
     2834    double *check = new double[coinFactorizationA_->numberRows()];
     2835    int returnCode = coinFactorizationA_->updateColumn(regionSparse,
     2836      regionSparse2,
     2837      noPermute);
     2838    networkBasis_->updateColumn(regionSparse, save, -1);
     2839    int i;
     2840    double *array = regionSparse2->denseVector();
     2841    int *indices = regionSparse2->getIndices();
     2842    int n = regionSparse2->getNumElements();
     2843    memset(check, 0, coinFactorizationA_->numberRows() * sizeof(double));
     2844    double *array2 = save->denseVector();
     2845    int *indices2 = save->getIndices();
     2846    int n2 = save->getNumElements();
     2847    assert(n == n2);
     2848    if (save->packedMode()) {
     2849      for (i = 0; i < n; i++) {
     2850        check[indices[i]] = array[i];
     2851      }
     2852      for (i = 0; i < n; i++) {
     2853        double value2 = array2[i];
     2854        assert(check[indices2[i]] == value2);
     2855      }
     2856    } else {
     2857      int numberRows = coinFactorizationA_->numberRows();
     2858      for (i = 0; i < numberRows; i++) {
     2859        double value1 = array[i];
     2860        double value2 = array2[i];
     2861        assert(value1 == value2);
     2862      }
     2863    }
     2864    delete save;
     2865    delete[] check;
     2866    return returnCode;
    29152867#else
    2916 networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
    2917 return 1;
    2918 #endif
    2919      }
     2868    networkBasis_->updateColumn(regionSparse, regionSparse2, -1);
     2869    return 1;
     2870#endif
     2871  }
    29202872#endif
    29212873}
     
    29252877   Also updates region3
    29262878   region1 starts as zero and is zero at end */
    2927 int
    2928 ClpFactorization::updateTwoColumnsFT ( CoinIndexedVector * regionSparse1,
    2929                                        CoinIndexedVector * regionSparse2,
    2930                                        CoinIndexedVector * regionSparse3,
    2931                                        bool noPermuteRegion3)
     2879int ClpFactorization::updateTwoColumnsFT(CoinIndexedVector *regionSparse1,
     2880  CoinIndexedVector *regionSparse2,
     2881  CoinIndexedVector *regionSparse3,
     2882  bool noPermuteRegion3)
    29322883{
    29332884#ifdef CLP_DEBUG
    2934      regionSparse1->checkClear();
    2935 #endif
    2936      if (!numberRows())
    2937           return 0;
    2938      int returnCode = 0;
    2939 #ifndef SLIM_CLP
    2940      if (!networkBasis_) {
    2941 #endif
    2942 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2943           factorization_instrument(-1);
    2944 #endif
    2945           if (coinFactorizationA_) {
    2946                coinFactorizationA_->setCollectStatistics(true);
    2947                if (coinFactorizationA_->spaceForForrestTomlin()) {
    2948                     assert (regionSparse2->packedMode());
    2949                     assert (!regionSparse3->packedMode());
    2950                     returnCode = coinFactorizationA_->updateTwoColumnsFT(regionSparse1,
    2951                                  regionSparse2,
    2952                                  regionSparse3,
    2953                                  noPermuteRegion3);
    2954                } else {
    2955                     returnCode = coinFactorizationA_->updateColumnFT(regionSparse1,
    2956                                  regionSparse2);
    2957                     coinFactorizationA_->updateColumn(regionSparse1,
    2958                                                       regionSparse3,
    2959                                                       noPermuteRegion3);
    2960                }
    2961                coinFactorizationA_->setCollectStatistics(false);
    2962           } else {
     2885  regionSparse1->checkClear();
     2886#endif
     2887  if (!numberRows())
     2888    return 0;
     2889  int returnCode = 0;
     2890#ifndef SLIM_CLP
     2891  if (!networkBasis_) {
     2892#endif
     2893#ifdef CLP_FACTORIZATION_INSTRUMENT
     2894    factorization_instrument(-1);
     2895#endif
     2896    if (coinFactorizationA_) {
     2897      coinFactorizationA_->setCollectStatistics(true);
     2898      if (coinFactorizationA_->spaceForForrestTomlin()) {
     2899        assert(regionSparse2->packedMode());
     2900        assert(!regionSparse3->packedMode());
     2901        returnCode = coinFactorizationA_->updateTwoColumnsFT(regionSparse1,
     2902          regionSparse2,
     2903          regionSparse3,
     2904          noPermuteRegion3);
     2905      } else {
     2906        returnCode = coinFactorizationA_->updateColumnFT(regionSparse1,
     2907          regionSparse2);
     2908        coinFactorizationA_->updateColumn(regionSparse1,
     2909          regionSparse3,
     2910          noPermuteRegion3);
     2911      }
     2912      coinFactorizationA_->setCollectStatistics(false);
     2913    } else {
    29632914#if 0
    29642915               CoinSimpFactorization * fact =
     
    29782929#else
    29792930#ifdef CLP_REUSE_ETAS
    2980                     int tempInfo[2];
    2981                     tempInfo[0] = model_->numberIterations();
    2982                     tempInfo[1] = model_->sequenceIn();
    2983                     coinFactorizationB_->setUsefulInformation(tempInfo, 3);
    2984 #endif
    2985                     returnCode =
    2986                       coinFactorizationB_->updateTwoColumnsFT(
    2987                                                               regionSparse1,
    2988                                                               regionSparse2,
    2989                                                               regionSparse3,
    2990                                                               noPermuteRegion3);
    2991 #endif
    2992           }
    2993 #ifdef CLP_FACTORIZATION_INSTRUMENT
    2994           factorization_instrument(9);
     2931      int tempInfo[2];
     2932      tempInfo[0] = model_->numberIterations();
     2933      tempInfo[1] = model_->sequenceIn();
     2934      coinFactorizationB_->setUsefulInformation(tempInfo, 3);
     2935#endif
     2936      returnCode = coinFactorizationB_->updateTwoColumnsFT(
     2937        regionSparse1,
     2938        regionSparse2,
     2939        regionSparse3,
     2940        noPermuteRegion3);
     2941#endif
     2942    }
     2943#ifdef CLP_FACTORIZATION_INSTRUMENT
     2944    factorization_instrument(9);
    29952945#endif
    29962946#ifdef PRINT_VECTOR
    2997           printf("UpdateTwoFT\n");
    2998           regionSparse2->print();
    2999           regionSparse3->print();
    3000 #endif
    3001           return returnCode;
    3002 #ifndef SLIM_CLP
    3003      } else {
    3004           returnCode = updateColumnFT(regionSparse1, regionSparse2);
    3005           updateColumn(regionSparse1, regionSparse3, noPermuteRegion3);
    3006      }
    3007 #endif
    3008      return returnCode;
     2947    printf("UpdateTwoFT\n");
     2948    regionSparse2->print();
     2949    regionSparse3->print();
     2950#endif
     2951    return returnCode;
     2952#ifndef SLIM_CLP
     2953  } else {
     2954    returnCode = updateColumnFT(regionSparse1, regionSparse2);
     2955    updateColumn(regionSparse1, regionSparse3, noPermuteRegion3);
     2956  }
     2957#endif
     2958  return returnCode;
    30092959}
    30102960/* Updates one column (FTRAN) from region2
    30112961   number returned is negative if no room
    30122962   region1 starts as zero and is zero at end */
    3013 int
    3014 ClpFactorization::updateColumnForDebug ( CoinIndexedVector * regionSparse,
    3015           CoinIndexedVector * regionSparse2,
    3016           bool noPermute) const
    3017 {
    3018      if (!noPermute)
    3019           regionSparse->checkClear();
    3020      if (!coinFactorizationA_->numberRows())
    3021           return 0;
    3022      coinFactorizationA_->setCollectStatistics(false);
    3023      // above doesn't work any more - do by hand
    3024      double save[15];
    3025      memcpy(save,&coinFactorizationA_->ftranCountInput_,sizeof(save));
    3026      int returnCode = coinFactorizationA_->updateColumn(regionSparse,
    3027                       regionSparse2,
    3028                       noPermute);
    3029      memcpy(&coinFactorizationA_->ftranCountInput_,save,sizeof(save));
    3030      return returnCode;
     2963int ClpFactorization::updateColumnForDebug(CoinIndexedVector *regionSparse,
     2964  CoinIndexedVector *regionSparse2,
     2965  bool noPermute) const
     2966{
     2967  if (!noPermute)
     2968    regionSparse->checkClear();
     2969  if (!coinFactorizationA_->numberRows())
     2970    return 0;
     2971  coinFactorizationA_->setCollectStatistics(false);
     2972  // above doesn't work any more - do by hand
     2973  double save[15];
     2974  memcpy(save, &coinFactorizationA_->ftranCountInput_, sizeof(save));
     2975  int returnCode = coinFactorizationA_->updateColumn(regionSparse,
     2976    regionSparse2,
     2977    noPermute);
     2978  memcpy(&coinFactorizationA_->ftranCountInput_, save, sizeof(save));
     2979  return returnCode;
    30312980}
    30322981/* Updates one column (BTRAN) from region2
    30332982   region1 starts as zero and is zero at end */
    3034 int
    3035 ClpFactorization::updateColumnTranspose ( CoinIndexedVector * regionSparse,
    3036           CoinIndexedVector * regionSparse2) const
    3037 {
    3038      if (!numberRows())
    3039           return 0;
    3040 #ifndef SLIM_CLP
    3041      if (!networkBasis_) {
    3042 #endif
    3043 #ifdef CLP_FACTORIZATION_INSTRUMENT
    3044           factorization_instrument(-1);
    3045 #endif
    3046           int returnCode;
     2983int ClpFactorization::updateColumnTranspose(CoinIndexedVector *regionSparse,
     2984  CoinIndexedVector *regionSparse2) const
     2985{
     2986  if (!numberRows())
     2987    return 0;
     2988#ifndef SLIM_CLP
     2989  if (!networkBasis_) {
     2990#endif
     2991#ifdef CLP_FACTORIZATION_INSTRUMENT
     2992    factorization_instrument(-1);
     2993#endif
     2994    int returnCode;
    30472995
    3048           if (coinFactorizationA_) {
    3049                coinFactorizationA_->setCollectStatistics(doStatistics_);
    3050                returnCode = coinFactorizationA_->updateColumnTranspose(regionSparse,
    3051                              regionSparse2);
    3052                coinFactorizationA_->setCollectStatistics(false);
    3053           } else {
    3054                returnCode = coinFactorizationB_->updateColumnTranspose(regionSparse,
    3055                             regionSparse2);
    3056           }
    3057 #ifdef CLP_FACTORIZATION_INSTRUMENT
    3058           factorization_instrument(6);
     2996    if (coinFactorizationA_) {
     2997      coinFactorizationA_->setCollectStatistics(doStatistics_);
     2998      returnCode = coinFactorizationA_->updateColumnTranspose(regionSparse,
     2999        regionSparse2);
     3000      coinFactorizationA_->setCollectStatistics(false);
     3001    } else {
     3002      returnCode = coinFactorizationB_->updateColumnTranspose(regionSparse,
     3003        regionSparse2);
     3004    }
     3005#ifdef CLP_FACTORIZATION_INSTRUMENT
     3006    factorization_instrument(6);
    30593007#endif
    30603008#ifdef PRINT_VECTOR
    3061           printf("UpdateTranspose\n");
    3062           regionSparse2->print();
    3063 #endif
    3064           return returnCode;
    3065 #ifndef SLIM_CLP
    3066      } else {
     3009    printf("UpdateTranspose\n");
     3010    regionSparse2->print();
     3011#endif
     3012    return returnCode;
     3013#ifndef SLIM_CLP
     3014  } else {
    30673015#ifdef CHECK_NETWORK
    3068           CoinIndexedVector * save = new CoinIndexedVector(*regionSparse2);
    3069           double * check = new double[coinFactorizationA_->numberRows()];
    3070           int returnCode = coinFactorizationA_->updateColumnTranspose(regionSparse,
    3071                            regionSparse2);
    3072           networkBasis_->updateColumnTranspose(regionSparse, save);
    3073           int i;
    3074           double * array = regionSparse2->denseVector();
    3075           int * indices = regionSparse2->getIndices();
    3076           int n = regionSparse2->getNumElements();
    3077           memset(check, 0, coinFactorizationA_->numberRows()*sizeof(double));
    3078           double * array2 = save->denseVector();
    3079           int * indices2 = save->getIndices();
    3080           int n2 = save->getNumElements();
    3081           assert (n == n2);
    3082           if (save->packedMode()) {
    3083                for (i = 0; i < n; i++) {
    3084                     check[indices[i]] = array[i];
    3085                }
    3086                for (i = 0; i < n; i++) {
    3087                     double value2 = array2[i];
    3088                     assert (check[indices2[i]] == value2);
    3089                }
    3090           } else {
    3091                int numberRows = coinFactorizationA_->numberRows();
    3092                for (i = 0; i < numberRows; i++) {
    3093                     double value1 = array[i];
    3094                     double value2 = array2[i];
    3095                     assert (value1 == value2);
    3096                }
    3097           }
    3098           delete save;
    3099           delete [] check;
    3100           return returnCode;
     3016    CoinIndexedVector *save = new CoinIndexedVector(*regionSparse2);
     3017    double *check = new double[coinFactorizationA_->numberRows()];
     3018    int returnCode = coinFactorizationA_->updateColumnTranspose(regionSparse,
     3019      regionSparse2);
     3020    networkBasis_->updateColumnTranspose(regionSparse, save);
     3021    int i;
     3022    double *array = regionSparse2->denseVector();
     3023    int *indices = regionSparse2->getIndices();
     3024    int n = regionSparse2->getNumElements();
     3025    memset(check, 0, coinFactorizationA_->numberRows() * sizeof(double));
     3026    double *array2 = save->denseVector();
     3027    int *indices2 = save->getIndices();
     3028    int n2 = save->getNumElements();
     3029    assert(n == n2);
     3030    if (save->packedMode()) {
     3031      for (i = 0; i < n; i++) {
     3032        check[indices[i]] = array[i];
     3033      }
     3034      for (i = 0; i < n; i++) {
     3035        double value2 = array2[i];
     3036        assert(check[indices2[i]] == value2);
     3037      }
     3038    } else {
     3039      int numberRows = coinFactorizationA_->numberRows();
     3040      for (i = 0; i < numberRows; i++) {
     3041        double value1 = array[i];
     3042        double value2 = array2[i];
     3043        assert(value1 == value2);
     3044      }
     3045    }
     3046    delete save;
     3047    delete[] check;
     3048    return returnCode;
    31013049#else
    3102 return networkBasis_->updateColumnTranspose(regionSparse, regionSparse2);
    3103 #endif
    3104      }
     3050    return networkBasis_->updateColumnTranspose(regionSparse, regionSparse2);
     3051#endif
     3052  }
    31053053#endif
    31063054}
     
    31093057   Note - if regionSparse2 packed on input - will be packed on output - same for 3
    31103058*/
    3111 void
    3112 ClpFactorization::updateTwoColumnsTranspose ( CoinIndexedVector * regionSparse,
    3113                                   CoinIndexedVector * regionSparse2,
    3114                                   CoinIndexedVector * regionSparse3) const
     3059void ClpFactorization::updateTwoColumnsTranspose(CoinIndexedVector *regionSparse,
     3060  CoinIndexedVector *regionSparse2,
     3061  CoinIndexedVector *regionSparse3) const
    31153062{
    31163063  if (!numberRows())
    3117           return;
    3118 #ifndef SLIM_CLP
    3119      if (!networkBasis_) {
    3120 #endif
    3121 #ifdef CLP_FACTORIZATION_INSTRUMENT
    3122           factorization_instrument(-1);
    3123 #endif
    3124           if (coinFactorizationA_) {
    3125                coinFactorizationA_->setCollectStatistics(doStatistics_);
     3064    return;
     3065#ifndef SLIM_CLP
     3066  if (!networkBasis_) {
     3067#endif
     3068#ifdef CLP_FACTORIZATION_INSTRUMENT
     3069    factorization_instrument(-1);
     3070#endif
     3071    if (coinFactorizationA_) {
     3072      coinFactorizationA_->setCollectStatistics(doStatistics_);
    31263073#if ABOCA_LITE_FACTORIZATION
    3127                coinFactorizationA_->updateTwoColumnsTranspose(regionSparse,
    3128                                                               regionSparse2,regionSparse3,abcState());
     3074      coinFactorizationA_->updateTwoColumnsTranspose(regionSparse,
     3075        regionSparse2, regionSparse3, abcState());
    31293076#else
    3130                coinFactorizationA_->updateTwoColumnsTranspose(regionSparse,
    3131               regionSparse2,regionSparse3,0);
    3132 #endif
    3133                coinFactorizationA_->setCollectStatistics(false);
    3134           } else {
    3135                coinFactorizationB_->updateColumnTranspose(regionSparse,
    3136                             regionSparse2);
    3137                coinFactorizationB_->updateColumnTranspose(regionSparse,
    3138                             regionSparse3);
    3139           }
    3140 #ifdef CLP_FACTORIZATION_INSTRUMENT
    3141           factorization_instrument(6);
    3142 #endif
    3143 #ifndef SLIM_CLP
    3144      } else {
    3145        updateColumnTranspose(regionSparse, regionSparse2);
    3146        updateColumnTranspose(regionSparse, regionSparse3);
    3147      }
     3077      coinFactorizationA_->updateTwoColumnsTranspose(regionSparse,
     3078        regionSparse2, regionSparse3, 0);
     3079#endif
     3080      coinFactorizationA_->setCollectStatistics(false);
     3081    } else {
     3082      coinFactorizationB_->updateColumnTranspose(regionSparse,
     3083        regionSparse2);
     3084      coinFactorizationB_->updateColumnTranspose(regionSparse,
     3085        regionSparse3);
     3086    }
     3087#ifdef CLP_FACTORIZATION_INSTRUMENT
     3088    factorization_instrument(6);
     3089#endif
     3090#ifndef SLIM_CLP
     3091  } else {
     3092    updateColumnTranspose(regionSparse, regionSparse2);
     3093    updateColumnTranspose(regionSparse, regionSparse3);
     3094  }
    31483095#endif
    31493096}
    31503097/* makes a row copy of L for speed and to allow very sparse problems */
    3151 void
    3152 ClpFactorization::goSparse()
    3153 {
    3154 #ifndef SLIM_CLP
    3155      if (!networkBasis_) {
    3156 #endif
    3157           if (coinFactorizationA_) {
    3158 #ifdef CLP_FACTORIZATION_INSTRUMENT
    3159                factorization_instrument(-1);
    3160 #endif
    3161                coinFactorizationA_->goSparse();
    3162 #ifdef CLP_FACTORIZATION_INSTRUMENT
    3163                factorization_instrument(7);
    3164 #endif
    3165           }
    3166      }
     3098void ClpFactorization::goSparse()
     3099{
     3100#ifndef SLIM_CLP
     3101  if (!networkBasis_) {
     3102#endif
     3103    if (coinFactorizationA_) {
     3104#ifdef CLP_FACTORIZATION_INSTRUMENT
     3105      factorization_instrument(-1);
     3106#endif
     3107      coinFactorizationA_->goSparse();
     3108#ifdef CLP_FACTORIZATION_INSTRUMENT
     3109      factorization_instrument(7);
     3110#endif
     3111    }
     3112  }
    31673113}
    31683114// Cleans up i.e. gets rid of network basis
    3169 void
    3170 ClpFactorization::cleanUp()
    3171 {
    3172 #ifndef SLIM_CLP
    3173      delete networkBasis_;
    3174      networkBasis_ = NULL;
    3175 #endif
    3176      if (coinFactorizationA_)
    3177           coinFactorizationA_->resetStatistics();
     3115void ClpFactorization::cleanUp()
     3116{
     3117#ifndef SLIM_CLP
     3118  delete networkBasis_;
     3119  networkBasis_ = NULL;
     3120#endif
     3121  if (coinFactorizationA_)
     3122    coinFactorizationA_->resetStatistics();
    31783123}
    31793124/// Says whether to redo pivot order
    3180 bool
    3181 ClpFactorization::needToReorder() const
     3125bool ClpFactorization::needToReorder() const
    31823126{
    31833127#ifdef CHECK_NETWORK
    3184      return true;
    3185 #endif
    3186 #ifndef SLIM_CLP
    3187      if (!networkBasis_)
    3188 #endif
    3189           return true;
    3190 #ifndef SLIM_CLP
    3191      else
    3192           return false;
     3128  return true;
     3129#endif
     3130#ifndef SLIM_CLP
     3131  if (!networkBasis_)
     3132#endif
     3133    return true;
     3134#ifndef SLIM_CLP
     3135  else
     3136    return false;
    31933137#endif
    31943138}
    31953139// Get weighted row list
    3196 void
    3197 ClpFactorization::getWeights(int * weights) const
    3198 {
    3199 #ifdef CLP_FACTORIZATION_INSTRUMENT
    3200      factorization_instrument(-1);
    3201 #endif
    3202 #ifndef SLIM_CLP
    3203      if (networkBasis_) {
    3204           // Network - just unit
    3205           int numberRows = coinFactorizationA_->numberRows();
    3206           for (int i = 0; i < numberRows; i++)
    3207                weights[i] = 1;
    3208           return;
    3209      }
    3210 #endif
    3211      int * numberInRow = coinFactorizationA_->numberInRow();
    3212      int * numberInColumn = coinFactorizationA_->numberInColumn();
    3213      int * permuteBack = coinFactorizationA_->pivotColumnBack();
    3214      int * indexRowU = coinFactorizationA_->indexRowU();
    3215      const int * startColumnU = coinFactorizationA_->startColumnU();
    3216      const int * startRowL = coinFactorizationA_->startRowL();
    3217      int numberRows = coinFactorizationA_->numberRows();
    3218      if (!startRowL || !coinFactorizationA_->numberInRow()) {
    3219           int * temp = new int[numberRows];
    3220           memset(temp, 0, numberRows * sizeof(int));
    3221           int i;
    3222           for (i = 0; i < numberRows; i++) {
    3223                // one for pivot
    3224                temp[i]++;
    3225                int j;
    3226                for (j = startColumnU[i]; j < startColumnU[i] + numberInColumn[i]; j++) {
    3227                     int iRow = indexRowU[j];
    3228                     temp[iRow]++;
    3229                }
    3230           }
    3231           int * startColumnL = coinFactorizationA_->startColumnL();
    3232           int * indexRowL = coinFactorizationA_->indexRowL();
    3233           int numberL = coinFactorizationA_->numberL();
    3234           int baseL = coinFactorizationA_->baseL();
    3235           for (i = baseL; i < baseL + numberL; i++) {
    3236                int j;
    3237                for (j = startColumnL[i]; j < startColumnL[i+1]; j++) {
    3238                     int iRow = indexRowL[j];