PARALUTION  1.0.0
PARALUTION
paralution::LocalMatrix< ValueType > Class Template Reference

#include <local_matrix.hpp>

Inheritance diagram for paralution::LocalMatrix< ValueType >:
paralution::Operator< ValueType > paralution::BaseParalution< ValueType > paralution::ParalutionObj

Public Member Functions

 LocalMatrix ()
 
virtual ~LocalMatrix ()
 
virtual void info (void) const
 Print the object information (properties, backends) More...
 
unsigned int get_format (void) const
 Return the matrix format id (see matrix_formats.hpp) More...
 
virtual int get_nrow (void) const
 Return the number of rows in the matrix/stencil. More...
 
virtual int get_ncol (void) const
 Return the number of columns in the matrix/stencil. More...
 
virtual int get_nnz (void) const
 Return the number of non-zeros in the matrix/stencil. More...
 
bool Check (void) const
 Return true if the matrix is ok (empty matrix is also ok) More...
 
void AllocateCSR (const std::string name, const int nnz, const int nrow, const int ncol)
 Allocate CSR Matrix. More...
 
void AllocateBCSR (void)
 Allocate BCSR Matrix - not implemented. More...
 
void AllocateMCSR (const std::string name, const int nnz, const int nrow, const int ncol)
 Allocate MCSR Matrix. More...
 
void AllocateCOO (const std::string name, const int nnz, const int nrow, const int ncol)
 Allocate COO Matrix. More...
 
void AllocateDIA (const std::string name, const int nnz, const int nrow, const int ncol, const int ndiag)
 Allocate DIA Matrix. More...
 
void AllocateELL (const std::string name, const int nnz, const int nrow, const int ncol, const int max_row)
 Allocate ELL Matrix. More...
 
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. More...
 
void AllocateDENSE (const std::string name, const int nrow, const int ncol)
 Allocate DENSE Matrix. More...
 
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. More...
 
void LeaveDataPtrCOO (int **row, int **col, ValueType **val)
 Leave a COO matrix to Host pointers. More...
 
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. More...
 
void LeaveDataPtrCSR (int **row_offset, int **col, ValueType **val)
 Leave a CSR matrix to Host pointers. More...
 
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. More...
 
void LeaveDataPtrELL (int **col, ValueType **val, int &max_row)
 Leave an ELL matrix to Host pointers. More...
 
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. More...
 
void LeaveDataPtrDIA (int **offset, ValueType **val, int &num_diag)
 Leave a DIA matrix to Host pointers. More...
 
void SetDataPtrDENSE (ValueType **val, std::string name, const int nrow, const int ncol)
 Initialize a DENSE matrix on the Host with externally allocated data. More...
 
void LeaveDataPtrDENSE (ValueType **val)
 Leave a DENSE matrix to Host pointers. More...
 
void Clear (void)
 Clear (free all data) the object. More...
 
void Zeros (void)
 Set all the values to zero. More...
 
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. More...
 
void AssembleUpdate (const ValueType *v)
 
void Scale (const ValueType alpha)
 Scale all the values in the matrix. More...
 
void ScaleDiagonal (const ValueType alpha)
 Scale the diagonal entries of the matrix with alpha, the diagonal elements must exist. More...
 
void ScaleOffDiagonal (const ValueType alpha)
 Scale the off-diagonal entries of the matrix with alpha, the diagonal elements must exist. More...
 
void AddScalar (const ValueType alpha)
 Add a scalar to all values. More...
 
void AddScalarDiagonal (const ValueType alpha)
 Add alpha to the diagonal entries of the matrix, the diagonal elements must exist. More...
 
void AddScalarOffDiagonal (const ValueType alpha)
 Add alpha to the off-diagonal entries of the matrix, the diagonal elements must exist. More...
 
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. More...
 
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; row/col_offset have sizes col/row_num_blocks+1, where [i+1]-[i] defines the i-th size of the sub-matrix) More...
 
void ExtractDiagonal (LocalVector< ValueType > *vec_diag) const
 Extract the diagonal values of the matrix into a LocalVector. More...
 
void ExtractInverseDiagonal (LocalVector< ValueType > *vec_inv_diag) const
 Extract the inverse (reciprocal) diagonal values of the matrix into a LocalVector. More...
 
void ExtractU (LocalMatrix< ValueType > *U, const bool diag) const
 Extract the upper triangular matrix. More...
 
void ExtractL (LocalMatrix< ValueType > *L, const bool diag) const
 Extract the lower triangular matrix. More...
 
void Permute (const LocalVector< int > &permutation)
 Perform (forward) permutation of the matrix. More...
 
void PermuteBackward (const LocalVector< int > &permutation)
 Perform (backward) permutation of the matrix. More...
 
void CMK (LocalVector< int > *permutation) const
 Create permutation vector for CMK reordering of the matrix. More...
 
void RCMK (LocalVector< int > *permutation) const
 Create permutation vector for reverse CMK reordering of the matrix. More...
 
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 (the array is allocated in the function) and the permutation. More...
 
void MaximalIndependentSet (int &size, LocalVector< int > *permutation) const
 Perform maximal independent set decomposition of the matrix; Returns the size of the maximal independent set and the corresponding permutation. More...
 
