source: trunk/Clp/src/AbcMatrix.hpp @ 2385

Last change on this file since 2385 was 2385, checked in by unxusr, 4 months ago

formatting

  • Property svn:keywords set to Id
File size: 21.0 KB
Line 
1/* $Id: AbcMatrix.hpp 2385 2019-01-06 19:43:06Z unxusr $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others, Copyright (C) 2012, FasterCoin.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6#ifndef AbcMatrix_H
7#define AbcMatrix_H
8
9#include "CoinPragma.hpp"
10
11#include "ClpMatrixBase.hpp"
12#include "AbcSimplex.hpp"
13#include "CoinAbcHelperFunctions.hpp"
14/** This implements a scaled version of CoinPackedMatrix
15    It may have THREE! copies
16    1) scaled CoinPackedMatrix without gaps
17    2) row copy non-basic,basic, fixed
18    3) vector copy
19*/
20class AbcMatrix2;
21class AbcMatrix3;
22class AbcMatrix {
23
24public:
25  /**@name Useful methods */
26  //@{
27  /// Return a complete CoinPackedMatrix
28  inline CoinPackedMatrix *getPackedMatrix() const
29  {
30    return matrix_;
31  }
32  /** Whether the packed matrix is column major ordered or not. */
33  inline bool isColOrdered() const
34  {
35    return true;
36  }
37  /** Number of entries in the packed matrix. */
38  inline CoinBigIndex getNumElements() const
39  {
40    return matrix_->getNumElements();
41  }
42  /** Number of columns. */
43  inline int getNumCols() const
44  {
45    assert(matrix_->getNumCols() == model_->numberColumns());
46    return matrix_->getNumCols();
47  }
48  /** Number of rows. */
49  inline int getNumRows() const
50  {
51    assert(matrix_->getNumRows() == model_->numberRows());
52    return matrix_->getNumRows();
53  }
54  /// Sets model
55  void setModel(AbcSimplex *model);
56  /// A vector containing the elements in the packed matrix.
57  inline const double *getElements() const
58  {
59    return matrix_->getElements();
60  }
61  /// Mutable elements
62  inline double *getMutableElements() const
63  {
64    return matrix_->getMutableElements();
65  }
66  /// A vector containing the minor indices of the elements in the packed matrix.
67  inline const int *getIndices() const
68  {
69    return matrix_->getIndices();
70  }
71  /// A vector containing the minor indices of the elements in the packed matrix.
72  inline int *getMutableIndices() const
73  {
74    return matrix_->getMutableIndices();
75  }
76  /// Starts
77  inline const CoinBigIndex *getVectorStarts() const
78  {
79    return matrix_->getVectorStarts();
80  }
81  inline CoinBigIndex *getMutableVectorStarts() const
82  {
83    return matrix_->getMutableVectorStarts();
84  }
85  /** The lengths of the major-dimension vectors. */
86  inline const int *getVectorLengths() const
87  {
88    return matrix_->getVectorLengths();
89  }
90  /** The lengths of the major-dimension vectors. */
91  inline int *getMutableVectorLengths() const
92  {
93    return matrix_->getMutableVectorLengths();
94  }
95  /// Row starts
96  CoinBigIndex *rowStart() const;
97  /// Row ends
98  CoinBigIndex *rowEnd() const;
99  /// Row elements
100  double *rowElements() const;
101  /// Row columns
102  CoinSimplexInt *rowColumns() const;
103  /** Returns a new matrix in reverse order without gaps */
104  CoinPackedMatrix *reverseOrderedCopy() const;
105  /// Returns number of elements in column part of basis
106  CoinBigIndex countBasis(const int *whichColumn,
107    int &numberColumnBasic);
108  /// Fills in column part of basis
109  void fillBasis(const int *whichColumn,
110    int &numberColumnBasic,
111    int *row, int *start,
112    int *rowCount, int *columnCount,
113    CoinSimplexDouble *element);
114  /// Fills in column part of basis
115  void fillBasis(const int *whichColumn,
116    int &numberColumnBasic,
117    int *row, int *start,
118    int *rowCount, int *columnCount,
119    long double *element);
120  /** Scales and creates row copy
121  */
122  void scale(int numberRowsAlreadyScaled);
123  /// Creates row copy
124  void createRowCopy();
125  /// Take out of useful
126  void takeOutOfUseful(int sequence, CoinIndexedVector &spare);
127  /// Put into useful
128  void putIntofUseful(int sequence, CoinIndexedVector &spare);
129  /// Put in and out for useful
130  void inOutUseful(int sequenceIn, int sequenceOut);
131  /// Make all useful
132  void makeAllUseful(CoinIndexedVector &spare);
133  /// Sort into useful
134  void sortUseful(CoinIndexedVector &spare);
135  /// Move largest in column to beginning (not used as doesn't help factorization)
136  void moveLargestToStart();
137
138  /** Unpacks a column into an CoinIndexedVector
139   */
140  void unpack(CoinIndexedVector &rowArray,
141    int column) const;
142  /** Adds multiple of a column (or slack) into an CoinIndexedvector
143      You can use quickAdd to add to vector */
144  void add(CoinIndexedVector &rowArray, int column, double multiplier) const;
145  //@}
146
147  /**@name Matrix times vector methods */
148  //@{
149  /** Return <code>y + A * scalar *x</code> in <code>y</code>.
150      @pre <code>x</code> must be of size <code>numColumns()</code>
151      @pre <code>y</code> must be of size <code>numRows()</code> */
152  void timesModifyExcludingSlacks(double scalar,
153    const double *x, double *y) const;
154  /** Return <code>y + A * scalar(+-1) *x</code> in <code>y</code>.
155      @pre <code>x</code> must be of size <code>numColumns()+numRows()</code>
156      @pre <code>y</code> must be of size <code>numRows()</code> */
157  void timesModifyIncludingSlacks(double scalar,
158    const double *x, double *y) const;
159  /** Return <code>A * scalar(+-1) *x</code> in <code>y</code>.
160      @pre <code>x</code> must be of size <code>numColumns()+numRows()</code>
161      @pre <code>y</code> must be of size <code>numRows()</code> */
162  void timesIncludingSlacks(double scalar,
163    const double *x, double *y) const;
164  /** Return A * scalar(+-1) *x + y</code> in <code>y</code>.
165      @pre <code>x</code> must be of size <code>numRows()</code>
166      @pre <code>y</code> must be of size <code>numRows()+numColumns()</code> */
167  void transposeTimesNonBasic(double scalar,
168    const double *x, double *y) const;
169  /** Return y - A * x</code> in <code>y</code>.
170      @pre <code>x</code> must be of size <code>numRows()</code>
171      @pre <code>y</code> must be of size <code>numRows()+numColumns()</code> */
172  void transposeTimesAll(const double *x, double *y) const;
173  /** Return y + A * scalar(+-1) *x</code> in <code>y</code>.
174      @pre <code>x</code> must be of size <code>numRows()</code>
175      @pre <code>y</code> must be of size <code>numRows()</code> */
176  void transposeTimesBasic(double scalar,
177    const double *x, double *y) const;
178  /** Return <code>x * scalar * A/code> in <code>z</code>.
179      Note - x unpacked mode - z packed mode including slacks
180      All these return atLo/atUp first then free/superbasic
181      number of first set returned
182      pivotVariable is extended to have that order
183      reversePivotVariable used to update that list
184      free/superbasic only stored in normal format
185      can use spare array to get this effect
186      may put djs alongside atLo/atUp
187      Squashes small elements and knows about AbcSimplex */
188  int transposeTimesNonBasic(double scalar,
189    const CoinIndexedVector &x,
190    CoinIndexedVector &z) const;
191  /// gets sorted tableau row and a possible value of theta
192  double dualColumn1(const CoinIndexedVector &update,
193    CoinPartitionedVector &tableauRow,
194    CoinPartitionedVector &candidateList) const;
195  /// gets sorted tableau row and a possible value of theta
196  double dualColumn1Row(int iBlock, double upperThetaSlack, int &freeSequence,
197    const CoinIndexedVector &update,
198    CoinPartitionedVector &tableauRow,
199    CoinPartitionedVector &candidateList) const;
200  /// gets sorted tableau row and a possible value of theta
201  double dualColumn1RowFew(int iBlock, double upperThetaSlack, int &freeSequence,
202    const CoinIndexedVector &update,
203    CoinPartitionedVector &tableauRow,
204    CoinPartitionedVector &candidateList) const;
205  /// gets sorted tableau row and a possible value of theta
206  double dualColumn1Row2(double upperThetaSlack, int &freeSequence,
207    const CoinIndexedVector &update,
208    CoinPartitionedVector &tableauRow,
209    CoinPartitionedVector &candidateList) const;
210  /// gets sorted tableau row and a possible value of theta
211  double dualColumn1Row1(double upperThetaSlack, int &freeSequence,
212    const CoinIndexedVector &update,
213    CoinPartitionedVector &tableauRow,
214    CoinPartitionedVector &candidateList) const;
215  /** gets sorted tableau row and a possible value of theta
216      On input first,last give what to scan
217      On output is number in tableauRow and candidateList */
218  void dualColumn1Part(int iBlock, int &sequenceIn, double &upperTheta,
219    const CoinIndexedVector &update,
220    CoinPartitionedVector &tableauRow,
221    CoinPartitionedVector &candidateList) const;
222  /// rebalance for parallel
223  void rebalance() const;
224  /// Get sequenceIn when Dantzig
225  int pivotColumnDantzig(const CoinIndexedVector &updates,
226    CoinPartitionedVector &spare) const;
227  /// Get sequenceIn when Dantzig (One block)
228  int pivotColumnDantzig(int iBlock, bool doByRow, const CoinIndexedVector &updates,
229    CoinPartitionedVector &spare,
230    double &bestValue) const;
231  /// gets tableau row - returns number of slacks in block
232  int primalColumnRow(int iBlock, bool doByRow, const CoinIndexedVector &update,
233    CoinPartitionedVector &tableauRow) const;
234  /// gets tableau row and dj row - returns number of slacks in block
235  int primalColumnRowAndDjs(int iBlock, const CoinIndexedVector &updateTableau,
236    const CoinIndexedVector &updateDjs,
237    CoinPartitionedVector &tableauRow) const;
238  /** Chooses best weighted dj
239   */
240  int chooseBestDj(int iBlock, const CoinIndexedVector &infeasibilities,
241    const double *weights) const;
242  /** does steepest edge double or triple update
243      If scaleFactor!=0 then use with tableau row to update djs
244      otherwise use updateForDjs
245      Returns best sequence
246   */
247  int primalColumnDouble(int iBlock, CoinPartitionedVector &updateForTableauRow,
248    CoinPartitionedVector &updateForDjs,
249    const CoinIndexedVector &updateForWeights,
250    CoinPartitionedVector &spareColumn1,
251    double *infeasibilities,
252    double referenceIn, double devex,
253    // Array for exact devex to say what is in reference framework
254    unsigned int *reference,
255    double *weights, double scaleFactor) const;
256  /** does steepest edge double or triple update
257      If scaleFactor!=0 then use with tableau row to update djs
258      otherwise use updateForDjs
259      Returns best sequence
260   */
261  int primalColumnSparseDouble(int iBlock, CoinPartitionedVector &updateForTableauRow,
262    CoinPartitionedVector &updateForDjs,
263    const CoinIndexedVector &updateForWeights,
264    CoinPartitionedVector &spareColumn1,
265    double *infeasibilities,
266    double referenceIn, double devex,
267    // Array for exact devex to say what is in reference framework
268    unsigned int *reference,
269    double *weights, double scaleFactor) const;
270  /** does steepest edge double or triple update
271      If scaleFactor!=0 then use with tableau row to update djs
272      otherwise use updateForDjs
273      Returns best sequence
274   */
275  int primalColumnDouble(CoinPartitionedVector &updateForTableauRow,
276    CoinPartitionedVector &updateForDjs,
277    const CoinIndexedVector &updateForWeights,
278    CoinPartitionedVector &spareColumn1,
279    CoinIndexedVector &infeasible,
280    double referenceIn, double devex,
281    // Array for exact devex to say what is in reference framework
282    unsigned int *reference,
283    double *weights, double scaleFactor) const;
284  /// gets subset updates
285  void primalColumnSubset(int iBlock, const CoinIndexedVector &update,
286    const CoinPartitionedVector &tableauRow,
287    CoinPartitionedVector &weights) const;
288  /// Partial pricing
289  void partialPricing(double startFraction, double endFraction,
290    int &bestSequence, int &numberWanted);
291  /** Return <code>x *A</code> in <code>z</code> but
292      just for indices Already in z.
293      Note - z always packed mode */
294  void subsetTransposeTimes(const CoinIndexedVector &x,
295    CoinIndexedVector &z) const;
296  /// Return <code>-x *A</code> in <code>z</code>
297  void transposeTimes(const CoinIndexedVector &x,
298    CoinIndexedVector &z) const;
299  //@}
300
301  /**@name Other */
302  //@{
303  /// Returns CoinPackedMatrix (non const)
304  inline CoinPackedMatrix *matrix() const
305  {
306    return matrix_;
307  }
308  /** Partial pricing tuning parameter - minimum number of "objects" to scan.
309      e.g. number of Gub sets but could be number of variables */
310  inline int minimumObjectsScan() const
311  {
312    return minimumObjectsScan_;
313  }
314  inline void setMinimumObjectsScan(int value)
315  {
316    minimumObjectsScan_ = value;
317  }
318  /// Partial pricing tuning parameter - minimum number of negative reduced costs to get
319  inline int minimumGoodReducedCosts() const
320  {
321    return minimumGoodReducedCosts_;
322  }
323  inline void setMinimumGoodReducedCosts(int value)
324  {
325    minimumGoodReducedCosts_ = value;
326  }
327  /// Current start of search space in matrix (as fraction)
328  inline double startFraction() const
329  {
330    return startFraction_;
331  }
332  inline void setStartFraction(double value)
333  {
334    startFraction_ = value;
335  }
336  /// Current end of search space in matrix (as fraction)
337  inline double endFraction() const
338  {
339    return endFraction_;
340  }
341  inline void setEndFraction(double value)
342  {
343    endFraction_ = value;
344  }
345  /// Current best reduced cost
346  inline double savedBestDj() const
347  {
348    return savedBestDj_;
349  }
350  inline void setSavedBestDj(double value)
351  {
352    savedBestDj_ = value;
353  }
354  /// Initial number of negative reduced costs wanted
355  inline int originalWanted() const
356  {
357    return originalWanted_;
358  }
359  inline void setOriginalWanted(int value)
360  {
361    originalWanted_ = value;
362  }
363  /// Current number of negative reduced costs which we still need
364  inline int currentWanted() const
365  {
366    return currentWanted_;
367  }
368  inline void setCurrentWanted(int value)
369  {
370    currentWanted_ = value;
371  }
372  /// Current best sequence
373  inline int savedBestSequence() const
374  {
375    return savedBestSequence_;
376  }
377  inline void setSavedBestSequence(int value)
378  {
379    savedBestSequence_ = value;
380  }
381  /// Start of each column block
382  inline int *startColumnBlock() const
383  {
384    return startColumnBlock_;
385  }
386  /// Start of each block (in stored)
387  inline const int *blockStart() const
388  {
389    return blockStart_;
390  }
391  inline bool gotRowCopy() const
392  {
393    return rowStart_ != 0;
394  }
395  /// Start of each block (in stored)
396  inline int blockStart(int block) const
397  {
398    return blockStart_[block];
399  }
400  /// Number of actual column blocks
401  inline int numberColumnBlocks() const
402  {
403    return numberColumnBlocks_;
404  }
405  /// Number of actual row blocks
406  inline int numberRowBlocks() const
407  {
408    return numberRowBlocks_;
409  }
410  //@}
411
412  /**@name Constructors, destructor */
413  //@{
414  /** Default constructor. */
415  AbcMatrix();
416  /** Destructor */
417  ~AbcMatrix();
418  //@}
419
420  /**@name Copy method */
421  //@{
422  /** The copy constructor. */
423  AbcMatrix(const AbcMatrix &);
424  /** The copy constructor from an CoinPackedMatrix. */
425  AbcMatrix(const CoinPackedMatrix &);
426  /** Subset constructor (without gaps).  Duplicates are allowed
427      and order is as given */
428  AbcMatrix(const AbcMatrix &wholeModel,
429    int numberRows, const int *whichRows,
430    int numberColumns, const int *whichColumns);
431  AbcMatrix(const CoinPackedMatrix &wholeModel,
432    int numberRows, const int *whichRows,
433    int numberColumns, const int *whichColumns);
434
435  AbcMatrix &operator=(const AbcMatrix &);
436  /// Copy contents - resizing if necessary - otherwise re-use memory
437  void copy(const AbcMatrix *from);
438  //@}
439private:
440protected:
441  /**@name Data members
442     The data members are protected to allow access for derived classes. */
443  //@{
444  /// Data
445  CoinPackedMatrix *matrix_;
446  /// Model
447  mutable AbcSimplex *model_;
448#if ABC_PARALLEL == 0
449#define NUMBER_ROW_BLOCKS 1
450#define NUMBER_COLUMN_BLOCKS 1
451#elif ABC_PARALLEL == 1
452#define NUMBER_ROW_BLOCKS 4
453#define NUMBER_COLUMN_BLOCKS 4
454#else
455#define NUMBER_ROW_BLOCKS 8
456#define NUMBER_COLUMN_BLOCKS 8
457#endif
458  /** Start of each row (per block) - last lot are useless
459      first all row starts for block 0, then for block2
460      so NUMBER_ROW_BLOCKS+2 times number rows */
461  CoinBigIndex *rowStart_;
462  /// Values by row
463  double *element_;
464  /// Columns
465  int *column_;
466  /// Start of each column block
467  mutable int startColumnBlock_[NUMBER_COLUMN_BLOCKS + 1];
468  /// Start of each block (in stored)
469  int blockStart_[NUMBER_ROW_BLOCKS + 1];
470  /// Number of actual column blocks
471  mutable int numberColumnBlocks_;
472  /// Number of actual row blocks
473  int numberRowBlocks_;
474  //#define COUNT_COPY
475#ifdef COUNT_COPY
476#define MAX_COUNT 13
477  /// Start in elements etc
478  CoinBigIndex countStart_[MAX_COUNT + 1];
479  /// First column
480  int countFirst_[MAX_COUNT + 1];
481  // later int countEndUseful_[MAX_COUNT+1];
482  int *countRealColumn_;
483  // later int * countInverseRealColumn_;
484  CoinBigIndex *countStartLarge_;
485  int *countRow_;
486  double *countElement_;
487  int smallestCount_;
488  int largestCount_;
489#endif
490  /// Special row copy
491  //AbcMatrix2 * rowCopy_;
492  /// Special column copy
493  //AbcMatrix3 * columnCopy_;
494  /// Current start of search space in matrix (as fraction)
495  double startFraction_;
496  /// Current end of search space in matrix (as fraction)
497  double endFraction_;
498  /// Best reduced cost so far
499  double savedBestDj_;
500  /// Initial number of negative reduced costs wanted
501  int originalWanted_;
502  /// Current number of negative reduced costs which we still need
503  int currentWanted_;
504  /// Saved best sequence in pricing
505  int savedBestSequence_;
506  /// Partial pricing tuning parameter - minimum number of "objects" to scan
507  int minimumObjectsScan_;
508  /// Partial pricing tuning parameter - minimum number of negative reduced costs to get
509  int minimumGoodReducedCosts_;
510  //@}
511};
512#ifdef THREAD
513#include <pthread.h>
514typedef struct {
515  double acceptablePivot;
516  const AbcSimplex *model;
517  double *spare;
518  int *spareIndex;
519  double *arrayTemp;
520  int *indexTemp;
521  int *numberInPtr;
522  double *bestPossiblePtr;
523  double *upperThetaPtr;
524  int *posFreePtr;
525  double *freePivotPtr;
526  int *numberOutPtr;
527  const unsigned short *count;
528  const double *pi;
529  const CoinBigIndex *rowStart;
530  const double *element;
531  const unsigned short *column;
532  int offset;
533  int numberInRowArray;
534  int numberLook;
535} dualColumn0Struct;
536#endif
537class AbcMatrix2 {
538
539public:
540  /**@name Useful methods */
541  //@{
542  /** Return <code>x * -1 * A in <code>z</code>.
543      Note - x packed and z will be packed mode
544      Squashes small elements and knows about AbcSimplex */
545  void transposeTimes(const AbcSimplex *model,
546    const CoinPackedMatrix *rowCopy,
547    const CoinIndexedVector &x,
548    CoinIndexedVector &spareArray,
549    CoinIndexedVector &z) const;
550  /// Returns true if copy has useful information
551  inline bool usefulInfo() const
552  {
553    return rowStart_ != NULL;
554  }
555  //@}
556
557  /**@name Constructors, destructor */
558  //@{
559  /** Default constructor. */
560  AbcMatrix2();
561  /** Constructor from copy. */
562  AbcMatrix2(AbcSimplex *model, const CoinPackedMatrix *rowCopy);
563  /** Destructor */
564  ~AbcMatrix2();
565  //@}
566
567  /**@name Copy method */
568  //@{
569  /** The copy constructor. */
570  AbcMatrix2(const AbcMatrix2 &);
571  AbcMatrix2 &operator=(const AbcMatrix2 &);
572  //@}
573
574protected:
575  /**@name Data members
576     The data members are protected to allow access for derived classes. */
577  //@{
578  /// Number of blocks
579  int numberBlocks_;
580  /// Number of rows
581  int numberRows_;
582  /// Column offset for each block (plus one at end)
583  int *offset_;
584  /// Counts of elements in each part of row
585  mutable unsigned short *count_;
586  /// Row starts
587  mutable CoinBigIndex *rowStart_;
588  /// columns within block
589  unsigned short *column_;
590  /// work arrays
591  double *work_;
592#ifdef THREAD
593  pthread_t *threadId_;
594  dualColumn0Struct *info_;
595#endif
596  //@}
597};
598typedef struct {
599  CoinBigIndex startElements_; // point to data
600  int startIndices_; // point to column_
601  int numberInBlock_;
602  int numberPrice_; // at beginning
603  int numberElements_; // number elements per column
604} blockStruct3;
605class AbcMatrix3 {
606
607public:
608  /**@name Useful methods */
609  //@{
610  /** Return <code>x * -1 * A in <code>z</code>.
611      Note - x packed and z will be packed mode
612      Squashes small elements and knows about AbcSimplex */
613  void transposeTimes(const AbcSimplex *model,
614    const double *pi,
615    CoinIndexedVector &output) const;
616  /// Updates two arrays for steepest
617  void transposeTimes2(const AbcSimplex *model,
618    const double *pi, CoinIndexedVector &dj1,
619    const double *piWeight,
620    double referenceIn, double devex,
621    // Array for exact devex to say what is in reference framework
622    unsigned int *reference,
623    double *weights, double scaleFactor);
624  //@}
625
626  /**@name Constructors, destructor */
627  //@{
628  /** Default constructor. */
629  AbcMatrix3();
630  /** Constructor from copy. */
631  AbcMatrix3(AbcSimplex *model, const CoinPackedMatrix *columnCopy);
632  /** Destructor */
633  ~AbcMatrix3();
634  //@}
635
636  /**@name Copy method */
637  //@{
638  /** The copy constructor. */
639  AbcMatrix3(const AbcMatrix3 &);
640  AbcMatrix3 &operator=(const AbcMatrix3 &);
641  //@}
642  /**@name Sort methods */
643  //@{
644  /** Sort blocks */
645  void sortBlocks(const AbcSimplex *model);
646  /// Swap one variable
647  void swapOne(const AbcSimplex *model, const AbcMatrix *matrix,
648    int iColumn);
649  //@}
650
651protected:
652  /**@name Data members
653     The data members are protected to allow access for derived classes. */
654  //@{
655  /// Number of blocks
656  int numberBlocks_;
657  /// Number of columns
658  int numberColumns_;
659  /// Column indices and reverse lookup (within block)
660  int *column_;
661  /// Starts for odd/long vectors
662  CoinBigIndex *start_;
663  /// Rows
664  int *row_;
665  /// Elements
666  double *element_;
667  /// Blocks (ordinary start at 0 and go to first block)
668  blockStruct *block_;
669  //@}
670};
671
672#endif
673
674/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
675*/
Note: See TracBrowser for help on using the repository browser.