Opm::gpuistl::GpuSparseMatrix< T > Class Template Reference

The GpuSparseMatrix class simple wrapper class for a CuSparse matrix. More...

#include <GpuSparseMatrix.hpp>

Inheritance diagram for Opm::gpuistl::GpuSparseMatrix< T >:
Inheritance graph

Public Types

using field_type = T
 

Public Member Functions

 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 &)
 
GpuSparseMatrixoperator= (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::GpuSparseMatrixDescriptiongetDescription ()
 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 FunctionType >
auto dispatchOnBlocksize (FunctionType function) const
 Dispatches a function based on the block size of the matrix. More...
 

Static Public Member Functions

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...
 

Static Public Attributes

static constexpr int max_block_size = 6
 Maximum block size supported by this implementation. More...
 

Detailed Description

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
Tthe 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.
This class uses the legacy cuSPARSE API, to be compatible with CuSparse's ilu0 preconditioner. However, this preconditioner is deprecated and will be removed in future versions of CuSparse. So we should migrate to the new cuSPARSE generic API in the future.
To also support block size 1, we use the GpuSparseMatrixGeneric class which uses the new cuSPARSE generic API. This is a temporary solution, and we should migrate to the new API for all block sizes in the future by replacing this class with GpuSparseMatrixGeneric.

Member Typedef Documentation

◆ field_type

template<typename T >
using Opm::gpuistl::GpuSparseMatrix< T >::field_type = T

Constructor & Destructor Documentation

◆ 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]nonZeroElementsthe non-zero values of the matrix
[in]rowIndicesthe row indices of the non-zero elements
[in]columnIndicesthe column indices of the non-zero elements
[in]numberOfNonzeroBlocksnumber of nonzero elements
[in]blockSizesize of each block matrix (typically 3)
[in]numberOfRowsthe 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]rowIndicesthe row indices of the non-zero elements
[in]columnIndicesthe column indices of the non-zero elements
[in]blockSizesize 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]

template<typename T >
Opm::gpuistl::GpuSparseMatrix< T >::GpuSparseMatrix ( const GpuSparseMatrix< T > &  )

◆ ~GpuSparseMatrix()

template<typename T >
virtual Opm::gpuistl::GpuSparseMatrix< T >::~GpuSparseMatrix ( )
virtual

Member Function Documentation

◆ blockSize()

template<typename T >
size_t Opm::gpuistl::GpuSparseMatrix< T >::blockSize ( ) const
inline

◆ dim()

template<typename T >
size_t Opm::gpuistl::GpuSparseMatrix< T >::dim ( ) const
inline

dim returns the dimension of the vector space on which this matrix acts

This is equivalent to matrix.N() * matrix.blockSize()

Returns
matrix.N() * matrix.blockSize()

References Opm::gpuistl::detail::to_size_t().

Referenced by Opm::gpuistl::GpuSeqILU0< M, X, Y, l >::GpuSeqILU0().

◆ dispatchOnBlocksize()

template<typename T >
template<class FunctionType >
auto Opm::gpuistl::GpuSparseMatrix< T >::dispatchOnBlocksize ( FunctionType  function) const
inline

Dispatches a function based on the block size of the matrix.

This method allows executing different code paths depending on the block size of the matrix, up to the maximum block size specified by max_block_size.

Use this function if you need the block size to be known at compile time.

Template Parameters
FunctionTypeType of the function to be dispatched
Parameters
functionThe function to be executed based on the block size
Returns
The result of the function execution

You can use this function as

matrix.dispatchOnBlocksize([](auto val) {
constexpr int blockSize = decltype(val)::value;
});
size_t blockSize() const
blockSize size of the blocks
Definition: GpuSparseMatrix.hpp:273

◆ fromMatrix()

template<typename T >
template<class MatrixType >
static GpuSparseMatrix< T > Opm::gpuistl::GpuSparseMatrix< T >::fromMatrix ( const MatrixType &  matrix,
bool  copyNonZeroElementsDirectly = false 
)
static

fromMatrix creates a new matrix with the same block size and values as the given matrix

Parameters
matrixthe matrix to copy from
copyNonZeroElementsDirectlyif 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
MatrixTypeis assumed to be a Dune::BCRSMatrix compatible matrix.

◆ getColumnIndices() [1/2]

template<typename T >
GpuVector< int > & Opm::gpuistl::GpuSparseMatrix< T >::getColumnIndices ( )
inline

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.

Referenced by Opm::gpuistl::AmgxInterface::updateAmgxMatrixFromGpuSparseMatrix().

◆ getColumnIndices() [2/2]

template<typename T >
const GpuVector< int > & Opm::gpuistl::GpuSparseMatrix< T >::getColumnIndices ( ) const
inline

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()

template<typename T >
detail::GpuSparseMatrixDescription & Opm::gpuistl::GpuSparseMatrix< T >::getDescription ( )
inline

getDescription the cusparse matrix description.

This description is needed for most calls to the CuSparse library

◆ getNonZeroValues() [1/2]

template<typename T >
GpuVector< T > & Opm::gpuistl::GpuSparseMatrix< T >::getNonZeroValues ( )
inline

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.

