My Project
FPGABILU0.hpp
1 /*
2  Copyright 2020 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 FPGA_BILU0_HEADER_INCLUDED
21 #define FPGA_BILU0_HEADER_INCLUDED
22 
23 #include <vector>
24 
25 #include <opm/simulators/linalg/bda/ILUReorder.hpp>
26 #include <opm/simulators/linalg/bda/BlockedMatrix.hpp>
27 
28 namespace Opm
29 {
30 namespace Accelerator
31 {
32 
33 /*
34  * This class implements a Blocked ILU0 preconditioner, with output data
35  * specifically formatted for the FPGA.
36  * The decomposition and reorders of the rows of the matrix are done on CPU.
37  */
38 
39 template <unsigned int block_size>
40 class FPGABILU0
41 {
42 
43 private:
44  int N; // number of rows of the matrix
45  int Nb; // number of blockrows of the matrix
46  int nnz; // number of nonzeroes of the matrix (scalar)
47  int nnzbs; // number of blocks of the matrix
48  std::unique_ptr<BlockedMatrix> LMat = nullptr, UMat = nullptr, LUMat = nullptr;
49  std::shared_ptr<BlockedMatrix> rMat = nullptr; // reordered mat
50  double *invDiagVals = nullptr;
51  std::vector<int> diagIndex;
52  std::vector<int> toOrder, fromOrder;
53  std::vector<int> rowsPerColor;
54  int numColors;
55  int verbosity;
56 
57  // sizes and arrays used during RDF generation
58  std::vector<std::vector<double> > nnzValues, LnnzValues, UnnzValues;
59  std::vector<short int> colIndices, LColIndices, UColIndices;
60  std::vector<unsigned char> NROffsets, LNROffsets, UNROffsets;
61  std::vector<int> PIndicesAddr, LPIndicesAddr, UPIndicesAddr;
62  std::vector<int> colorSizes, LColorSizes, UColorSizes;
63  std::vector<int> nnzValsSizes, LnnzValsSizes, UnnzValsSizes;
64  std::vector<std::vector<int> > colIndicesInColor, LColIndicesInColor, UColIndicesInColor;
65 
66  int rowSize, valSize;
67  int LRowSize, LValSize, LNumColors;
68  int URowSize, UValSize, UNumColors;
69  std::vector<double> blockDiag;
70  ILUReorder opencl_ilu_reorder;
71  bool level_scheduling = false, graph_coloring = false;
72  int numResultPointers = 21;
73  std::vector<void *> resultPointers;
74  int numResultSizes = 18;
75  std::vector<int> resultSizes;
76  int maxRowsPerColor, maxColsPerColor, maxNNZsPerRow, maxNumColors; // are set via the constructor
77 
78 public:
79 
80  FPGABILU0(ILUReorder opencl_ilu_reorder, int verbosity, int maxRowsPerColor, int maxColsPerColor, int maxNNZsPerRow, int maxNumColors);
81 
82  ~FPGABILU0();
83 
84  // analysis (optional)
85  bool init(BlockedMatrix *mat);
86 
87  // ilu_decomposition
88  bool create_preconditioner(BlockedMatrix *mat);
89 
90  int* getToOrder()
91  {
92  return toOrder.data();
93  }
94 
95  int* getFromOrder()
96  {
97  return fromOrder.data();
98  }
99 
100  BlockedMatrix* getRMat()
101  {
102  return rMat.get();
103  }
104 
105  void **getResultPointers()
106  {
107  return resultPointers.data();
108  }
109 
110  int *getResultSizes()
111  {
112  return resultSizes.data();
113  }
114 
115 };
116 
117 } // namespace Accelerator
118 } // namespace Opm
119 
120 #endif // FPGA_BILU0_HEADER_INCLUDED
This struct resembles a blocked csr matrix, like Dune::BCRSMatrix.
Definition: BlockedMatrix.hpp:37
Definition: FPGABILU0.hpp:41
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27