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 Opm
32 {
33 namespace Accelerator
34 {
35 
37 template <unsigned int block_size>
38 class cusparseSolverBackend : public BdaSolver<block_size> {
39 
41 
42  using Base::N;
43  using Base::Nb;
44  using Base::nnz;
45  using Base::nnzb;
46  using Base::verbosity;
47  using Base::deviceID;
48  using Base::maxit;
49  using Base::tolerance;
50  using Base::initialized;
51 
52 private:
53 
54  cublasHandle_t cublasHandle;
55  cusparseHandle_t cusparseHandle;
56  cudaStream_t stream;
57  cusparseMatDescr_t descr_B, descr_M, descr_L, descr_U;
58  bsrilu02Info_t info_M;
59  bsrsv2Info_t info_L, info_U;
60  // b: bsr matrix, m: preconditioner
61  double *d_bVals, *d_mVals;
62  int *d_bCols, *d_mCols;
63  int *d_bRows, *d_mRows;
64  double *d_x, *d_b, *d_r, *d_rw, *d_p; // vectors, used during linear solve
65  double *d_pw, *d_s, *d_t, *d_v;
66  void *d_buffer;
67  double *vals_contiguous; // only used if COPY_ROW_BY_ROW is true in cusparseSolverBackend.cpp
68 
69  bool analysis_done = false;
70 
71 
75  void gpu_pbicgstab(WellContributions& wellContribs, BdaResult& res);
76 
81  void initialize(int N, int nnz, int dim);
82 
84  void finalize();
85 
91  void copy_system_to_gpu(double *vals, int *rows, int *cols, double *b);
92 
93  // Update linear system on GPU, don't copy rowpointers and colindices, they stay the same
97  void update_system_on_gpu(double *vals, int *rows, double *b);
98 
100  void reset_prec_on_gpu();
101 
104  bool analyse_matrix();
105 
108  bool create_preconditioner();
109 
113  void solve_system(WellContributions& wellContribs, BdaResult &res);
114 
115 public:
116 
117 
123  cusparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int deviceID);
124 
127 
139  SolverStatus solve_system(int N, int nnz, int dim, double *vals, int *rows, int *cols, double *b, WellContributions& wellContribs, BdaResult &res) override;
140 
143  void get_result(double *x) override;
144 
145 }; // end class cusparseSolverBackend
146 
147 } // namespace Accelerator
148 } // namespace Opm
149 
150 #endif
151 
This class is based on InverseOperatorResult struct from dune/istl/solver.hh It is needed to prevent ...
Definition: BdaResult.hpp:31
This class serves to simplify choosing between different backend solvers, such as cusparseSolver and ...
Definition: BdaSolver.hpp:44
This class implements a cusparse-based ilu0-bicgstab solver on GPU.
Definition: cusparseSolverBackend.hpp:38
~cusparseSolverBackend()
Destroy a cusparseSolver, and free memory.
cusparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int deviceID)
Construct a cusparseSolver.
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.
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:53
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27