My Project
Loading...
Searching...
No Matches
Connection.hpp
1/*
2 Copyright 2013 Statoil 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
21#ifndef COMPLETION_HPP_
22#define COMPLETION_HPP_
23
24#include <array>
25#include <cstddef>
26#include <map>
27#include <memory>
28#include <string>
29#include <vector>
30#include <optional>
31#include <limits>
32
33#include <opm/input/eclipse/Schedule/Well/FilterCake.hpp>
34#include <opm/input/eclipse/Schedule/Well/WINJMULT.hpp>
35
36namespace Opm {
37
38namespace RestartIO {
39 struct RstConnection;
40}
41
42 class DeckKeyword;
43 class DeckRecord;
44 class ScheduleGrid;
45 class FieldPropsManager;
46
47 class Connection {
48 public:
49 enum class State {
50 OPEN = 1,
51 SHUT = 2,
52 AUTO = 3 // Seems like the AUTO state can not be serialized to restart files.
53 };
54
55 static const std::string State2String( State enumValue );
56 static State StateFromString( const std::string& stringValue );
57
58
59 enum class Direction{
60 X = 1,
61 Y = 2,
62 Z = 3
63 };
64
65 static std::string Direction2String(const Direction enumValue);
66 static Direction DirectionFromString(const std::string& stringValue);
67
68
69 enum class Order {
70 DEPTH,
71 INPUT,
72 TRACK
73 };
74
75 static const std::string Order2String( Order enumValue );
76 static Order OrderFromString(const std::string& comporderStringValue);
77
78 enum class CTFKind {
80 Defaulted,
81 };
82
83 Connection();
84 Connection(int i, int j , int k ,
85 std::size_t global_index,
86 int complnum,
87 double depth,
88 State state,
89 double CF,
90 double Kh,
91 double rw,
92 double r0,
93 double re,
94 double connection_length,
95 double skin_factor,
96 const int satTableId,
97 const Direction direction,
98 const CTFKind ctf_kind,
99 const std::size_t sort_value,
100 const bool defaultSatTabId);
101
102 Connection(const RestartIO::RstConnection& rst_connection, const ScheduleGrid& grid, const FieldPropsManager& fp);
103
104 static Connection serializationTestObject();
105
106 bool attachedToSegment() const;
107 bool sameCoordinate(const int i, const int j, const int k) const;
108 int getI() const;
109 int getJ() const;
110 int getK() const;
111 std::size_t global_index() const;
112 State state() const;
113 Direction dir() const;
114 double depth() const;
115 int satTableId() const;
116 int complnum() const;
117 int segment() const;
118 double CF() const;
119 double wpimult() const;
120 double Kh() const;
121 double rw() const;
122 double r0() const;
123 double re() const;
124 double connectionLength() const;
125 double skinFactor() const;
126 CTFKind kind() const;
127 const InjMult& injmult() const;
128 bool activeInjMult() const;
129 void setInjMult(const InjMult& inj_mult);
130 void setFilterCake(const FilterCake& filter_cake);
131 const FilterCake& getFilterCake() const;
132 bool filterCakeActive() const;
133 double getFilterCakeRadius() const;
134 double getFilterCakeArea() const;
135
136 void setState(State state);
137 void setComplnum(int compnum);
138 void setSkinFactor(double skin_factor);
139 void setCF(double CF);
140 void scaleWellPi(double wellPi);
141 bool prepareWellPIScaling();
142 bool applyWellPIScaling(const double scaleFactor);
143 void updateSegmentRST(int segment_number_arg,
144 double center_depth_arg);
145 void updateSegment(int segment_number_arg,
146 double center_depth_arg,
147 std::size_t compseg_insert_index,
148 const std::pair<double,double>& perf_range);
149 std::size_t sort_value() const;
150 const bool& getDefaultSatTabId() const;
151 void setDefaultSatTabId(bool id);
152 const std::optional<std::pair<double, double>>& perf_range() const;
153 std::string str() const;
154 bool ctfAssignedFromInput() const
155 {
156 return this->m_ctfkind == CTFKind::DeckValue;
157 }
158
159 bool operator==( const Connection& ) const;
160 bool operator!=( const Connection& ) const;
161
162 template<class Serializer>
163 void serializeOp(Serializer& serializer)
164 {
165 serializer(direction);
166 serializer(center_depth);
167 serializer(open_state);
168 serializer(sat_tableId);
169 serializer(m_complnum);
170 serializer(m_CF);
171 serializer(m_Kh);
172 serializer(m_rw);
173 serializer(m_r0);
174 serializer(m_re);
175 serializer(m_connection_length);
176 serializer(m_skin_factor);
177 serializer(ijk);
178 serializer(m_global_index);
179 serializer(m_ctfkind);
180 serializer(m_injmult);
181 serializer(m_sort_value);
182 serializer(m_perf_range);
183 serializer(m_defaultSatTabId);
184 serializer(segment_number);
185 serializer(m_subject_to_welpi);
186 serializer(m_filter_cake);
187 }
188
189 private:
190 Direction direction;
191 double center_depth;
192 State open_state;
193 int sat_tableId;
194 int m_complnum;
195 double m_CF;
196 double m_Kh;
197 double m_rw;
198 double m_r0;
199 double m_re;
200 double m_connection_length;
201 double m_skin_factor;
202
203 std::array<int,3> ijk;
204 CTFKind m_ctfkind;
205 std::optional<InjMult> m_injmult;
206 std::size_t m_global_index;
207 /*
208 The sort_value member is a peculiar quantity. The connections are
209 assembled in the WellConnections class. During the lifetime of the
210 connections there are three different sort orders which are all
211 relevant:
212
213 input: This is the ordering implied be the order of the
214 connections in the input deck.
215
216 simulation: This is the ordering the connections have in
217 WellConnections container during the simulation and RFT output.
218
219 restart: This is the ordering the connections have when they are
220 written out to a restart file.
221
222 Exactly what consitutes input, simulation and restart ordering, and
223 how the connections transition between the three during application
224 lifetime is different from MSW and normal wells.
225
226 normal wells: For normal wells the simulation order is given by the
227 COMPORD keyword, and then when the connections are serialized to the
228 restart file they are written in input order; i.e. we have:
229
230 input == restart and simulation given COMPORD
231
232 To recover the input order when creating the restart files the
233 sort_value member corresponds to the insert index for normal wells.
234
235 MSW wells: For MSW wells the wells simulator order[*] is given by
236 COMPSEGS keyword, the COMPORD keyword is ignored. The connections are
237 sorted in WellConnections::order() and then retain that order for all
238 eternity i.e.
239
240 input and simulation == restart
241
242 Now the important point is that the COMPSEGS detail used to perform
243 this sorting is not available when loading from a restart file, but
244 then the connections are already sorted correctly. I.e. *after* a
245 restart we will have:
246
247 input(from restart) == simulation == restart
248
249 The sort_value member is used to sort the connections into restart
250 ordering. In the case of normal wells this corresponds to recovering
251 the input order, whereas for MSW wells this is equivalent to the
252 simulation order.
253
254 [*]: For MSW wells the topology is given by the segments and entered
255 explicitly, so the truth is probably that the storage order
256 during simulation makes no difference?
257 */
258
259 std::size_t m_sort_value;
260 std::optional<std::pair<double,double>> m_perf_range;
261 bool m_defaultSatTabId;
262
263 // related segment number
264 // 0 means the completion is not related to segment
265 int segment_number = 0;
266
267 // Whether or not this Connection is subject to WELPI scaling.
268 bool m_subject_to_welpi = false;
269
270 // For applying last known WPIMULT to when calculating connection transmissibilty factor in CSKIN
271 double m_wpimult = 1.0;
272
273 std::optional<FilterCake> m_filter_cake;
274
275 static std::string CTFKindToString(const CTFKind);
276 };
277}
278
279#endif /* COMPLETION_HPP_ */
280
Definition Connection.hpp:47
Definition DeckValue.hpp:30
Definition FieldPropsManager.hpp:41
Definition ScheduleGrid.hpp:29
Class for (de-)serializing.
Definition Serializer.hpp:84
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition FilterCake.hpp:31
Definition WINJMULT.hpp:30
Definition connection.hpp:33