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 <unordered_map>
29#include <utility>
30#include <vector>
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
89
90
108 void compress(Offset maxNumVertices,
109 bool expandExistingIdxMap = false);
110
114 Offset numVertices() const;
115
118 Offset numEdges() const;
119
121 const Start& startPointers() const
122 {
123 return this->csr_.startPointers();
124 }
125
129 {
130 return this->csr_.columnIndices();
131 }
132
140 template <typename Ret = const Start&>
141 std::enable_if_t<TrackCompressedIdx, Ret> compressedIndexMap() const
142 {
143 return this->csr_.compressedIndexMap();
144 }
145
146 // MessageBufferType API should be similar to Dune::MessageBufferIF
147 template <class MessageBufferType>
148 void write(MessageBufferType& buffer) const
149 {
150 this->csr_.write(buffer);
151 }
152
153 // MessageBufferType API should be similar to Dune::MessageBufferIF
154 template <class MessageBufferType>
155 void read(MessageBufferType& buffer)
156 {
157 auto other = CSR{};
158 other.read(buffer);
159
160 this->uncompressed_
161 .add(other.maxRowID(),
162 other.maxColID(),
163 other.coordinateFormatRowIndices(),
164 other.columnIndices());
165 }
166
171 void addVertexGroup(const std::vector<VertexID>& vertices);
172
178 VertexID getFinalVertexID(VertexID originalVertexID) const;
179
180 private:
183 class Connections
184 {
185 public:
191 void add(VertexID v1, VertexID v2);
192
206 void add(VertexID maxRowIdx,
207 VertexID maxColIdx,
208 const Neighbours& rows,
209 const Neighbours& cols);
210
212 void clear();
213
217 bool empty() const;
218
221 bool isValid() const;
222
224 std::optional<BaseVertexID> maxRow() const;
225
227 std::optional<BaseVertexID> maxCol() const;
228
230 typename Neighbours::size_type numContributions() const;
231
233 const Neighbours& rowIndices() const;
234
236 const Neighbours& columnIndices() const;
237
239 VertexID findMergedVertexID(VertexID v, const std::unordered_map<VertexID, VertexID>& vertex_merges) const;
240
242 std::unordered_map<VertexID, VertexID> applyVertexMerges(const std::unordered_map<VertexID, VertexID>& vertex_merges);
243
244 private:
246 Neighbours i_{};
247
249 Neighbours j_{};
250
252 std::optional<VertexID> max_i_{};
253
255 std::optional<VertexID> max_j_{};
256 };
257
262 class CSR
263 {
264 public:
280 void merge(const Connections& conns,
281 const Offset maxNumVertices,
282 const bool expandExistingIdxMap);
283
285 Offset numRows() const;
286
288 BaseVertexID maxRowID() const;
289
291 BaseVertexID maxColID() const;
292
294 const Start& startPointers() const;
295
298 const Neighbours& columnIndices() const;
299
302 Neighbours coordinateFormatRowIndices() const;
303
304 template <typename Ret = const Start&>
305 std::enable_if_t<TrackCompressedIdx, Ret> compressedIndexMap() const
306 {
307 return this->compressedIdx_;
308 }
309
310 // MessageBufferType API should be similar to Dune::MessageBufferIF
311 template <class MessageBufferType>
312 void write(MessageBufferType& buffer) const
313 {
314 this->writeVector(this->ia_, buffer);
315 this->writeVector(this->ja_, buffer);
316
317 if constexpr (TrackCompressedIdx) {
318 this->writeVector(this->compressedIdx_, buffer);
319 }
320
321 buffer.write(this->numRows_);
322 buffer.write(this->numCols_);
323 }
324
325 // MessageBufferType API should be similar to Dune::MessageBufferIF
326 template <class MessageBufferType>
327 void read(MessageBufferType& buffer)
328 {
329 this->readVector(buffer, this->ia_);
330 this->readVector(buffer, this->ja_);
331
332 if constexpr (TrackCompressedIdx) {
333 this->readVector(buffer, this->compressedIdx_);
334 }
335
336 buffer.read(this->numRows_);
337 buffer.read(this->numCols_);
338 }
339
341 void clear();
342
343 private:
344 struct EmptyPlaceHolder {};
345
347 Start ia_{};
348
351 Neighbours ja_{};
352
358 std::conditional_t<TrackCompressedIdx, Start, EmptyPlaceHolder> compressedIdx_{};
359
361 BaseVertexID numRows_{};
362
365 BaseVertexID numCols_{};
366
367 // ---------------------------------------------------------
368 // Implementation of read()/write()
369 // ---------------------------------------------------------
370
371 template <typename T, class A, class MessageBufferType>
372 void writeVector(const std::vector<T,A>& vec,
373 MessageBufferType& buffer) const
374 {
375 const auto n = vec.size();
376 buffer.write(n);
377
378 for (const auto& x : vec) {
379 buffer.write(x);
380 }
381 }
382
383 template <class MessageBufferType, typename T, class A>
384 void readVector(MessageBufferType& buffer,
385 std::vector<T,A>& vec)
386 {
387 auto n = 0 * vec.size();
388 buffer.read(n);
389
390 vec.resize(n);
391
392 for (auto& x : vec) {
393 buffer.read(x);
394 }
395 }
396
397 // ---------------------------------------------------------
398 // Implementation of merge()
399 // ---------------------------------------------------------
400
429 void assemble(const Neighbours& rows,
430 const Neighbours& cols,
431 BaseVertexID maxRowID,
432 BaseVertexID maxColID,
433 bool expandExistingIdxMap);
434
445 void compress(const Offset maxNumVertices);
446
453 void sortColumnIndicesPerRow();
454
464 void condenseDuplicates();
465
466 // ---------------------------------------------------------
467 // Implementation of assemble()
468 // ---------------------------------------------------------
469
483 void preparePushbackRowGrouping(const int numRows,
484 const Neighbours& rowIdx);
485
497 void groupAndTrackColumnIndicesByRow(const Neighbours& rowIdx,
498 const Neighbours& colIdx);
499
500 // ---------------------------------------------------------
501 // General utilities
502 // ---------------------------------------------------------
503
508 void transpose();
509
526 void condenseAndTrackUniqueColumnsForSingleRow(typename Neighbours::const_iterator begin,
527 typename Neighbours::const_iterator end);
528
538 void remapCompressedIndex(Start&& compressedIdx,
539 std::optional<typename Start::size_type> numOrigNNZ = std::nullopt);
540
541 };
542
545 Connections uncompressed_;
546
548 CSR csr_;
549
551 std::unordered_map<VertexID, VertexID> parent_{};
552
554 std::unordered_map<VertexID, VertexID> vertex_mapping_{};
555
557 VertexID find(VertexID v);
558
560 void unionSets(VertexID a, VertexID b);
561 };
562
563}} // namespace Opm::utility
564
565// Actual implementation of member functions in _impl.hpp file.
566#include <opm/common/utility/CSRGraphFromCoordinates_impl.hpp>
567
568#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:121
Offset numEdges() const
Retrieve number of edges (non-zero matrix elements) in input graph.
Definition CSRGraphFromCoordinates_impl.hpp:810
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:141
void clear()
Clear all internal buffers, but preserve allocated capacity.
Definition CSRGraphFromCoordinates_impl.hpp:644
void addConnection(VertexID v1, VertexID v2)
Add flow rate connection between regions.
Definition CSRGraphFromCoordinates_impl.hpp:653
const Neighbours & columnIndices() const
Read-only access to compressed structure's column indices, ascendingly sorted per row.
Definition CSRGraphFromCoordinates.hpp:128
VertexID getFinalVertexID(VertexID originalVertexID) const
Get the final vertex ID after all merges and renumbering for a given original vertex ID.
Definition CSRGraphFromCoordinates_impl.hpp:791
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:803
void addVertexGroup(const std::vector< VertexID > &vertices)
Add a group of vertices that should be merged together.
Definition CSRGraphFromCoordinates_impl.hpp:676
Offset applyVertexMerges()
Apply vertex merges to all vertex groups.
Definition CSRGraphFromCoordinates_impl.hpp:741
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:770
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