void ZeroBlockPermutation (int &size, LocalVector< int > *permutation) const
 Return a permutation for saddle-point problems (zero diagonal entries), where all zero diagonal elements are mapped to the last block; the return size is the size of the first block. More...
 
void ILU0Factorize (void)
 Perform ILU(0) factorization. More...
 
void LUFactorize (void)
 Perform LU factorization. More...
 
void ILUTFactorize (const double t, const int maxrow)
 Perform ILU(t,m) factorization based on threshold and maximum number of elements per row. More...
 
void ILUpFactorize (const int p, const bool level=true)
 Perform ILU(p) factorization based on power (see power(q)-pattern method, D. Lukarski "Parallel Sparse Linear Algebra for Multi-core and Many-core Platforms - Parallel Solvers and Preconditioners", PhD Thesis, 2012, KIT); level==true build the structure based on levels; level==false build the structure only based on the power(p+1) More...
 
void LUAnalyse (void)
 Analyse the structure (level-scheduling) More...
 
void LUAnalyseClear (void)
 Delete the analysed data (see LUAnalyse) More...
 
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 in parallel. More...
 
void ICFactorize (LocalVector< ValueType > *inv_diag)
 Perform IC(0) factorization. More...
 
void LLAnalyse (void)
 Analyse the structure (level-scheduling) More...
 
void LLAnalyseClear (void)
 Delete the analysed data (see LLAnalyse) More...
 
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 in parallel. More...
 
void LLSolve (const LocalVector< ValueType > &in, const LocalVector< ValueType > &inv_diag, LocalVector< ValueType > *out) const
 
void LAnalyse (const bool diag_unit=false)
 Analyse the structure (level-scheduling) L-part diag_unit == true the diag is 1; diag_unit == false the diag is 0;. More...
 
void LAnalyseClear (void)
 Delete the analysed data (see LAnalyse) L-party. More...
 
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 parallel. More...
 
void UAnalyse (const bool diag_unit=false)
 Analyse the structure (level-scheduling) U-part; diag_unit == true the diag is 1; diag_unit == false the diag is 0;. More...
 
void UAnalyseClear (void)
 Delete the analysed data (see UAnalyse) U-party. More...
 
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 parallel. More...
 
void Householder (const int idx, ValueType &beta, LocalVector< ValueType > *vec) const
 Compute Householder vector. More...
 
void QRDecompose (void)
 QR Decomposition. More...
 
void QRSolve (const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
 Solve QR out = in. More...
 
void Invert (void)
 Matrix inversion using QR decomposition. More...
 
void ReadFileMTX (const std::string filename)
 Read matrix from MTX (Matrix Market Format) file. More...
 
void WriteFileMTX (const std::string filename) const
 Write matrix to MTX (Matrix Market Format) file. More...
 
void ReadFileCSR (const std::string filename)
 Read matrix from CSR (PARALUTION binary format) file. More...
 
void WriteFileCSR (const std::string filename) const
 Write matrix to CSR (PARALUTION binary format) file. More...
 
virtual void MoveToAccelerator (void)
 Move the object to the Accelerator backend. More...
 
virtual void MoveToAcceleratorAsync (void)
 Move the object to the Accelerator backend with async move. More...
 
virtual void MoveToHost (void)
 Move the object to the Host backend. More...
 
virtual void MoveToHostAsync (void)
 Move the object to the Host backend with async move. More...
 
virtual void Sync (void)
 
void CopyFrom (const LocalMatrix< ValueType > &src)
 Copy matrix (values and structure) from another LocalMatrix. More...
 
void CopyFromAsync (const LocalMatrix< ValueType > &src)
 Async copy matrix (values and structure) from another LocalMatrix. More...
 
void CloneFrom (const LocalMatrix< ValueType > &src)
 Clone the entire matrix (values,structure+backend descr) from another LocalMatrix. More...
 
void UpdateValuesCSR (ValueType *val)
 Update CSR matrix entries only, structure will remain the same. More...
 
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 be allocated (call AllocateCSR first) More...
 
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 to be allocated. More...
 
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 allocated (call AllocateCOO first) More...
 
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 be allocated. More...
 
void CreateFromMap (const LocalVector< int > &map, const int n, const int m)
 Create a restriction matrix operator based on an int vector map. More...
 
void ConvertToCSR (void)
 Convert the matrix to CSR structure. More...
 
void ConvertToMCSR (void)
 Convert the matrix to MCSR structure. More...
 
void ConvertToBCSR (void)
 Convert the matrix to BCSR structure. More...
 
void ConvertToCOO (void)
 Convert the matrix to COO structure. More...
 
void ConvertToELL (void)
 Convert the matrix to ELL structure. More...
 
void ConvertToDIA (void)
 Convert the matrix to DIA structure. More...
 
void ConvertToHYB (void)
 Convert the matrix to HYB structure. More...
 
void ConvertToDENSE (void)
 Convert the matrix to DENSE structure. More...
 
void ConvertTo (const unsigned int matrix_format)
 Convert the matrix to specified matrix ID format. More...
 
virtual void Apply (const LocalVector< ValueType > &in, LocalVector< ValueType > *out) const
 Apply the operator, out = Operator(in), where in, out are local vectors. More...
 
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. More...
 
void SymbolicPower (const int p)
 Perform symbolic computation (structure only) of |this|^p. More...
 
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 the matrix is not changed; if structure==true a new sparsity pattern is computed. More...
 
void MatrixMult (const LocalMatrix< ValueType > &A, const LocalMatrix< ValueType > &B)
 Multiply two matrices, this = A * B. More...
 
void DiagonalMatrixMult (const LocalVector< ValueType > &diag)
 Multiply the matrix with diagonal matrix (stored in LocalVector), as DiagonalMatrixMultR() More...
 
void DiagonalMatrixMultL (const LocalVector< ValueType > &diag)
 Multiply the matrix with diagonal matrix (stored in LocalVector), this=diag*this. More...
 
void DiagonalMatrixMultR (const LocalVector< ValueType > &diag)
 Multiply the matrix with diagonal matrix (stored in LocalVector), this=this*diag. More...
 
void Gershgorin (ValueType &lambda_min, ValueType &lambda_max) const
 Compute the spectrum approximation with Gershgorin circles theorem. More...
 
void Compress (const double drop_off)
 Delete all entries in the matrix which abs(a_ij) <= drop_off; the diagonal elements are never deleted;. More...
 
void Transpose (void)
 Transpose the matrix. More...
 
void Sort (void)
 Sort the matrix indices. More...
 
void ReplaceColumnVector (const int idx, const LocalVector< ValueType > &vec)
 Replace a column vector of a matrix. More...
 
void ReplaceRowVector (const int idx, const LocalVector< ValueType > &vec)
 Replace a row vector of a matrix. More...
 
void ExtractColumnVector (const int idx, LocalVector< ValueType > *vec) const
 Extract values from a column of a matrix to a vector. More...
 
void ExtractRowVector (const int idx, LocalVector< ValueType > *vec) const
 Extract values from a row of a matrix to a vector. More...
 
void AMGConnect (const ValueType eps, LocalVector< int > *connections) const
 Strong couplings for aggregation-based AMG. More...
 
void AMGAggregate (const LocalVector< int > &connections, LocalVector< int > *aggregates) const
 Plain aggregation - Modification of a greedy aggregation scheme from Vanek (1996) More...
 
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) More...
 
