CoinFactorization.hpp
Go to the documentation of this file.
1 /* $Id: CoinFactorization.hpp 1277 2010-04-21 20:19:32Z forrest $ */
2 // Copyright (C) 2002, International Business Machines
3 // Corporation and others. All Rights Reserved.
4 
5 /*
6  Authors
7 
8  John Forrest
9 
10  */
11 #ifndef CoinFactorization_H
12 #define CoinFactorization_H
13 //#define COIN_ONE_ETA_COPY 100
14 
15 #include <iostream>
16 #include <string>
17 #include <cassert>
18 #include <cstdio>
19 #include "CoinFinite.hpp"
20 #include "CoinIndexedVector.hpp"
21 class CoinPackedMatrix;
48  friend void CoinFactorizationUnitTest( const std::string & mpsDir );
49 
50 public:
51 
54 
57  CoinFactorization ( const CoinFactorization &other);
58 
62  void almostDestructor();
64  void show_self ( ) const;
66  int saveFactorization (const char * file ) const;
70  int restoreFactorization (const char * file , bool factor=false) ;
72  void sort ( ) const;
76 
86  int factorize ( const CoinPackedMatrix & matrix,
87  int rowIsBasic[], int columnIsBasic[] ,
88  double areaFactor = 0.0 );
99  int factorize ( int numberRows,
100  int numberColumns,
102  CoinBigIndex maximumL,
103  CoinBigIndex maximumU,
104  const int indicesRow[],
105  const int indicesColumn[], const double elements[] ,
106  int permutation[],
107  double areaFactor = 0.0);
112  int factorizePart1 ( int numberRows,
113  int numberColumns,
114  CoinBigIndex estimateNumberElements,
115  int * indicesRow[],
116  int * indicesColumn[],
117  CoinFactorizationDouble * elements[],
118  double areaFactor = 0.0);
125  int factorizePart2 (int permutation[],int exactNumberElements);
127  double conditionNumber() const;
128 
130 
133 
134  inline int status ( ) const {
135  return status_;
136  }
138  inline void setStatus ( int value)
139  { status_=value; }
141  inline int pivots ( ) const {
142  return numberPivots_;
143  }
145  inline void setPivots ( int value )
146  { numberPivots_=value; }
148  inline int *permute ( ) const {
149  return permute_.array();
150  }
152  inline int *pivotColumn ( ) const {
153  return pivotColumn_.array();
154  }
157  return pivotRegion_.array();
158  }
160  inline int *permuteBack ( ) const {
161  return permuteBack_.array();
162  }
165  inline int *pivotColumnBack ( ) const {
166  //return firstCount_.array();
167  return pivotColumnBack_.array();
168  }
170  inline CoinBigIndex * startRowL() const
171  { return startRowL_.array();}
172 
174  inline CoinBigIndex * startColumnL() const
175  { return startColumnL_.array();}
176 
178  inline int * indexColumnL() const
179  { return indexColumnL_.array();}
180 
182  inline int * indexRowL() const
183  { return indexRowL_.array();}
184 
187  { return elementByRowL_.array();}
188 
190  inline int numberRowsExtra ( ) const {
191  return numberRowsExtra_;
192  }
194  inline void setNumberRows(int value)
195  { numberRows_ = value; }
197  inline int numberRows ( ) const {
198  return numberRows_;
199  }
201  inline CoinBigIndex numberL() const
202  { return numberL_;}
203 
205  inline CoinBigIndex baseL() const
206  { return baseL_;}
208  inline int maximumRowsExtra ( ) const {
209  return maximumRowsExtra_;
210  }
212  inline int numberColumns ( ) const {
213  return numberColumns_;
214  }
216  inline int numberElements ( ) const {
217  return totalElements_;
218  }
220  inline int numberForrestTomlin ( ) const {
222  }
224  inline int numberGoodColumns ( ) const {
225  return numberGoodU_;
226  }
228  inline double areaFactor ( ) const {
229  return areaFactor_;
230  }
231  inline void areaFactor ( double value ) {
232  areaFactor_=value;
233  }
235  double adjustedAreaFactor() const;
237  inline void relaxAccuracyCheck(double value)
238  { relaxCheck_ = value;}
239  inline double getAccuracyCheck() const
240  { return relaxCheck_;}
242  inline int messageLevel ( ) const {
243  return messageLevel_ ;
244  }
245  void messageLevel ( int value );
247  inline int maximumPivots ( ) const {
248  return maximumPivots_ ;
249  }
250  void maximumPivots ( int value );
251 
253  inline int denseThreshold() const
254  { return denseThreshold_;}
256  inline void setDenseThreshold(int value)
257  { denseThreshold_ = value;}
259  inline double pivotTolerance ( ) const {
260  return pivotTolerance_ ;
261  }
262  void pivotTolerance ( double value );
264  inline double zeroTolerance ( ) const {
265  return zeroTolerance_ ;
266  }
267  void zeroTolerance ( double value );
268 #ifndef COIN_FAST_CODE
269 
270  inline double slackValue ( ) const {
271  return slackValue_ ;
272  }
273  void slackValue ( double value );
274 #endif
275 
276  double maximumCoefficient() const;
278  inline bool forrestTomlin() const
279  { return doForrestTomlin_;}
280  inline void setForrestTomlin(bool value)
281  { doForrestTomlin_=value;}
283  inline bool spaceForForrestTomlin() const
284  {
286  CoinBigIndex space = lengthAreaU_ - ( start + numberRowsExtra_ );
287  return (space>=0)&&doForrestTomlin_;
288  }
290 
293 
295  inline int numberDense() const
296  { return numberDense_;}
297 
299  inline CoinBigIndex numberElementsU ( ) const {
300  return lengthU_;
301  }
303  inline void setNumberElementsU(CoinBigIndex value)
304  { lengthU_ = value; }
306  inline CoinBigIndex lengthAreaU ( ) const {
307  return lengthAreaU_;
308  }
310  inline CoinBigIndex numberElementsL ( ) const {
311  return lengthL_;
312  }
314  inline CoinBigIndex lengthAreaL ( ) const {
315  return lengthAreaL_;
316  }
318  inline CoinBigIndex numberElementsR ( ) const {
319  return lengthR_;
320  }
323  { return numberCompressions_;}
325  inline int * numberInRow() const
326  { return numberInRow_.array();}
328  inline int * numberInColumn() const
329  { return numberInColumn_.array();}
332  { return elementU_.array();}
334  inline int * indexRowU() const
335  { return indexRowU_.array();}
337  inline CoinBigIndex * startColumnU() const
338  { return startColumnU_.array();}
340  inline int maximumColumnsExtra()
341  { return maximumColumnsExtra_;}
345  inline int biasLU() const
346  { return biasLU_;}
347  inline void setBiasLU(int value)
348  { biasLU_=value;}
354  inline int persistenceFlag() const
355  { return persistenceFlag_;}
356  void setPersistenceFlag(int value);
358 
361 
369  int replaceColumn ( CoinIndexedVector * regionSparse,
370  int pivotRow,
371  double pivotCheck ,
372  bool checkBeforeModifying=false,
373  double acceptablePivot=1.0e-8);
378  void replaceColumnU ( CoinIndexedVector * regionSparse,
379  CoinBigIndex * deleted,
380  int internalPivotRow);
382 
392  int updateColumnFT ( CoinIndexedVector * regionSparse,
393  CoinIndexedVector * regionSparse2);
396  int updateColumn ( CoinIndexedVector * regionSparse,
397  CoinIndexedVector * regionSparse2,
398  bool noPermute=false) const;
404  int updateTwoColumnsFT ( CoinIndexedVector * regionSparse1,
405  CoinIndexedVector * regionSparse2,
406  CoinIndexedVector * regionSparse3,
407  bool noPermuteRegion3=false) ;
412  int updateColumnTranspose ( CoinIndexedVector * regionSparse,
413  CoinIndexedVector * regionSparse2) const;
415  void goSparse();
417  inline int sparseThreshold ( ) const
418  { return sparseThreshold_;}
420  void sparseThreshold ( int value );
422 
423 
427 
428  inline void clearArrays()
429  { gutsOfDestructor();}
431 
434 
437  int add ( CoinBigIndex numberElements,
438  int indicesRow[],
439  int indicesColumn[], double elements[] );
440 
443  int addColumn ( CoinBigIndex numberElements,
444  int indicesRow[], double elements[] );
445 
448  int addRow ( CoinBigIndex numberElements,
449  int indicesColumn[], double elements[] );
450 
452  int deleteColumn ( int Row );
454  int deleteRow ( int Row );
455 
459  int replaceRow ( int whichRow, int numberElements,
460  const int indicesColumn[], const double elements[] );
462  void emptyRows(int numberToEmpty, const int which[]);
464 
465 
466  void checkSparse();
468  inline bool collectStatistics() const
469  { return collectStatistics_;}
471  inline void setCollectStatistics(bool onOff) const
472  { collectStatistics_ = onOff;}
474  void gutsOfDestructor(int type=1);
476  void gutsOfInitialize(int type);
477  void gutsOfCopy(const CoinFactorization &other);
478 
480  void resetStatistics();
481 
482 
484 
486 
487  void getAreas ( int numberRows,
488  int numberColumns,
489  CoinBigIndex maximumL,
490  CoinBigIndex maximumU );
491 
494  void preProcess ( int state,
495  int possibleDuplicates = -1 );
497  int factor ( );
498 protected:
501  int factorSparse ( );
504  int factorSparseSmall ( );
507  int factorSparseLarge ( );
510  int factorDense ( );
511 
513  bool pivotOneOtherRow ( int pivotRow,
514  int pivotColumn );
516  bool pivotRowSingleton ( int pivotRow,
517  int pivotColumn );
519  bool pivotColumnSingleton ( int pivotRow,
520  int pivotColumn );
521 
526  bool getColumnSpace ( int iColumn,
527  int extraNeeded );
528 
531  bool reorderU();
535  bool getColumnSpaceIterateR ( int iColumn, double value,
536  int iRow);
542  CoinBigIndex getColumnSpaceIterate ( int iColumn, double value,
543  int iRow);
547  bool getRowSpace ( int iRow, int extraNeeded );
548 
552  bool getRowSpaceIterate ( int iRow,
553  int extraNeeded );
555  void checkConsistency ( );
557  inline void addLink ( int index, int count ) {
558  int *nextCount = nextCount_.array();
559  int *firstCount = firstCount_.array();
560  int *lastCount = lastCount_.array();
561  int next = firstCount[count];
562  lastCount[index] = -2 - count;
563  if ( next < 0 ) {
564  //first with that count
565  firstCount[count] = index;
566  nextCount[index] = -1;
567  } else {
568  firstCount[count] = index;
569  nextCount[index] = next;
570  lastCount[next] = index;
571  }}
573  inline void deleteLink ( int index ) {
574  int *nextCount = nextCount_.array();
575  int *firstCount = firstCount_.array();
576  int *lastCount = lastCount_.array();
577  int next = nextCount[index];
578  int last = lastCount[index];
579  if ( last >= 0 ) {
580  nextCount[last] = next;
581  } else {
582  int count = -last - 2;
583 
584  firstCount[count] = next;
585  }
586  if ( next >= 0 ) {
587  lastCount[next] = last;
588  }
589  nextCount[index] = -2;
590  lastCount[index] = -2;
591  return;
592  }
594  void separateLinks(int count,bool rowsFirst);
596  void cleanup ( );
597 
599  void updateColumnL ( CoinIndexedVector * region, int * indexIn ) const;
601  void updateColumnLDensish ( CoinIndexedVector * region, int * indexIn ) const;
603  void updateColumnLSparse ( CoinIndexedVector * region, int * indexIn ) const;
605  void updateColumnLSparsish ( CoinIndexedVector * region, int * indexIn ) const;
606 
608  void updateColumnR ( CoinIndexedVector * region ) const;
611  void updateColumnRFT ( CoinIndexedVector * region, int * indexIn );
612 
614  void updateColumnU ( CoinIndexedVector * region, int * indexIn) const;
615 
617  void updateColumnUSparse ( CoinIndexedVector * regionSparse,
618  int * indexIn) const;
620  void updateColumnUSparsish ( CoinIndexedVector * regionSparse,
621  int * indexIn) const;
623  int updateColumnUDensish ( double * COIN_RESTRICT region,
624  int * COIN_RESTRICT regionIndex) const;
627  int & numberNonZero1,
628  double * COIN_RESTRICT region1,
629  int * COIN_RESTRICT index1,
630  int & numberNonZero2,
631  double * COIN_RESTRICT region2,
632  int * COIN_RESTRICT index2) const;
634  void updateColumnPFI ( CoinIndexedVector * regionSparse) const;
636  void permuteBack ( CoinIndexedVector * regionSparse,
637  CoinIndexedVector * outVector) const;
638 
640  void updateColumnTransposePFI ( CoinIndexedVector * region) const;
644  int smallestIndex) const;
648  int smallestIndex) const;
652  int smallestIndex) const;
655  void updateColumnTransposeUSparse ( CoinIndexedVector * region) const;
659  int smallestIndex) const;
660 
662  void updateColumnTransposeR ( CoinIndexedVector * region ) const;
664  void updateColumnTransposeRDensish ( CoinIndexedVector * region ) const;
666  void updateColumnTransposeRSparse ( CoinIndexedVector * region ) const;
667 
669  void updateColumnTransposeL ( CoinIndexedVector * region ) const;
671  void updateColumnTransposeLDensish ( CoinIndexedVector * region ) const;
673  void updateColumnTransposeLByRow ( CoinIndexedVector * region ) const;
675  void updateColumnTransposeLSparsish ( CoinIndexedVector * region ) const;
677  void updateColumnTransposeLSparse ( CoinIndexedVector * region ) const;
678 public:
683  int replaceColumnPFI ( CoinIndexedVector * regionSparse,
684  int pivotRow, double alpha);
685 protected:
688  int checkPivot(double saveFromU, double oldPivot) const;
689  /********************************* START LARGE TEMPLATE ********/
690 #ifdef INT_IS_8
691 #define COINFACTORIZATION_BITS_PER_INT 64
692 #define COINFACTORIZATION_SHIFT_PER_INT 6
693 #define COINFACTORIZATION_MASK_PER_INT 0x3f
694 #else
695 #define COINFACTORIZATION_BITS_PER_INT 32
696 #define COINFACTORIZATION_SHIFT_PER_INT 5
697 #define COINFACTORIZATION_MASK_PER_INT 0x1f
698 #endif
699  template <class T> inline bool
700  pivot ( int pivotRow,
701  int pivotColumn,
702  CoinBigIndex pivotRowPosition,
703  CoinBigIndex pivotColumnPosition,
705  unsigned int workArea2[],
706  int increment2,
707  T markRow[] ,
708  int largeInteger)
709 {
710  int *indexColumnU = indexColumnU_.array();
714  int *indexRowU = indexRowU_.array();
715  CoinBigIndex *startRowU = startRowU_.array();
716  int *numberInRow = numberInRow_.array();
718  int *indexRowL = indexRowL_.array();
719  int *saveColumn = saveColumn_.array();
720  int *nextRow = nextRow_.array();
721  int *lastRow = lastRow_.array() ;
722 
723  //store pivot columns (so can easily compress)
724  int numberInPivotRow = numberInRow[pivotRow] - 1;
725  CoinBigIndex startColumn = startColumnU[pivotColumn];
726  int numberInPivotColumn = numberInColumn[pivotColumn] - 1;
727  CoinBigIndex endColumn = startColumn + numberInPivotColumn + 1;
728  int put = 0;
729  CoinBigIndex startRow = startRowU[pivotRow];
730  CoinBigIndex endRow = startRow + numberInPivotRow + 1;
731 
732  if ( pivotColumnPosition < 0 ) {
733  for ( pivotColumnPosition = startRow; pivotColumnPosition < endRow; pivotColumnPosition++ ) {
734  int iColumn = indexColumnU[pivotColumnPosition];
735  if ( iColumn != pivotColumn ) {
736  saveColumn[put++] = iColumn;
737  } else {
738  break;
739  }
740  }
741  } else {
742  for (CoinBigIndex i = startRow ; i < pivotColumnPosition ; i++ ) {
743  saveColumn[put++] = indexColumnU[i];
744  }
745  }
746  assert (pivotColumnPosition<endRow);
747  assert (indexColumnU[pivotColumnPosition]==pivotColumn);
748  pivotColumnPosition++;
749  for ( ; pivotColumnPosition < endRow; pivotColumnPosition++ ) {
750  saveColumn[put++] = indexColumnU[pivotColumnPosition];
751  }
752  //take out this bit of indexColumnU
753  int next = nextRow[pivotRow];
754  int last = lastRow[pivotRow];
755 
756  nextRow[last] = next;
757  lastRow[next] = last;
758  nextRow[pivotRow] = numberGoodU_; //use for permute
759  lastRow[pivotRow] = -2;
760  numberInRow[pivotRow] = 0;
761  //store column in L, compress in U and take column out
763 
764  if ( l + numberInPivotColumn > lengthAreaL_ ) {
765  //need more memory
766  if ((messageLevel_&4)!=0)
767  printf("more memory needed in middle of invert\n");
768  return false;
769  }
770  //l+=currentAreaL_->elementByColumn-elementL;
771  CoinBigIndex lSave = l;
772 
774  startColumnL[numberGoodL_] = l; //for luck and first time
775  numberGoodL_++;
776  startColumnL[numberGoodL_] = l + numberInPivotColumn;
777  lengthL_ += numberInPivotColumn;
778  if ( pivotRowPosition < 0 ) {
779  for ( pivotRowPosition = startColumn; pivotRowPosition < endColumn; pivotRowPosition++ ) {
780  int iRow = indexRowU[pivotRowPosition];
781  if ( iRow != pivotRow ) {
782  indexRowL[l] = iRow;
783  elementL[l] = elementU[pivotRowPosition];
784  markRow[iRow] = static_cast<T>(l - lSave);
785  l++;
786  //take out of row list
787  CoinBigIndex start = startRowU[iRow];
788  CoinBigIndex end = start + numberInRow[iRow];
789  CoinBigIndex where = start;
790 
791  while ( indexColumnU[where] != pivotColumn ) {
792  where++;
793  } /* endwhile */
794 #if DEBUG_COIN
795  if ( where >= end ) {
796  abort ( );
797  }
798 #endif
799  indexColumnU[where] = indexColumnU[end - 1];
800  numberInRow[iRow]--;
801  } else {
802  break;
803  }
804  }
805  } else {
806  CoinBigIndex i;
807 
808  for ( i = startColumn; i < pivotRowPosition; i++ ) {
809  int iRow = indexRowU[i];
810 
811  markRow[iRow] = static_cast<T>(l - lSave);
812  indexRowL[l] = iRow;
813  elementL[l] = elementU[i];
814  l++;
815  //take out of row list
816  CoinBigIndex start = startRowU[iRow];
817  CoinBigIndex end = start + numberInRow[iRow];
818  CoinBigIndex where = start;
819 
820  while ( indexColumnU[where] != pivotColumn ) {
821  where++;
822  } /* endwhile */
823 #if DEBUG_COIN
824  if ( where >= end ) {
825  abort ( );
826  }
827 #endif
828  indexColumnU[where] = indexColumnU[end - 1];
829  numberInRow[iRow]--;
830  assert (numberInRow[iRow]>=0);
831  }
832  }
833  assert (pivotRowPosition<endColumn);
834  assert (indexRowU[pivotRowPosition]==pivotRow);
835  CoinFactorizationDouble pivotElement = elementU[pivotRowPosition];
836  CoinFactorizationDouble pivotMultiplier = 1.0 / pivotElement;
837 
838  pivotRegion_.array()[numberGoodU_] = pivotMultiplier;
839  pivotRowPosition++;
840  for ( ; pivotRowPosition < endColumn; pivotRowPosition++ ) {
841  int iRow = indexRowU[pivotRowPosition];
842 
843  markRow[iRow] = static_cast<T>(l - lSave);
844  indexRowL[l] = iRow;
845  elementL[l] = elementU[pivotRowPosition];
846  l++;
847  //take out of row list
848  CoinBigIndex start = startRowU[iRow];
849  CoinBigIndex end = start + numberInRow[iRow];
850  CoinBigIndex where = start;
851 
852  while ( indexColumnU[where] != pivotColumn ) {
853  where++;
854  } /* endwhile */
855 #if DEBUG_COIN
856  if ( where >= end ) {
857  abort ( );
858  }
859 #endif
860  indexColumnU[where] = indexColumnU[end - 1];
861  numberInRow[iRow]--;
862  assert (numberInRow[iRow]>=0);
863  }
864  markRow[pivotRow] = static_cast<T>(largeInteger);
865  //compress pivot column (move pivot to front including saved)
866  numberInColumn[pivotColumn] = 0;
867  //use end of L for temporary space
868  int *indexL = &indexRowL[lSave];
869  CoinFactorizationDouble *multipliersL = &elementL[lSave];
870 
871  //adjust
872  int j;
873 
874  for ( j = 0; j < numberInPivotColumn; j++ ) {
875  multipliersL[j] *= pivotMultiplier;
876  }
877  //zero out fill
878  CoinBigIndex iErase;
879  for ( iErase = 0; iErase < increment2 * numberInPivotRow;
880  iErase++ ) {
881  workArea2[iErase] = 0;
882  }
883  CoinBigIndex added = numberInPivotRow * numberInPivotColumn;
884  unsigned int *temp2 = workArea2;
885  int * nextColumn = nextColumn_.array();
886 
887  //pack down and move to work
888  int jColumn;
889  for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
890  int iColumn = saveColumn[jColumn];
891  CoinBigIndex startColumn = startColumnU[iColumn];
892  CoinBigIndex endColumn = startColumn + numberInColumn[iColumn];
893  int iRow = indexRowU[startColumn];
894  CoinFactorizationDouble value = elementU[startColumn];
895  double largest;
896  CoinBigIndex put = startColumn;
897  CoinBigIndex positionLargest = -1;
898  CoinFactorizationDouble thisPivotValue = 0.0;
899 
900  //compress column and find largest not updated
901  bool checkLargest;
902  int mark = markRow[iRow];
903 
904  if ( mark == largeInteger+1 ) {
905  largest = fabs ( value );
906  positionLargest = put;
907  put++;
908  checkLargest = false;
909  } else {
910  //need to find largest
911  largest = 0.0;
912  checkLargest = true;
913  if ( mark != largeInteger ) {
914  //will be updated
915  work[mark] = value;
916  int word = mark >> COINFACTORIZATION_SHIFT_PER_INT;
917  int bit = mark & COINFACTORIZATION_MASK_PER_INT;
918 
919  temp2[word] = temp2[word] | ( 1 << bit ); //say already in counts
920  added--;
921  } else {
922  thisPivotValue = value;
923  }
924  }
925  CoinBigIndex i;
926  for ( i = startColumn + 1; i < endColumn; i++ ) {
927  iRow = indexRowU[i];
928  value = elementU[i];
929  int mark = markRow[iRow];
930 
931  if ( mark == largeInteger+1 ) {
932  //keep
933  indexRowU[put] = iRow;
934  elementU[put] = value;
935  if ( checkLargest ) {
936  double absValue = fabs ( value );
937 
938  if ( absValue > largest ) {
939  largest = absValue;
940  positionLargest = put;
941  }
942  }
943  put++;
944  } else if ( mark != largeInteger ) {
945  //will be updated
946  work[mark] = value;
947  int word = mark >> COINFACTORIZATION_SHIFT_PER_INT;
948  int bit = mark & COINFACTORIZATION_MASK_PER_INT;
949 
950  temp2[word] = temp2[word] | ( 1 << bit ); //say already in counts
951  added--;
952  } else {
953  thisPivotValue = value;
954  }
955  }
956  //slot in pivot
957  elementU[put] = elementU[startColumn];
958  indexRowU[put] = indexRowU[startColumn];
959  if ( positionLargest == startColumn ) {
960  positionLargest = put; //follow if was largest
961  }
962  put++;
963  elementU[startColumn] = thisPivotValue;
964  indexRowU[startColumn] = pivotRow;
965  //clean up counts
966  startColumn++;
967  numberInColumn[iColumn] = put - startColumn;
968  int * numberInColumnPlus = numberInColumnPlus_.array();
969  numberInColumnPlus[iColumn]++;
970  startColumnU[iColumn]++;
971  //how much space have we got
972  int next = nextColumn[iColumn];
973  CoinBigIndex space;
974 
975  space = startColumnU[next] - put - numberInColumnPlus[next];
976  //assume no zero elements
977  if ( numberInPivotColumn > space ) {
978  //getColumnSpace also moves fixed part
979  if ( !getColumnSpace ( iColumn, numberInPivotColumn ) ) {
980  return false;
981  }
982  //redo starts
983  positionLargest = positionLargest + startColumnU[iColumn] - startColumn;
984  startColumn = startColumnU[iColumn];
985  put = startColumn + numberInColumn[iColumn];
986  }
987  double tolerance = zeroTolerance_;
988 
989  int *nextCount = nextCount_.array();
990  for ( j = 0; j < numberInPivotColumn; j++ ) {
991  value = work[j] - thisPivotValue * multipliersL[j];
992  double absValue = fabs ( value );
993 
994  if ( absValue > tolerance ) {
995  work[j] = 0.0;
996  assert (put<lengthAreaU_);
997  elementU[put] = value;
998  indexRowU[put] = indexL[j];
999  if ( absValue > largest ) {
1000  largest = absValue;
1001  positionLargest = put;
1002  }
1003  put++;
1004  } else {
1005  work[j] = 0.0;
1006  added--;
1007  int word = j >> COINFACTORIZATION_SHIFT_PER_INT;
1008  int bit = j & COINFACTORIZATION_MASK_PER_INT;
1009 
1010  if ( temp2[word] & ( 1 << bit ) ) {
1011  //take out of row list
1012  iRow = indexL[j];
1013  CoinBigIndex start = startRowU[iRow];
1014  CoinBigIndex end = start + numberInRow[iRow];
1015  CoinBigIndex where = start;
1016 
1017  while ( indexColumnU[where] != iColumn ) {
1018  where++;
1019  } /* endwhile */
1020 #if DEBUG_COIN
1021  if ( where >= end ) {
1022  abort ( );
1023  }
1024 #endif
1025  indexColumnU[where] = indexColumnU[end - 1];
1026  numberInRow[iRow]--;
1027  } else {
1028  //make sure won't be added
1029  int word = j >> COINFACTORIZATION_SHIFT_PER_INT;
1030  int bit = j & COINFACTORIZATION_MASK_PER_INT;
1031 
1032  temp2[word] = temp2[word] | ( 1 << bit ); //say already in counts
1033  }
1034  }
1035  }
1036  numberInColumn[iColumn] = put - startColumn;
1037  //move largest
1038  if ( positionLargest >= 0 ) {
1039  value = elementU[positionLargest];
1040  iRow = indexRowU[positionLargest];
1041  elementU[positionLargest] = elementU[startColumn];
1042  indexRowU[positionLargest] = indexRowU[startColumn];
1043  elementU[startColumn] = value;
1044  indexRowU[startColumn] = iRow;
1045  }
1046  //linked list for column
1047  if ( nextCount[iColumn + numberRows_] != -2 ) {
1048  //modify linked list
1049  deleteLink ( iColumn + numberRows_ );
1050  addLink ( iColumn + numberRows_, numberInColumn[iColumn] );
1051  }
1052  temp2 += increment2;
1053  }
1054  //get space for row list
1055  unsigned int *putBase = workArea2;
1056  int bigLoops = numberInPivotColumn >> COINFACTORIZATION_SHIFT_PER_INT;
1057  int i = 0;
1058 
1059  // do linked lists and update counts
1060  while ( bigLoops ) {
1061  bigLoops--;
1062  int bit;
1063  for ( bit = 0; bit < COINFACTORIZATION_BITS_PER_INT; i++, bit++ ) {
1064  unsigned int *putThis = putBase;
1065  int iRow = indexL[i];
1066 
1067  //get space
1068  int number = 0;
1069  int jColumn;
1070 
1071  for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
1072  unsigned int test = *putThis;
1073 
1074  putThis += increment2;
1075  test = 1 - ( ( test >> bit ) & 1 );
1076  number += test;
1077  }
1078  int next = nextRow[iRow];
1079  CoinBigIndex space;
1080 
1081  space = startRowU[next] - startRowU[iRow];
1082  number += numberInRow[iRow];
1083  if ( space < number ) {
1084  if ( !getRowSpace ( iRow, number ) ) {
1085  return false;
1086  }
1087  }
1088  // now do
1089  putThis = putBase;
1090  next = nextRow[iRow];
1091  number = numberInRow[iRow];
1092  CoinBigIndex end = startRowU[iRow] + number;
1093  int saveIndex = indexColumnU[startRowU[next]];
1094 
1095  //add in
1096  for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
1097  unsigned int test = *putThis;
1098 
1099  putThis += increment2;
1100  test = 1 - ( ( test >> bit ) & 1 );
1101  indexColumnU[end] = saveColumn[jColumn];
1102  end += test;
1103  }
1104  //put back next one in case zapped
1105  indexColumnU[startRowU[next]] = saveIndex;
1106  markRow[iRow] = static_cast<T>(largeInteger+1);
1107  number = end - startRowU[iRow];
1108  numberInRow[iRow] = number;
1109  deleteLink ( iRow );
1110  addLink ( iRow, number );
1111  }
1112  putBase++;
1113  } /* endwhile */
1114  int bit;
1115 
1116  for ( bit = 0; i < numberInPivotColumn; i++, bit++ ) {
1117  unsigned int *putThis = putBase;
1118  int iRow = indexL[i];
1119 
1120  //get space
1121  int number = 0;
1122  int jColumn;
1123 
1124  for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
1125  unsigned int test = *putThis;
1126 
1127  putThis += increment2;
1128  test = 1 - ( ( test >> bit ) & 1 );
1129  number += test;
1130  }
1131  int next = nextRow[iRow];
1132  CoinBigIndex space;
1133 
1134  space = startRowU[next] - startRowU[iRow];
1135  number += numberInRow[iRow];
1136  if ( space < number ) {
1137  if ( !getRowSpace ( iRow, number ) ) {
1138  return false;
1139  }
1140  }
1141  // now do
1142  putThis = putBase;
1143  next = nextRow[iRow];
1144  number = numberInRow[iRow];
1145  CoinBigIndex end = startRowU[iRow] + number;
1146  int saveIndex;
1147 
1148  saveIndex = indexColumnU[startRowU[next]];
1149 
1150  //add in
1151  for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
1152  unsigned int test = *putThis;
1153 
1154  putThis += increment2;
1155  test = 1 - ( ( test >> bit ) & 1 );
1156 
1157  indexColumnU[end] = saveColumn[jColumn];
1158  end += test;
1159  }
1160  indexColumnU[startRowU[next]] = saveIndex;
1161  markRow[iRow] = static_cast<T>(largeInteger+1);
1162  number = end - startRowU[iRow];
1163  numberInRow[iRow] = number;
1164  deleteLink ( iRow );
1165  addLink ( iRow, number );
1166  }
1167  markRow[pivotRow] = static_cast<T>(largeInteger+1);
1168  //modify linked list for pivots
1169  deleteLink ( pivotRow );
1170  deleteLink ( pivotColumn + numberRows_ );
1171  totalElements_ += added;
1172  return true;
1173 }
1174 
1175  /********************************* END LARGE TEMPLATE ********/
1177 
1178 protected:
1179 
1182 
1186 #ifndef COIN_FAST_CODE
1187 
1188  double slackValue_;
1189 #else
1190 #ifndef slackValue_
1191 #define slackValue_ -1.0
1192 #endif
1193 #endif
1194 
1195  double areaFactor_;
1197  double relaxCheck_;
1232  int status_;
1233 
1238  //int increasingRows_;
1239 
1244 
1247 
1250 
1253 
1257 
1260 
1263 
1266 
1269 
1272 
1275 
1278 
1281 
1284 
1287 
1290 
1293 
1296 
1299 
1302 
1305 
1307  //int baseU_;
1308 
1311 
1314 
1317 
1320 
1323 
1326 
1329 
1332 
1335 
1338 
1341 
1344 
1347 
1350 
1353 
1356 
1359 
1362 
1365 
1368 
1370  double * denseArea_;
1371 
1374 
1377 
1380 
1383 
1386 
1389 
1391  mutable double ftranCountInput_;
1392  mutable double ftranCountAfterL_;
1393  mutable double ftranCountAfterR_;
1394  mutable double ftranCountAfterU_;
1395  mutable double btranCountInput_;
1396  mutable double btranCountAfterU_;
1397  mutable double btranCountAfterR_;
1398  mutable double btranCountAfterL_;
1399 
1401  mutable int numberFtranCounts_;
1402  mutable int numberBtranCounts_;
1403 
1411 
1413  mutable bool collectStatistics_;
1414 
1417 
1420 
1423 
1426 
1429 
1435  int biasLU_;
1443 };
1444 // Dense coding
1445 #ifdef COIN_HAS_LAPACK
1446 #define DENSE_CODE 1
1447 /* Type of Fortran integer translated into C */
1448 #ifndef ipfint
1449 //typedef ipfint FORTRAN_INTEGER_TYPE ;
1450 typedef int ipfint;
1451 typedef const int cipfint;
1452 #endif
1453 #endif
1454 #endif
1455 // Extra for ugly include
1456 #ifdef UGLY_COIN_FACTOR_CODING
1457 #define FAC_UNSET (FAC_SET+1)
1458 {
1459  goodPivot=false;
1460  //store pivot columns (so can easily compress)
1461  CoinBigIndex startColumnThis = startColumn[iPivotColumn];
1462  CoinBigIndex endColumn = startColumnThis + numberDoColumn + 1;
1463  int put = 0;
1464  CoinBigIndex startRowThis = startRow[iPivotRow];
1465  CoinBigIndex endRow = startRowThis + numberDoRow + 1;
1466  if ( pivotColumnPosition < 0 ) {
1467  for ( pivotColumnPosition = startRowThis; pivotColumnPosition < endRow; pivotColumnPosition++ ) {
1468  int iColumn = indexColumn[pivotColumnPosition];
1469  if ( iColumn != iPivotColumn ) {
1470  saveColumn[put++] = iColumn;
1471  } else {
1472  break;
1473  }
1474  }
1475  } else {
1476  for (CoinBigIndex i = startRowThis ; i < pivotColumnPosition ; i++ ) {
1477  saveColumn[put++] = indexColumn[i];
1478  }
1479  }
1480  assert (pivotColumnPosition<endRow);
1481  assert (indexColumn[pivotColumnPosition]==iPivotColumn);
1482  pivotColumnPosition++;
1483  for ( ; pivotColumnPosition < endRow; pivotColumnPosition++ ) {
1484  saveColumn[put++] = indexColumn[pivotColumnPosition];
1485  }
1486  //take out this bit of indexColumn
1487  int next = nextRow[iPivotRow];
1488  int last = lastRow[iPivotRow];
1489 
1490  nextRow[last] = next;
1491  lastRow[next] = last;
1492  nextRow[iPivotRow] = numberGoodU_; //use for permute
1493  lastRow[iPivotRow] = -2;
1494  numberInRow[iPivotRow] = 0;
1495  //store column in L, compress in U and take column out
1496  CoinBigIndex l = lengthL_;
1497  // **** HORRID coding coming up but a goto seems best!
1498  {
1499  if ( l + numberDoColumn > lengthAreaL_ ) {
1500  //need more memory
1501  if ((messageLevel_&4)!=0)
1502  printf("more memory needed in middle of invert\n");
1503  goto BAD_PIVOT;
1504  }
1505  //l+=currentAreaL_->elementByColumn-elementL;
1506  CoinBigIndex lSave = l;
1507 
1508  CoinBigIndex * startColumnL = startColumnL_.array();
1509  startColumnL[numberGoodL_] = l; //for luck and first time
1510  numberGoodL_++;
1511  startColumnL[numberGoodL_] = l + numberDoColumn;
1512  lengthL_ += numberDoColumn;
1513  if ( pivotRowPosition < 0 ) {
1514  for ( pivotRowPosition = startColumnThis; pivotRowPosition < endColumn; pivotRowPosition++ ) {
1515  int iRow = indexRow[pivotRowPosition];
1516  if ( iRow != iPivotRow ) {
1517  indexRowL[l] = iRow;
1518  elementL[l] = element[pivotRowPosition];
1519  markRow[iRow] = l - lSave;
1520  l++;
1521  //take out of row list
1522  CoinBigIndex start = startRow[iRow];
1523  CoinBigIndex end = start + numberInRow[iRow];
1524  CoinBigIndex where = start;
1525 
1526  while ( indexColumn[where] != iPivotColumn ) {
1527  where++;
1528  } /* endwhile */
1529 #if DEBUG_COIN
1530  if ( where >= end ) {
1531  abort ( );
1532  }
1533 #endif
1534  indexColumn[where] = indexColumn[end - 1];
1535  numberInRow[iRow]--;
1536  } else {
1537  break;
1538  }
1539  }
1540  } else {
1541  CoinBigIndex i;
1542 
1543  for ( i = startColumnThis; i < pivotRowPosition; i++ ) {
1544  int iRow = indexRow[i];
1545 
1546  markRow[iRow] = l - lSave;
1547  indexRowL[l] = iRow;
1548  elementL[l] = element[i];
1549  l++;
1550  //take out of row list
1551  CoinBigIndex start = startRow[iRow];
1552  CoinBigIndex end = start + numberInRow[iRow];
1553  CoinBigIndex where = start;
1554 
1555  while ( indexColumn[where] != iPivotColumn ) {
1556  where++;
1557  } /* endwhile */
1558 #if DEBUG_COIN
1559  if ( where >= end ) {
1560  abort ( );
1561  }
1562 #endif
1563  indexColumn[where] = indexColumn[end - 1];
1564  numberInRow[iRow]--;
1565  assert (numberInRow[iRow]>=0);
1566  }
1567  }
1568  assert (pivotRowPosition<endColumn);
1569  assert (indexRow[pivotRowPosition]==iPivotRow);
1570  CoinFactorizationDouble pivotElement = element[pivotRowPosition];
1571  CoinFactorizationDouble pivotMultiplier = 1.0 / pivotElement;
1572 
1573  pivotRegion_.array()[numberGoodU_] = pivotMultiplier;
1574  pivotRowPosition++;
1575  for ( ; pivotRowPosition < endColumn; pivotRowPosition++ ) {
1576  int iRow = indexRow[pivotRowPosition];
1577 
1578  markRow[iRow] = l - lSave;
1579  indexRowL[l] = iRow;
1580  elementL[l] = element[pivotRowPosition];
1581  l++;
1582  //take out of row list
1583  CoinBigIndex start = startRow[iRow];
1584  CoinBigIndex end = start + numberInRow[iRow];
1585  CoinBigIndex where = start;
1586 
1587  while ( indexColumn[where] != iPivotColumn ) {
1588  where++;
1589  } /* endwhile */
1590 #if DEBUG_COIN
1591  if ( where >= end ) {
1592  abort ( );
1593  }
1594 #endif
1595  indexColumn[where] = indexColumn[end - 1];
1596  numberInRow[iRow]--;
1597  assert (numberInRow[iRow]>=0);
1598  }
1599  markRow[iPivotRow] = FAC_SET;
1600  //compress pivot column (move pivot to front including saved)
1601  numberInColumn[iPivotColumn] = 0;
1602  //use end of L for temporary space
1603  int *indexL = &indexRowL[lSave];
1604  CoinFactorizationDouble *multipliersL = &elementL[lSave];
1605 
1606  //adjust
1607  int j;
1608 
1609  for ( j = 0; j < numberDoColumn; j++ ) {
1610  multipliersL[j] *= pivotMultiplier;
1611  }
1612  //zero out fill
1613  CoinBigIndex iErase;
1614  for ( iErase = 0; iErase < increment2 * numberDoRow;
1615  iErase++ ) {
1616  workArea2[iErase] = 0;
1617  }
1618  CoinBigIndex added = numberDoRow * numberDoColumn;
1619  unsigned int *temp2 = workArea2;
1620  int * nextColumn = nextColumn_.array();
1621 
1622  //pack down and move to work
1623  int jColumn;
1624  for ( jColumn = 0; jColumn < numberDoRow; jColumn++ ) {
1625  int iColumn = saveColumn[jColumn];
1626  CoinBigIndex startColumnThis = startColumn[iColumn];
1627  CoinBigIndex endColumn = startColumnThis + numberInColumn[iColumn];
1628  int iRow = indexRow[startColumnThis];
1629  CoinFactorizationDouble value = element[startColumnThis];
1630  double largest;
1631  CoinBigIndex put = startColumnThis;
1632  CoinBigIndex positionLargest = -1;
1633  CoinFactorizationDouble thisPivotValue = 0.0;
1634 
1635  //compress column and find largest not updated
1636  bool checkLargest;
1637  int mark = markRow[iRow];
1638 
1639  if ( mark == FAC_UNSET ) {
1640  largest = fabs ( value );
1641  positionLargest = put;
1642  put++;
1643  checkLargest = false;
1644  } else {
1645  //need to find largest
1646  largest = 0.0;
1647  checkLargest = true;
1648  if ( mark != FAC_SET ) {
1649  //will be updated
1650  workArea[mark] = value;
1651  int word = mark >> COINFACTORIZATION_SHIFT_PER_INT;
1652  int bit = mark & COINFACTORIZATION_MASK_PER_INT;
1653 
1654  temp2[word] = temp2[word] | ( 1 << bit ); //say already in counts
1655  added--;
1656  } else {
1657  thisPivotValue = value;
1658  }
1659  }
1660  CoinBigIndex i;
1661  for ( i = startColumnThis + 1; i < endColumn; i++ ) {
1662  iRow = indexRow[i];
1663  value = element[i];
1664  int mark = markRow[iRow];
1665 
1666  if ( mark == FAC_UNSET ) {
1667  //keep
1668  indexRow[put] = iRow;
1669  element[put] = value;
1670  if ( checkLargest ) {
1671  double absValue = fabs ( value );
1672 
1673  if ( absValue > largest ) {
1674  largest = absValue;
1675  positionLargest = put;
1676  }
1677  }
1678  put++;
1679  } else if ( mark != FAC_SET ) {
1680  //will be updated
1681  workArea[mark] = value;
1682  int word = mark >> COINFACTORIZATION_SHIFT_PER_INT;
1683  int bit = mark & COINFACTORIZATION_MASK_PER_INT;
1684 
1685  temp2[word] = temp2[word] | ( 1 << bit ); //say already in counts
1686  added--;
1687  } else {
1688  thisPivotValue = value;
1689  }
1690  }
1691  //slot in pivot
1692  element[put] = element[startColumnThis];
1693  indexRow[put] = indexRow[startColumnThis];
1694  if ( positionLargest == startColumnThis ) {
1695  positionLargest = put; //follow if was largest
1696  }
1697  put++;
1698  element[startColumnThis] = thisPivotValue;
1699  indexRow[startColumnThis] = iPivotRow;
1700  //clean up counts
1701  startColumnThis++;
1702  numberInColumn[iColumn] = put - startColumnThis;
1703  int * numberInColumnPlus = numberInColumnPlus_.array();
1704  numberInColumnPlus[iColumn]++;
1705  startColumn[iColumn]++;
1706  //how much space have we got
1707  int next = nextColumn[iColumn];
1708  CoinBigIndex space;
1709 
1710  space = startColumn[next] - put - numberInColumnPlus[next];
1711  //assume no zero elements
1712  if ( numberDoColumn > space ) {
1713  //getColumnSpace also moves fixed part
1714  if ( !getColumnSpace ( iColumn, numberDoColumn ) ) {
1715  goto BAD_PIVOT;
1716  }
1717  //redo starts
1718  positionLargest = positionLargest + startColumn[iColumn] - startColumnThis;
1719  startColumnThis = startColumn[iColumn];
1720  put = startColumnThis + numberInColumn[iColumn];
1721  }
1722  double tolerance = zeroTolerance_;
1723 
1724  int *nextCount = nextCount_.array();
1725  for ( j = 0; j < numberDoColumn; j++ ) {
1726  value = workArea[j] - thisPivotValue * multipliersL[j];
1727  double absValue = fabs ( value );
1728 
1729  if ( absValue > tolerance ) {
1730  workArea[j] = 0.0;
1731  element[put] = value;
1732  indexRow[put] = indexL[j];
1733  if ( absValue > largest ) {
1734  largest = absValue;
1735  positionLargest = put;
1736  }
1737  put++;
1738  } else {
1739  workArea[j] = 0.0;
1740  added--;
1741  int word = j >> COINFACTORIZATION_SHIFT_PER_INT;
1742  int bit = j & COINFACTORIZATION_MASK_PER_INT;
1743 
1744  if ( temp2[word] & ( 1 << bit ) ) {
1745  //take out of row list
1746  iRow = indexL[j];
1747  CoinBigIndex start = startRow[iRow];
1748  CoinBigIndex end = start + numberInRow[iRow];
1749  CoinBigIndex where = start;
1750 
1751  while ( indexColumn[where] != iColumn ) {
1752  where++;
1753  } /* endwhile */
1754 #if DEBUG_COIN
1755  if ( where >= end ) {
1756  abort ( );
1757  }
1758 #endif
1759  indexColumn[where] = indexColumn[end - 1];
1760  numberInRow[iRow]--;
1761  } else {
1762  //make sure won't be added
1763  int word = j >> COINFACTORIZATION_SHIFT_PER_INT;
1764  int bit = j & COINFACTORIZATION_MASK_PER_INT;
1765 
1766  temp2[word] = temp2[word] | ( 1 << bit ); //say already in counts
1767  }
1768  }
1769  }
1770  numberInColumn[iColumn] = put - startColumnThis;
1771  //move largest
1772  if ( positionLargest >= 0 ) {
1773  value = element[positionLargest];
1774  iRow = indexRow[positionLargest];
1775  element[positionLargest] = element[startColumnThis];
1776  indexRow[positionLargest] = indexRow[startColumnThis];
1777  element[startColumnThis] = value;
1778  indexRow[startColumnThis] = iRow;
1779  }
1780  //linked list for column
1781  if ( nextCount[iColumn + numberRows_] != -2 ) {
1782  //modify linked list
1783  deleteLink ( iColumn + numberRows_ );
1784  addLink ( iColumn + numberRows_, numberInColumn[iColumn] );
1785  }
1786  temp2 += increment2;
1787  }
1788  //get space for row list
1789  unsigned int *putBase = workArea2;
1790  int bigLoops = numberDoColumn >> COINFACTORIZATION_SHIFT_PER_INT;
1791  int i = 0;
1792 
1793  // do linked lists and update counts
1794  while ( bigLoops ) {
1795  bigLoops--;
1796  int bit;
1797  for ( bit = 0; bit < COINFACTORIZATION_BITS_PER_INT; i++, bit++ ) {
1798  unsigned int *putThis = putBase;
1799  int iRow = indexL[i];
1800 
1801  //get space
1802  int number = 0;
1803  int jColumn;
1804 
1805  for ( jColumn = 0; jColumn < numberDoRow; jColumn++ ) {
1806  unsigned int test = *putThis;
1807 
1808  putThis += increment2;
1809  test = 1 - ( ( test >> bit ) & 1 );
1810  number += test;
1811  }
1812  int next = nextRow[iRow];
1813  CoinBigIndex space;
1814 
1815  space = startRow[next] - startRow[iRow];
1816  number += numberInRow[iRow];
1817  if ( space < number ) {
1818  if ( !getRowSpace ( iRow, number ) ) {
1819  goto BAD_PIVOT;
1820  }
1821  }
1822  // now do
1823  putThis = putBase;
1824  next = nextRow[iRow];
1825  number = numberInRow[iRow];
1826  CoinBigIndex end = startRow[iRow] + number;
1827  int saveIndex = indexColumn[startRow[next]];
1828 
1829  //add in
1830  for ( jColumn = 0; jColumn < numberDoRow; jColumn++ ) {
1831  unsigned int test = *putThis;
1832 
1833  putThis += increment2;
1834  test = 1 - ( ( test >> bit ) & 1 );
1835  indexColumn[end] = saveColumn[jColumn];
1836  end += test;
1837  }
1838  //put back next one in case zapped
1839  indexColumn[startRow[next]] = saveIndex;
1840  markRow[iRow] = FAC_UNSET;
1841  number = end - startRow[iRow];
1842  numberInRow[iRow] = number;
1843  deleteLink ( iRow );
1844  addLink ( iRow, number );
1845  }
1846  putBase++;
1847  } /* endwhile */
1848  int bit;
1849 
1850  for ( bit = 0; i < numberDoColumn; i++, bit++ ) {
1851  unsigned int *putThis = putBase;
1852  int iRow = indexL[i];
1853 
1854  //get space
1855  int number = 0;
1856  int jColumn;
1857 
1858  for ( jColumn = 0; jColumn < numberDoRow; jColumn++ ) {
1859  unsigned int test = *putThis;
1860 
1861  putThis += increment2;
1862  test = 1 - ( ( test >> bit ) & 1 );
1863  number += test;
1864  }
1865  int next = nextRow[iRow];
1866  CoinBigIndex space;
1867 
1868  space = startRow[next] - startRow[iRow];
1869  number += numberInRow[iRow];
1870  if ( space < number ) {
1871  if ( !getRowSpace ( iRow, number ) ) {
1872  goto BAD_PIVOT;
1873  }
1874  }
1875  // now do
1876  putThis = putBase;
1877  next = nextRow[iRow];
1878  number = numberInRow[iRow];
1879  CoinBigIndex end = startRow[iRow] + number;
1880  int saveIndex;
1881 
1882  saveIndex = indexColumn[startRow[next]];
1883 
1884  //add in
1885  for ( jColumn = 0; jColumn < numberDoRow; jColumn++ ) {
1886  unsigned int test = *putThis;
1887 
1888  putThis += increment2;
1889  test = 1 - ( ( test >> bit ) & 1 );
1890 
1891  indexColumn[end] = saveColumn[jColumn];
1892  end += test;
1893  }
1894  indexColumn[startRow[next]] = saveIndex;
1895  markRow[iRow] = FAC_UNSET;
1896  number = end - startRow[iRow];
1897  numberInRow[iRow] = number;
1898  deleteLink ( iRow );
1899  addLink ( iRow, number );
1900  }
1901  markRow[iPivotRow] = FAC_UNSET;
1902  //modify linked list for pivots
1903  deleteLink ( iPivotRow );
1904  deleteLink ( iPivotColumn + numberRows_ );
1905  totalElements_ += added;
1906  goodPivot= true;
1907  // **** UGLY UGLY UGLY
1908  }
1909  BAD_PIVOT:
1910 
1911  ;
1912 }
1913 #undef FAC_UNSET
1914 #endif