MFEM  v3.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
sparsemat.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_SPARSEMAT
13 #define MFEM_SPARSEMAT
14 
15 // Data types for sparse matrix
16 
17 #include "../general/mem_alloc.hpp"
18 #include "../general/table.hpp"
19 #include "densemat.hpp"
20 #include <iostream>
21 
22 namespace mfem
23 {
24 
25 class
26 #if defined(__alignas_is_defined)
27  alignas(double)
28 #endif
29  RowNode
30 {
31 public:
32  double Value;
33  RowNode *Prev;
34  int Column;
35 };
36 
39 {
40 private:
41 
45  int *I, *J;
46 
48  double *A;
49 
51  RowNode **Rows;
52 
53  mutable int current_row;
54  mutable int* ColPtrJ;
55  mutable RowNode ** ColPtrNode;
56 
57 #ifdef MFEM_USE_MEMALLOC
59  RowNodeAlloc * NodesMem;
60 #endif
61 
63  bool ownGraph;
65  bool ownData;
66 
68  bool isSorted;
69 
70  void Destroy(); // Delete all owned data
71  void SetEmpty(); // Init all entries with empty values
72 
73 public:
75  SparseMatrix() { SetEmpty(); }
76 
81  explicit SparseMatrix(int nrows, int ncols = 0);
82 
85  SparseMatrix(int *i, int *j, double *data, int m, int n);
86 
89  SparseMatrix(int *i, int *j, double *data, int m, int n, bool ownij,
90  bool owna, bool issorted);
91 
95  SparseMatrix(const SparseMatrix &mat, bool copy_graph = true);
96 
100  void MakeRef(const SparseMatrix &master);
101 
103  int Size() const { return Height(); }
104 
106  void Clear() { Destroy(); SetEmpty(); }
107 
109  bool Empty() { return (A == NULL) && (Rows == NULL); }
110 
112  inline int *GetI() const { return I; }
114  inline int *GetJ() const { return J; }
116  inline double *GetData() const { return A; }
118  int RowSize(const int i) const;
120  int MaxRowSize() const;
122  int *GetRowColumns(const int row);
123  const int *GetRowColumns(const int row) const;
125  double *GetRowEntries(const int row);
126  const double *GetRowEntries(const int row) const;
127 
129 
136  void SetWidth(int width_ = -1);
137 
139 
140  int ActualWidth();
141 
143  void SortColumnIndices();
144 
146  virtual double &Elem(int i, int j);
147 
149  virtual const double &Elem(int i, int j) const;
150 
152  double &operator()(int i, int j);
153 
155  const double &operator()(int i, int j) const;
156 
158  void GetDiag(Vector & d) const;
159 
161  virtual void Mult(const Vector &x, Vector &y) const;
162 
164  void AddMult(const Vector &x, Vector &y, const double a = 1.0) const;
165 
167  void MultTranspose(const Vector &x, Vector &y) const;
168 
170  void AddMultTranspose(const Vector &x, Vector &y,
171  const double a = 1.0) const;
172 
173  void PartMult(const Array<int> &rows, const Vector &x, Vector &y) const;
174  void PartAddMult(const Array<int> &rows, const Vector &x, Vector &y,
175  const double a=1.0) const;
176 
178  void BooleanMult(const Array<int> &x, Array<int> &y) const;
180  void BooleanMultTranspose(const Array<int> &x, Array<int> &y) const;
181 
183  double InnerProduct(const Vector &x, const Vector &y) const;
184 
186  void GetRowSums(Vector &x) const;
188  double GetRowNorml1(int irow) const;
189 
191  virtual MatrixInverse *Inverse() const;
192 
194  void EliminateRow(int row, const double sol, Vector &rhs);
196 
202  void EliminateRow(int row, int setOneDiagonal = 0);
203  void EliminateCol(int col);
205  void EliminateCols(Array<int> &cols, Vector *x = NULL, Vector *b = NULL);
206 
211  void EliminateRowCol(int rc, const double sol, Vector &rhs, int d = 0);
212 
215  void EliminateRowColMultipleRHS(int rc, const Vector &sol,
216  DenseMatrix &rhs, int d = 0);
217 
218  void EliminateRowCol(int rc, int d = 0);
220  void EliminateRowColDiag(int rc, double value);
221  // Same as above + save the eliminated entries in Ae so that
222  // (*this) + Ae is the original matrix
223  void EliminateRowCol(int rc, SparseMatrix &Ae, int d = 0);
224 
226  void SetDiagIdentity();
228  void EliminateZeroRows();
229 
231  void Gauss_Seidel_forw(const Vector &x, Vector &y) const;
232  void Gauss_Seidel_back(const Vector &x, Vector &y) const;
233 
235  double GetJacobiScaling() const;
238  void Jacobi(const Vector &b, const Vector &x0, Vector &x1, double sc) const;
239 
240  void DiagScale(const Vector &b, Vector &x, double sc = 1.0) const;
241 
243  void Jacobi2(const Vector &b, const Vector &x0, Vector &x1,
244  double sc = 1.0) const;
245 
247  void Jacobi3(const Vector &b, const Vector &x0, Vector &x1,
248  double sc = 1.0) const;
249 
254  virtual void Finalize(int skip_zeros = 1);
255 
256  bool Finalized() const { return (A != NULL); }
257  bool areColumnsSorted() const { return isSorted; }
258 
262  void GetBlocks(Array2D<SparseMatrix *> &blocks) const;
263 
264  void GetSubMatrix(const Array<int> &rows, const Array<int> &cols,
265  DenseMatrix &subm);
266 
267  inline void SetColPtr(const int row) const;
268  inline void ClearColPtr() const;
269  inline double &SearchRow(const int col);
270  inline void _Add_(const int col, const double a)
271  { SearchRow(col) += a; }
272  inline void _Set_(const int col, const double a)
273  { SearchRow(col) = a; }
274  inline double _Get_(const int col) const;
275 
276  inline double &SearchRow(const int row, const int col);
277  inline void _Add_(const int row, const int col, const double a)
278  { SearchRow(row, col) += a; }
279  inline void _Set_(const int row, const int col, const double a)
280  { SearchRow(row, col) = a; }
281 
282  void Set(const int i, const int j, const double a);
283  void Add(const int i, const int j, const double a);
284 
285  void SetSubMatrix(const Array<int> &rows, const Array<int> &cols,
286  const DenseMatrix &subm, int skip_zeros = 1);
287 
288  void SetSubMatrixTranspose(const Array<int> &rows, const Array<int> &cols,
289  const DenseMatrix &subm, int skip_zeros = 1);
290 
291  void AddSubMatrix(const Array<int> &rows, const Array<int> &cols,
292  const DenseMatrix &subm, int skip_zeros = 1);
293 
294  bool RowIsEmpty(const int row) const;
295 
304  virtual int GetRow(const int row, Array<int> &cols, Vector &srow) const;
305 
306  void SetRow(const int row, const Array<int> &cols, const Vector &srow);
307  void AddRow(const int row, const Array<int> &cols, const Vector &srow);
308 
309  void ScaleRow(const int row, const double scale);
310  // this = diag(sl) * this;
311  void ScaleRows(const Vector & sl);
312  // this = this * diag(sr);
313  void ScaleColumns(const Vector & sr);
314 
318 
321  void Add(const double a, const SparseMatrix &B);
322 
323  SparseMatrix &operator=(double a);
324 
325  SparseMatrix &operator*=(double a);
326 
328  void Print(std::ostream &out = std::cout, int width_ = 4) const;
329 
331  void PrintMatlab(std::ostream &out = std::cout) const;
332 
334  void PrintMM(std::ostream &out = std::cout) const;
335 
337  void PrintCSR(std::ostream &out) const;
338 
340  void PrintCSR2(std::ostream &out) const;
341 
343  int Walk(int &i, int &j, double &a);
344 
346  double IsSymmetric() const;
347 
349  void Symmetrize();
350 
352  virtual int NumNonZeroElems() const;
353 
354  double MaxNorm() const;
355 
357  int CountSmallElems(double tol) const;
358 
360  void SetGraphOwner(bool ownij) { ownGraph = ownij; }
362  void SetDataOwner(bool owna) { ownData = owna; }
364  bool OwnsGraph() const { return ownGraph; }
366  bool OwnsData() const { return ownData; }
368  void LoseData() { ownGraph = ownData = false; }
369 
370  void Swap(SparseMatrix &other);
371 
373  virtual ~SparseMatrix() { Destroy(); }
374 };
375 
377 void SparseMatrixFunction(SparseMatrix &S, double (*f)(double));
378 
379 
381 SparseMatrix *Transpose(const SparseMatrix &A);
383 SparseMatrix *TransposeAbstractSparseMatrix (const AbstractSparseMatrix &A,
384  int useActualWidth);
385 
392 SparseMatrix *Mult(const SparseMatrix &A, const SparseMatrix &B,
393  SparseMatrix *OAB = NULL);
394 
396 SparseMatrix *MultAbstractSparseMatrix (const AbstractSparseMatrix &A,
397  const AbstractSparseMatrix &B);
398 
399 
402 SparseMatrix *RAP(const SparseMatrix &A, const SparseMatrix &R,
403  SparseMatrix *ORAP = NULL);
404 
406 SparseMatrix *RAP(const SparseMatrix &Rt, const SparseMatrix &A,
407  const SparseMatrix &P);
408 
410 SparseMatrix *Mult_AtDA(const SparseMatrix &A, const Vector &D,
411  SparseMatrix *OAtDA = NULL);
412 
413 
415 SparseMatrix * Add(const SparseMatrix & A, const SparseMatrix & B);
417 SparseMatrix * Add(double a, const SparseMatrix & A, double b,
418  const SparseMatrix & B);
420 SparseMatrix * Add(Array<SparseMatrix *> & Ai);
421 
422 
423 // Inline methods
424 
425 inline void SparseMatrix::SetColPtr(const int row) const
426 {
427  if (Rows)
428  {
429  if (ColPtrNode == NULL)
430  {
431  ColPtrNode = new RowNode *[width];
432  for (int i = 0; i < width; i++)
433  {
434  ColPtrNode[i] = NULL;
435  }
436  }
437  for (RowNode *node_p = Rows[row]; node_p != NULL; node_p = node_p->Prev)
438  {
439  ColPtrNode[node_p->Column] = node_p;
440  }
441  }
442  else
443  {
444  if (ColPtrJ == NULL)
445  {
446  ColPtrJ = new int[width];
447  for (int i = 0; i < width; i++)
448  {
449  ColPtrJ[i] = -1;
450  }
451  }
452  for (int j = I[row], end = I[row+1]; j < end; j++)
453  {
454  ColPtrJ[J[j]] = j;
455  }
456  }
457  current_row = row;
458 }
459 
460 inline void SparseMatrix::ClearColPtr() const
461 {
462  if (Rows)
463  for (RowNode *node_p = Rows[current_row]; node_p != NULL;
464  node_p = node_p->Prev)
465  {
466  ColPtrNode[node_p->Column] = NULL;
467  }
468  else
469  for (int j = I[current_row], end = I[current_row+1]; j < end; j++)
470  {
471  ColPtrJ[J[j]] = -1;
472  }
473 }
474 
475 inline double &SparseMatrix::SearchRow(const int col)
476 {
477  if (Rows)
478  {
479  RowNode *node_p = ColPtrNode[col];
480  if (node_p == NULL)
481  {
482 #ifdef MFEM_USE_MEMALLOC
483  node_p = NodesMem->Alloc();
484 #else
485  node_p = new RowNode;
486 #endif
487  node_p->Prev = Rows[current_row];
488  node_p->Column = col;
489  node_p->Value = 0.0;
490  Rows[current_row] = ColPtrNode[col] = node_p;
491  }
492  return node_p->Value;
493  }
494  else
495  {
496  const int j = ColPtrJ[col];
497  MFEM_VERIFY(j != -1, "Entry for column " << col << " is not allocated.");
498  return A[j];
499  }
500 }
501 
502 inline double SparseMatrix::_Get_(const int col) const
503 {
504  if (Rows)
505  {
506  RowNode *node_p = ColPtrNode[col];
507  return (node_p == NULL) ? 0.0 : node_p->Value;
508  }
509  else
510  {
511  const int j = ColPtrJ[col];
512  return (j == -1) ? 0.0 : A[j];
513  }
514 }
515 
516 inline double &SparseMatrix::SearchRow(const int row, const int col)
517 {
518  if (Rows)
519  {
520  RowNode *node_p;
521 
522  for (node_p = Rows[row]; 1; node_p = node_p->Prev)
523  {
524  if (node_p == NULL)
525  {
526 #ifdef MFEM_USE_MEMALLOC
527  node_p = NodesMem->Alloc();
528 #else
529  node_p = new RowNode;
530 #endif
531  node_p->Prev = Rows[row];
532  node_p->Column = col;
533  node_p->Value = 0.0;
534  Rows[row] = node_p;
535  break;
536  }
537  else if (node_p->Column == col)
538  {
539  break;
540  }
541  }
542  return node_p->Value;
543  }
544  else
545  {
546  int *Ip = I+row, *Jp = J;
547  for (int k = Ip[0], end = Ip[1]; k < end; k++)
548  {
549  if (Jp[k] == col)
550  {
551  return A[k];
552  }
553  }
554  MFEM_ABORT("Could not find entry for row = " << row << ", col = " << col);
555  }
556  return A[0];
557 }
558 
560 template<> inline void Swap<SparseMatrix>(SparseMatrix &a, SparseMatrix &b)
561 {
562  a.Swap(b);
563 }
564 
565 }
566 
567 #endif
Elem * Alloc()
Definition: mem_alloc.hpp:124
void _Add_(const int row, const int col, const double a)
Definition: sparsemat.hpp:277
double GetJacobiScaling() const
Determine appropriate scaling for Jacobi iteration.
Definition: sparsemat.cpp:1589
SparseMatrix * TransposeAbstractSparseMatrix(const AbstractSparseMatrix &A, int useActualWidth)
Transpose of a sparse matrix. A does not need to be a CSR matrix.
Definition: sparsemat.cpp:2472
int RowSize(const int i) const
Returns the number of elements in row i.
Definition: sparsemat.cpp:190
virtual int NumNonZeroElems() const
Returns the number of the nonzero elements in the matrix.
Definition: sparsemat.cpp:879
void Add(const int i, const int j, const double a)
Definition: sparsemat.cpp:1790
void _Add_(const int col, const double a)
Definition: sparsemat.hpp:270
void Clear()
Clear the contents of the SparseMatrix.
Definition: sparsemat.hpp:106
void _Set_(const int row, const int col, const double a)
Definition: sparsemat.hpp:279
void Print(std::ostream &out=std::cout, int width_=4) const
Prints matrix to stream out.
Definition: sparsemat.cpp:2201
void DiagScale(const Vector &b, Vector &x, double sc=1.0) const
Definition: sparsemat.cpp:1653
void MakeRef(const SparseMatrix &master)
Definition: sparsemat.cpp:164
double & SearchRow(const int col)
Definition: sparsemat.hpp:475
void PrintMM(std::ostream &out=std::cout) const
Prints matrix in Matrix Market sparse format.
Definition: sparsemat.cpp:2263
int * GetRowColumns(const int row)
Return a pointer to the column indices in a row.
Definition: sparsemat.cpp:237
void BooleanMult(const Array< int > &x, Array< int > &y) const
y = A * x, but treat all elements as booleans (zero=false, nonzero=true).
Definition: sparsemat.cpp:573
void Mult(const Table &A, const Table &B, Table &C)
C = A * B (as boolean matrices)
Definition: table.cpp:445
Abstract data type for sparse matrices.
Definition: matrix.hpp:69
void EliminateRowCol(int rc, const double sol, Vector &rhs, int d=0)
Definition: sparsemat.cpp:1045
Data type dense matrix using column-major storage.
Definition: densemat.hpp:22
void AddMult(const Vector &x, Vector &y, const double a=1.0) const
y += A * x (default) or y += a * A * x
Definition: sparsemat.cpp:423
double MaxNorm() const
Definition: sparsemat.cpp:899
virtual ~SparseMatrix()
Destroys sparse matrix.
Definition: sparsemat.hpp:373
void AddMultTranspose(const Vector &x, Vector &y, const double a=1.0) const
y += At * x (default) or y += a * At * x
Definition: sparsemat.cpp:498
Abstract data type for matrix inverse.
Definition: matrix.hpp:58
void EliminateRowColMultipleRHS(int rc, const Vector &sol, DenseMatrix &rhs, int d=0)
Definition: sparsemat.cpp:1125
virtual double & Elem(int i, int j)
Returns reference to a_{ij}.
Definition: sparsemat.cpp:332
void LoseData()
Lose the ownership of the graph (I, J) and data (A) arrays.
Definition: sparsemat.hpp:368
void AddRow(const int row, const Array< int > &cols, const Vector &srow)
Definition: sparsemat.cpp:1993
double * GetRowEntries(const int row)
Return a pointer to the entries in a row.
Definition: sparsemat.cpp:251
void SortColumnIndices()
Sort the column indices corresponding to each row.
Definition: sparsemat.cpp:303
void _Set_(const int col, const double a)
Definition: sparsemat.hpp:272
void BooleanMultTranspose(const Array< int > &x, Array< int > &y) const
y = At * x, but treat all elements as booleans (zero=false, nonzero=true).
Definition: sparsemat.cpp:596
SparseMatrix * Mult_AtDA(const SparseMatrix &A, const Vector &D, SparseMatrix *OAtDA)
Matrix multiplication A^t D A. All matrices must be finalized.
Definition: sparsemat.cpp:2794
void PartAddMult(const Array< int > &rows, const Vector &x, Vector &y, const double a=1.0) const
Definition: sparsemat.cpp:555
bool RowIsEmpty(const int row) const
Definition: sparsemat.cpp:1905
void GetBlocks(Array2D< SparseMatrix * > &blocks) const
Definition: sparsemat.cpp:752
void ScaleRow(const int row, const double scale)
Definition: sparsemat.cpp:2025
int ActualWidth()
Returns the actual Width of the matrix.
Definition: sparsemat.cpp:2380
void Symmetrize()
(*this) = 1/2 ((*this) + (*this)^t)
Definition: sparsemat.cpp:864
void SetGraphOwner(bool ownij)
Set the graph ownership flag (I and J arrays).
Definition: sparsemat.hpp:360
void Add(const DenseMatrix &A, const DenseMatrix &B, double alpha, DenseMatrix &C)
C = A + alpha*B.
Definition: densemat.cpp:2769
bool OwnsData() const
Get the data ownership flag (A array).
Definition: sparsemat.hpp:366
void ClearColPtr() const
Definition: sparsemat.hpp:460
void EliminateCol(int col)
Definition: sparsemat.cpp:1000
void ScaleColumns(const Vector &sr)
Definition: sparsemat.cpp:2084
void PrintCSR2(std::ostream &out) const
Prints a sparse matrix to stream out in CSR format.
Definition: sparsemat.cpp:2307
virtual int GetRow(const int row, Array< int > &cols, Vector &srow) const
Definition: sparsemat.cpp:1926
SparseMatrix()
Create an empty SparseMatrix.
Definition: sparsemat.hpp:75
Data type sparse matrix.
Definition: sparsemat.hpp:38
int Height() const
Get the height (size of output) of the Operator. Synonym with NumRows.
Definition: operator.hpp:35
void EliminateCols(Array< int > &cols, Vector *x=NULL, Vector *b=NULL)
Eliminate all columns &#39;i&#39; for which cols[i] != 0.
Definition: sparsemat.cpp:1014
bool OwnsGraph() const
Get the graph ownership flag (I and J arrays).
Definition: sparsemat.hpp:364
void SetWidth(int width_=-1)
Change the width of a SparseMatrix.
Definition: sparsemat.cpp:265
bool areColumnsSorted() const
Definition: sparsemat.hpp:257
HypreParMatrix * RAP(HypreParMatrix *A, HypreParMatrix *P)
Returns the matrix P^t * A * P.
Definition: hypre.cpp:1318
double IsSymmetric() const
Returns max_{i,j} |(i,j)-(j,i)| for a finalized matrix.
Definition: sparsemat.cpp:842
void SetDiagIdentity()
If a row contains only one diag entry of zero, set it to 1.
Definition: sparsemat.cpp:1403
void Transpose(const Table &A, Table &At, int _ncols_A)
Transpose a Table.
Definition: table.cpp:385
void ScaleRows(const Vector &sl)
Definition: sparsemat.cpp:2053
void Swap< SparseMatrix >(SparseMatrix &a, SparseMatrix &b)
Specialization of the template function Swap&lt;&gt; for class SparseMatrix.
Definition: sparsemat.hpp:560
virtual void Finalize(int skip_zeros=1)
Definition: sparsemat.cpp:680
void EliminateZeroRows()
If a row contains only zeros, set its diagonal to 1.
Definition: sparsemat.cpp:1412
double * GetData() const
Return element data.
Definition: sparsemat.hpp:116
void AddSubMatrix(const Array< int > &rows, const Array< int > &cols, const DenseMatrix &subm, int skip_zeros=1)
Definition: sparsemat.cpp:1733
int * GetI() const
Return the array I.
Definition: sparsemat.hpp:112
void Jacobi(const Vector &b, const Vector &x0, Vector &x1, double sc) const
Definition: sparsemat.cpp:1622
void SetSubMatrix(const Array< int > &rows, const Array< int > &cols, const DenseMatrix &subm, int skip_zeros=1)
Definition: sparsemat.cpp:1809
bool Finalized() const
Definition: sparsemat.hpp:256
void SparseMatrixFunction(SparseMatrix &S, double(*f)(double))
Applies f() to each element of the matrix (after it is finalized).
Definition: sparsemat.cpp:2407
virtual MatrixInverse * Inverse() const
Returns a pointer to approximation of the matrix inverse.
Definition: sparsemat.cpp:952
void EliminateRowColDiag(int rc, double value)
Perform elimination and set the diagonal entry to the given value.
Definition: sparsemat.cpp:1275
virtual void Mult(const Vector &x, Vector &y) const
Matrix vector multiplication.
Definition: sparsemat.cpp:417
void SetColPtr(const int row) const
Definition: sparsemat.hpp:425
SparseMatrix & operator=(double a)
Definition: sparsemat.cpp:2165
void SetSubMatrixTranspose(const Array< int > &rows, const Array< int > &cols, const DenseMatrix &subm, int skip_zeros=1)
Definition: sparsemat.cpp:1842
void Set(const int i, const int j, const double a)
Definition: sparsemat.cpp:1771
void PrintMatlab(std::ostream &out=std::cout) const
Prints matrix in matlab format.
Definition: sparsemat.cpp:2245
int Size() const
For backward compatibility define Size to be synonym of Height()
Definition: sparsemat.hpp:103
void Gauss_Seidel_back(const Vector &x, Vector &y) const
Definition: sparsemat.cpp:1514
void Gauss_Seidel_forw(const Vector &x, Vector &y) const
Gauss-Seidel forward and backward iterations over a vector x.
Definition: sparsemat.cpp:1439
void EliminateRow(int row, const double sol, Vector &rhs)
Eliminates a column from the transpose matrix.
Definition: sparsemat.cpp:957
void GetDiag(Vector &d) const
Returns the Diagonal of A.
Definition: sparsemat.cpp:389
double GetRowNorml1(int irow) const
For i = irow compute .
Definition: sparsemat.cpp:660
double & operator()(int i, int j)
Returns reference to A[i][j].
Definition: sparsemat.cpp:342
void MultTranspose(const Vector &x, Vector &y) const
Multiply a vector with the transposed matrix. y = At * x.
Definition: sparsemat.cpp:492
SparseMatrix & operator*=(double a)
Definition: sparsemat.cpp:2183
void SetDataOwner(bool owna)
Set the data ownership flag (A array).
Definition: sparsemat.hpp:362
bool Empty()
Check if the SparseMatrix is empty.
Definition: sparsemat.hpp:109
void SetRow(const int row, const Array< int > &cols, const Vector &srow)
Definition: sparsemat.cpp:1965
Vector data type.
Definition: vector.hpp:33
int Walk(int &i, int &j, double &a)
Walks the sparse matrix.
void GetSubMatrix(const Array< int > &rows, const Array< int > &cols, DenseMatrix &subm)
Definition: sparsemat.cpp:1877
void Jacobi2(const Vector &b, const Vector &x0, Vector &x1, double sc=1.0) const
Definition: sparsemat.cpp:1685
int MaxRowSize() const
Returns the maximum number of elements among all rows.
Definition: sparsemat.cpp:213
void Jacobi3(const Vector &b, const Vector &x0, Vector &x1, double sc=1.0) const
Definition: sparsemat.cpp:1709
SparseMatrix * MultAbstractSparseMatrix(const AbstractSparseMatrix &A, const AbstractSparseMatrix &B)
Matrix product of sparse matrices. A and B do not need to be CSR matrices.
Definition: sparsemat.cpp:2675
double _Get_(const int col) const
Definition: sparsemat.hpp:502
SparseMatrix & operator+=(SparseMatrix &B)
Definition: sparsemat.cpp:2112
void Swap(SparseMatrix &other)
Definition: sparsemat.cpp:2922
int * GetJ() const
Return the array J.
Definition: sparsemat.hpp:114
void PrintCSR(std::ostream &out) const
Prints matrix to stream out in hypre_CSRMatrix format.
Definition: sparsemat.cpp:2283
double InnerProduct(const Vector &x, const Vector &y) const
Compute y^t A x.
Definition: sparsemat.cpp:619
int CountSmallElems(double tol) const
Count the number of entries with |a_ij| &lt; tol.
Definition: sparsemat.cpp:922
void PartMult(const Array< int > &rows, const Vector &x, Vector &y) const
Definition: sparsemat.cpp:537
void GetRowSums(Vector &x) const
For all i compute .
Definition: sparsemat.cpp:641