void AMGAggregation (const LocalVector< int > &aggregates, LocalMatrix< ValueType > *prolong, LocalMatrix< ValueType > *restrict) const
 Aggregation-based interpolation scheme. More...
 
void FSAI (const int power, const LocalMatrix< ValueType > *pattern)
 Factorized Sparse Approximate Inverse assembly for given system matrix power pattern or external sparsity pattern. More...
 
void SPAI (void)
 SParse Approximate Inverse assembly for given system matrix pattern. More...
 
virtual int get_local_nrow (void) const
 Return the number of rows in the local matrix/stencil. More...
 
virtual int get_local_ncol (void) const
 Return the number of columns in the local matrix/stencil. More...
 
virtual int get_local_nnz (void) const
 Return the number of non-zeros in the local matrix/stencil. More...
 
virtual int get_ghost_nrow (void) const
 Return the number of rows in the ghost matrix/stencil. More...
 
virtual int get_ghost_ncol (void) const
 Return the number of columns in the ghost matrix/stencil. More...
 
virtual int get_ghost_nnz (void) const
 Return the number of non-zeros in the ghost matrix/stencil. More...
 
void CloneBackend (const BaseParalution< ValueType > &src)
 Clone the Backend descriptor from another object. More...
 
template<typename ValueType2 >
void CloneBackend (const BaseParalution< ValueType2 > &src)
 Clone the Backend descriptor from another object with different template ValueType. More...
 

Protected Member Functions

virtual bool is_host (void) const
 Return true if the object is on the host. More...
 
virtual bool is_accel (void) const
 Return true if the object is on the accelerator. More...
 

Protected Attributes

std::string object_name_
 Name of the object. More...
 
Paralution_Backend_Descriptor local_backend_
 Backend descriptor. More...
 
bool asyncf
 
size_t global_obj_id
 

Private Member Functions

void free_assembly_data (void)
 

Private Attributes

BaseMatrix< ValueType > * matrix_
 Pointer from the base matrix class to the current allocated matrix (host_ or accel_) More...
 
HostMatrix< ValueType > * matrix_host_
 Host Matrix. More...
 
AcceleratorMatrix< ValueType > * matrix_accel_
 Accelerator Matrix. More...
 
intassembly_rank
 
intassembly_irank
 
intassembly_loop_start
 
intassembly_loop_end
 
int assembly_threads
 

Friends

class LocalVector< ValueType >
 

Constructor & Destructor Documentation

template<typename ValueType >
paralution::LocalMatrix< ValueType >::LocalMatrix ( )
template<typename ValueType >
paralution::LocalMatrix< ValueType >::~LocalMatrix ( )
virtual

Member Function Documentation

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::AddScalar ( const ValueType  alpha)

Add a scalar to all values.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::AddScalarDiagonal ( const ValueType  alpha)

