20#ifndef OPM_ROCSPARSESOLVER_BACKEND_HEADER_INCLUDED
21#define OPM_ROCSPARSESOLVER_BACKEND_HEADER_INCLUDED
25#include <opm/simulators/linalg/bda/opencl/opencl.hpp>
26#include <opm/simulators/linalg/bda/BdaResult.hpp>
27#include <opm/simulators/linalg/bda/BdaSolver.hpp>
28#include <opm/simulators/linalg/bda/WellContributions.hpp>
30#include <rocblas/rocblas.h>
31#include <rocsparse/rocsparse.h>
33#include <hip/hip_version.h>
41template <
unsigned int block_size>
50 using Base::verbosity;
51 using Base::platformID;
54 using Base::tolerance;
55 using Base::initialized;
59 bool useJacMatrix =
false;
61 bool analysis_done =
false;
62 std::shared_ptr<BlockedMatrix> mat =
nullptr;
63 std::shared_ptr<BlockedMatrix> jacMat =
nullptr;
66 rocsparse_direction dir = rocsparse_direction_row;
67 rocsparse_operation operation = rocsparse_operation_none;
68 rocsparse_handle handle;
69 rocblas_handle blas_handle;
70 rocsparse_mat_descr descr_A, descr_M, descr_L, descr_U;
71 rocsparse_mat_info ilu_info;
72#if HIP_VERSION >= 50400000
73 rocsparse_mat_info spmv_info;
77 rocsparse_int *d_Arows, *d_Mrows;
78 rocsparse_int *d_Acols, *d_Mcols;
79 double *d_Avals, *d_Mvals;
80 double *d_x, *d_b, *d_r, *d_rw, *d_p;
81 double *d_pw, *d_s, *d_t, *d_v;
95 void initialize(std::shared_ptr<BlockedMatrix> matrix, std::shared_ptr<BlockedMatrix> jacMatrix);
99 void copy_system_to_gpu(
double *b);
103 void update_system_on_gpu(
double *b);
107 bool analyze_matrix();
111 bool create_preconditioner();
125 rocsparseSolverBackend(
int linear_solver_verbosity,
int maxit,
double tolerance,
unsigned int platformID,
unsigned int deviceID);
140 SolverStatus solve_system(std::shared_ptr<BlockedMatrix> matrix,
double *b,
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:46
This class implements a rocsparse-based ilu0-bicgstab solver on GPU.
Definition: rocsparseSolverBackend.hpp:43
void get_result(double *x) override
Solve scalar linear system, for example a coarse system of an AMG preconditioner Data is already on t...
Definition: rocsparseSolverBackend.cpp:554
~rocsparseSolverBackend()
For the CPR coarse solver.
Definition: rocsparseSolverBackend.cpp:110
rocsparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int platformID, unsigned int deviceID)
Construct a openclSolver.
Definition: rocsparseSolverBackend.cpp:85
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:52
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27