My Project
Loading...
Searching...
No Matches
GpuSparseMatrix.hpp
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>
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>
29#include <vector>
30
31namespace Opm::gpuistl
32{
33
46template <typename T>
48{
49public:
63 GpuSparseMatrix(const T* nonZeroElements,
64 const int* rowIndices,
65 const int* columnIndices,
66 size_t numberOfNonzeroBlocks,
67 size_t blockSize,
68 size_t numberOfRows);
69
80 GpuSparseMatrix(const GpuVector<int>& rowIndices,
81 const GpuVector<int>& columnIndices,
82 size_t blockSize);
83
88
93
94 virtual ~GpuSparseMatrix();
95
105 template <class MatrixType>
106 static GpuSparseMatrix<T> fromMatrix(const MatrixType& matrix, bool copyNonZeroElementsDirectly = false);
107
111 void setUpperTriangular();
112
116 void setLowerTriangular();
117
121 void setUnitDiagonal();
122
123
127 void setNonUnitDiagonal();
128
132 size_t N() const
133 {
134 // Technically this safe conversion is not needed since we enforce these to be
135 // non-negative in the constructor, but keeping them for added sanity for now.
136 //
137 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
138 // but should that be false, they can be removed.
139 return detail::to_size_t(m_numberOfRows);
140 }
141
146 size_t nonzeroes() const
147 {
148 // Technically this safe conversion is not needed since we enforce these to be
149 // non-negative in the constructor, but keeping them for added sanity for now.
150 //
151 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
152 // but should that be false, they can be removed.
153 return detail::to_size_t(m_numberOfNonzeroBlocks);
154 }
155
161 GpuVector<T>& getNonZeroValues()
162 {
163 return m_nonZeroElements;
164 }
165
171 const GpuVector<T>& getNonZeroValues() const
172 {
173 return m_nonZeroElements;
174 }
175
181 GpuVector<int>& getRowIndices()
182 {
183 return m_rowIndices;
184 }
185
191 const GpuVector<int>& getRowIndices() const
192 {
193 return m_rowIndices;
194 }
195
201 GpuVector<int>& getColumnIndices()
202 {
203 return m_columnIndices;
204 }
205
211 const GpuVector<int>& getColumnIndices() const
212 {
213 return m_columnIndices;
214 }
215
222 size_t dim() const
223 {
224 // Technically this safe conversion is not needed since we enforce these to be
225 // non-negative in the constructor, but keeping them for added sanity for now.
226 //
227 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
228 // but should that be false, they can be removed.
229 return detail::to_size_t(m_blockSize) * detail::to_size_t(m_numberOfRows);
230 }
231
235 size_t blockSize() const
236 {
237 // Technically this safe conversion is not needed since we enforce these to be
238 // non-negative in the constructor, but keeping them for added sanity for now.
239 //
240 // We don't believe this will yield any performance penality (it's used too far away from the inner loop),
241 // but should that be false, they can be removed.
242 return detail::to_size_t(m_blockSize);
243 }
244
251 {
252 return *m_matrixDescription;
253 }
254
262 virtual void mv(const GpuVector<T>& x, GpuVector<T>& y) const;
263
271 virtual void umv(const GpuVector<T>& x, GpuVector<T>& y) const;
272
273
281 virtual void usmv(T alpha, const GpuVector<T>& x, GpuVector<T>& y) const;
282
293 template <class MatrixType>
294 void updateNonzeroValues(const MatrixType& matrix, bool copyNonZeroElementsDirectly = false);
295
296private:
297 GpuVector<T> m_nonZeroElements;
298 GpuVector<int> m_columnIndices;
299 GpuVector<int> m_rowIndices;
300
301 // Notice that we store these three as int to make sure we are cusparse compatible.
302 //
303 // This gives the added benefit of checking the size constraints at construction of the matrix
304 // rather than in some call to cusparse.
305 const int m_numberOfNonzeroBlocks;
306 const int m_numberOfRows;
307 const int m_blockSize;
308
309 detail::GpuSparseMatrixDescriptionPtr m_matrixDescription;
310 detail::CuSparseHandle& m_cusparseHandle;
311
312 template <class VectorType>
313 void assertSameSize(const VectorType& vector) const;
314};
315} // namespace Opm::gpuistl
316#endif
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 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.
Definition GpuSparseMatrix.cpp:64
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