My Project
Loading...
Searching...
No Matches
CSRGraphFromCoordinates.hpp
Go to the documentation of this file.
1/*
2 Copyright 2016 SINTEF ICT, Applied Mathematics.
3 Copyright 2016 Statoil ASA.
4 Copyright 2022 Equinor ASA
5
6 This file is part of the Open Porous Media Project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#ifndef OPM_UTILITY_CSRGRAPHFROMCOORDINATES_HPP
23#define OPM_UTILITY_CSRGRAPHFROMCOORDINATES_HPP
24
25#include <cstddef>
26#include <optional>
27#include <type_traits>
28#include <utility>
29#include <vector>
30
36
37namespace Opm { namespace utility {
38
54 template <typename VertexID = int, bool TrackCompressedIdx = false, bool PermitSelfConnections = false>
56 {
57 private:
58 using BaseVertexID = std::remove_cv_t<VertexID>;
59
60 static_assert(std::is_integral_v<BaseVertexID>,
61 "The VertexID must be an integral type");
62
63 public:
65 using Neighbours = std::vector<BaseVertexID>;
66
68 using Offset = typename Neighbours::size_type;
69
71 using Start = std::vector<Offset>;
72
74 void clear();
75
85 void addConnection(VertexID v1, VertexID v2);
86
104 void compress(Offset maxNumVertices,
105 bool expandExistingIdxMap = false);
106
110 Offset numVertices() const;
111
114 Offset numEdges() const;
115
117 const Start& startPointers() const
118 {
119 return this->csr_.startPointers();
120 }
121
125 {
126 return this->csr_.columnIndices();
127 }
128
136 template <typename Ret = const Start&>
137 std::enable_if_t<TrackCompressedIdx, Ret> compressedIndexMap() const
138 {
139 return this->csr_.compressedIndexMap();
140 }
141
142 // MessageBufferType API should be similar to Dune::MessageBufferIF
143 template <class MessageBufferType>
144 void write(MessageBufferType& buffer) const
145 {
146 this->csr_.write(buffer);
147 }
148
149 // MessageBufferType API should be similar to Dune::MessageBufferIF
150 template <class MessageBufferType>
151 void read(MessageBufferType& buffer)
152 {
153 auto other = CSR{};
154 other.read(buffer);
155
156 this->uncompressed_
157 .add(other.maxRowID(),
158 other.maxColID(),
159 other.coordinateFormatRowIndices(),
160 other.columnIndices());
161 }
162
163 private:
166 class Connections
167 {
168 public:
174 void add(VertexID v1, VertexID v2);
175
189 void add(VertexID maxRowIdx,
190 VertexID maxColIdx,
191 const Neighbours& rows,
192 const Neighbours& cols);
193
195 void clear();
196
200 bool empty() const;
201
204 bool isValid() const;
205
207 std::optional<BaseVertexID> maxRow() const;
208
210 std::optional<BaseVertexID> maxCol() const;
211
213 typename Neighbours::size_type numContributions() const;
214
216 const Neighbours& rowIndices() const;
217
219 const Neighbours& columnIndices() const;
220
221 private:
223 Neighbours i_{};
224
226 Neighbours j_{};
227
229 std::optional<VertexID> max_i_{};
230
232 std::optional<VertexID> max_j_{};
233 };
234
239 class CSR
240 {
241 public:
257 void merge(const Connections& conns,
258 const Offset maxNumVertices,
259 const bool expandExistingIdxMap);
260
262 Offset numRows() const;
263
265 BaseVertexID maxRowID() const;
266
268 BaseVertexID maxColID() const;
269
271 const Start& startPointers() const;
272
275 const Neighbours& columnIndices() const;
276
279 Neighbours coordinateFormatRowIndices() const;
280
281 template <typename Ret = const Start&>
282 std::enable_if_t<TrackCompressedIdx, Ret> compressedIndexMap() const
283 {
284 return this->compressedIdx_;
285 }
286
287 // MessageBufferType API should be similar to Dune::MessageBufferIF
288 template <class MessageBufferType>
289 void write(MessageBufferType& buffer) const
290 {
291 this->writeVector(this->ia_, buffer);
292 this->writeVector(this->ja_, buffer);
293
294 if constexpr (TrackCompressedIdx) {
295 this->writeVector(this->compressedIdx_, buffer);
296 }
297
298 buffer.write(this->numRows_);
299 buffer.write(this->numCols_);
300 }
301
302 // MessageBufferType API should be similar to Dune::MessageBufferIF
303 template <class MessageBufferType>
304 void read(MessageBufferType& buffer)
305 {
306 this->readVector(buffer, this->ia_);
307 this->readVector(buffer, this->ja_);
308
309 if constexpr (TrackCompressedIdx) {
310 this->readVector(buffer, this->compressedIdx_);
311 }
312
313 buffer.read(this->numRows_);
314 buffer.read(this->numCols_);
315 }
316
318 void clear();
319
320 private:
321 struct EmptyPlaceHolder {};
322
324 Start ia_{};
325
328 Neighbours ja_{};
329
335 std::conditional_t<TrackCompressedIdx, Start, EmptyPlaceHolder> compressedIdx_{};
336
338 BaseVertexID numRows_{};
339
342 BaseVertexID numCols_{};
343
344 // ---------------------------------------------------------
345 // Implementation of read()/write()
346 // ---------------------------------------------------------
347
348 template <typename T, class A, class MessageBufferType>
349 void writeVector(const std::vector<T,A>& vec,
350 MessageBufferType& buffer) const
351 {
352 const auto n = vec.size();
353 buffer.write(n);
354
355 for (const auto& x : vec) {
356 buffer.write(x);
357 }
358 }
359
360 template <class MessageBufferType, typename T, class A>
361 void readVector(MessageBufferType& buffer,
362 std::vector<T,A>& vec)
363 {
364 auto n = 0 * vec.size();
365 buffer.read(n);
366
367 vec.resize(n);
368
369 for (auto& x : vec) {
370 buffer.read(x);
371 }
372 }
373
374 // ---------------------------------------------------------
375 // Implementation of merge()
376 // ---------------------------------------------------------
377
406 void assemble(const Neighbours& rows,
407 const Neighbours& cols,
408 BaseVertexID maxRowID,
409 BaseVertexID maxColID,
410 bool expandExistingIdxMap);
411
422 void compress(const Offset maxNumVertices);
423
430 void sortColumnIndicesPerRow();
431
441 void condenseDuplicates();
442
443 // ---------------------------------------------------------
444 // Implementation of assemble()
445 // ---------------------------------------------------------
446
460 void preparePushbackRowGrouping(const int numRows,
461 const Neighbours& rowIdx);
462
474 void groupAndTrackColumnIndicesByRow(const Neighbours& rowIdx,
475 const Neighbours& colIdx);
476
477 // ---------------------------------------------------------
478 // General utilities
479 // ---------------------------------------------------------
480
485 void transpose();
486
503 void condenseAndTrackUniqueColumnsForSingleRow(typename Neighbours::const_iterator begin,
504 typename Neighbours::const_iterator end);
505
515 void remapCompressedIndex(Start&& compressedIdx,
516 std::optional<typename Start::size_type> numOrigNNZ = std::nullopt);
517 };
518
521 Connections uncompressed_;
522
524 CSR csr_;
525 };
526
527}} // namespace Opm::utility
528
529// Actual implementation of member functions in _impl.hpp file.
530#include <opm/common/utility/CSRGraphFromCoordinates_impl.hpp>
531
532#endif // OPM_UTILITY_CSRGRAPHFROMCOORDINATES_HPP
Form CSR adjacency matrix representation of unstructured graphs.
Definition CSRGraphFromCoordinates.hpp:56
const Start & startPointers() const
Read-only access to compressed structure's start pointers.
Definition CSRGraphFromCoordinates.hpp:117
Offset numEdges() const
Retrieve number of edges (non-zero matrix elements) in input graph.
Definition CSRGraphFromCoordinates_impl.hpp:605
std::enable_if_t< TrackCompressedIdx, Ret > compressedIndexMap() const
Read-only access to mapping from input order vertex pairs to compressed structure's edge index (locat...
Definition CSRGraphFromCoordinates.hpp:137
void clear()
Clear all internal buffers, but preserve allocated capacity.
Definition CSRGraphFromCoordinates_impl.hpp:551
void addConnection(VertexID v1, VertexID v2)
Add flow rate connection between regions.
Definition CSRGraphFromCoordinates_impl.hpp:560
const Neighbours & columnIndices() const
Read-only access to compressed structure's column indices, ascendingly sorted per row.
Definition CSRGraphFromCoordinates.hpp:124
typename Neighbours::size_type Offset
Offset into neighbour array.
Definition CSRGraphFromCoordinates.hpp:68
Offset numVertices() const
Retrieve number of rows (source entities) in input graph.
Definition CSRGraphFromCoordinates_impl.hpp:598
void compress(Offset maxNumVertices, bool expandExistingIdxMap=false)
Form CSR adjacency matrix representation of input graph from connections established in previous call...
Definition CSRGraphFromCoordinates_impl.hpp:583
std::vector< BaseVertexID > Neighbours
Representation of neighbouring regions.
Definition CSRGraphFromCoordinates.hpp:65
std::vector< Offset > Start
CSR start pointers.
Definition CSRGraphFromCoordinates.hpp:71
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30