Add alpha to the diagonal entries of the matrix, the diagonal elements must exist.

Referenced by main().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::AddScalarOffDiagonal ( const ValueType  alpha)

Add alpha to the off-diagonal entries of the matrix, the diagonal elements must exist.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::AllocateBCSR ( void  )
inline

Allocate BCSR Matrix - not implemented.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::AllocateCOO ( const std::string  name,
const int  nnz,
const int  nrow,
const int  ncol 
)

Allocate COO Matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::AllocateCSR ( const std::string  name,
const int  nnz,
const int  nrow,
const int  ncol 
)

Allocate CSR Matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::AllocateDENSE ( const std::string  name,
const int  nrow,
const int  ncol 
)

Allocate DENSE Matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::AllocateDIA ( const std::string  name,
const int  nnz,
const int  nrow,
const int  ncol,
const int  ndiag 
)

Allocate DIA Matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::AllocateELL ( const std::string  name,
const int  nnz,
const int  nrow,
const int  ncol,
const int  max_row 
)

Allocate ELL Matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::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.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::AllocateMCSR ( const std::string  name,
const int  nnz,
const int  nrow,
const int  ncol 
)

Allocate MCSR Matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::AMGAggregate ( const LocalVector< int > &  connections,
LocalVector< int > *  aggregates 
) const

Plain aggregation - Modification of a greedy aggregation scheme from Vanek (1996)

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::AMGAggregation ( const LocalVector< int > &  aggregates,
LocalMatrix< ValueType > *  prolong,
LocalMatrix< ValueType > *  restrict 
) const

Aggregation-based interpolation scheme.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::AMGConnect ( const ValueType  eps,
LocalVector< int > *  connections 
) const

Strong couplings for aggregation-based AMG.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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)

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::Apply ( const LocalVector< ValueType > &  in,
LocalVector< ValueType > *  out 
) const
virtual

Apply the operator, out = Operator(in), where in, out are local vectors.

Reimplemented from paralution::Operator< ValueType >.

Referenced by main().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ApplyAdd ( const LocalVector< ValueType > &  in,
const ValueType  scalar,
LocalVector< ValueType > *  out 
) const
virtual

Apply and add the operator, out = out + scalar*Operator(in), where in, out are local vectors.

Reimplemented from paralution::Operator< ValueType >.

Referenced by main().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::Assemble ( const int i,
const int j,
const ValueType *  v,
int  size,
std::string  name,
const int  n = 0,
const int  m = 0 
)

Assembling.

Referenced by main().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::AssembleUpdate ( const ValueType *  v)

Referenced by main().

template<typename ValueType >
bool paralution::LocalMatrix< ValueType >::Check ( void  ) const

Return true if the matrix is ok (empty matrix is also ok)

Referenced by main().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Clear ( void  )
virtual
template<typename ValueType >
template<typename ValueType2 >
void paralution::BaseParalution< ValueType >::CloneBackend ( const BaseParalution< ValueType2 > &  src)
inherited

Clone the Backend descriptor from another object with different template ValueType.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::CloneFrom ( const LocalMatrix< ValueType > &  src)

Clone the entire matrix (values,structure+backend descr) from another LocalMatrix.

Referenced by paralution::LocalMatrix< ValueType >::ILUpFactorize().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::CMK ( LocalVector< int > *  permutation) const

Create permutation vector for CMK reordering of the matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Compress ( const double  drop_off)

Delete all entries in the matrix which abs(a_ij) <= drop_off; the diagonal elements are never deleted;.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertTo ( const unsigned int  matrix_format)

Convert the matrix to specified matrix ID format.

Referenced by paralution::LocalMatrix< ValueType >::AMGAggregate(), paralution::LocalMatrix< ValueType >::AMGAggregation(), paralution::LocalMatrix< ValueType >::AMGConnect(), paralution::LocalMatrix< ValueType >::AMGSmoothedAggregation(), paralution::LocalMatrix< ValueType >::Apply(), paralution::LocalMatrix< ValueType >::ApplyAdd(), paralution::LocalMatrix< ValueType >::Check(), paralution::LocalMatrix< ValueType >::CMK(), paralution::LocalMatrix< ValueType >::ExtractColumnVector(), paralution::LocalMatrix< ValueType >::ExtractDiagonal(), paralution::LocalMatrix< ValueType >::ExtractInverseDiagonal(), paralution::LocalMatrix< ValueType >::ExtractL(), paralution::LocalMatrix< ValueType >::ExtractRowVector(), paralution::LocalMatrix< ValueType >::ExtractSubMatrix(), paralution::LocalMatrix< ValueType >::ExtractU(), paralution::LocalMatrix< ValueType >::Gershgorin(), paralution::LocalMatrix< ValueType >::Householder(), paralution::LocalMatrix< ValueType >::LLSolve(), paralution::LocalMatrix< ValueType >::LSolve(), paralution::LocalMatrix< ValueType >::LUSolve(), paralution::LocalMatrix< ValueType >::MatrixAdd(), paralution::LocalMatrix< ValueType >::MatrixMult(), paralution::LocalMatrix< ValueType >::MaximalIndependentSet(), paralution::LocalMatrix< ValueType >::MultiColoring(), paralution::LocalMatrix< ValueType >::QRSolve(), paralution::LocalMatrix< ValueType >::RCMK(), paralution::LocalMatrix< ValueType >::USolve(), paralution::LocalMatrix< ValueType >::WriteFileCSR(), paralution::LocalMatrix< ValueType >::WriteFileMTX(), and paralution::LocalMatrix< ValueType >::ZeroBlockPermutation().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertToBCSR ( void  )

