My Project
UDQState.hpp
1 /*
2  Copyright 2020 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 UDQSTATE_HPP_
21 #define UDQSTATE_HPP_
22 
23 #include <string>
24 #include <unordered_map>
25 #include <vector>
26 
27 #include <opm/input/eclipse/Schedule/UDQ/UDQSet.hpp>
28 
29 
30 namespace Opm {
31 
32 namespace RestartIO {
33  struct RstState;
34 }
35 
36 class UDQState {
37 public:
38  UDQState() = default;
39  UDQState(double undefined);
40 
41  bool has(const std::string& key) const;
42  void load_rst(const RestartIO::RstState& rst_state);
43 
44  bool has_well_var(const std::string& well, const std::string& key) const;
45  bool has_group_var(const std::string& group, const std::string& key) const;
46 
47  double get(const std::string& key) const;
48  double get_group_var(const std::string& well, const std::string& var) const;
49  double get_well_var(const std::string& well, const std::string& var) const;
50  void add_define(std::size_t report_step, const std::string& udq_key, const UDQSet& result);
51  void add_assign(std::size_t report_step, const std::string& udq_key, const UDQSet& result);
52  bool assign(std::size_t report_step, const std::string& udq_key) const;
53  bool define(const std::string& udq_key, std::pair<UDQUpdate, std::size_t> update_status) const;
54  double undefined_value() const;
55 
56  std::vector<char> serialize() const;
57  void deserialize(const std::vector<char>& buffer);
58  bool operator==(const UDQState& other) const;
59 
60  static UDQState serializeObject() {
61  UDQState st;
62  st.undef_value = 78;
63  st.scalar_values = {{"FU1", 100}, {"FU2", 200}};
64  st.assignments = {{"GU1", 99}, {"GU2", 199}};
65  st.defines = {{"DU1", 299}, {"DU2", 399}};
66 
67  st.well_values.emplace("W1", std::unordered_map<std::string, double>{{"U1", 100}, {"U2", 200}});
68  st.well_values.emplace("W2", std::unordered_map<std::string, double>{{"U1", 700}, {"32", 600}});
69 
70  st.group_values.emplace("G1", std::unordered_map<std::string, double>{{"U1", 100}, {"U2", 200}});
71  st.group_values.emplace("G2", std::unordered_map<std::string, double>{{"U1", 700}, {"32", 600}});
72  return st;
73  }
74 
75 
76  template<class Serializer>
77  void pack_unpack_wgmap(Serializer& serializer, std::unordered_map<std::string, std::unordered_map<std::string, double>>& wgmap) {
78  std::size_t map_size = wgmap.size();
79  serializer(map_size);
80  if (serializer.isSerializing()) {
81  for (auto& [udq_key, values] : wgmap) {
82  serializer(udq_key);
83  serializer.template map<std::unordered_map<std::string, double>, false>(values);
84  }
85  } else {
86  for (std::size_t index=0; index < map_size; index++) {
87  std::string udq_key;
88  std::unordered_map<std::string, double> inner_map;
89  serializer(udq_key);
90  serializer.template map<std::unordered_map<std::string, double>, false>(inner_map);
91 
92  wgmap.emplace(udq_key, inner_map);
93  }
94  }
95  }
96 
97  template<class Serializer>
98  void serializeOp(Serializer& serializer)
99  {
100  serializer(this->undef_value);
101  serializer.template map<std::unordered_map<std::string, double>, false>(this->scalar_values);
102  serializer.template map<std::unordered_map<std::string, std::size_t>, false>(this->assignments);
103  serializer.template map<std::unordered_map<std::string, std::size_t>, false>(this->defines);
104 
105  pack_unpack_wgmap(serializer, this->well_values);
106  pack_unpack_wgmap(serializer, this->group_values);
107  }
108 
109 
110 private:
111  void add(const std::string& udq_key, const UDQSet& result);
112  double get_wg_var(const std::string& well, const std::string& key, UDQVarType var_type) const;
113  double undef_value;
114  std::unordered_map<std::string, double> scalar_values;
115  std::unordered_map<std::string, std::unordered_map<std::string, double>> well_values;
116  std::unordered_map<std::string, std::unordered_map<std::string, double>> group_values;
117  std::unordered_map<std::string, std::size_t> assignments;
118  std::unordered_map<std::string, std::size_t> defines;
119 };
120 }
121 
122 #endif
Definition: Serializer.hpp:38
Definition: UDQSet.hpp:66
Definition: UDQState.hpp:36
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29
Definition: state.hpp:51