My Project
Loading...
Searching...
No Matches
GuideRate.hpp
1/*
2 Copyright 2019, 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 GUIDE_RATE_HPP
21#define GUIDE_RATE_HPP
22
23#include <opm/input/eclipse/Schedule/Group/Group.hpp>
24#include <opm/input/eclipse/Schedule/Group/GuideRateModel.hpp>
25
26#include <cstddef>
27#include <ctime>
28#include <limits>
29#include <memory>
30#include <string>
31#include <unordered_map>
32#include <utility>
33
34namespace Opm {
35
36class Schedule;
37enum class WellGuideRateTarget;
38
39} // namespace Opm
40
41namespace Opm {
42
44{
45public:
46 // used for potentials and well rates
47 struct RateVector {
48 RateVector() = default;
49 RateVector(const double orat, const double grat, const double wrat)
50 : oil_rat(orat)
51 , gas_rat(grat)
52 , wat_rat(wrat)
53 {}
54
55 static RateVector serializationTestObject()
56 {
57 return RateVector{1.0, 2.0, 3.0};
58 }
59
60 double eval(const WellGuideRateTarget target) const;
61 double eval(const Group::GuideRateProdTarget target) const;
62 double eval(const GuideRateModel::Target target) const;
63
64 template<class Serializer>
65 void serializeOp(Serializer& serializer)
66 {
67 serializer(oil_rat);
68 serializer(gas_rat);
69 serializer(wat_rat);
70 }
71
72 double oil_rat{0.0};
73 double gas_rat{0.0};
74 double wat_rat{0.0};
75 };
76
78 GuideRateValue() = default;
79 GuideRateValue(const double t, const double v, const GuideRateModel::Target tg)
80 : sim_time(t)
81 , value (v)
82 , target (tg)
83 {}
84
85 static GuideRateValue serializationTestObject()
86 {
87 return GuideRateValue{1.0, 2.0, GuideRateModel::Target::LIQ};
88 }
89
90 bool operator==(const GuideRateValue& other) const
91 {
92 return (this->sim_time == other.sim_time)
93 && (this->value == other.value);
94 }
95
96 bool operator!=(const GuideRateValue& other) const
97 {
98 return !(*this == other);
99 }
100
101 template<class Serializer>
102 void serializeOp(Serializer& serializer)
103 {
104 serializer(sim_time);
105 serializer(value);
106 serializer(target);
107 }
108
109 double sim_time { std::numeric_limits<double>::lowest() };
110 double value { std::numeric_limits<double>::lowest() };
111 GuideRateModel::Target target { GuideRateModel::Target::NONE };
112 };
113
114 explicit GuideRate(const Schedule& schedule);
115
116 void setSerializationTestData();
117
118 void compute(const std::string& wgname,
119 const std::size_t report_step,
120 const double sim_time,
121 const double oil_pot,
122 const double gas_pot,
123 const double wat_pot);
124
125 void compute(const std::string& wgname,
126 const Phase& phase,
127 const std::size_t report_step,
128 const double guide_rate);
129
130 bool has(const std::string& name) const;
131 bool hasPotentials(const std::string& name) const;
132 bool has(const std::string& name, const Phase& phase) const;
133
134 double get(const std::string& well, const WellGuideRateTarget target, const RateVector& rates) const;
135 double get(const std::string& group, const Group::GuideRateProdTarget target, const RateVector& rates) const;
136 double get(const std::string& name, const GuideRateModel::Target model_target, const RateVector& rates) const;
137 double get(const std::string& group, const Phase& phase) const;
138
139 double getSI(const std::string& well, const WellGuideRateTarget target, const RateVector& rates) const;
140 double getSI(const std::string& group, const Group::GuideRateProdTarget target, const RateVector& rates) const;
141 double getSI(const std::string& wgname, const GuideRateModel::Target target, const RateVector& rates) const;
142 double getSI(const std::string& group, const Phase& phase) const;
143
144 void init_grvalue(const std::size_t report_step, const std::string& wgname, GuideRateValue value);
145 void init_grvalue_SI(const std::size_t report_step, const std::string& wgname, GuideRateValue value);
146
147 void updateGuideRateExpiration(const double sim_time,
148 const std::size_t report_step);
149
150 template<class Serializer>
151 void serializeOp(Serializer& serializer)
152 {
153 if (serializer.isSerializing()) {
154 serializer(values.size());
155 for (const auto& [key, value] : values) {
156 serializer(key);
157 serializer(*value);
158 }
159 } else {
160 std::size_t size = 0;
161 serializer(size);
162 for (size_t i = 0; i < size; ++i) {
163 std::string key;
164 serializer(key);
165 auto val = values.emplace(key, std::make_unique<GRValState>());
166 if (val.first != values.end())
167 serializer(*val.first->second);
168 }
169 }
170 serializer(injection_group_values);
171 serializer(potentials);
172 serializer(guide_rates_expired);
173 }
174
175private:
176 struct GRValState
177 {
178 GuideRateValue curr{};
179 GuideRateValue prev{};
180
181 static GRValState serializationTestObject()
182 {
183 return GRValState{GuideRateValue::serializationTestObject(),
184 GuideRateValue::serializationTestObject()};
185 }
186
187 template<class Serializer>
188 void serializeOp(Serializer& serializer)
189 {
190 serializer(curr);
191 serializer(prev);
192 }
193 };
194
195 struct pair_hash
196 {
197 template <class T1, class T2>
198 std::size_t operator()(const std::pair<T1, T2>& pair) const
199 {
200 return std::hash<T1>()(pair.first) ^ std::hash<T2>()(pair.second);
201 }
202 };
203
204 using GRValPtr = std::unique_ptr<GRValState>;
205 using pair = std::pair<Phase, std::string>;
206
207 void well_compute(const std::string& wgname,
208 const std::size_t report_step,
209 const double sim_time,
210 const double oil_pot,
211 const double gas_pot,
212 const double wat_pot);
213
214 void group_compute(const std::string& wgname,
215 const std::size_t report_step,
216 const double sim_time,
217 const double oil_pot,
218 const double gas_pot,
219 const double wat_pot);
220
221 double eval_form(const GuideRateModel& model,
222 const double oil_pot,
223 const double gas_pot,
224 const double wat_pot) const;
225 double eval_group_pot() const;
226 double eval_group_resvinj() const;
227
228 void assign_grvalue(const std::string& wgname,
229 const GuideRateModel& model,
230 GuideRateValue&& value);
231 double get_grvalue_result(const GRValState& gr) const;
232
233 const Schedule& schedule;
234
235 std::unordered_map<std::string, GRValPtr> values{};
236 std::unordered_map<pair, double, pair_hash> injection_group_values{};
237 std::unordered_map<std::string, RateVector> potentials{};
238 bool guide_rates_expired {false};
239};
240
241} // namespace Opm
242
243#endif
Definition GuideRate.hpp:44
Definition Schedule.hpp:88
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
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition GuideRate.hpp:77
Definition GuideRate.hpp:47