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

formatting

File:
1 edited

Legend:

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

    r2271 r2385  
    33// Corporation and others.  All Rights Reserved.
    44// This code is licensed under the terms of the Eclipse Public License (EPL).
    5 
    65
    76#include <cstdio>
     
    1817#include "ClpMessage.hpp"
    1918#ifdef CLP_PLUS_ONE_MATRIX
    20 static int oneitcount[13]={0,0,0,0,0,0,0,0,0,0,0,0,0};
     19static int oneitcount[13] = { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 };
    2120static void oneit(int i)
    2221{
    2322  if (!oneitcount[i])
    24     printf("Plus ones for call %d\n",i);
     23    printf("Plus ones for call %d\n", i);
    2524  oneitcount[i]++;
    2625  oneitcount[12]++;
    27   if ((oneitcount[12]%1000)==0) {
     26  if ((oneitcount[12] % 1000) == 0) {
    2827    printf("Plus counts");
    29     for (int j=0;j<12;j++)
    30       printf(" %d",oneitcount[j]);
     28    for (int j = 0; j < 12; j++)
     29      printf(" %d", oneitcount[j]);
    3130    printf("\n");
    32   }   
     31  }
    3332}
    3433#endif
     
    4039// Default Constructor
    4140//-------------------------------------------------------------------
    42 ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix ()
    43      : ClpMatrixBase()
    44 {
    45      setType(12);
    46      matrix_ = NULL;
    47      startPositive_ = NULL;
    48      startNegative_ = NULL;
    49      lengths_ = NULL;
    50      indices_ = NULL;
    51      numberRows_ = 0;
    52      numberColumns_ = 0;
     41ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix()
     42  : ClpMatrixBase()
     43{
     44  setType(12);
     45  matrix_ = NULL;
     46  startPositive_ = NULL;
     47  startNegative_ = NULL;
     48  lengths_ = NULL;
     49  indices_ = NULL;
     50  numberRows_ = 0;
     51  numberColumns_ = 0;
    5352#ifdef CLP_PLUS_ONE_MATRIX
    54      otherFlags_ = 0;
    55 #endif
    56      columnOrdered_ = true;
     53  otherFlags_ = 0;
     54#endif
     55  columnOrdered_ = true;
    5756}
    5857
     
    6059// Copy constructor
    6160//-------------------------------------------------------------------
    62 ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix (const ClpPlusMinusOneMatrix & rhs)
    63      : ClpMatrixBase(rhs)
    64 {
    65      matrix_ = NULL;
    66      startPositive_ = NULL;
    67      startNegative_ = NULL;
    68      lengths_ = NULL;
    69      indices_ = NULL;
    70      numberRows_ = rhs.numberRows_;
    71      numberColumns_ = rhs.numberColumns_;
     61ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix(const ClpPlusMinusOneMatrix &rhs)
     62  : ClpMatrixBase(rhs)
     63{
     64  matrix_ = NULL;
     65  startPositive_ = NULL;
     66  startNegative_ = NULL;
     67  lengths_ = NULL;
     68  indices_ = NULL;
     69  numberRows_ = rhs.numberRows_;
     70  numberColumns_ = rhs.numberColumns_;
    7271#ifdef CLP_PLUS_ONE_MATRIX
    73      otherFlags_ = rhs.otherFlags_;
    74 #endif
    75      columnOrdered_ = rhs.columnOrdered_;
    76      if (numberColumns_) {
    77           CoinBigIndex numberElements = rhs.startPositive_[numberColumns_];
    78           indices_ = new int [ numberElements];
    79           CoinMemcpyN(rhs.indices_, numberElements, indices_);
    80           startPositive_ = new CoinBigIndex [ numberColumns_+1];
    81           CoinMemcpyN(rhs.startPositive_, (numberColumns_ + 1), startPositive_);
    82           startNegative_ = new CoinBigIndex [ numberColumns_];
    83           CoinMemcpyN(rhs.startNegative_, numberColumns_, startNegative_);
    84      }
    85      int numberRows = getNumRows();
    86      if (rhs.rhsOffset_ && numberRows) {
    87           rhsOffset_ = ClpCopyOfArray(rhs.rhsOffset_, numberRows);
    88      } else {
    89           rhsOffset_ = NULL;
    90      }
     72  otherFlags_ = rhs.otherFlags_;
     73#endif
     74  columnOrdered_ = rhs.columnOrdered_;
     75  if (numberColumns_) {
     76    CoinBigIndex numberElements = rhs.startPositive_[numberColumns_];
     77    indices_ = new int[numberElements];
     78    CoinMemcpyN(rhs.indices_, numberElements, indices_);
     79    startPositive_ = new CoinBigIndex[numberColumns_ + 1];
     80    CoinMemcpyN(rhs.startPositive_, (numberColumns_ + 1), startPositive_);
     81    startNegative_ = new CoinBigIndex[numberColumns_];
     82    CoinMemcpyN(rhs.startNegative_, numberColumns_, startNegative_);
     83  }
     84  int numberRows = getNumRows();
     85  if (rhs.rhsOffset_ && numberRows) {
     86    rhsOffset_ = ClpCopyOfArray(rhs.rhsOffset_, numberRows);
     87  } else {
     88    rhsOffset_ = NULL;
     89  }
    9190}
    9291// Constructor from arrays
    9392ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix(int numberRows, int numberColumns,
    94           bool columnOrdered, const int * indices,
    95           const CoinBigIndex * startPositive,
    96           const CoinBigIndex * startNegative)
    97      : ClpMatrixBase()
    98 {
    99      setType(12);
    100      matrix_ = NULL;
    101      lengths_ = NULL;
    102      numberRows_ = numberRows;
    103      numberColumns_ = numberColumns;
    104      columnOrdered_ = columnOrdered;
     93  bool columnOrdered, const int *indices,
     94  const CoinBigIndex *startPositive,
     95  const CoinBigIndex *startNegative)
     96  : ClpMatrixBase()
     97{
     98  setType(12);
     99  matrix_ = NULL;
     100  lengths_ = NULL;
     101  numberRows_ = numberRows;
     102  numberColumns_ = numberColumns;
     103  columnOrdered_ = columnOrdered;
    105104#ifdef CLP_PLUS_ONE_MATRIX
    106      otherFlags_ = 0;
    107 #endif
    108      int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    109      CoinBigIndex numberElements = startPositive[numberMajor];
    110      startPositive_ = ClpCopyOfArray(startPositive, numberMajor + 1);
    111      startNegative_ = ClpCopyOfArray(startNegative, numberMajor);
    112      indices_ = ClpCopyOfArray(indices, numberElements);
    113      // Check valid
    114      checkValid(false);
     105  otherFlags_ = 0;
     106#endif
     107  int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     108  CoinBigIndex numberElements = startPositive[numberMajor];
     109  startPositive_ = ClpCopyOfArray(startPositive, numberMajor + 1);
     110  startNegative_ = ClpCopyOfArray(startNegative, numberMajor);
     111  indices_ = ClpCopyOfArray(indices, numberElements);
     112  // Check valid
     113  checkValid(false);
    115114}
    116115
    117 ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix (const CoinPackedMatrix & rhs)
    118      : ClpMatrixBase()
    119 {
    120      setType(12);
    121      matrix_ = NULL;
    122      startPositive_ = NULL;
    123      startNegative_ = NULL;
    124      lengths_ = NULL;
    125      indices_ = NULL;
    126      int iColumn;
    127      assert (rhs.isColOrdered());
    128      // get matrix data pointers
    129      const int * row = rhs.getIndices();
    130      const CoinBigIndex * columnStart = rhs.getVectorStarts();
    131      const int * columnLength = rhs.getVectorLengths();
    132      const double * elementByColumn = rhs.getElements();
    133      numberColumns_ = rhs.getNumCols();
    134      numberRows_ = -1;
    135      indices_ = new int[rhs.getNumElements()];
    136      startPositive_ = new CoinBigIndex [numberColumns_+1];
    137      startNegative_ = new CoinBigIndex [numberColumns_];
    138      int * temp = new int [rhs.getNumRows()];
    139      CoinBigIndex j = 0;
    140      CoinBigIndex numberGoodP = 0;
    141      CoinBigIndex numberGoodM = 0;
    142      CoinBigIndex numberBad = 0;
    143      for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    144           CoinBigIndex k;
    145           int iNeg = 0;
    146           startPositive_[iColumn] = j;
    147           for (k = columnStart[iColumn]; k < columnStart[iColumn] + columnLength[iColumn];
    148                     k++) {
    149                int iRow;
    150                if (fabs(elementByColumn[k] - 1.0) < 1.0e-10) {
    151                     iRow = row[k];
    152                     numberRows_ = CoinMax(numberRows_, iRow);
    153                     indices_[j++] = iRow;
    154                     numberGoodP++;
    155                } else if (fabs(elementByColumn[k] + 1.0) < 1.0e-10) {
    156                     iRow = row[k];
    157                     numberRows_ = CoinMax(numberRows_, iRow);
    158                     temp[iNeg++] = iRow;
    159                     numberGoodM++;
    160                } else {
    161                     numberBad++;
    162                }
    163           }
    164           // move negative
    165           startNegative_[iColumn] = j;
    166           for (k = 0; k < iNeg; k++) {
    167                indices_[j++] = temp[k];
    168           }
    169      }
    170      startPositive_[numberColumns_] = j;
    171      delete [] temp;
    172      if (numberBad) {
    173           delete [] indices_;
    174           indices_ = NULL;
    175           numberRows_ = 0;
    176           numberColumns_ = 0;
    177           delete [] startPositive_;
    178           delete [] startNegative_;
    179           // Put in statistics
    180           startPositive_ = new CoinBigIndex [3];
    181           startPositive_[0] = numberGoodP;
    182           startPositive_[1] = numberGoodM;
    183           startPositive_[2] = numberBad;
    184           startNegative_ = NULL;
    185      } else {
    186           numberRows_ ++; //  correct
    187           // but number should be same as rhs
    188           assert (numberRows_ <= rhs.getNumRows());
    189           numberRows_ = rhs.getNumRows();
    190           columnOrdered_ = true;
    191      }
    192      // Check valid
    193      if (!numberBad)
    194           checkValid(false);
     116ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix(const CoinPackedMatrix &rhs)
     117  : ClpMatrixBase()
     118{
     119  setType(12);
     120  matrix_ = NULL;
     121  startPositive_ = NULL;
     122  startNegative_ = NULL;
     123  lengths_ = NULL;
     124  indices_ = NULL;
     125  int iColumn;
     126  assert(rhs.isColOrdered());
     127  // get matrix data pointers
     128  const int *row = rhs.getIndices();
     129  const CoinBigIndex *columnStart = rhs.getVectorStarts();
     130  const int *columnLength = rhs.getVectorLengths();
     131  const double *elementByColumn = rhs.getElements();
     132  numberColumns_ = rhs.getNumCols();
     133  numberRows_ = -1;
     134  indices_ = new int[rhs.getNumElements()];
     135  startPositive_ = new CoinBigIndex[numberColumns_ + 1];
     136  startNegative_ = new CoinBigIndex[numberColumns_];
     137  int *temp = new int[rhs.getNumRows()];
     138  CoinBigIndex j = 0;
     139  CoinBigIndex numberGoodP = 0;
     140  CoinBigIndex numberGoodM = 0;
     141  CoinBigIndex numberBad = 0;
     142  for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     143    CoinBigIndex k;
     144    int iNeg = 0;
     145    startPositive_[iColumn] = j;
     146    for (k = columnStart[iColumn]; k < columnStart[iColumn] + columnLength[iColumn];
     147         k++) {
     148      int iRow;
     149      if (fabs(elementByColumn[k] - 1.0) < 1.0e-10) {
     150        iRow = row[k];
     151        numberRows_ = CoinMax(numberRows_, iRow);
     152        indices_[j++] = iRow;
     153        numberGoodP++;
     154      } else if (fabs(elementByColumn[k] + 1.0) < 1.0e-10) {
     155        iRow = row[k];
     156        numberRows_ = CoinMax(numberRows_, iRow);
     157        temp[iNeg++] = iRow;
     158        numberGoodM++;
     159      } else {
     160        numberBad++;
     161      }
     162    }
     163    // move negative
     164    startNegative_[iColumn] = j;
     165    for (k = 0; k < iNeg; k++) {
     166      indices_[j++] = temp[k];
     167    }
     168  }
     169  startPositive_[numberColumns_] = j;
     170  delete[] temp;
     171  if (numberBad) {
     172    delete[] indices_;
     173    indices_ = NULL;
     174    numberRows_ = 0;
     175    numberColumns_ = 0;
     176    delete[] startPositive_;
     177    delete[] startNegative_;
     178    // Put in statistics
     179    startPositive_ = new CoinBigIndex[3];
     180    startPositive_[0] = numberGoodP;
     181    startPositive_[1] = numberGoodM;
     182    startPositive_[2] = numberBad;
     183    startNegative_ = NULL;
     184  } else {
     185    numberRows_++; //  correct
     186    // but number should be same as rhs
     187    assert(numberRows_ <= rhs.getNumRows());
     188    numberRows_ = rhs.getNumRows();
     189    columnOrdered_ = true;
     190  }
     191  // Check valid
     192  if (!numberBad)
     193    checkValid(false);
    195194}
    196195
     
    198197// Destructor
    199198//-------------------------------------------------------------------
    200 ClpPlusMinusOneMatrix::~ClpPlusMinusOneMatrix ()
    201 {
    202      delete matrix_;
    203      delete [] startPositive_;
    204      delete [] startNegative_;
    205      delete [] lengths_;
    206      delete [] indices_;
     199ClpPlusMinusOneMatrix::~ClpPlusMinusOneMatrix()
     200{
     201  delete matrix_;
     202  delete[] startPositive_;
     203  delete[] startNegative_;
     204  delete[] lengths_;
     205  delete[] indices_;
    207206}
    208207
     
    211210//-------------------------------------------------------------------
    212211ClpPlusMinusOneMatrix &
    213 ClpPlusMinusOneMatrix::operator=(const ClpPlusMinusOneMatrix& rhs)
    214 {
    215      if (this != &rhs) {
    216           ClpMatrixBase::operator=(rhs);
    217           delete matrix_;
    218           delete [] startPositive_;
    219           delete [] startNegative_;
    220           delete [] lengths_;
    221           delete [] indices_;
    222           matrix_ = NULL;
    223           startPositive_ = NULL;
    224           lengths_ = NULL;
    225           indices_ = NULL;
    226           numberRows_ = rhs.numberRows_;
    227           numberColumns_ = rhs.numberColumns_;
    228           columnOrdered_ = rhs.columnOrdered_;
     212ClpPlusMinusOneMatrix::operator=(const ClpPlusMinusOneMatrix &rhs)
     213{
     214  if (this != &rhs) {
     215    ClpMatrixBase::operator=(rhs);
     216    delete matrix_;
     217    delete[] startPositive_;
     218    delete[] startNegative_;
     219    delete[] lengths_;
     220    delete[] indices_;
     221    matrix_ = NULL;
     222    startPositive_ = NULL;
     223    lengths_ = NULL;
     224    indices_ = NULL;
     225    numberRows_ = rhs.numberRows_;
     226    numberColumns_ = rhs.numberColumns_;
     227    columnOrdered_ = rhs.columnOrdered_;
    229228#ifdef CLP_PLUS_ONE_MATRIX
    230           otherFlags_ = rhs.otherFlags_;
    231 #endif
    232           if (numberColumns_) {
    233                CoinBigIndex numberElements = rhs.startPositive_[numberColumns_];
    234                indices_ = new int [ numberElements];
    235                CoinMemcpyN(rhs.indices_, numberElements, indices_);
    236                startPositive_ = new CoinBigIndex [ numberColumns_+1];
    237                CoinMemcpyN(rhs.startPositive_, (numberColumns_ + 1), startPositive_);
    238                startNegative_ = new CoinBigIndex [ numberColumns_];
    239                CoinMemcpyN(rhs.startNegative_, numberColumns_, startNegative_);
    240           }
    241      }
    242      return *this;
     229    otherFlags_ = rhs.otherFlags_;
     230#endif
     231    if (numberColumns_) {
     232      CoinBigIndex numberElements = rhs.startPositive_[numberColumns_];
     233      indices_ = new int[numberElements];
     234      CoinMemcpyN(rhs.indices_, numberElements, indices_);
     235      startPositive_ = new CoinBigIndex[numberColumns_ + 1];
     236      CoinMemcpyN(rhs.startPositive_, (numberColumns_ + 1), startPositive_);
     237      startNegative_ = new CoinBigIndex[numberColumns_];
     238      CoinMemcpyN(rhs.startNegative_, numberColumns_, startNegative_);
     239    }
     240  }
     241  return *this;
    243242}
    244243//-------------------------------------------------------------------
    245244// Clone
    246245//-------------------------------------------------------------------
    247 ClpMatrixBase * ClpPlusMinusOneMatrix::clone() const
    248 {
    249      return new ClpPlusMinusOneMatrix(*this);
     246ClpMatrixBase *ClpPlusMinusOneMatrix::clone() const
     247{
     248  return new ClpPlusMinusOneMatrix(*this);
    250249}
    251250/* Subset clone (without gaps).  Duplicates are allowed
    252251   and order is as given */
    253252ClpMatrixBase *
    254 ClpPlusMinusOneMatrix::subsetClone (int numberRows, const int * whichRows,
    255                                     int numberColumns,
    256                                     const int * whichColumns) const
    257 {
    258      return new ClpPlusMinusOneMatrix(*this, numberRows, whichRows,
    259                                       numberColumns, whichColumns);
     253ClpPlusMinusOneMatrix::subsetClone(int numberRows, const int *whichRows,
     254  int numberColumns,
     255  const int *whichColumns) const
     256{
     257  return new ClpPlusMinusOneMatrix(*this, numberRows, whichRows,
     258    numberColumns, whichColumns);
    260259}
    261260/* Subset constructor (without gaps).  Duplicates are allowed
    262261   and order is as given */
    263 ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix (
    264      const ClpPlusMinusOneMatrix & rhs,
    265      int numberRows, const int * whichRow,
    266      int numberColumns, const int * whichColumn)
    267      : ClpMatrixBase(rhs)
    268 {
    269      matrix_ = NULL;
    270      startPositive_ = NULL;
    271      startNegative_ = NULL;
    272      lengths_ = NULL;
    273      indices_ = NULL;
    274      numberRows_ = 0;
    275      numberColumns_ = 0;
    276      columnOrdered_ = rhs.columnOrdered_;
     262ClpPlusMinusOneMatrix::ClpPlusMinusOneMatrix(
     263  const ClpPlusMinusOneMatrix &rhs,
     264  int numberRows, const int *whichRow,
     265  int numberColumns, const int *whichColumn)
     266  : ClpMatrixBase(rhs)
     267{
     268  matrix_ = NULL;
     269  startPositive_ = NULL;
     270  startNegative_ = NULL;
     271  lengths_ = NULL;
     272  indices_ = NULL;
     273  numberRows_ = 0;
     274  numberColumns_ = 0;
     275  columnOrdered_ = rhs.columnOrdered_;
    277276#ifdef CLP_PLUS_ONE_MATRIX
    278      otherFlags_ = rhs.otherFlags_;
    279 #endif
    280      if (numberRows <= 0 || numberColumns <= 0) {
    281           startPositive_ = new CoinBigIndex [1];
    282           startPositive_[0] = 0;
    283      } else {
    284           numberColumns_ = numberColumns;
    285           numberRows_ = numberRows;
    286           const int * index1 = rhs.indices_;
    287           CoinBigIndex * startPositive1 = rhs.startPositive_;
     277  otherFlags_ = rhs.otherFlags_;
     278#endif
     279  if (numberRows <= 0 || numberColumns <= 0) {
     280    startPositive_ = new CoinBigIndex[1];
     281    startPositive_[0] = 0;
     282  } else {
     283    numberColumns_ = numberColumns;
     284    numberRows_ = numberRows;
     285    const int *index1 = rhs.indices_;
     286    CoinBigIndex *startPositive1 = rhs.startPositive_;
    288287
    289           int numberMinor = (!columnOrdered_) ? numberColumns_ : numberRows_;
    290           int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    291           int numberMinor1 = (!columnOrdered_) ? rhs.numberColumns_ : rhs.numberRows_;
    292           int numberMajor1 = (columnOrdered_) ? rhs.numberColumns_ : rhs.numberRows_;
    293           // Also swap incoming if not column ordered
    294           if (!columnOrdered_) {
    295                int temp1 = numberRows;
    296                numberRows = numberColumns;
    297                numberColumns = temp1;
    298                const int * temp2;
    299                temp2 = whichRow;
    300                whichRow = whichColumn;
    301                whichColumn = temp2;
     288    int numberMinor = (!columnOrdered_) ? numberColumns_ : numberRows_;
     289    int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     290    int numberMinor1 = (!columnOrdered_) ? rhs.numberColumns_ : rhs.numberRows_;
     291    int numberMajor1 = (columnOrdered_) ? rhs.numberColumns_ : rhs.numberRows_;
     292    // Also swap incoming if not column ordered
     293    if (!columnOrdered_) {
     294      int temp1 = numberRows;
     295      numberRows = numberColumns;
     296      numberColumns = temp1;
     297      const int *temp2;
     298      temp2 = whichRow;
     299      whichRow = whichColumn;
     300      whichColumn = temp2;
     301    }
     302    // Throw exception if rhs empty
     303    if (numberMajor1 <= 0 || numberMinor1 <= 0)
     304      throw CoinError("empty rhs", "subset constructor", "ClpPlusMinusOneMatrix");
     305    // Array to say if an old row is in new copy
     306    int *newRow = new int[numberMinor1];
     307    int iRow;
     308    for (iRow = 0; iRow < numberMinor1; iRow++)
     309      newRow[iRow] = -1;
     310    // and array for duplicating rows
     311    int *duplicateRow = new int[numberMinor];
     312    int numberBad = 0;
     313    for (iRow = 0; iRow < numberMinor; iRow++) {
     314      duplicateRow[iRow] = -1;
     315      int kRow = whichRow[iRow];
     316      if (kRow >= 0 && kRow < numberMinor1) {
     317        if (newRow[kRow] < 0) {
     318          // first time
     319          newRow[kRow] = iRow;
     320        } else {
     321          // duplicate
     322          int lastRow = newRow[kRow];
     323          newRow[kRow] = iRow;
     324          duplicateRow[iRow] = lastRow;
     325        }
     326      } else {
     327        // bad row
     328        numberBad++;
     329      }
     330    }
     331
     332    if (numberBad)
     333      throw CoinError("bad minor entries",
     334        "subset constructor", "ClpPlusMinusOneMatrix");
     335    // now get size and check columns
     336    CoinBigIndex size = 0;
     337    int iColumn;
     338    numberBad = 0;
     339    for (iColumn = 0; iColumn < numberMajor; iColumn++) {
     340      int kColumn = whichColumn[iColumn];
     341      if (kColumn >= 0 && kColumn < numberMajor1) {
     342        CoinBigIndex i;
     343        for (i = startPositive1[kColumn]; i < startPositive1[kColumn + 1]; i++) {
     344          int kRow = index1[i];
     345          kRow = newRow[kRow];
     346          while (kRow >= 0) {
     347            size++;
     348            kRow = duplicateRow[kRow];
    302349          }
    303           // Throw exception if rhs empty
    304           if (numberMajor1 <= 0 || numberMinor1 <= 0)
    305                throw CoinError("empty rhs", "subset constructor", "ClpPlusMinusOneMatrix");
    306           // Array to say if an old row is in new copy
    307           int * newRow = new int [numberMinor1];
    308           int iRow;
    309           for (iRow = 0; iRow < numberMinor1; iRow++)
    310                newRow[iRow] = -1;
    311           // and array for duplicating rows
    312           int * duplicateRow = new int [numberMinor];
    313           int numberBad = 0;
    314           for (iRow = 0; iRow < numberMinor; iRow++) {
    315                duplicateRow[iRow] = -1;
    316                int kRow = whichRow[iRow];
    317                if (kRow >= 0  && kRow < numberMinor1) {
    318                     if (newRow[kRow] < 0) {
    319                          // first time
    320                          newRow[kRow] = iRow;
    321                     } else {
    322                          // duplicate
    323                          int lastRow = newRow[kRow];
    324                          newRow[kRow] = iRow;
    325                          duplicateRow[iRow] = lastRow;
    326                     }
    327                } else {
    328                     // bad row
    329                     numberBad++;
    330                }
    331           }
    332 
    333           if (numberBad)
    334                throw CoinError("bad minor entries",
    335                                "subset constructor", "ClpPlusMinusOneMatrix");
    336           // now get size and check columns
    337           CoinBigIndex size = 0;
    338           int iColumn;
    339           numberBad = 0;
    340           for (iColumn = 0; iColumn < numberMajor; iColumn++) {
    341                int kColumn = whichColumn[iColumn];
    342                if (kColumn >= 0  && kColumn < numberMajor1) {
    343                     CoinBigIndex i;
    344                     for (i = startPositive1[kColumn]; i < startPositive1[kColumn+1]; i++) {
    345                          int kRow = index1[i];
    346                          kRow = newRow[kRow];
    347                          while (kRow >= 0) {
    348                               size++;
    349                               kRow = duplicateRow[kRow];
    350                          }
    351                     }
    352                } else {
    353                     // bad column
    354                     numberBad++;
    355                     printf("%d %d %d %d\n", iColumn, numberMajor, numberMajor1, kColumn);
    356                }
    357           }
    358           if (numberBad)
    359                throw CoinError("bad major entries",
    360                                "subset constructor", "ClpPlusMinusOneMatrix");
    361           // now create arrays
    362           startPositive_ = new CoinBigIndex [numberMajor+1];
    363           startNegative_ = new CoinBigIndex [numberMajor];
    364           indices_ = new int[size];
    365           // and fill them
    366           size = 0;
    367           startPositive_[0] = 0;
    368           CoinBigIndex * startNegative1 = rhs.startNegative_;
    369           for (iColumn = 0; iColumn < numberMajor; iColumn++) {
    370                int kColumn = whichColumn[iColumn];
    371                CoinBigIndex i;
    372                for (i = startPositive1[kColumn]; i < startNegative1[kColumn]; i++) {
    373                     int kRow = index1[i];
    374                     kRow = newRow[kRow];
    375                     while (kRow >= 0) {
    376                          indices_[size++] = kRow;
    377                          kRow = duplicateRow[kRow];
    378                     }
    379                }
    380                startNegative_[iColumn] = size;
    381                for (; i < startPositive1[kColumn+1]; i++) {
    382                     int kRow = index1[i];
    383                     kRow = newRow[kRow];
    384                     while (kRow >= 0) {
    385                          indices_[size++] = kRow;
    386                          kRow = duplicateRow[kRow];
    387                     }
    388                }
    389                startPositive_[iColumn+1] = size;
    390           }
    391           delete [] newRow;
    392           delete [] duplicateRow;
    393      }
    394      // Check valid
    395      checkValid(false);
    396 }
    397 
     350        }
     351      } else {
     352        // bad column
     353        numberBad++;
     354        printf("%d %d %d %d\n", iColumn, numberMajor, numberMajor1, kColumn);
     355      }
     356    }
     357    if (numberBad)
     358      throw CoinError("bad major entries",
     359        "subset constructor", "ClpPlusMinusOneMatrix");
     360    // now create arrays
     361    startPositive_ = new CoinBigIndex[numberMajor + 1];
     362    startNegative_ = new CoinBigIndex[numberMajor];
     363    indices_ = new int[size];
     364    // and fill them
     365    size = 0;
     366    startPositive_[0] = 0;
     367    CoinBigIndex *startNegative1 = rhs.startNegative_;
     368    for (iColumn = 0; iColumn < numberMajor; iColumn++) {
     369      int kColumn = whichColumn[iColumn];
     370      CoinBigIndex i;
     371      for (i = startPositive1[kColumn]; i < startNegative1[kColumn]; i++) {
     372        int kRow = index1[i];
     373        kRow = newRow[kRow];
     374        while (kRow >= 0) {
     375          indices_[size++] = kRow;
     376          kRow = duplicateRow[kRow];
     377        }
     378      }
     379      startNegative_[iColumn] = size;
     380      for (; i < startPositive1[kColumn + 1]; i++) {
     381        int kRow = index1[i];
     382        kRow = newRow[kRow];
     383        while (kRow >= 0) {
     384          indices_[size++] = kRow;
     385          kRow = duplicateRow[kRow];
     386        }
     387      }
     388      startPositive_[iColumn + 1] = size;
     389    }
     390    delete[] newRow;
     391    delete[] duplicateRow;
     392  }
     393  // Check valid
     394  checkValid(false);
     395}
    398396
    399397/* Returns a new matrix in reverse order without gaps */
     
    401399ClpPlusMinusOneMatrix::reverseOrderedCopy() const
    402400{
    403      int numberMinor = (!columnOrdered_) ? numberColumns_ : numberRows_;
    404      int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    405      // count number in each row/column
    406      CoinBigIndex * tempP = new CoinBigIndex [numberMinor];
    407      CoinBigIndex * tempN = new CoinBigIndex [numberMinor];
    408      memset(tempP, 0, numberMinor * sizeof(CoinBigIndex));
    409      memset(tempN, 0, numberMinor * sizeof(CoinBigIndex));
    410      CoinBigIndex j = 0;
    411      int i;
    412      for (i = 0; i < numberMajor; i++) {
    413           for (; j < startNegative_[i]; j++) {
    414                int iRow = indices_[j];
    415                tempP[iRow]++;
    416           }
    417           for (; j < startPositive_[i+1]; j++) {
    418                int iRow = indices_[j];
    419                tempN[iRow]++;
    420           }
    421      }
    422      int * newIndices = new int [startPositive_[numberMajor]];
    423      CoinBigIndex * newP = new CoinBigIndex [numberMinor+1];
    424      CoinBigIndex * newN = new CoinBigIndex[numberMinor];
    425      int iRow;
    426      j = 0;
    427      // do starts
    428      for (iRow = 0; iRow < numberMinor; iRow++) {
    429           newP[iRow] = j;
    430           j += tempP[iRow];
    431           tempP[iRow] = newP[iRow];
    432           newN[iRow] = j;
    433           j += tempN[iRow];
    434           tempN[iRow] = newN[iRow];
    435      }
    436      newP[numberMinor] = j;
    437      j = 0;
    438      for (i = 0; i < numberMajor; i++) {
    439           for (; j < startNegative_[i]; j++) {
    440                int iRow = indices_[j];
    441                CoinBigIndex put = tempP[iRow];
    442                newIndices[put++] = i;
    443                tempP[iRow] = put;
    444           }
    445           for (; j < startPositive_[i+1]; j++) {
    446                int iRow = indices_[j];
    447                CoinBigIndex put = tempN[iRow];
    448                newIndices[put++] = i;
    449                tempN[iRow] = put;
    450           }
    451      }
    452      delete [] tempP;
    453      delete [] tempN;
    454      ClpPlusMinusOneMatrix * newCopy = new ClpPlusMinusOneMatrix();
    455      newCopy->passInCopy(numberMinor, numberMajor,
    456                          !columnOrdered_, newIndices, newP, newN);
    457      return newCopy;
     401  int numberMinor = (!columnOrdered_) ? numberColumns_ : numberRows_;
     402  int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     403  // count number in each row/column
     404  CoinBigIndex *tempP = new CoinBigIndex[numberMinor];
     405  CoinBigIndex *tempN = new CoinBigIndex[numberMinor];
     406  memset(tempP, 0, numberMinor * sizeof(CoinBigIndex));
     407  memset(tempN, 0, numberMinor * sizeof(CoinBigIndex));
     408  CoinBigIndex j = 0;
     409  int i;
     410  for (i = 0; i < numberMajor; i++) {
     411    for (; j < startNegative_[i]; j++) {
     412      int iRow = indices_[j];
     413      tempP[iRow]++;
     414    }
     415    for (; j < startPositive_[i + 1]; j++) {
     416      int iRow = indices_[j];
     417      tempN[iRow]++;
     418    }
     419  }
     420  int *newIndices = new int[startPositive_[numberMajor]];
     421  CoinBigIndex *newP = new CoinBigIndex[numberMinor + 1];
     422  CoinBigIndex *newN = new CoinBigIndex[numberMinor];
     423  int iRow;
     424  j = 0;
     425  // do starts
     426  for (iRow = 0; iRow < numberMinor; iRow++) {
     427    newP[iRow] = j;
     428    j += tempP[iRow];
     429    tempP[iRow] = newP[iRow];
     430    newN[iRow] = j;
     431    j += tempN[iRow];
     432    tempN[iRow] = newN[iRow];
     433  }
     434  newP[numberMinor] = j;
     435  j = 0;
     436  for (i = 0; i < numberMajor; i++) {
     437    for (; j < startNegative_[i]; j++) {
     438      int iRow = indices_[j];
     439      CoinBigIndex put = tempP[iRow];
     440      newIndices[put++] = i;
     441      tempP[iRow] = put;
     442    }
     443    for (; j < startPositive_[i + 1]; j++) {
     444      int iRow = indices_[j];
     445      CoinBigIndex put = tempN[iRow];
     446      newIndices[put++] = i;
     447      tempN[iRow] = put;
     448    }
     449  }
     450  delete[] tempP;
     451  delete[] tempN;
     452  ClpPlusMinusOneMatrix *newCopy = new ClpPlusMinusOneMatrix();
     453  newCopy->passInCopy(numberMinor, numberMajor,
     454    !columnOrdered_, newIndices, newP, newN);
     455  return newCopy;
    458456}
    459457//static bool doPlusOnes=true;
    460458//unscaled versions
    461 void
    462 ClpPlusMinusOneMatrix::times(double scalar,
    463                              const double * COIN_RESTRICT x, double * COIN_RESTRICT y) const
    464 {
    465      int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    466      int i;
    467      CoinBigIndex j;
    468      assert (columnOrdered_);
     459void ClpPlusMinusOneMatrix::times(double scalar,
     460  const double *COIN_RESTRICT x, double *COIN_RESTRICT y) const
     461{
     462  int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     463  int i;
     464  CoinBigIndex j;
     465  assert(columnOrdered_);
    469466#ifdef CLP_PLUS_ONE_MATRIX
    470      if ((otherFlags_&1)==0||!doPlusOnes) {
    471 #endif
    472        for (i = 0; i < numberMajor; i++) {
    473           double value = scalar * x[i];
    474           if (value) {
    475                for (j = startPositive_[i]; j < startNegative_[i]; j++) {
    476                     int iRow = indices_[j];
    477                     y[iRow] += value;
    478                }
    479                for (; j < startPositive_[i+1]; j++) {
    480                     int iRow = indices_[j];
    481                     y[iRow] -= value;
    482                }
    483           }
    484        }
     467  if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     468#endif
     469    for (i = 0; i < numberMajor; i++) {
     470      double value = scalar * x[i];
     471      if (value) {
     472        for (j = startPositive_[i]; j < startNegative_[i]; j++) {
     473          int iRow = indices_[j];
     474          y[iRow] += value;
     475        }
     476        for (; j < startPositive_[i + 1]; j++) {
     477          int iRow = indices_[j];
     478          y[iRow] -= value;
     479        }
     480      }
     481    }
    485482#ifdef CLP_PLUS_ONE_MATRIX
    486      } else {
    487        // plus one
    488        oneit(0);
    489        for (i = 0; i < numberMajor; i++) {
    490           double value = scalar * x[i];
    491           if (value) {
    492                for (j = startPositive_[i]; j < startPositive_[i+1]; j++) {
    493                     int iRow = indices_[j];
    494                     y[iRow] += value;
    495                }
    496           }
    497        }
    498      }
    499 #endif
    500 }
    501 void
    502 ClpPlusMinusOneMatrix::transposeTimes(double scalar,
    503                                       const double * COIN_RESTRICT x, double * COIN_RESTRICT y) const
    504 {
    505      int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    506      int i;
    507      CoinBigIndex j = 0;
    508      assert (columnOrdered_);
     483  } else {
     484    // plus one
     485    oneit(0);
     486    for (i = 0; i < numberMajor; i++) {
     487      double value = scalar * x[i];
     488      if (value) {
     489        for (j = startPositive_[i]; j < startPositive_[i + 1]; j++) {
     490          int iRow = indices_[j];
     491          y[iRow] += value;
     492        }
     493      }
     494    }
     495  }
     496#endif
     497}
     498void ClpPlusMinusOneMatrix::transposeTimes(double scalar,
     499  const double *COIN_RESTRICT x, double *COIN_RESTRICT y) const
     500{
     501  int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     502  int i;
     503  CoinBigIndex j = 0;
     504  assert(columnOrdered_);
    509505#ifdef CLP_PLUS_ONE_MATRIX
    510      if ((otherFlags_&1)==0||!doPlusOnes) {
    511 #endif
    512        for (i = 0; i < numberMajor; i++) {
    513           double value = 0.0;
    514           for (; j < startNegative_[i]; j++) {
    515                int iRow = indices_[j];
    516                value += x[iRow];
    517           }
    518           for (; j < startPositive_[i+1]; j++) {
    519                int iRow = indices_[j];
    520                value -= x[iRow];
    521           }
    522           y[i] += scalar * value;
    523        }
     506  if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     507#endif
     508    for (i = 0; i < numberMajor; i++) {
     509      double value = 0.0;
     510      for (; j < startNegative_[i]; j++) {
     511        int iRow = indices_[j];
     512        value += x[iRow];
     513      }
     514      for (; j < startPositive_[i + 1]; j++) {
     515        int iRow = indices_[j];
     516        value -= x[iRow];
     517      }
     518      y[i] += scalar * value;
     519    }
    524520#ifdef CLP_PLUS_ONE_MATRIX
    525      } else {
    526        // plus one
    527        oneit(1);
    528        for (i = 0; i < numberMajor; i++) {
    529           double value = 0.0;
    530           for (; j < startPositive_[i+1]; j++) {
    531                int iRow = indices_[j];
    532                value += x[iRow];
    533           }
    534           y[i] += scalar * value;
    535        }
    536      }
    537 #endif
    538 }
    539 void
    540 ClpPlusMinusOneMatrix::times(double scalar,
    541                              const double * COIN_RESTRICT x, double * COIN_RESTRICT y,
    542                              const double * /*rowScale*/,
    543                              const double * /*columnScale*/) const
    544 {
    545      // we know it is not scaled
    546      times(scalar, x, y);
    547 }
    548 void
    549 ClpPlusMinusOneMatrix::transposeTimes( double scalar,
    550                                        const double * COIN_RESTRICT x, double * COIN_RESTRICT y,
    551                                        const double * /*rowScale*/,
    552                                        const double * /*columnScale*/,
    553                                        double * /*spare*/) const
    554 {
    555      // we know it is not scaled
    556      transposeTimes(scalar, x, y);
     521  } else {
     522    // plus one
     523    oneit(1);
     524    for (i = 0; i < numberMajor; i++) {
     525      double value = 0.0;
     526      for (; j < startPositive_[i + 1]; j++) {
     527        int iRow = indices_[j];
     528        value += x[iRow];
     529      }
     530      y[i] += scalar * value;
     531    }
     532  }
     533#endif
     534}
     535void ClpPlusMinusOneMatrix::times(double scalar,
     536  const double *COIN_RESTRICT x, double *COIN_RESTRICT y,
     537  const double * /*rowScale*/,
     538  const double * /*columnScale*/) const
     539{
     540  // we know it is not scaled
     541  times(scalar, x, y);
     542}
     543void ClpPlusMinusOneMatrix::transposeTimes(double scalar,
     544  const double *COIN_RESTRICT x, double *COIN_RESTRICT y,
     545  const double * /*rowScale*/,
     546  const double * /*columnScale*/,
     547  double * /*spare*/) const
     548{
     549  // we know it is not scaled
     550  transposeTimes(scalar, x, y);
    557551}
    558552/* Return <code>x * A + y</code> in <code>z</code>.
    559553        Squashes small elements and knows about ClpSimplex */
    560 void
    561 ClpPlusMinusOneMatrix::transposeTimes(const ClpSimplex * model, double scalar,
    562                                       const CoinIndexedVector * rowArray,
    563                                       CoinIndexedVector * y,
    564                                       CoinIndexedVector * columnArray) const
    565 {
    566      // we know it is not scaled
    567      columnArray->clear();
    568      double * COIN_RESTRICT pi = rowArray->denseVector();
    569      int numberNonZero = 0;
    570      int * COIN_RESTRICT index = columnArray->getIndices();
    571      double * COIN_RESTRICT array = columnArray->denseVector();
    572      int numberInRowArray = rowArray->getNumElements();
    573      // maybe I need one in OsiSimplex
    574      double zeroTolerance = model->zeroTolerance();
    575      int numberRows = model->numberRows();
    576      bool packed = rowArray->packedMode();
     554void ClpPlusMinusOneMatrix::transposeTimes(const ClpSimplex *model, double scalar,
     555  const CoinIndexedVector *rowArray,
     556  CoinIndexedVector *y,
     557  CoinIndexedVector *columnArray) const
     558{
     559  // we know it is not scaled
     560  columnArray->clear();
     561  double *COIN_RESTRICT pi = rowArray->denseVector();
     562  int numberNonZero = 0;
     563  int *COIN_RESTRICT index = columnArray->getIndices();
     564  double *COIN_RESTRICT array = columnArray->denseVector();
     565  int numberInRowArray = rowArray->getNumElements();
     566  // maybe I need one in OsiSimplex
     567  double zeroTolerance = model->zeroTolerance();
     568  int numberRows = model->numberRows();
     569  bool packed = rowArray->packedMode();
    577570#ifndef NO_RTTI
    578      ClpPlusMinusOneMatrix* rowCopy =
    579           dynamic_cast< ClpPlusMinusOneMatrix*>(model->rowCopy());
     571  ClpPlusMinusOneMatrix *rowCopy = dynamic_cast< ClpPlusMinusOneMatrix * >(model->rowCopy());
    580572#else
    581      ClpPlusMinusOneMatrix* rowCopy =
    582           static_cast< ClpPlusMinusOneMatrix*>(model->rowCopy());
    583 #endif
    584      double factor = 0.3;
    585      // We may not want to do by row if there may be cache problems
    586      int numberColumns = model->numberColumns();
    587      // It would be nice to find L2 cache size - for moment 512K
    588      // Be slightly optimistic
    589      if (numberColumns * sizeof(double) > 1000000) {
    590           if (numberRows * 10 < numberColumns)
    591                factor = 0.1;
    592           else if (numberRows * 4 < numberColumns)
    593                factor = 0.15;
    594           else if (numberRows * 2 < numberColumns)
    595                factor = 0.2;
    596      }
    597      if (numberInRowArray > factor * numberRows || !rowCopy) {
    598           assert (!y->getNumElements());
    599           // do by column
    600           // Need to expand if packed mode
    601           int iColumn;
    602           CoinBigIndex j = 0;
    603           assert (columnOrdered_);
    604           if (packed) {
    605                // need to expand pi into y
    606                assert(y->capacity() >= numberRows);
    607                double * COIN_RESTRICT piOld = pi;
    608                pi = y->denseVector();
    609                const int * COIN_RESTRICT whichRow = rowArray->getIndices();
    610                int i;
    611                // modify pi so can collapse to one loop
    612                for (i = 0; i < numberInRowArray; i++) {
    613                     int iRow = whichRow[i];
    614                     pi[iRow] = scalar * piOld[i];
    615                }
     573  ClpPlusMinusOneMatrix *rowCopy = static_cast< ClpPlusMinusOneMatrix * >(model->rowCopy());
     574#endif
     575  double factor = 0.3;
     576  // We may not want to do by row if there may be cache problems
     577  int numberColumns = model->numberColumns();
     578  // It would be nice to find L2 cache size - for moment 512K
     579  // Be slightly optimistic
     580  if (numberColumns * sizeof(double) > 1000000) {
     581    if (numberRows * 10 < numberColumns)
     582      factor = 0.1;
     583    else if (numberRows * 4 < numberColumns)
     584      factor = 0.15;
     585    else if (numberRows * 2 < numberColumns)
     586      factor = 0.2;
     587  }
     588  if (numberInRowArray > factor * numberRows || !rowCopy) {
     589    assert(!y->getNumElements());
     590    // do by column
     591    // Need to expand if packed mode
     592    int iColumn;
     593    CoinBigIndex j = 0;
     594    assert(columnOrdered_);
     595    if (packed) {
     596      // need to expand pi into y
     597      assert(y->capacity() >= numberRows);
     598      double *COIN_RESTRICT piOld = pi;
     599      pi = y->denseVector();
     600      const int *COIN_RESTRICT whichRow = rowArray->getIndices();
     601      int i;
     602      // modify pi so can collapse to one loop
     603      for (i = 0; i < numberInRowArray; i++) {
     604        int iRow = whichRow[i];
     605        pi[iRow] = scalar * piOld[i];
     606      }
    616607#ifdef CLP_PLUS_ONE_MATRIX
    617                if ((otherFlags_&1)==0||!doPlusOnes) {
    618 #endif
    619                 for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    620                     double value = 0.0;
    621                     for (; j < startNegative_[iColumn]; j++) {
    622                          int iRow = indices_[j];
    623                          value += pi[iRow];
    624                     }
    625                     for (; j < startPositive_[iColumn+1]; j++) {
    626                          int iRow = indices_[j];
    627                          value -= pi[iRow];
    628                     }
    629                     if (fabs(value) > zeroTolerance) {
    630                          array[numberNonZero] = value;
    631                          index[numberNonZero++] = iColumn;
    632                     }
    633                 }
     608      if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     609#endif
     610        for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     611          double value = 0.0;
     612          for (; j < startNegative_[iColumn]; j++) {
     613            int iRow = indices_[j];
     614            value += pi[iRow];
     615          }
     616          for (; j < startPositive_[iColumn + 1]; j++) {
     617            int iRow = indices_[j];
     618            value -= pi[iRow];
     619          }
     620          if (fabs(value) > zeroTolerance) {
     621            array[numberNonZero] = value;
     622            index[numberNonZero++] = iColumn;
     623          }
     624        }
    634625#ifdef CLP_PLUS_ONE_MATRIX
    635                } else {
    636                  // plus one
    637                  oneit(2);
    638                  for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    639                     double value = 0.0;
    640                     for (; j < startPositive_[iColumn+1]; j++) {
    641                          int iRow = indices_[j];
    642                          value += pi[iRow];
    643                     }
    644                     if (fabs(value) > zeroTolerance) {
    645                          array[numberNonZero] = value;
    646                          index[numberNonZero++] = iColumn;
    647                     }
    648                  }
    649                }
    650 #endif
    651                for (i = 0; i < numberInRowArray; i++) {
    652                     int iRow = whichRow[i];
    653                     pi[iRow] = 0.0;
    654                }
    655           } else {
    656                  for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    657                     double value = 0.0;
    658                     for (; j < startNegative_[iColumn]; j++) {
    659                          int iRow = indices_[j];
    660                          value += pi[iRow];
    661                     }
    662                     for (; j < startPositive_[iColumn+1]; j++) {
    663                          int iRow = indices_[j];
    664                          value -= pi[iRow];
    665                     }
    666                     value *= scalar;
    667                     if (fabs(value) > zeroTolerance) {
    668                          index[numberNonZero++] = iColumn;
    669                          array[iColumn] = value;
    670                     }
    671                  }
     626      } else {
     627        // plus one
     628        oneit(2);
     629        for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     630          double value = 0.0;
     631          for (; j < startPositive_[iColumn + 1]; j++) {
     632            int iRow = indices_[j];
     633            value += pi[iRow];
    672634          }
    673           columnArray->setNumElements(numberNonZero);
    674      } else {
    675           // do by row
    676           rowCopy->transposeTimesByRow(model, scalar, rowArray, y, columnArray);
    677      }
     635          if (fabs(value) > zeroTolerance) {
     636            array[numberNonZero] = value;
     637            index[numberNonZero++] = iColumn;
     638          }
     639        }
     640      }
     641#endif
     642      for (i = 0; i < numberInRowArray; i++) {
     643        int iRow = whichRow[i];
     644        pi[iRow] = 0.0;
     645      }
     646    } else {
     647      for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     648        double value = 0.0;
     649        for (; j < startNegative_[iColumn]; j++) {
     650          int iRow = indices_[j];
     651          value += pi[iRow];
     652        }
     653        for (; j < startPositive_[iColumn + 1]; j++) {
     654          int iRow = indices_[j];
     655          value -= pi[iRow];
     656        }
     657        value *= scalar;
     658        if (fabs(value) > zeroTolerance) {
     659          index[numberNonZero++] = iColumn;
     660          array[iColumn] = value;
     661        }
     662      }
     663    }
     664    columnArray->setNumElements(numberNonZero);
     665  } else {
     666    // do by row
     667    rowCopy->transposeTimesByRow(model, scalar, rowArray, y, columnArray);
     668  }
    678669}
    679670/* Return <code>x * A + y</code> in <code>z</code>.
    680671        Squashes small elements and knows about ClpSimplex */
    681 void
    682 ClpPlusMinusOneMatrix::transposeTimesByRow(const ClpSimplex * model, double scalar,
    683           const CoinIndexedVector * rowArray,
    684           CoinIndexedVector * y,
    685           CoinIndexedVector * columnArray) const
    686 {
    687      columnArray->clear();
    688      double * COIN_RESTRICT pi = rowArray->denseVector();
    689      int numberNonZero = 0;
    690      int * COIN_RESTRICT index = columnArray->getIndices();
    691      double * COIN_RESTRICT array = columnArray->denseVector();
    692      int numberInRowArray = rowArray->getNumElements();
    693      // maybe I need one in OsiSimplex
    694      double zeroTolerance = model->zeroTolerance();
    695      const int * COIN_RESTRICT column = indices_;
    696      const CoinBigIndex * COIN_RESTRICT startPositive = startPositive_;
    697      const CoinBigIndex * COIN_RESTRICT startNegative = startNegative_;
    698      const int * COIN_RESTRICT whichRow = rowArray->getIndices();
    699      bool packed = rowArray->packedMode();
    700      if (numberInRowArray > 2) {
    701           // do by rows
    702           int iRow;
    703           double * COIN_RESTRICT markVector = y->denseVector(); // probably empty .. but
    704           int numberOriginal = 0;
    705           int i;
    706           if (packed) {
    707                CoinBigIndex numberCovered=0;
    708                int numberColumns = getNumCols();
    709                bool sparse=true;
    710                int target=1*numberColumns;
    711                for (int i = 0; i < numberInRowArray; i++) {
    712                  int iRow = whichRow[i];
    713                  numberCovered += startPositive[iRow+1] - startPositive[iRow];
    714                  if (numberCovered>target) {
    715                    sparse=false;
    716                    break;
    717                  }
    718                }
    719                numberNonZero = 0;
    720                if (sparse) {
    721                  // and set up mark as char array
    722                  char * COIN_RESTRICT marked = reinterpret_cast<char *> (index + columnArray->capacity());
    723                  double * COIN_RESTRICT array2 = y->denseVector();
     672void ClpPlusMinusOneMatrix::transposeTimesByRow(const ClpSimplex *model, double scalar,
     673  const CoinIndexedVector *rowArray,
     674  CoinIndexedVector *y,
     675  CoinIndexedVector *columnArray) const
     676{
     677  columnArray->clear();
     678  double *COIN_RESTRICT pi = rowArray->denseVector();
     679  int numberNonZero = 0;
     680  int *COIN_RESTRICT index = columnArray->getIndices();
     681  double *COIN_RESTRICT array = columnArray->denseVector();
     682  int numberInRowArray = rowArray->getNumElements();
     683  // maybe I need one in OsiSimplex
     684  double zeroTolerance = model->zeroTolerance();
     685  const int *COIN_RESTRICT column = indices_;
     686  const CoinBigIndex *COIN_RESTRICT startPositive = startPositive_;
     687  const CoinBigIndex *COIN_RESTRICT startNegative = startNegative_;
     688  const int *COIN_RESTRICT whichRow = rowArray->getIndices();
     689  bool packed = rowArray->packedMode();
     690  if (numberInRowArray > 2) {
     691    // do by rows
     692    int iRow;
     693    double *COIN_RESTRICT markVector = y->denseVector(); // probably empty .. but
     694    int numberOriginal = 0;
     695    int i;
     696    if (packed) {
     697      CoinBigIndex numberCovered = 0;
     698      int numberColumns = getNumCols();
     699      bool sparse = true;
     700      int target = 1 * numberColumns;
     701      for (int i = 0; i < numberInRowArray; i++) {
     702        int iRow = whichRow[i];
     703        numberCovered += startPositive[iRow + 1] - startPositive[iRow];
     704        if (numberCovered > target) {
     705          sparse = false;
     706          break;
     707        }
     708      }
     709      numberNonZero = 0;
     710      if (sparse) {
     711        // and set up mark as char array
     712        char *COIN_RESTRICT marked = reinterpret_cast< char * >(index + columnArray->capacity());
     713        double *COIN_RESTRICT array2 = y->denseVector();
    724714#ifdef CLP_DEBUG
    725                  int numberColumns = model->numberColumns();
    726                  for (int i = 0; i < numberColumns; i++) {
    727                     assert(!marked[i]);
    728                     assert(!array2[i]);
    729                  }
     715        int numberColumns = model->numberColumns();
     716        for (int i = 0; i < numberColumns; i++) {
     717          assert(!marked[i]);
     718          assert(!array2[i]);
     719        }
    730720#endif
    731721#ifdef CLP_PLUS_ONE_MATRIX
    732                  if ((otherFlags_&1)==0||!doPlusOnes) {
    733 #endif
    734                    for (int i = 0; i < numberInRowArray; i++) {
    735                      iRow = whichRow[i];
    736                      double value = pi[i] * scalar;
    737                      CoinBigIndex j;
    738                      for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
    739                          int iColumn = column[j];
    740                          if (!marked[iColumn]) {
    741                               marked[iColumn] = 1;
    742                               index[numberNonZero++] = iColumn;
    743                          }
    744                          array2[iColumn] += value;
    745                       }
    746                       for (j = startNegative[iRow]; j < startPositive[iRow+1]; j++) {
    747                          int iColumn = column[j];
    748                          if (!marked[iColumn]) {
    749                               marked[iColumn] = 1;
    750                               index[numberNonZero++] = iColumn;
    751                          }
    752                          array2[iColumn] -= value;
    753                       }
    754                     }
     722        if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     723#endif
     724          for (int i = 0; i < numberInRowArray; i++) {
     725            iRow = whichRow[i];
     726            double value = pi[i] * scalar;
     727            CoinBigIndex j;
     728            for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
     729              int iColumn = column[j];
     730              if (!marked[iColumn]) {
     731                marked[iColumn] = 1;
     732                index[numberNonZero++] = iColumn;
     733              }
     734              array2[iColumn] += value;
     735            }
     736            for (j = startNegative[iRow]; j < startPositive[iRow + 1]; j++) {
     737              int iColumn = column[j];
     738              if (!marked[iColumn]) {
     739                marked[iColumn] = 1;
     740                index[numberNonZero++] = iColumn;
     741              }
     742              array2[iColumn] -= value;
     743            }
     744          }
    755745#ifdef CLP_PLUS_ONE_MATRIX
    756                  } else {
    757                     // plus one
    758                     oneit(4);
    759                    for (int i = 0; i < numberInRowArray; i++) {
    760                      iRow = whichRow[i];
    761                      double value = pi[i] * scalar;
    762                      CoinBigIndex j;
    763                      for (j = startPositive[iRow]; j < startPositive[iRow+1]; j++) {
    764                          int iColumn = column[j];
    765                          if (!marked[iColumn]) {
    766                               marked[iColumn] = 1;
    767                               index[numberNonZero++] = iColumn;
    768                          }
    769                          array2[iColumn] += value;
    770                       }
    771                     }
    772                 }
    773 #endif
    774                 // get rid of tiny values and zero out marked
    775                  numberOriginal = numberNonZero;
    776                  numberNonZero = 0;
    777                  for (int i = 0; i < numberOriginal; i++) {
    778                     int iColumn = index[i];
    779                     if (marked[iColumn]) {
    780                          double value = array2[iColumn];
    781                          array2[iColumn] = 0.0;
    782                          marked[iColumn] = 0;
    783                          if (fabs(value) > zeroTolerance) {
    784                               array[numberNonZero] = value;
    785                               index[numberNonZero++] = iColumn;
    786                          }
    787                     }
    788                  }
    789                } else {
    790                  // not sparse
     746        } else {
     747          // plus one
     748          oneit(4);
     749          for (int i = 0; i < numberInRowArray; i++) {
     750            iRow = whichRow[i];
     751            double value = pi[i] * scalar;
     752            CoinBigIndex j;
     753            for (j = startPositive[iRow]; j < startPositive[iRow + 1]; j++) {
     754              int iColumn = column[j];
     755              if (!marked[iColumn]) {
     756                marked[iColumn] = 1;
     757                index[numberNonZero++] = iColumn;
     758              }
     759              array2[iColumn] += value;
     760            }
     761          }
     762        }
     763#endif
     764        // get rid of tiny values and zero out marked
     765        numberOriginal = numberNonZero;
     766        numberNonZero = 0;
     767        for (int i = 0; i < numberOriginal; i++) {
     768          int iColumn = index[i];
     769          if (marked[iColumn]) {
     770            double value = array2[iColumn];
     771            array2[iColumn] = 0.0;
     772            marked[iColumn] = 0;
     773            if (fabs(value) > zeroTolerance) {
     774              array[numberNonZero] = value;
     775              index[numberNonZero++] = iColumn;
     776            }
     777          }
     778        }
     779      } else {
     780        // not sparse
    791781#ifdef CLP_PLUS_ONE_MATRIX
    792                  if ((otherFlags_&1)==0||!doPlusOnes) {
    793 #endif
    794                    for (int i = 0; i < numberInRowArray; i++) {
    795                      iRow = whichRow[i];
    796                      double value = pi[i] * scalar;
    797                      CoinBigIndex j;
    798                      for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
    799                          int iColumn = column[j];
    800                          array[iColumn] += value;
    801                      }
    802                      for (j = startNegative[iRow]; j < startPositive[iRow+1]; j++) {
    803                          int iColumn = column[j];
    804                          array[iColumn] -= value;
    805                      }
    806                    }
     782        if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     783#endif
     784          for (int i = 0; i < numberInRowArray; i++) {
     785            iRow = whichRow[i];
     786            double value = pi[i] * scalar;
     787            CoinBigIndex j;
     788            for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
     789              int iColumn = column[j];
     790              array[iColumn] += value;
     791            }
     792            for (j = startNegative[iRow]; j < startPositive[iRow + 1]; j++) {
     793              int iColumn = column[j];
     794              array[iColumn] -= value;
     795            }
     796          }
    807797#ifdef CLP_PLUS_ONE_MATRIX
    808                 } else {
    809                   // plus one
    810                   oneit(5);
    811                    for (int i = 0; i < numberInRowArray; i++) {
    812                      iRow = whichRow[i];
    813                      double value = pi[i] * scalar;
    814                      CoinBigIndex j;
    815                      for (j = startPositive[iRow]; j < startPositive[iRow+1]; j++) {
    816                          int iColumn = column[j];
    817                          array[iColumn] += value;
    818                      }
    819                    }
    820                  }
    821 #endif
    822                  // get rid of tiny values and count
    823                  for (int i = 0; i < numberColumns; i++) {
    824                    double value = array[i];
    825                    if (value) {
    826                      array[i] = 0.0;
    827                      if (fabs(value) > zeroTolerance) {
    828                        array[numberNonZero] = value;
    829                        index[numberNonZero++] = i;
    830                      }
    831                    }
    832                  }
    833                }
     798        } else {
     799          // plus one
     800          oneit(5);
     801          for (int i = 0; i < numberInRowArray; i++) {
     802            iRow = whichRow[i];
     803            double value = pi[i] * scalar;
     804            CoinBigIndex j;
     805            for (j = startPositive[iRow]; j < startPositive[iRow + 1]; j++) {
     806              int iColumn = column[j];
     807              array[iColumn] += value;
     808            }
     809          }
     810        }
     811#endif
     812        // get rid of tiny values and count
     813        for (int i = 0; i < numberColumns; i++) {
     814          double value = array[i];
     815          if (value) {
     816            array[i] = 0.0;
     817            if (fabs(value) > zeroTolerance) {
     818              array[numberNonZero] = value;
     819              index[numberNonZero++] = i;
     820            }
     821          }
     822        }
     823      }
     824    } else {
     825      numberNonZero = 0;
     826      // and set up mark as char array
     827      char *COIN_RESTRICT marked = reinterpret_cast< char * >(markVector);
     828      for (i = 0; i < numberOriginal; i++) {
     829        int iColumn = index[i];
     830        marked[iColumn] = 0;
     831      }
     832      for (i = 0; i < numberInRowArray; i++) {
     833        iRow = whichRow[i];
     834        double value = pi[iRow] * scalar;
     835        CoinBigIndex j;
     836        for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
     837          int iColumn = column[j];
     838          if (!marked[iColumn]) {
     839            marked[iColumn] = 1;
     840            index[numberNonZero++] = iColumn;
     841          }
     842          array[iColumn] += value;
     843        }
     844        for (j = startNegative[iRow]; j < startPositive[iRow + 1]; j++) {
     845          int iColumn = column[j];
     846          if (!marked[iColumn]) {
     847            marked[iColumn] = 1;
     848            index[numberNonZero++] = iColumn;
     849          }
     850          array[iColumn] -= value;
     851        }
     852      }
     853      // get rid of tiny values and zero out marked
     854      numberOriginal = numberNonZero;
     855      numberNonZero = 0;
     856      for (i = 0; i < numberOriginal; i++) {
     857        int iColumn = index[i];
     858        marked[iColumn] = 0;
     859        if (fabs(array[iColumn]) > zeroTolerance) {
     860          index[numberNonZero++] = iColumn;
     861        } else {
     862          array[iColumn] = 0.0;
     863        }
     864      }
     865    }
     866  } else if (numberInRowArray == 2) {
     867    /* do by rows when two rows (do longer first when not packed
     868             and shorter first if packed */
     869    int iRow0 = whichRow[0];
     870    int iRow1 = whichRow[1];
     871    CoinBigIndex j;
     872    if (packed) {
     873      double pi0 = pi[0];
     874      double pi1 = pi[1];
     875      if (startPositive[iRow0 + 1] - startPositive[iRow0] > startPositive[iRow1 + 1] - startPositive[iRow1]) {
     876        int temp = iRow0;
     877        iRow0 = iRow1;
     878        iRow1 = temp;
     879        pi0 = pi1;
     880        pi1 = pi[0];
     881      }
     882      // and set up mark as char array
     883      char *COIN_RESTRICT marked = reinterpret_cast< char * >(index + columnArray->capacity());
     884      int *COIN_RESTRICT lookup = y->getIndices();
     885      double value = pi0 * scalar;
     886      int numberOriginal;
     887#ifdef CLP_PLUS_ONE_MATRIX
     888      if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     889#endif
     890        for (j = startPositive[iRow0]; j < startNegative[iRow0]; j++) {
     891          int iColumn = column[j];
     892          array[numberNonZero] = value;
     893          marked[iColumn] = 1;
     894          lookup[iColumn] = numberNonZero;
     895          index[numberNonZero++] = iColumn;
     896        }
     897        for (j = startNegative[iRow0]; j < startPositive[iRow0 + 1]; j++) {
     898          int iColumn = column[j];
     899          array[numberNonZero] = -value;
     900          marked[iColumn] = 1;
     901          lookup[iColumn] = numberNonZero;
     902          index[numberNonZero++] = iColumn;
     903        }
     904        numberOriginal = numberNonZero;
     905        value = pi1 * scalar;
     906        for (j = startPositive[iRow1]; j < startNegative[iRow1]; j++) {
     907          int iColumn = column[j];
     908          if (marked[iColumn]) {
     909            int iLookup = lookup[iColumn];
     910            array[iLookup] += value;
    834911          } else {
    835                numberNonZero = 0;
    836                // and set up mark as char array
    837                char * COIN_RESTRICT marked = reinterpret_cast<char *> (markVector);
    838                for (i = 0; i < numberOriginal; i++) {
    839                     int iColumn = index[i];
    840                     marked[iColumn] = 0;
    841                }
    842                for (i = 0; i < numberInRowArray; i++) {
    843                  iRow = whichRow[i];
    844                  double value = pi[iRow] * scalar;
    845                  CoinBigIndex j;
    846                  for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
    847                    int iColumn = column[j];
    848                    if (!marked[iColumn]) {
    849                      marked[iColumn] = 1;
    850                      index[numberNonZero++] = iColumn;
    851                    }
    852                    array[iColumn] += value;
    853                  }
    854                  for (j = startNegative[iRow]; j < startPositive[iRow+1]; j++) {
    855                    int iColumn = column[j];
    856                    if (!marked[iColumn]) {
    857                      marked[iColumn] = 1;
    858                      index[numberNonZero++] = iColumn;
    859                    }
    860                    array[iColumn] -= value;
    861                  }
    862                }
    863                // get rid of tiny values and zero out marked
    864                numberOriginal = numberNonZero;
    865                numberNonZero = 0;
    866                for (i = 0; i < numberOriginal; i++) {
    867                     int iColumn = index[i];
    868                     marked[iColumn] = 0;
    869                     if (fabs(array[iColumn]) > zeroTolerance) {
    870                          index[numberNonZero++] = iColumn;
    871                     } else {
    872                          array[iColumn] = 0.0;
    873                     }
    874                }
     912            if (fabs(value) > zeroTolerance) {
     913              array[numberNonZero] = value;
     914              index[numberNonZero++] = iColumn;
     915            }
    875916          }
    876      } else if (numberInRowArray == 2) {
    877           /* do by rows when two rows (do longer first when not packed
    878              and shorter first if packed */
    879           int iRow0 = whichRow[0];
    880           int iRow1 = whichRow[1];
    881           CoinBigIndex j;
    882           if (packed) {
    883                double pi0 = pi[0];
    884                double pi1 = pi[1];
    885                if (startPositive[iRow0+1] - startPositive[iRow0] >
    886                          startPositive[iRow1+1] - startPositive[iRow1]) {
    887                     int temp = iRow0;
    888                     iRow0 = iRow1;
    889                     iRow1 = temp;
    890                     pi0 = pi1;
    891                     pi1 = pi[0];
    892                }
    893                // and set up mark as char array
    894                char * COIN_RESTRICT marked = reinterpret_cast<char *> (index + columnArray->capacity());
    895                int * COIN_RESTRICT lookup = y->getIndices();
    896                double value = pi0 * scalar;
    897                int numberOriginal;
     917        }
     918        for (j = startNegative[iRow1]; j < startPositive[iRow1 + 1]; j++) {
     919          int iColumn = column[j];
     920          if (marked[iColumn]) {
     921            int iLookup = lookup[iColumn];
     922            array[iLookup] -= value;
     923          } else {
     924            if (fabs(value) > zeroTolerance) {
     925              array[numberNonZero] = -value;
     926              index[numberNonZero++] = iColumn;
     927            }
     928          }
     929        }
    898930#ifdef CLP_PLUS_ONE_MATRIX
    899                if ((otherFlags_&1)==0||!doPlusOnes) {
    900 #endif
    901                  for (j = startPositive[iRow0]; j < startNegative[iRow0]; j++) {
    902                     int iColumn = column[j];
    903                     array[numberNonZero] = value;
    904                     marked[iColumn] = 1;
    905                     lookup[iColumn] = numberNonZero;
    906                     index[numberNonZero++] = iColumn;
    907                  }
    908                  for (j = startNegative[iRow0]; j < startPositive[iRow0+1]; j++) {
    909                     int iColumn = column[j];
    910                     array[numberNonZero] = -value;
    911                     marked[iColumn] = 1;
    912                     lookup[iColumn] = numberNonZero;
    913                     index[numberNonZero++] = iColumn;
    914                  }
    915                  numberOriginal = numberNonZero;
    916                  value = pi1 * scalar;
    917                  for (j = startPositive[iRow1]; j < startNegative[iRow1]; j++) {
    918                     int iColumn = column[j];
    919                     if (marked[iColumn]) {
    920                          int iLookup = lookup[iColumn];
    921                          array[iLookup] += value;
    922                     } else {
    923                          if (fabs(value) > zeroTolerance) {
    924                               array[numberNonZero] = value;
    925                               index[numberNonZero++] = iColumn;
    926                          }
    927                     }
    928                  }
    929                  for (j = startNegative[iRow1]; j < startPositive[iRow1+1]; j++) {
    930                     int iColumn = column[j];
    931                     if (marked[iColumn]) {
    932                          int iLookup = lookup[iColumn];
    933                          array[iLookup] -= value;
    934                     } else {
    935                          if (fabs(value) > zeroTolerance) {
    936                               array[numberNonZero] = -value;
    937                               index[numberNonZero++] = iColumn;
    938                          }
    939                     }
    940                  }
     931      } else {
     932        // plus one
     933        oneit(7);
     934        for (j = startPositive[iRow0]; j < startPositive[iRow0 + 1]; j++) {
     935          int iColumn = column[j];
     936          array[numberNonZero] = value;
     937          marked[iColumn] = 1;
     938          lookup[iColumn] = numberNonZero;
     939          index[numberNonZero++] = iColumn;
     940        }
     941        numberOriginal = numberNonZero;
     942        value = pi1 * scalar;
     943        for (j = startPositive[iRow1]; j < startPositive[iRow1 + 1]; j++) {
     944          int iColumn = column[j];
     945          if (marked[iColumn]) {
     946            int iLookup = lookup[iColumn];
     947            array[iLookup] += value;
     948          } else {
     949            if (fabs(value) > zeroTolerance) {
     950              array[numberNonZero] = value;
     951              index[numberNonZero++] = iColumn;
     952            }
     953          }
     954        }
     955      }
     956#endif
     957      // get rid of tiny values and zero out marked
     958      int nDelete = 0;
     959      for (j = 0; j < numberOriginal; j++) {
     960        int iColumn = index[j];
     961        marked[iColumn] = 0;
     962        if (fabs(array[j]) <= zeroTolerance)
     963          nDelete++;
     964      }
     965      if (nDelete) {
     966        numberOriginal = numberNonZero;
     967        numberNonZero = 0;
     968        for (j = 0; j < numberOriginal; j++) {
     969          int iColumn = index[j];
     970          double value = array[j];
     971          array[j] = 0.0;
     972          if (fabs(value) > zeroTolerance) {
     973            array[numberNonZero] = value;
     974            index[numberNonZero++] = iColumn;
     975          }
     976        }
     977      }
     978    } else {
     979      if (startPositive[iRow0 + 1] - startPositive[iRow0] < startPositive[iRow1 + 1] - startPositive[iRow1]) {
     980        int temp = iRow0;
     981        iRow0 = iRow1;
     982        iRow1 = temp;
     983      }
     984      int numberOriginal;
     985      int i;
     986      numberNonZero = 0;
     987      double value;
     988      value = pi[iRow0] * scalar;
     989      CoinBigIndex j;
     990      for (j = startPositive[iRow0]; j < startNegative[iRow0]; j++) {
     991        int iColumn = column[j];
     992        index[numberNonZero++] = iColumn;
     993        array[iColumn] = value;
     994      }
     995      for (j = startNegative[iRow0]; j < startPositive[iRow0 + 1]; j++) {
     996        int iColumn = column[j];
     997        index[numberNonZero++] = iColumn;
     998        array[iColumn] = -value;
     999      }
     1000      value = pi[iRow1] * scalar;
     1001      for (j = startPositive[iRow1]; j < startNegative[iRow1]; j++) {
     1002        int iColumn = column[j];
     1003        double value2 = array[iColumn];
     1004        if (value2) {
     1005          value2 += value;
     1006        } else {
     1007          value2 = value;
     1008          index[numberNonZero++] = iColumn;
     1009        }
     1010        array[iColumn] = value2;
     1011      }
     1012      for (j = startNegative[iRow1]; j < startPositive[iRow1 + 1]; j++) {
     1013        int iColumn = column[j];
     1014        double value2 = array[iColumn];
     1015        if (value2) {
     1016          value2 -= value;
     1017        } else {
     1018          value2 = -value;
     1019          index[numberNonZero++] = iColumn;
     1020        }
     1021        array[iColumn] = value2;
     1022      }
     1023      // get rid of tiny values and zero out marked
     1024      numberOriginal = numberNonZero;
     1025      numberNonZero = 0;
     1026      for (i = 0; i < numberOriginal; i++) {
     1027        int iColumn = index[i];
     1028        if (fabs(array[iColumn]) > zeroTolerance) {
     1029          index[numberNonZero++] = iColumn;
     1030        } else {
     1031          array[iColumn] = 0.0;
     1032        }
     1033      }
     1034    }
     1035  } else if (numberInRowArray == 1) {
     1036    // Just one row
     1037    int iRow = rowArray->getIndices()[0];
     1038    numberNonZero = 0;
     1039    double value;
     1040    iRow = whichRow[0];
     1041    CoinBigIndex j;
     1042    if (packed) {
     1043      value = pi[0] * scalar;
     1044      if (fabs(value) > zeroTolerance) {
    9411045#ifdef CLP_PLUS_ONE_MATRIX
    942                } else {
    943                  // plus one
    944                  oneit(7);
    945                  for (j = startPositive[iRow0]; j < startPositive[iRow0+1]; j++) {
    946                     int iColumn = column[j];
    947                     array[numberNonZero] = value;
    948                     marked[iColumn] = 1;
    949                     lookup[iColumn] = numberNonZero;
    950                     index[numberNonZero++] = iColumn;
    951                  }
    952                  numberOriginal = numberNonZero;
    953                  value = pi1 * scalar;
    954                  for (j = startPositive[iRow1]; j < startPositive[iRow1+1]; j++) {
    955                     int iColumn = column[j];
    956                     if (marked[iColumn]) {
    957                          int iLookup = lookup[iColumn];
    958                          array[iLookup] += value;
    959                     } else {
    960                          if (fabs(value) > zeroTolerance) {
    961                               array[numberNonZero] = value;
    962                               index[numberNonZero++] = iColumn;
    963                          }
    964                     }
    965                  }
    966                }
    967 #endif
    968                // get rid of tiny values and zero out marked
    969                int nDelete = 0;
    970                for (j = 0; j < numberOriginal; j++) {
    971                     int iColumn = index[j];
    972                     marked[iColumn] = 0;
    973                     if (fabs(array[j]) <= zeroTolerance)
    974                          nDelete++;
    975                }
    976                if (nDelete) {
    977                     numberOriginal = numberNonZero;
    978                     numberNonZero = 0;
    979                     for (j = 0; j < numberOriginal; j++) {
    980                          int iColumn = index[j];
    981                          double value = array[j];
    982                          array[j] = 0.0;
    983                          if (fabs(value) > zeroTolerance) {
    984                               array[numberNonZero] = value;
    985                               index[numberNonZero++] = iColumn;
    986                          }
    987                     }
    988                }
    989           } else {
    990                if (startPositive[iRow0+1] - startPositive[iRow0] <
    991                          startPositive[iRow1+1] - startPositive[iRow1]) {
    992                     int temp = iRow0;
    993                     iRow0 = iRow1;
    994                     iRow1 = temp;
    995                }
    996                int numberOriginal;
    997                int i;
    998                numberNonZero = 0;
    999                double value;
    1000                value = pi[iRow0] * scalar;
    1001                CoinBigIndex j;
    1002                for (j = startPositive[iRow0]; j < startNegative[iRow0]; j++) {
    1003                  int iColumn = column[j];
    1004                  index[numberNonZero++] = iColumn;
    1005                  array[iColumn] = value;
    1006                }
    1007                for (j = startNegative[iRow0]; j < startPositive[iRow0+1]; j++) {
    1008                  int iColumn = column[j];
    1009                  index[numberNonZero++] = iColumn;
    1010                  array[iColumn] = -value;
    1011                }
    1012                value = pi[iRow1] * scalar;
    1013                for (j = startPositive[iRow1]; j < startNegative[iRow1]; j++) {
    1014                  int iColumn = column[j];
    1015                  double value2 = array[iColumn];
    1016                  if (value2) {
    1017                    value2 += value;
    1018                  } else {
    1019                    value2 = value;
    1020                    index[numberNonZero++] = iColumn;
    1021                  }
    1022                  array[iColumn] = value2;
    1023                }
    1024                for (j = startNegative[iRow1]; j < startPositive[iRow1+1]; j++) {
    1025                  int iColumn = column[j];
    1026                  double value2 = array[iColumn];
    1027                  if (value2) {
    1028                    value2 -= value;
    1029                  } else {
    1030                    value2 = -value;
    1031                    index[numberNonZero++] = iColumn;
    1032                  }
    1033                  array[iColumn] = value2;
    1034                }
    1035                // get rid of tiny values and zero out marked
    1036                numberOriginal = numberNonZero;
    1037                numberNonZero = 0;
    1038                for (i = 0; i < numberOriginal; i++) {
    1039                     int iColumn = index[i];
    1040                     if (fabs(array[iColumn]) > zeroTolerance) {
    1041                          index[numberNonZero++] = iColumn;
    1042                     } else {
    1043                          array[iColumn] = 0.0;
    1044                     }
    1045                }
     1046        if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     1047#endif
     1048          for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
     1049            int iColumn = column[j];
     1050            array[numberNonZero] = value;
     1051            index[numberNonZero++] = iColumn;
    10461052          }
    1047      } else if (numberInRowArray == 1) {
    1048           // Just one row
    1049           int iRow = rowArray->getIndices()[0];
    1050           numberNonZero = 0;
    1051           double value;
    1052           iRow = whichRow[0];
    1053           CoinBigIndex j;
    1054           if (packed) {
    1055                value = pi[0] * scalar;
    1056                if (fabs(value) > zeroTolerance) {
     1053          for (j = startNegative[iRow]; j < startPositive[iRow + 1]; j++) {
     1054            int iColumn = column[j];
     1055            array[numberNonZero] = -value;
     1056            index[numberNonZero++] = iColumn;
     1057          }
    10571058#ifdef CLP_PLUS_ONE_MATRIX
    1058                  if ((otherFlags_&1)==0||!doPlusOnes) {
    1059 #endif
    1060                     for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
    1061                          int iColumn = column[j];
    1062                          array[numberNonZero] = value;
    1063                          index[numberNonZero++] = iColumn;
    1064                     }
    1065                     for (j = startNegative[iRow]; j < startPositive[iRow+1]; j++) {
    1066                          int iColumn = column[j];
    1067                          array[numberNonZero] = -value;
    1068                          index[numberNonZero++] = iColumn;
    1069                     }
    1070 #ifdef CLP_PLUS_ONE_MATRIX
    1071                  } else {
    1072                    // plus one
    1073                    oneit(9);
    1074                     for (j = startPositive[iRow]; j < startPositive[iRow+1]; j++) {
    1075                          int iColumn = column[j];
    1076                          array[numberNonZero] = value;
    1077                          index[numberNonZero++] = iColumn;
    1078                     }
    1079                  }
    1080 #endif
    1081                }
    1082           } else {
    1083                value = pi[iRow] * scalar;
    1084                if (fabs(value) > zeroTolerance) {
    1085                   for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
    1086                          int iColumn = column[j];
    1087                          array[iColumn] = value;
    1088                          index[numberNonZero++] = iColumn;
    1089                     }
    1090                     for (j = startNegative[iRow]; j < startPositive[iRow+1]; j++) {
    1091                          int iColumn = column[j];
    1092                          array[iColumn] = -value;
    1093                          index[numberNonZero++] = iColumn;
    1094                     }
    1095                }
     1059        } else {
     1060          // plus one
     1061          oneit(9);
     1062          for (j = startPositive[iRow]; j < startPositive[iRow + 1]; j++) {
     1063            int iColumn = column[j];
     1064            array[numberNonZero] = value;
     1065            index[numberNonZero++] = iColumn;
    10961066          }
    1097      }
    1098      columnArray->setNumElements(numberNonZero);
    1099      if (packed)
    1100           columnArray->setPacked();
    1101      y->setNumElements(0);
     1067        }
     1068#endif
     1069      }
     1070    } else {
     1071      value = pi[iRow] * scalar;
     1072      if (fabs(value) > zeroTolerance) {
     1073        for (j = startPositive[iRow]; j < startNegative[iRow]; j++) {
     1074          int iColumn = column[j];
     1075          array[iColumn] = value;
     1076          index[numberNonZero++] = iColumn;
     1077        }
     1078        for (j = startNegative[iRow]; j < startPositive[iRow + 1]; j++) {
     1079          int iColumn = column[j];
     1080          array[iColumn] = -value;
     1081          index[numberNonZero++] = iColumn;
     1082        }
     1083      }
     1084    }
     1085  }
     1086  columnArray->setNumElements(numberNonZero);
     1087  if (packed)
     1088    columnArray->setPacked();
     1089  y->setNumElements(0);
    11021090}
    11031091/* Return <code>x *A in <code>z</code> but
    11041092   just for indices in y. */
    1105 void
    1106 ClpPlusMinusOneMatrix::subsetTransposeTimes(const ClpSimplex * ,
    1107           const CoinIndexedVector * rowArray,
    1108           const CoinIndexedVector * y,
    1109           CoinIndexedVector * columnArray) const
    1110 {
    1111      columnArray->clear();
    1112      double * COIN_RESTRICT pi = rowArray->denseVector();
    1113      double * COIN_RESTRICT array = columnArray->denseVector();
    1114      int jColumn;
    1115      int numberToDo = y->getNumElements();
    1116      const int * COIN_RESTRICT which = y->getIndices();
    1117      assert (!rowArray->packedMode());
    1118      columnArray->setPacked();
     1093void ClpPlusMinusOneMatrix::subsetTransposeTimes(const ClpSimplex *,
     1094  const CoinIndexedVector *rowArray,
     1095  const CoinIndexedVector *y,
     1096  CoinIndexedVector *columnArray) const
     1097{
     1098  columnArray->clear();
     1099  double *COIN_RESTRICT pi = rowArray->denseVector();
     1100  double *COIN_RESTRICT array = columnArray->denseVector();
     1101  int jColumn;
     1102  int numberToDo = y->getNumElements();
     1103  const int *COIN_RESTRICT which = y->getIndices();
     1104  assert(!rowArray->packedMode());
     1105  columnArray->setPacked();
    11191106#ifdef CLP_PLUS_ONE_MATRIX
    1120      if ((otherFlags_&1)==0||!doPlusOnes) {
    1121 #endif
    1122        for (jColumn = 0; jColumn < numberToDo; jColumn++) {
    1123           int iColumn = which[jColumn];
    1124           double value = 0.0;
    1125           CoinBigIndex j = startPositive_[iColumn];
    1126           for (; j < startNegative_[iColumn]; j++) {
    1127                int iRow = indices_[j];
    1128                value += pi[iRow];
    1129           }
    1130           for (; j < startPositive_[iColumn+1]; j++) {
    1131                int iRow = indices_[j];
    1132                value -= pi[iRow];
    1133           }
    1134           array[jColumn] = value;
    1135        }
     1107  if ((otherFlags_ & 1) == 0 || !doPlusOnes) {
     1108#endif
     1109    for (jColumn = 0; jColumn < numberToDo; jColumn++) {
     1110      int iColumn = which[jColumn];
     1111      double value = 0.0;
     1112      CoinBigIndex j = startPositive_[iColumn];
     1113      for (; j < startNegative_[iColumn]; j++) {
     1114        int iRow = indices_[j];
     1115        value += pi[iRow];
     1116      }
     1117      for (; j < startPositive_[iColumn + 1]; j++) {
     1118        int iRow = indices_[j];
     1119        value -= pi[iRow];
     1120      }
     1121      array[jColumn] = value;
     1122    }
    11361123#ifdef CLP_PLUS_ONE_MATRIX
    1137      } else {
    1138        // plus one
    1139        oneit(11);
    1140        for (jColumn = 0; jColumn < numberToDo; jColumn++) {
    1141           int iColumn = which[jColumn];
    1142           double value = 0.0;
    1143           CoinBigIndex j = startPositive_[iColumn];
    1144           for (; j < startPositive_[iColumn+1]; j++) {
    1145                int iRow = indices_[j];
    1146                value += pi[iRow];
    1147           }
    1148           array[jColumn] = value;
    1149        }
    1150      }
     1124  } else {
     1125    // plus one
     1126    oneit(11);
     1127    for (jColumn = 0; jColumn < numberToDo; jColumn++) {
     1128      int iColumn = which[jColumn];
     1129      double value = 0.0;
     1130      CoinBigIndex j = startPositive_[iColumn];
     1131      for (; j < startPositive_[iColumn + 1]; j++) {
     1132        int iRow = indices_[j];
     1133        value += pi[iRow];
     1134      }
     1135      array[jColumn] = value;
     1136    }
     1137  }
    11511138#endif
    11521139}
    11531140/// returns number of elements in column part of basis,
    1154 int
    1155 ClpPlusMinusOneMatrix::countBasis(const int * whichColumn,
    1156                                   int & numberColumnBasic)
    1157 {
    1158      int i;
    1159      CoinBigIndex numberElements = 0;
    1160      for (i = 0; i < numberColumnBasic; i++) {
    1161           int iColumn = whichColumn[i];
    1162           numberElements += startPositive_[iColumn+1] - startPositive_[iColumn];
    1163      }
     1141int ClpPlusMinusOneMatrix::countBasis(const int *whichColumn,
     1142  int &numberColumnBasic)
     1143{
     1144  int i;
     1145  CoinBigIndex numberElements = 0;
     1146  for (i = 0; i < numberColumnBasic; i++) {
     1147    int iColumn = whichColumn[i];
     1148    numberElements += startPositive_[iColumn + 1] - startPositive_[iColumn];
     1149  }
    11641150#if COIN_BIG_INDEX
    1165      if (numberElements>COIN_INT_MAX) {
    1166        printf("Factorization too large\n");
    1167        abort();
    1168      }
    1169 #endif
    1170      return static_cast<int>(numberElements);
    1171 }
    1172 void
    1173 ClpPlusMinusOneMatrix::fillBasis(ClpSimplex * ,
    1174                                  const int * whichColumn,
    1175                                  int & numberColumnBasic,
    1176                                  int * indexRowU, int * start,
    1177                                  int * rowCount, int * columnCount,
    1178                                  CoinFactorizationDouble * elementU)
    1179 {
    1180      int i;
    1181      CoinBigIndex numberElements = start[0];
    1182      assert (columnOrdered_);
    1183      for (i = 0; i < numberColumnBasic; i++) {
    1184           int iColumn = whichColumn[i];
    1185           CoinBigIndex j = startPositive_[iColumn];
    1186           for (; j < startNegative_[iColumn]; j++) {
    1187                int iRow = indices_[j];
    1188                indexRowU[numberElements] = iRow;
    1189                rowCount[iRow]++;
    1190                elementU[numberElements++] = 1.0;
    1191           }
    1192           for (; j < startPositive_[iColumn+1]; j++) {
    1193                int iRow = indices_[j];
    1194                indexRowU[numberElements] = iRow;
    1195                rowCount[iRow]++;
    1196                elementU[numberElements++] = -1.0;
    1197           }
    1198           start[i+1] = static_cast<int>(numberElements);
    1199           columnCount[i] = static_cast<int>(numberElements - start[i]);
    1200      }
     1151  if (numberElements > COIN_INT_MAX) {
     1152    printf("Factorization too large\n");
     1153    abort();
     1154  }
     1155#endif
     1156  return static_cast< int >(numberElements);
     1157}
     1158void ClpPlusMinusOneMatrix::fillBasis(ClpSimplex *,
     1159  const int *whichColumn,
     1160  int &numberColumnBasic,
     1161  int *indexRowU, int *start,
     1162  int *rowCount, int *columnCount,
     1163  CoinFactorizationDouble *elementU)
     1164{
     1165  int i;
     1166  CoinBigIndex numberElements = start[0];
     1167  assert(columnOrdered_);
     1168  for (i = 0; i < numberColumnBasic; i++) {
     1169    int iColumn = whichColumn[i];
     1170    CoinBigIndex j = startPositive_[iColumn];
     1171    for (; j < startNegative_[iColumn]; j++) {
     1172      int iRow = indices_[j];
     1173      indexRowU[numberElements] = iRow;
     1174      rowCount[iRow]++;
     1175      elementU[numberElements++] = 1.0;
     1176    }
     1177    for (; j < startPositive_[iColumn + 1]; j++) {
     1178      int iRow = indices_[j];
     1179      indexRowU[numberElements] = iRow;
     1180      rowCount[iRow]++;
     1181      elementU[numberElements++] = -1.0;
     1182    }
     1183    start[i + 1] = static_cast< int >(numberElements);
     1184    columnCount[i] = static_cast< int >(numberElements - start[i]);
     1185  }
    12011186#if COIN_BIG_INDEX
    1202      if (numberElements>COIN_INT_MAX)
    1203        abort();
     1187  if (numberElements > COIN_INT_MAX)
     1188    abort();
    12041189#endif
    12051190}
    12061191/* Unpacks a column into an CoinIndexedvector
    12071192 */
    1208 void
    1209 ClpPlusMinusOneMatrix::unpack(const ClpSimplex * ,
    1210                               CoinIndexedVector * rowArray,
    1211                               int iColumn) const
    1212 {
    1213      CoinBigIndex j = startPositive_[iColumn];
    1214      for (; j < startNegative_[iColumn]; j++) {
    1215           int iRow = indices_[j];
    1216           rowArray->add(iRow, 1.0);
    1217      }
    1218      for (; j < startPositive_[iColumn+1]; j++) {
    1219           int iRow = indices_[j];
    1220           rowArray->add(iRow, -1.0);
    1221      }
     1193void ClpPlusMinusOneMatrix::unpack(const ClpSimplex *,
     1194  CoinIndexedVector *rowArray,
     1195  int iColumn) const
     1196{
     1197  CoinBigIndex j = startPositive_[iColumn];
     1198  for (; j < startNegative_[iColumn]; j++) {
     1199    int iRow = indices_[j];
     1200    rowArray->add(iRow, 1.0);
     1201  }
     1202  for (; j < startPositive_[iColumn + 1]; j++) {
     1203    int iRow = indices_[j];
     1204    rowArray->add(iRow, -1.0);
     1205  }
    12221206}
    12231207/* Unpacks a column into an CoinIndexedvector
     
    12251209Note that model is NOT const.  Bounds and objective could
    12261210be modified if doing column generation (just for this variable) */
    1227 void
    1228 ClpPlusMinusOneMatrix::unpackPacked(ClpSimplex * ,
    1229                                     CoinIndexedVector * rowArray,
    1230                                     int iColumn) const
    1231 {
    1232      int * COIN_RESTRICT index = rowArray->getIndices();
    1233      double * COIN_RESTRICT array = rowArray->denseVector();
    1234      int number = 0;
    1235      CoinBigIndex j = startPositive_[iColumn];
    1236      for (; j < startNegative_[iColumn]; j++) {
    1237           int iRow = indices_[j];
    1238           array[number] = 1.0;
    1239           index[number++] = iRow;
    1240      }
    1241      for (; j < startPositive_[iColumn+1]; j++) {
    1242           int iRow = indices_[j];
    1243           array[number] = -1.0;
    1244           index[number++] = iRow;
    1245      }
    1246      rowArray->setNumElements(number);
    1247      rowArray->setPackedMode(true);
     1211void ClpPlusMinusOneMatrix::unpackPacked(ClpSimplex *,
     1212  CoinIndexedVector *rowArray,
     1213  int iColumn) const
     1214{
     1215  int *COIN_RESTRICT index = rowArray->getIndices();
     1216  double *COIN_RESTRICT array = rowArray->denseVector();
     1217  int number = 0;
     1218  CoinBigIndex j = startPositive_[iColumn];
     1219  for (; j < startNegative_[iColumn]; j++) {
     1220    int iRow = indices_[j];
     1221    array[number] = 1.0;
     1222    index[number++] = iRow;
     1223  }
     1224  for (; j < startPositive_[iColumn + 1]; j++) {
     1225    int iRow = indices_[j];
     1226    array[number] = -1.0;
     1227    index[number++] = iRow;
     1228  }
     1229  rowArray->setNumElements(number);
     1230  rowArray->setPackedMode(true);
    12481231}
    12491232/* Adds multiple of a column into an CoinIndexedvector
    12501233      You can use quickAdd to add to vector */
    1251 void
    1252 ClpPlusMinusOneMatrix::add(const ClpSimplex * , CoinIndexedVector * rowArray,
    1253                            int iColumn, double multiplier) const
    1254 {
    1255      CoinBigIndex j = startPositive_[iColumn];
    1256      for (; j < startNegative_[iColumn]; j++) {
    1257           int iRow = indices_[j];
    1258           rowArray->quickAdd(iRow, multiplier);
    1259      }
    1260      for (; j < startPositive_[iColumn+1]; j++) {
    1261           int iRow = indices_[j];
    1262           rowArray->quickAdd(iRow, -multiplier);
    1263      }
     1234void ClpPlusMinusOneMatrix::add(const ClpSimplex *, CoinIndexedVector *rowArray,
     1235  int iColumn, double multiplier) const
     1236{
     1237  CoinBigIndex j = startPositive_[iColumn];
     1238  for (; j < startNegative_[iColumn]; j++) {
     1239    int iRow = indices_[j];
     1240    rowArray->quickAdd(iRow, multiplier);
     1241  }
     1242  for (; j < startPositive_[iColumn + 1]; j++) {
     1243    int iRow = indices_[j];
     1244    rowArray->quickAdd(iRow, -multiplier);
     1245  }
    12641246}
    12651247/* Adds multiple of a column into an array */
    1266 void
    1267 ClpPlusMinusOneMatrix::add(const ClpSimplex * , double * array,
    1268                            int iColumn, double multiplier) const
    1269 {
    1270      CoinBigIndex j = startPositive_[iColumn];
    1271      for (; j < startNegative_[iColumn]; j++) {
    1272           int iRow = indices_[j];
    1273           array[iRow] += multiplier;
    1274      }
    1275      for (; j < startPositive_[iColumn+1]; j++) {
    1276           int iRow = indices_[j];
    1277           array[iRow] -= multiplier;
    1278      }
     1248void ClpPlusMinusOneMatrix::add(const ClpSimplex *, double *array,
     1249  int iColumn, double multiplier) const
     1250{
     1251  CoinBigIndex j = startPositive_[iColumn];
     1252  for (; j < startNegative_[iColumn]; j++) {
     1253    int iRow = indices_[j];
     1254    array[iRow] += multiplier;
     1255  }
     1256  for (; j < startPositive_[iColumn + 1]; j++) {
     1257    int iRow = indices_[j];
     1258    array[iRow] -= multiplier;
     1259  }
    12791260}
    12801261
     
    12831264ClpPlusMinusOneMatrix::getPackedMatrix() const
    12841265{
    1285      if (!matrix_) {
    1286           int numberMinor = (!columnOrdered_) ? numberColumns_ : numberRows_;
    1287           int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    1288           CoinBigIndex numberElements = startPositive_[numberMajor];
    1289           double * elements = new double [numberElements];
    1290           CoinBigIndex j = 0;
    1291           int i;
    1292           for (i = 0; i < numberMajor; i++) {
    1293                for (; j < startNegative_[i]; j++) {
    1294                     elements[j] = 1.0;
    1295                }
    1296                for (; j < startPositive_[i+1]; j++) {
    1297                     elements[j] = -1.0;
    1298                }
    1299           }
    1300           matrix_ = new CoinPackedMatrix(columnOrdered_, numberMinor, numberMajor,
    1301                                           getNumElements(),
    1302                                           elements, indices_,
    1303                                           startPositive_, getVectorLengths());
    1304           delete [] elements;
    1305           delete [] lengths_;
    1306           lengths_ = NULL;
    1307      }
    1308      return matrix_;
     1266  if (!matrix_) {
     1267    int numberMinor = (!columnOrdered_) ? numberColumns_ : numberRows_;
     1268    int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     1269    CoinBigIndex numberElements = startPositive_[numberMajor];
     1270    double *elements = new double[numberElements];
     1271    CoinBigIndex j = 0;
     1272    int i;
     1273    for (i = 0; i < numberMajor; i++) {
     1274      for (; j < startNegative_[i]; j++) {
     1275        elements[j] = 1.0;
     1276      }
     1277      for (; j < startPositive_[i + 1]; j++) {
     1278        elements[j] = -1.0;
     1279      }
     1280    }
     1281    matrix_ = new CoinPackedMatrix(columnOrdered_, numberMinor, numberMajor,
     1282      getNumElements(),
     1283      elements, indices_,
     1284      startPositive_, getVectorLengths());
     1285    delete[] elements;
     1286    delete[] lengths_;
     1287    lengths_ = NULL;
     1288  }
     1289  return matrix_;
    13091290}
    13101291/* A vector containing the elements in the packed matrix. Note that there
     
    13151296ClpPlusMinusOneMatrix::getElements() const
    13161297{
    1317      if (!matrix_)
    1318           getPackedMatrix();
    1319      return matrix_->getElements();
     1298  if (!matrix_)
     1299    getPackedMatrix();
     1300  return matrix_->getElements();
    13201301}
    13211302
     
    13231304ClpPlusMinusOneMatrix::getVectorStarts() const
    13241305{
    1325      return startPositive_;
     1306  return startPositive_;
    13261307}
    13271308/* The lengths of the major-dimension vectors. */
     
    13291310ClpPlusMinusOneMatrix::getVectorLengths() const
    13301311{
    1331      if (!lengths_) {
    1332           int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    1333           lengths_ = new int [numberMajor];
    1334           int i;
    1335           for (i = 0; i < numberMajor; i++) {
    1336             lengths_[i] = static_cast<int>(startPositive_[i+1] - startPositive_[i]);
    1337           }
    1338      }
    1339      return lengths_;
     1312  if (!lengths_) {
     1313    int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     1314    lengths_ = new int[numberMajor];
     1315    int i;
     1316    for (i = 0; i < numberMajor; i++) {
     1317      lengths_[i] = static_cast< int >(startPositive_[i + 1] - startPositive_[i]);
     1318    }
     1319  }
     1320  return lengths_;
    13401321}
    13411322/* Delete the columns whose indices are listed in <code>indDel</code>. */
    1342 void
    1343 ClpPlusMinusOneMatrix::deleteCols(const int numDel, const int * indDel)
    1344 {
    1345      int iColumn;
    1346      CoinBigIndex newSize = startPositive_[numberColumns_];;
    1347      int numberBad = 0;
    1348      // Use array to make sure we can have duplicates
    1349      int * which = new int[numberColumns_];
    1350      memset(which, 0, numberColumns_ * sizeof(int));
    1351      int nDuplicate = 0;
    1352      for (iColumn = 0; iColumn < numDel; iColumn++) {
    1353           int jColumn = indDel[iColumn];
    1354           if (jColumn < 0 || jColumn >= numberColumns_) {
    1355                numberBad++;
    1356           } else {
    1357                newSize -= startPositive_[jColumn+1] - startPositive_[jColumn];
    1358                if (which[jColumn])
    1359                     nDuplicate++;
    1360                else
    1361                     which[jColumn] = 1;
    1362           }
    1363      }
    1364      if (numberBad)
    1365           throw CoinError("Indices out of range", "deleteCols", "ClpPlusMinusOneMatrix");
    1366      int newNumber = numberColumns_ - numDel + nDuplicate;
    1367      // Get rid of temporary arrays
    1368      delete [] lengths_;
    1369      lengths_ = NULL;
    1370      delete matrix_;
    1371      matrix_ = NULL;
    1372      CoinBigIndex * newPositive = new CoinBigIndex [newNumber+1];
    1373      CoinBigIndex * newNegative = new CoinBigIndex [newNumber];
    1374      int * newIndices = new int [newSize];
    1375      newNumber = 0;
    1376      newSize = 0;
    1377      for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    1378           if (!which[iColumn]) {
    1379                CoinBigIndex start, end;
    1380                CoinBigIndex i;
    1381                start = startPositive_[iColumn];
    1382                end = startNegative_[iColumn];
    1383                newPositive[newNumber] = newSize;
    1384                for (i = start; i < end; i++)
    1385                     newIndices[newSize++] = indices_[i];
    1386                start = startNegative_[iColumn];
    1387                end = startPositive_[iColumn+1];
    1388                newNegative[newNumber++] = newSize;
    1389                for (i = start; i < end; i++)
    1390                     newIndices[newSize++] = indices_[i];
    1391           }
    1392      }
    1393      newPositive[newNumber] = newSize;
    1394      delete [] which;
    1395      delete [] startPositive_;
    1396      startPositive_ = newPositive;
    1397      delete [] startNegative_;
    1398      startNegative_ = newNegative;
    1399      delete [] indices_;
    1400      indices_ = newIndices;
    1401      numberColumns_ = newNumber;
     1323void ClpPlusMinusOneMatrix::deleteCols(const int numDel, const int *indDel)
     1324{
     1325  int iColumn;
     1326  CoinBigIndex newSize = startPositive_[numberColumns_];
     1327  ;
     1328  int numberBad = 0;
     1329  // Use array to make sure we can have duplicates
     1330  int *which = new int[numberColumns_];
     1331  memset(which, 0, numberColumns_ * sizeof(int));
     1332  int nDuplicate = 0;
     1333  for (iColumn = 0; iColumn < numDel; iColumn++) {
     1334    int jColumn = indDel[iColumn];
     1335    if (jColumn < 0 || jColumn >= numberColumns_) {
     1336      numberBad++;
     1337    } else {
     1338      newSize -= startPositive_[jColumn + 1] - startPositive_[jColumn];
     1339      if (which[jColumn])
     1340        nDuplicate++;
     1341      else
     1342        which[jColumn] = 1;
     1343    }
     1344  }
     1345  if (numberBad)
     1346    throw CoinError("Indices out of range", "deleteCols", "ClpPlusMinusOneMatrix");
     1347  int newNumber = numberColumns_ - numDel + nDuplicate;
     1348  // Get rid of temporary arrays
     1349  delete[] lengths_;
     1350  lengths_ = NULL;
     1351  delete matrix_;
     1352  matrix_ = NULL;
     1353  CoinBigIndex *newPositive = new CoinBigIndex[newNumber + 1];
     1354  CoinBigIndex *newNegative = new CoinBigIndex[newNumber];
     1355  int *newIndices = new int[newSize];
     1356  newNumber = 0;
     1357  newSize = 0;
     1358  for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     1359    if (!which[iColumn]) {
     1360      CoinBigIndex start, end;
     1361      CoinBigIndex i;
     1362      start = startPositive_[iColumn];
     1363      end = startNegative_[iColumn];
     1364      newPositive[newNumber] = newSize;
     1365      for (i = start; i < end; i++)
     1366        newIndices[newSize++] = indices_[i];
     1367      start = startNegative_[iColumn];
     1368      end = startPositive_[iColumn + 1];
     1369      newNegative[newNumber++] = newSize;
     1370      for (i = start; i < end; i++)
     1371        newIndices[newSize++] = indices_[i];
     1372    }
     1373  }
     1374  newPositive[newNumber] = newSize;
     1375  delete[] which;
     1376  delete[] startPositive_;
     1377  startPositive_ = newPositive;
     1378  delete[] startNegative_;
     1379  startNegative_ = newNegative;
     1380  delete[] indices_;
     1381  indices_ = newIndices;
     1382  numberColumns_ = newNumber;
    14021383}
    14031384/* Delete the rows whose indices are listed in <code>indDel</code>. */
    1404 void
    1405 ClpPlusMinusOneMatrix::deleteRows(const int numDel, const int * indDel)
    1406 {
    1407      int iRow;
    1408      int numberBad = 0;
    1409      // Use array to make sure we can have duplicates
    1410      int * which = new int[numberRows_];
    1411      memset(which, 0, numberRows_ * sizeof(int));
    1412      int nDuplicate = 0;
    1413      for (iRow = 0; iRow < numDel; iRow++) {
    1414           int jRow = indDel[iRow];
    1415           if (jRow < 0 || jRow >= numberRows_) {
    1416                numberBad++;
    1417           } else {
    1418                if (which[jRow])
    1419                     nDuplicate++;
    1420                else
    1421                     which[jRow] = 1;
    1422           }
    1423      }
    1424      if (numberBad)
    1425           throw CoinError("Indices out of range", "deleteRows", "ClpPlusMinusOneMatrix");
    1426      CoinBigIndex iElement;
    1427      CoinBigIndex numberElements = startPositive_[numberColumns_];
    1428      CoinBigIndex newSize = 0;
    1429      for (iElement = 0; iElement < numberElements; iElement++) {
    1430           iRow = indices_[iElement];
    1431           if (!which[iRow])
    1432                newSize++;
    1433      }
    1434      int newNumber = numberRows_ - numDel + nDuplicate;
    1435      // Get rid of temporary arrays
    1436      delete [] lengths_;
    1437      lengths_ = NULL;
    1438      delete matrix_;
    1439      matrix_ = NULL;
    1440      // redo which
    1441      int numberRows=0;
    1442      for (iRow = 0; iRow < numberRows_; iRow++) {
    1443        if (which[iRow]) {
    1444          which[iRow]=-1; // not wanted
    1445        } else {
    1446          which[iRow]=numberRows;
    1447          numberRows++;
    1448        }
    1449      }
    1450      int * newIndices = new int [newSize];
    1451      newSize = 0;
    1452      int iColumn;
    1453      for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    1454           CoinBigIndex start, end;
    1455           CoinBigIndex i;
    1456           start = startPositive_[iColumn];
    1457           end = startNegative_[iColumn];
    1458           startPositive_[newNumber] = newSize;
    1459           for (i = start; i < end; i++) {
    1460                iRow = which[indices_[i]];
    1461                if (iRow>=0)
    1462                     newIndices[newSize++] = iRow;
    1463           }
    1464           start = startNegative_[iColumn];
    1465           end = startPositive_[iColumn+1];
    1466           startNegative_[newNumber] = newSize;
    1467           for (i = start; i < end; i++) {
    1468                iRow = which[indices_[i]];
    1469                if (iRow>=0)
    1470                     newIndices[newSize++] = iRow;
    1471           }
    1472      }
    1473      startPositive_[numberColumns_] = newSize;
    1474      delete [] which;
    1475      delete [] indices_;
    1476      indices_ = newIndices;
    1477      numberRows_ = newNumber;
    1478 }
    1479 bool
    1480 ClpPlusMinusOneMatrix::isColOrdered() const
    1481 {
    1482      return columnOrdered_;
     1385void ClpPlusMinusOneMatrix::deleteRows(const int numDel, const int *indDel)
     1386{
     1387  int iRow;
     1388  int numberBad = 0;
     1389  // Use array to make sure we can have duplicates
     1390  int *which = new int[numberRows_];
     1391  memset(which, 0, numberRows_ * sizeof(int));
     1392  int nDuplicate = 0;
     1393  for (iRow = 0; iRow < numDel; iRow++) {
     1394    int jRow = indDel[iRow];
     1395    if (jRow < 0 || jRow >= numberRows_) {
     1396      numberBad++;
     1397    } else {
     1398      if (which[jRow])
     1399        nDuplicate++;
     1400      else
     1401        which[jRow] = 1;
     1402    }
     1403  }
     1404  if (numberBad)
     1405    throw CoinError("Indices out of range", "deleteRows", "ClpPlusMinusOneMatrix");
     1406  CoinBigIndex iElement;
     1407  CoinBigIndex numberElements = startPositive_[numberColumns_];
     1408  CoinBigIndex newSize = 0;
     1409  for (iElement = 0; iElement < numberElements; iElement++) {
     1410    iRow = indices_[iElement];
     1411    if (!which[iRow])
     1412      newSize++;
     1413  }
     1414  int newNumber = numberRows_ - numDel + nDuplicate;
     1415  // Get rid of temporary arrays
     1416  delete[] lengths_;
     1417  lengths_ = NULL;
     1418  delete matrix_;
     1419  matrix_ = NULL;
     1420  // redo which
     1421  int numberRows = 0;
     1422  for (iRow = 0; iRow < numberRows_; iRow++) {
     1423    if (which[iRow]) {
     1424      which[iRow] = -1; // not wanted
     1425    } else {
     1426      which[iRow] = numberRows;
     1427      numberRows++;
     1428    }
     1429  }
     1430  int *newIndices = new int[newSize];
     1431  newSize = 0;
     1432  int iColumn;
     1433  for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     1434    CoinBigIndex start, end;
     1435    CoinBigIndex i;
     1436    start = startPositive_[iColumn];
     1437    end = startNegative_[iColumn];
     1438    startPositive_[newNumber] = newSize;
     1439    for (i = start; i < end; i++) {
     1440      iRow = which[indices_[i]];
     1441      if (iRow >= 0)
     1442        newIndices[newSize++] = iRow;
     1443    }
     1444    start = startNegative_[iColumn];
     1445    end = startPositive_[iColumn + 1];
     1446    startNegative_[newNumber] = newSize;
     1447    for (i = start; i < end; i++) {
     1448      iRow = which[indices_[i]];
     1449      if (iRow >= 0)
     1450        newIndices[newSize++] = iRow;
     1451    }
     1452  }
     1453  startPositive_[numberColumns_] = newSize;
     1454  delete[] which;
     1455  delete[] indices_;
     1456  indices_ = newIndices;
     1457  numberRows_ = newNumber;
     1458}
     1459bool ClpPlusMinusOneMatrix::isColOrdered() const
     1460{
     1461  return columnOrdered_;
    14831462}
    14841463/* Number of entries in the packed matrix. */
     
    14861465ClpPlusMinusOneMatrix::getNumElements() const
    14871466{
    1488      int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
    1489      if (startPositive_)
    1490           return startPositive_[numberMajor];
    1491      else
    1492           return 0;
     1467  int numberMajor = (columnOrdered_) ? numberColumns_ : numberRows_;
     1468  if (startPositive_)
     1469    return startPositive_[numberMajor];
     1470  else
     1471    return 0;
    14931472}
    14941473// pass in copy (object takes ownership)
    1495 void
    1496 ClpPlusMinusOneMatrix::passInCopy(int numberRows, int numberColumns,
    1497                                   bool columnOrdered, int * indices,
    1498                                   CoinBigIndex * startPositive, CoinBigIndex * startNegative)
    1499 {
    1500      columnOrdered_ = columnOrdered;
     1474void ClpPlusMinusOneMatrix::passInCopy(int numberRows, int numberColumns,
     1475  bool columnOrdered, int *indices,
     1476  CoinBigIndex *startPositive, CoinBigIndex *startNegative)
     1477{
     1478  columnOrdered_ = columnOrdered;
    15011479#ifdef CLP_PLUS_ONE_MATRIX
    1502      otherFlags_ = 0;
    1503 #endif
    1504      startPositive_ = startPositive;
    1505      startNegative_ = startNegative;
    1506      indices_ = indices;
    1507      numberRows_ = numberRows;
    1508      numberColumns_ = numberColumns;
    1509      // Check valid
    1510      checkValid(false);
     1480  otherFlags_ = 0;
     1481#endif
     1482  startPositive_ = startPositive;
     1483  startNegative_ = startNegative;
     1484  indices_ = indices;
     1485  numberRows_ = numberRows;
     1486  numberColumns_ = numberColumns;
     1487  // Check valid
     1488  checkValid(false);
    15111489}
    15121490// Just checks matrix valid - will say if dimensions not quite right if detail
    1513 void
    1514 ClpPlusMinusOneMatrix::checkValid(bool detail) const
    1515 {
    1516      int maxIndex = -1;
    1517      int minIndex = columnOrdered_ ? numberRows_ : numberColumns_;
    1518      int number = !columnOrdered_ ? numberRows_ : numberColumns_;
    1519      CoinBigIndex numberElements = getNumElements();
    1520      CoinBigIndex last = -1;
    1521      int bad = 0;
    1522      // say if all +1
     1491void ClpPlusMinusOneMatrix::checkValid(bool detail) const
     1492{
     1493  int maxIndex = -1;
     1494  int minIndex = columnOrdered_ ? numberRows_ : numberColumns_;
     1495  int number = !columnOrdered_ ? numberRows_ : numberColumns_;
     1496  CoinBigIndex numberElements = getNumElements();
     1497  CoinBigIndex last = -1;
     1498  int bad = 0;
     1499  // say if all +1
    15231500#ifdef CLP_PLUS_ONE_MATRIX
    1524      otherFlags_ = 1;
    1525 #endif
    1526      for (int i = 0; i < number; i++) {
    1527           if(startPositive_[i] < last)
    1528                bad++;
    1529           else
    1530                last = startPositive_[i];
     1501  otherFlags_ = 1;
     1502#endif
     1503  for (int i = 0; i < number; i++) {
     1504    if (startPositive_[i] < last)
     1505      bad++;
     1506    else
     1507      last = startPositive_[i];
    15311508#ifdef CLP_PLUS_ONE_MATRIX
    1532           if (startNegative_[i] < startPositive_[i+1])
    1533             otherFlags_ &= ~1; // not all +1
    1534 #endif
    1535           if(startNegative_[i] < last)
    1536                bad++;
    1537           else
    1538                last = startNegative_[i];
    1539      }
    1540      if(startPositive_[number] < last)
    1541           bad++;
    1542      CoinAssertHint(!bad, "starts are not monotonic");
    1543      for (CoinBigIndex cbi = 0; cbi < numberElements; cbi++) {
    1544           maxIndex = CoinMax(indices_[cbi], maxIndex);
    1545           minIndex = CoinMin(indices_[cbi], minIndex);
    1546      }
    1547      CoinAssert(maxIndex < (columnOrdered_ ? numberRows_ : numberColumns_));
    1548      CoinAssert(minIndex >= 0);
    1549      if (detail) {
    1550           if (minIndex > 0 || maxIndex + 1 < (columnOrdered_ ? numberRows_ : numberColumns_))
    1551                printf("Not full range of indices - %d to %d\n", minIndex, maxIndex);
    1552      }
     1509    if (startNegative_[i] < startPositive_[i + 1])
     1510      otherFlags_ &= ~1; // not all +1
     1511#endif
     1512    if (startNegative_[i] < last)
     1513      bad++;
     1514    else
     1515      last = startNegative_[i];
     1516  }
     1517  if (startPositive_[number] < last)
     1518    bad++;
     1519  CoinAssertHint(!bad, "starts are not monotonic");
     1520  for (CoinBigIndex cbi = 0; cbi < numberElements; cbi++) {
     1521    maxIndex = CoinMax(indices_[cbi], maxIndex);
     1522    minIndex = CoinMin(indices_[cbi], minIndex);
     1523  }
     1524  CoinAssert(maxIndex < (columnOrdered_ ? numberRows_ : numberColumns_));
     1525  CoinAssert(minIndex >= 0);
     1526  if (detail) {
     1527    if (minIndex > 0 || maxIndex + 1 < (columnOrdered_ ? numberRows_ : numberColumns_))
     1528      printf("Not full range of indices - %d to %d\n", minIndex, maxIndex);
     1529  }
    15531530}
    15541531/* Given positive integer weights for each row fills in sum of weights
     
    15571534*/
    15581535CoinBigIndex *
    1559 ClpPlusMinusOneMatrix::dubiousWeights(const ClpSimplex * model, int * inputWeights) const
    1560 {
    1561      int numberRows = model->numberRows();
    1562      int numberColumns = model->numberColumns();
    1563      int number = numberRows + numberColumns;
    1564      CoinBigIndex * weights = new CoinBigIndex[number];
    1565      int i;
    1566      for (i = 0; i < numberColumns; i++) {
    1567           CoinBigIndex j;
    1568           CoinBigIndex count = 0;
    1569           for (j = startPositive_[i]; j < startPositive_[i+1]; j++) {
    1570                int iRow = indices_[j];
    1571                count += inputWeights[iRow];
    1572           }
    1573           weights[i] = count;
    1574      }
    1575      for (i = 0; i < numberRows; i++) {
    1576           weights[i+numberColumns] = inputWeights[i];
    1577      }
    1578      return weights;
     1536ClpPlusMinusOneMatrix::dubiousWeights(const ClpSimplex *model, int *inputWeights) const
     1537{
     1538  int numberRows = model->numberRows();
     1539  int numberColumns = model->numberColumns();
     1540  int number = numberRows + numberColumns;
     1541  CoinBigIndex *weights = new CoinBigIndex[number];
     1542  int i;
     1543  for (i = 0; i < numberColumns; i++) {
     1544    CoinBigIndex j;
     1545    CoinBigIndex count = 0;
     1546    for (j = startPositive_[i]; j < startPositive_[i + 1]; j++) {
     1547      int iRow = indices_[j];
     1548      count += inputWeights[iRow];
     1549    }
     1550    weights[i] = count;
     1551  }
     1552  for (i = 0; i < numberRows; i++) {
     1553    weights[i + numberColumns] = inputWeights[i];
     1554  }
     1555  return weights;
    15791556}
    15801557// Append Columns
    1581 void
    1582 ClpPlusMinusOneMatrix::appendCols(int number, const CoinPackedVectorBase * const * columns)
    1583 {
    1584      int iColumn;
    1585      CoinBigIndex size = 0;
    1586      int numberBad = 0;
    1587      // say if all +1
     1558void ClpPlusMinusOneMatrix::appendCols(int number, const CoinPackedVectorBase *const *columns)
     1559{
     1560  int iColumn;
     1561  CoinBigIndex size = 0;
     1562  int numberBad = 0;
     1563  // say if all +1
    15881564#ifdef CLP_PLUS_ONE_MATRIX
    1589      otherFlags_ = 0;
    1590 #endif
    1591      for (iColumn = 0; iColumn < number; iColumn++) {
    1592           int n = columns[iColumn]->getNumElements();
    1593           const double * element = columns[iColumn]->getElements();
    1594           size += n;
    1595           int i;
    1596           for (i = 0; i < n; i++) {
    1597                if (fabs(element[i]) != 1.0)
    1598                     numberBad++;
    1599           }
    1600      }
    1601      if (numberBad)
    1602           throw CoinError("Not +- 1", "appendCols", "ClpPlusMinusOneMatrix");
    1603      // Get rid of temporary arrays
    1604      delete [] lengths_;
    1605      lengths_ = NULL;
    1606      delete matrix_;
    1607      matrix_ = NULL;
    1608      CoinBigIndex numberNow = startPositive_[numberColumns_];
    1609      CoinBigIndex * temp;
    1610      temp = new CoinBigIndex [numberColumns_+1+number];
    1611      CoinMemcpyN(startPositive_, (numberColumns_ + 1), temp);
    1612      delete [] startPositive_;
    1613      startPositive_ = temp;
    1614      temp = new CoinBigIndex [numberColumns_+number];
    1615      CoinMemcpyN(startNegative_, numberColumns_, temp);
    1616      delete [] startNegative_;
    1617      startNegative_ = temp;
    1618      int * temp2 = new int [numberNow+size];
    1619      CoinMemcpyN(indices_, numberNow, temp2);
    1620      delete [] indices_;
    1621      indices_ = temp2;
    1622      // now add
    1623      size = numberNow;
    1624      for (iColumn = 0; iColumn < number; iColumn++) {
    1625           int n = columns[iColumn]->getNumElements();
    1626           const int * row = columns[iColumn]->getIndices();
    1627           const double * element = columns[iColumn]->getElements();
    1628           int i;
    1629           for (i = 0; i < n; i++) {
    1630                if (element[i] == 1.0)
    1631                     indices_[size++] = row[i];
    1632           }
    1633           startNegative_[iColumn+numberColumns_] = size;
    1634           for (i = 0; i < n; i++) {
    1635                if (element[i] == -1.0)
    1636                     indices_[size++] = row[i];
    1637           }
    1638           startPositive_[iColumn+numberColumns_+1] = size;
    1639      }
     1565  otherFlags_ = 0;
     1566#endif
     1567  for (iColumn = 0; iColumn < number; iColumn++) {
     1568    int n = columns[iColumn]->getNumElements();
     1569    const double *element = columns[iColumn]->getElements();
     1570    size += n;
     1571    int i;
     1572    for (i = 0; i < n; i++) {
     1573      if (fabs(element[i]) != 1.0)
     1574        numberBad++;
     1575    }
     1576  }
     1577  if (numberBad)
     1578    throw CoinError("Not +- 1", "appendCols", "ClpPlusMinusOneMatrix");
     1579  // Get rid of temporary arrays
     1580  delete[] lengths_;
     1581  lengths_ = NULL;
     1582  delete matrix_;
     1583  matrix_ = NULL;
     1584  CoinBigIndex numberNow = startPositive_[numberColumns_];
     1585  CoinBigIndex *temp;
     1586  temp = new CoinBigIndex[numberColumns_ + 1 + number];
     1587  CoinMemcpyN(startPositive_, (numberColumns_ + 1), temp);
     1588  delete[] startPositive_;
     1589  startPositive_ = temp;
     1590  temp = new CoinBigIndex[numberColumns_ + number];
     1591  CoinMemcpyN(startNegative_, numberColumns_, temp);
     1592  delete[] startNegative_;
     1593  startNegative_ = temp;
     1594  int *temp2 = new int[numberNow + size];
     1595  CoinMemcpyN(indices_, numberNow, temp2);
     1596  delete[] indices_;
     1597  indices_ = temp2;
     1598  // now add
     1599  size = numberNow;
     1600  for (iColumn = 0; iColumn < number; iColumn++) {
     1601    int n = columns[iColumn]->getNumElements();
     1602    const int *row = columns[iColumn]->getIndices();
     1603    const double *element = columns[iColumn]->getElements();
     1604    int i;
     1605    for (i = 0; i < n; i++) {
     1606      if (element[i] == 1.0)
     1607        indices_[size++] = row[i];
     1608    }
     1609    startNegative_[iColumn + numberColumns_] = size;
     1610    for (i = 0; i < n; i++) {
     1611      if (element[i] == -1.0)
     1612        indices_[size++] = row[i];
     1613    }
     1614    startPositive_[iColumn + numberColumns_ + 1] = size;
     1615  }
    16401616
    1641      numberColumns_ += number;
     1617  numberColumns_ += number;
    16421618}
    16431619// Append Rows
    1644 void
    1645 ClpPlusMinusOneMatrix::appendRows(int number, const CoinPackedVectorBase * const * rows)
    1646 {
    1647      // Allocate arrays to use for counting
    1648      int * countPositive = new int [numberColumns_+1];
    1649      memset(countPositive, 0, numberColumns_ * sizeof(int));
    1650      int * countNegative = new int [numberColumns_];
    1651      memset(countNegative, 0, numberColumns_ * sizeof(int));
    1652      int iRow;
    1653      CoinBigIndex size = 0;
    1654      int numberBad = 0;
    1655      // say if all +1
     1620void ClpPlusMinusOneMatrix::appendRows(int number, const CoinPackedVectorBase *const *rows)
     1621{
     1622  // Allocate arrays to use for counting
     1623  int *countPositive = new int[numberColumns_ + 1];
     1624  memset(countPositive, 0, numberColumns_ * sizeof(int));
     1625  int *countNegative = new int[numberColumns_];
     1626  memset(countNegative, 0, numberColumns_ * sizeof(int));
     1627  int iRow;
     1628  CoinBigIndex size = 0;
     1629  int numberBad = 0;
     1630  // say if all +1
    16561631#ifdef CLP_PLUS_ONE_MATRIX
    1657      otherFlags_ = 0;
    1658 #endif
    1659      for (iRow = 0; iRow < number; iRow++) {
    1660           int n = rows[iRow]->getNumElements();
    1661           const int * column = rows[iRow]->getIndices();
    1662           const double * element = rows[iRow]->getElements();
    1663           size += n;
    1664           int i;
    1665           for (i = 0; i < n; i++) {
    1666                int iColumn = column[i];
    1667                if (element[i] == 1.0)
    1668                     countPositive[iColumn]++;
    1669                else if (element[i] == -1.0)
    1670                     countNegative[iColumn]++;
    1671                else
    1672                     numberBad++;
    1673           }
    1674      }
    1675      if (numberBad)
    1676           throw CoinError("Not +- 1", "appendRows", "ClpPlusMinusOneMatrix");
    1677      // Get rid of temporary arrays
    1678      delete [] lengths_;
    1679      lengths_ = NULL;
    1680      delete matrix_;
    1681      matrix_ = NULL;
    1682      CoinBigIndex numberNow = startPositive_[numberColumns_];
    1683      int * newIndices = new int [numberNow+size];
    1684      // Update starts and turn counts into positions
    1685      // also move current indices
    1686      int iColumn;
    1687      CoinBigIndex numberAdded = 0;
    1688      for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    1689           int n, move;
    1690           CoinBigIndex now;
    1691           now = startPositive_[iColumn];
    1692           move = static_cast<int>(startNegative_[iColumn] - now);
    1693           n = countPositive[iColumn];
    1694           startPositive_[iColumn] += numberAdded;
    1695           CoinMemcpyN(indices_+now, move, newIndices + startPositive_[iColumn]);
    1696           countPositive[iColumn] = static_cast<int>(startNegative_[iColumn] + numberAdded);
    1697           numberAdded += n;
    1698           now = startNegative_[iColumn];
    1699           move = static_cast<int>(startPositive_[iColumn+1] - now);
    1700           n = countNegative[iColumn];
    1701           startNegative_[iColumn] += numberAdded;
    1702           CoinMemcpyN(indices_+now, move, newIndices + startNegative_[iColumn]);
    1703           countNegative[iColumn] = static_cast<int>(startPositive_[iColumn+1] + numberAdded);
    1704           numberAdded += n;
    1705      }
    1706      delete [] indices_;
    1707      indices_ = newIndices;
    1708      startPositive_[numberColumns_] += numberAdded;
    1709      // Now put in
    1710      for (iRow = 0; iRow < number; iRow++) {
    1711           int newRow = numberRows_ + iRow;
    1712           int n = rows[iRow]->getNumElements();
    1713           const int * column = rows[iRow]->getIndices();
    1714           const double * element = rows[iRow]->getElements();
    1715           int i;
    1716           for (i = 0; i < n; i++) {
    1717                int iColumn = column[i];
    1718                int put;
    1719                if (element[i] == 1.0) {
    1720                     put = countPositive[iColumn];
    1721                     countPositive[iColumn] = put + 1;
    1722                } else {
    1723                     put = countNegative[iColumn];
    1724                     countNegative[iColumn] = put + 1;
    1725                }
    1726                indices_[put] = newRow;
    1727           }
    1728      }
    1729      delete [] countPositive;
    1730      delete [] countNegative;
    1731      numberRows_ += number;
     1632  otherFlags_ = 0;
     1633#endif
     1634  for (iRow = 0; iRow < number; iRow++) {
     1635    int n = rows[iRow]->getNumElements();
     1636    const int *column = rows[iRow]->getIndices();
     1637    const double *element = rows[iRow]->getElements();
     1638    size += n;
     1639    int i;
     1640    for (i = 0; i < n; i++) {
     1641      int iColumn = column[i];
     1642      if (element[i] == 1.0)
     1643        countPositive[iColumn]++;
     1644      else if (element[i] == -1.0)
     1645        countNegative[iColumn]++;
     1646      else
     1647        numberBad++;
     1648    }
     1649  }
     1650  if (numberBad)
     1651    throw CoinError("Not +- 1", "appendRows", "ClpPlusMinusOneMatrix");
     1652  // Get rid of temporary arrays
     1653  delete[] lengths_;
     1654  lengths_ = NULL;
     1655  delete matrix_;
     1656  matrix_ = NULL;
     1657  CoinBigIndex numberNow = startPositive_[numberColumns_];
     1658  int *newIndices = new int[numberNow + size];
     1659  // Update starts and turn counts into positions
     1660  // also move current indices
     1661  int iColumn;
     1662  CoinBigIndex numberAdded = 0;
     1663  for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     1664    int n, move;
     1665    CoinBigIndex now;
     1666    now = startPositive_[iColumn];
     1667    move = static_cast< int >(startNegative_[iColumn] - now);
     1668    n = countPositive[iColumn];
     1669    startPositive_[iColumn] += numberAdded;
     1670    CoinMemcpyN(indices_ + now, move, newIndices + startPositive_[iColumn]);
     1671    countPositive[iColumn] = static_cast< int >(startNegative_[iColumn] + numberAdded);
     1672    numberAdded += n;
     1673    now = startNegative_[iColumn];
     1674    move = static_cast< int >(startPositive_[iColumn + 1] - now);
     1675    n = countNegative[iColumn];
     1676    startNegative_[iColumn] += numberAdded;
     1677    CoinMemcpyN(indices_ + now, move, newIndices + startNegative_[iColumn]);
     1678    countNegative[iColumn] = static_cast< int >(startPositive_[iColumn + 1] + numberAdded);
     1679    numberAdded += n;
     1680  }
     1681  delete[] indices_;
     1682  indices_ = newIndices;
     1683  startPositive_[numberColumns_] += numberAdded;
     1684  // Now put in
     1685  for (iRow = 0; iRow < number; iRow++) {
     1686    int newRow = numberRows_ + iRow;
     1687    int n = rows[iRow]->getNumElements();
     1688    const int *column = rows[iRow]->getIndices();
     1689    const double *element = rows[iRow]->getElements();
     1690    int i;
     1691    for (i = 0; i < n; i++) {
     1692      int iColumn = column[i];
     1693      int put;
     1694      if (element[i] == 1.0) {
     1695        put = countPositive[iColumn];
     1696        countPositive[iColumn] = put + 1;
     1697      } else {
     1698        put = countNegative[iColumn];
     1699        countNegative[iColumn] = put + 1;
     1700      }
     1701      indices_[put] = newRow;
     1702    }
     1703  }
     1704  delete[] countPositive;
     1705  delete[] countNegative;
     1706  numberRows_ += number;
    17321707}
    17331708/* Returns largest and smallest elements of both signs.
    17341709   Largest refers to largest absolute value.
    17351710*/
    1736 void
    1737 ClpPlusMinusOneMatrix::rangeOfElements(double & smallestNegative, double & largestNegative,
    1738                                        double & smallestPositive, double & largestPositive)
    1739 {
    1740      int iColumn;
    1741      bool plusOne = false;
    1742      bool minusOne = false;
    1743      for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    1744           if (startNegative_[iColumn] > startPositive_[iColumn])
    1745                plusOne = true;
    1746           if (startPositive_[iColumn+1] > startNegative_[iColumn])
    1747                minusOne = true;
    1748      }
    1749      if (minusOne) {
    1750           smallestNegative = -1.0;
    1751           largestNegative = -1.0;
    1752      } else {
    1753           smallestNegative = 0.0;
    1754           largestNegative = 0.0;
    1755      }
    1756      if (plusOne) {
    1757           smallestPositive = 1.0;
    1758           largestPositive = 1.0;
    1759      } else {
    1760           smallestPositive = 0.0;
    1761           largestPositive = 0.0;
    1762      }
     1711void ClpPlusMinusOneMatrix::rangeOfElements(double &smallestNegative, double &largestNegative,
     1712  double &smallestPositive, double &largestPositive)
     1713{
     1714  int iColumn;
     1715  bool plusOne = false;
     1716  bool minusOne = false;
     1717  for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     1718    if (startNegative_[iColumn] > startPositive_[iColumn])
     1719      plusOne = true;
     1720    if (startPositive_[iColumn + 1] > startNegative_[iColumn])
     1721      minusOne = true;
     1722  }
     1723  if (minusOne) {
     1724    smallestNegative = -1.0;
     1725    largestNegative = -1.0;
     1726  } else {
     1727    smallestNegative = 0.0;
     1728    largestNegative = 0.0;
     1729  }
     1730  if (plusOne) {
     1731    smallestPositive = 1.0;
     1732    largestPositive = 1.0;
     1733  } else {
     1734    smallestPositive = 0.0;
     1735    largestPositive = 0.0;
     1736  }
    17631737}
    17641738// Says whether it can do partial pricing
    1765 bool
    1766 ClpPlusMinusOneMatrix::canDoPartialPricing() const
    1767 {
    1768      return true;
     1739bool ClpPlusMinusOneMatrix::canDoPartialPricing() const
     1740{
     1741  return true;
    17691742}
    17701743// Partial pricing
    1771 void
    1772 ClpPlusMinusOneMatrix::partialPricing(ClpSimplex * model, double startFraction, double endFraction,
    1773                                       int & bestSequence, int & numberWanted)
    1774 {
    1775      numberWanted = currentWanted_;
    1776      int start = static_cast<int> (startFraction * numberColumns_);
    1777      int end = CoinMin(static_cast<int> (endFraction * numberColumns_ + 1), numberColumns_);
    1778      CoinBigIndex j;
    1779      double tolerance = model->currentDualTolerance();
    1780      double * COIN_RESTRICT reducedCost = model->djRegion();
    1781      const double * COIN_RESTRICT duals = model->dualRowSolution();
    1782      const double * COIN_RESTRICT cost = model->costRegion();
    1783      double bestDj;
    1784      if (bestSequence >= 0)
    1785           bestDj = fabs(reducedCost[bestSequence]);
    1786      else
    1787           bestDj = tolerance;
    1788      int sequenceOut = model->sequenceOut();
    1789      int saveSequence = bestSequence;
    1790      int iSequence;
    1791      for (iSequence = start; iSequence < end; iSequence++) {
    1792           if (iSequence != sequenceOut) {
    1793                double value;
    1794                ClpSimplex::Status status = model->getStatus(iSequence);
     1744void ClpPlusMinusOneMatrix::partialPricing(ClpSimplex *model, double startFraction, double endFraction,
     1745  int &bestSequence, int &numberWanted)
     1746{
     1747  numberWanted = currentWanted_;
     1748  int start = static_cast< int >(startFraction * numberColumns_);
     1749  int end = CoinMin(static_cast< int >(endFraction * numberColumns_ + 1), numberColumns_);
     1750  CoinBigIndex j;
     1751  double tolerance = model->currentDualTolerance();
     1752  double *COIN_RESTRICT reducedCost = model->djRegion();
     1753  const double *COIN_RESTRICT duals = model->dualRowSolution();
     1754  const double *COIN_RESTRICT cost = model->costRegion();
     1755  double bestDj;
     1756  if (bestSequence >= 0)
     1757    bestDj = fabs(reducedCost[bestSequence]);
     1758  else
     1759    bestDj = tolerance;
     1760  int sequenceOut = model->sequenceOut();
     1761  int saveSequence = bestSequence;
     1762  int iSequence;
     1763  for (iSequence = start; iSequence < end; iSequence++) {
     1764    if (iSequence != sequenceOut) {
     1765      double value;
     1766      ClpSimplex::Status status = model->getStatus(iSequence);
    17951767
    1796                switch(status) {
     1768      switch (status) {
    17971769
    1798                case ClpSimplex::basic:
    1799                case ClpSimplex::isFixed:
    1800                     break;
    1801                case ClpSimplex::isFree:
    1802                case ClpSimplex::superBasic:
    1803                     value = cost[iSequence];
    1804                     j = startPositive_[iSequence];
    1805                     for (; j < startNegative_[iSequence]; j++) {
    1806                          int iRow = indices_[j];
    1807                          value -= duals[iRow];
    1808                     }
    1809                     for (; j < startPositive_[iSequence+1]; j++) {
    1810                          int iRow = indices_[j];
    1811                          value += duals[iRow];
    1812                     }
    1813                     value = fabs(value);
    1814                     if (value > FREE_ACCEPT * tolerance) {
    1815                          numberWanted--;
    1816                          // we are going to bias towards free (but only if reasonable)
    1817                          value *= FREE_BIAS;
    1818                          if (value > bestDj) {
    1819                               // check flagged variable and correct dj
    1820                               if (!model->flagged(iSequence)) {
    1821                                    bestDj = value;
    1822                                    bestSequence = iSequence;
    1823                               } else {
    1824                                    // just to make sure we don't exit before got something
    1825                                    numberWanted++;
    1826                               }
    1827                          }
    1828                     }
    1829                     break;
    1830                case ClpSimplex::atUpperBound:
    1831                     value = cost[iSequence];
    1832                     j = startPositive_[iSequence];
    1833                     for (; j < startNegative_[iSequence]; j++) {
    1834                          int iRow = indices_[j];
    1835                          value -= duals[iRow];
    1836                     }
    1837                     for (; j < startPositive_[iSequence+1]; j++) {
    1838                          int iRow = indices_[j];
    1839                          value += duals[iRow];
    1840                     }
    1841                     if (value > tolerance) {
    1842                          numberWanted--;
    1843                          if (value > bestDj) {
    1844                               // check flagged variable and correct dj
    1845                               if (!model->flagged(iSequence)) {
    1846                                    bestDj = value;
    1847                                    bestSequence = iSequence;
    1848                               } else {
    1849                                    // just to make sure we don't exit before got something
    1850                                    numberWanted++;
    1851                               }
    1852                          }
    1853                     }
    1854                     break;
    1855                case ClpSimplex::atLowerBound:
    1856                     value = cost[iSequence];
    1857                     j = startPositive_[iSequence];
    1858                     for (; j < startNegative_[iSequence]; j++) {
    1859                          int iRow = indices_[j];
    1860                          value -= duals[iRow];
    1861                     }
    1862                     for (; j < startPositive_[iSequence+1]; j++) {
    1863                          int iRow = indices_[j];
    1864                          value += duals[iRow];
    1865                     }
    1866                     value = -value;
    1867                     if (value > tolerance) {
    1868                          numberWanted--;
    1869                          if (value > bestDj) {
    1870                               // check flagged variable and correct dj
    1871                               if (!model->flagged(iSequence)) {
    1872                                    bestDj = value;
    1873                                    bestSequence = iSequence;
    1874                               } else {
    1875                                    // just to make sure we don't exit before got something
    1876                                    numberWanted++;
    1877                               }
    1878                          }
    1879                     }
    1880                     break;
    1881                }
     1770      case ClpSimplex::basic:
     1771      case ClpSimplex::isFixed:
     1772        break;
     1773      case ClpSimplex::isFree:
     1774      case ClpSimplex::superBasic:
     1775        value = cost[iSequence];
     1776        j = startPositive_[iSequence];
     1777        for (; j < startNegative_[iSequence]; j++) {
     1778          int iRow = indices_[j];
     1779          value -= duals[iRow];
     1780        }
     1781        for (; j < startPositive_[iSequence + 1]; j++) {
     1782          int iRow = indices_[j];
     1783          value += duals[iRow];
     1784        }
     1785        value = fabs(value);
     1786        if (value > FREE_ACCEPT * tolerance) {
     1787          numberWanted--;
     1788          // we are going to bias towards free (but only if reasonable)
     1789          value *= FREE_BIAS;
     1790          if (value > bestDj) {
     1791            // check flagged variable and correct dj
     1792            if (!model->flagged(iSequence)) {
     1793              bestDj = value;
     1794              bestSequence = iSequence;
     1795            } else {
     1796              // just to make sure we don't exit before got something
     1797              numberWanted++;
     1798            }
    18821799          }
    1883           if (!numberWanted)
    1884                break;
    1885      }
    1886      if (bestSequence != saveSequence) {
    1887           // recompute dj
    1888           double value = cost[bestSequence];
    1889           j = startPositive_[bestSequence];
    1890           for (; j < startNegative_[bestSequence]; j++) {
    1891                int iRow = indices_[j];
    1892                value -= duals[iRow];
     1800        }
     1801        break;
     1802      case ClpSimplex::atUpperBound:
     1803        value = cost[iSequence];
     1804        j = startPositive_[iSequence];
     1805        for (; j < startNegative_[iSequence]; j++) {
     1806          int iRow = indices_[j];
     1807          value -= duals[iRow];
     1808        }
     1809        for (; j < startPositive_[iSequence + 1]; j++) {
     1810          int iRow = indices_[j];
     1811          value += duals[iRow];
     1812        }
     1813        if (value > tolerance) {
     1814          numberWanted--;
     1815          if (value > bestDj) {
     1816            // check flagged variable and correct dj
     1817            if (!model->flagged(iSequence)) {
     1818              bestDj = value;
     1819              bestSequence = iSequence;
     1820            } else {
     1821              // just to make sure we don't exit before got something
     1822              numberWanted++;
     1823            }
    18931824          }
    1894           for (; j < startPositive_[bestSequence+1]; j++) {
    1895                int iRow = indices_[j];
    1896                value += duals[iRow];
     1825        }
     1826        break;
     1827      case ClpSimplex::atLowerBound:
     1828        value = cost[iSequence];
     1829        j = startPositive_[iSequence];
     1830        for (; j < startNegative_[iSequence]; j++) {
     1831          int iRow = indices_[j];
     1832          value -= duals[iRow];
     1833        }
     1834        for (; j < startPositive_[iSequence + 1]; j++) {
     1835          int iRow = indices_[j];
     1836          value += duals[iRow];
     1837        }
     1838        value = -value;
     1839        if (value > tolerance) {
     1840          numberWanted--;
     1841          if (value > bestDj) {
     1842            // check flagged variable and correct dj
     1843            if (!model->flagged(iSequence)) {
     1844              bestDj = value;
     1845              bestSequence = iSequence;
     1846            } else {
     1847              // just to make sure we don't exit before got something
     1848              numberWanted++;
     1849            }
    18971850          }
    1898           reducedCost[bestSequence] = value;
    1899           savedBestSequence_ = bestSequence;
    1900           savedBestDj_ = reducedCost[savedBestSequence_];
    1901      }
    1902      currentWanted_ = numberWanted;
     1851        }
     1852        break;
     1853      }
     1854    }
     1855    if (!numberWanted)
     1856      break;
     1857  }
     1858  if (bestSequence != saveSequence) {
     1859    // recompute dj
     1860    double value = cost[bestSequence];
     1861    j = startPositive_[bestSequence];
     1862    for (; j < startNegative_[bestSequence]; j++) {
     1863      int iRow = indices_[j];
     1864      value -= duals[iRow];
     1865    }
     1866    for (; j < startPositive_[bestSequence + 1]; j++) {
     1867      int iRow = indices_[j];
     1868      value += duals[iRow];
     1869    }
     1870    reducedCost[bestSequence] = value;
     1871    savedBestSequence_ = bestSequence;
     1872    savedBestDj_ = reducedCost[savedBestSequence_];
     1873  }
     1874  currentWanted_ = numberWanted;
    19031875}
    19041876// Allow any parts of a created CoinMatrix to be deleted
    1905 void
    1906 ClpPlusMinusOneMatrix::releasePackedMatrix() const
    1907 {
    1908      delete matrix_;
    1909      delete [] lengths_;
    1910      matrix_ = NULL;
    1911      lengths_ = NULL;
    1912      // say if all +1
     1877void ClpPlusMinusOneMatrix::releasePackedMatrix() const
     1878{
     1879  delete matrix_;
     1880  delete[] lengths_;
     1881  matrix_ = NULL;
     1882  lengths_ = NULL;
     1883  // say if all +1
    19131884#ifdef CLP_PLUS_ONE_MATRIX
    1914      otherFlags_ = 0;
     1885  otherFlags_ = 0;
    19151886#endif
    19161887}
    19171888/* Returns true if can combine transposeTimes and subsetTransposeTimes
    19181889   and if it would be faster */
    1919 bool
    1920 ClpPlusMinusOneMatrix::canCombine(const ClpSimplex * model,
    1921                                   const CoinIndexedVector * pi) const
    1922 {
    1923      int numberInRowArray = pi->getNumElements();
    1924      int numberRows = model->numberRows();
    1925      bool packed = pi->packedMode();
    1926      // factor should be smaller if doing both with two pi vectors
    1927      double factor = 0.27;
    1928      // We may not want to do by row if there may be cache problems
    1929      // It would be nice to find L2 cache size - for moment 512K
    1930      // Be slightly optimistic
    1931      if (numberColumns_ * sizeof(double) > 1000000) {
    1932           if (numberRows * 10 < numberColumns_)
    1933                factor *= 0.333333333;
    1934           else if (numberRows * 4 < numberColumns_)
    1935                factor *= 0.5;
    1936           else if (numberRows * 2 < numberColumns_)
    1937                factor *= 0.66666666667;
    1938           //if (model->numberIterations()%50==0)
    1939           //printf("%d nonzero\n",numberInRowArray);
    1940      }
    1941      // if not packed then bias a bit more towards by column
    1942      if (!packed)
    1943           factor *= 0.9;
    1944      return (numberInRowArray > factor * numberRows || !model->rowCopy());
     1890bool ClpPlusMinusOneMatrix::canCombine(const ClpSimplex *model,
     1891  const CoinIndexedVector *pi) const
     1892{
     1893  int numberInRowArray = pi->getNumElements();
     1894  int numberRows = model->numberRows();
     1895  bool packed = pi->packedMode();
     1896  // factor should be smaller if doing both with two pi vectors
     1897  double factor = 0.27;
     1898  // We may not want to do by row if there may be cache problems
     1899  // It would be nice to find L2 cache size - for moment 512K
     1900  // Be slightly optimistic
     1901  if (numberColumns_ * sizeof(double) > 1000000) {
     1902    if (numberRows * 10 < numberColumns_)
     1903      factor *= 0.333333333;
     1904    else if (numberRows * 4 < numberColumns_)
     1905      factor *= 0.5;
     1906    else if (numberRows * 2 < numberColumns_)
     1907      factor *= 0.66666666667;
     1908    //if (model->numberIterations()%50==0)
     1909    //printf("%d nonzero\n",numberInRowArray);
     1910  }
     1911  // if not packed then bias a bit more towards by column
     1912  if (!packed)
     1913    factor *= 0.9;
     1914  return (numberInRowArray > factor * numberRows || !model->rowCopy());
    19451915}
    19461916// These have to match ClpPrimalColumnSteepest version
    1947 #define reference(i)  (((reference[i>>5]>>(i&31))&1)!=0)
     1917#define reference(i) (((reference[i >> 5] >> (i & 31)) & 1) != 0)
    19481918/* Updates two arrays for steepest and does devex weights
    19491919   Returns nonzero if updates reduced cost and infeas -
    19501920   new infeas in dj1  - This does not at present*/
    1951 int
    1952 ClpPlusMinusOneMatrix::transposeTimes2(const ClpSimplex * model,
    1953                                        const CoinIndexedVector * pi1, CoinIndexedVector * dj1,
    1954                                        const CoinIndexedVector * pi2,
    1955                                        CoinIndexedVector * spare,
    1956                                        double * , double * ,
    1957                                        double referenceIn, double devex,
    1958                                        // Array for exact devex to say what is in reference framework
    1959                                        unsigned int * COIN_RESTRICT reference,
    1960                                        double * COIN_RESTRICT weights, double scaleFactor)
    1961 {
    1962      // put row of tableau in dj1
    1963      double * COIN_RESTRICT pi = pi1->denseVector();
    1964      int numberNonZero = 0;
    1965      int * COIN_RESTRICT index = dj1->getIndices();
    1966      double * COIN_RESTRICT array = dj1->denseVector();
    1967      int numberInRowArray = pi1->getNumElements();
    1968      double zeroTolerance = model->zeroTolerance();
    1969      bool packed = pi1->packedMode();
    1970      // do by column
    1971      int iColumn;
    1972      assert (!spare->getNumElements());
    1973      double * COIN_RESTRICT piWeight = pi2->denseVector();
    1974      assert (!pi2->packedMode());
    1975      bool killDjs = (scaleFactor == 0.0);
    1976      if (!scaleFactor)
    1977           scaleFactor = 1.0;
    1978      // Note scale factor was -1.0
    1979      if (packed) {
    1980           // need to expand pi into y
    1981           assert(spare->capacity() >= model->numberRows());
    1982           double * COIN_RESTRICT piOld = pi;
    1983           pi = spare->denseVector();
    1984           const int * COIN_RESTRICT whichRow = pi1->getIndices();
    1985           int i;
    1986           // modify pi so can collapse to one loop
    1987           for (i = 0; i < numberInRowArray; i++) {
    1988                int iRow = whichRow[i];
    1989                pi[iRow] = piOld[i];
     1921int ClpPlusMinusOneMatrix::transposeTimes2(const ClpSimplex *model,
     1922  const CoinIndexedVector *pi1, CoinIndexedVector *dj1,
     1923  const CoinIndexedVector *pi2,
     1924  CoinIndexedVector *spare,
     1925  double *, double *,
     1926  double referenceIn, double devex,
     1927  // Array for exact devex to say what is in reference framework
     1928  unsigned int *COIN_RESTRICT reference,
     1929  double *COIN_RESTRICT weights, double scaleFactor)
     1930{
     1931  // put row of tableau in dj1
     1932  double *COIN_RESTRICT pi = pi1->denseVector();
     1933  int numberNonZero = 0;
     1934  int *COIN_RESTRICT index = dj1->getIndices();
     1935  double *COIN_RESTRICT array = dj1->denseVector();
     1936  int numberInRowArray = pi1->getNumElements();
     1937  double zeroTolerance = model->zeroTolerance();
     1938  bool packed = pi1->packedMode();
     1939  // do by column
     1940  int iColumn;
     1941  assert(!spare->getNumElements());
     1942  double *COIN_RESTRICT piWeight = pi2->denseVector();
     1943  assert(!pi2->packedMode());
     1944  bool killDjs = (scaleFactor == 0.0);
     1945  if (!scaleFactor)
     1946    scaleFactor = 1.0;
     1947  // Note scale factor was -1.0
     1948  if (packed) {
     1949    // need to expand pi into y
     1950    assert(spare->capacity() >= model->numberRows());
     1951    double *COIN_RESTRICT piOld = pi;
     1952    pi = spare->denseVector();
     1953    const int *COIN_RESTRICT whichRow = pi1->getIndices();
     1954    int i;
     1955    // modify pi so can collapse to one loop
     1956    for (i = 0; i < numberInRowArray; i++) {
     1957      int iRow = whichRow[i];
     1958      pi[iRow] = piOld[i];
     1959    }
     1960    CoinBigIndex j;
     1961    for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     1962      ClpSimplex::Status status = model->getStatus(iColumn);
     1963      if (status == ClpSimplex::basic || status == ClpSimplex::isFixed)
     1964        continue;
     1965      double value = 0.0;
     1966      for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
     1967        int iRow = indices_[j];
     1968        value -= pi[iRow];
     1969      }
     1970      for (; j < startPositive_[iColumn + 1]; j++) {
     1971        int iRow = indices_[j];
     1972        value += pi[iRow];
     1973      }
     1974      if (fabs(value) > zeroTolerance) {
     1975        // and do other array
     1976        double modification = 0.0;
     1977        for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
     1978          int iRow = indices_[j];
     1979          modification += piWeight[iRow];
     1980        }
     1981        for (; j < startPositive_[iColumn + 1]; j++) {
     1982          int iRow = indices_[j];
     1983          modification -= piWeight[iRow];
     1984        }
     1985        double thisWeight = weights[iColumn];
     1986        double pivot = value * scaleFactor;
     1987        double pivotSquared = pivot * pivot;
     1988        thisWeight += pivotSquared * devex + pivot * modification;
     1989        if (thisWeight < DEVEX_TRY_NORM) {
     1990          if (referenceIn < 0.0) {
     1991            // steepest
     1992            thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
     1993          } else {
     1994            // exact
     1995            thisWeight = referenceIn * pivotSquared;
     1996            if (reference(iColumn))
     1997              thisWeight += 1.0;
     1998            thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
    19901999          }
    1991           CoinBigIndex j;
    1992           for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    1993                ClpSimplex::Status status = model->getStatus(iColumn);
    1994                if (status == ClpSimplex::basic || status == ClpSimplex::isFixed) continue;
    1995                double value = 0.0;
    1996                for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
    1997                     int iRow = indices_[j];
    1998                     value -= pi[iRow];
    1999                }
    2000                for (; j < startPositive_[iColumn+1]; j++) {
    2001                     int iRow = indices_[j];
    2002                     value += pi[iRow];
    2003                }
    2004                if (fabs(value) > zeroTolerance) {
    2005                     // and do other array
    2006                     double modification = 0.0;
    2007                     for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
    2008                          int iRow = indices_[j];
    2009                          modification += piWeight[iRow];
    2010                     }
    2011                     for (; j < startPositive_[iColumn+1]; j++) {
    2012                          int iRow = indices_[j];
    2013                          modification -= piWeight[iRow];
    2014                     }
    2015                     double thisWeight = weights[iColumn];
    2016                     double pivot = value * scaleFactor;
    2017                     double pivotSquared = pivot * pivot;
    2018                     thisWeight += pivotSquared * devex + pivot * modification;
    2019                     if (thisWeight < DEVEX_TRY_NORM) {
    2020                          if (referenceIn < 0.0) {
    2021                               // steepest
    2022                               thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
    2023                          } else {
    2024                               // exact
    2025                               thisWeight = referenceIn * pivotSquared;
    2026                               if (reference(iColumn))
    2027                                    thisWeight += 1.0;
    2028                               thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
    2029                          }
    2030                     }
    2031                     weights[iColumn] = thisWeight;
    2032                     if (!killDjs) {
    2033                          array[numberNonZero] = value;
    2034                          index[numberNonZero++] = iColumn;
    2035                     }
    2036                }
     2000        }
     2001        weights[iColumn] = thisWeight;
     2002        if (!killDjs) {
     2003          array[numberNonZero] = value;
     2004          index[numberNonZero++] = iColumn;
     2005        }
     2006      }
     2007    }
     2008    // zero out
     2009    for (i = 0; i < numberInRowArray; i++) {
     2010      int iRow = whichRow[i];
     2011      pi[iRow] = 0.0;
     2012    }
     2013  } else {
     2014    CoinBigIndex j;
     2015    for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     2016      ClpSimplex::Status status = model->getStatus(iColumn);
     2017      if (status == ClpSimplex::basic || status == ClpSimplex::isFixed)
     2018        continue;
     2019      double value = 0.0;
     2020      for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
     2021        int iRow = indices_[j];
     2022        value -= pi[iRow];
     2023      }
     2024      for (; j < startPositive_[iColumn + 1]; j++) {
     2025        int iRow = indices_[j];
     2026        value += pi[iRow];
     2027      }
     2028      if (fabs(value) > zeroTolerance) {
     2029        // and do other array
     2030        double modification = 0.0;
     2031        for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
     2032          int iRow = indices_[j];
     2033          modification += piWeight[iRow];
     2034        }
     2035        for (; j < startPositive_[iColumn + 1]; j++) {
     2036          int iRow = indices_[j];
     2037          modification -= piWeight[iRow];
     2038        }
     2039        double thisWeight = weights[iColumn];
     2040        double pivot = value * scaleFactor;
     2041        double pivotSquared = pivot * pivot;
     2042        thisWeight += pivotSquared * devex + pivot * modification;
     2043        if (thisWeight < DEVEX_TRY_NORM) {
     2044          if (referenceIn < 0.0) {
     2045            // steepest
     2046            thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
     2047          } else {
     2048            // exact
     2049            thisWeight = referenceIn * pivotSquared;
     2050            if (reference(iColumn))
     2051              thisWeight += 1.0;
     2052            thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
    20372053          }
    2038           // zero out
    2039           for (i = 0; i < numberInRowArray; i++) {
    2040                int iRow = whichRow[i];
    2041                pi[iRow] = 0.0;
    2042           }
    2043      } else {
    2044           CoinBigIndex j;
    2045           for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    2046                ClpSimplex::Status status = model->getStatus(iColumn);
    2047                if (status == ClpSimplex::basic || status == ClpSimplex::isFixed) continue;
    2048                double value = 0.0;
    2049                for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
    2050                     int iRow = indices_[j];
    2051                     value -= pi[iRow];
    2052                }
    2053                for (; j < startPositive_[iColumn+1]; j++) {
    2054                     int iRow = indices_[j];
    2055                     value += pi[iRow];
    2056                }
    2057                if (fabs(value) > zeroTolerance) {
    2058                     // and do other array
    2059                     double modification = 0.0;
    2060                     for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
    2061                          int iRow = indices_[j];
    2062                          modification += piWeight[iRow];
    2063                     }
    2064                     for (; j < startPositive_[iColumn+1]; j++) {
    2065                          int iRow = indices_[j];
    2066                          modification -= piWeight[iRow];
    2067                     }
    2068                     double thisWeight = weights[iColumn];
    2069                     double pivot = value * scaleFactor;
    2070                     double pivotSquared = pivot * pivot;
    2071                     thisWeight += pivotSquared * devex + pivot * modification;
    2072                     if (thisWeight < DEVEX_TRY_NORM) {
    2073                          if (referenceIn < 0.0) {
    2074                               // steepest
    2075                               thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
    2076                          } else {
    2077                               // exact
    2078                               thisWeight = referenceIn * pivotSquared;
    2079                               if (reference(iColumn))
    2080                                    thisWeight += 1.0;
    2081                               thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
    2082                          }
    2083                     }
    2084                     weights[iColumn] = thisWeight;
    2085                     if (!killDjs) {
    2086                          array[iColumn] = value;
    2087                          index[numberNonZero++] = iColumn;
    2088                     }
    2089                }
    2090           }
    2091      }
    2092      dj1->setNumElements(numberNonZero);
    2093      spare->setNumElements(0);
    2094      if (packed)
    2095           dj1->setPackedMode(true);
    2096      return 0;
     2054        }
     2055        weights[iColumn] = thisWeight;
     2056        if (!killDjs) {
     2057          array[iColumn] = value;
     2058          index[numberNonZero++] = iColumn;
     2059        }
     2060      }
     2061    }
     2062  }
     2063  dj1->setNumElements(numberNonZero);
     2064  spare->setNumElements(0);
     2065  if (packed)
     2066    dj1->setPackedMode(true);
     2067  return 0;
    20972068}
    20982069// Updates second array for steepest and does devex weights
    2099 void
    2100 ClpPlusMinusOneMatrix::subsetTimes2(const ClpSimplex * ,
    2101                                     CoinIndexedVector * dj1,
    2102                                     const CoinIndexedVector * pi2, CoinIndexedVector *,
    2103                                     double referenceIn, double devex,
    2104                                     // Array for exact devex to say what is in reference framework
    2105                                     unsigned int * COIN_RESTRICT reference,
    2106                                     double * COIN_RESTRICT weights, double scaleFactor)
    2107 {
    2108      int number = dj1->getNumElements();
    2109      const int * COIN_RESTRICT index = dj1->getIndices();
    2110      double * COIN_RESTRICT array = dj1->denseVector();
    2111      assert( dj1->packedMode());
     2070void ClpPlusMinusOneMatrix::subsetTimes2(const ClpSimplex *,
     2071  CoinIndexedVector *dj1,
     2072  const CoinIndexedVector *pi2, CoinIndexedVector *,
     2073  double referenceIn, double devex,
     2074  // Array for exact devex to say what is in reference framework
     2075  unsigned int *COIN_RESTRICT reference,
     2076  double *COIN_RESTRICT weights, double scaleFactor)
     2077{
     2078  int number = dj1->getNumElements();
     2079  const int *COIN_RESTRICT index = dj1->getIndices();
     2080  double *COIN_RESTRICT array = dj1->denseVector();
     2081  assert(dj1->packedMode());
    21122082
    2113      double * COIN_RESTRICT piWeight = pi2->denseVector();
    2114      bool killDjs = (scaleFactor == 0.0);
    2115      if (!scaleFactor)
    2116           scaleFactor = 1.0;
    2117      for (int k = 0; k < number; k++) {
    2118           int iColumn = index[k];
    2119           double pivot = array[k] * scaleFactor;
    2120           if (killDjs)
    2121                array[k] = 0.0;
    2122           // and do other array
    2123           double modification = 0.0;
    2124           CoinBigIndex j;
    2125           for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
    2126                int iRow = indices_[j];
    2127                modification += piWeight[iRow];
    2128           }
    2129           for (; j < startPositive_[iColumn+1]; j++) {
    2130                int iRow = indices_[j];
    2131                modification -= piWeight[iRow];
    2132           }
    2133           double thisWeight = weights[iColumn];
    2134           double pivotSquared = pivot * pivot;
    2135           thisWeight += pivotSquared * devex + pivot * modification;
    2136           if (thisWeight < DEVEX_TRY_NORM) {
    2137                if (referenceIn < 0.0) {
    2138                     // steepest
    2139                     thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
    2140                } else {
    2141                     // exact
    2142                     thisWeight = referenceIn * pivotSquared;
    2143                     if (reference(iColumn))
    2144                          thisWeight += 1.0;
    2145                     thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
    2146                }
    2147           }
    2148           weights[iColumn] = thisWeight;
    2149      }
     2083  double *COIN_RESTRICT piWeight = pi2->denseVector();
     2084  bool killDjs = (scaleFactor == 0.0);
     2085  if (!scaleFactor)
     2086    scaleFactor = 1.0;
     2087  for (int k = 0; k < number; k++) {
     2088    int iColumn = index[k];
     2089    double pivot = array[k] * scaleFactor;
     2090    if (killDjs)
     2091      array[k] = 0.0;
     2092    // and do other array
     2093    double modification = 0.0;
     2094    CoinBigIndex j;
     2095    for (j = startPositive_[iColumn]; j < startNegative_[iColumn]; j++) {
     2096      int iRow = indices_[j];
     2097      modification += piWeight[iRow];
     2098    }
     2099    for (; j < startPositive_[iColumn + 1]; j++) {
     2100      int iRow = indices_[j];
     2101      modification -= piWeight[iRow];
     2102    }
     2103    double thisWeight = weights[iColumn];
     2104    double pivotSquared = pivot * pivot;
     2105    thisWeight += pivotSquared * devex + pivot * modification;
     2106    if (thisWeight < DEVEX_TRY_NORM) {
     2107      if (referenceIn < 0.0) {
     2108        // steepest
     2109        thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
     2110      } else {
     2111        // exact
     2112        thisWeight = referenceIn * pivotSquared;
     2113        if (reference(iColumn))
     2114          thisWeight += 1.0;
     2115        thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
     2116      }
     2117    }
     2118    weights[iColumn] = thisWeight;
     2119  }
    21502120}
    21512121/* Set the dimensions of the matrix. In effect, append new empty
     
    21542124   MUST be at least as large as the current ones otherwise an exception
    21552125   is thrown. */
    2156 void
    2157 ClpPlusMinusOneMatrix::setDimensions(int newnumrows, int newnumcols)
    2158 {
    2159      if (newnumrows < 0)
    2160           newnumrows = numberRows_;
    2161      if (newnumrows < numberRows_)
    2162           throw CoinError("Bad new rownum (less than current)",
    2163                           "setDimensions", "CoinPackedMatrix");
     2126void ClpPlusMinusOneMatrix::setDimensions(int newnumrows, int newnumcols)
     2127{
     2128  if (newnumrows < 0)
     2129    newnumrows = numberRows_;
     2130  if (newnumrows < numberRows_)
     2131    throw CoinError("Bad new rownum (less than current)",
     2132      "setDimensions", "CoinPackedMatrix");
    21642133
    2165      if (newnumcols < 0)
    2166           newnumcols = numberColumns_;
    2167      if (newnumcols < numberColumns_)
    2168           throw CoinError("Bad new colnum (less than current)",
    2169                           "setDimensions", "CoinPackedMatrix");
     2134  if (newnumcols < 0)
     2135    newnumcols = numberColumns_;
     2136  if (newnumcols < numberColumns_)
     2137    throw CoinError("Bad new colnum (less than current)",
     2138      "setDimensions", "CoinPackedMatrix");
    21702139
    2171      int number = 0;
    2172      int length = 0;
    2173      if (columnOrdered_) {
    2174           length = numberColumns_;
    2175           numberColumns_ = newnumcols;
    2176           number = numberColumns_;
     2140  int number = 0;
     2141  int length = 0;
     2142  if (columnOrdered_) {
     2143    length = numberColumns_;
     2144    numberColumns_ = newnumcols;
     2145    number = numberColumns_;
    21772146
    2178      } else {
    2179           length = numberRows_;
    2180           numberRows_ = newnumrows;
    2181           number = numberRows_;
    2182      }
    2183      if (number > length) {
    2184           CoinBigIndex * temp;
    2185           int i;
    2186           CoinBigIndex end = startPositive_[length];
    2187           temp = new CoinBigIndex [number+1];
    2188           CoinMemcpyN(startPositive_, (length + 1), temp);
    2189           delete [] startPositive_;
    2190           for (i = length + 1; i < number + 1; i++)
    2191                temp[i] = end;
    2192           startPositive_ = temp;
    2193           temp = new CoinBigIndex [number];
    2194           CoinMemcpyN(startNegative_, length, temp);
    2195           delete [] startNegative_;
    2196           for (i = length; i < number; i++)
    2197                temp[i] = end;
    2198           startNegative_ = temp;
    2199      }
     2147  } else {
     2148    length = numberRows_;
     2149    numberRows_ = newnumrows;
     2150    number = numberRows_;
     2151  }
     2152  if (number > length) {
     2153    CoinBigIndex *temp;
     2154    int i;
     2155    CoinBigIndex end = startPositive_[length];
     2156    temp = new CoinBigIndex[number + 1];
     2157    CoinMemcpyN(startPositive_, (length + 1), temp);
     2158    delete[] startPositive_;
     2159    for (i = length + 1; i < number + 1; i++)
     2160      temp[i] = end;
     2161    startPositive_ = temp;
     2162    temp = new CoinBigIndex[number];
     2163    CoinMemcpyN(startNegative_, length, temp);
     2164    delete[] startNegative_;
     2165    for (i = length; i < number; i++)
     2166      temp[i] = end;
     2167    startNegative_ = temp;
     2168  }
    22002169}
    22012170#ifndef SLIM_CLP
     
    22042173   number of columns-1/rows-1 (if numberOther>0) or duplicates
    22052174   If 0 then rows, 1 if columns */
    2206 int
    2207 ClpPlusMinusOneMatrix::appendMatrix(int number, int type,
    2208                                     const CoinBigIndex * starts, const int * index,
    2209                                     const double * element, int /*numberOther*/)
    2210 {
    2211      int numberErrors = 0;
    2212      // say if all +1
     2175int ClpPlusMinusOneMatrix::appendMatrix(int number, int type,
     2176  const CoinBigIndex *starts, const int *index,
     2177  const double *element, int /*numberOther*/)
     2178{
     2179  int numberErrors = 0;
     2180  // say if all +1
    22132181#ifdef CLP_PLUS_ONE_MATRIX
    2214      otherFlags_ = 0;
    2215 #endif
    2216      // make into CoinPackedVector
    2217      CoinPackedVectorBase ** vectors =
    2218           new CoinPackedVectorBase * [number];
    2219      int iVector;
    2220      for (iVector = 0; iVector < number; iVector++) {
    2221           CoinBigIndex iStart = starts[iVector];
    2222           vectors[iVector] =
    2223             new CoinPackedVector(static_cast<int>(starts[iVector+1] - iStart),
    2224                                     index + iStart, element + iStart);
    2225      }
    2226      if (type == 0) {
    2227           // rows
    2228           appendRows(number, vectors);
    2229      } else {
    2230           // columns
    2231           appendCols(number, vectors);
    2232      }
    2233      for (iVector = 0; iVector < number; iVector++)
    2234           delete vectors[iVector];
    2235      delete [] vectors;
    2236      return numberErrors;
     2182  otherFlags_ = 0;
     2183#endif
     2184  // make into CoinPackedVector
     2185  CoinPackedVectorBase **vectors = new CoinPackedVectorBase *[number];
     2186  int iVector;
     2187  for (iVector = 0; iVector < number; iVector++) {
     2188    CoinBigIndex iStart = starts[iVector];
     2189    vectors[iVector] = new CoinPackedVector(static_cast< int >(starts[iVector + 1] - iStart),
     2190      index + iStart, element + iStart);
     2191  }
     2192  if (type == 0) {
     2193    // rows
     2194    appendRows(number, vectors);
     2195  } else {
     2196    // columns
     2197    appendCols(number, vectors);
     2198  }
     2199  for (iVector = 0; iVector < number; iVector++)
     2200    delete vectors[iVector];
     2201  delete[] vectors;
     2202  return numberErrors;
    22372203}
    22382204#endif
     
    22452211// Default Constructor
    22462212//-------------------------------------------------------------------
    2247 ClpPoolMatrix::ClpPoolMatrix ()
     2213ClpPoolMatrix::ClpPoolMatrix()
    22482214  : ClpMatrixBase()
    22492215{
     
    22622228// Copy constructor
    22632229//-------------------------------------------------------------------
    2264 ClpPoolMatrix::ClpPoolMatrix (const ClpPoolMatrix & rhs)
     2230ClpPoolMatrix::ClpPoolMatrix(const ClpPoolMatrix &rhs)
    22652231  : ClpMatrixBase(rhs)
    22662232{
     
    22752241  numberDifferent_ = rhs.numberDifferent_;
    22762242  if (numberColumns_) {
    2277     columnStart_ = CoinCopyOfArray(rhs.columnStart_,numberColumns_+1);
     2243    columnStart_ = CoinCopyOfArray(rhs.columnStart_, numberColumns_ + 1);
    22782244    CoinBigIndex numberElements = columnStart_[numberColumns_];
    2279     stuff_ = CoinCopyOfArray(rhs.stuff_,numberElements);
    2280     elements_ = CoinCopyOfArray(rhs.elements_,numberDifferent_);
     2245    stuff_ = CoinCopyOfArray(rhs.stuff_, numberElements);
     2246    elements_ = CoinCopyOfArray(rhs.elements_, numberDifferent_);
    22812247  }
    22822248}
    22832249static const int mmult[] = {
    2284   262139, 259459, 256889, 254291, 251701, 249133, 246709, 244247};
     2250  262139, 259459, 256889, 254291, 251701, 249133, 246709, 244247
     2251};
    22852252inline int hashit(double value)
    22862253{
    2287   const unsigned char * chars = reinterpret_cast<unsigned char *>(&value);
     2254  const unsigned char *chars = reinterpret_cast< unsigned char * >(&value);
    22882255  int n = 0;
    2289   for (int j=0;j<8;j++)
    2290     n += mmult[j]*chars[j];
    2291   n = n % (1<<CLP_POOL_SIZE);
     2256  for (int j = 0; j < 8; j++)
     2257    n += mmult[j] * chars[j];
     2258  n = n % (1 << CLP_POOL_SIZE);
    22922259  return n;
    22932260}
    2294 ClpPoolMatrix::ClpPoolMatrix (const CoinPackedMatrix & rhs)
     2261ClpPoolMatrix::ClpPoolMatrix(const CoinPackedMatrix &rhs)
    22952262  : ClpMatrixBase()
    22962263{
     
    22992266  lengths_ = NULL;
    23002267  numberDifferent_ = 0;
    2301   assert (rhs.isColOrdered());
     2268  assert(rhs.isColOrdered());
    23022269  // get matrix data pointers
    2303   const int * row = rhs.getIndices();
    2304   const CoinBigIndex * columnStart = rhs.getVectorStarts();
    2305   const int * columnLength = rhs.getVectorLengths();
    2306   const double * elementByColumn = rhs.getElements();
     2270  const int *row = rhs.getIndices();
     2271  const CoinBigIndex *columnStart = rhs.getVectorStarts();
     2272  const int *columnLength = rhs.getVectorLengths();
     2273  const double *elementByColumn = rhs.getElements();
    23072274  numberColumns_ = rhs.getNumCols();
    23082275  numberRows_ = rhs.getNumRows();
    2309   assert (numberRows_<(1<<CLP_POOL_MATRIX));
     2276  assert(numberRows_ < (1 << CLP_POOL_MATRIX));
    23102277  elements_ = NULL;
    2311   columnStart_ = new CoinBigIndex [numberColumns_+1];
    2312   stuff_ = new poolInfo [rhs.getNumElements()];
    2313   int maxPool = 1<<CLP_POOL_SIZE;
    2314   double * tempDifferent  = new double [maxPool];
     2278  columnStart_ = new CoinBigIndex[numberColumns_ + 1];
     2279  stuff_ = new poolInfo[rhs.getNumElements()];
     2280  int maxPool = 1 << CLP_POOL_SIZE;
     2281  double *tempDifferent = new double[maxPool];
    23152282#define HASH 2
    23162283#if HASH > 1
     
    23192286    int index, next;
    23202287  } CoinHashLink;
    2321   CoinHashLink * hashThis = new CoinHashLink[4*maxPool];
    2322   for ( int i = 0; i < 4*maxPool; i++ ) {
     2288  CoinHashLink *hashThis = new CoinHashLink[4 * maxPool];
     2289  for (int i = 0; i < 4 * maxPool; i++) {
    23232290    hashThis[i].index = -1;
    23242291    hashThis[i].next = -1;
    23252292  }
    23262293#endif
    2327   int numberElements=0;
    2328   int hashDifferent=0;
     2294  int numberElements = 0;
     2295  int hashDifferent = 0;
    23292296  for (int iColumn = 0; iColumn < numberColumns_; iColumn++) {
    23302297    CoinBigIndex k;
    2331     columnStart_[iColumn]=numberElements;
     2298    columnStart_[iColumn] = numberElements;
    23322299    for (k = columnStart[iColumn]; k < columnStart[iColumn] + columnLength[iColumn];
    2333         k++) {
    2334       int iRow=row[k];
     2300        k++) {
     2301      int iRow = row[k];
    23352302      double value = elementByColumn[k];
    23362303      poolInfo stuff;
    2337       stuff.row_=iRow;
    2338 #if HASH==1 || HASH ==3
     2304      stuff.row_ = iRow;
     2305#if HASH == 1 || HASH == 3
    23392306      int j;
    2340       for (j=0;j<numberDifferent_;j++) {
    2341         if (value==tempDifferent[j])
    2342           break;
    2343       }
    2344 #endif
    2345 #if HASH==2 || HASH ==3
     2307      for (j = 0; j < numberDifferent_; j++) {
     2308        if (value == tempDifferent[j])
     2309          break;
     2310      }
     2311#endif
     2312#if HASH == 2 || HASH == 3
    23462313      int ipos = hashit(value);
    23472314      int j1;
    2348       while ( true ) {
    2349         j1 = hashThis[ipos].index;
    2350         if (j1 == -1) {
    2351           hashThis[ipos].index=numberDifferent_;
    2352           j1=numberDifferent_;
    2353           break;
    2354         } else if (value==tempDifferent[j1]) {
    2355           break;
    2356         } else {
    2357           int k = hashThis[ipos].next;
    2358           if ( k == -1 ) {
    2359             j1 = numberDifferent_;
    2360             while ( true ) {
    2361               ++hashDifferent;
    2362               if ( hashThis[hashDifferent].index == -1 ) {
    2363                 break;
    2364               }
    2365             }
    2366             hashThis[ipos].next = hashDifferent;
    2367             hashThis[hashDifferent].index = numberDifferent_;
    2368             break;
    2369           } else {
    2370             ipos = k;
    2371             /* nothing worked - try it again */
    2372           }
    2373         }
    2374       }
    2375 #if HASH ==2
    2376       int j=j1;
    2377 #endif
    2378 #endif
    2379 #if HASH ==3
    2380       assert (j==j1);
    2381 #endif
    2382       if (j==numberDifferent_) {
    2383         if (j==maxPool)
    2384           break; // too many
    2385         tempDifferent[j]=value;
    2386         numberDifferent_++;
    2387       }
    2388       stuff.pool_=j;
    2389       stuff_[numberElements++]=stuff;
    2390     }
    2391   }
    2392   columnStart_[numberColumns_]=numberElements;
    2393 #if HASH>1
    2394   delete [] hashThis;
    2395 #endif
    2396   elements_ = new double [numberDifferent_];
    2397   memcpy(elements_,tempDifferent,numberDifferent_*sizeof(double));
    2398   delete [] tempDifferent;
    2399   if (numberDifferent_==maxPool) {
    2400     delete [] stuff_;
    2401     delete [] elements_;
    2402     delete [] columnStart_;
     2315      while (true) {
     2316        j1 = hashThis[ipos].index;
     2317        if (j1 == -1) {
     2318          hashThis[ipos].index = numberDifferent_;
     2319          j1 = numberDifferent_;
     2320          break;
     2321        } else if (value == tempDifferent[j1]) {
     2322          break;
     2323        } else {
     2324          int k = hashThis[ipos].next;
     2325          if (k == -1) {
     2326            j1 = numberDifferent_;
     2327            while (true) {
     2328              ++hashDifferent;
     2329              if (hashThis[hashDifferent].index == -1) {
     2330                break;
     2331              }
     2332            }
     2333            hashThis[ipos].next = hashDifferent;
     2334            hashThis[hashDifferent].index = numberDifferent_;
     2335            break;
     2336          } else {
     2337            ipos = k;
     2338            /* nothing worked - try it again */
     2339          }
     2340        }
     2341      }
     2342#if HASH == 2
     2343      int j = j1;
     2344#endif
     2345#endif
     2346#if HASH == 3
     2347      assert(j == j1);
     2348#endif
     2349      if (j == numberDifferent_) {
     2350        if (j == maxPool)
     2351          break; // too many
     2352        tempDifferent[j] = value;
     2353        numberDifferent_++;
     2354      }
     2355      stuff.pool_ = j;
     2356      stuff_[numberElements++] = stuff;
     2357    }
     2358  }
     2359  columnStart_[numberColumns_] = numberElements;
     2360#if HASH > 1
     2361  delete[] hashThis;
     2362#endif
     2363  elements_ = new double[numberDifferent_];
     2364  memcpy(elements_, tempDifferent, numberDifferent_ * sizeof(double));
     2365  delete[] tempDifferent;
     2366  if (numberDifferent_ == maxPool) {
     2367    delete[] stuff_;
     2368    delete[] elements_;
     2369    delete[] columnStart_;
    24032370    elements_ = NULL;
    24042371    columnStart_ = NULL;
     
    24132380// Destructor
    24142381//-------------------------------------------------------------------
    2415 ClpPoolMatrix::~ClpPoolMatrix ()
     2382ClpPoolMatrix::~ClpPoolMatrix()
    24162383{
    24172384  delete matrix_;
    2418   delete [] elements_;
    2419   delete [] columnStart_;
    2420   delete [] lengths_;
    2421   delete [] stuff_;
     2385  delete[] elements_;
     2386  delete[] columnStart_;
     2387  delete[] lengths_;
     2388  delete[] stuff_;
    24222389}
    24232390
     
    24262393//-------------------------------------------------------------------
    24272394ClpPoolMatrix &
    2428 ClpPoolMatrix::operator=(const ClpPoolMatrix& rhs)
     2395ClpPoolMatrix::operator=(const ClpPoolMatrix &rhs)
    24292396{
    24302397  if (this != &rhs) {
    24312398    ClpMatrixBase::operator=(rhs);
    24322399    delete matrix_;
    2433     delete [] elements_;
    2434     delete [] columnStart_;
    2435     delete [] lengths_;
    2436     delete [] stuff_;
     2400    delete[] elements_;
     2401    delete[] columnStart_;
     2402    delete[] lengths_;
     2403    delete[] stuff_;
    24372404    matrix_ = NULL;
    24382405    lengths_ = NULL;
     
    24402407    numberColumns_ = rhs.numberColumns_;
    24412408    if (numberColumns_) {
    2442       columnStart_ = CoinCopyOfArray(rhs.columnStart_,numberColumns_+1);
     2409      columnStart_ = CoinCopyOfArray(rhs.columnStart_, numberColumns_ + 1);
    24432410      CoinBigIndex numberElements = columnStart_[numberColumns_];
    2444       stuff_ = CoinCopyOfArray(rhs.stuff_,numberElements);
    2445       elements_ = CoinCopyOfArray(rhs.elements_,numberDifferent_);
     2411      stuff_ = CoinCopyOfArray(rhs.stuff_, numberElements);
     2412      elements_ = CoinCopyOfArray(rhs.elements_, numberDifferent_);
    24462413    }
    24472414  }
     
    24492416}
    24502417// Constructor from arrays - handing over ownership
    2451 ClpPoolMatrix::ClpPoolMatrix(int numberColumns,CoinBigIndex * columnStart,
    2452                              poolInfo * stuff, double * elements)
     2418ClpPoolMatrix::ClpPoolMatrix(int numberColumns, CoinBigIndex *columnStart,
     2419  poolInfo *stuff, double *elements)
    24532420  : ClpMatrixBase()
    24542421{
     
    24642431  for (int iColumn = 0; iColumn < numberColumns_; iColumn++) {
    24652432    CoinBigIndex k;
    2466     for (k = columnStart_[iColumn]; k < columnStart_[iColumn+1];
    2467         k++) {
    2468       int iRow=stuff_[k].row_;
    2469       int iPool=stuff_[k].pool_;
    2470       numberDifferent_ = CoinMax(numberDifferent_,iPool);
    2471       numberRows_ = CoinMax(numberRows_,iRow);
     2433    for (k = columnStart_[iColumn]; k < columnStart_[iColumn + 1];
     2434        k++) {
     2435      int iRow = stuff_[k].row_;
     2436      int iPool = stuff_[k].pool_;
     2437      numberDifferent_ = CoinMax(numberDifferent_, iPool);
     2438      numberRows_ = CoinMax(numberRows_, iRow);
    24722439    }
    24732440  }
     
    24792446// Clone
    24802447//-------------------------------------------------------------------
    2481 ClpMatrixBase * ClpPoolMatrix::clone() const
     2448ClpMatrixBase *ClpPoolMatrix::clone() const
    24822449{
    24832450  return new ClpPoolMatrix(*this);
     
    24862453   and order is as given */
    24872454ClpMatrixBase *
    2488 ClpPoolMatrix::subsetClone (int numberRows, const int * whichRows,
    2489                             int numberColumns,
    2490                             const int * whichColumns) const
     2455ClpPoolMatrix::subsetClone(int numberRows, const int *whichRows,
     2456  int numberColumns,
     2457  const int *whichColumns) const
    24912458{
    24922459  return new ClpPoolMatrix(*this, numberRows, whichRows,
    2493                            numberColumns, whichColumns);
     2460    numberColumns, whichColumns);
    24942461}
    24952462/* Subset constructor (without gaps).  Duplicates are allowed
    24962463   and order is as given */
    2497 ClpPoolMatrix::ClpPoolMatrix (
    2498                               const ClpPoolMatrix & rhs,
    2499                               int numberRows, const int * whichRow,
    2500                               int numberColumns, const int * whichColumn)
     2464ClpPoolMatrix::ClpPoolMatrix(
     2465  const ClpPoolMatrix &rhs,
     2466  int numberRows, const int *whichRow,
     2467  int numberColumns, const int *whichColumn)
    25012468  : ClpMatrixBase(rhs)
    25022469{
     
    25102477  numberColumns_ = numberColumns;
    25112478  numberDifferent_ = 0;
    2512   if (!numberRows_||!numberColumns_)
     2479  if (!numberRows_ || !numberColumns_)
    25132480    return;
    25142481  numberDifferent_ = rhs.numberDifferent_;
    2515   elements_ = CoinCopyOfArray(rhs.elements_,numberDifferent_);
     2482  elements_ = CoinCopyOfArray(rhs.elements_, numberDifferent_);
    25162483  int numberRowsOther = rhs.numberRows_;
    25172484  int numberColumnsOther = rhs.numberColumns_;
    2518   int * newRow = new int [numberRowsOther+numberRows_];
    2519   int * duplicateRow = newRow+numberRowsOther;
     2485  int *newRow = new int[numberRowsOther + numberRows_];
     2486  int *duplicateRow = newRow + numberRowsOther;
    25202487  for (int iRow = 0; iRow < numberRowsOther; iRow++)
    25212488    newRow[iRow] = -1;
     
    25252492    duplicateRow[iRow] = -1;
    25262493    int kRow = whichRow[iRow];
    2527     if (kRow >= 0  && kRow < numberRowsOther) {
     2494    if (kRow >= 0 && kRow < numberRowsOther) {
    25282495      if (newRow[kRow] < 0) {
    2529         // first time
    2530         newRow[kRow] = iRow;
     2496        // first time
     2497        newRow[kRow] = iRow;
    25312498      } else {
    2532         // duplicate
    2533         int lastRow = newRow[kRow];
    2534         newRow[kRow] = iRow;
    2535         duplicateRow[iRow] = lastRow;
     2499        // duplicate
     2500        int lastRow = newRow[kRow];
     2501        newRow[kRow] = iRow;
     2502        duplicateRow[iRow] = lastRow;
    25362503      }
    25372504    } else {
     
    25402507    }
    25412508  }
    2542  
     2509
    25432510  if (numberBad)
    25442511    throw CoinError("bad row entries",
    2545                     "subset constructor", "ClpPoolMatrix");
     2512      "subset constructor", "ClpPoolMatrix");
    25462513  // now get size and check columns
    25472514  CoinBigIndex size = 0;
    25482515  numberBad = 0;
    2549   const CoinBigIndex * starts = rhs.columnStart_;
    2550   const poolInfo * stuff = rhs.stuff_;
     2516  const CoinBigIndex *starts = rhs.columnStart_;
     2517  const poolInfo *stuff = rhs.stuff_;
    25512518  for (int iColumn = 0; iColumn < numberColumns_; iColumn++) {
    25522519    int kColumn = whichColumn[iColumn];
    2553     if (kColumn >= 0  && kColumn < numberColumnsOther) {
     2520    if (kColumn >= 0 && kColumn < numberColumnsOther) {
    25542521      CoinBigIndex i;
    2555       for (i = starts[kColumn]; i < starts[kColumn+1]; i++) {
    2556         int kRow = stuff[i].row_;
    2557         kRow = newRow[kRow];
    2558         while (kRow >= 0) {
    2559           size++;
    2560           kRow = duplicateRow[kRow];
    2561         }
     2522      for (i = starts[kColumn]; i < starts[kColumn + 1]; i++) {
     2523        int kRow = stuff[i].row_;
     2524        kRow = newRow[kRow];
     2525        while (kRow >= 0) {
     2526          size++;
     2527          kRow = duplicateRow[kRow];
     2528        }
    25622529      }
    25632530    } else {
     
    25682535  if (numberBad)
    25692536    throw CoinError("bad major entries",
    2570                     "subset constructor", "ClpPoolMatrix");
     2537      "subset constructor", "ClpPoolMatrix");
    25712538  // now create arrays
    25722539  stuff_ = new poolInfo[size];
    2573   columnStart_ = new CoinBigIndex [numberColumns_+1];
    2574   size=0;
     2540  columnStart_ = new CoinBigIndex[numberColumns_ + 1];
     2541  size = 0;
    25752542  for (int iColumn = 0; iColumn < numberColumns_; iColumn++) {
    25762543    int kColumn = whichColumn[iColumn];
    2577     columnStart_[iColumn]=size;
     2544    columnStart_[iColumn] = size;
    25782545    CoinBigIndex i;
    2579     for (i = starts[kColumn]; i < starts[kColumn+1]; i++) {
     2546    for (i = starts[kColumn]; i < starts[kColumn + 1]; i++) {
    25802547      int kRow = stuff[i].row_;
    25812548      int iPool = stuff[i].pool_;
    25822549      kRow = newRow[kRow];
    25832550      while (kRow >= 0) {
    2584         stuff_[size].row_=kRow;
    2585         stuff_[size++].pool_=iPool;
    2586         kRow = duplicateRow[kRow];
    2587       }
    2588     }
    2589   }
    2590   columnStart_[numberColumns_]=size;
    2591   delete [] newRow;
     2551        stuff_[size].row_ = kRow;
     2552        stuff_[size++].pool_ = iPool;
     2553        kRow = duplicateRow[kRow];
     2554      }
     2555    }
     2556  }
     2557  columnStart_[numberColumns_] = size;
     2558  delete[] newRow;
    25922559}
    25932560// Create matrix_
    2594 ClpPackedMatrix * 
     2561ClpPackedMatrix *
    25952562ClpPoolMatrix::createMatrix() const
    25962563{
    25972564  if (!matrix_) {
    25982565    CoinBigIndex numberElements = columnStart_[numberColumns_];
    2599     double * elements = new double [numberElements];
    2600     int * rows = new int [numberElements];
    2601     bool lengthsExist = lengths_!=NULL;
    2602     numberElements=0;
     2566    double *elements = new double[numberElements];
     2567    int *rows = new int[numberElements];
     2568    bool lengthsExist = lengths_ != NULL;
     2569    numberElements = 0;
    26032570    if (!lengthsExist)
    26042571      lengths_ = new int[numberColumns_];
    26052572    for (int i = 0; i < numberColumns_; i++) {
    2606       lengths_[i]=static_cast<int>(columnStart_[i+1]-columnStart_[i]);
    2607       for (CoinBigIndex j = columnStart_[i]; j < columnStart_[i+1]; j++) {
    2608         elements[numberElements] = elements_[stuff_[j].pool_];
    2609         rows[numberElements++] = stuff_[j].row_;
    2610       }
    2611     }
    2612     CoinPackedMatrix * matrix = new CoinPackedMatrix(true, numberRows_, numberColumns_,
    2613                                                       numberElements,
    2614                                                       elements, rows,
    2615                                                       columnStart_, lengths_);
    2616     delete [] elements;
    2617     delete [] rows;
     2573      lengths_[i] = static_cast< int >(columnStart_[i + 1] - columnStart_[i]);
     2574      for (CoinBigIndex j = columnStart_[i]; j < columnStart_[i + 1]; j++) {
     2575        elements[numberElements] = elements_[stuff_[j].pool_];
     2576        rows[numberElements++] = stuff_[j].row_;
     2577      }
     2578    }
     2579    CoinPackedMatrix *matrix = new CoinPackedMatrix(true, numberRows_, numberColumns_,
     2580      numberElements,
     2581      elements, rows,
     2582      columnStart_, lengths_);
     2583    delete[] elements;
     2584    delete[] rows;
    26182585    matrix_ = new ClpPackedMatrix(matrix);
    26192586  }
    26202587  return matrix_;
    26212588}
    2622 
    26232589
    26242590/* Returns a new matrix in reverse order without gaps */
     
    26292595#define BIG_MATRIX
    26302596#ifndef BIG_MATRIX
    2631   printf("XcreateMatrix at file %s line %d\n",__FILE__,__LINE__);
     2597  printf("XcreateMatrix at file %s line %d\n", __FILE__, __LINE__);
    26322598  createMatrix();
    26332599  return matrix_->reverseOrderedCopy();
     
    26382604#undef ABOCA_LITE
    26392605//unscaled versions
    2640 void
    2641 ClpPoolMatrix::times(double scalar,
    2642                      const double * COIN_RESTRICT x, double * COIN_RESTRICT y) const
     2606void ClpPoolMatrix::times(double scalar,
     2607  const double *COIN_RESTRICT x, double *COIN_RESTRICT y) const
    26432608{
    26442609#if 0
     
    26512616    if (value) {
    26522617      CoinBigIndex start = columnStart_[iColumn];
    2653       CoinBigIndex end = columnStart_[iColumn+1];
     2618      CoinBigIndex end = columnStart_[iColumn + 1];
    26542619      value *= scalar;
    26552620      for (j = start; j < end; j++) {
    2656         int iRow = stuff_[j].row_;
    2657         y[iRow] += value * elements_[stuff_[j].pool_];
    2658       }
    2659     }
    2660   }
    2661 #endif
    2662 }
    2663 void
    2664 ClpPoolMatrix::transposeTimes(double scalar,
    2665                               const double * COIN_RESTRICT x, double * COIN_RESTRICT y) const
     2621        int iRow = stuff_[j].row_;
     2622        y[iRow] += value * elements_[stuff_[j].pool_];
     2623      }
     2624    }
     2625  }
     2626#endif
     2627}
     2628void ClpPoolMatrix::transposeTimes(double scalar,
     2629  const double *COIN_RESTRICT x, double *COIN_RESTRICT y) const
    26662630{
    26672631#if 0
     
    26742638    for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    26752639      CoinBigIndex j;
    2676       CoinBigIndex next = columnStart_[iColumn+1];
     2640      CoinBigIndex next = columnStart_[iColumn + 1];
    26772641      double value = 0.0;
    26782642      // scaled
    26792643      for (j = start; j < next; j++) {
    2680         int jRow = stuff_[j].row_;
    2681         value += x[jRow] * elements_[stuff_[j].pool_];
     2644        int jRow = stuff_[j].row_;
     2645        value += x[jRow] * elements_[stuff_[j].pool_];
    26822646      }
    26832647      start = next;
     
    26872651    for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    26882652      CoinBigIndex j;
    2689       CoinBigIndex next = columnStart_[iColumn+1];
     2653      CoinBigIndex next = columnStart_[iColumn + 1];
    26902654      double value = 0.0;
    26912655      // scaled
    26922656      for (j = start; j < next; j++) {
    2693         int jRow = stuff_[j].row_;
    2694         value += x[jRow] * elements_[stuff_[j].pool_];
     2657        int jRow = stuff_[j].row_;
     2658        value += x[jRow] * elements_[stuff_[j].pool_];
    26952659      }
    26962660      start = next;
     
    27002664#endif
    27012665}
    2702 void
    2703 ClpPoolMatrix::times(double scalar,
    2704                      const double * COIN_RESTRICT x, double * COIN_RESTRICT y,
    2705                      const double * rowScale,
    2706                      const double * columnScale) const
    2707 {
    2708   printf("createMatrix at file %s line %d\n",__FILE__,__LINE__);
    2709   createMatrix()->times(scalar,x,y);
    2710 }
    2711 void
    2712 ClpPoolMatrix::transposeTimes( double scalar,
    2713                                const double * COIN_RESTRICT x, double * COIN_RESTRICT y,
    2714                                const double * rowScale,
    2715                                const double * columnScale,
    2716                                double * spare) const
     2666void ClpPoolMatrix::times(double scalar,
     2667  const double *COIN_RESTRICT x, double *COIN_RESTRICT y,
     2668  const double *rowScale,
     2669  const double *columnScale) const
     2670{
     2671  printf("createMatrix at file %s line %d\n", __FILE__, __LINE__);
     2672  createMatrix()->times(scalar, x, y);
     2673}
     2674void ClpPoolMatrix::transposeTimes(double scalar,
     2675  const double *COIN_RESTRICT x, double *COIN_RESTRICT y,
     2676  const double *rowScale,
     2677  const double *columnScale,
     2678  double *spare) const
    27172679{
    27182680#if 0
     
    27252687      CoinBigIndex start = columnStart_[0];
    27262688      if (scalar == -1.0) {
    2727         for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    2728           CoinBigIndex j;
    2729           CoinBigIndex next = columnStart_[iColumn+1];
    2730           double value = 0.0;
    2731           // scaled
    2732           for (j = start; j < next; j++) {
    2733             int jRow = stuff_[j].row_;
    2734             value += x[jRow] * elements_[stuff_[j].pool_] * rowScale[jRow];
    2735           }
    2736           start = next;
    2737           y[iColumn] -= value * columnScale[iColumn];
    2738         }
     2689        for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     2690          CoinBigIndex j;
     2691          CoinBigIndex next = columnStart_[iColumn + 1];
     2692          double value = 0.0;
     2693          // scaled
     2694          for (j = start; j < next; j++) {
     2695            int jRow = stuff_[j].row_;
     2696            value += x[jRow] * elements_[stuff_[j].pool_] * rowScale[jRow];
     2697          }
     2698          start = next;
     2699          y[iColumn] -= value * columnScale[iColumn];
     2700        }
    27392701      } else {
    2740         for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    2741           CoinBigIndex j;
    2742           CoinBigIndex next = columnStart_[iColumn+1];
    2743           double value = 0.0;
    2744           // scaled
    2745           for (j = start; j < next; j++) {
    2746             int jRow = stuff_[j].row_;
    2747             value += x[jRow] * elements_[stuff_[j].pool_] * rowScale[jRow];
    2748           }
    2749           start = next;
    2750           y[iColumn] += value * scalar * columnScale[iColumn];
    2751         }
     2702        for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
     2703          CoinBigIndex j;
     2704          CoinBigIndex next = columnStart_[iColumn + 1];
     2705          double value = 0.0;
     2706          // scaled
     2707          for (j = start; j < next; j++) {
     2708            int jRow = stuff_[j].row_;
     2709            value += x[jRow] * elements_[stuff_[j].pool_] * rowScale[jRow];
     2710          }
     2711          start = next;
     2712          y[iColumn] += value * scalar * columnScale[iColumn];
     2713        }
    27522714      }
    27532715    } else {
     
    27562718      int numberRows = matrix_->getNumRows();
    27572719      for (iRow = 0; iRow < numberRows; iRow++) {
    2758         double value = x[iRow];
    2759         if (value)
    2760           spare[iRow] = value * rowScale[iRow];
    2761         else
    2762           spare[iRow] = 0.0;
     2720        double value = x[iRow];
     2721        if (value)
     2722          spare[iRow] = value * rowScale[iRow];
     2723        else
     2724          spare[iRow] = 0.0;
    27632725      }
    27642726      CoinBigIndex start = columnStart_[0];
    27652727      for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    2766         CoinBigIndex j;
    2767         CoinBigIndex next = columnStart_[iColumn+1];
    2768         double value = 0.0;
    2769         // scaled
    2770         for (j = start; j < next; j++) {
    2771           int jRow = stuff_[j].row_;
    2772           value += spare[jRow] * elements_[stuff_[j].pool_];
    2773         }
    2774         start = next;
    2775         y[iColumn] += value * scalar * columnScale[iColumn];
     2728        CoinBigIndex j;
     2729        CoinBigIndex next = columnStart_[iColumn + 1];
     2730        double value = 0.0;
     2731        // scaled
     2732        for (j = start; j < next; j++) {
     2733          int jRow = stuff_[j].row_;
     2734          value += spare[jRow] * elements_[stuff_[j].pool_];
     2735        }
     2736        start = next;
     2737        y[iColumn] += value * scalar * columnScale[iColumn];
    27762738      }
    27772739    }
     
    27832745/* Return <code>x * A + y</code> in <code>z</code>.
    27842746   Squashes small elements and knows about ClpSimplex */
    2785 void
    2786 ClpPoolMatrix::transposeTimes(const ClpSimplex * model, double scalar,
    2787                               const CoinIndexedVector * rowArray,
    2788                               CoinIndexedVector * y,
    2789                               CoinIndexedVector * columnArray) const
     2747void ClpPoolMatrix::transposeTimes(const ClpSimplex *model, double scalar,
     2748  const CoinIndexedVector *rowArray,
     2749  CoinIndexedVector *y,
     2750  CoinIndexedVector *columnArray) const
    27902751{
    27912752#if 0
     
    27972758  columnArray->clear();
    27982759  // do by column
    2799   double * COIN_RESTRICT pi = rowArray->denseVector();
     2760  double *COIN_RESTRICT pi = rowArray->denseVector();
    28002761  int numberNonZero = 0;
    2801   int * COIN_RESTRICT index = columnArray->getIndices();
    2802   double * COIN_RESTRICT array = columnArray->denseVector();
     2762  int *COIN_RESTRICT index = columnArray->getIndices();
     2763  double *COIN_RESTRICT array = columnArray->denseVector();
    28032764  int numberInRowArray = rowArray->getNumElements();
    28042765  // maybe I need one in OsiSimplex
    28052766  double zeroTolerance = model->zeroTolerance();
    28062767  bool packed = rowArray->packedMode();
    2807   assert (packed);
     2768  assert(packed);
    28082769  // do by column
    2809   const double * COIN_RESTRICT rowScale = model->rowScale();
    2810   assert (!y->getNumElements());
     2770  const double *COIN_RESTRICT rowScale = model->rowScale();
     2771  assert(!y->getNumElements());
    28112772  // need to expand pi into y
    28122773  assert(y->capacity() >= model->numberRows());
    2813   double * COIN_RESTRICT piOld = pi;
     2774  double *COIN_RESTRICT piOld = pi;
    28142775  pi = y->denseVector();
    2815   const int * COIN_RESTRICT whichRow = rowArray->getIndices();
     2776  const int *COIN_RESTRICT whichRow = rowArray->getIndices();
    28162777  int i;
    28172778  // later combine phases
     
    28322793    for (iColumn = 0; iColumn < numberColumns_ - 1; iColumn++) {
    28332794      CoinBigIndex start = end;
    2834       end = columnStart_[iColumn+2];
     2795      end = columnStart_[iColumn + 2];
    28352796      if (fabs(value) > zeroTolerance) {
    2836         array[numberNonZero] = value;
    2837         index[numberNonZero++] = iColumn;
     2797        array[numberNonZero] = value;
     2798        index[numberNonZero++] = iColumn;
    28382799      }
    28392800      value = 0.0;
    28402801      for (j = start; j < end; j++) {
    2841         int iRow = stuff_[j].row_;
    2842         value += pi[iRow] * elements_[stuff_[j].pool_];
     2802        int iRow = stuff_[j].row_;
     2803        value += pi[iRow] * elements_[stuff_[j].pool_];
    28432804      }
    28442805    }
     
    28542815      pi[iRow] = scalar * piOld[i] * rowScale[iRow];
    28552816    }
    2856     const double * columnScale = model->columnScale();
     2817    const double *columnScale = model->columnScale();
    28572818    double value = 0.0;
    28582819    double scale = columnScale[0];
     
    28672828      value *= scale;
    28682829      CoinBigIndex start = end;
    2869       scale = columnScale[iColumn+1];
    2870       end = columnStart_[iColumn+2];
     2830      scale = columnScale[iColumn + 1];
     2831      end = columnStart_[iColumn + 2];
    28712832      if (fabs(value) > zeroTolerance) {
    2872         array[numberNonZero] = value;
    2873         index[numberNonZero++] = iColumn;
     2833        array[numberNonZero] = value;
     2834        index[numberNonZero++] = iColumn;
    28742835      }
    28752836      value = 0.0;
    28762837      for (j = start; j < end; j++) {
    2877         int iRow = stuff_[j].row_;
    2878         value += pi[iRow] * elements_[stuff_[j].pool_];
     2838        int iRow = stuff_[j].row_;
     2839        value += pi[iRow] * elements_[stuff_[j].pool_];
    28792840      }
    28802841    }
     
    29022863/* Return <code>x * A + y</code> in <code>z</code>.
    29032864   Squashes small elements and knows about ClpSimplex */
    2904 void
    2905 ClpPoolMatrix::transposeTimesByRow(const ClpSimplex * model, double scalar,
    2906                                    const CoinIndexedVector * rowArray,
    2907                                    CoinIndexedVector * y,
    2908                                    CoinIndexedVector * columnArray) const
    2909 {
    2910   printf("createMatrix at file %s line %d\n",__FILE__,__LINE__);
    2911   createMatrix()->transposeTimesByRow(model,scalar,rowArray,y,columnArray);
     2865void ClpPoolMatrix::transposeTimesByRow(const ClpSimplex *model, double scalar,
     2866  const CoinIndexedVector *rowArray,
     2867  CoinIndexedVector *y,
     2868  CoinIndexedVector *columnArray) const
     2869{
     2870  printf("createMatrix at file %s line %d\n", __FILE__, __LINE__);
     2871  createMatrix()->transposeTimesByRow(model, scalar, rowArray, y, columnArray);
    29122872}
    29132873/* Return <code>x *A in <code>z</code> but
    29142874   just for indices in y. */
    2915 void
    2916 ClpPoolMatrix::subsetTransposeTimes(const ClpSimplex * model,
    2917                                     const CoinIndexedVector * rowArray,
    2918                                     const CoinIndexedVector * y,
    2919                                     CoinIndexedVector * columnArray) const
     2875void ClpPoolMatrix::subsetTransposeTimes(const ClpSimplex *model,
     2876  const CoinIndexedVector *rowArray,
     2877  const CoinIndexedVector *y,
     2878  CoinIndexedVector *columnArray) const
    29202879{
    29212880#if 0
     
    29242883#else
    29252884  columnArray->clear();
    2926   double * COIN_RESTRICT pi = rowArray->denseVector();
    2927   double * COIN_RESTRICT array = columnArray->denseVector();
     2885  double *COIN_RESTRICT pi = rowArray->denseVector();
     2886  double *COIN_RESTRICT array = columnArray->denseVector();
    29282887  int jColumn;
    29292888  // get matrix data pointers
    2930   const double * COIN_RESTRICT rowScale = model->rowScale();
     2889  const double *COIN_RESTRICT rowScale = model->rowScale();
    29312890  int numberToDo = y->getNumElements();
    2932   const int * COIN_RESTRICT which = y->getIndices();
    2933   assert (!rowArray->packedMode());
     2891  const int *COIN_RESTRICT which = y->getIndices();
     2892  assert(!rowArray->packedMode());
    29342893  columnArray->setPacked();
    29352894  if (!rowScale) {
    29362895    for (jColumn = 0; jColumn < numberToDo; jColumn++) {
    2937       int iColumn=which[jColumn];
     2896      int iColumn = which[jColumn];
    29382897      double value = 0.0;
    2939       CoinBigIndex start=columnStart_[iColumn];
    2940       CoinBigIndex end=columnStart_[iColumn+1];
     2898      CoinBigIndex start = columnStart_[iColumn];
     2899      CoinBigIndex end = columnStart_[iColumn + 1];
    29412900      array[jColumn] = value;
    29422901      CoinBigIndex j;
    29432902      for (j = start; j < end; j++) {
    2944         int iRow = stuff_[j].row_;
    2945         value += pi[iRow] * elements_[stuff_[j].pool_];
     2903        int iRow = stuff_[j].row_;
     2904        value += pi[iRow] * elements_[stuff_[j].pool_];
    29462905      }
    29472906      array[jColumn] = value;
     
    29492908  } else {
    29502909    // scaled
    2951     const double * columnScale = model->columnScale();
     2910    const double *columnScale = model->columnScale();
    29522911    for (jColumn = 0; jColumn < numberToDo; jColumn++) {
    2953       int iColumn=which[jColumn];
     2912      int iColumn = which[jColumn];
    29542913      double value = 0.0;
    2955       CoinBigIndex start=columnStart_[iColumn];
    2956       CoinBigIndex end=columnStart_[iColumn+1];
     2914      CoinBigIndex start = columnStart_[iColumn];
     2915      CoinBigIndex end = columnStart_[iColumn + 1];
    29572916      CoinBigIndex j;
    29582917      for (j = start; j < end; j++) {
    2959         int iRow = stuff_[j].row_;
    2960         value += pi[iRow] * elements_[stuff_[j].pool_] * rowScale[iRow];
    2961       }
    2962       array[jColumn] = value*columnScale[iColumn];
     2918        int iRow = stuff_[j].row_;
     2919        value += pi[iRow] * elements_[stuff_[j].pool_] * rowScale[iRow];
     2920      }
     2921      array[jColumn] = value * columnScale[iColumn];
    29632922    }
    29642923  }
     
    29662925}
    29672926/// returns number of elements in column part of basis,
    2968 int
    2969 ClpPoolMatrix::countBasis(const int * whichColumn,
    2970                           int & numberColumnBasic)
     2927int ClpPoolMatrix::countBasis(const int *whichColumn,
     2928  int &numberColumnBasic)
    29712929{
    29722930  int i;
     
    29742932  for (i = 0; i < numberColumnBasic; i++) {
    29752933    int iColumn = whichColumn[i];
    2976     numberElements += columnStart_[iColumn+1] - columnStart_[iColumn];
     2934    numberElements += columnStart_[iColumn + 1] - columnStart_[iColumn];
    29772935  }
    29782936#if COIN_BIG_INDEX
    2979   if (numberElements>COIN_INT_MAX) {
     2937  if (numberElements > COIN_INT_MAX) {
    29802938    printf("Factorization too large\n");
    29812939    abort();
    29822940  }
    29832941#endif
    2984   return static_cast<int>(numberElements);
    2985 }
    2986 void
    2987 ClpPoolMatrix::fillBasis(ClpSimplex * ,
    2988                          const int * whichColumn,
    2989                          int & numberColumnBasic,
    2990                          int * indexRowU, int * start,
    2991                          int * rowCount, int * columnCount,
    2992                          CoinFactorizationDouble * elementU)
     2942  return static_cast< int >(numberElements);
     2943}
     2944void ClpPoolMatrix::fillBasis(ClpSimplex *,
     2945  const int *whichColumn,
     2946  int &numberColumnBasic,
     2947  int *indexRowU, int *start,
     2948  int *rowCount, int *columnCount,
     2949  CoinFactorizationDouble *elementU)
    29932950{
    29942951  CoinBigIndex numberElements = start[0];
     
    29962953    int iColumn = whichColumn[i];
    29972954    CoinBigIndex j = columnStart_[iColumn];
    2998     for (; j < columnStart_[iColumn+1]; j++) {
     2955    for (; j < columnStart_[iColumn + 1]; j++) {
    29992956      int iRow = stuff_[j].row_;
    30002957      indexRowU[numberElements] = iRow;
    30012958      rowCount[iRow]++;
    3002       elementU[numberElements++] = elements_[stuff_[j].pool_];;
    3003     }
    3004     start[i+1] = static_cast<int>(numberElements);
    3005     columnCount[i] = static_cast<int>(numberElements - start[i]);
     2959      elementU[numberElements++] = elements_[stuff_[j].pool_];
     2960      ;
     2961    }
     2962    start[i + 1] = static_cast< int >(numberElements);
     2963    columnCount[i] = static_cast< int >(numberElements - start[i]);
    30062964  }
    30072965#if COIN_BIG_INDEX
    3008   if (numberElements>COIN_INT_MAX)
     2966  if (numberElements > COIN_INT_MAX)
    30092967    abort();
    30102968#endif
     
    30122970/* Unpacks a column into an CoinIndexedvector
    30132971 */
    3014 void
    3015 ClpPoolMatrix::unpack(const ClpSimplex * ,
    3016                       CoinIndexedVector * rowArray,
    3017                       int iColumn) const
     2972void ClpPoolMatrix::unpack(const ClpSimplex *,
     2973  CoinIndexedVector *rowArray,
     2974  int iColumn) const
    30182975{
    30192976  CoinBigIndex j = columnStart_[iColumn];
    3020   for (; j < columnStart_[iColumn+1]; j++) {
     2977  for (; j < columnStart_[iColumn + 1]; j++) {
    30212978    int iRow = stuff_[j].row_;
    30222979    rowArray->add(iRow, elements_[stuff_[j].pool_]);
     
    30272984Note that model is NOT const.  Bounds and objective could
    30282985be modified if doing column generation (just for this variable) */
    3029 void
    3030 ClpPoolMatrix::unpackPacked(ClpSimplex * ,
    3031                             CoinIndexedVector * rowArray,
    3032                             int iColumn) const
    3033 {
    3034   int * COIN_RESTRICT index = rowArray->getIndices();
    3035   double * COIN_RESTRICT array = rowArray->denseVector();
     2986void ClpPoolMatrix::unpackPacked(ClpSimplex *,
     2987  CoinIndexedVector *rowArray,
     2988  int iColumn) const
     2989{
     2990  int *COIN_RESTRICT index = rowArray->getIndices();
     2991  double *COIN_RESTRICT array = rowArray->denseVector();
    30362992  int number = 0;
    30372993  CoinBigIndex j = columnStart_[iColumn];
    3038   for (; j < columnStart_[iColumn+1]; j++) {
     2994  for (; j < columnStart_[iColumn + 1]; j++) {
    30392995    int iRow = stuff_[j].row_;
    30402996    array[number] = elements_[stuff_[j].pool_];
     
    30463002/* Adds multiple of a column into an CoinIndexedvector
    30473003   You can use quickAdd to add to vector */
    3048 void
    3049 ClpPoolMatrix::add(const ClpSimplex * , CoinIndexedVector * rowArray,
    3050                    int iColumn, double multiplier) const
     3004void ClpPoolMatrix::add(const ClpSimplex *, CoinIndexedVector *rowArray,
     3005  int iColumn, double multiplier) const
    30513006{
    30523007  CoinBigIndex j = columnStart_[iColumn];
    3053   for (; j < columnStart_[iColumn+1]; j++) {
     3008  for (; j < columnStart_[iColumn + 1]; j++) {
    30543009    int iRow = stuff_[j].row_;
    3055     rowArray->quickAdd(iRow, multiplier*elements_[stuff_[j].pool_]);
     3010    rowArray->quickAdd(iRow, multiplier * elements_[stuff_[j].pool_]);
    30563011  }
    30573012}
    30583013/* Adds multiple of a column into an array */
    3059 void
    3060 ClpPoolMatrix::add(const ClpSimplex * , double * array,
    3061                    int iColumn, double multiplier) const
     3014void ClpPoolMatrix::add(const ClpSimplex *, double *array,
     3015  int iColumn, double multiplier) const
    30623016{
    30633017  CoinBigIndex j = columnStart_[iColumn];
    3064   for (; j < columnStart_[iColumn+1]; j++) {
     3018  for (; j < columnStart_[iColumn + 1]; j++) {
    30653019    int iRow = stuff_[j].row_;
    3066     array[iRow] += multiplier*elements_[stuff_[j].pool_];
     3020    array[iRow] += multiplier * elements_[stuff_[j].pool_];
    30673021  }
    30683022}
     
    30823036ClpPoolMatrix::getElements() const
    30833037{
    3084   printf("createMatrix at file %s line %d\n",__FILE__,__LINE__);
     3038  printf("createMatrix at file %s line %d\n", __FILE__, __LINE__);
    30853039  return createMatrix()->getElements();
    30863040}
     
    30943048ClpPoolMatrix::getIndices() const
    30953049{
    3096   printf("createMatrix at file %s line %d\n",__FILE__,__LINE__);
     3050  printf("createMatrix at file %s line %d\n", __FILE__, __LINE__);
    30973051  return createMatrix()->getIndices();
    30983052}
     
    31003054ClpPoolMatrix::getVectorStarts() const
    31013055{
    3102   printf("createMatrix at file %s line %d\n",__FILE__,__LINE__);
     3056  printf("createMatrix at file %s line %d\n", __FILE__, __LINE__);
    31033057  return createMatrix()->getVectorStarts();
    31043058}
     
    31133067  if (!lengths_) {
    31143068    lengths_ = new int[numberColumns_];
    3115     for (int i = 0; i < numberColumns_; i++) 
    3116       lengths_[i]=static_cast<int>(columnStart_[i+1]-columnStart_[i]);
     3069    for (int i = 0; i < numberColumns_; i++)
     3070      lengths_[i] = static_cast< int >(columnStart_[i + 1] - columnStart_[i]);
    31173071  }
    31183072  return lengths_;
     
    31203074}
    31213075/* The length of a major-dimension vector. */
    3122 int
    3123 ClpPoolMatrix::getVectorLength(int index) const
    3124 {
    3125   return static_cast<int>(columnStart_[index+1]-columnStart_[index]);
     3076int ClpPoolMatrix::getVectorLength(int index) const
     3077{
     3078  return static_cast< int >(columnStart_[index + 1] - columnStart_[index]);
    31263079}
    31273080/* Delete the columns whose indices are listed in <code>indDel</code>. */
    3128 void
    3129 ClpPoolMatrix::deleteCols(const int numDel, const int * indDel)
     3081void ClpPoolMatrix::deleteCols(const int numDel, const int *indDel)
    31303082{
    31313083  int iColumn;
    3132   CoinBigIndex newSize = columnStart_[numberColumns_];;
     3084  CoinBigIndex newSize = columnStart_[numberColumns_];
     3085  ;
    31333086  int numberBad = 0;
    31343087  // Use array to make sure we can have duplicates
    3135   int * which = new int[numberColumns_];
     3088  int *which = new int[numberColumns_];
    31363089  memset(which, 0, numberColumns_ * sizeof(int));
    31373090  int nDuplicate = 0;
     
    31413094      numberBad++;
    31423095    } else {
    3143       newSize -= columnStart_[jColumn+1] - columnStart_[jColumn];
     3096      newSize -= columnStart_[jColumn + 1] - columnStart_[jColumn];
    31443097      if (which[jColumn])
    3145         nDuplicate++;
     3098        nDuplicate++;
    31463099      else
    3147         which[jColumn] = 1;
     3100        which[jColumn] = 1;
    31483101    }
    31493102  }
     
    31543107  delete matrix_;
    31553108  matrix_ = NULL;
    3156   CoinBigIndex * columnStart = new CoinBigIndex [newNumber+1];
    3157   poolInfo * stuff = new poolInfo [newSize];
     3109  CoinBigIndex *columnStart = new CoinBigIndex[newNumber + 1];
     3110  poolInfo *stuff = new poolInfo[newSize];
    31583111  newNumber = 0;
    31593112  newSize = 0;
     
    31633116      CoinBigIndex i;
    31643117      start = columnStart_[iColumn];
    3165       end = columnStart_[iColumn+1];
     3118      end = columnStart_[iColumn + 1];
    31663119      columnStart[newNumber] = newSize;
    31673120      for (i = start; i < end; i++)
    3168         stuff[newSize++] = stuff_[i];
     3121        stuff[newSize++] = stuff_[i];
    31693122    }
    31703123  }
    31713124  columnStart[newNumber] = newSize;
    3172   delete [] which;
    3173   delete [] columnStart_;
     3125  delete[] which;
     3126  delete[] columnStart_;
    31743127  columnStart_ = columnStart;
    3175   delete [] stuff_;
     3128  delete[] stuff_;
    31763129  stuff_ = stuff;
    31773130  numberColumns_ = newNumber;
    31783131}
    31793132/* Delete the rows whose indices are listed in <code>indDel</code>. */
    3180 void
    3181 ClpPoolMatrix::deleteRows(const int numDel, const int * indDel)
     3133void ClpPoolMatrix::deleteRows(const int numDel, const int *indDel)
    31823134{
    31833135  int iRow;
    31843136  int numberBad = 0;
    31853137  // Use array to make sure we can have duplicates
    3186   int * which = new int[numberRows_];
     3138  int *which = new int[numberRows_];
    31873139  memset(which, 0, numberRows_ * sizeof(int));
    31883140  int nDuplicate = 0;
     
    31933145    } else {
    31943146      if (which[jRow])
    3195         nDuplicate++;
     3147        nDuplicate++;
    31963148      else
    3197         which[jRow] = 1;
     3149        which[jRow] = 1;
    31983150    }
    31993151  }
     
    32123164  delete matrix_;
    32133165  matrix_ = NULL;
    3214   poolInfo * stuff = new poolInfo [newSize];
     3166  poolInfo *stuff = new poolInfo[newSize];
    32153167  newSize = 0;
    32163168  // redo which
    3217   int numberRows=0;
     3169  int numberRows = 0;
    32183170  for (iRow = 0; iRow < numberRows_; iRow++) {
    32193171    if (which[iRow]) {
    3220       which[iRow]=-1; // not wanted
     3172      which[iRow] = -1; // not wanted
    32213173    } else {
    3222       which[iRow]=numberRows;
     3174      which[iRow] = numberRows;
    32233175      numberRows++;
    32243176    }
     
    32293181    CoinBigIndex i;
    32303182    start = columnStart_[iColumn];
    3231     end = columnStart_[iColumn+1];
     3183    end = columnStart_[iColumn + 1];
    32323184    columnStart_[newNumber] = newSize;
    32333185    for (i = start; i < end; i++) {
    3234       poolInfo newStuff=stuff_[i];
     3186      poolInfo newStuff = stuff_[i];
    32353187      iRow = which[newStuff.row_];
    3236       if (iRow>=0) {
    3237         newStuff.row_=iRow;
    3238         stuff[newSize++] = newStuff;
     3188      if (iRow >= 0) {
     3189        newStuff.row_ = iRow;
     3190        stuff[newSize++] = newStuff;
    32393191      }
    32403192    }
    32413193  }
    32423194  columnStart_[numberColumns_] = newSize;
    3243   delete [] which;
    3244   delete [] stuff_;
     3195  delete[] which;
     3196  delete[] stuff_;
    32453197  stuff_ = stuff;
    32463198  numberRows_ = newNumber;
    32473199}
    3248 bool
    3249 ClpPoolMatrix::isColOrdered() const
     3200bool ClpPoolMatrix::isColOrdered() const
    32503201{
    32513202  return true;
     
    32633214   Largest refers to largest absolute value.
    32643215*/
    3265 void
    3266 ClpPoolMatrix::rangeOfElements(double & smallestNegative, double & largestNegative,
    3267                                double & smallestPositive, double & largestPositive)
     3216void ClpPoolMatrix::rangeOfElements(double &smallestNegative, double &largestNegative,
     3217  double &smallestPositive, double &largestPositive)
    32683218{
    32693219  smallestNegative = -COIN_DBL_MAX;
     
    32833233}
    32843234// Says whether it can do partial pricing
    3285 bool
    3286 ClpPoolMatrix::canDoPartialPricing() const
     3235bool ClpPoolMatrix::canDoPartialPricing() const
    32873236{
    32883237  return true;
    32893238}
    32903239// Partial pricing
    3291 void
    3292 ClpPoolMatrix::partialPricing(ClpSimplex * model, double startFraction, double endFraction,
    3293                               int & bestSequence, int & numberWanted)
    3294 {
    3295   printf("createMatrix at file %s line %d\n",__FILE__,__LINE__);
    3296   createMatrix()->partialPricing(model,startFraction,
    3297                                  endFraction,bestSequence,
    3298                                  numberWanted);
     3240void ClpPoolMatrix::partialPricing(ClpSimplex *model, double startFraction, double endFraction,
     3241  int &bestSequence, int &numberWanted)
     3242{
     3243  printf("createMatrix at file %s line %d\n", __FILE__, __LINE__);
     3244  createMatrix()->partialPricing(model, startFraction,
     3245    endFraction, bestSequence,
     3246    numberWanted);
    32993247}
    33003248// Allow any parts of a created CoinMatrix to be deleted
    3301 void
    3302 ClpPoolMatrix::releasePackedMatrix() const
     3249void ClpPoolMatrix::releasePackedMatrix() const
    33033250{
    33043251  delete matrix_;
     
    33063253}
    33073254// These have to match ClpPrimalColumnSteepest version
    3308 #define reference(i)  (((reference[i>>5]>>(i&31))&1)!=0)
     3255#define reference(i) (((reference[i >> 5] >> (i & 31)) & 1) != 0)
    33093256/* Updates two arrays for steepest and does devex weights
    33103257   Returns nonzero if updates reduced cost and infeas -
    33113258   new infeas in dj1  - This does not at present*/
    3312 int
    3313 ClpPoolMatrix::transposeTimes2(const ClpSimplex * model,
    3314                                const CoinIndexedVector * pi1, CoinIndexedVector * dj1,
    3315                                const CoinIndexedVector * pi2,
    3316                                CoinIndexedVector * spare,
    3317                                double * infeas, double * reducedCost,
    3318                                double referenceIn, double devex,
    3319                                // Array for exact devex to say what is in reference framework
    3320                                unsigned int * COIN_RESTRICT reference,
    3321                                double * COIN_RESTRICT weights, double scaleFactor)
     3259int ClpPoolMatrix::transposeTimes2(const ClpSimplex *model,
     3260  const CoinIndexedVector *pi1, CoinIndexedVector *dj1,
     3261  const CoinIndexedVector *pi2,
     3262  CoinIndexedVector *spare,
     3263  double *infeas, double *reducedCost,
     3264  double referenceIn, double devex,
     3265  // Array for exact devex to say what is in reference framework
     3266  unsigned int *COIN_RESTRICT reference,
     3267  double *COIN_RESTRICT weights, double scaleFactor)
    33223268{
    33233269#if 0
     
    33273273                                         reference,weights,scaleFactor);
    33283274#else
    3329   int returnCode=0;
     3275  int returnCode = 0;
    33303276  // put row of tableau in dj1
    3331   double * COIN_RESTRICT pi = pi1->denseVector();
     3277  double *COIN_RESTRICT pi = pi1->denseVector();
    33323278  int numberNonZero = 0;
    3333   int * COIN_RESTRICT index = dj1->getIndices();
    3334   double * COIN_RESTRICT array = dj1->denseVector();
     3279  int *COIN_RESTRICT index = dj1->getIndices();
     3280  double *COIN_RESTRICT array = dj1->denseVector();
    33353281  int numberInRowArray = pi1->getNumElements();
    33363282  double zeroTolerance = model->zeroTolerance();
     
    33403286  double error = CoinMin(1.0e-2, model->largestDualError());
    33413287  // allow tolerance at least slightly bigger than standard
    3342   dualTolerance = dualTolerance +  error;
     3288  dualTolerance = dualTolerance + error;
    33433289  bool packed = pi1->packedMode();
    3344   assert (packed);
     3290  assert(packed);
    33453291  // do by column
    33463292  int iColumn;
    3347   const double * COIN_RESTRICT rowScale = model->rowScale();
    3348   assert (!spare->getNumElements());
    3349   assert (numberColumns_ > 0);
    3350   double * COIN_RESTRICT piWeight = pi2->denseVector();
    3351   assert (!pi2->packedMode());
     3293  const double *COIN_RESTRICT rowScale = model->rowScale();
     3294  assert(!spare->getNumElements());
     3295  assert(numberColumns_ > 0);
     3296  double *COIN_RESTRICT piWeight = pi2->denseVector();
     3297  assert(!pi2->packedMode());
    33523298  bool killDjs = (scaleFactor == 0.0);
    33533299  if (!scaleFactor)
     
    33553301  // need to expand pi into y
    33563302  assert(spare->capacity() >= model->numberRows());
    3357   double * COIN_RESTRICT piOld = pi;
     3303  double *COIN_RESTRICT piOld = pi;
    33583304  pi = spare->denseVector();
    3359   const int * COIN_RESTRICT whichRow = pi1->getIndices();
     3305  const int *COIN_RESTRICT whichRow = pi1->getIndices();
    33603306  int i;
    33613307  if (!rowScale) {
     
    33663312    }
    33673313    if (infeas)
    3368       returnCode=1;
     3314      returnCode = 1;
    33693315    CoinBigIndex j;
    33703316    CoinBigIndex end = columnStart_[0];
    33713317    for (iColumn = 0; iColumn < numberColumns_; iColumn++) {
    33723318      CoinBigIndex start = end;
    3373       end = columnStart_[iColumn+1];
     3319      end = columnStart_[iColumn + 1];
    33743320      ClpSimplex::Status status = model->getStatus(iColumn);
    3375       if (status == ClpSimplex::basic || status == ClpSimplex::isFixed) continue;
     3321      if (status == ClpSimplex::basic || status == ClpSimplex::isFixed)
     3322        continue;
    33763323      double value = 0.0;
    33773324      for (j = start; j < end; j++) {
    3378         int iRow = stuff_[j].row_;
    3379         value -= pi[iRow]*elements_[stuff_[j].pool_];
     3325        int iRow = stuff_[j].row_;
     3326        value -= pi[iRow] * elements_[stuff_[j].pool_];
    33803327      }
    33813328      if (fabs(value) > zeroTolerance) {
    3382         // and do other array
    3383         double modification = 0.0;
    3384         for (j = start; j < end; j++) {
    3385           int iRow = stuff_[j].row_;
    3386           modification += piWeight[iRow]*elements_[stuff_[j].pool_];
    3387         }
    3388         double thisWeight = weights[iColumn];
    3389         double pivot = value * scaleFactor;
    3390         double pivotSquared = pivot * pivot;
    3391         thisWeight += pivotSquared * devex + pivot * modification;
    3392         if (thisWeight < DEVEX_TRY_NORM) {
    3393           if (referenceIn < 0.0) {
    3394             // steepest
    3395             thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
    3396           } else {
    3397             // exact
    3398             thisWeight = referenceIn * pivotSquared;
    3399             if (reference(iColumn))
    3400               thisWeight += 1.0;
    3401             thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
    3402           }
    3403         }
    3404         weights[iColumn] = thisWeight;
    3405         if (!killDjs) {
    3406           value = reducedCost[iColumn]-value;
    3407           reducedCost[iColumn] = value;
    3408           // simplify status
    3409           ClpSimplex::Status status = model->getStatus(iColumn);
    3410          
    3411           switch(status) {
    3412            
    3413           case ClpSimplex::basic:
    3414           case ClpSimplex::isFixed:
    3415             break;
    3416           case ClpSimplex::isFree:
    3417           case ClpSimplex::superBasic:
    3418             if (fabs(value) > FREE_ACCEPT * dualTolerance) {
    3419               // we are going to bias towards free (but only if reasonable)
    3420               value *= FREE_BIAS;
    3421               value *= value;
    3422               // store square in list
    3423               if (infeas[iColumn]) {
    3424                 infeas[iColumn] = value; // already there
    3425               } else {
    3426                 array[numberNonZero]=value;
    3427                 index[numberNonZero++]=iColumn;
    3428               }
    3429             } else {
    3430               array[numberNonZero]=0.0;
    3431               index[numberNonZero++]=iColumn;
    3432             }
    3433             break;
    3434           case ClpSimplex::atUpperBound:
    3435             if (value > dualTolerance) {
    3436               value *= value;
    3437               // store square in list
    3438               if (infeas[iColumn]) {
    3439                 infeas[iColumn] = value; // already there
    3440               } else {
    3441                 array[numberNonZero]=value;
    3442                 index[numberNonZero++]=iColumn;
    3443               }
    3444             } else {
    3445               array[numberNonZero]=0.0;
    3446               index[numberNonZero++]=iColumn;
    3447             }
    3448             break;
    3449           case ClpSimplex::atLowerBound:
    3450             if (value < -dualTolerance) {
    3451               value *= value;
    3452               // store square in list
    3453               if (infeas[iColumn]) {
    3454                 infeas[iColumn] = value; // already there
    3455               } else {
    3456                 array[numberNonZero]=value;
    3457                 index[numberNonZero++]=iColumn;
    3458               }
    3459             } else {
    3460               array[numberNonZero]=0.0;
    3461               index[numberNonZero++]=iColumn;
    3462             }
    3463           }
    3464         }
     3329        // and do other array
     3330        double modification = 0.0;
     3331        for (j = start; j < end; j++) {
     3332          int iRow = stuff_[j].row_;
     3333          modification += piWeight[iRow] * elements_[stuff_[j].pool_];
     3334        }
     3335        double thisWeight = weights[iColumn];
     3336        double pivot = value * scaleFactor;
     3337        double pivotSquared = pivot * pivot;
     3338        thisWeight += pivotSquared * devex + pivot * modification;
     3339        if (thisWeight < DEVEX_TRY_NORM) {
     3340          if (referenceIn < 0.0) {
     3341            // steepest
     3342            thisWeight = CoinMax(DEVEX_TRY_NORM, DEVEX_ADD_ONE + pivotSquared);
     3343          } else {
     3344            // exact
     3345            thisWeight = referenceIn * pivotSquared;
     3346            if (reference(iColumn))
     3347              thisWeight += 1.0;
     3348            thisWeight = CoinMax(thisWeight, DEVEX_TRY_NORM);
     3349          }
     3350        }
     3351        weights[iColumn] = thisWeight;
     3352        if (!killDjs) {
     3353          value = reducedCost[iColumn] - value;
     3354          reducedCost[iColumn] = value;
     3355          // simplify status
     3356          ClpSimplex::Status status = model->getStatus(iColumn);
     3357
     3358          switch (status) {
     3359
     3360          case ClpSimplex::basic:
     3361          case ClpSimplex::isFixed:
     3362            break;
     3363          case ClpSimplex::isFree:
     3364          case ClpSimplex::superBasic:
     3365            if (fabs(value) > FREE_ACCEPT * dualTolerance) {
     3366              // we are going to bias towards free (but only if reasonable)
     3367              value *= FREE_BIAS;
     3368              value *= value;
     3369              // store square in list
     3370              if (infeas[iColumn]) {
     3371                infeas[iColumn] = value; // already there
     3372              } else {
     3373                array[numberNonZero] = value;
     3374                index[numberNonZero++] = iColumn;
     3375              }
     3376            } else {
     3377              array[numberNonZero] = 0.0;
     3378              index[numberNonZero++] = iColumn;
     3379            }
     3380            break;
     3381          case ClpSimplex::atUpperBound:
     3382            if (value > dualTolerance) {
     3383              value *= value;
     3384              // store square in list
     3385              if (infeas[iColumn]) {
     3386                infeas[iColumn] = value; // already there
     3387              } else {
     3388                array[numberNonZero] = value;
     3389                index[numberNonZero++] = iColumn;
     3390              }
     3391            } else {
     3392              array[numberNonZero] = 0.0;
     3393              index[numberNonZero++] = iColumn;
     3394            }
     3395            break;
     3396          case ClpSimplex::atLowerBound:
     3397            if (value < -dualTolerance) {
     3398              value *= value;
     3399              // store square in list
     3400              if (infeas[iColumn]) {
     3401