My Project
Loading...
Searching...
No Matches
OpmGpuILU0.hpp
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>
24#include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
25#include <opm/simulators/linalg/gpuistl/GpuSparseMatrix.hpp>
26#include <opm/simulators/linalg/gpuistl/GpuVector.hpp>
27#include <opm/simulators/linalg/gpuistl/gpu_resources.hpp>
28#include <opm/simulators/linalg/gpuistl/detail/kernel_enums.hpp>
29#include <optional>
30#include <type_traits>
31#include <vector>
32
33
34namespace Opm::gpuistl
35{
46template <class M, class X, class Y, int l = 1>
48{
49public:
51 using matrix_type = typename std::remove_const<M>::type;
53 using domain_type = X;
55 using range_type = Y;
57 using field_type = typename X::field_type;
62
69 explicit OpmGpuILU0(const M& A, bool splitMatrix, bool tuneKernels, int mixedPrecisionScheme);
70
73 void pre(X& x, Y& b) override;
74
76 void apply(X& v, const Y& d) override;
77
80 void post(X& x) override;
81
83 Dune::SolverCategory::Category category() const override;
84
86 void update() final;
87
90
93
96
99 {
100 return false;
101 }
102
104 static constexpr bool shouldCallPost()
105 {
106 return false;
107 }
108
109 virtual bool hasPerfectUpdate() const override {
110 return true;
111 }
112
113
114private:
116 void apply(X& v, const Y& d, int lowerSolveThreadBlockSize, int upperSolveThreadBlockSize);
120 const M& m_cpuMatrix;
122 static constexpr const size_t blocksize_ = matrix_type::block_type::cols;
124 Opm::SparseTable<size_t> m_levelSets;
126 std::vector<int> m_reorderedToNatural;
128 std::vector<int> m_naturalToReordered;
130 GpuMat m_gpuMatrix;
131 std::unique_ptr<GpuMat> m_gpuReorderedLU;
133 std::unique_ptr<GpuMat> m_gpuMatrixReorderedLower;
134 std::unique_ptr<GpuMat> m_gpuMatrixReorderedUpper;
136 std::unique_ptr<FloatMat> m_gpuMatrixReorderedLowerFloat;
137 std::unique_ptr<FloatMat> m_gpuMatrixReorderedUpperFloat;
138 std::optional<GpuVector<float>> m_gpuMatrixReorderedDiagFloat;
140 std::optional<GpuVector<field_type>> m_gpuMatrixReorderedDiag;
142 GpuVector<int> m_gpuNaturalToReorder;
144 GpuVector<int> m_gpuReorderToNatural;
146 GpuVector<field_type> m_gpuDInv;
148 bool m_splitMatrix;
150 bool m_tuneThreadBlockSizes;
152 const MatrixStorageMPScheme m_mixedPrecisionScheme;
155 int m_upperSolveThreadBlockSize = -1;
156 int m_lowerSolveThreadBlockSize = -1;
157 int m_moveThreadBlockSize = -1;
158 int m_ILU0FactorizationThreadBlockSize = -1;
159
160 // Graphs for Apply
161 std::map<std::pair<field_type*, const field_type*>, GPUGraph> m_apply_graphs;
162 std::map<std::pair<field_type*, const field_type*>, GPUGraphExec> m_executableGraphs;
163
164 // Stream for the DILU operations on the GPU
165 GPUStream m_stream{};
166 // Events for synchronization with main stream
167 GPUEvent m_before{};
168 GPUEvent m_after{};
169};
170} // end namespace Opm::gpuistl
171
172#endif
Interface class adding the update() method to the preconditioner interface.
Definition PreconditionerWithUpdate.hpp:32
Definition BlackoilWellModel.hpp:83
ILU0 preconditioner on the GPU.
Definition OpmGpuILU0.hpp:48
GpuSparseMatrix< field_type > GpuMat
The GPU matrix type.
Definition OpmGpuILU0.hpp:59
static constexpr bool shouldCallPost()
Definition OpmGpuILU0.hpp:104
void LUFactorizeMatrix(int factorizationThreadBlockSize)
Compute LU factorization, and update the data of the reordered matrix.
Definition OpmGpuILU0.cpp:329
Y range_type
The range type of the preconditioner.
Definition OpmGpuILU0.hpp:55
void pre(X &x, Y &b) override
Prepare the preconditioner.
Definition OpmGpuILU0.cpp:111
void update() final
Updates the matrix data.
Definition OpmGpuILU0.cpp:272
void apply(X &v, const Y &d) override
Apply the preconditoner.
Definition OpmGpuILU0.cpp:117
void reorderAndSplitMatrix(int moveThreadBlockSize)
perform matrix splitting and reordering
Definition OpmGpuILU0.cpp:301
static constexpr bool shouldCallPre()
Definition OpmGpuILU0.hpp:98
typename X::field_type field_type
The field type of the preconditioner.
Definition OpmGpuILU0.hpp:57
Dune::SolverCategory::Category category() const override
Category of the preconditioner (see SolverCategory::Category)
Definition OpmGpuILU0.cpp:265
void post(X &x) override
Post processing.
Definition OpmGpuILU0.cpp:259
void tuneThreadBlockSizes()
function that will experimentally tune the thread block sizes of the important cuda kernels
Definition OpmGpuILU0.cpp:407
typename std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition OpmGpuILU0.hpp:51
X domain_type
The domain type of the preconditioner.
Definition OpmGpuILU0.hpp:53
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242