20#ifndef OPM_ROCSPARSEBILU0_HPP
21#define OPM_ROCSPARSEBILU0_HPP
23#include <opm/simulators/linalg/bda/BlockedMatrix.hpp>
25#include <opm/simulators/linalg/bda/rocm/rocsparsePreconditioner.hpp>
27#include <rocblas/rocblas.h>
28#include <rocsparse/rocsparse.h>
30#include <hip/hip_version.h>
32namespace Opm::Accelerator {
37template <
class Scalar,
unsigned int block_size>
46 using Base::verbosity;
50 rocsparse_mat_descr descr_M, descr_L, descr_U;
51 rocsparse_mat_info ilu_info;
52#if HIP_VERSION >= 50400000
53 rocsparse_mat_info spmv_info;
56 rocsparse_int *d_Mrows, *d_Mcols;
57 Scalar *d_Mvals, *d_t;
69 rocsparse_int *d_Arows,
70 rocsparse_int *d_Acols)
override;
102 void apply(
const cl::Buffer&, cl::Buffer&)
override {}
This struct resembles a blocked csr matrix, like Dune::BCRSMatrix.
Definition rocsparsePreconditioner.hpp:29
This class implements a Blocked ILU0 preconditioner The decomposition is done on GPU,...
Definition rocsparseBILU0.hpp:39
void apply(Scalar &y, Scalar &x) override
Apply preconditioner, x = prec(y) via Lz = y and Ux = z.
Definition rocsparseBILU0.cpp:329
bool create_preconditioner(BlockedMatrix< Scalar > *mat) override
ILU decomposition.
Definition rocsparseBILU0.cpp:222
void update_system_on_gpu(Scalar *b) override
Reassign pointers, in case the addresses of the Dune variables have changed --> TODO: check when/if w...
Definition rocsparseBILU0.cpp:300
bool initialize(std::shared_ptr< BlockedMatrix< Scalar > > matrix, std::shared_ptr< BlockedMatrix< Scalar > > jacMatrix, rocsparse_int *d_Arows, rocsparse_int *d_Acols) override
Initialize GPU and allocate memory.
Definition rocsparseBILU0.cpp:55
bool analyze_matrix(BlockedMatrix< Scalar > *mat) override
Analysis, extract parallelism if specified.
Definition rocsparseBILU0.cpp:89
void copy_system_to_gpu(Scalar *mVals) override
Copy matrix A values to GPU.
Definition rocsparseBILU0.cpp:268
Definition rocsparsePreconditioner.hpp:33