Referenced by Opm::gpuistl::AmgxInterface::updateAmgxMatrixCoefficientsFromGpuSparseMatrix(), Opm::gpuistl::AmgxInterface::updateAmgxMatrixFromGpuSparseMatrix(), and Opm::gpuistl::AmgxInterface::updateGpuSparseMatrixFromAmgxMatrix().

◆ getNonZeroValues() [2/2]

template<typename T >
const GpuVector< T > & Opm::gpuistl::GpuSparseMatrix< T >::getNonZeroValues ( ) const
inline

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]

template<typename T >
GpuVector< int > & Opm::gpuistl::GpuSparseMatrix< T >::getRowIndices ( )
inline

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.

Referenced by Opm::gpuistl::AmgxInterface::updateAmgxMatrixFromGpuSparseMatrix().

◆ getRowIndices() [2/2]

template<typename T >
const GpuVector< int > & Opm::gpuistl::GpuSparseMatrix< T >::getRowIndices ( ) const
inline

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()

template<typename T >
virtual void Opm::gpuistl::GpuSparseMatrix< T >::mv ( const GpuVector< T > &  x,
GpuVector< T > &  y 
) const
virtual

mv performs matrix vector multiply y = Ax

Parameters
[in]xthe vector to multiply the matrix with
[out]ythe output vector

◆ N()

◆ nonzeroes()

template<typename T >
size_t Opm::gpuistl::GpuSparseMatrix< T >::nonzeroes ( ) const
inline

nonzeroes behaves as the Dune::BCRSMatrix::nonzeros() function and returns the number of non zero blocks

Returns
number of non zero blocks.

References Opm::gpuistl::detail::to_size_t().

Referenced by Opm::gpuistl::GpuSeqILU0< M, X, Y, l >::GpuSeqILU0(), Opm::gpuistl::AmgxInterface::updateAmgxMatrixCoefficientsFromGpuSparseMatrix(), and Opm::gpuistl::AmgxInterface::updateAmgxMatrixFromGpuSparseMatrix().

◆ operator=()

template<typename T >
GpuSparseMatrix & Opm::gpuistl::GpuSparseMatrix< T >::operator= ( const GpuSparseMatrix< T > &  )
delete

◆ setLowerTriangular()

template<typename T >
void Opm::gpuistl::GpuSparseMatrix< T >::setLowerTriangular ( )

setLowerTriangular sets the CuSparse flag that this is an lower diagonal (with non-unit diagonal) matrix.

◆ setNonUnitDiagonal()

template<typename T >
void Opm::gpuistl::GpuSparseMatrix< T >::setNonUnitDiagonal ( )

setNonUnitDiagonal sets the CuSparse flag that this has non-unit diagional.

◆ setUnitDiagonal()

template<typename T >
void Opm::gpuistl::GpuSparseMatrix< T >::setUnitDiagonal ( )

setUnitDiagonal sets the CuSparse flag that this has unit diagional.

◆ setUpperTriangular()

template<typename T >
void Opm::gpuistl::GpuSparseMatrix< T >::setUpperTriangular ( )

setUpperTriangular sets the CuSparse flag that this is an upper diagonal (with unit diagonal) matrix.

◆ umv()

template<typename T >
virtual void Opm::gpuistl::GpuSparseMatrix< T >::umv ( const GpuVector< T > &  x,
GpuVector< T > &  y 
) const
virtual

umv computes y=Ax+y

Parameters
[in]xthe vector to multiply with A
[in,out]ythe vector to add and store the output in

◆ updateNonzeroValues() [1/2]

template<typename T >
void Opm::gpuistl::GpuSparseMatrix< T >::updateNonzeroValues ( const GpuSparseMatrix< T > &  matrix)

updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix

Parameters
matrixthe 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 >
void Opm::gpuistl::GpuSparseMatrix< T >::updateNonzeroValues ( const MatrixType &  matrix,
bool  copyNonZeroElementsDirectly = false 
)

updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix

Parameters
matrixthe matrix to extract the non-zero values from
copyNonZeroElementsDirectlyif 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 matrix has been initialized. If unsure, leave it as false.
Note
This assumes the given matrix has the same sparsity pattern.
Template Parameters
MatrixTypeis assumed to be a Dune::BCRSMatrix compatible matrix.

Referenced by Opm::gpuistl::SolverAdapter< Operator, UnderlyingSolver, X >::apply(), and Opm::gpuistl::PreconditionerCPUMatrixToGPUMatrix< X, Y, CudaPreconditionerType, CPUMatrixType >::update().

◆ usmv()

template<typename T >
virtual void Opm::gpuistl::GpuSparseMatrix< T >::usmv ( alpha,
const GpuVector< T > &  x,
GpuVector< T > &  y 
) const
virtual

umv computes y=alpha * Ax + y

Parameters
[in]alphaThe scaling factor for the matrix-vector product
[in]xthe vector to multiply with A
[in,out]ythe vector to add and store the output in

Member Data Documentation

◆ max_block_size

template<typename T >
constexpr int Opm::gpuistl::GpuSparseMatrix< T >::max_block_size = 6
staticconstexpr

Maximum block size supported by this implementation.

This constant defines an upper bound on the block size to ensure reasonable compilation times. While this class itself could support larger values, functions that call dispatchOnBlocksize() might have limitations. This value can be increased if needed, but will increase compilation time due to template instantiations.


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