My Project
rocsparseSolverBackend.hpp
1/*
2 Copyright 2023 Equinor ASA
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20#ifndef OPM_ROCSPARSESOLVER_BACKEND_HEADER_INCLUDED
21#define OPM_ROCSPARSESOLVER_BACKEND_HEADER_INCLUDED
22
23#include <memory>
24
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>
29
30#include <rocblas/rocblas.h>
31#include <rocsparse/rocsparse.h>
32
33#include <hip/hip_version.h>
34
35namespace Opm
36{
37namespace Accelerator
38{
39
41template <unsigned int block_size>
42class rocsparseSolverBackend : public BdaSolver<block_size>
43{
45
46 using Base::N;
47 using Base::Nb;
48 using Base::nnz;
49 using Base::nnzb;
50 using Base::verbosity;
51 using Base::platformID;
52 using Base::deviceID;
53 using Base::maxit;
54 using Base::tolerance;
55 using Base::initialized;
56
57private:
58
59 bool useJacMatrix = false;
60
61 bool analysis_done = false;
62 std::shared_ptr<BlockedMatrix> mat = nullptr; // original matrix
63 std::shared_ptr<BlockedMatrix> jacMat = nullptr; // matrix for preconditioner
64 int nnzbs_prec = 0; // number of nnz blocks in preconditioner matrix M
65
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;
74#endif
75 hipStream_t stream;
76
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; // vectors, used during linear solve
81 double *d_pw, *d_s, *d_t, *d_v;
82 void *d_buffer; // buffer space, used by rocsparse ilu0 analysis
83 int ver;
84 char rev[64];
85
86
90 void gpu_pbicgstab(WellContributions& wellContribs, BdaResult& res);
91
95 void initialize(std::shared_ptr<BlockedMatrix> matrix, std::shared_ptr<BlockedMatrix> jacMatrix);
96
99 void copy_system_to_gpu(double *b);
100
103 void update_system_on_gpu(double *b);
104
107 bool analyze_matrix();
108
111 bool create_preconditioner();
112
116 void solve_system(WellContributions &wellContribs, BdaResult &res);
117
118public:
125 rocsparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int platformID, unsigned int deviceID);
126
128 // rocsparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, ILUReorder opencl_ilu_reorder);
129
132
140 SolverStatus solve_system(std::shared_ptr<BlockedMatrix> matrix, double *b,
141 std::shared_ptr<BlockedMatrix> jacMatrix, WellContributions& wellContribs, BdaResult &res) override;
142
145 // SolverStatus solve_system(BdaResult &res);
146
149 void get_result(double *x) override;
150
151}; // end class rocsparseSolverBackend
152
153} // namespace Accelerator
154} // namespace Opm
155
156#endif
157
158
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