My Project
Loading...
Searching...
No Matches
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/BdaResult.hpp>
26#include <opm/simulators/linalg/bda/BdaSolver.hpp>
27#include <opm/simulators/linalg/bda/WellContributions.hpp>
28
29#include <rocblas/rocblas.h>
30#include <rocsparse/rocsparse.h>
31
32#include <hip/hip_version.h>
33
34namespace Opm
35{
36namespace Accelerator
37{
38
40template <unsigned int block_size>
41class rocsparseSolverBackend : public BdaSolver<block_size>
42{
44
45 using Base::N;
46 using Base::Nb;
47 using Base::nnz;
48 using Base::nnzb;
49 using Base::verbosity;
50 using Base::platformID;
51 using Base::deviceID;
52 using Base::maxit;
53 using Base::tolerance;
54 using Base::initialized;
55
56private:
57
58 bool useJacMatrix = false;
59
60 bool analysis_done = false;
61 std::shared_ptr<BlockedMatrix> mat = nullptr; // original matrix
62 std::shared_ptr<BlockedMatrix> jacMat = nullptr; // matrix for preconditioner
63 int nnzbs_prec = 0; // number of nnz blocks in preconditioner matrix M
64
65 rocsparse_direction dir = rocsparse_direction_row;
66 rocsparse_operation operation = rocsparse_operation_none;
67 rocsparse_handle handle;
68 rocblas_handle blas_handle;
69 rocsparse_mat_descr descr_A, descr_M, descr_L, descr_U;
70 rocsparse_mat_info ilu_info;
71#if HIP_VERSION >= 50400000
72 rocsparse_mat_info spmv_info;
73#endif
74 hipStream_t stream;
75
76 rocsparse_int *d_Arows, *d_Mrows;
77 rocsparse_int *d_Acols, *d_Mcols;
78 double *d_Avals, *d_Mvals;
79 double *d_x, *d_b, *d_r, *d_rw, *d_p; // vectors, used during linear solve
80 double *d_pw, *d_s, *d_t, *d_v;
81 void *d_buffer; // buffer space, used by rocsparse ilu0 analysis
82 int ver;
83 char rev[64];
84
85
89 void gpu_pbicgstab(WellContributions& wellContribs, BdaResult& res);
90
94 void initialize(std::shared_ptr<BlockedMatrix> matrix, std::shared_ptr<BlockedMatrix> jacMatrix);
95
98 void copy_system_to_gpu(double *b);
99
102 void update_system_on_gpu(double *b);
103
106 bool analyze_matrix();
107
110 bool create_preconditioner();
111
115 void solve_system(WellContributions &wellContribs, BdaResult &res);
116
117public:
124 rocsparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, unsigned int platformID, unsigned int deviceID);
125
127 // rocsparseSolverBackend(int linear_solver_verbosity, int maxit, double tolerance, ILUReorder opencl_ilu_reorder);
128
131
139 SolverStatus solve_system(std::shared_ptr<BlockedMatrix> matrix, double *b,
140 std::shared_ptr<BlockedMatrix> jacMatrix, WellContributions& wellContribs, BdaResult &res) override;
141
144 // SolverStatus solve_system(BdaResult &res);
145
148 void get_result(double *x) override;
149
150}; // end class rocsparseSolverBackend
151
152} // namespace Accelerator
153} // namespace Opm
154
155#endif
156
157
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:42
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:600
~rocsparseSolverBackend()
For the CPR coarse solver.
Definition rocsparseSolverBackend.cpp:125
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