HDK
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
UT_SparseMatrix.h
Go to the documentation of this file.
1 /*
2  * PROPRIETARY INFORMATION. This software is proprietary to
3  * Side Effects Software Inc., and is not to be reproduced,
4  * transmitted, or disclosed in any way without written permission.
5  *
6  * NAME: Sparse Matrix Class (C++)
7  *
8  * COMMENTS:
9  *
10  * Storage format: Co-ordinate (COO). See
11  * Y. Saad. SPARSKIT: A basic toolkit for sparse matrix computations.
12  * Research Institute for Advanced Computer Science, RIACS-90-20,
13  * 1990.
14  * for details.
15  */
16 
17 #ifndef __UT_SparseMatrix_H__
18 #define __UT_SparseMatrix_H__
19 
20 #include "UT_API.h"
21 
22 #include "UT_IntArray.h"
23 #include "UT_ThreadedAlgorithm.h"
24 #include "UT_Vector.h"
25 #include "UT_Interrupt.h"
26 #include "UT_ParallelUtil.h"
27 
28 #include "UT_VectorTypes.h"
29 
30 #include <iterator>
31 
32 class UT_IStream;
33 class UT_BitArray;
34 
35 template <typename T> class UT_API UT_SparseMatrixRowT;
36 
37 template <typename T, bool IsPaged>
39 {
40  class ut_MatrixCell
41  {
42  public:
43  int myRow;
44  int myCol;
45  T myValue;
46 
47  inline bool operator<(const ut_MatrixCell &o) const
48  {
49  if (myRow < o.myRow)
50  return true;
51  if (myRow == o.myRow && myCol < o.myCol)
52  return true;
53  return false;
54  }
55  };
56 
57  // For paged matrices, store the size of each page
58  // We always allocate full pages.
59  static const int CELL_PAGESIZE = 1024;
60  static const int CELL_PAGEMASK = 1023;
61  static const int CELL_PAGEBITS = 10;
62 
63  class ut_CellIterator : public std::iterator<std::random_access_iterator_tag, ut_MatrixCell>
64  {
65  public:
66  ut_CellIterator(const UT_SparseMatrixT<T, IsPaged> *matrix, int idx)
67  { myMatrix = matrix; myPos = idx; repage(); }
68 
69  ut_CellIterator operator+(ptrdiff_t n) { return ut_CellIterator(myMatrix, myPos+n); }
70  ut_CellIterator operator-(ptrdiff_t n) { return ut_CellIterator(myMatrix, myPos-n); }
71  ut_CellIterator &operator+=(ptrdiff_t n) { myPos += n; repage(); return *this; }
72  ut_CellIterator &operator-=(ptrdiff_t n) { myPos -= n; repage(); return *this; }
73  ut_CellIterator &operator++() { myPos++; myOffset++; myData++; if (myOffset >= CELL_PAGESIZE) repage(); return *this; }
74  ut_CellIterator &operator--() { myPos--; myOffset--; myData--; if (myOffset < 0) repage(); return *this; }
75  ut_CellIterator operator++(int) { ut_CellIterator result = *this; myPos++; myOffset++; myData++; if (myOffset >= CELL_PAGESIZE) repage(); return result; }
76  ut_CellIterator operator--(int) { ut_CellIterator result = *this; myPos--; myOffset--; myData--; if (myOffset < 0) repage(); return result; }
77  int operator-(ut_CellIterator b) { return myPos - b.myPos; }
78  ut_MatrixCell &operator[](ptrdiff_t idx) { return myMatrix->getCell(idx + myPos); }
79 
80  bool operator<(ut_CellIterator b) const { return myPos < b.myPos; }
81  bool operator>(ut_CellIterator b) const { return myPos > b.myPos; }
82  bool operator>=(ut_CellIterator b) const { return myPos >= b.myPos; }
83  bool operator==(ut_CellIterator b) const { return myPos == b.myPos; }
84  bool operator!=(ut_CellIterator b) const { return myPos != b.myPos; }
85 
86  ut_MatrixCell &operator*() const { return *myData; }
87  ut_MatrixCell *operator->() const { return myData; }
88 
89  protected:
90 
91  void repage()
92  { myPage = myPos >> CELL_PAGEBITS; myOffset = myPos & CELL_PAGEMASK;
93  myData = 0;
94  if (myPage < myMatrix->myCellPages.entries())
95  myData = &myMatrix->myCellPages(myPage)[myOffset];
96  }
97 
98  const UT_SparseMatrixT<T, IsPaged> *myMatrix;
99  ut_MatrixCell *myData;
100  ptrdiff_t myPos;
101  int myPage, myOffset;
102  };
103 
104  // For compact compiled matrix cells, specify the number of values to
105  // pack into a single cell.
106  static const int CELLBITS = 2;
107  static const int CELLSIZE = 1 << CELLBITS;
108  static const int CELLMASK = CELLSIZE-1;
109 
110  class ut_4MatrixCell
111  {
112  public:
113  T myValue[CELLSIZE];
114  int myRow;
115  int myCol;
116  } SYS_ALIGN16;
117 
118 public:
120 
121  // creates a valid 'rows' x 'cols' zero matrix
122  UT_SparseMatrixT(int rows, int cols);
124  ~UT_SparseMatrixT();
125 
126  int getNumRows() const { return myNumRows; }
127  int getNumCols() const { return myNumCols; }
128 
129  /// Return the amount of memory used by this array.
130  int64 getMemoryUsage() const;
131  int64 getIdealMemoryUsage() const;
132 
133  /// Ensures at least this number of cells are available, useful
134  /// if you can predict your size.
135  void reserve(int numcells);
136 
137  /// Shrinks the array to exactly fit the number of cells present
138  void shrinkToFit();
139 
140  // sets equal to a 'rows' x 'cols' zero matrix
141  void init(int rows, int cols);
142 
143  // Sets all entries of the matrix to zero.
144  void zero();
145 
146  // Determines if we are large enough to justify mulithreading.
147  // The 5000 was determined from experiments with UT_Vector.
148  bool shouldMultiThread() const
149  {
150 #ifdef CELLBE
151  return false;
152 #else
153  return getNumRows() > 5000;
154 #endif
155  }
156 
157  // does M(row,col) += value
158  bool addToElement(int row, int col, T value);
159 
160  // Returns the index of the cell for which cell->row >= row is true
161  // first. This could be myCount if row is too high.
162  // myCells[result].myRow may be strictly greater than row if it isn't
163  // present.
164  int findCellFromRow(int row) const;
165 
166  // computes M * v
167  inline void multVec(const UT_VectorT<T> &v, UT_VectorT<T> &result) const
168  {
169  compile();
170  multVecInternal(v, result);
171  }
172 
173  /// Do result = result - (M * v)
174  inline void subtractMultVec(const UT_VectorT<T> &v, UT_VectorT<T> &result) const
175  {
176  compile();
177  subtractMultVecInternal(v, result);
178  }
179  inline void subtractMultVecNoThread(const UT_VectorT<T> &v, UT_VectorT<T> &result) const
180  {
181  compile();
182  subtractMultVecInternalNoThread(v, result);
183  }
184 
185  // computes M^T * v
186  void transposeMultVec(const UT_VectorT<T> &v, UT_VectorT<T> &result) const;
187 
188  /// Square of L2-norm of all columns
189  /// Use to form the diagonal part of A^tA of an overdetermined system A.
190  /// Needed for both Jacobi and Gauss-Seidel iterations.
191  void allColNorm2(UT_VectorT<T> &result) const;
192 
193  /// Initializes out to a matrix where the (i,j) elements are
194  /// (rowstart+i, colstart+j) elements of this. The rowend and
195  /// colend represent the exclusive end - [rowstart,rowend) will
196  /// be extracted for numrows = rowend - rowstart
197  /// This function will compile *this, if it wasn't before the call.
198  /// The resulting matrix will be compiled (by construction).
199  void extractSubMatrix(UT_SparseMatrixT<T, IsPaged> &out,
200  int rowstart, int rowend,
201  int colstart, int colend) const;
202  /// Extract submatrix without compiling *this.
203  /// The resulting submatrix won't be compiled either.
204  void extractSubMatrixUncompiled(UT_SparseMatrixT<T, IsPaged> &out,
205  int rowstart, int rowend,
206  int colstart, int colend) const;
207 
208  // Get fixed-size square submatrix
209  // This will trigger compilation
210  void getSmallSquareSubMatrix(T submatrix[],
211  const exint indices[], const int num_indices) const;
212 
213  // For each row, column, and value call visitor(row, column, value)
214  // visitor may be called with the same (row, column) pair multiple times
215  // the sum of all the associated values is the component in the matrix for (row, column)
216  template<typename Visitor>
217  void accept(Visitor& visitor) const
218  {
219  for(exint ci = 0; ci < myCount; ++ci)
220  {
221  ut_MatrixCell *cell = &getCell(ci);
222  visitor(cell->myRow, cell->myCol, cell->myValue);
223  }
224  }
225 
226  /// Clear the rows and columns specified by the given bit field to a row
227  /// from the identity matrix
228  void clearRowsAndColumns(const UT_BitArray &toclear);
229 
230  /// Extract the diagonal vector from a matrix. This does not
231  /// require the matrix to be compiled, nor does it try to compile it.
232  /// For entries (i, j), this method will extract the entries where
233  /// (i-j == idx).
234  void extractDiagonal(UT_VectorT<T> &out, int idx = 0) const;
235 
236  /// Extract everything but the diagonal:
237  /// out(i, j) = 0, if i == j,
238  // *this(i, j), otherwise
239  void extractNondiagonal(UT_SparseMatrixT<T, IsPaged> &out) const;
240 
241  /// Incomplete Cholesky Factorization.
242  /// Does a Cholesky Factorization but does not make any elements non-zero
243  /// The result of this operation is an UPPER triangular matrix G
244  /// such that A = Gt G.
245  /// The input must be a symmetric matrix.
246  /// Note that this factorization is not always stable.
247  /// Returns 0 on success, -1 if semi-definite (diagonal contained zeros
248  /// within the given tolerance) and -2 if not positive definite (diagonal
249  /// contained negative numbers)
250  int incompleteCholeskyFactorization(T tol=1e-5);
251 
252  /// Modified Incomplete Cholesky
253  /// Same as incomplte cholesky, except attempt to account for the
254  /// discarded entries by adjusting the diagonal.
255  /// tau is a tuning constant.
256  int modifiedIncompleteCholesky(T tau = 0.97,
257  T mindiagratio = 0.25,
258  T tol=1e-5);
259 
260  /// Assumes this is a lower triangular matrix. Solves the equation
261  /// A x = b
262  /// If the diagonal of A is zero within tolerance, the corresponding
263  /// x coordinate is zero.
264  /// Returned is the number of artifical zeros places into x. 0 means
265  /// the solution encountered no singularities, 10 would mean 10
266  /// singularities.
267  int solveLowerTriangular(UT_VectorT<T> &x, const UT_VectorT<T> &b,
268  T tol=1e-5) const;
269  int solveUpperTriangular(UT_VectorT<T> &x, const UT_VectorT<T> &b,
270  T tol=1e-5) const;
271  /// Given an upper triangular matrix, solves the lower triangular
272  /// transposed of it and negates the result.
273  int solveLowerTriangularTransposeNegate(UT_VectorT<T> &x,
274  const UT_VectorT<T> &b,
275  T tol=1e-5) const;
276 
277  // Use the conjugate gradient method to solve Ax=b.
278  // NB: 'x' should be initialized with the initial estimate
279  bool solveConjugateGradient(UT_VectorT<T> &x, const UT_VectorT<T> &b,
281  bool (*callback_func)(void *) = 0,
282  void *callback_data = 0, T tol2=1e-5,
283  int max_iters = -1) const;
284 
285  /// Transposes this matrix.
286  THREADED_METHOD(UT_SparseMatrixT, shouldMultiThread(), transpose)
287  void transposePartial(const UT_JobInfo &info);
288 
289  /// Makes this a transposed copy of source. By not working
290  /// in place we can avoid a sort if the source is compiled.
291  /// If source is not compiled, we compile it.
292  void transposeCompiled(const UT_SparseMatrixT<T, IsPaged> &src);
293 
294  /// *this = -*this.
295  THREADED_METHOD(UT_SparseMatrixT, shouldMultiThread(), negate)
296  void negatePartial(const UT_JobInfo &info);
297 
298  UT_SparseMatrixT<T, IsPaged> &operator=(const UT_SparseMatrixT<T, IsPaged> &m);
299  UT_SparseMatrixT<T, IsPaged> &operator*=(T scalar);
300  UT_SparseMatrixT<T, IsPaged> &operator+=(const UT_SparseMatrixT<T, IsPaged> &m);
301 
302  // Writes out the full matrix as an Row X Column array. Likely
303  // not too useful with big matrices...
304  void printFull(std::ostream &os) const;
305 
306  // Writes out all of the non-zero elements.
307  void printSparse(std::ostream &os) const;
308 
309  // Writes out all of the non-zero elements as a Matlab "sparse" expression.
310  void printSparseMatlab(std::ostream &os,
311  const UT_String &varname) const;
312 
313  // Binary save and load of a sparse matrix
314  void save(std::ostream &os) const;
315  void load(UT_IStream &is);
316 
317  /// Reorders the cell array to be sorted by row than column.
318  /// Consolidates any duplicate entries by adding them together.
319  /// Culls any zero entries.
320  /// compile will be a no-op if the matrix's structure has not changed
321  /// since the last invocation.
322  /// While compilation is expensive, it also performs the collapse/cull
323  /// operation, so you may want to invoke it explicitly if you are
324  /// repeatedly gathering similar matrices, say with +=.
325  void compile() const;
326  bool isCompiled() const { return myCompiledFlag; }
327  bool isStillSorted() const { return myStillSortedFlag; }
328 
329  bool hasNan() const;
330  void testForNan() const;
331 
332 private:
333  /// Always recompiles the matrix regardless of the compile flag setting.
334  void forceCompile() const;
335 
336  /// Extract a cell by its index.
337  inline ut_MatrixCell &getCell(int idx) const
338  {
339  if (IsPaged)
340  return myCellPages(idx >> CELL_PAGEBITS)[idx & CELL_PAGEMASK];
341  else
342  return myCells[idx];
343  }
344 
345  THREADED_METHOD2_CONST(UT_SparseMatrixT, shouldMultiThread(),
346  multVecInternal, const UT_VectorT<T> &, v,
347  UT_VectorT<T> &, result)
348  void multVecInternalPartial(const UT_VectorT<T> &v, UT_VectorT<T> &result,
349  const UT_JobInfo &info) const;
350 
351  THREADED_METHOD2_CONST(UT_SparseMatrixT, shouldMultiThread(),
352  subtractMultVecInternal,
353  const UT_VectorT<T>&, v,
354  UT_VectorT<T>&, result)
355  void subtractMultVecInternalPartial(const UT_VectorT<T>& v, UT_VectorT<T>& result,
356  const UT_JobInfo& info) const;
357 
358  /// Compile into 4-value cell data structure
359  void compile4() const;
360 
361  int myNumRows;
362  int myNumCols;
363  mutable UT_ValArray<ut_MatrixCell *> myCellPages;
364  mutable ut_MatrixCell *myCells;
365  mutable int myCount;
366  mutable int *my4RowOffsets;
367  mutable ut_4MatrixCell *my4Cells;
368  mutable int my4Count;
369  int myMaxSize;
370  mutable bool myCompiledFlag;
371  mutable bool myStillSortedFlag;
372  mutable bool my4CompiledFlag;
373 
374  friend class UT_SparseMatrixRowT<T>;
375 };
376 
377 typedef UT_SparseMatrixT<fpreal32, false> UT_SparseMatrixF;
378 typedef UT_SparseMatrixT<fpreal64, false> UT_SparseMatrixD;
379 typedef UT_SparseMatrixT<fpreal64, false> UT_SparseMatrix;
380 
381 
382 ///
383 /// This is a highly specialized varient of the SparseMatrix
384 /// which does not let you change the number or position of cells
385 /// It has, however, compiled fixed row offsets to make it fast
386 /// to look up specific rows in the matrix.
387 ///
388 template <typename T>
390 {
391  class ut_MatrixCell
392  {
393  public:
394  int myCol;
395  T myValue;
396  };
397 
398 public:
400  ~UT_SparseMatrixRowT();
401 
402  // Can't mix paged and non paged, and this is a non-paged algorithm
403  // Destroys the source matrix in the process as we reuse its
404  // myCells array.
405  void buildFrom(UT_SparseMatrixT<T, false> &m,
406  bool invertdiag = false,
407  T tol=1e-5f);
408 
409  int getNumRows() const { return myNumRows; }
410  int getNumCols() const { return myNumCols; }
411 
412  /// Return the amount of memory used by this array.
413  int64 getMemoryUsage() const;
414 
415  // The 5000 was determined from experiments with UT_Vector.
416  bool shouldMultiThread() const
417  {
418  return getNumRows() > 5000;
419  }
420 
421  // Returns the index of the cell for which cell->row >= row is true
422  // first. This could be myCount if row is too high.
423  // myCells[result].myRow may be strictly greater than row if it isn't
424  // present.
425  int findCellFromRow(int row) const
426  { return myRowOffsets(row); }
427 
428  // Get the inverted diagonal.
429  void getDiagonalInv(UT_VectorT<T> &out) const;
430 
431  // computes M * v
432  THREADED_METHOD2_CONST(UT_SparseMatrixRowT, shouldMultiThread(), multVec,
433  const UT_VectorT<T> &, v,
434  UT_VectorT<T> &, result)
435  void multVecPartial(const UT_VectorT<T> &v, UT_VectorT<T> &result,
436  const UT_JobInfo &info) const;
437 
438  // computes M * v and accumulates the dot product of v & result
439  void multVecAndDot(
440  const UT_VectorT<T> &v,
441  UT_VectorT<T> &result,
442  fpreal64 *dotpq) const;
443 
444  /// Assumes this is a lower triangular matrix. Solves the equation
445  /// A x = b
446  /// If the diagonal of A is zero within tolerance, the corresponding
447  /// x coordinate is zero.
448  /// Returned is the number of artifical zeros places into x. 0 means
449  /// the solution encountered no singularities, 10 would mean 10
450  /// singularities.
451  int solveUpperTriangular(UT_VectorT<T> &x, const UT_VectorT<T> &b,
452  T tol=1e-5) const;
453  /// Given an upper triangular matrix, solves the lower triangular
454  /// transposed of it and negates the result.
455  int solveLowerTriangularTransposeNegate(UT_VectorT<T> &x,
456  const UT_VectorT<T> &b,
457  T tol=1e-5) const;
458 
459  /// Solves conjugate gradient using our specialized functions.
460  /// These allow us to perform some normal and dot operations
461  /// while the cache is still hot.
462  /// This matrix is the matrix to solve. The provided GT matrix
463  /// is the upper triangular result of cholesky factoriziation.
464  /// Norm is hardcoded to 2.
465  /// If the GT matrix is null, no preconditioner will be used.
466  float solveConjugateGradient(UT_VectorT<T> &x, const UT_VectorT<T> &b,
467  const UT_SparseMatrixRowT<T> *GT, T tol2=1e-5,
468  int max_iters = -1, int *iterout = NULL) const;
469 private:
470  THREADED_METHOD3_CONST(UT_SparseMatrixRowT, shouldMultiThread(),
471  multVecAndDotInternal,
472  const UT_VectorT<T> &, v,
473  UT_VectorT<T> &, result,
474  fpreal64 *, dotpq)
475  void multVecAndDotInternalPartial(const UT_VectorT<T> &v, UT_VectorT<T> &result,
476  fpreal64 *dotpq,
477  const UT_JobInfo &info) const;
478 
479  static const exint PARALLEL_BLOCK_SIZE = 1024;
480 
481  /// Extract a cell by its index.
482  inline ut_MatrixCell &getCell(int idx) const
483  {
484  return myCells[idx];
485  }
486 
487  int myNumRows;
488  int myNumCols;
489  ut_MatrixCell *myCells;
490  int myCount;
491  UT_IntArray myRowOffsets;
492  UT_VectorT<T> myDiagonal;
493 };
494 
497 
498 
499 
500 /// This is another highly specialized varient of the SparseMatrix
501 /// which allows multithreaded addition of elements to individual rows and
502 /// supports a conjugate gradient solve. This can be useful for solvers which
503 /// need to create and solve large linear systems in parallel.
504 /// It uses the ELL format which is appropriate for matrices with a limited
505 /// but fairly consistent numer of non-zeros per row, for example the result
506 /// of discretizations across regular grids.
507 template <typename T, bool colmajor=false>
509 {
510 public:
512 
513  UT_SparseMatrixELLT(int rows, int nzeros);
514 
515  void init(int rows, int nzeros);
516 
517  void zero();
518 
519  int getNumRows() const { return myRows; }
520 
521  int getNonZerosPerRow() const { return myNonZeros; }
522 
523  const UT_ValArray<T> &getRowValues() const { return myRowVals; }
524  UT_ValArray<T> &getRowValues() { return myRowVals; }
525 
526  const UT_IntArray &getColumns() const { return myColumns; }
527  UT_IntArray &getColumns() { return myColumns; }
528 
529  size_t index(int row, int nz) const
530  {
531  if (colmajor)
532  return nz * myRows + row;
533  return row * myNonZeros + nz;
534  }
535 
536  bool appendRowElement(int row, int col, T val, int &rowidx)
537  {
538  UT_ASSERT(rowidx < myNonZeros);
539  size_t idx = index(row, rowidx);
540  myRowVals(idx) = val;
541  myColumns(idx) = col;
542  rowidx++;
543  return true;
544  }
545 
546  // Sort rows by column index.
547  THREADED_METHOD(UT_SparseMatrixELLT, shouldMultiThread(),
548  sortRows)
549  void sortRowsPartial(const UT_JobInfo &info);
550 
551  /// Every reference to given rows is removed from the diagonal.
552  /// Because searching the whole matrix is expensive, assumes it
553  /// is a standard fluid matrix, so for each specified row, it
554  /// looks at all columns that are not itself and removes those.
555  void removeEffectOfRows(const UT_ExintArray &rows);
556 
557  bool shouldMultiThread() const
558  {
559  return getNumRows() > 5000;
560  }
561 
562  /// Copy to the supplied SparseMatrix.
563  void copyToSparse(UT_SparseMatrixT<T, false> &A) const;
564 
565  // Extract the inverted diagonal
566  THREADED_METHOD1_CONST(UT_SparseMatrixELLT, shouldMultiThread(),
567  getDiagonalInv, UT_VectorT<T> &, out)
568 
569  /// Get the inverse of the diagonal. Required for the conjugate gradient
570  /// solve when using Jacobi pre-conditioning.
571  void getDiagonalInvPartial(UT_VectorT<T> &out,
572  const UT_JobInfo &info) const;
573 
574  // computes M * v
575  THREADED_METHOD2_CONST(UT_SparseMatrixELLT, shouldMultiThread(),
576  multVec, const UT_VectorT<T> &, v,
577  UT_VectorT<T> &, result)
578 
579  void multVecPartial(const UT_VectorT<T> &v, UT_VectorT<T> &result,
580  const UT_JobInfo &info) const;
581 
582  // computes M * v and accumulates the dot product of v & result
583  void multVecAndDot(
584  const UT_VectorT<T> &v,
585  UT_VectorT<T> &result,
586  fpreal64 *dotpq) const;
587 
588  // computes M * v and accumulates the dot product of v & result
589  void multVecAndDotUpTo(
590  const UT_VectorT<T> &v,
591  UT_VectorT<T> &result,
592  fpreal64 *dotpq,
593  exint solverbase) const;
594 
595  /// Solves conjugate gradient using our specialized functions.
596  /// These allow us to perform some normal and dot operations
597  /// while the cache is still hot.
598  /// This matrix is the matrix to solve. The provided GT matrix
599  /// is the upper triangular result of cholesky factoriziation.
600  /// Norm is hardcoded to 2.
601  /// If the GT matrix is null, no preconditioner will be used.
602  float solveConjugateGradient(UT_VectorT<T> &x, const UT_VectorT<T> &b,
603  const UT_SparseMatrixRowT<T> *GT, T tol2=1e-5,
604  int max_iters = -1, int *iterout = NULL) const;
605 
606  // Writes out all of the non-zero elements as a Matlab "sparse" expression.
607  void printSparseMatlab(std::ostream &os,
608  const UT_String &varname) const;
609 
610  // Writes out all of the non-zero elements in Matrix Market format.
611  void printSparseMatrixMarket(std::ostream &os) const;
612 private:
613  THREADED_METHOD3_CONST(UT_SparseMatrixELLT, shouldMultiThread(),
614  multVecAndDotInternal,
615  const UT_VectorT<T> &, v,
616  UT_VectorT<T> &, result,
617  fpreal64 *, dotpq)
618  void multVecAndDotInternalPartial(const UT_VectorT<T> &v, UT_VectorT<T> &result,
619  fpreal64 *dotpq,
620  const UT_JobInfo &info) const;
621  THREADED_METHOD4_CONST(UT_SparseMatrixELLT, shouldMultiThread(),
622  multVecAndDotUpToInternal,
623  const UT_VectorT<T> &, v,
624  UT_VectorT<T> &, result,
625  fpreal64 *, dotpq,
626  exint, solverbase)
627  void multVecAndDotUpToInternalPartial(const UT_VectorT<T> &v, UT_VectorT<T> &result,
628  fpreal64 *dotpq, exint solverbase,
629  const UT_JobInfo &info) const;
630 
631  static const exint PARALLEL_BLOCK_SIZE = 1024;
632 
633 private:
634  UT_ValArray<T> myRowVals;
635  UT_IntArray myColumns;
636  int myRows;
637  int myNonZeros;
638 };
639 
640 typedef UT_SparseMatrixELLT<fpreal32> UT_SparseMatrixELLF;
641 typedef UT_SparseMatrixELLT<fpreal64> UT_SparseMatrixELLD;
642 
643 
644 template <typename T> class UT_API UT_SparseMatrixBuilderT;
645 
646 /// Sparse matrix class that uses the compressed sparse row storage scheme.
647 /// This is a general-purpose sparse matrix class that can contain an arbitrary
648 /// number of non-zero elements per row.
649 template <typename T>
651 {
652 public:
653  /// Triplet to construct matrix. Is used only for setValues and addValues
654  /// to construct the matrix.
655  struct Triplet {
656  exint myRow;
657  exint myCol;
659 
660  /// Construct a zero triplet.
661  Triplet() : myRow(0), myCol(0), myValue(0) {}
662 
663  /// Construct a triplet for the given row and column with value v.
664  Triplet(exint r, exint c, T v) : myRow(r), myCol(c), myValue(v) {}
665 
666  /// Returns true if the triplet's coordinate comes before the given
667  /// triplet's coordinate.
668  inline bool operator<(const Triplet &t) const {
669  return (myRow < t.myRow || (myRow == t.myRow && myCol < t.myCol));
670  }
671  };
672 
673  /// Defines one cell in the matrix. Can also be used to directly construct
674  /// the matrix with setValuesFromRows, which is easier to multithread.
675  struct ColumnValue {
676  exint myCol;
678 
679  /// Construct a zero column-value.
680  ColumnValue() : myCol(0), myValue(0) {}
681 
682  /// Construct a column value for the given column with value v.
683  ColumnValue(exint c, T v) : myCol(c), myValue(v) {}
684 
685  /// Returns true if the column-value's column comes before the given
686  /// column-value's column.
687  inline bool operator<(const ColumnValue &cv) const {
688  return myCol < cv.myCol;
689  }
690  };
691 
692  /// Constructs an uninitialized sparse matrix
693  UT_SparseMatrixCSRT();
694 
695  /// Constructs a sparse matrix initialized with the given rows and columns.
696  UT_SparseMatrixCSRT(exint rows, exint columns);
697 
698  /// Constructs a sparse matrix that is a deep copy of the given sparse
699  /// matrix.
700  explicit UT_SparseMatrixCSRT(const UT_SparseMatrixCSRT &matrix);
701 
702  /// Sets the size of the matrix with the given rows and columns and clears
703  /// the matrix.
704  void init(exint rows, exint columns);
705 
706  /// Clears all data in the matrix
707  void zero();
708 
709  /// Swaps size and data elements with the given sparse matrix.
710  void swap(UT_SparseMatrixCSRT &other);
711 
712  /// The number of rows in the matrix.
713  exint getNumRows() const { return myNumRows; }
714 
715  /// The number of columns in the matrix.
716  exint getNumCols() const { return myNumCols; }
717 
718  /// The total number of non-zero elements in the matrix.
719  exint getNumNonZeros() const { return myCells.size(); }
720 
721  /// The number of non-zero elements in the given row.
722  exint getNumCellsInRow(exint row) const { return rowEnd(row) - rowBegin(row); }
723 
724  /// Clears then fills the matrix with the given array of triplets.
725  /// This sorts the triplets array.
726  /// Triplets with the same indices will be summed together.
727  /// A cell whose value is smaller than the tolerance will be deleted.
728  void setValues(UT_Array<Triplet> &triplets, T tolerance = 0.0);
729 
730  /// Clears then fills the matrix with the given array of ColumnValue arrays.
731  /// rowValues must have the same size as the number of rows in this matrix.
732  /// This data layout makes it easier to multithread the construction of the
733  /// sparse matrix entries.
734  void setValuesFromRows(UT_Array<UT_Array<ColumnValue>> &rowValues,
735  T tolerance = 0.0);
736 
737  /// Adds the given cells to the matrix. Triplets with the same coordinates
738  /// as cells already present in the matrix will be added together.
739  /// The array of triplets will be sorted unless alreadySorted is true.
740  /// Cells whose magnitude is below the tolerance value after the addition
741  /// will be discarded
742  void addValues(UT_Array<Triplet> &triplets, bool alreadySorted = false,
743  T tolerance=1e-5);
744 
745  /// Writes the diagonal values of the matrix into the given vector
746  /// If invert is true, the values will be inverted (i.e. 1/x).
747  /// If a row diagonal's value <= tolerance, then defaultValue is placed
748  /// in its result row.
749  void extractDiagonal(UT_VectorT<T> &result, bool invert = false,
750  T tolerance = 0, T defaultValue = 0) const;
751 
752  /// Matrix-vector multiplication. Computes result = M * v.
753  void multVec(const UT_VectorT<T> &v, UT_VectorT<T> &result) const;
754 
755  /// Matrix-vector multiplication.Computes result = M * v and returns the dot
756  /// product of v and the result.
757  UT_PreciseT<T> multVecAndDot(const UT_VectorT<T> &v, UT_VectorT<T> &result) const;
758 
759  /// Matrix-vector multiplication. Computes result = y - Ax. result can refer
760  /// to the same vector as y.
761  void subtractMultVec(const UT_VectorT<T> &x, const UT_VectorT<T> &y,
762  UT_VectorT<T> &result) const;
763 
764  /// Matrix-matrix multiplication. Computes result = M * B. Resultant
765  /// elements whose values are smaller than the given tolerance are deleted.
766  void multMatrix(const UT_SparseMatrixCSRT<T> &B, UT_SparseMatrixCSRT &result,
767  T tolerance = 0.0) const;
768 
769  /// Computes the transpose of this matrix. result can refer to this matrix.
770  void transpose(UT_SparseMatrixCSRT<T> &result) const;
771 
772  /// Element-wise scalar multiplication. Every element in this matrix is
773  /// multiplied by the given scalar.
774  void scale(const T& scalar);
775 
776  /// Element-wise scalar multiplication. Every element in row i of the matrix
777  /// is multiplied by scalars[i]
778  void scaleRows(const UT_VectorT<T> &scalars);
779 
780  /// Multiples each row such that the row diagonal becomes 1. If the row
781  /// diagonal is below the tolerance, the row is unchanged.
782  void normalizeRows(T tolerance = 1e-5);
783 
784  /// Calculates the norms of each row. result must have the same size as the
785  /// number of rows in this matrix. p determines the exponent of each column.
786  /// i.e. result[i] = SUM( |Aij|^p ) for j in row(i)
787  /// If invert is true, then inverts each sum
788  /// i.e. result[i] = 1.0 / SUM( |Aij|^p ) for j in row(i)
789  /// If the row norm is smaller than the tolerance, defaultValue is
790  /// placed in its result entry instead.
791  void rowNorms(UT_VectorT<T> &result, int p = 2, bool invert = false,
792  T tolerance = 0, T defaultValue = 0) const;
793 
794  /// Analyzes the matrix and constructs a directed acyclic graph as a level
795  /// set, where each level contains rows that are independent of each other.
796  /// More specifically, level i of the levelSet contains rows that are only
797  /// dependent on rows present in levels [0, i-1].
798  /// The size of the levelSet is the number of levels.
799  /// If buildForUpper is true, this constructs a level set for the upper
800  /// triangular part of the matrix. Otherwise, this constructs a level set
801  /// for the lower triangular part of the matrix.
803  void buildDependencyLevelSet(LevelSet &levelSet, bool buildForUpper = false) const;
804 
805  /// Assumes this is a lower triangular matrix. Solves the equation Ax = b.
806  /// Values may be present in the upper triangular part of the matrix
807  /// but they will be ignored.
808  /// If unitDiagonal is true, then 1 will be used instead of the
809  /// current diagonal values.
810  /// If a diagonal of A is zero within the tolerance, 0 will be placed
811  /// in the corresponding row for x.
812  /// If negate is true, then x is negated
813  /// x can refer to the same vector as b to solve in-place.
814  void solveLowerTriangular(UT_VectorT<T> &x, const UT_VectorT<T> &b,
815  bool unitDiagonal = false, T tolerance = 1e-5,
816  bool negate = false) const;
817 
818  /// Assumes this is an upper triangular matrix. Solves the equation Ax = b.
819  /// Values may be present in the lower triangular part of the matrix
820  /// but they will be ignored.
821  /// If unitDiagonal is true, then 1 will be used instead of the
822  /// current diagonal values.
823  /// If a diagonal of A is zero within the tolerance, 0 will be placed
824  /// in the corresponding row for x.
825  /// If negate is true, then x is negated
826  /// x can refer to the same vector as b to solve in-place.
827  void solveUpperTriangular(UT_VectorT<T> &x, const UT_VectorT<T> &b,
828  bool unitDiagonal = false, T tolerance = 1e-5,
829  bool negate = false) const;
830 
831  /// Solves the lower triangular system A^T x = b.
832  /// A is assumed to be an upper triangular matrix.
833  /// Values may be present in the lower triangular part of the matrix
834  /// but they will be ignored.
835  /// If unitDiagonal is true, then 1 will be used instead of the
836  /// current diagonal values.
837  /// If a diagonal of A is zero within the tolerance, 0 will be placed
838  /// in the corresponding row for x.
839  /// If negate is true, then x is negated.
840  void solveLowerTriangularTranspose(UT_VectorT<T> &x, const UT_VectorT<T> &b,
841  bool unitDiagonal = false, T tolerance = 1e-5,
842  bool negate = false) const;
843 
844  /// Solves the upper triangular system A^T x = b.
845  /// A is assumed to be a lower triangular matrix.
846  /// Values may be present in the upper triangular part of the matrix
847  /// but they will be ignored.
848  /// If unitDiagonal is true, then 1 will be used instead of the
849  /// current diagonal values.
850  /// If a diagonal of A is zero within the tolerance, 0 will be placed
851  /// in the corresponding row for x.
852  /// If negate is true, then x is negated.
853  void solveUpperTriangularTranspose(UT_VectorT<T> &x, const UT_VectorT<T> &b,
854  bool unitDiagonal = false, T tolerance = 1e-5,
855  bool negate = false) const;
856 
857  /// Factorizes this matrix into a triangular matrix M such that M*M^T ~= A
858  /// This keeps the non-zero structure of the matrix. The matrix must be
859  /// symmetric.
860  /// This is unstable depending on the diagonal; Returns -1 if the diagonal
861  /// contains zeros within the given tolerance, and -2 if not positive
862  /// definite (diagonal has negative numbers).
863  /// Returns 0 on success.
864  /// If useModified = true, then this uses the Modified Incomplete Cholesky
865  /// Factorization, where discarded entries are accounted for by adjusting
866  /// the diagonal. tau is a tuning constant that scales the adjustment from
867  /// the missing entries. mindiagratio is the relative tolerance at which
868  /// the original diagonal will be used if the matrix becomes negative
869  /// definite as a result of the computation.
870  int incompleteCholeskyFactorization(T tol=1e-5, bool useModified = false,
871  T tau = 0.97, T mindiagratio = 0.25);
872 
873  /// Factorizes this matrix into L and U where LU approximates this matrix.
874  /// The non-zero structure of LU matches this matrix; that is, this produces
875  /// ILU(0) of this matrix. A lower triangular level set is needed for
876  /// multithreading.
877  /// e.g.
878  /// LevelSet levelSet;
879  /// A.buildDependencyLevelSet(levelSet);
880  /// A.incompleteLUFactorization(L, U, 1e-5, levelSet);
881  void incompleteLUFactorization(UT_SparseMatrixCSRT<T> &L,
882  UT_SparseMatrixCSRT<T> &U, T tolerance,
883  const LevelSet &workUnits, UT_Interrupt *boss = NULL) const;
884 
885  /// Initializes out to a matrix with size (rowend-rowstart, colend-colstart)
886  /// and fills it with the elements defined by the submatrix in rows
887  /// [rowstart, rowend), [colstart, colend).
888  void extractSubMatrix(UT_SparseMatrixCSRT<T> &out,
889  exint rowstart, exint rowend, exint colstart, exint colend) const;
890 
891  /// Calls inspect(row, col, value) on each cell. Each (row, col) is visited
892  /// exactly once. This is done in parallel over the rows by default, but can
893  /// be done serially by setting serial = true. inspect can modify the cell
894  /// values if it has the signature (exint, exint, &T).
895  template <typename CellFunction>
896  void visit(const CellFunction& inspect, bool serial = false)
897  {
898  parallelExecute(UT_BlockedRange<exint>(0, myNumRows), [this, &inspect]
899  (const UT_BlockedRange<exint> &r)
900  {
901  for (exint row = r.begin(); row < r.end(); ++row)
902  for (exint i = rowBegin(row), n = rowEnd(row); i < n; ++i)
903  inspect(row, cellColumn(i), cellValue(i));
904  }, serial, myLightRowOpSubRatio, myLightRowOpGrainSize);
905  }
906 
907  /// Condition is called for each cell exactly once. If the condition returns
908  /// true, then the cell is removed from the matrix. The function signature
909  /// of the condition should be (exint row, exint col, T value) -> bool.
910  /// condition can modify the cell vakyes if it has the signature
911  /// (exint, exint, &T)
912  template <typename CellFunction>
913  void removeIf(const CellFunction& condition)
914  {
915  exint idx = -1, i = rowBegin(0);
916  for (exint row = 0; row < myNumRows; ++row)
917  {
918  for (exint n = rowEnd(row); i < n; ++i)
919  if (!condition(row, cellColumn(i), cellValue(i)) && ++idx != i)
920  myCells[idx] = myCells[i];
921  myRowOffsets[row+1] = idx+1;
922  }
923  myCells.setSize(idx+1);
924  }
925 
926  /// Returns true if any non-zero element in this matrix has value NaN.
927  bool hasNan() const;
928 
929  /// Asserts that the matrix contains no element with value NaN.
930  void testForNan() const;
931 
932  /// Writes the full matrix, including zeros, to the given output stream.
933  void printFull(std::ostream &os) const;
934 
935  /// Serializes the matrix to the given output stream.
936  void save(std::ostream &os) const;
937 
938  /// Deserializes the matrix from the given input stream.
939  void load(UT_IStream &is);
940 
941  /// Deep copies all data from the given matrix.
942  UT_SparseMatrixCSRT<T> &operator=(const UT_SparseMatrixCSRT<T> &m);
943 
944  /// Equivalent to A.scale(scalar)
945  UT_SparseMatrixCSRT<T> &operator*=(const T& scalar) { scale(scalar); return *this; }
946 
947 private:
948 
949  // Helper functions to get value and column
950  inline T& cellValue(exint offset) { return myCells[offset].myValue; }
951  inline const T &cellValue(exint offset) const { return myCells[offset].myValue; }
952  inline const exint &cellColumn(exint offset) const { return myCells[offset].myCol; }
953 
954  // Helper functions to return the beginning and end indexes into myValues
955  // for the given row. This can be used to iterate over the non-zero values
956  // of a particular row.
957  // e.g. for (exint i = rowBegin(row); i < rowEnd(row); ++i)
958  // ... Do something with cellValue(i) and cellColumn(i)
959  inline const exint &rowBegin(exint row) const { return myRowOffsets[row]; }
960  inline const exint &rowEnd(exint row) const { return myRowOffsets[row + 1]; }
961 
962  // Multithreading constants
963  // If the number of rows is less than this number, all operations run
964  // serially
965  static const int myParallelThreshold = 1000;
966 
967  // The subratio and grainsize used for row operations with cheap work units
968  static const int myLightRowOpSubRatio = 0;
969  static const int myLightRowOpGrainSize = 64;
970 
971  // Executes the given body over the given range. Can force this method to
972  // run serially by setting serial = true
973  template <typename Body>
974  inline void parallelExecute(const UT_BlockedRange<exint> &range,
975  const Body &body, bool serial = false, int subratio = 0,
976  int grainsize = 32, exint parallelThreshold = -1) const {
977  // Use internal parallel threshold by default
978  if (parallelThreshold < 0)
979  parallelThreshold = myParallelThreshold;
980  // Perform serially if specified or if number of rows is below threshold
981  if (myNumRows < parallelThreshold || serial)
982  UTserialFor(range, body);
983  else
984  UTparallelFor(range, body, subratio, grainsize);
985  };
986 
987 private:
988  exint myNumRows;
989  exint myNumCols;
990 
991  // Stores the non-zero cells of the matrix, sorted by row and column.
992  // Has size Nz, where Nz is the number of non-zero values.
993  UT_Array<ColumnValue> myCells;
994 
995  // Stores the offset of myValues at which each row starts.
996  // Has size myNumRows + 1.
997  // myRowOffsets[r] is the offset of myValues where row r begins.
998  // myRowOffsets[r + 1] - myRowOffsets[r] is the number of
999  // non-zero values in row r.
1000  // myRowOffsets[myNumRows] has the sentinel value Nz
1001  UT_ExintArray myRowOffsets;
1002 
1003  friend class UT_SparseMatrixBuilderT<T>;
1004 };
1005 
1008 
1009 /// Simple helper class to construct a SparseMatrixCSRT. This operates on the
1010 /// matrix data in-place and is generally much more performant than constructing
1011 /// via triplets, but the builder has the limitation of only being able to add
1012 /// to rows sequentially; i.e. each call to startRow must start the construction
1013 /// of a row greater than the last, while addToColumn can be called with any
1014 /// column ordering.
1015 /// 1. UT_SparseMatrixCSRD matrix(M, N);
1016 /// 2. UT_SparseMatrixBuilderT builder(matrix);
1017 /// 3. for (exint r = startRow; r < endRow; ++r)
1018 /// 4. builder.startRow(r);
1019 /// 5. for (any column)
1020 /// 6. builder.addToColumn(column, value);
1021 /// 7. builder.finish();
1022 template <typename T>
1023 class UT_API UT_SparseMatrixBuilderT
1024 {
1025 public:
1027 
1028  /// Construct an uninitialized sparse matrix builder.
1029  UT_SparseMatrixBuilderT() : myMatrix(NULL) {}
1030 
1031  /// Construct a sparse matrix builder for the given matrix.
1032  UT_SparseMatrixBuilderT(SparseMatrixType &matrix) { init(matrix); }
1033 
1034  /// Initializes the builder to start building the given matrix
1035  void init(SparseMatrixType &matrix);
1036 
1037  /// Sets the minimum capacity of non-zero elements for the given matrix.
1038  /// This is used to preallocate memory for the matrix to avoid excessive
1039  /// reallocation as elements are added.
1040  void setCapacity(exint capacity);
1041 
1042  /// Advances the current row to the given row. The row must be greater than
1043  /// the current row.
1044  void startRow(exint row);
1045 
1046  /// Adds the given value at the given column to the current row. If a cell
1047  /// is already present for the given column, the given value is added to it.
1048  void addToColumn(exint col, T value);
1049 
1050  /// Completes the remaining rows and finalizes the matrix.
1051  void finish();
1052 
1053 private:
1054  exint myCurrentRow;
1055  SparseMatrixType *myMatrix;
1056 };
1057 
1060 
1061 #endif
int getNumRows() const
int getNumRows() const
#define THREADED_METHOD4_CONST(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2, PARMTYPE3, PARMNAME3, PARMTYPE4, PARMNAME4)
Color4< T > operator*(S a, const Color4< T > &v)
Definition: ImathColor.h:728
exint getNumCols() const
The number of columns in the matrix.
ColumnValue()
Construct a zero column-value.
GLsizei GLenum const void * indices
Definition: glcorearb.h:405
GLenum GLint * range
Definition: glcorearb.h:1924
bool operator<(const ColumnValue &cv) const
int getNumRows() const
size_t index(int row, int nz) const
UT_SparseMatrixRowT< fpreal64 > UT_SparseMatrixRowD
GLboolean invert
Definition: glcorearb.h:548
exint getNumNonZeros() const
The total number of non-zero elements in the matrix.
void swap(UT::ArraySet< Key, MULTI, MAX_LOAD_FACTOR_256, Clearer, Hash, KeyEqual > &a, UT::ArraySet< Key, MULTI, MAX_LOAD_FACTOR_256, Clearer, Hash, KeyEqual > &b)
Definition: UT_ArraySet.h:1561
const GLdouble * v
Definition: glcorearb.h:836
bool shouldMultiThread() const
typedef void(APIENTRYP PFNGLCULLFACEPROC)(GLenum mode)
UT_SparseMatrixCSRT< T > SparseMatrixType
Plane3< T > operator-(const Plane3< T > &plane)
Definition: ImathPlane.h:249
const UT_ValArray< T > & getRowValues() const
ColumnValue(exint c, T v)
Construct a column value for the given column with value v.
class UT_API UT_SparseMatrixRowT
#define UT_API
Definition: UT_API.h:12
GLint y
Definition: glcorearb.h:102
bool isStillSorted() const
png_uint_32 i
Definition: png.h:2877
void subtractMultVecNoThread(const UT_VectorT< T > &v, UT_VectorT< T > &result) const
int getNumCols() const
UT_SparseMatrixRowT< fpreal32 > UT_SparseMatrixRowF
bool operator<(const Triplet &t) const
UT_SparseMatrixCSRT< T > & operator*=(const T &scalar)
Equivalent to A.scale(scalar)
GLsizei GLboolean transpose
Definition: glcorearb.h:831
void multVec(const UT_VectorT< T > &v, UT_VectorT< T > &result) const
bool operator==(const BaseDimensions< T > &a, const BaseDimensions< Y > &b)
Definition: Dimensions.h:137
long long int64
Definition: SYS_Types.h:100
GA_API const UT_StringHolder scale
GLdouble n
Definition: glcorearb.h:2007
GLfloat f
Definition: glcorearb.h:1925
exint getNumRows() const
The number of rows in the matrix.
#define UT_ASSERT(ZZ)
Definition: UT_Assert.h:102
void removeIf(const CellFunction &condition)
int64 exint
Definition: SYS_Types.h:109
UT_SparseMatrixBuilderT(SparseMatrixType &matrix)
Construct a sparse matrix builder for the given matrix.
double fpreal64
Definition: SYS_Types.h:185
int getNonZerosPerRow() const
GLintptr offset
Definition: glcorearb.h:664
typename UT_TypePromoteT< T >::PreciseType UT_PreciseT
bool operator>=(const UT_Vector2T< T > &v1, const UT_Vector2T< T > &v2)
Definition: UT_Vector2.h:386
bool operator<(const GU_TetrahedronFacet &a, const GU_TetrahedronFacet &b)
#define SYS_ALIGN16
Definition: SYS_Align.h:117
UT_SparseMatrixBuilderT< fpreal32 > UT_SparseMatrixBuilderF
void visit(const CellFunction &inspect, bool serial=false)
png_bytepp row
Definition: png.h:1836
GLboolean GLboolean GLboolean b
Definition: glcorearb.h:1221
Triplet()
Construct a zero triplet.
void subtractMultVec(const UT_VectorT< T > &v, UT_VectorT< T > &result) const
Do result = result - (M * v)
bool operator>(const UT_Vector2T< T > &v1, const UT_Vector2T< T > &v2)
Definition: UT_Vector2.h:381
#define THREADED_METHOD2_CONST(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2)
#define THREADED_METHOD1_CONST(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1)
const UT_IntArray & getColumns() const
GLsizei const GLfloat * value
Definition: glcorearb.h:823
Triplet(exint r, exint c, T v)
Construct a triplet for the given row and column with value v.
exint getNumCellsInRow(exint row) const
The number of non-zero elements in the given row.
void accept(Visitor &visitor) const
UT_SparseMatrixBuilderT< fpreal64 > UT_SparseMatrixBuilderD
int findCellFromRow(int row) const
GLuint index
Definition: glcorearb.h:785
UT_IntArray & getColumns()
GLint GLenum GLint x
Definition: glcorearb.h:408
GLuint GLfloat * val
Definition: glcorearb.h:1607
UT_ValArray< T > & getRowValues()
GLboolean r
Definition: glcorearb.h:1221
bool operator!=(const BaseDimensions< T > &a, const BaseDimensions< Y > &b)
Definition: Dimensions.h:165
#define const
Definition: zconf.h:214
UT_SparseMatrixCSRT< fpreal32 > UT_SparseMatrixCSRF
#define THREADED_METHOD3_CONST(CLASSNAME, DOMULTI, METHOD, PARMTYPE1, PARMNAME1, PARMTYPE2, PARMNAME2, PARMTYPE3, PARMNAME3)
UT_SparseMatrixBuilderT()
Construct an uninitialized sparse matrix builder.
Declare prior to use.
float fpreal32
Definition: SYS_Types.h:184
void UTparallelFor(const Range &range, const Body &body, const int subscribe_ratio=2, const int min_grain_size=1)
#define THREADED_METHOD(CLASSNAME, DOMULTI, METHOD)
bool appendRowElement(int row, int col, T val, int &rowidx)
Quat< T > operator+(const Quat< T > &q1, const Quat< T > &q2)
Definition: ImathQuat.h:919
void UTserialFor(const Range &range, const Body &body)
bool shouldMultiThread() const
int getNumCols() const
UT_Array< UT_ExintArray > LevelSet
GLenum src
Definition: glcorearb.h:1792
UT_SparseMatrixCSRT< fpreal64 > UT_SparseMatrixCSRD