My Project
UDQConfig.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 UDQINPUT_HPP_
21#define UDQINPUT_HPP_
22
23#include <opm/input/eclipse/Schedule/UDQ/UDQAssign.hpp>
24#include <opm/input/eclipse/Schedule/UDQ/UDQDefine.hpp>
25#include <opm/input/eclipse/Schedule/UDQ/UDQEnums.hpp>
26#include <opm/input/eclipse/Schedule/UDQ/UDQFunctionTable.hpp>
27#include <opm/input/eclipse/Schedule/UDQ/UDQInput.hpp>
28#include <opm/input/eclipse/Schedule/UDQ/UDQParams.hpp>
29
30#include <opm/input/eclipse/EclipseState/Util/OrderedMap.hpp>
31#include <opm/input/eclipse/EclipseState/Util/IOrderSet.hpp>
32
33#include <cstddef>
34#include <map>
35#include <string>
36#include <unordered_map>
37#include <unordered_set>
38#include <vector>
39
40namespace Opm {
41
42 class DeckRecord;
43 class SummaryState;
44 class UDQState;
45 class KeywordLocation;
46 class WellMatcher;
47
48} // namespace Opm
49
50namespace Opm { namespace RestartIO {
51 struct RstState;
52}} // namespace Opm::RestartIO
53
54namespace Opm {
55
57 {
58 public:
59 UDQConfig() = default;
60 explicit UDQConfig(const UDQParams& params);
61 UDQConfig(const UDQParams& params, const RestartIO::RstState& rst_state);
62
63 static UDQConfig serializationTestObject();
64
65 const std::string& unit(const std::string& key) const;
66 bool has_unit(const std::string& keyword) const;
67 bool has_keyword(const std::string& keyword) const;
68 void add_record(const DeckRecord& record, const KeywordLocation& location, std::size_t report_step);
69
70 void add_unit(const std::string& keyword, const std::string& unit);
71 void add_update(const std::string& keyword, std::size_t report_step, const KeywordLocation& location, const std::vector<std::string>& data);
72 void add_assign(const std::string& quantity, const std::vector<std::string>& selector, double value, std::size_t report_step);
73 void add_assign(const std::string& quantity, const std::unordered_set<std::string>& selector, double value, std::size_t report_step);
74 void add_define(const std::string& quantity, const KeywordLocation& location, const std::vector<std::string>& expression, std::size_t report_step);
75
76 void eval_assign(std::size_t report_step, const WellMatcher& wm, SummaryState& st, UDQState& udq_state) const;
77 void eval(std::size_t report_step, const WellMatcher& wm, SummaryState& st, UDQState& udq_state) const;
78 const UDQDefine& define(const std::string& key) const;
79 const UDQAssign& assign(const std::string& key) const;
80 std::vector<UDQDefine> definitions() const;
81 std::vector<UDQDefine> definitions(UDQVarType var_type) const;
82 std::vector<UDQInput> input() const;
83
84 // The size() method will return the number of active DEFINE and ASSIGN
85 // statements; this will correspond to the length of the vector returned
86 // from input().
87 size_t size() const;
88
89 UDQInput operator[](const std::string& keyword) const;
90 UDQInput operator[](std::size_t insert_index) const;
91
92 std::vector<UDQAssign> assignments() const;
93 std::vector<UDQAssign> assignments(UDQVarType var_type) const;
94 const UDQParams& params() const;
95 const UDQFunctionTable& function_table() const;
96
97 bool operator==(const UDQConfig& config) const;
98 void required_summary(std::unordered_set<std::string>& summary_keys) const;
99
100 template<class Serializer>
101 void serializeOp(Serializer& serializer)
102 {
103 serializer(udq_params);
104 serializer(m_definitions);
105 serializer(m_assignments);
106 serializer(units);
107 serializer(input_index);
108 serializer(type_count);
109 // The UDQFunction table is constant up to udq_params.
110 // So we can just construct a new instance here.
111 if (!serializer.isSerializing())
112 udqft = UDQFunctionTable(udq_params);
113 }
114
115 private:
116 void add_node(const std::string& quantity, UDQAction action);
117 UDQAction action_type(const std::string& udq_key) const;
118 void eval_assign(std::size_t report_step, SummaryState& st, UDQState& udq_state, UDQContext& context) const;
119 void eval_define(std::size_t report_step, UDQState& udq_state, UDQContext& context) const;
120
121 UDQParams udq_params;
122 UDQFunctionTable udqft;
123
124 // The choices of data structures are strongly motivated by the
125 // constraints imposed by the ECLIPSE formatted restart files; for
126 // writing restart files it is essential to keep meticulous control
127 // over the ordering of the keywords. In this class the ordering is
128 // mainly maintained by the input_index map which keeps track of the
129 // insert order of each keyword, and whether the keyword is
130 // currently DEFINE'ed or ASSIGN'ed.
131 std::unordered_map<std::string, UDQDefine> m_definitions;
132 std::unordered_map<std::string, UDQAssign> m_assignments;
133 std::unordered_map<std::string, std::string> units;
134
135 IOrderSet<std::string> define_order;
136 OrderedMap<UDQIndex> input_index;
137 std::map<UDQVarType, std::size_t> type_count;
138 };
139}
140
141
142
143#endif
Definition: DeckRecord.hpp:32
Definition: KeywordLocation.hpp:27
A map with iteration in the order of insertion.
Definition: OrderedMap.hpp:114
Class for (de-)serializing.
Definition: Serializer.hpp:84
bool isSerializing() const
Returns true if we are currently doing a serialization operation.
Definition: Serializer.hpp:183
Definition: SummaryState.hpp:68
Definition: UDQAssign.hpp:35
Definition: UDQConfig.hpp:57
Definition: UDQContext.hpp:40
Definition: UDQDefine.hpp:51
Definition: UDQFunctionTable.hpp:32
Definition: UDQInput.hpp:84
Definition: UDQParams.hpp:31
Definition: UDQState.hpp:38
Definition: WellMatcher.hpp:32
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30
Definition: state.hpp:54