20#ifndef OPM_CUSPARSESOLVER_BACKEND_HEADER_INCLUDED
21#define OPM_CUSPARSESOLVER_BACKEND_HEADER_INCLUDED
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>
37template <
unsigned int block_size>
46 using Base::verbosity;
49 using Base::tolerance;
50 using Base::initialized;
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;
71 bool useJacMatrix =
false;
83 void initialize(std::shared_ptr<BlockedMatrix> matrix, std::shared_ptr<BlockedMatrix> jacMatrix);
93 void copy_system_to_gpu(std::shared_ptr<BlockedMatrix> matrix,
double *
b, std::shared_ptr<BlockedMatrix> jacMatrix);
100 void update_system_on_gpu(std::shared_ptr<BlockedMatrix> matrix,
double *
b, std::shared_ptr<BlockedMatrix> jacMatrix);
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:46
This class implements a cusparse-based ilu0-bicgstab solver on GPU.
Definition cusparseSolverBackend.hpp:38
~cusparseSolverBackend()
Destroy a cusparseSolver, and free memory.
SolverStatus solve_system(std::shared_ptr< BlockedMatrix > matrix, double *b, std::shared_ptr< BlockedMatrix > jacMatrix, WellContributions &wellContribs, BdaResult &res) override
Solve linear system, A*x = b, matrix A must be in blocked-CSR format.
cusparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int deviceID)
Construct a cusparseSolver.
void get_result(double *x) override
Get resulting vector x after linear solve, also includes post processing if necessary.
Definition AquiferInterface.hpp:35
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