Convert the matrix to BCSR structure.

Referenced by paralution_fortran_solve().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertToCOO ( void  )

Convert the matrix to COO structure.

Referenced by main(), paralution_fortran_solve(), and paralution::LocalMatrix< ValueType >::WriteFileMTX().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertToCSR ( void  )

Convert the matrix to CSR structure.

Referenced by paralution::LocalMatrix< ValueType >::AMGAggregate(), paralution::LocalMatrix< ValueType >::AMGAggregation(), paralution::LocalMatrix< ValueType >::AMGConnect(), paralution::LocalMatrix< ValueType >::AMGSmoothedAggregation(), paralution::LocalMatrix< ValueType >::Check(), paralution::LocalMatrix< ValueType >::CMK(), paralution::LocalMatrix< ValueType >::ExtractColumnVector(), paralution::LocalMatrix< ValueType >::ExtractDiagonal(), paralution::LocalMatrix< ValueType >::ExtractInverseDiagonal(), paralution::LocalMatrix< ValueType >::ExtractL(), paralution::LocalMatrix< ValueType >::ExtractRowVector(), paralution::LocalMatrix< ValueType >::ExtractSubMatrix(), paralution::LocalMatrix< ValueType >::ExtractU(), paralution::LocalMatrix< ValueType >::Gershgorin(), paralution::LocalMatrix< ValueType >::ILUpFactorize(), import_dealii_matrix(), paralution::LocalMatrix< ValueType >::LLSolve(), paralution::LocalMatrix< ValueType >::LSolve(), paralution::LocalMatrix< ValueType >::LUSolve(), main(), paralution::LocalMatrix< ValueType >::MatrixAdd(), paralution::LocalMatrix< ValueType >::MatrixMult(), paralution::LocalMatrix< ValueType >::MaximalIndependentSet(), paralution::LocalMatrix< ValueType >::MultiColoring(), paralution_fortran_solve(), paralution_fortran_solve_coo(), paralution::LocalMatrix< ValueType >::RCMK(), paralution::LocalMatrix< ValueType >::USolve(), paralution::LocalMatrix< ValueType >::WriteFileCSR(), and paralution::LocalMatrix< ValueType >::ZeroBlockPermutation().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertToDIA ( void  )

Convert the matrix to DIA structure.

Referenced by main(), and paralution_fortran_solve().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertToELL ( void  )

Convert the matrix to ELL structure.

Referenced by main(), and paralution_fortran_solve().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertToHYB ( void  )

Convert the matrix to HYB structure.

Referenced by main(), and paralution_fortran_solve().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ConvertToMCSR ( void  )

Convert the matrix to MCSR structure.

Referenced by main(), and paralution_fortran_solve().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::CopyFrom ( const LocalMatrix< ValueType > &  src)

Copy matrix (values and structure) from another LocalMatrix.

Referenced by paralution::LocalMatrix< ValueType >::AMGAggregate(), paralution::LocalMatrix< ValueType >::AMGAggregation(), paralution::LocalMatrix< ValueType >::AMGConnect(), paralution::LocalMatrix< ValueType >::AMGSmoothedAggregation(), paralution::LocalMatrix< ValueType >::Apply(), paralution::LocalMatrix< ValueType >::ApplyAdd(), paralution::LocalMatrix< ValueType >::Check(), paralution::LocalMatrix< ValueType >::CMK(), paralution::LocalMatrix< ValueType >::ExtractColumnVector(), paralution::LocalMatrix< ValueType >::ExtractDiagonal(), paralution::LocalMatrix< ValueType >::ExtractInverseDiagonal(), paralution::LocalMatrix< ValueType >::ExtractL(), paralution::LocalMatrix< ValueType >::ExtractRowVector(), paralution::LocalMatrix< ValueType >::ExtractSubMatrix(), paralution::LocalMatrix< ValueType >::ExtractU(), paralution::LocalMatrix< ValueType >::FSAI(), paralution::LocalMatrix< ValueType >::Gershgorin(), paralution::LocalMatrix< ValueType >::Householder(), paralution::LocalMatrix< ValueType >::LLSolve(), paralution::LocalMatrix< ValueType >::LSolve(), paralution::LocalMatrix< ValueType >::LUSolve(), paralution::LocalMatrix< ValueType >::MatrixAdd(), paralution::LocalMatrix< ValueType >::MatrixMult(), paralution::LocalMatrix< ValueType >::MaximalIndependentSet(), paralution::LocalMatrix< ValueType >::MultiColoring(), paralution::LocalMatrix< ValueType >::QRSolve(), paralution::LocalMatrix< ValueType >::RCMK(), paralution::LocalMatrix< ValueType >::USolve(), paralution::LocalMatrix< ValueType >::WriteFileCSR(), paralution::LocalMatrix< ValueType >::WriteFileMTX(), and paralution::LocalMatrix< ValueType >::ZeroBlockPermutation().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::CopyFromAsync ( const LocalMatrix< ValueType > &  src)

