20 #ifndef OPM_CUSPARSESOLVER_BACKEND_HEADER_INCLUDED
21 #define OPM_CUSPARSESOLVER_BACKEND_HEADER_INCLUDED
24 #include "cublas_v2.h"
25 #include "cusparse_v2.h"
27 #include <opm/simulators/linalg/bda/BdaResult.hpp>
28 #include <opm/simulators/linalg/bda/BdaSolver.hpp>
29 #include <opm/simulators/linalg/bda/WellContributions.hpp>
37 template <
unsigned int block_size>
46 using Base::verbosity;
49 using Base::tolerance;
50 using Base::initialized;
54 cublasHandle_t cublasHandle;
55 cusparseHandle_t cusparseHandle;
57 cusparseMatDescr_t descr_B, descr_M, descr_L, descr_U;
58 bsrilu02Info_t info_M;
59 bsrsv2Info_t info_L, info_U;
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;
65 double *d_pw, *d_s, *d_t, *d_v;
67 double *vals_contiguous;
69 bool analysis_done =
false;
81 void initialize(
int N,
int nnz,
int dim);
91 void copy_system_to_gpu(
double *vals,
int *rows,
int *cols,
double *b);
97 void update_system_on_gpu(
double *vals,
int *rows,
double *b);
100 void reset_prec_on_gpu();
104 bool analyse_matrix();
108 bool create_preconditioner();
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