The GpuSparseMatrix class simple wrapper class for a CuSparse matrix.
More...
#include <GpuSparseMatrix.hpp>
|
| GpuSparseMatrix (const T *nonZeroElements, const int *rowIndices, const int *columnIndices, size_t numberOfNonzeroBlocks, size_t blockSize, size_t numberOfRows) |
|
| GpuSparseMatrix (const GpuVector< int > &rowIndices, const GpuVector< int > &columnIndices, size_t blockSize) |
|
| GpuSparseMatrix (const GpuSparseMatrix &) |
|
GpuSparseMatrix & | operator= (const GpuSparseMatrix &)=delete |
|
virtual | ~GpuSparseMatrix () |
|
void | setUpperTriangular () |
| setUpperTriangular sets the CuSparse flag that this is an upper diagonal (with unit diagonal) matrix. More...
|
|
void | setLowerTriangular () |
| setLowerTriangular sets the CuSparse flag that this is an lower diagonal (with non-unit diagonal) matrix. More...
|
|
void | setUnitDiagonal () |
| setUnitDiagonal sets the CuSparse flag that this has unit diagional. More...
|
|
void | setNonUnitDiagonal () |
| setNonUnitDiagonal sets the CuSparse flag that this has non-unit diagional. More...
|
|
size_t | N () const |
| N returns the number of rows (which is equal to the number of columns) More...
|
|
size_t | nonzeroes () const |
| nonzeroes behaves as the Dune::BCRSMatrix::nonzeros() function and returns the number of non zero blocks More...
|
|
GpuVector< T > & | getNonZeroValues () |
| getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block) More...
|
|
const GpuVector< T > & | getNonZeroValues () const |
| getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block) More...
|
|
GpuVector< int > & | getRowIndices () |
| getRowIndices returns the row indices used to represent the BSR structure. More...
|
|
const GpuVector< int > & | getRowIndices () const |
| getRowIndices returns the row indices used to represent the BSR structure. More...
|
|
GpuVector< int > & | getColumnIndices () |
| getColumnIndices returns the column indices used to represent the BSR structure. More...
|
|
const GpuVector< int > & | getColumnIndices () const |
| getColumnIndices returns the column indices used to represent the BSR structure. More...
|
|
size_t | dim () const |
| dim returns the dimension of the vector space on which this matrix acts More...
|
|
size_t | blockSize () const |
| blockSize size of the blocks More...
|
|
detail::GpuSparseMatrixDescription & | getDescription () |
| getDescription the cusparse matrix description. More...
|
|
virtual void | mv (const GpuVector< T > &x, GpuVector< T > &y) const |
| mv performs matrix vector multiply y = Ax More...
|
|
virtual void | umv (const GpuVector< T > &x, GpuVector< T > &y) const |
| umv computes y=Ax+y More...
|
|
virtual void | usmv (T alpha, const GpuVector< T > &x, GpuVector< T > &y) const |
| umv computes y=alpha * Ax + y More...
|
|
template<class MatrixType > |
void | updateNonzeroValues (const MatrixType &matrix, bool copyNonZeroElementsDirectly=false) |
| updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix More...
|
|
void | updateNonzeroValues (const GpuSparseMatrix< T > &matrix) |
| updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix More...
|
|
|
template<class MatrixType > |
static GpuSparseMatrix< T > | fromMatrix (const MatrixType &matrix, bool copyNonZeroElementsDirectly=false) |
| fromMatrix creates a new matrix with the same block size and values as the given matrix More...
|
|
template<typename T>
class Opm::gpuistl::GpuSparseMatrix< T > The GpuSparseMatrix class simple wrapper class for a CuSparse matrix.
- Note
- we currently only support simple raw primitives for T (double and float). Block size is handled through the block size parameter
- Template Parameters
-
T | the type to store. Can be either float, double or int. |
- Note
- we only support square matrices.
-
We only support Block Compressed Sparse Row Format (BSR) for now.
◆ field_type
◆ GpuSparseMatrix() [1/3]
template<typename T >
Opm::gpuistl::GpuSparseMatrix< T >::GpuSparseMatrix |
( |
const T * |
nonZeroElements, |
|
|
const int * |
rowIndices, |
|
|
const int * |
columnIndices, |
|
|
size_t |
numberOfNonzeroBlocks, |
|
|
size_t |
blockSize, |
|
|
size_t |
numberOfRows |
|
) |
| |
Create the sparse matrix specified by the raw data.
- Note
- Prefer to use the constructor taking a const reference to a matrix instead.
- Parameters
-
[in] | nonZeroElements | the non-zero values of the matrix |
[in] | rowIndices | the row indices of the non-zero elements |
[in] | columnIndices | the column indices of the non-zero elements |
[in] | numberOfNonzeroElements | number of nonzero elements |
[in] | blockSize | size of each block matrix (typically 3) |
[in] | numberOfRows | the number of rows |
- Note
- We assume numberOfNonzeroBlocks, blockSize and numberOfRows all are representable as int due to restrictions in the current version of cusparse. This might change in future versions.
◆ GpuSparseMatrix() [2/3]
template<typename T >
Opm::gpuistl::GpuSparseMatrix< T >::GpuSparseMatrix |
( |
const GpuVector< int > & |
rowIndices, |
|
|
const GpuVector< int > & |
columnIndices, |
|
|
size_t |
blockSize |
|
) |
| |
Create a sparse matrix by copying the sparsity structure of another matrix, not filling in the values
- Note
- Prefer to use the constructor taking a const reference to a matrix instead.
- Parameters
-
[in] | rowIndices | the row indices of the non-zero elements |
[in] | columnIndices | the column indices of the non-zero elements |
[in] | blockSize | size of each block matrix (typically 3) |
- Note
- We assume numberOfNonzeroBlocks, blockSize and numberOfRows all are representable as int due to restrictions in the current version of cusparse. This might change in future versions.
◆ GpuSparseMatrix() [3/3]
◆ ~GpuSparseMatrix()
◆ blockSize()
◆ dim()
◆ fromMatrix()
template<typename T >
template<class MatrixType >
fromMatrix creates a new matrix with the same block size and values as the given matrix
- Parameters
-
matrix | the matrix to copy from |
copyNonZeroElementsDirectly | if true will do a memcpy from matrix[0][0][0][0], otherwise will build up the non-zero elements by looping over the matrix. Note that setting this to true will yield a performance increase, but might not always yield correct results depending on how the matrix has been initialized. If unsure, leave it as false. |
- Template Parameters
-
MatrixType | is assumed to be a Dune::BCRSMatrix compatible matrix. |
◆ getColumnIndices() [1/2]
getColumnIndices returns the column indices used to represent the BSR structure.
- Returns
- Read the CuSPARSE documentation on Block Compressed Sparse Row Format (BSR) for the exact ordering.
◆ getColumnIndices() [2/2]
getColumnIndices returns the column indices used to represent the BSR structure.
- Returns
- Read the CuSPARSE documentation on Block Compressed Sparse Row Format (BSR) for the exact ordering.
◆ getDescription()
getDescription the cusparse matrix description.
This description is needed for most calls to the CuSparse library
◆ getNonZeroValues() [1/2]
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
- Note
- Read the CuSPARSE documentation on Block Compressed Sparse Row Format (BSR) for the exact ordering.
◆ getNonZeroValues() [2/2]
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
- Note
- Read the CuSPARSE documentation on Block Compressed Sparse Row Format (BSR) for the exact ordering.
◆ getRowIndices() [1/2]
getRowIndices returns the row indices used to represent the BSR structure.
- Note
- Read the CuSPARSE documentation on Block Compressed Sparse Row Format (BSR) for the exact ordering.
◆ getRowIndices() [2/2]
getRowIndices returns the row indices used to represent the BSR structure.
- Note
- Read the CuSPARSE documentation on Block Compressed Sparse Row Format (BSR) for the exact ordering.
◆ mv()
mv performs matrix vector multiply y = Ax
- Parameters
-
[in] | x | the vector to multiply the matrix with |
[out] | y | the output vector |
- Note
- Due to limitations of CuSparse, this is only supported for block sizes greater than 1.
◆ N()
◆ nonzeroes()
◆ operator=()
◆ setLowerTriangular()
setLowerTriangular sets the CuSparse flag that this is an lower diagonal (with non-unit diagonal) matrix.
◆ setNonUnitDiagonal()
setNonUnitDiagonal sets the CuSparse flag that this has non-unit diagional.
◆ setUnitDiagonal()
setUnitDiagonal sets the CuSparse flag that this has unit diagional.
◆ setUpperTriangular()
setUpperTriangular sets the CuSparse flag that this is an upper diagonal (with unit diagonal) matrix.
◆ umv()
umv computes y=Ax+y
- Parameters
-
[in] | x | the vector to multiply with A |
[in,out] | y | the vector to add and store the output in |
- Note
- Due to limitations of CuSparse, this is only supported for block sizes greater than 1.
◆ updateNonzeroValues() [1/2]
updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix
- Parameters
-
matrix | the matrix to extract the non-zero values from |
- Note
- This assumes the given matrix has the same sparsity pattern.
◆ updateNonzeroValues() [2/2]
template<typename T >
template<class MatrixType >
◆ usmv()
umv computes y=alpha * Ax + y
- Parameters
-
[in] | x | the vector to multiply with A |
[in,out] | y | the vector to add and store the output in |
- Note
- Due to limitations of CuSparse, this is only supported for block sizes greater than 1.
The documentation for this class was generated from the following file:
|