19#ifndef OPM_GPUSPARSEMATRIX_HPP
20#define OPM_GPUSPARSEMATRIX_HPP
24#include <opm/common/ErrorMacros.hpp>
25#include <opm/simulators/linalg/gpuistl/GpuVector.hpp>
26#include <opm/simulators/linalg/gpuistl/detail/CuMatrixDescription.hpp>
27#include <opm/simulators/linalg/gpuistl/detail/CuSparseHandle.hpp>
28#include <opm/simulators/linalg/gpuistl/detail/safe_conversion.hpp>
64 const int* rowIndices,
105 template <
class MatrixType>
163 return m_nonZeroElements;
173 return m_nonZeroElements;
203 return m_columnIndices;
213 return m_columnIndices;
252 return *m_matrixDescription;
293 template <
class MatrixType>
305 const int m_numberOfNonzeroBlocks;
306 const int m_numberOfRows;
307 const int m_blockSize;
312 template <
class VectorType>
313 void assertSameSize(
const VectorType&
vector)
const;
The GpuSparseMatrix class simple wrapper class for a CuSparse matrix.
Definition GpuSparseMatrix.hpp:48
GpuSparseMatrix & operator=(const GpuSparseMatrix &)=delete
We don't want to be able to copy this for now (too much hassle in copying the cusparse resources)
const GpuVector< int > & getRowIndices() const
getRowIndices returns the row indices used to represent the BSR structure.
Definition GpuSparseMatrix.hpp:191
const GpuVector< int > & getColumnIndices() const
getColumnIndices returns the column indices used to represent the BSR structure.
Definition GpuSparseMatrix.hpp:211
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
Definition GpuSparseMatrix.cpp:108
void setNonUnitDiagonal()
setNonUnitDiagonal sets the CuSparse flag that this has non-unit diagional.
Definition GpuSparseMatrix.cpp:209
void updateNonzeroValues(const MatrixType &matrix, bool copyNonZeroElementsDirectly=false)
updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix
Definition GpuSparseMatrix.cpp:170
virtual void umv(const GpuVector< T > &x, GpuVector< T > &y) const
umv computes y=Ax+y
Definition GpuSparseMatrix.cpp:250
detail::GpuSparseMatrixDescription & getDescription()
getDescription the cusparse matrix description.
Definition GpuSparseMatrix.hpp:250
GpuVector< int > & getRowIndices()
getRowIndices returns the row indices used to represent the BSR structure.
Definition GpuSparseMatrix.hpp:181
const GpuVector< T > & getNonZeroValues() const
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
Definition GpuSparseMatrix.hpp:171
size_t nonzeroes() const
nonzeroes behaves as the Dune::BCRSMatrix::nonzeros() function and returns the number of non zero blo...
Definition GpuSparseMatrix.hpp:146
size_t blockSize() const
blockSize size of the blocks
Definition GpuSparseMatrix.hpp:235
void setUpperTriangular()
setUpperTriangular sets the CuSparse flag that this is an upper diagonal (with unit diagonal) matrix.
Definition GpuSparseMatrix.cpp:188
size_t N() const
N returns the number of rows (which is equal to the number of columns)
Definition GpuSparseMatrix.hpp:132
void setLowerTriangular()
setLowerTriangular sets the CuSparse flag that this is an lower diagonal (with non-unit diagonal) mat...
Definition GpuSparseMatrix.cpp:195
virtual void mv(const GpuVector< T > &x, GpuVector< T > &y) const
mv performs matrix vector multiply y = Ax
Definition GpuSparseMatrix.cpp:216
virtual void usmv(T alpha, const GpuVector< T > &x, GpuVector< T > &y) const
umv computes y=alpha * Ax + y
Definition GpuSparseMatrix.cpp:285
void setUnitDiagonal()
setUnitDiagonal sets the CuSparse flag that this has unit diagional.
Definition GpuSparseMatrix.cpp:202
GpuVector< int > & getColumnIndices()
getColumnIndices returns the column indices used to represent the BSR structure.
Definition GpuSparseMatrix.hpp:201
GpuSparseMatrix(const GpuSparseMatrix &)=delete
We don't want to be able to copy this for now (too much hassle in copying the cusparse resources)
GpuVector< T > & getNonZeroValues()
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
Definition GpuSparseMatrix.hpp:161
size_t dim() const
dim returns the dimension of the vector space on which this matrix acts
Definition GpuSparseMatrix.hpp:222
The CuSparseHandle class provides a singleton for the simulator universal cuSparseHandle.
Definition CuSparseHandle.hpp:41
The CuSparseResource class wraps a CuSparse resource in a proper RAII pattern.
Definition CuSparseResource.hpp:55
std::shared_ptr< CuSparseResource< cusparseMatDescr_t > > GpuSparseMatrixDescriptionPtr
Pointer to GpuSparseMatrixDescription holder.
Definition CuMatrixDescription.hpp:35
__host__ __device__ std::size_t to_size_t(int i)
to_size_t converts a (on most relevant platforms) a 32 bit signed int to a 64 bits unsigned int
Definition safe_conversion.hpp:86
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242