1 #ifndef PARALUTION_LOCAL_MATRIX_HPP_
2 #define PARALUTION_LOCAL_MATRIX_HPP_
10 template <
typename ValueType>
14 template <
typename ValueType>
22 virtual void info(
void)
const;
29 virtual int get_nnz(
void)
const;
33 bool Check(
void)
const;
36 void AllocateCSR(
const std::string name,
const int nnz,
const int nrow,
const int ncol);
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);
55 std::string name,
const int nnz,
const int nrow,
const int ncol);
60 void SetDataPtrCSR(
int **row_offset,
int **col, ValueType **val,
61 std::string name,
const int nnz,
const int nrow,
const int ncol);
68 const int nnz,
const int nrow,
const int ncol,
const int max_row);
75 const int nnz,
const int nrow,
const int ncol,
const int num_diag);
80 void SetDataPtrDENSE(ValueType **val, std::string name,
const int nrow,
const int ncol);
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);
96 void Scale(
const ValueType alpha);
115 const int col_offset,
124 const int col_num_blocks,
125 const int *row_offset,
126 const int *col_offset,
209 void LAnalyse(
const bool diag_unit=
false);
219 void UAnalyse(
const bool diag_unit=
false);
250 virtual void Sync(
void);
267 void CopyFromCSR(
const int *row_offsets,
const int *col,
const ValueType *val);
271 void CopyToCSR(
int *row_offsets,
int *col, ValueType *val)
const;
275 void CopyFromCOO(
const int *row,
const int *col,
const ValueType *val);
279 void CopyToCOO(
int *row,
int *col, ValueType *val)
const;
301 void ConvertTo(
const unsigned int matrix_format);
314 const ValueType beta=ValueType(1.0),
const bool structure=
false);
333 ValueType &lambda_max)
const;
337 void Compress(
const double drop_off);
381 virtual bool is_host(
void)
const;
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