OpmGpuILU0.hpp
Go to the documentation of this file.
1/*
2 Copyright 2024 SINTEF AS
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#ifndef OPM_GPUILU0_OPM_Impl_HPP
20#define OPM_GPUILU0_OPM_Impl_HPP
21
22#include <memory>
23#include <opm/grid/utility/SparseTable.hpp>
29#include <optional>
30#include <type_traits>
31#include <vector>
32
33
34namespace Opm::gpuistl
35{
49template <class CPUMatrixT, class X, class Y, int l = 1>
51{
52public:
54 using domain_type = X;
56 using range_type = Y;
58 using field_type = typename X::field_type;
63
66
67
74 explicit OpmGpuILU0(const GpuMatrix& gpuMatrix, const CPUMatrixT& cpuMatrix, bool splitMatrix, bool tuneKernels, int mixedPrecisionScheme);
75
78 void pre(X& x, Y& b) override;
79
81 void apply(X& v, const Y& d) override;
82
85 void post(X& x) override;
86
88 Dune::SolverCategory::Category category() const override;
89
91 void update() final;
92
94 void reorderAndSplitMatrix(int moveThreadBlockSize);
95
97 void LUFactorizeMatrix(int factorizationThreadBlockSize);
98
101
103 static constexpr bool shouldCallPre()
104 {
105 return false;
106 }
107
109 static constexpr bool shouldCallPost()
110 {
111 return false;
112 }
113
114 virtual bool hasPerfectUpdate() const override {
115 return true;
116 }
117
118
119private:
121 void apply(X& v, const Y& d, int lowerSolveThreadBlockSize, int upperSolveThreadBlockSize);
123 void update(int moveThreadBlockSize, int factorizationThreadBlockSize);
124
126 static constexpr const size_t blocksize_ = CPUMatrixT::block_type::cols;
128 Opm::SparseTable<size_t> m_levelSets;
130 std::vector<int> m_reorderedToNatural;
132 std::vector<int> m_naturalToReordered;
134 const GpuMatrix& m_gpuMatrix;
135 std::unique_ptr<GpuMatrix> m_gpuReorderedLU;
137 std::unique_ptr<GpuMatrix> m_gpuMatrixReorderedLower;
138 std::unique_ptr<GpuMatrix> m_gpuMatrixReorderedUpper;
140 std::unique_ptr<FloatMat> m_gpuMatrixReorderedLowerFloat;
141 std::unique_ptr<FloatMat> m_gpuMatrixReorderedUpperFloat;
142 std::optional<GpuVector<float>> m_gpuMatrixReorderedDiagFloat;
144 std::optional<GpuVector<field_type>> m_gpuMatrixReorderedDiag;
146 GpuVector<int> m_gpuNaturalToReorder;
148 GpuVector<int> m_gpuReorderToNatural;
150 GpuVector<field_type> m_gpuDInv;
152 bool m_splitMatrix;
154 bool m_tuneThreadBlockSizes;
156 const MatrixStorageMPScheme m_mixedPrecisionScheme;
159 int m_upperSolveThreadBlockSize = -1;
160 int m_lowerSolveThreadBlockSize = -1;
161 int m_moveThreadBlockSize = -1;
162 int m_ILU0FactorizationThreadBlockSize = -1;
163
164 // Graphs for Apply
165 std::map<std::pair<field_type*, const field_type*>, GPUGraph> m_apply_graphs;
166 std::map<std::pair<field_type*, const field_type*>, GPUGraphExec> m_executableGraphs;
167
168 // Stream for the DILU operations on the GPU
169 GPUStream m_stream{};
170 // Events for synchronization with main stream
171 GPUEvent m_before{};
172 GPUEvent m_after{};
173};
174} // end namespace Opm::gpuistl
175
176#endif
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:103
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:114
void reorderAndSplitMatrix(int moveThreadBlockSize)
perform matrix splitting and reordering
static constexpr bool shouldCallPost()
Definition: OpmGpuILU0.hpp:109
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: autotuner.hpp:30
MatrixStorageMPScheme
Definition: kernel_enums.hpp:31