19#ifndef OPM_GPUILU0_OPM_Impl_HPP
20#define OPM_GPUILU0_OPM_Impl_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/GpuVector.hpp>
44template <
class M,
class X,
class Y,
int l = 1>
67 explicit OpmGpuILU0(
const M& A,
bool splitMatrix,
bool tuneKernels,
bool storeFactorizationAsFloat);
71 void pre(X& x, Y& b)
override;
74 void apply(X& v,
const Y& d)
override;
78 void post(X& x)
override;
81 Dune::SolverCategory::Category
category()
const override;
104 virtual bool hasPerfectUpdate()
const override {
111 void apply(X& v,
const Y& d,
int lowerSolveThreadBlockSize,
int upperSolveThreadBlockSize);
113 void update(
int moveThreadBlockSize,
int factorizationThreadBlockSize);
115 const M& m_cpuMatrix;
117 static constexpr const size_t blocksize_ = matrix_type::block_type::cols;
119 Opm::SparseTable<size_t> m_levelSets;
121 std::vector<int> m_reorderedToNatural;
123 std::vector<int> m_naturalToReordered;
126 std::unique_ptr<GpuMat> m_gpuReorderedLU;
128 std::unique_ptr<GpuMat> m_gpuMatrixReorderedLower;
129 std::unique_ptr<GpuMat> m_gpuMatrixReorderedUpper;
131 std::unique_ptr<FloatMat> m_gpuMatrixReorderedLowerFloat;
132 std::unique_ptr<FloatMat> m_gpuMatrixReorderedUpperFloat;
133 std::optional<GpuVector<float>> m_gpuMatrixReorderedDiagFloat;
135 std::optional<GpuVector<field_type>> m_gpuMatrixReorderedDiag;
137 GpuVector<int> m_gpuNaturalToReorder;
139 GpuVector<int> m_gpuReorderToNatural;
141 GpuVector<field_type> m_gpuDInv;
145 bool m_tuneThreadBlockSizes;
148 bool m_storeFactorizationAsFloat;
151 int m_upperSolveThreadBlockSize = -1;
152 int m_lowerSolveThreadBlockSize = -1;
153 int m_moveThreadBlockSize = -1;
154 int m_ILU0FactorizationThreadBlockSize = -1;
Interface class adding the update() method to the preconditioner interface.
Definition PreconditionerWithUpdate.hpp:32
ILU0 preconditioner on the GPU.
Definition OpmGpuILU0.hpp:46
GpuSparseMatrix< field_type > GpuMat
The GPU matrix type.
Definition OpmGpuILU0.hpp:57
static constexpr bool shouldCallPost()
Definition OpmGpuILU0.hpp:99
Y range_type
The range type of the preconditioner.
Definition OpmGpuILU0.hpp:53
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition OpmGpuILU0.cpp:106
void update() final
Updates the matrix data.
Definition OpmGpuILU0.cpp:226
void apply(X &v, const Y &d) override
Apply the preconditoner.
Definition OpmGpuILU0.cpp:112
static constexpr bool shouldCallPre()
Definition OpmGpuILU0.hpp:93
typename X::field_type field_type
The field type of the preconditioner.
Definition OpmGpuILU0.hpp:55
void LUFactorizeAndMoveData(int moveThreadBlockSize, int factorizationThreadBlockSize)
Compute LU factorization, and update the data of the reordered matrix.
Definition OpmGpuILU0.cpp:246
Dune::SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition OpmGpuILU0.cpp:219
void post(X &x) override
Post processing.
Definition OpmGpuILU0.cpp:213
void tuneThreadBlockSizes()
function that will experimentally tune the thread block sizes of the important cuda kernels
Definition OpmGpuILU0.cpp:328
typename std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition OpmGpuILU0.hpp:49
OpmGpuILU0(const M &A, bool splitMatrix, bool tuneKernels, bool storeFactorizationAsFloat)
Constructor.
Definition OpmGpuILU0.cpp:44
X domain_type
The domain type of the preconditioner.
Definition OpmGpuILU0.hpp:51