My Project
Loading...
Searching...
No Matches
udq.hpp
1/*
2 Copyright 2021 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 RST_UDQ
21#define RST_UDQ
22
23#include <opm/input/eclipse/EclipseState/Phase.hpp>
24
26
27#include <opm/input/eclipse/Schedule/UDQ/UDQEnums.hpp>
28
29#include <cstddef>
30#include <iterator>
31#include <optional>
32#include <string>
33#include <string_view>
34#include <utility>
35#include <variant>
36#include <vector>
37
38namespace Opm::RestartIO {
39
62class RstUDQ
63{
64public:
66 class ValueRange
67 {
68 public:
71 class Iterator
72 {
73 public:
75 using iterator_category = std::forward_iterator_tag;
76
82 using value_type = std::pair<std::size_t, double>;
83
85 using difference_type = int;
86
88 using pointer = const value_type*;
89
91 using reference = const value_type&;
92
96 Iterator& operator++()
97 {
98 ++this->ix_;
99
100 return *this;
101 }
102
106 Iterator operator++(int)
107 {
108 auto iter = *this;
109
110 ++(*this);
111
112 return iter;
113 }
114
118 reference operator*()
119 {
120 this->setValue();
121 return this->deref_value_;
122 }
123
127 pointer operator->()
128 {
129 this->setValue();
130 return &this->deref_value_;
131 }
132
138 bool operator==(const Iterator& that) const
139 {
140 return (this->ix_ == that.ix_)
141 && (this->i_ == that.i_)
142 && (this->value_ == that.value_);
143 }
144
150 bool operator!=(const Iterator& that) const
151 {
152 return ! (*this == that);
153 }
154
155 friend class ValueRange;
156
157 private:
163 Iterator(std::size_t ix,
164 const int* i,
165 const std::variant<double, const double*>& value)
166 : ix_{ix}, i_{i}, value_{value}
167 {}
168
170 std::size_t ix_;
171 const int* i_;
172 std::variant<double, const double*> value_;
173
174 value_type deref_value_{};
175
176 void setValue();
177 };
178
180 Iterator begin() const
181 { return this->makeIterator(this->begin_); }
182
184 Iterator end() const
185 { return this->makeIterator(this->end_); }
186
187 friend class RstUDQ;
188
189 private:
203 ValueRange(std::size_t begin_arg,
204 std::size_t end_arg,
205 const int* i,
206 const double value)
207 : begin_{begin_arg}, end_{end_arg}, i_{i}, value_{value}
208 {}
209
223 ValueRange(std::size_t begin_arg,
224 std::size_t end_arg,
225 const int* i,
226 const double* value)
227 : begin_{begin_arg}, end_{end_arg}, i_{i}, value_{value}
228 {}
229
231 std::size_t begin_{};
232
234 std::size_t end_{};
235
237 const int* i_{};
238
243 std::variant<double, const double*> value_{};
244
250 Iterator makeIterator(const std::size_t ix) const
251 {
252 return { ix, this->i_, this->value_ };
253 }
254 };
255
259 std::string name{};
260
262 std::string unit{};
263
268 UDQVarType category{UDQVarType::NONE};
269
283 RstUDQ(const std::string& name_arg,
284 const std::string& unit_arg,
285 const std::string& define_arg,
286 UDQUpdate status_arg);
287
295 RstUDQ(const std::string& name_arg,
296 const std::string& unit_arg);
297
302 void prepareValues();
303
318 void addValue(const int entity, const int subEntity, const double value);
319
325 void commitValues();
326
335 void assignScalarValue(const double value);
336
340 void addEntityName(std::string_view wgname);
341
349 double scalarValue() const;
350
352 std::size_t numEntities() const;
353
360 ValueRange operator[](const std::size_t i) const;
361
366 const std::vector<std::string>& entityNames() const
367 {
368 return this->wgnames_;
369 }
370
379 const std::vector<int>& nameIndex() const;
380
385 std::string_view definingExpression() const;
386
388 UDQUpdate currentUpdateStatus() const;
389
394 bool isDefine() const
395 {
396 return this->definition_.has_value();
397 }
398
401 bool isScalar() const
402 {
403 return std::holds_alternative<double>(this->sa_);
404 }
405
406private:
412 using Graph = utility::CSRGraphFromCoordinates<int, true, true>;
413
415 struct Definition
416 {
424 Definition(const std::string& expression_arg,
425 const UDQUpdate status_arg)
426 : expression { expression_arg }
427 , status { status_arg }
428 {}
429
431 std::string expression{};
432
434 UDQUpdate status{};
435 };
436
448 Graph entityMap_{};
449
456 std::variant<std::monostate, double, std::vector<double>> sa_;
457
461 std::vector<std::string> wgnames_{};
462
466 std::optional<int> maxEntityIdx_{};
467
469 mutable std::optional<std::vector<int>> wgNameIdx_{};
470
475 std::optional<Definition> definition_{};
476
478 bool isUDQSet() const
479 {
480 return std::holds_alternative<std::vector<double>>(this->sa_);
481 }
482
489 ValueRange scalarRange(const std::size_t i) const;
490
498 ValueRange udqSetRange(const std::size_t i) const;
499
500 void ensureValidNameIndex() const;
501};
502
504struct RstUDQActive
505{
507 struct RstRecord
508 {
526 RstRecord(UDAControl c,
527 std::size_t i,
528 std::size_t numIuap,
529 std::size_t u1,
530 std::size_t u2);
531
533 UDAControl control;
534
536 std::size_t input_index;
537
541 std::size_t use_count;
542
544 std::size_t wg_offset;
545
547 std::size_t num_wg_elems;
548 };
549
553 RstUDQActive() = default;
554
566 RstUDQActive(const std::vector<int>& iuad,
567 const std::vector<int>& iuap,
568 const std::vector<int>& igph);
569
571 std::vector<int> wg_index{};
572
574 std::vector<RstRecord> iuad{};
575
577 std::vector<Phase> ig_phase{};
578};
579
580} // namespace Opm::RestartIO
581
582#endif // RST_UDQ
Facility for converting collection of region ID pairs into a sparse (CSR) adjacency matrix representa...