PARALUTION  1.0.0
PARALUTION
local_matrix.hpp
Go to the documentation of this file.
1 #ifndef PARALUTION_LOCAL_MATRIX_HPP_
2 #define PARALUTION_LOCAL_MATRIX_HPP_
3 
4 #include "operator.hpp"
5 #include "base_matrix.hpp"
6 #include "backend_manager.hpp"
7 
8 namespace paralution {
9 
10 template <typename ValueType>
11 class LocalVector;
12 
13 // Local Matrix
14 template <typename ValueType>
15 class LocalMatrix : public Operator<ValueType> {
16 
17 public:
18 
19  LocalMatrix();
20  virtual ~LocalMatrix();
21 
22  virtual void info(void) const;
23 
25  unsigned int get_format(void) const;
26 
27  virtual int get_nrow(void) const;
28  virtual int get_ncol(void) const;
29  virtual int get_nnz(void) const;
30 
32  // and false if there is something wrong with the strcture or some of values are NaN
33  bool Check(void) const;
34 
36  void AllocateCSR(const std::string name, const int nnz, const int nrow, const int ncol);
38  void AllocateBCSR(void) {};
40  void AllocateMCSR(const std::string name, const int nnz, const int nrow, const int ncol);
42  void AllocateCOO(const std::string name, const int nnz, const int nrow, const int ncol);
44  void AllocateDIA(const std::string name, const int nnz, const int nrow, const int ncol, const int ndiag);
46  void AllocateELL(const std::string name, const int nnz, const int nrow, const int ncol, const int max_row);
48  void AllocateHYB(const std::string name, const int ell_nnz, const int coo_nnz, const int ell_max_row,
49  const int nrow, const int ncol);
51  void AllocateDENSE(const std::string name, const int nrow, const int ncol);
52 
54  void SetDataPtrCOO(int **row, int **col, ValueType **val,
55  std::string name, const int nnz, const int nrow, const int ncol);
57  void LeaveDataPtrCOO(int **row, int **col, ValueType **val);
58 
60  void SetDataPtrCSR(int **row_offset, int **col, ValueType **val,
61  std::string name, const int nnz, const int nrow, const int ncol);
63  void LeaveDataPtrCSR(int **row_offset, int **col, ValueType **val);
64 
66  void SetDataPtrELL(int **col, ValueType **val,
67  std::string name,
68  const int nnz, const int nrow, const int ncol, const int max_row);
70  void LeaveDataPtrELL(int **col, ValueType **val, int &max_row);
71 
73  void SetDataPtrDIA(int **offset, ValueType **val,
74  std::string name,
75  const int nnz, const int nrow, const int ncol, const int num_diag);
77  void LeaveDataPtrDIA(int **offset, ValueType **val, int &num_diag);
78 
80  void SetDataPtrDENSE(ValueType **val, std::string name, const int nrow, const int ncol);
82  void LeaveDataPtrDENSE(ValueType **val);
83 
84  void Clear(void);
85 
87  void Zeros(void);
88 
90  void Assemble(const int *i, const int *j, const ValueType *v, int size,
91  std::string name, const int n=0, const int m=0);
92 
93  void AssembleUpdate(const ValueType *v);
94 
96  void Scale(const ValueType alpha);
99  void ScaleDiagonal(const ValueType alpha);
102  void ScaleOffDiagonal(const ValueType alpha);
103 
105  void AddScalar(const ValueType alpha);
108  void AddScalarDiagonal(const ValueType alpha);
111  void AddScalarOffDiagonal(const ValueType alpha);
112 
114  void ExtractSubMatrix(const int row_offset,
115  const int col_offset,
116  const int row_size,
117  const int col_size,
118  LocalMatrix<ValueType> *mat) const;
119 
123  void ExtractSubMatrices(const int row_num_blocks,
124  const int col_num_blocks,
125  const int *row_offset,
126  const int *col_offset,
127  LocalMatrix<ValueType> ***mat) const;
128 
130  void ExtractDiagonal(LocalVector<ValueType> *vec_diag) const;
131 
133  void ExtractInverseDiagonal(LocalVector<ValueType> *vec_inv_diag) const;
134 
136  void ExtractU(LocalMatrix<ValueType> *U, const bool diag) const;
138  void ExtractL(LocalMatrix<ValueType> *L, const bool diag) const;
139 
141  void Permute(const LocalVector<int> &permutation);
142 
144  void PermuteBackward(const LocalVector<int> &permutation);
145 
147  void CMK(LocalVector<int> *permutation) const;
149  void RCMK(LocalVector<int> *permutation) const;
150 
154  void MultiColoring(int &num_colors,
155  int **size_colors,
156  LocalVector<int> *permutation) const;
157 
161  void MaximalIndependentSet(int &size,
162  LocalVector<int> *permutation) const;
163 
167  void ZeroBlockPermutation(int &size,
168  LocalVector<int> *permutation) const;
169 
171  void ILU0Factorize(void);
173  void LUFactorize(void);
174 
177  void ILUTFactorize(const double t, const int maxrow);
178 
184  void ILUpFactorize(const int p, const bool level=true);
186  void LUAnalyse(void);
188  void LUAnalyseClear(void);
191  void LUSolve(const LocalVector<ValueType> &in, LocalVector<ValueType> *out) const;
192 
194  void ICFactorize(LocalVector<ValueType> *inv_diag);
195 
197  void LLAnalyse(void);
199  void LLAnalyseClear(void);
202  void LLSolve(const LocalVector<ValueType> &in, LocalVector<ValueType> *out) const;
203  void LLSolve(const LocalVector<ValueType> &in, const LocalVector<ValueType> &inv_diag,
204  LocalVector<ValueType> *out) const;
205 
209  void LAnalyse(const bool diag_unit=false);
211  void LAnalyseClear(void);
214  void LSolve(const LocalVector<ValueType> &in, LocalVector<ValueType> *out) const;
215 
219  void UAnalyse(const bool diag_unit=false);
221  void UAnalyseClear(void);
224  void USolve(const LocalVector<ValueType> &in, LocalVector<ValueType> *out) const;
225 
227  void Householder(const int idx, ValueType &beta, LocalVector<ValueType> *vec) const;
229  void QRDecompose(void);
231  void QRSolve(const LocalVector<ValueType> &in, LocalVector<ValueType> *out) const;
232 
234  void Invert(void);
235 
237  void ReadFileMTX(const std::string filename);
239  void WriteFileMTX(const std::string filename) const;
240 
242  void ReadFileCSR(const std::string filename);
244  void WriteFileCSR(const std::string filename) const;
245 
246  virtual void MoveToAccelerator(void);
247  virtual void MoveToAcceleratorAsync(void);
248  virtual void MoveToHost(void);
249  virtual void MoveToHostAsync(void);
250  virtual void Sync(void);
251 
253  void CopyFrom(const LocalMatrix<ValueType> &src);
254 
256  void CopyFromAsync(const LocalMatrix<ValueType> &src);
257 
260  void CloneFrom(const LocalMatrix<ValueType> &src);
261 
263  void UpdateValuesCSR(ValueType *val);
264 
267  void CopyFromCSR(const int *row_offsets, const int *col, const ValueType *val);
268 
271  void CopyToCSR(int *row_offsets, int *col, ValueType *val) const;
272 
275  void CopyFromCOO(const int *row, const int *col, const ValueType *val);
276 
279  void CopyToCOO(int *row, int *col, ValueType *val) const;
280 
282  void CreateFromMap(const LocalVector<int> &map, const int n, const int m);
283 
285  void ConvertToCSR(void);
287  void ConvertToMCSR(void);
289  void ConvertToBCSR(void);
291  void ConvertToCOO(void);
293  void ConvertToELL(void);
295  void ConvertToDIA(void);
297  void ConvertToHYB(void);
299  void ConvertToDENSE(void);
301  void ConvertTo(const unsigned int matrix_format);
302 
303  virtual void Apply(const LocalVector<ValueType> &in, LocalVector<ValueType> *out) const;
304  virtual void ApplyAdd(const LocalVector<ValueType> &in, const ValueType scalar,
305  LocalVector<ValueType> *out) const;
306 
308  void SymbolicPower(const int p);
309 
313  void MatrixAdd(const LocalMatrix<ValueType> &mat, const ValueType alpha=ValueType(1.0),
314  const ValueType beta=ValueType(1.0), const bool structure=false);
315 
318 
321  void DiagonalMatrixMult(const LocalVector<ValueType> &diag);
322 
326 
330 
332  void Gershgorin(ValueType &lambda_min,
333  ValueType &lambda_max) const;
334 
337  void Compress(const double drop_off);
338 
340  void Transpose(void);
341 
343  void Sort(void);
344 
346  void ReplaceColumnVector(const int idx, const LocalVector<ValueType> &vec);
347 
349  void ReplaceRowVector(const int idx, const LocalVector<ValueType> &vec);
350 
352  void ExtractColumnVector(const int idx, LocalVector<ValueType> *vec) const;
353 
355  void ExtractRowVector(const int idx, LocalVector<ValueType> *vec) const;
356 
358  void AMGConnect(const ValueType eps, LocalVector<int> *connections) const;
360  void AMGAggregate(const LocalVector<int> &connections, LocalVector<int> *aggregates) const;
362  void AMGSmoothedAggregation(const ValueType relax,
363  const LocalVector<int> &aggregates,
364  const LocalVector<int> &connections,
365  LocalMatrix<ValueType> *prolong,
366  LocalMatrix<ValueType> *restrict) const;
368  void AMGAggregation(const LocalVector<int> &aggregates,
369  LocalMatrix<ValueType> *prolong,
370  LocalMatrix<ValueType> *restrict) const;
371 
374  void FSAI(const int power, const LocalMatrix<ValueType> *pattern);
375 
377  void SPAI(void);
378 
379 protected:
380 
381  virtual bool is_host(void) const;
382  virtual bool is_accel(void) const;
383 
384 private:
385 
389 
392 
395 
401 
402  void free_assembly_data(void);
403 
404  friend class LocalVector<ValueType>;
405 
406 };
407 
408 
409 }
410 
411 #endif // PARALUTION_LOCAL_MATRIX_HPP_
A
Definition: pcg_example.m:10
void LeaveDataPtrCSR(int **row_offset, int **col, ValueType **val)
Leave a CSR matrix to Host pointers.
Definition: local_matrix.cpp:572
void SetDataPtrCSR(int **row_offset, int **col, ValueType **val, std::string name, const int nnz, const int nrow, const int ncol)
Initialize a CSR matrix on the Host with externally allocated data.
Definition: local_matrix.cpp:542
void CopyFromAsync(const LocalMatrix< ValueType > &src)
Async copy matrix (values and structure) from another LocalMatrix.
Definition: local_matrix.cpp:1099
void Gershgorin(ValueType &lambda_min, ValueType &lambda_max) const
Compute the spectrum approximation with Gershgorin circles theorem.
Definition: local_matrix.cpp:3415
void PermuteBackward(const LocalVector< int > &permutation)
Perform (backward) permutation of the matrix.
Definition: local_matrix.cpp:3112
IndexType i
Definition: cuda_kernels_coo.hpp:195
void Sort(void)
Sort the matrix indices.
Definition: local_matrix.cpp:4112
void ConvertToMCSR(void)
Convert the matrix to MCSR structure.
Definition: local_matrix.cpp:1394
void LLAnalyseClear(void)
Delete the analysed data (see LLAnalyse)
Definition: local_matrix.cpp:2089
void CopyToCSR(int *row_offsets, int *col, ValueType *val) const
Copy (export) CSR matrix described in three arrays (offsets, columns, values). The output arrays have...
Definition: local_matrix.cpp:753
void AssembleUpdate(const ValueType *v)
Definition: local_matrix.cpp:867
virtual void Sync(void)
Definition: local_matrix.cpp:1335
void ILUTFactorize(const double t, const int maxrow)
Perform ILU(t,m) factorization based on threshold and maximum number of elements per row...
Definition: local_matrix.cpp:2447
void ConvertToDIA(void)
Convert the matrix to DIA structure.
Definition: local_matrix.cpp:1422
void ExtractL(LocalMatrix< ValueType > *L, const bool diag) const
Extract the lower triangular matrix.
Definition: local_matrix.cpp:1915
const IndexType idx
Definition: cuda_kernels_coo.hpp:115
void WriteFileCSR(const std::string filename) const
Write matrix to CSR (PARALUTION binary format) file.
Definition: local_matrix.cpp:1052
void AMGAggregation(const LocalVector< int > &aggregates, LocalMatrix< ValueType > *prolong, LocalMatrix< ValueType > *restrict) const
Aggregation-based interpolation scheme.
Definition: local_matrix.cpp:4371
void ConvertToCSR(void)
Convert the matrix to CSR structure.
Definition: local_matrix.cpp:1387
void Compress(const double drop_off)
Delete all entries in the matrix which abs(a_ij) <= drop_off; the diagonal elements are never deleted...
Definition: local_matrix.cpp:3998
void AMGAggregate(const LocalVector< int > &connections, LocalVector< int > *aggregates) const
Plain aggregation - Modification of a greedy aggregation scheme from Vanek (1996) ...
Definition: local_matrix.cpp:4224
void MultiColoring(int &num_colors, int **size_colors, LocalVector< int > *permutation) const
Perform multi-coloring decomposition of the matrix; Returns number of colors, the corresponding sizes...
Definition: local_matrix.cpp:2694
void free_assembly_data(void)
Definition: local_matrix.cpp:4764
void FSAI(const int power, const LocalMatrix< ValueType > *pattern)
Factorized Sparse Approximate Inverse assembly for given system matrix power pattern or external spar...
Definition: local_matrix.cpp:4567
void LLAnalyse(void)
Analyse the structure (level-scheduling)
Definition: local_matrix.cpp:2077
void WriteFileMTX(const std::string filename) const
Write matrix to MTX (Matrix Market Format) file.
Definition: local_matrix.cpp:973
void ConvertToHYB(void)
Convert the matrix to HYB structure.
Definition: local_matrix.cpp:1429
virtual void Apply(const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
Apply the operator, out = Operator(in), where in, out are local vectors.
Definition: local_matrix.cpp:1522
void Permute(const LocalVector< int > &permutation)
Perform (forward) permutation of the matrix.
Definition: local_matrix.cpp:3046
void Householder(const int idx, ValueType &beta, LocalVector< ValueType > *vec) const
Compute Householder vector.
Definition: local_matrix.cpp:2878
void LUAnalyse(void)
Analyse the structure (level-scheduling)
Definition: local_matrix.cpp:1988
void ExtractInverseDiagonal(LocalVector< ValueType > *vec_inv_diag) const
Extract the inverse (reciprocal) diagonal values of the matrix into a LocalVector.
Definition: local_matrix.cpp:1656
void LAnalyse(const bool diag_unit=false)
Analyse the structure (level-scheduling) L-part diag_unit == true the diag is 1; diag_unit == false t...
Definition: local_matrix.cpp:2227
void ConvertToELL(void)
Convert the matrix to ELL structure.
Definition: local_matrix.cpp:1415
void Invert(void)
Matrix inversion using QR decomposition.
Definition: local_matrix.cpp:4708
void SetDataPtrELL(int **col, ValueType **val, std::string name, const int nnz, const int nrow, const int ncol, const int max_row)
Initialize an ELL matrix on the Host with externally allocated data.
Definition: local_matrix.cpp:592
void ILU0Factorize(void)
Perform ILU(0) factorization.
Definition: local_matrix.cpp:2391
int assembly_threads
Definition: local_matrix.hpp:400
void ILUpFactorize(const int p, const bool level=true)
Perform ILU(p) factorization based on power (see power(q)-pattern method, D. Lukarski "Parallel Spars...
Definition: local_matrix.cpp:2503
virtual bool is_accel(void) const
Return true if the object is on the accelerator.
Definition: local_matrix.cpp:1200
nnz
Definition: pcg_example.m:8
void ConvertToCOO(void)
Convert the matrix to COO structure.
Definition: local_matrix.cpp:1408
void DiagonalMatrixMultL(const LocalVector< ValueType > &diag)
Multiply the matrix with diagonal matrix (stored in LocalVector), this=diag*this. ...
Definition: local_matrix.cpp:3935
unsigned int get_format(void) const
Return the matrix format id (see matrix_formats.hpp)
Definition: local_matrix.cpp:81
void ExtractColumnVector(const int idx, LocalVector< ValueType > *vec) const
Extract values from a column of a matrix to a vector.
Definition: local_matrix.cpp:4854
void ConvertToBCSR(void)
Convert the matrix to BCSR structure.
Definition: local_matrix.cpp:1401
end if j
Definition: pcg_example.m:22
void Zeros(void)
Set all the values to zero.
Definition: local_matrix.cpp:99
virtual void MoveToHost(void)
Move the object to the Host backend.
Definition: local_matrix.cpp:1269
virtual void ApplyAdd(const LocalVector< ValueType > &in, const ValueType scalar, LocalVector< ValueType > *out) const
Apply and add the operator, out = out + scalar*Operator(in), where in, out are local vectors...
Definition: local_matrix.cpp:1560
const IndexType const IndexType const IndexType const ValueType const ValueType scalar
Definition: cuda_kernels_coo.hpp:91
void ExtractDiagonal(LocalVector< ValueType > *vec_diag) const
Extract the diagonal values of the matrix into a LocalVector.
Definition: local_matrix.cpp:1599
void MaximalIndependentSet(int &size, LocalVector< int > *permutation) const
Perform maximal independent set decomposition of the matrix; Returns the size of the maximal independ...
Definition: local_matrix.cpp:2757
void UAnalyseClear(void)
Delete the analysed data (see UAnalyse) U-party.
Definition: local_matrix.cpp:2321
void QRSolve(const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
Solve QR out = in.
Definition: local_matrix.cpp:2986
void ICFactorize(LocalVector< ValueType > *inv_diag)
Perform IC(0) factorization.
Definition: local_matrix.cpp:2631
void ExtractU(LocalMatrix< ValueType > *U, const bool diag) const
Extract the upper triangular matrix.
Definition: local_matrix.cpp:1842
void CreateFromMap(const LocalVector< int > &map, const int n, const int m)
Create a restriction matrix operator based on an int vector map.
Definition: local_matrix.cpp:4448
void LeaveDataPtrDIA(int **offset, ValueType **val, int &num_diag)
Leave a DIA matrix to Host pointers.
Definition: local_matrix.cpp:674
void CopyFromCSR(const int *row_offsets, const int *col, const ValueType *val)
Copy (import) CSR matrix described in three arrays (offsets, columns, values). The object data has to...
Definition: local_matrix.cpp:735
void AllocateBCSR(void)
Allocate BCSR Matrix - not implemented.
Definition: local_matrix.hpp:38
AcceleratorMatrix< ValueType > * matrix_accel_
Accelerator Matrix.
Definition: local_matrix.hpp:394
void CopyFromCOO(const int *row, const int *col, const ValueType *val)
Copy (import) COO matrix described in three arrays (rows, columns, values). The object data has to be...
Definition: local_matrix.cpp:769
void ExtractSubMatrix(const int row_offset, const int col_offset, const int row_size, const int col_size, LocalMatrix< ValueType > *mat) const
Extract a sub-matrix with row/col_offset and row/col_size.
Definition: local_matrix.cpp:1713
void AllocateHYB(const std::string name, const int ell_nnz, const int coo_nnz, const int ell_max_row, const int nrow, const int ncol)
Allocate HYB Matrix.
Definition: local_matrix.cpp:359
bool Check(void) const
Return true if the matrix is ok (empty matrix is also ok)
Definition: local_matrix.cpp:441
void Assemble(const int *i, const int *j, const ValueType *v, int size, std::string name, const int n=0, const int m=0)
Assembling.
Definition: local_matrix.cpp:803
int * assembly_loop_start
Definition: local_matrix.hpp:398
void CloneFrom(const LocalMatrix< ValueType > &src)
Clone the entire matrix (values,structure+backend descr) from another LocalMatrix.
Definition: local_matrix.cpp:1115
void AllocateMCSR(const std::string name, const int nnz, const int nrow, const int ncol)
Allocate MCSR Matrix.
Definition: local_matrix.cpp:277
void Scale(const ValueType alpha)
Scale all the values in the matrix.
Definition: local_matrix.cpp:3458
void power(const int mic_dev, const int size, const double val, ValueType *vec)
Definition: mic_vector_kernel.cpp:241
void ConvertTo(const unsigned int matrix_format)
Convert the matrix to specified matrix ID format.
Definition: local_matrix.cpp:1443
void SetDataPtrDENSE(ValueType **val, std::string name, const int nrow, const int ncol)
Initialize a DENSE matrix on the Host with externally allocated data.
Definition: local_matrix.cpp:693
virtual int get_nnz(void) const
Return the number of non-zeros in the matrix/stencil.
Definition: local_matrix.cpp:74
Definition: backend_manager.hpp:12
void ConvertToDENSE(void)
Convert the matrix to DENSE structure.
Definition: local_matrix.cpp:1436
void SymbolicPower(const int p)
Perform symbolic computation (structure only) of |this|^p.
Definition: local_matrix.cpp:3296
void RCMK(LocalVector< int > *permutation) const
Create permutation vector for reverse CMK reordering of the matrix.
Definition: local_matrix.cpp:3237
void DiagonalMatrixMult(const LocalVector< ValueType > &diag)
Multiply the matrix with diagonal matrix (stored in LocalVector), as DiagonalMatrixMultR() ...
Definition: local_matrix.cpp:3928
void ScaleDiagonal(const ValueType alpha)
Scale the diagonal entries of the matrix with alpha, the diagonal elements must exist.
Definition: local_matrix.cpp:3512
void SPAI(void)
SParse Approximate Inverse assembly for given system matrix pattern.
Definition: local_matrix.cpp:4652
virtual void MoveToHostAsync(void)
Move the object to the Host backend with async move.
Definition: local_matrix.cpp:1315
void MatrixMult(const LocalMatrix< ValueType > &A, const LocalMatrix< ValueType > &B)
Multiply two matrices, this = A * B.
Definition: local_matrix.cpp:3782
void AddScalarOffDiagonal(const ValueType alpha)
Add alpha to the off-diagonal entries of the matrix, the diagonal elements must exist.
Definition: local_matrix.cpp:3728
void LLSolve(const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
Solve LL^T out = in; if level-scheduling algorithm is provided then the graph traversing is performed...
Definition: local_matrix.cpp:2101
int * assembly_loop_end
Definition: local_matrix.hpp:399
void AllocateDENSE(const std::string name, const int nrow, const int ncol)
Allocate DENSE Matrix.
Definition: local_matrix.cpp:403
void CopyToCOO(int *row, int *col, ValueType *val) const
Copy (export) COO matrix described in three arrays (rows, columns, values). The output arrays have to...
Definition: local_matrix.cpp:787
void DiagonalMatrixMultR(const LocalVector< ValueType > &diag)
Multiply the matrix with diagonal matrix (stored in LocalVector), this=this*diag. ...
Definition: local_matrix.cpp:3865
virtual int get_nrow(void) const
Return the number of rows in the matrix/stencil.
Definition: local_matrix.cpp:60
void LeaveDataPtrELL(int **col, ValueType **val, int &max_row)
Leave an ELL matrix to Host pointers.
Definition: local_matrix.cpp:621
Base class for all host/accelerator matrices.
Definition: base_matrix.hpp:92
void ReplaceColumnVector(const int idx, const LocalVector< ValueType > &vec)
Replace a column vector of a matrix.
Definition: local_matrix.cpp:4783
int * assembly_rank
Definition: local_matrix.hpp:396
void LUAnalyseClear(void)
Delete the analysed data (see LUAnalyse)
Definition: local_matrix.cpp:2000
void USolve(const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
Solve U out = in; if level-scheduling algorithm is provided then the graph traversing is performed in...
Definition: local_matrix.cpp:2333
void LeaveDataPtrCOO(int **row, int **col, ValueType **val)
Leave a COO matrix to Host pointers.
Definition: local_matrix.cpp:522
void ScaleOffDiagonal(const ValueType alpha)
Scale the off-diagonal entries of the matrix with alpha, the diagonal elements must exist...
Definition: local_matrix.cpp:3566
virtual void MoveToAccelerator(void)
Move the object to the Accelerator backend.
Definition: local_matrix.cpp:1243
void LUSolve(const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
Solve LU out = in; if level-scheduling algorithm is provided then the graph traversing is performed i...
Definition: local_matrix.cpp:2012
void AMGConnect(const ValueType eps, LocalVector< int > *connections) const
Strong couplings for aggregation-based AMG.
Definition: local_matrix.cpp:4166
void ReplaceRowVector(const int idx, const LocalVector< ValueType > &vec)
Replace a row vector of a matrix.
Definition: local_matrix.cpp:4910
HostMatrix< ValueType > * matrix_host_
Host Matrix.
Definition: local_matrix.hpp:391
Operator class defines the generic interface for applying an operator (e.g. matrix, stencil) from/to global and local vectors.
Definition: operator.hpp:19
void AllocateCSR(const std::string name, const int nnz, const int nrow, const int ncol)
Allocate CSR Matrix.
Definition: local_matrix.cpp:153
Definition: backend_manager.cpp:43
void MatrixAdd(const LocalMatrix< ValueType > &mat, const ValueType alpha=ValueType(1.0), const ValueType beta=ValueType(1.0), const bool structure=false)
Perform matrix addition, this = alpha*this + beta*mat; if structure==false the sparsity pattern of th...
Definition: local_matrix.cpp:3354
int * assembly_irank
Definition: local_matrix.hpp:397
virtual ~LocalMatrix()
Definition: local_matrix.cpp:49
void AMGSmoothedAggregation(const ValueType relax, const LocalVector< int > &aggregates, const LocalVector< int > &connections, LocalMatrix< ValueType > *prolong, LocalMatrix< ValueType > *restrict) const
Interpolation scheme based on smoothed aggregation from Vanek (1996)
Definition: local_matrix.cpp:4284
Definition: backend_manager.hpp:14
void SetDataPtrCOO(int **row, int **col, ValueType **val, std::string name, const int nnz, const int nrow, const int ncol)
Initialize a COO matrix on the Host with externally allocated data.
Definition: local_matrix.cpp:492
virtual void MoveToAcceleratorAsync(void)
Move the object to the Accelerator backend with async move.
Definition: local_matrix.cpp:1292
virtual int get_ncol(void) const
Return the number of columns in the matrix/stencil.
Definition: local_matrix.cpp:67
void ReadFileCSR(const std::string filename)
Read matrix from CSR (PARALUTION binary format) file.
Definition: local_matrix.cpp:1007
void UpdateValuesCSR(ValueType *val)
Update CSR matrix entries only, structure will remain the same.
Definition: local_matrix.cpp:1151
Definition: local_matrix.hpp:15
void CopyFrom(const LocalMatrix< ValueType > &src)
Copy matrix (values and structure) from another LocalMatrix.
Definition: local_matrix.cpp:1086
void LUFactorize(void)
Perform LU factorization.
Definition: local_matrix.cpp:4511
void LSolve(const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
Solve L out = in; if level-scheduling algorithm is provided then the graph traversing is performed in...
Definition: local_matrix.cpp:2251
void Clear(void)
Clear (free all data) the object.
Definition: local_matrix.cpp:88
virtual void info(void) const
Print the object information (properties, backends)
Definition: local_matrix.cpp:1207
void ReadFileMTX(const std::string filename)
Read matrix from MTX (Matrix Market Format) file.
Definition: local_matrix.cpp:929
void ExtractRowVector(const int idx, LocalVector< ValueType > *vec) const
Extract values from a row of a matrix to a vector.
Definition: local_matrix.cpp:4981
void CMK(LocalVector< int > *permutation) const
Create permutation vector for CMK reordering of the matrix.
Definition: local_matrix.cpp:3178
void Transpose(void)
Transpose the matrix.
Definition: local_matrix.cpp:4056
void LeaveDataPtrDENSE(ValueType **val)
Leave a DENSE matrix to Host pointers.
Definition: local_matrix.cpp:717
void QRDecompose(void)
QR Decomposition.
Definition: local_matrix.cpp:2930
void ExtractSubMatrices(const int row_num_blocks, const int col_num_blocks, const int *row_offset, const int *col_offset, LocalMatrix< ValueType > ***mat) const
Extract array of non-overlapping sub-matrices (row/col_num_blocks define the blocks for rows/columns;...
Definition: local_matrix.cpp:1808
void AllocateELL(const std::string name, const int nnz, const int nrow, const int ncol, const int max_row)
Allocate ELL Matrix.
Definition: local_matrix.cpp:318
Definition: host_vector.hpp:13
void SetDataPtrDIA(int **offset, ValueType **val, std::string name, const int nnz, const int nrow, const int ncol, const int num_diag)
Initialize a DIA matrix on the Host with externally allocated data.
Definition: local_matrix.cpp:640
void AllocateCOO(const std::string name, const int nnz, const int nrow, const int ncol)
Allocate COO Matrix.
Definition: local_matrix.cpp:194
void UAnalyse(const bool diag_unit=false)
Analyse the structure (level-scheduling) U-part; diag_unit == true the diag is 1; diag_unit == false ...
Definition: local_matrix.cpp:2309
LocalMatrix()
Definition: local_matrix.cpp:25
void AddScalarDiagonal(const ValueType alpha)
Add alpha to the diagonal entries of the matrix, the diagonal elements must exist.
Definition: local_matrix.cpp:3674
BaseMatrix< ValueType > * matrix_
Pointer from the base matrix class to the current allocated matrix (host_ or accel_) ...
Definition: local_matrix.hpp:388
void ZeroBlockPermutation(int &size, LocalVector< int > *permutation) const
Return a permutation for saddle-point problems (zero diagonal entries), where all zero diagonal eleme...
Definition: local_matrix.cpp:2818
void AllocateDIA(const std::string name, const int nnz, const int nrow, const int ncol, const int ndiag)
Allocate DIA Matrix.
Definition: local_matrix.cpp:236
virtual bool is_host(void) const
Return true if the object is on the host.
Definition: local_matrix.cpp:1193
void LAnalyseClear(void)
Delete the analysed data (see LAnalyse) L-party.
Definition: local_matrix.cpp:2239
void AddScalar(const ValueType alpha)
Add a scalar to all values.
Definition: local_matrix.cpp:3620