20 #ifndef OPM_OPENCLSOLVER_BACKEND_HEADER_INCLUDED
21 #define OPM_OPENCLSOLVER_BACKEND_HEADER_INCLUDED
23 #include <opm/simulators/linalg/bda/opencl.hpp>
24 #include <opm/simulators/linalg/bda/openclKernels.hpp>
25 #include <opm/simulators/linalg/bda/BdaResult.hpp>
26 #include <opm/simulators/linalg/bda/BdaSolver.hpp>
27 #include <opm/simulators/linalg/bda/ILUReorder.hpp>
28 #include <opm/simulators/linalg/bda/WellContributions.hpp>
29 #include <opm/simulators/linalg/bda/BILU0.hpp>
37 template <
unsigned int block_size>
47 using Base::verbosity;
48 using Base::platformID;
51 using Base::tolerance;
52 using Base::initialized;
56 double *vals_contiguous =
nullptr;
59 cl::Buffer d_Avals, d_Acols, d_Arows;
60 cl::Buffer d_x, d_b, d_rb, d_r, d_rw, d_p;
61 cl::Buffer d_pw, d_s, d_t, d_v;
64 double *tmp =
nullptr;
67 std::vector<cl::Device> devices;
68 std::unique_ptr<cl::KernelFunctor<cl::Buffer&, cl::Buffer&, cl::Buffer&, const unsigned int, cl::LocalSpaceArg> > dot_k;
69 std::unique_ptr<cl::KernelFunctor<cl::Buffer&, cl::Buffer&, const unsigned int, cl::LocalSpaceArg> > norm_k;
70 std::unique_ptr<cl::KernelFunctor<cl::Buffer&, const double, cl::Buffer&, const unsigned int> > axpy_k;
71 std::unique_ptr<cl::KernelFunctor<cl::Buffer&, const double, const unsigned int> > scale_k;
72 std::unique_ptr<cl::KernelFunctor<cl::Buffer&, cl::Buffer&, cl::Buffer&, const double, const double, const unsigned int> > custom_k;
73 std::unique_ptr<spmv_kernel_type> spmv_blocked_k;
74 std::shared_ptr<ilu_apply1_kernel_type> ILU_apply1_k;
75 std::shared_ptr<ilu_apply2_kernel_type> ILU_apply2_k;
76 std::shared_ptr<stdwell_apply_kernel_type> stdwell_apply_k;
77 std::shared_ptr<stdwell_apply_no_reorder_kernel_type> stdwell_apply_no_reorder_k;
78 std::shared_ptr<ilu_decomp_kernel_type> ilu_decomp_k;
81 int *toOrder =
nullptr, *fromOrder =
nullptr;
82 bool analysis_done =
false;
83 std::unique_ptr<BlockedMatrix<block_size> > mat =
nullptr;
85 ILUReorder opencl_ilu_reorder;
86 std::vector<cl::Event> events;
93 unsigned int ceilDivision(
const unsigned int A,
const unsigned int B);
100 double dot_w(cl::Buffer in1, cl::Buffer in2, cl::Buffer out);
107 double norm_w(cl::Buffer in, cl::Buffer out);
113 void axpy_w(cl::Buffer in,
const double a, cl::Buffer out);
118 void scale_w(cl::Buffer vec,
const double a);
127 void custom_w(cl::Buffer p, cl::Buffer v, cl::Buffer r,
const double omega,
const double beta);
137 void spmv_blocked_w(cl::Buffer vals, cl::Buffer cols, cl::Buffer rows, cl::Buffer x, cl::Buffer b);
151 void initialize(
int N,
int nnz,
int dim,
double *vals,
int *rows,
int *cols);
154 void get_opencl_kernels();
160 void copy_system_to_gpu();
169 void update_system_on_gpu();
173 bool analyse_matrix();
177 bool create_preconditioner();
185 std::shared_ptr<cl::Context> context;
186 std::shared_ptr<cl::CommandQueue> queue;
195 openclSolverBackend(
int linear_solver_verbosity,
int maxit,
double tolerance,
unsigned int platformID,
unsigned int deviceID, ILUReorder opencl_ilu_reorder);
215 SolverStatus solve_system(
int N,
int nnz,
int dim,
double *vals,
int *rows,
int *cols,
double *b,
WellContributions& wellContribs,
BdaResult &res)
override;
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:61
This class implementa a Blocked ILU0 preconditioner The decomposition is done on CPU,...
Definition: BILU0.hpp:51
This class is based on InverseOperatorResult struct from dune/istl/solver.hh It is needed to prevent ...
Definition: BdaResult.hpp:29
This class serves to simplify choosing between different backend solvers, such as cusparseSolver and ...
Definition: BdaSolver.hpp:43
This struct resembles a blocked csr matrix, like Dune::BCRSMatrix.
Definition: BlockedMatrix.hpp:36
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:191
openclSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int platformID, unsigned int deviceID, ILUReorder opencl_ilu_reorder)
Construct a openclSolver.
Definition: openclSolverBackend.cpp:45
void get_result(double *x) override
Get result after linear solve, and peform postprocessing if necessary.
Definition: openclSolverBackend.cpp:767