19#ifndef OPM_CUSPARSEMATRIX_HPP
20#define OPM_CUSPARSEMATRIX_HPP
24#include <opm/common/ErrorMacros.hpp>
25#include <opm/simulators/linalg/cuistl/CuVector.hpp>
26#include <opm/simulators/linalg/cuistl/detail/CuMatrixDescription.hpp>
27#include <opm/simulators/linalg/cuistl/detail/CuSparseHandle.hpp>
28#include <opm/simulators/linalg/cuistl/detail/safe_conversion.hpp>
64 const int* rowIndices,
91 template <
class MatrixType>
149 return m_nonZeroElements;
159 return m_nonZeroElements;
189 return m_columnIndices;
199 return m_columnIndices;
238 return *m_matrixDescription;
279 template <
class MatrixType>
291 const int m_numberOfNonzeroBlocks;
292 const int m_numberOfRows;
293 const int m_blockSize;
298 template <
class VectorType>
299 void assertSameSize(
const VectorType&
vector)
const;
Definition AquiferInterface.hpp:35
The CuSparseMatrix class simple wrapper class for a CuSparse matrix.
Definition CuSparseMatrix.hpp:48
const CuVector< T > & getNonZeroValues() const
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
Definition CuSparseMatrix.hpp:157
size_t nonzeroes() const
nonzeroes behaves as the Dune::BCRSMatrix::nonzeros() function and returns the number of non zero blo...
Definition CuSparseMatrix.hpp:132
detail::CuSparseMatrixDescription & getDescription()
getDescription the cusparse matrix description.
Definition CuSparseMatrix.hpp:236
size_t N() const
N returns the number of rows (which is equal to the number of columns)
Definition CuSparseMatrix.hpp:118
CuVector< int > & getColumnIndices()
getColumnIndices returns the column indices used to represent the BSR structure.
Definition CuSparseMatrix.hpp:187
const CuVector< int > & getRowIndices() const
getRowIndices returns the row indices used to represent the BSR structure.
Definition CuSparseMatrix.hpp:177
void setNonUnitDiagonal()
setNonUnitDiagonal sets the CuSparse flag that this has non-unit diagional.
Definition CuSparseMatrix.cpp:178
virtual void usmv(T alpha, const CuVector< T > &x, CuVector< T > &y) const
umv computes y=alpha * Ax + y
Definition CuSparseMatrix.cpp:254
CuSparseMatrix(const CuSparseMatrix &)=delete
We don't want to be able to copy this for now (too much hassle in copying the cusparse resources)
CuSparseMatrix & operator=(const CuSparseMatrix &)=delete
We don't want to be able to copy this for now (too much hassle in copying the cusparse resources)
size_t dim() const
dim returns the dimension of the vector space on which this matrix acts
Definition CuSparseMatrix.hpp:208
void setLowerTriangular()
setLowerTriangular sets the CuSparse flag that this is an lower diagonal (with non-unit diagonal) mat...
Definition CuSparseMatrix.cpp:164
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 CuSparseMatrix.cpp:139
void setUnitDiagonal()
setUnitDiagonal sets the CuSparse flag that this has unit diagional.
Definition CuSparseMatrix.cpp:171
size_t blockSize() const
blockSize size of the blocks
Definition CuSparseMatrix.hpp:221
const CuVector< int > & getColumnIndices() const
getColumnIndices returns the column indices used to represent the BSR structure.
Definition CuSparseMatrix.hpp:197
virtual void mv(const CuVector< T > &x, CuVector< T > &y) const
mv performs matrix vector multiply y = Ax
Definition CuSparseMatrix.cpp:185
static CuSparseMatrix< 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 CuSparseMatrix.cpp:92
CuVector< T > & getNonZeroValues()
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
Definition CuSparseMatrix.hpp:147
void setUpperTriangular()
setUpperTriangular sets the CuSparse flag that this is an upper diagonal (with unit diagonal) matrix.
Definition CuSparseMatrix.cpp:157
virtual void umv(const CuVector< T > &x, CuVector< T > &y) const
umv computes y=Ax+y
Definition CuSparseMatrix.cpp:219
CuVector< int > & getRowIndices()
getRowIndices returns the row indices used to represent the BSR structure.
Definition CuSparseMatrix.hpp:167
The CuSparseHandle class provides a singleton for the simulator universal cuSparseHandle.
Definition CuSparseHandle.hpp:41
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:85
std::shared_ptr< CuSparseResource< cusparseMatDescr_t > > CuSparseMatrixDescriptionPtr
Pointer to CuSparseMatrixDescription holder.
Definition CuMatrixDescription.hpp:35