My Project
BISAI.hpp
1 /*
2  Copyright 2022 Equinor ASA
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 
20 #ifndef BISAI_HPP
21 #define BISAI_HPP
22 
23 #include <mutex>
24 
25 #include <opm/simulators/linalg/bda/opencl/opencl.hpp>
26 #include <opm/simulators/linalg/bda/opencl/BILU0.hpp>
27 #include <opm/simulators/linalg/bda/opencl/Preconditioner.hpp>
28 
29 namespace Opm
30 {
31 namespace Accelerator
32 {
33 
34 class BlockedMatrix;
35 
38 template <unsigned int block_size>
39 class BISAI : public Preconditioner<block_size>
40 {
42 
43  using Base::N;
44  using Base::Nb;
45  using Base::nnz;
46  using Base::nnzb;
47  using Base::verbosity;
48  using Base::context;
49  using Base::queue;
50  using Base::events;
51  using Base::err;
52 
53 private:
54  std::once_flag initialize;
55 
56  std::vector<int> colPointers;
57  std::vector<int> rowIndices;
58  std::vector<int> diagIndex;
59  std::vector<int> csrToCscOffsetMap;
60  std::vector<double> invLvals;
61  std::vector<double> invUvals;
62 
63  cl::Buffer d_colPointers;
64  cl::Buffer d_rowIndices;
65  cl::Buffer d_csrToCscOffsetMap;
66  cl::Buffer d_diagIndex;
67  cl::Buffer d_LUvals;
68  cl::Buffer d_invDiagVals;
69  cl::Buffer d_invLvals;
70  cl::Buffer d_invUvals;
71  cl::Buffer d_invL_x;
72 
73  ILUReorder opencl_ilu_reorder;
74  std::unique_ptr<BILU0<block_size> > bilu0;
75 
77  typedef struct{
80  std::vector<int> subsystemPointers;
85  std::vector<int> nzIndices;
88  std::vector<int> knownRhsIndices;
90  std::vector<int> unknownRhsIndices;
91  } subsystemStructure;
92 
94  typedef struct{
95  cl::Buffer subsystemPointers;
96  cl::Buffer nzIndices;
97  cl::Buffer knownRhsIndices;
98  cl::Buffer unknownRhsIndices;
99  } subsystemStructureGPU;
100 
101  subsystemStructure lower, upper;
102  subsystemStructureGPU d_lower, d_upper;
103 
106  void buildLowerSubsystemsStructures();
107 
110  void buildUpperSubsystemsStructures();
111 
112 public:
113  BISAI(ILUReorder opencl_ilu_reorder, int verbosity);
114 
115  // set own Opencl variables, but also that of the bilu0 preconditioner
116  void setOpencl(std::shared_ptr<cl::Context>& context, std::shared_ptr<cl::CommandQueue>& queue) override;
117 
118  // analysis, find reordering if specified
119  bool analyze_matrix(BlockedMatrix *mat) override;
120 
121  // ilu_decomposition
122  bool create_preconditioner(BlockedMatrix *mat) override;
123 
124  // apply preconditioner, x = prec(y)
125  void apply(const cl::Buffer& y, cl::Buffer& x) override;
126 
127  int* getToOrder() override
128  {
129  return bilu0->getToOrder();
130  }
131 
132  int* getFromOrder() override
133  {
134  return bilu0->getFromOrder();
135  }
136 
137  BlockedMatrix* getRMat() override
138  {
139  return bilu0->getRMat();
140  }
141 };
142 
146 std::vector<int> buildCsrToCscOffsetMap(std::vector<int> colPointers, std::vector<int> rowIndices);
147 
148 } // namespace Accelerator
149 } // namespace Opm
150 
151 #endif
This class implements a Blocked version of the Incomplete Sparse Approximate Inverse (ISAI) precondit...
Definition: BISAI.hpp:40
This struct resembles a blocked csr matrix, like Dune::BCRSMatrix.
Definition: BlockedMatrix.hpp:37
Definition: Preconditioner.hpp:36
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27