Async copy matrix (values and structure) from another LocalMatrix.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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 allocated (call AllocateCOO first)

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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 be allocated (call AllocateCSR first)

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::CopyToCOO ( int row,
int col,
ValueType *  val 
) const

Copy (export) COO matrix described in three arrays (rows, columns, values). The output arrays have to be allocated.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::CopyToCSR ( int row_offsets,
int col,
ValueType *  val 
) const

Copy (export) CSR matrix described in three arrays (offsets, columns, values). The output arrays have to be allocated.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::CreateFromMap ( const LocalVector< int > &  map,
const int  n,
const int  m 
)

Create a restriction matrix operator based on an int vector map.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::DiagonalMatrixMult ( const LocalVector< ValueType > &  diag)

Multiply the matrix with diagonal matrix (stored in LocalVector), as DiagonalMatrixMultR()

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::DiagonalMatrixMultL ( const LocalVector< ValueType > &  diag)

Multiply the matrix with diagonal matrix (stored in LocalVector), this=diag*this.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::DiagonalMatrixMultR ( const LocalVector< ValueType > &  diag)

Multiply the matrix with diagonal matrix (stored in LocalVector), this=this*diag.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ExtractColumnVector ( const int  idx,
LocalVector< ValueType > *  vec 
) const

Extract values from a column of a matrix to a vector.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ExtractDiagonal ( LocalVector< ValueType > *  vec_diag) const

Extract the diagonal values of the matrix into a LocalVector.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ExtractInverseDiagonal ( LocalVector< ValueType > *  vec_inv_diag) const

Extract the inverse (reciprocal) diagonal values of the matrix into a LocalVector.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ExtractL ( LocalMatrix< ValueType > *  L,
const bool  diag 
) const

Extract the lower triangular matrix.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ExtractRowVector ( const int  idx,
LocalVector< ValueType > *  vec 
) const

Extract values from a row of a matrix to a vector.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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; row/col_offset have sizes col/row_num_blocks+1, where [i+1]-[i] defines the i-th size of the sub-matrix)

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ExtractU ( LocalMatrix< ValueType > *  U,
const bool  diag 
) const

Extract the upper triangular matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::free_assembly_data ( void  )
private
template<typename ValueType>
void paralution::LocalMatrix< ValueType >::FSAI ( const int  power,
const LocalMatrix< ValueType > *  pattern 
)

Factorized Sparse Approximate Inverse assembly for given system matrix power pattern or external sparsity pattern.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::Gershgorin ( ValueType &  lambda_min,
ValueType &  lambda_max 
) const

Compute the spectrum approximation with Gershgorin circles theorem.

Referenced by main().

template<typename ValueType >
unsigned int paralution::LocalMatrix< ValueType >::get_format ( void  ) const
template<typename ValueType >
int paralution::Operator< ValueType >::get_ghost_ncol ( void  ) const
virtualinherited

Return the number of columns in the ghost matrix/stencil.

template<typename ValueType >
int paralution::Operator< ValueType >::get_ghost_nnz ( void  ) const
virtualinherited

Return the number of non-zeros in the ghost matrix/stencil.

template<typename ValueType >
int paralution::Operator< ValueType >::get_ghost_nrow ( void  ) const
virtualinherited

Return the number of rows in the ghost matrix/stencil.

template<typename ValueType >
int paralution::Operator< ValueType >::get_local_ncol ( void  ) const
virtualinherited

Return the number of columns in the local matrix/stencil.

template<typename ValueType >
int paralution::Operator< ValueType >::get_local_nnz ( void  ) const
virtualinherited

Return the number of non-zeros in the local matrix/stencil.

template<typename ValueType >
int paralution::Operator< ValueType >::get_local_nrow ( void  ) const
virtualinherited

Return the number of rows in the local matrix/stencil.

template<typename ValueType >
int paralution::LocalMatrix< ValueType >::get_ncol ( void  ) const
virtual

Return the number of columns in the matrix/stencil.

Implements paralution::Operator< ValueType >.

Referenced by main(), and paralution::LocalMatrix< ValueType >::MatrixMult().

template<typename ValueType >
int paralution::LocalMatrix< ValueType >::get_nnz ( void  ) const
virtual

Return the number of non-zeros in the matrix/stencil.

Implements paralution::Operator< ValueType >.

Referenced by main(), and paralution::LocalMatrix< ValueType >::MatrixMult().

template<typename ValueType >
int paralution::LocalMatrix< ValueType >::get_nrow ( void  ) const
virtual

Return the number of rows in the matrix/stencil.

Implements paralution::Operator< ValueType >.

Referenced by main(), and paralution::LocalMatrix< ValueType >::MatrixMult().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::Householder ( const int  idx,
ValueType &  beta,
LocalVector< ValueType > *  vec 
) const

