19#ifndef OPM_GPUDILU_HPP
20#define OPM_GPUDILU_HPP
23#include <opm/grid/utility/SparseTable.hpp>
24#include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
25#include <opm/simulators/linalg/gpuistl/GpuSparseMatrix.hpp>
26#include <opm/simulators/linalg/gpuistl/detail/kernel_enums.hpp>
27#include <opm/simulators/linalg/gpuistl/gpu_resources.hpp>
45template <
class M,
class X,
class Y,
int l = 1>
72 void pre(X& x, Y&
b)
override;
75 void apply(X&
v,
const Y&
d)
override;
79 void post(X& x)
override;
82 Dune::SolverCategory::Category
category()
const override;
109 virtual bool hasPerfectUpdate()
const override {
120 const M& m_cpuMatrix;
122 static constexpr const size_t blocksize_ = matrix_type::block_type::cols;
126 std::vector<int> m_reorderedToNatural;
128 std::vector<int> m_naturalToReordered;
132 std::unique_ptr<CuMat> m_gpuMatrixReordered;
134 std::unique_ptr<CuMat> m_gpuMatrixReorderedLower;
135 std::unique_ptr<CuMat> m_gpuMatrixReorderedUpper;
137 std::unique_ptr<GpuVector<field_type>> m_gpuMatrixReorderedDiag;
139 std::unique_ptr<FloatMat> m_gpuMatrixReorderedLowerFloat;
140 std::unique_ptr<FloatMat> m_gpuMatrixReorderedUpperFloat;
141 std::unique_ptr<FloatVec> m_gpuMatrixReorderedDiagFloat;
142 std::unique_ptr<FloatVec> m_gpuDInvFloat;
152 bool m_tuneThreadBlockSizes;
154 const MatrixStorageMPScheme m_mixedPrecisionScheme;
157 int m_upperSolveThreadBlockSize = -1;
158 int m_lowerSolveThreadBlockSize = -1;
159 int m_moveThreadBlockSize = -1;
160 int m_DILUFactorizationThreadBlockSize = -1;
163 std::map<std::pair<field_type*, const field_type*>,
GPUGraph> m_apply_graphs;
164 std::map<std::pair<field_type*, const field_type*>,
GPUGraphExec> m_executableGraphs;
Interface class adding the update() method to the preconditioner interface.
Definition PreconditionerWithUpdate.hpp:32
Definition BlackoilWellModel.hpp:83
DILU preconditioner on the GPU.
Definition GpuDILU.hpp:47
void apply(X &v, const Y &d) override
Apply the preconditoner.
Definition GpuDILU.cpp:113
Y range_type
The range type of the preconditioner.
Definition GpuDILU.hpp:54
void computeDiagonal(int factorizationThreadBlockSize)
Compute the diagonal of the DILU, and update the data of the reordered matrix.
Definition GpuDILU.cpp:339
GpuSparseMatrix< field_type > CuMat
The GPU matrix type.
Definition GpuDILU.hpp:58
typename std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition GpuDILU.hpp:50
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition GpuDILU.cpp:107
void post(X &x) override
Post processing.
Definition GpuDILU.cpp:269
typename X::field_type field_type
The field type of the preconditioner.
Definition GpuDILU.hpp:56
static constexpr bool shouldCallPre()
Definition GpuDILU.hpp:98
Dune::SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition GpuDILU.cpp:275
void reorderAndSplitMatrix(int moveThreadBlockSize)
perform matrix splitting and reordering
Definition GpuDILU.cpp:311
void tuneThreadBlockSizes()
function that will experimentally tune the thread block sizes of the important cuda kernels
Definition GpuDILU.cpp:419
X domain_type
The domain type of the preconditioner.
Definition GpuDILU.hpp:52
static constexpr bool shouldCallPost()
Definition GpuDILU.hpp:104
void update() final
Updates the matrix data.
Definition GpuDILU.cpp:282
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242