GpuSparseMatrix.hpp
Go to the documentation of this file.
1/*
2 Copyright 2022-2023 SINTEF AS
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#ifndef OPM_GPUSPARSEMATRIX_HPP
20#define OPM_GPUSPARSEMATRIX_HPP
21#include <cusparse.h>
22#include <iostream>
23#include <memory>
24#include <opm/common/ErrorMacros.hpp>
29#include <vector>
30
31namespace Opm::gpuistl
32{
33
46template <typename T>
48{
49public:
50 using field_type = T;
51
65 GpuSparseMatrix(const T* nonZeroElements,
66 const int* rowIndices,
67 const int* columnIndices,
68 size_t numberOfNonzeroBlocks,
69 size_t blockSize,
70 size_t numberOfRows);
71
82 GpuSparseMatrix(const GpuVector<int>& rowIndices,
83 const GpuVector<int>& columnIndices,
84 size_t blockSize);
85
87
88 // We want to have this as non-mutable as possible, that is we do not want
89 // to deal with changing matrix sizes and sparsity patterns.
91
93
103 template <class MatrixType>
104 static GpuSparseMatrix<T> fromMatrix(const MatrixType& matrix, bool copyNonZeroElementsDirectly = false);
105
110
115
120
121
126
130 size_t N() const
131 {
132 // Technically this safe conversion is not needed since we enforce these to be
133 // non-negative in the constructor, but keeping them for added sanity for now.
134 //
135 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
136 // but should that be false, they can be removed.
137 return detail::to_size_t(m_numberOfRows);
138 }
139
144 size_t nonzeroes() const
145 {
146 // Technically this safe conversion is not needed since we enforce these to be
147 // non-negative in the constructor, but keeping them for added sanity for now.
148 //
149 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
150 // but should that be false, they can be removed.
151 return detail::to_size_t(m_numberOfNonzeroBlocks);
152 }
153
159 GpuVector<T>& getNonZeroValues()
160 {
161 return m_nonZeroElements;
162 }
163
169 const GpuVector<T>& getNonZeroValues() const
170 {
171 return m_nonZeroElements;
172 }
173
179 GpuVector<int>& getRowIndices()
180 {
181 return m_rowIndices;
182 }
183
189 const GpuVector<int>& getRowIndices() const
190 {
191 return m_rowIndices;
192 }
193
199 GpuVector<int>& getColumnIndices()
200 {
201 return m_columnIndices;
202 }
203
209 const GpuVector<int>& getColumnIndices() const
210 {
211 return m_columnIndices;
212 }
213
220 size_t dim() const
221 {
222 // Technically this safe conversion is not needed since we enforce these to be
223 // non-negative in the constructor, but keeping them for added sanity for now.
224 //
225 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
226 // but should that be false, they can be removed.
227 return detail::to_size_t(m_blockSize) * detail::to_size_t(m_numberOfRows);
228 }
229
233 size_t blockSize() const
234 {
235 // Technically this safe conversion is not needed since we enforce these to be
236 // non-negative in the constructor, but keeping them for added sanity for now.
237 //
238 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
239 // but should that be false, they can be removed.
240 return detail::to_size_t(m_blockSize);
241 }
242
249 {
250 return *m_matrixDescription;
251 }
252
260 virtual void mv(const GpuVector<T>& x, GpuVector<T>& y) const;
261
269 virtual void umv(const GpuVector<T>& x, GpuVector<T>& y) const;
270
271
279 virtual void usmv(T alpha, const GpuVector<T>& x, GpuVector<T>& y) const;
280
291 template <class MatrixType>
292 void updateNonzeroValues(const MatrixType& matrix, bool copyNonZeroElementsDirectly = false);
293
300
301private:
302 GpuVector<T> m_nonZeroElements;
303 GpuVector<int> m_columnIndices;
304 GpuVector<int> m_rowIndices;
305
306 // Notice that we store these three as int to make sure we are cusparse compatible.
307 //
308 // This gives the added benefit of checking the size constraints at construction of the matrix
309 // rather than in some call to cusparse.
310 const int m_numberOfNonzeroBlocks;
311 const int m_numberOfRows;
312 const int m_blockSize;
313
314 detail::GpuSparseMatrixDescriptionPtr m_matrixDescription;
315 detail::CuSparseHandle& m_cusparseHandle;
316
317 template <class VectorType>
318 void assertSameSize(const VectorType& vector) const;
319};
320} // namespace Opm::gpuistl
321#endif
The GpuSparseMatrix class simple wrapper class for a CuSparse matrix.
Definition: GpuSparseMatrix.hpp:48
GpuSparseMatrix & operator=(const GpuSparseMatrix &)=delete
virtual void mv(const GpuVector< T > &x, GpuVector< T > &y) const
mv performs matrix vector multiply y = Ax
const GpuVector< int > & getRowIndices() const
getRowIndices returns the row indices used to represent the BSR structure.
Definition: GpuSparseMatrix.hpp:189
GpuSparseMatrix(const GpuVector< int > &rowIndices, const GpuVector< int > &columnIndices, size_t blockSize)
const GpuVector< int > & getColumnIndices() const
getColumnIndices returns the column indices used to represent the BSR structure.
Definition: GpuSparseMatrix.hpp:209
void setNonUnitDiagonal()
setNonUnitDiagonal sets the CuSparse flag that this has non-unit diagional.
T field_type
Definition: GpuSparseMatrix.hpp:50
void updateNonzeroValues(const MatrixType &matrix, bool copyNonZeroElementsDirectly=false)
updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix
void updateNonzeroValues(const GpuSparseMatrix< T > &matrix)
updateNonzeroValues updates the non-zero values by using the non-zero values of the supplied matrix
detail::GpuSparseMatrixDescription & getDescription()
getDescription the cusparse matrix description.
Definition: GpuSparseMatrix.hpp:248
GpuSparseMatrix(const GpuSparseMatrix &)
GpuVector< int > & getRowIndices()
getRowIndices returns the row indices used to represent the BSR structure.
Definition: GpuSparseMatrix.hpp:179
const GpuVector< T > & getNonZeroValues() const
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
Definition: GpuSparseMatrix.hpp:169
size_t nonzeroes() const
nonzeroes behaves as the Dune::BCRSMatrix::nonzeros() function and returns the number of non zero blo...
Definition: GpuSparseMatrix.hpp:144
size_t blockSize() const
blockSize size of the blocks
Definition: GpuSparseMatrix.hpp:233
void setUpperTriangular()
setUpperTriangular sets the CuSparse flag that this is an upper diagonal (with unit diagonal) matrix.
size_t N() const
N returns the number of rows (which is equal to the number of columns)
Definition: GpuSparseMatrix.hpp:130
virtual void umv(const GpuVector< T > &x, GpuVector< T > &y) const
umv computes y=Ax+y
void setLowerTriangular()
setLowerTriangular sets the CuSparse flag that this is an lower diagonal (with non-unit diagonal) mat...
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
void setUnitDiagonal()
setUnitDiagonal sets the CuSparse flag that this has unit diagional.
virtual void usmv(T alpha, const GpuVector< T > &x, GpuVector< T > &y) const
umv computes y=alpha * Ax + y
GpuVector< int > & getColumnIndices()
getColumnIndices returns the column indices used to represent the BSR structure.
Definition: GpuSparseMatrix.hpp:199
GpuSparseMatrix(const T *nonZeroElements, const int *rowIndices, const int *columnIndices, size_t numberOfNonzeroBlocks, size_t blockSize, size_t numberOfRows)
GpuVector< T > & getNonZeroValues()
getNonZeroValues returns the GPU vector containing the non-zero values (ordered by block)
Definition: GpuSparseMatrix.hpp:159
size_t dim() const
dim returns the dimension of the vector space on which this matrix acts
Definition: GpuSparseMatrix.hpp:220
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
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
Definition: autotuner.hpp:30