My Project
cusparseSolverBackend.hpp
1 /*
2  Copyright 2019 Equinor ASA
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 
20 #ifndef OPM_CUSPARSESOLVER_BACKEND_HEADER_INCLUDED
21 #define OPM_CUSPARSESOLVER_BACKEND_HEADER_INCLUDED
22 
23 
24 #include "cublas_v2.h"
25 #include "cusparse_v2.h"
26 
27 #include <opm/simulators/linalg/bda/BdaResult.hpp>
28 #include <opm/simulators/linalg/bda/BdaSolver.hpp>
29 #include <opm/simulators/linalg/bda/WellContributions.hpp>
30 
31 namespace bda
32 {
33 
35 template <unsigned int block_size>
36 class cusparseSolverBackend : public BdaSolver<block_size> {
37 
39 
40  using Base::N;
41  using Base::Nb;
42  using Base::nnz;
43  using Base::nnzb;
44  using Base::verbosity;
45  using Base::deviceID;
46  using Base::maxit;
47  using Base::tolerance;
48  using Base::initialized;
49 
50 private:
51 
52  cublasHandle_t cublasHandle;
53  cusparseHandle_t cusparseHandle;
54  cudaStream_t stream;
55  cusparseMatDescr_t descr_B, descr_M, descr_L, descr_U;
56  bsrilu02Info_t info_M;
57  bsrsv2Info_t info_L, info_U;
58  // b: bsr matrix, m: preconditioner
59  double *d_bVals, *d_mVals;
60  int *d_bCols, *d_mCols;
61  int *d_bRows, *d_mRows;
62  double *d_x, *d_b, *d_r, *d_rw, *d_p; // vectors, used during linear solve
63  double *d_pw, *d_s, *d_t, *d_v;
64  void *d_buffer;
65  double *vals_contiguous; // only used if COPY_ROW_BY_ROW is true in cusparseSolverBackend.cpp
66 
67  bool analysis_done = false;
68 
69 
73  void gpu_pbicgstab(WellContributions& wellContribs, BdaResult& res);
74 
79  void initialize(int N, int nnz, int dim);
80 
82  void finalize();
83 
89  void copy_system_to_gpu(double *vals, int *rows, int *cols, double *b);
90 
91  // Update linear system on GPU, don't copy rowpointers and colindices, they stay the same
95  void update_system_on_gpu(double *vals, int *rows, double *b);
96 
98  void reset_prec_on_gpu();
99 
102  bool analyse_matrix();
103 
106  bool create_preconditioner();
107 
111  void solve_system(WellContributions& wellContribs, BdaResult &res);
112 
113 public:
114 
115 
121  cusparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int deviceID);
122 
125 
137  SolverStatus solve_system(int N, int nnz, int dim, double *vals, int *rows, int *cols, double *b, WellContributions& wellContribs, BdaResult &res) override;
138 
141  void get_result(double *x) override;
142 
143 }; // end class cusparseSolverBackend
144 
145 } // namespace bda
146 
147 #endif
148 
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:61
This class is based on InverseOperatorResult struct from dune/istl/solver.hh It is needed to prevent ...
Definition: BdaResult.hpp:29
This class serves to simplify choosing between different backend solvers, such as cusparseSolver and ...
Definition: BdaSolver.hpp:43
This class implements a cusparse-based ilu0-bicgstab solver on GPU.
Definition: cusparseSolverBackend.hpp:36
SolverStatus solve_system(int N, int nnz, int dim, double *vals, int *rows, int *cols, double *b, WellContributions &wellContribs, BdaResult &res) override
Solve linear system, A*x = b, matrix A must be in blocked-CSR format.
void get_result(double *x) override
Get resulting vector x after linear solve, also includes post processing if necessary.
cusparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int deviceID)
Construct a cusparseSolver.
~cusparseSolverBackend()
Destroy a cusparseSolver, and free memory.