Compute Householder vector.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ICFactorize ( LocalVector< ValueType > *  inv_diag)

Perform IC(0) factorization.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ILU0Factorize ( void  )

Perform ILU(0) factorization.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ILUpFactorize ( const int  p,
const bool  level = true 
)

Perform ILU(p) factorization based on power (see power(q)-pattern method, D. Lukarski "Parallel Sparse Linear Algebra for Multi-core and Many-core Platforms - Parallel Solvers and Preconditioners", PhD Thesis, 2012, KIT); level==true build the structure based on levels; level==false build the structure only based on the power(p+1)

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ILUTFactorize ( const double  t,
const int  maxrow 
)

Perform ILU(t,m) factorization based on threshold and maximum number of elements per row.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Invert ( void  )

Matrix inversion using QR decomposition.

template<typename ValueType >
bool paralution::LocalMatrix< ValueType >::is_accel ( void  ) const
protectedvirtual
template<typename ValueType >
bool paralution::LocalMatrix< ValueType >::is_host ( void  ) const
protectedvirtual

Return true if the object is on the host.

Implements paralution::BaseParalution< ValueType >.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::LAnalyse ( const bool  diag_unit = false)

Analyse the structure (level-scheduling) L-part diag_unit == true the diag is 1; diag_unit == false the diag is 0;.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::LAnalyseClear ( void  )

Delete the analysed data (see LAnalyse) L-party.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::LeaveDataPtrCOO ( int **  row,
int **  col,
ValueType **  val 
)

Leave a COO matrix to Host pointers.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::LeaveDataPtrCSR ( int **  row_offset,
int **  col,
ValueType **  val 
)

Leave a CSR matrix to Host pointers.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::LeaveDataPtrDENSE ( ValueType **  val)

Leave a DENSE matrix to Host pointers.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::LeaveDataPtrDIA ( int **  offset,
ValueType **  val,
int num_diag 
)

Leave a DIA matrix to Host pointers.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::LeaveDataPtrELL ( int **  col,
ValueType **  val,
int max_row 
)

Leave an ELL matrix to Host pointers.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::LLAnalyse ( void  )

Analyse the structure (level-scheduling)

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::LLAnalyseClear ( void  )

Delete the analysed data (see LLAnalyse)

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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 in parallel.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::LLSolve ( const LocalVector< ValueType > &  in,
const LocalVector< ValueType > &  inv_diag,
LocalVector< ValueType > *  out 
) const
template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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 parallel.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::LUAnalyse ( void  )

Analyse the structure (level-scheduling)

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::LUAnalyseClear ( void  )

Delete the analysed data (see LUAnalyse)

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::LUFactorize ( void  )

Perform LU factorization.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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 in parallel.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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 the matrix is not changed; if structure==true a new sparsity pattern is computed.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::MatrixMult ( const LocalMatrix< ValueType > &  A,
const LocalMatrix< ValueType > &  B 
)

Multiply two matrices, this = A * B.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::MaximalIndependentSet ( int size,
LocalVector< int > *  permutation 
) const

Perform maximal independent set decomposition of the matrix; Returns the size of the maximal independent set and the corresponding permutation.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::MoveToAcceleratorAsync ( void  )
virtual

Move the object to the Accelerator backend with async move.

Reimplemented from paralution::BaseParalution< ValueType >.

Referenced by main().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::MoveToHostAsync ( void  )
virtual

Move the object to the Host backend with async move.

Reimplemented from paralution::BaseParalution< ValueType >.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::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 (the array is allocated in the function) and the permutation.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Permute ( const LocalVector< int > &  permutation)

Perform (forward) permutation of the matrix.

Referenced by main().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::PermuteBackward ( const LocalVector< int > &  permutation)

Perform (backward) permutation of the matrix.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::QRSolve ( const LocalVector< ValueType > &  in,
LocalVector< ValueType > *  out 
) const
template<typename ValueType >
void paralution::LocalMatrix< ValueType >::RCMK ( LocalVector< int > *  permutation) const

Create permutation vector for reverse CMK reordering of the matrix.

Referenced by main().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ReadFileCSR ( const std::string  filename)

Read matrix from CSR (PARALUTION binary format) file.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ReadFileMTX ( const std::string  filename)

Read matrix from MTX (Matrix Market Format) file.

Referenced by main().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ReplaceColumnVector ( const int  idx,
const LocalVector< ValueType > &  vec 
)

Replace a column vector of a matrix.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ReplaceRowVector ( const int  idx,
const LocalVector< ValueType > &  vec 
)

Replace a row vector of a matrix.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::Scale ( const ValueType  alpha)

Scale all the values in the matrix.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ScaleDiagonal ( const ValueType  alpha)

Scale the diagonal entries of the matrix with alpha, the diagonal elements must exist.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::ScaleOffDiagonal ( const ValueType  alpha)

Scale the off-diagonal entries of the matrix with alpha, the diagonal elements must exist.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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.

Referenced by import_dealii_matrix(), and paralution_fortran_solve_coo().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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.

Referenced by mexFunction(), paralution_fortran_solve_csr(), and paralution_solve().

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::SetDataPtrDENSE ( ValueType **  val,
std::string  name,
const int  nrow,
const int  ncol 
)
template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Sort ( void  )

