20 #ifndef OPM_OPENCLSOLVER_BACKEND_HEADER_INCLUDED
21 #define OPM_OPENCLSOLVER_BACKEND_HEADER_INCLUDED
23 #include <opm/simulators/linalg/bda/opencl/opencl.hpp>
24 #include <opm/simulators/linalg/bda/BdaResult.hpp>
25 #include <opm/simulators/linalg/bda/BdaSolver.hpp>
26 #include <opm/simulators/linalg/bda/ILUReorder.hpp>
27 #include <opm/simulators/linalg/bda/WellContributions.hpp>
29 #include <opm/simulators/linalg/bda/opencl/Preconditioner.hpp>
37 template <
unsigned int block_size>
46 using Base::verbosity;
47 using Base::platformID;
50 using Base::tolerance;
51 using Base::initialized;
55 std::vector<double> vals_contiguous;
58 cl::Buffer d_Avals, d_Acols, d_Arows;
59 cl::Buffer d_x, d_b, d_rb, d_r, d_rw, d_p;
60 cl::Buffer d_pw, d_s, d_t, d_v;
64 std::vector<cl::Device> devices;
66 bool useJacMatrix =
false;
68 std::unique_ptr<Preconditioner<block_size> > prec;
71 int *toOrder =
nullptr, *fromOrder =
nullptr;
72 bool analysis_done =
false;
73 std::shared_ptr<BlockedMatrix> mat =
nullptr;
74 std::shared_ptr<BlockedMatrix> jacMat =
nullptr;
76 ILUReorder opencl_ilu_reorder;
77 std::vector<cl::Event> events;
84 unsigned int ceilDivision(
const unsigned int A,
const unsigned int B);
91 double dot_w(cl::Buffer in1, cl::Buffer in2, cl::Buffer out);
98 double norm_w(cl::Buffer in, cl::Buffer out);
104 void axpy_w(cl::Buffer in,
const double a, cl::Buffer out);
109 void scale_w(cl::Buffer vec,
const double a);
118 void custom_w(cl::Buffer p, cl::Buffer v, cl::Buffer r,
const double omega,
const double beta);
128 void spmv_blocked_w(cl::Buffer vals, cl::Buffer cols, cl::Buffer rows, cl::Buffer x, cl::Buffer b);
138 void initialize(std::shared_ptr<BlockedMatrix> matrix, std::shared_ptr<BlockedMatrix> jacMatrix);
144 void copy_system_to_gpu();
153 void update_system_on_gpu();
157 bool analyze_matrix();
161 bool create_preconditioner();
169 std::shared_ptr<cl::Context> context;
170 std::shared_ptr<cl::CommandQueue> queue;
181 openclSolverBackend(
int linear_solver_verbosity,
int maxit,
double tolerance,
unsigned int platformID,
unsigned int deviceID,
182 ILUReorder opencl_ilu_reorder, std::string linsolver);
185 openclSolverBackend(
int linear_solver_verbosity,
int maxit,
double tolerance, ILUReorder opencl_ilu_reorder);
197 SolverStatus solve_system(std::shared_ptr<BlockedMatrix> matrix,
double *b,
212 void setOpencl(std::shared_ptr<cl::Context>& context, std::shared_ptr<cl::CommandQueue>& queue);
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:45
This struct resembles a blocked csr matrix, like Dune::BCRSMatrix.
Definition: BlockedMatrix.hpp:37
This class implements a opencl-based ilu0-bicgstab solver on GPU.
Definition: openclSolverBackend.hpp:39
~openclSolverBackend()
Destroy a openclSolver, and free memory.
Definition: openclSolverBackend.cpp:236
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: openclSolverBackend.cpp:670
openclSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int platformID, unsigned int deviceID, ILUReorder opencl_ilu_reorder, std::string linsolver)
Construct a openclSolver.
Definition: openclSolverBackend.cpp:50
void setOpencl(std::shared_ptr< cl::Context > &context, std::shared_ptr< cl::CommandQueue > &queue)
Set OpenCL objects This class either creates them based on platformID and deviceID or receives them t...
Definition: openclSolverBackend.cpp:230
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