My Project
BlockedMatrix.hpp
1 /*
2  Copyright 2019 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 OPM_BLOCKED_MATRIX_HPP
21 #define OPM_BLOCKED_MATRIX_HPP
22 
23 #if HAVE_FPGA
24 #include <vector>
25 #include <opm/simulators/linalg/bda/Matrix.hpp>
26 #endif
27 
28 
29 namespace Opm
30 {
31 namespace Accelerator
32 {
33 
37 {
38 
39 public:
40 
45  BlockedMatrix(int Nb_, int nnzbs_, unsigned int block_size_)
46  : nnzValues(new double[nnzbs_*block_size_*block_size_]),
47  colIndices(new int[nnzbs_*block_size_*block_size_]),
48  rowPointers(new int[Nb_+1]),
49  Nb(Nb_),
50  nnzbs(nnzbs_),
51  block_size(block_size_),
52  deleteNnzs(true),
53  deleteSparsity(true)
54  {}
55 
59  : nnzValues(new double[M.nnzbs*M.block_size*M.block_size]),
60  colIndices(M.colIndices),
61  rowPointers(M.rowPointers),
62  Nb(M.Nb),
63  nnzbs(M.nnzbs),
64  block_size(M.block_size),
65  deleteNnzs(true),
66  deleteSparsity(false)
67  {}
68 
76  BlockedMatrix(int Nb_, int nnzbs_, unsigned int block_size_, double *nnzValues_, int *colIndices_, int *rowPointers_)
77  : nnzValues(nnzValues_),
78  colIndices(colIndices_),
79  rowPointers(rowPointers_),
80  Nb(Nb_),
81  nnzbs(nnzbs_),
82  block_size(block_size_),
83  deleteNnzs(false),
84  deleteSparsity(false)
85  {}
86 
87  ~BlockedMatrix(){
88  if (deleteNnzs) {
89  delete[] nnzValues;
90  }
91  if (deleteSparsity) {
92  delete[] colIndices;
93  delete[] rowPointers;
94  }
95  }
96 
97 #if HAVE_FPGA
98  constexpr static double nnzThreshold = 1e-80; // for unblocking, a nonzero must be bigger than this threshold to be considered a nonzero
99 
100  int countUnblockedNnzs();
101 
102  void unblock(Matrix *mat, bool isUMatrix);
103 
106  int toRDF(int numColors, int *nodesPerColor, bool isUMatrix,
107  std::vector<std::vector<int> >& colIndicesInColor, int nnzsPerRowLimit, int *nnzValsSizes,
108  std::vector<std::vector<double> >& nnzValues, short int *colIndices, unsigned char *NROffsets, int *colorSizes, int *valSize);
109 
112  int findPartitionColumns(int numColors, int *nodesPerColor,
113  int rowsPerColorLimit, int columnsPerColorLimit,
114  std::vector<std::vector<int> >& colIndicesInColor, int *PIndicesAddr, int *colorSizes,
115  std::vector<std::vector<int> >& LColIndicesInColor, int *LPIndicesAddr, int *LColorSizes,
116  std::vector<std::vector<int> >& UColIndicesInColor, int *UPIndicesAddr, int *UColorSizes);
117 #endif
118 
119  double *nnzValues;
120  int *colIndices;
121  int *rowPointers;
122  int Nb;
123  int nnzbs;
124  unsigned int block_size;
125  bool deleteNnzs;
126  bool deleteSparsity;
127 };
128 
129 
136 void sortBlockedRow(int *colIndices, double *data, int left, int right, unsigned block_size);
137 
144 void blockMultSub(double *a, double *b, double *c, unsigned int block_size);
145 
152 void blockMult(double *mat1, double *mat2, double *resMat, unsigned int block_size);
153 
154 
155 #if HAVE_FPGA
162 void blockSub(double *mat1, double *mat2, double *resMat, unsigned int block_size);
163 
173 void blockVectMult(double *mat, double *vect, double scale, double *resVect, bool resetRes, unsigned int block_size);
174 
185 void blockedDiagtoRDF(double *blockedDiagVals, int rowSize, int numColors, std::vector<int>& rowsPerColor, double *RDFDiag);
186 #endif
187 
188 } // namespace Accelerator
189 } // namespace Opm
190 
191 #endif
This struct resembles a blocked csr matrix, like Dune::BCRSMatrix.
Definition: BlockedMatrix.hpp:37
BlockedMatrix(const BlockedMatrix &M)
Allocate BlockedMatrix, but copy sparsity pattern instead of allocating new memory.
Definition: BlockedMatrix.hpp:58
BlockedMatrix(int Nb_, int nnzbs_, unsigned int block_size_, double *nnzValues_, int *colIndices_, int *rowPointers_)
Allocate BlockedMatrix, but let data arrays point to existing arrays.
Definition: BlockedMatrix.hpp:76
BlockedMatrix(int Nb_, int nnzbs_, unsigned int block_size_)
Allocate BlockedMatrix and data arrays with given sizes.
Definition: BlockedMatrix.hpp:45
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27