Sort the matrix indices.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::SPAI ( void  )

SParse Approximate Inverse assembly for given system matrix pattern.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::SymbolicPower ( const int  p)

Perform symbolic computation (structure only) of |this|^p.

Referenced by paralution::LocalMatrix< ValueType >::ILUpFactorize().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Sync ( void  )
virtual

Reimplemented from paralution::BaseParalution< ValueType >.

Referenced by main().

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Transpose ( void  )

Transpose the matrix.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::UAnalyse ( const bool  diag_unit = false)

Analyse the structure (level-scheduling) U-part; diag_unit == true the diag is 1; diag_unit == false the diag is 0;.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::UAnalyseClear ( void  )

Delete the analysed data (see UAnalyse) U-party.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::UpdateValuesCSR ( ValueType *  val)

Update CSR matrix entries only, structure will remain the same.

template<typename ValueType>
void paralution::LocalMatrix< ValueType >::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 parallel.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::WriteFileCSR ( const std::string  filename) const

Write matrix to CSR (PARALUTION binary format) file.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::WriteFileMTX ( const std::string  filename) const

Write matrix to MTX (Matrix Market Format) file.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::ZeroBlockPermutation ( int size,
LocalVector< int > *  permutation 
) const

Return a permutation for saddle-point problems (zero diagonal entries), where all zero diagonal elements are mapped to the last block; the return size is the size of the first block.

template<typename ValueType >
void paralution::LocalMatrix< ValueType >::Zeros ( void  )

Set all the values to zero.

Referenced by main().

Friends And Related Function Documentation

template<typename ValueType>
friend class LocalVector< ValueType >
friend

Field Documentation

template<typename ValueType>
int* paralution::LocalMatrix< ValueType >::assembly_irank
private
template<typename ValueType>
int* paralution::LocalMatrix< ValueType >::assembly_loop_end
private
template<typename ValueType>
int* paralution::LocalMatrix< ValueType >::assembly_loop_start
private
template<typename ValueType>
int* paralution::LocalMatrix< ValueType >::assembly_rank
private
template<typename ValueType>
int paralution::LocalMatrix< ValueType >::assembly_threads
private
template<typename ValueType>
bool paralution::BaseParalution< ValueType >::asyncf
protectedinherited
size_t paralution::ParalutionObj::global_obj_id
protectedinherited
template<typename ValueType>
Paralution_Backend_Descriptor paralution::BaseParalution< ValueType >::local_backend_
protectedinherited
template<typename ValueType>
BaseMatrix<ValueType>* paralution::LocalMatrix< ValueType >::matrix_
private

Pointer from the base matrix class to the current allocated matrix (host_ or accel_)

Referenced by paralution::LocalMatrix< ValueType >::AMGAggregate(), paralution::LocalMatrix< ValueType >::AMGAggregation(), paralution::LocalMatrix< ValueType >::AMGConnect(), paralution::LocalMatrix< ValueType >::AMGSmoothedAggregation(), paralution::LocalMatrix< ValueType >::Apply(), paralution::LocalMatrix< ValueType >::ApplyAdd(), paralution::LocalMatrix< ValueType >::Check(), paralution::LocalMatrix< ValueType >::CloneFrom(), paralution::LocalMatrix< ValueType >::CMK(), paralution::LocalMatrix< ValueType >::CopyFrom(), paralution::LocalMatrix< ValueType >::CopyFromAsync(), paralution::LocalMatrix< ValueType >::ExtractColumnVector(), paralution::LocalMatrix< ValueType >::ExtractDiagonal(), paralution::LocalMatrix< ValueType >::ExtractInverseDiagonal(), paralution::LocalMatrix< ValueType >::ExtractL(), paralution::LocalMatrix< ValueType >::ExtractRowVector(), paralution::LocalMatrix< ValueType >::ExtractSubMatrix(), paralution::LocalMatrix< ValueType >::ExtractU(), paralution::LocalMatrix< ValueType >::FSAI(), paralution::LocalMatrix< ValueType >::Gershgorin(), paralution::LocalMatrix< ValueType >::Householder(), paralution::LocalMatrix< ValueType >::ILUpFactorize(), paralution::LocalMatrix< ValueType >::LLSolve(), paralution::LocalMatrix< ValueType >::LSolve(), paralution::LocalMatrix< ValueType >::LUSolve(), paralution::LocalMatrix< ValueType >::MatrixAdd(), paralution::LocalMatrix< ValueType >::MatrixMult(), paralution::LocalMatrix< ValueType >::MaximalIndependentSet(), paralution::LocalMatrix< ValueType >::MultiColoring(), paralution::LocalMatrix< ValueType >::QRSolve(), paralution::LocalMatrix< ValueType >::RCMK(), paralution::LocalMatrix< ValueType >::USolve(), paralution::LocalMatrix< ValueType >::WriteFileCSR(), paralution::LocalMatrix< ValueType >::WriteFileMTX(), and paralution::LocalMatrix< ValueType >::ZeroBlockPermutation().


The documentation for this class was generated from the following files: