00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #ifndef CoinFactorization_H
00011 #define CoinFactorization_H
00012
00013 #include <iostream>
00014 #include <string>
00015 #include <cassert>
00016 #include "CoinFinite.hpp"
00017 #include "CoinIndexedVector.hpp"
00018 class CoinPackedMatrix;
00045 class CoinFactorization {
00046 friend void CoinFactorizationUnitTest( const std::string & mpsDir );
00047
00048 public:
00049
00052
00053 CoinFactorization ( );
00055 CoinFactorization ( const CoinFactorization &other);
00056
00058 ~CoinFactorization ( );
00060 void almostDestructor();
00062 void show_self ( ) const;
00064 int saveFactorization (const char * file ) const;
00068 int restoreFactorization (const char * file , bool factor=false) ;
00070 void sort ( ) const;
00072 CoinFactorization & operator = ( const CoinFactorization & other );
00074
00084 int factorize ( const CoinPackedMatrix & matrix,
00085 int rowIsBasic[], int columnIsBasic[] ,
00086 double areaFactor = 0.0 );
00097 int factorize ( int numberRows,
00098 int numberColumns,
00099 CoinBigIndex numberElements,
00100 CoinBigIndex maximumL,
00101 CoinBigIndex maximumU,
00102 const int indicesRow[],
00103 const int indicesColumn[], const double elements[] ,
00104 int permutation[],
00105 double areaFactor = 0.0);
00110 int factorizePart1 ( int numberRows,
00111 int numberColumns,
00112 CoinBigIndex estimateNumberElements,
00113 int * indicesRow[],
00114 int * indicesColumn[],
00115 double * elements[],
00116 double areaFactor = 0.0);
00123 int factorizePart2 (int permutation[],int exactNumberElements);
00125
00128
00129 inline int status ( ) const {
00130 return status_;
00131 };
00133 inline int pivots ( ) const {
00134 return numberPivots_;
00135 };
00137 inline int *permute ( ) const {
00138 return permute_.array();
00139 };
00141 inline int *pivotColumn ( ) const {
00142 return pivotColumn_.array();
00143 };
00145 inline int *permuteBack ( ) const {
00146 return permuteBack_.array();
00147 };
00149 inline int *pivotColumnBack ( ) const {
00150 return pivotColumnBack_.array();
00151 };
00153 inline int numberRowsExtra ( ) const {
00154 return numberRowsExtra_;
00155 };
00157 inline int numberRows ( ) const {
00158 return numberRows_;
00159 };
00161 inline int maximumRowsExtra ( ) const {
00162 return maximumRowsExtra_;
00163 };
00165 inline int numberColumns ( ) const {
00166 return numberColumns_;
00167 };
00169 inline int numberElements ( ) const {
00170 return totalElements_;
00171 };
00173 inline int numberForrestTomlin ( ) const {
00174 return numberInColumn_.array()[numberColumnsExtra_];
00175 };
00177 inline int numberGoodColumns ( ) const {
00178 return numberGoodU_;
00179 };
00181 inline double areaFactor ( ) const {
00182 return areaFactor_;
00183 };
00184 inline void areaFactor ( double value ) {
00185 areaFactor_=value;
00186 };
00188 double adjustedAreaFactor() const;
00190 inline void relaxAccuracyCheck(double value)
00191 { relaxCheck_ = value;};
00192 inline double getAccuracyCheck() const
00193 { return relaxCheck_;};
00195 inline bool increasingRows ( ) const
00196 { return true; };
00202 inline void increasingRows ( int value ) {};
00204 inline int messageLevel ( ) const {
00205 return messageLevel_ ;
00206 };
00207 void messageLevel ( int value );
00209 inline int maximumPivots ( ) const {
00210 return maximumPivots_ ;
00211 };
00212 void maximumPivots ( int value );
00213
00215 inline int denseThreshold() const
00216 { return denseThreshold_;};
00218 inline void setDenseThreshold(int value)
00219 { denseThreshold_ = value;};
00221 inline double pivotTolerance ( ) const {
00222 return pivotTolerance_ ;
00223 };
00224 void pivotTolerance ( double value );
00226 inline double zeroTolerance ( ) const {
00227 return zeroTolerance_ ;
00228 };
00229 void zeroTolerance ( double value );
00231 inline double slackValue ( ) const {
00232 return slackValue_ ;
00233 };
00234 void slackValue ( double value );
00236 double maximumCoefficient() const;
00238 inline bool forrestTomlin() const
00239 { return doForrestTomlin_;};
00240 inline void setForrestTomlin(bool value)
00241 { doForrestTomlin_=value;};
00243
00246
00248 inline int numberDense() const
00249 { return numberDense_;};
00250
00252 inline CoinBigIndex numberElementsU ( ) const {
00253 return lengthU_;
00254 };
00256 inline CoinBigIndex lengthAreaU ( ) const {
00257 return lengthAreaU_;
00258 };
00260 inline CoinBigIndex numberElementsL ( ) const {
00261 return lengthL_;
00262 };
00264 inline CoinBigIndex lengthAreaL ( ) const {
00265 return lengthAreaL_;
00266 };
00268 inline CoinBigIndex numberElementsR ( ) const {
00269 return lengthR_;
00270 };
00272 inline CoinBigIndex numberCompressions() const
00273 { return numberCompressions_;};
00277 inline int biasLU() const
00278 { return biasLU_;};
00279 inline void setBiasLU(int value)
00280 { biasLU_=value;};
00286 inline int persistenceFlag() const
00287 { return persistenceFlag_;};
00288 void setPersistenceFlag(int value);
00290
00293
00301 int replaceColumn ( CoinIndexedVector * regionSparse,
00302 int pivotRow,
00303 double pivotCheck ,
00304 bool checkBeforeModifying=false);
00306
00316 int updateColumnFT ( CoinIndexedVector * regionSparse,
00317 CoinIndexedVector * regionSparse2);
00320 int updateColumn ( CoinIndexedVector * regionSparse,
00321 CoinIndexedVector * regionSparse2,
00322 bool noPermute=false) const;
00327 int updateColumnTranspose ( CoinIndexedVector * regionSparse,
00328 CoinIndexedVector * regionSparse2) const;
00330 void goSparse();
00332 int sparseThreshold ( ) const;
00334 void sparseThreshold ( int value );
00336
00337
00341
00342 inline void clearArrays()
00343 { gutsOfDestructor();};
00345
00348
00351 int add ( CoinBigIndex numberElements,
00352 int indicesRow[],
00353 int indicesColumn[], double elements[] );
00354
00357 int addColumn ( CoinBigIndex numberElements,
00358 int indicesRow[], double elements[] );
00359
00362 int addRow ( CoinBigIndex numberElements,
00363 int indicesColumn[], double elements[] );
00364
00366 int deleteColumn ( int Row );
00368 int deleteRow ( int Row );
00369
00373 int replaceRow ( int whichRow, int numberElements,
00374 const int indicesColumn[], const double elements[] );
00376 void emptyRows(int numberToEmpty, const int which[]);
00378 protected:
00379
00381
00382 void getAreas ( int numberRows,
00383 int numberColumns,
00384 CoinBigIndex maximumL,
00385 CoinBigIndex maximumU );
00386
00389 void preProcess ( int state,
00390 int possibleDuplicates = -1 );
00392 int factor ( );
00395 int factorSparse ( );
00398 int factorDense ( );
00399
00401 bool pivotOneOtherRow ( int pivotRow,
00402 int pivotColumn );
00404 bool pivotRowSingleton ( int pivotRow,
00405 int pivotColumn );
00407 bool pivotColumnSingleton ( int pivotRow,
00408 int pivotColumn );
00409
00414 bool getColumnSpace ( int iColumn,
00415 int extraNeeded );
00416
00420 bool getColumnSpaceIterateR ( int iColumn, double value,
00421 int iRow);
00427 CoinBigIndex getColumnSpaceIterate ( int iColumn, double value,
00428 int iRow);
00432 bool getRowSpace ( int iRow, int extraNeeded );
00433
00437 bool getRowSpaceIterate ( int iRow,
00438 int extraNeeded );
00440 void checkConsistency ( );
00442 inline void addLink ( int index, int count ) {
00443 int *nextCount = nextCount_.array();
00444 int *firstCount = firstCount_.array();
00445 int *lastCount = lastCount_.array();
00446 int next = firstCount[count];
00447 lastCount[index] = -2 - count;
00448 if ( next < 0 ) {
00449
00450 firstCount[count] = index;
00451 nextCount[index] = -1;
00452 } else {
00453 firstCount[count] = index;
00454 nextCount[index] = next;
00455 lastCount[next] = index;
00456 }};
00458 inline void deleteLink ( int index ) {
00459 int *nextCount = nextCount_.array();
00460 int *firstCount = firstCount_.array();
00461 int *lastCount = lastCount_.array();
00462 int next = nextCount[index];
00463 int last = lastCount[index];
00464 if ( last >= 0 ) {
00465 nextCount[last] = next;
00466 } else {
00467 int count = -last - 2;
00468
00469 firstCount[count] = next;
00470 }
00471 if ( next >= 0 ) {
00472 lastCount[next] = last;
00473 }
00474 nextCount[index] = -2;
00475 lastCount[index] = -2;
00476 return;
00477 };
00479 void separateLinks(int count,bool rowsFirst);
00481 void cleanup ( );
00482
00484 void updateColumnL ( CoinIndexedVector * region, int * indexIn ) const;
00486 void updateColumnLDensish ( CoinIndexedVector * region, int * indexIn ) const;
00488 void updateColumnLSparse ( CoinIndexedVector * region, int * indexIn ) const;
00490 void updateColumnLSparsish ( CoinIndexedVector * region, int * indexIn ) const;
00491
00493 void updateColumnR ( CoinIndexedVector * region ) const;
00496 void updateColumnRFT ( CoinIndexedVector * region, int * indexIn );
00497
00499 void updateColumnU ( CoinIndexedVector * region, int * indexIn) const;
00500
00502 void updateColumnUSparse ( CoinIndexedVector * regionSparse,
00503 int * indexIn) const;
00505 void updateColumnUSparsish ( CoinIndexedVector * regionSparse,
00506 int * indexIn) const;
00508 void updateColumnUDensish ( CoinIndexedVector * regionSparse,
00509 int * indexIn) const;
00511 void updateColumnPFI ( CoinIndexedVector * regionSparse) const;
00513 void permuteBack ( CoinIndexedVector * regionSparse,
00514 CoinIndexedVector * outVector) const;
00515
00517 void updateColumnTransposePFI ( CoinIndexedVector * region) const;
00520 void updateColumnTransposeU ( CoinIndexedVector * region,
00521 int smallestIndex) const;
00524 void updateColumnTransposeUSparsish ( CoinIndexedVector * region,
00525 int smallestIndex) const;
00528 void updateColumnTransposeUDensish ( CoinIndexedVector * region,
00529 int smallestIndex) const;
00532 void updateColumnTransposeUSparse ( CoinIndexedVector * region) const;
00533
00535 void updateColumnTransposeR ( CoinIndexedVector * region ) const;
00537 void updateColumnTransposeRDensish ( CoinIndexedVector * region ) const;
00539 void updateColumnTransposeRSparse ( CoinIndexedVector * region ) const;
00540
00542 void updateColumnTransposeL ( CoinIndexedVector * region ) const;
00544 void updateColumnTransposeLDensish ( CoinIndexedVector * region ) const;
00546 void updateColumnTransposeLByRow ( CoinIndexedVector * region ) const;
00548 void updateColumnTransposeLSparsish ( CoinIndexedVector * region ) const;
00550 void updateColumnTransposeLSparse ( CoinIndexedVector * region ) const;
00555 int replaceColumnPFI ( CoinIndexedVector * regionSparse,
00556 int pivotRow, double alpha);
00557
00560 int checkPivot(double saveFromU, double oldPivot) const;
00562 void gutsOfDestructor(int type=1);
00564 void gutsOfInitialize(int type);
00565 void gutsOfCopy(const CoinFactorization &other);
00566
00568 void resetStatistics();
00569
00570
00571 #ifdef INT_IS_8
00572 #define COINFACTORIZATION_BITS_PER_INT 64
00573 #define COINFACTORIZATION_SHIFT_PER_INT 6
00574 #define COINFACTORIZATION_MASK_PER_INT 0x3f
00575 #else
00576 #define COINFACTORIZATION_BITS_PER_INT 32
00577 #define COINFACTORIZATION_SHIFT_PER_INT 5
00578 #define COINFACTORIZATION_MASK_PER_INT 0x1f
00579 #endif
00580 template <class T> inline bool
00581 pivot ( int pivotRow,
00582 int pivotColumn,
00583 CoinBigIndex pivotRowPosition,
00584 CoinBigIndex pivotColumnPosition,
00585 double work[],
00586 unsigned int workArea2[],
00587 int increment,
00588 int increment2,
00589 T markRow[] ,
00590 int largeInteger)
00591 {
00592 int *indexColumnU = indexColumnU_.array();
00593 CoinBigIndex *startColumnU = startColumnU_.array();
00594 int *numberInColumn = numberInColumn_.array();
00595 double *elementU = elementU_.array();
00596 int *indexRowU = indexRowU_.array();
00597 CoinBigIndex *startRowU = startRowU_.array();
00598 int *numberInRow = numberInRow_.array();
00599 double *elementL = elementL_.array();
00600 int *indexRowL = indexRowL_.array();
00601 int *saveColumn = saveColumn_.array();
00602 int *nextRow = nextRow_.array();
00603 int *lastRow = lastRow_.array() ;
00604
00605
00606 int numberInPivotRow = numberInRow[pivotRow] - 1;
00607 CoinBigIndex startColumn = startColumnU[pivotColumn];
00608 int numberInPivotColumn = numberInColumn[pivotColumn] - 1;
00609 CoinBigIndex endColumn = startColumn + numberInPivotColumn + 1;
00610 int put = 0;
00611 CoinBigIndex startRow = startRowU[pivotRow];
00612 CoinBigIndex endRow = startRow + numberInPivotRow + 1;
00613
00614 if ( pivotColumnPosition < 0 ) {
00615 for ( pivotColumnPosition = startRow; pivotColumnPosition < endRow; pivotColumnPosition++ ) {
00616 int iColumn = indexColumnU[pivotColumnPosition];
00617 if ( iColumn != pivotColumn ) {
00618 saveColumn[put++] = iColumn;
00619 } else {
00620 break;
00621 }
00622 }
00623 } else {
00624 for (CoinBigIndex i = startRow ; i < pivotColumnPosition ; i++ ) {
00625 saveColumn[put++] = indexColumnU[i];
00626 }
00627 }
00628 assert (pivotColumnPosition<endRow);
00629 assert (indexColumnU[pivotColumnPosition]==pivotColumn);
00630 pivotColumnPosition++;
00631 for ( ; pivotColumnPosition < endRow; pivotColumnPosition++ ) {
00632 saveColumn[put++] = indexColumnU[pivotColumnPosition];
00633 }
00634
00635 int next = nextRow[pivotRow];
00636 int last = lastRow[pivotRow];
00637
00638 nextRow[last] = next;
00639 lastRow[next] = last;
00640 nextRow[pivotRow] = numberGoodU_;
00641 lastRow[pivotRow] = -2;
00642 numberInRow[pivotRow] = 0;
00643
00644 CoinBigIndex l = lengthL_;
00645
00646 if ( l + numberInPivotColumn > lengthAreaL_ ) {
00647
00648 printf("more memory needed in middle of invert\n");
00649 return false;
00650 }
00651
00652 CoinBigIndex lSave = l;
00653
00654 pivotRowL_.array()[numberGoodL_] = pivotRow;
00655 CoinBigIndex * startColumnL = startColumnL_.array();
00656 startColumnL[numberGoodL_] = l;
00657 numberGoodL_++;
00658 startColumnL[numberGoodL_] = l + numberInPivotColumn;
00659 lengthL_ += numberInPivotColumn;
00660 if ( pivotRowPosition < 0 ) {
00661 for ( pivotRowPosition = startColumn; pivotRowPosition < endColumn; pivotRowPosition++ ) {
00662 int iRow = indexRowU[pivotRowPosition];
00663 if ( iRow != pivotRow ) {
00664 indexRowL[l] = iRow;
00665 elementL[l] = elementU[pivotRowPosition];
00666 markRow[iRow] = l - lSave;
00667 l++;
00668
00669 CoinBigIndex start = startRowU[iRow];
00670 CoinBigIndex end = start + numberInRow[iRow];
00671 CoinBigIndex where = start;
00672
00673 while ( indexColumnU[where] != pivotColumn ) {
00674 where++;
00675 }
00676 #if DEBUG_COIN
00677 if ( where >= end ) {
00678 abort ( );
00679 }
00680 #endif
00681 indexColumnU[where] = indexColumnU[end - 1];
00682 numberInRow[iRow]--;
00683 } else {
00684 break;
00685 }
00686 }
00687 } else {
00688 CoinBigIndex i;
00689
00690 for ( i = startColumn; i < pivotRowPosition; i++ ) {
00691 int iRow = indexRowU[i];
00692
00693 markRow[iRow] = l - lSave;
00694 indexRowL[l] = iRow;
00695 elementL[l] = elementU[i];
00696 l++;
00697
00698 CoinBigIndex start = startRowU[iRow];
00699 CoinBigIndex end = start + numberInRow[iRow];
00700 CoinBigIndex where = start;
00701
00702 while ( indexColumnU[where] != pivotColumn ) {
00703 where++;
00704 }
00705 #if DEBUG_COIN
00706 if ( where >= end ) {
00707 abort ( );
00708 }
00709 #endif
00710 indexColumnU[where] = indexColumnU[end - 1];
00711 numberInRow[iRow]--;
00712 assert (numberInRow[iRow]>=0);
00713 }
00714 }
00715 assert (pivotRowPosition<endColumn);
00716 assert (indexRowU[pivotRowPosition]==pivotRow);
00717 double pivotElement = elementU[pivotRowPosition];
00718 double pivotMultiplier = 1.0 / pivotElement;
00719
00720 pivotRegion_.array()[numberGoodU_] = pivotMultiplier;
00721 pivotRowPosition++;
00722 for ( ; pivotRowPosition < endColumn; pivotRowPosition++ ) {
00723 int iRow = indexRowU[pivotRowPosition];
00724
00725 markRow[iRow] = l - lSave;
00726 indexRowL[l] = iRow;
00727 elementL[l] = elementU[pivotRowPosition];
00728 l++;
00729
00730 CoinBigIndex start = startRowU[iRow];
00731 CoinBigIndex end = start + numberInRow[iRow];
00732 CoinBigIndex where = start;
00733
00734 while ( indexColumnU[where] != pivotColumn ) {
00735 where++;
00736 }
00737 #if DEBUG_COIN
00738 if ( where >= end ) {
00739 abort ( );
00740 }
00741 #endif
00742 indexColumnU[where] = indexColumnU[end - 1];
00743 numberInRow[iRow]--;
00744 assert (numberInRow[iRow]>=0);
00745 }
00746 markRow[pivotRow] = largeInteger;
00747
00748 numberInColumn[pivotColumn] = 0;
00749
00750 int *indexL = &indexRowL[lSave];
00751 double *multipliersL = &elementL[lSave];
00752
00753
00754 int j;
00755
00756 for ( j = 0; j < numberInPivotColumn; j++ ) {
00757 multipliersL[j] *= pivotMultiplier;
00758 }
00759
00760 CoinBigIndex iErase;
00761 for ( iErase = 0; iErase < increment2 * numberInPivotRow;
00762 iErase++ ) {
00763 workArea2[iErase] = 0;
00764 }
00765 CoinBigIndex added = numberInPivotRow * numberInPivotColumn;
00766 unsigned int *temp2 = workArea2;
00767 int * nextColumn = nextColumn_.array();
00768
00769
00770 int jColumn;
00771 for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
00772 int iColumn = saveColumn[jColumn];
00773 CoinBigIndex startColumn = startColumnU[iColumn];
00774 CoinBigIndex endColumn = startColumn + numberInColumn[iColumn];
00775 int iRow = indexRowU[startColumn];
00776 double value = elementU[startColumn];
00777 double largest;
00778 CoinBigIndex put = startColumn;
00779 CoinBigIndex positionLargest = -1;
00780 double thisPivotValue = 0.0;
00781
00782
00783 bool checkLargest;
00784 int mark = markRow[iRow];
00785
00786 if ( mark < 0 ) {
00787 largest = fabs ( value );
00788 positionLargest = put;
00789 put++;
00790 checkLargest = false;
00791 } else {
00792
00793 largest = 0.0;
00794 checkLargest = true;
00795 if ( mark != largeInteger ) {
00796
00797 work[mark] = value;
00798 int word = mark >> COINFACTORIZATION_SHIFT_PER_INT;
00799 int bit = mark & COINFACTORIZATION_MASK_PER_INT;
00800
00801 temp2[word] = temp2[word] | ( 1 << bit );
00802 added--;
00803 } else {
00804 thisPivotValue = value;
00805 }
00806 }
00807 CoinBigIndex i;
00808 for ( i = startColumn + 1; i < endColumn; i++ ) {
00809 iRow = indexRowU[i];
00810 value = elementU[i];
00811 int mark = markRow[iRow];
00812
00813 if ( mark < 0 ) {
00814
00815 indexRowU[put] = iRow;
00816 elementU[put] = value;;
00817 if ( checkLargest ) {
00818 double absValue = fabs ( value );
00819
00820 if ( absValue > largest ) {
00821 largest = absValue;
00822 positionLargest = put;
00823 }
00824 }
00825 put++;
00826 } else if ( mark != largeInteger ) {
00827
00828 work[mark] = value;;
00829 int word = mark >> COINFACTORIZATION_SHIFT_PER_INT;
00830 int bit = mark & COINFACTORIZATION_MASK_PER_INT;
00831
00832 temp2[word] = temp2[word] | ( 1 << bit );
00833 added--;
00834 } else {
00835 thisPivotValue = value;
00836 }
00837 }
00838
00839 elementU[put] = elementU[startColumn];
00840 indexRowU[put] = indexRowU[startColumn];
00841 if ( positionLargest == startColumn ) {
00842 positionLargest = put;
00843 }
00844 put++;
00845 elementU[startColumn] = thisPivotValue;
00846 indexRowU[startColumn] = pivotRow;
00847
00848 startColumn++;
00849 numberInColumn[iColumn] = put - startColumn;
00850 int * numberInColumnPlus = numberInColumnPlus_.array();
00851 numberInColumnPlus[iColumn]++;
00852 startColumnU[iColumn]++;
00853
00854 int next = nextColumn[iColumn];
00855 CoinBigIndex space;
00856
00857 space = startColumnU[next] - put - numberInColumnPlus[next];
00858
00859 if ( numberInPivotColumn > space ) {
00860
00861 if ( !getColumnSpace ( iColumn, numberInPivotColumn ) ) {
00862 return false;
00863 }
00864
00865 positionLargest = positionLargest + startColumnU[iColumn] - startColumn;
00866 startColumn = startColumnU[iColumn];
00867 put = startColumn + numberInColumn[iColumn];
00868 }
00869 double tolerance = zeroTolerance_;
00870
00871 int *nextCount = nextCount_.array();
00872 for ( j = 0; j < numberInPivotColumn; j++ ) {
00873 value = work[j] - thisPivotValue * multipliersL[j];
00874 double absValue = fabs ( value );
00875
00876 if ( absValue > tolerance ) {
00877 work[j] = 0.0;
00878 elementU[put] = value;
00879 indexRowU[put] = indexL[j];
00880 if ( absValue > largest ) {
00881 largest = absValue;
00882 positionLargest = put;
00883 }
00884 put++;
00885 } else {
00886 work[j] = 0.0;
00887 added--;
00888 int word = j >> COINFACTORIZATION_SHIFT_PER_INT;
00889 int bit = j & COINFACTORIZATION_MASK_PER_INT;
00890
00891 if ( temp2[word] & ( 1 << bit ) ) {
00892
00893 iRow = indexL[j];
00894 CoinBigIndex start = startRowU[iRow];
00895 CoinBigIndex end = start + numberInRow[iRow];
00896 CoinBigIndex where = start;
00897
00898 while ( indexColumnU[where] != iColumn ) {
00899 where++;
00900 }
00901 #if DEBUG_COIN
00902 if ( where >= end ) {
00903 abort ( );
00904 }
00905 #endif
00906 indexColumnU[where] = indexColumnU[end - 1];
00907 numberInRow[iRow]--;
00908 } else {
00909
00910 int word = j >> COINFACTORIZATION_SHIFT_PER_INT;
00911 int bit = j & COINFACTORIZATION_MASK_PER_INT;
00912
00913 temp2[word] = temp2[word] | ( 1 << bit );
00914 }
00915 }
00916 }
00917 numberInColumn[iColumn] = put - startColumn;
00918
00919 if ( positionLargest >= 0 ) {
00920 value = elementU[positionLargest];
00921 iRow = indexRowU[positionLargest];
00922 elementU[positionLargest] = elementU[startColumn];
00923 indexRowU[positionLargest] = indexRowU[startColumn];
00924 elementU[startColumn] = value;
00925 indexRowU[startColumn] = iRow;
00926 }
00927
00928 if ( nextCount[iColumn + numberRows_] != -2 ) {
00929
00930 deleteLink ( iColumn + numberRows_ );
00931 addLink ( iColumn + numberRows_, numberInColumn[iColumn] );
00932 }
00933 temp2 += increment2;
00934 }
00935
00936 unsigned int *putBase = workArea2;
00937 int bigLoops = numberInPivotColumn >> COINFACTORIZATION_SHIFT_PER_INT;
00938 int i = 0;
00939
00940
00941 while ( bigLoops ) {
00942 bigLoops--;
00943 int bit;
00944 for ( bit = 0; bit < COINFACTORIZATION_BITS_PER_INT; i++, bit++ ) {
00945 unsigned int *putThis = putBase;
00946 int iRow = indexL[i];
00947
00948
00949 int number = 0;
00950 int jColumn;
00951
00952 for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
00953 unsigned int test = *putThis;
00954
00955 putThis += increment2;
00956 test = 1 - ( ( test >> bit ) & 1 );
00957 number += test;
00958 }
00959 int next = nextRow[iRow];
00960 CoinBigIndex space;
00961
00962 space = startRowU[next] - startRowU[iRow];
00963 number += numberInRow[iRow];
00964 if ( space < number ) {
00965 if ( !getRowSpace ( iRow, number ) ) {
00966 return false;
00967 }
00968 }
00969
00970 putThis = putBase;
00971 next = nextRow[iRow];
00972 number = numberInRow[iRow];
00973 CoinBigIndex end = startRowU[iRow] + number;
00974 int saveIndex = indexColumnU[startRowU[next]];
00975
00976
00977 for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
00978 unsigned int test = *putThis;
00979
00980 putThis += increment2;
00981 test = 1 - ( ( test >> bit ) & 1 );
00982 indexColumnU[end] = saveColumn[jColumn];
00983 end += test;
00984 }
00985
00986 indexColumnU[startRowU[next]] = saveIndex;
00987 markRow[iRow] = -1;
00988 number = end - startRowU[iRow];
00989 numberInRow[iRow] = number;
00990 deleteLink ( iRow );
00991 addLink ( iRow, number );
00992 }
00993 putBase++;
00994 }
00995 int bit;
00996
00997 for ( bit = 0; i < numberInPivotColumn; i++, bit++ ) {
00998 unsigned int *putThis = putBase;
00999 int iRow = indexL[i];
01000
01001
01002 int number = 0;
01003 int jColumn;
01004
01005 for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
01006 unsigned int test = *putThis;
01007
01008 putThis += increment2;
01009 test = 1 - ( ( test >> bit ) & 1 );
01010 number += test;
01011 }
01012 int next = nextRow[iRow];
01013 CoinBigIndex space;
01014
01015 space = startRowU[next] - startRowU[iRow];
01016 number += numberInRow[iRow];
01017 if ( space < number ) {
01018 if ( !getRowSpace ( iRow, number ) ) {
01019 return false;
01020 }
01021 }
01022
01023 putThis = putBase;
01024 next = nextRow[iRow];
01025 number = numberInRow[iRow];
01026 CoinBigIndex end = startRowU[iRow] + number;
01027 int saveIndex;
01028
01029 saveIndex = indexColumnU[startRowU[next]];
01030
01031
01032 for ( jColumn = 0; jColumn < numberInPivotRow; jColumn++ ) {
01033 unsigned int test = *putThis;
01034
01035 putThis += increment2;
01036 test = 1 - ( ( test >> bit ) & 1 );
01037
01038 indexColumnU[end] = saveColumn[jColumn];
01039 end += test;
01040 }
01041 indexColumnU[startRowU[next]] = saveIndex;
01042 markRow[iRow] = -1;
01043 number = end - startRowU[iRow];
01044 numberInRow[iRow] = number;
01045 deleteLink ( iRow );
01046 addLink ( iRow, number );
01047 }
01048 markRow[pivotRow] = -1;
01049
01050 deleteLink ( pivotRow );
01051 deleteLink ( pivotColumn + numberRows_ );
01052 totalElements_ += added;
01053 return true;
01054 }
01055
01056
01058
01059 protected:
01060
01063
01064 double pivotTolerance_;
01066 double zeroTolerance_;
01068 double slackValue_;
01070 double areaFactor_;
01072 double relaxCheck_;
01074 int numberRows_;
01076 int numberRowsExtra_;
01078 int maximumRowsExtra_;
01080 int numberColumns_;
01082 int numberColumnsExtra_;
01084 int maximumColumnsExtra_;
01086 int numberGoodU_;
01088 int numberGoodL_;
01090 int maximumPivots_;
01092 int numberPivots_;
01095 CoinBigIndex totalElements_;
01097 CoinBigIndex factorElements_;
01099 CoinIntArrayWithLength pivotColumn_;
01101 CoinIntArrayWithLength permute_;
01103 CoinIntArrayWithLength permuteBack_;
01105 CoinIntArrayWithLength pivotColumnBack_;
01107 int status_;
01108
01113
01114
01116 int numberTrials_;
01118 CoinBigIndexArrayWithLength startRowU_;
01119
01121 CoinIntArrayWithLength numberInRow_;
01122
01124 CoinIntArrayWithLength numberInColumn_;
01125
01127 CoinIntArrayWithLength numberInColumnPlus_;
01128
01131 CoinIntArrayWithLength firstCount_;
01132
01134 CoinIntArrayWithLength nextCount_;
01135
01137 CoinIntArrayWithLength lastCount_;
01138
01140 CoinIntArrayWithLength nextColumn_;
01141
01143 CoinIntArrayWithLength lastColumn_;
01144
01146 CoinIntArrayWithLength nextRow_;
01147
01149 CoinIntArrayWithLength lastRow_;
01150
01152 CoinIntArrayWithLength saveColumn_;
01153
01155 CoinIntArrayWithLength markRow_;
01156
01158 int messageLevel_;
01159
01161 int biggerDimension_;
01162
01164 CoinIntArrayWithLength indexColumnU_;
01165
01167 CoinIntArrayWithLength pivotRowL_;
01168
01170 CoinDoubleArrayWithLength pivotRegion_;
01171
01173 int numberSlacks_;
01174
01176 int numberU_;
01177
01179 CoinBigIndex maximumU_;
01180
01182
01183
01185 CoinBigIndex lengthU_;
01186
01188 CoinBigIndex lengthAreaU_;
01189
01191 CoinDoubleArrayWithLength elementU_;
01192
01194 CoinIntArrayWithLength indexRowU_;
01195
01197 CoinBigIndexArrayWithLength startColumnU_;
01198
01200 CoinBigIndexArrayWithLength convertRowToColumnU_;
01201
01203 int numberL_;
01204
01206 int baseL_;
01207
01209 CoinBigIndex lengthL_;
01210
01212 CoinBigIndex lengthAreaL_;
01213
01215 CoinDoubleArrayWithLength elementL_;
01216
01218 CoinIntArrayWithLength indexRowL_;
01219
01221 CoinBigIndexArrayWithLength startColumnL_;
01222
01224 bool doForrestTomlin_;
01225
01227 int numberR_;
01228
01230 CoinBigIndex lengthR_;
01231
01233 CoinBigIndex lengthAreaR_;
01234
01236 double *elementR_;
01237
01239 int *indexRowR_;
01240
01242 CoinBigIndexArrayWithLength startColumnR_;
01243
01245 double * denseArea_;
01246
01248 int * densePermute_;
01249
01251 int numberDense_;
01252
01254 int denseThreshold_;
01255
01257 CoinDoubleArrayWithLength workArea_;
01258
01260 CoinUnsignedIntArrayWithLength workArea2_;
01261
01263 CoinBigIndex numberCompressions_;
01264
01266 mutable double ftranCountInput_;
01267 mutable double ftranCountAfterL_;
01268 mutable double ftranCountAfterR_;
01269 mutable double ftranCountAfterU_;
01270 mutable double btranCountInput_;
01271 mutable double btranCountAfterU_;
01272 mutable double btranCountAfterR_;
01273 mutable double btranCountAfterL_;
01274
01276 mutable int numberFtranCounts_;
01277 mutable int numberBtranCounts_;
01278
01280 double ftranAverageAfterL_;
01281 double ftranAverageAfterR_;
01282 double ftranAverageAfterU_;
01283 double btranAverageAfterU_;
01284 double btranAverageAfterR_;
01285 double btranAverageAfterL_;
01286
01288 mutable bool collectStatistics_;
01289
01291 int sparseThreshold_;
01292
01294 int sparseThreshold2_;
01295
01297 CoinBigIndexArrayWithLength startRowL_;
01298
01300 CoinIntArrayWithLength indexColumnL_;
01301
01303 CoinDoubleArrayWithLength elementByRowL_;
01304
01306 mutable CoinIntArrayWithLength sparse_;
01310 int biasLU_;
01316 int persistenceFlag_;
01318 };
01319
01320 #ifdef COIN_HAS_LAPACK
01321 #define DENSE_CODE 1
01322
01323 #ifndef ipfint
01324
01325 typedef int ipfint;
01326 typedef const int cipfint;
01327 #endif
01328 #endif
01329 #endif