19#ifndef OPM_GPUILU0_OPM_Impl_HPP
20#define OPM_GPUILU0_OPM_Impl_HPP
23#include <opm/grid/utility/SparseTable.hpp>
49template <
class CPUMatrixT,
class X,
class Y,
int l = 1>
78 const CPUMatrixT& cpuMatrix,
81 int mixedPrecisionScheme);
85 void pre(X& x, Y& b)
override;
88 void apply(X& v,
const Y& d)
override;
95 Dune::SolverCategory::Category
category()
const override;
128 void apply(X& v,
const Y& d,
int lowerSolveThreadBlockSize,
int upperSolveThreadBlockSize);
130 void update(
int moveThreadBlockSize,
int factorizationThreadBlockSize);
133 static constexpr const size_t blocksize_ = CPUMatrixT::block_type::cols;
137 std::vector<int> m_reorderedToNatural;
139 std::vector<int> m_naturalToReordered;
142 std::unique_ptr<GpuMatrix> m_gpuReorderedLU;
144 std::unique_ptr<GpuMatrix> m_gpuMatrixReorderedLower;
145 std::unique_ptr<GpuMatrix> m_gpuMatrixReorderedUpper;
147 std::unique_ptr<FloatMat> m_gpuMatrixReorderedLowerFloat;
148 std::unique_ptr<FloatMat> m_gpuMatrixReorderedUpperFloat;
149 std::optional<GpuVector<float>> m_gpuMatrixReorderedDiagFloat;
151 std::optional<GpuVector<field_type>> m_gpuMatrixReorderedDiag;
161 bool m_tuneThreadBlockSizes;
166 int m_upperSolveThreadBlockSize = -1;
167 int m_lowerSolveThreadBlockSize = -1;
168 int m_moveThreadBlockSize = -1;
169 int m_ILU0FactorizationThreadBlockSize = -1;
172 std::map<std::pair<field_type*, const field_type*>, GPUGraph> m_apply_graphs;
173 std::map<std::pair<field_type*, const field_type*>, GPUGraphExec> m_executableGraphs;
176 GPUStream m_stream{};
Interface class adding the update() method to the preconditioner interface.
Definition: PreconditionerWithUpdate.hpp:32
ILU0 preconditioner on the GPU.
Definition: OpmGpuILU0.hpp:51
static constexpr bool shouldCallPre()
Definition: OpmGpuILU0.hpp:110
void apply(X &v, const Y &d) override
Apply the preconditoner.
OpmGpuILU0(const GpuMatrix &gpuMatrix, const CPUMatrixT &cpuMatrix, bool splitMatrix, bool tuneKernels, int mixedPrecisionScheme)
Constructor.
Y range_type
The range type of the preconditioner.
Definition: OpmGpuILU0.hpp:56
void pre(X &x, Y &b) override
Prepare the preconditioner.
void LUFactorizeMatrix(int factorizationThreadBlockSize)
Compute LU factorization, and update the data of the reordered matrix.
void post(X &x) override
Post processing.
typename X::field_type field_type
The field type of the preconditioner.
Definition: OpmGpuILU0.hpp:58
Dune::SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
void tuneThreadBlockSizes()
function that will experimentally tune the thread block sizes of the important cuda kernels
void update() final
Updates the matrix data.
virtual bool hasPerfectUpdate() const override
Definition: OpmGpuILU0.hpp:121
void reorderAndSplitMatrix(int moveThreadBlockSize)
perform matrix splitting and reordering
static constexpr bool shouldCallPost()
Definition: OpmGpuILU0.hpp:116
X domain_type
The domain type of the preconditioner.
Definition: OpmGpuILU0.hpp:54
GpuSparseMatrix< field_type > GpuMatrix
The GPU matrix type.
Definition: OpmGpuILU0.hpp:60
Definition: AmgxInterface.hpp:38
MatrixStorageMPScheme
Definition: kernel_enums.hpp:31