My Project
Loading...
Searching...
No Matches
Segment.hpp
1/*
2 Copyright 2015 SINTEF ICT, Applied Mathematics.
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 SEGMENT_HPP_HEADER_INCLUDED
21#define SEGMENT_HPP_HEADER_INCLUDED
22
23#include <opm/input/eclipse/Schedule/MSW/AICD.hpp>
24#include <opm/input/eclipse/Schedule/MSW/SICD.hpp>
25#include <opm/input/eclipse/Schedule/MSW/Valve.hpp>
26
27#include <optional>
28#include <variant>
29#include <vector>
30
31#include <cmath>
32
33namespace Opm { namespace RestartIO {
34 struct RstSegment;
35}} // namespace Opm::RestartIO
36
37namespace Opm {
38
39 class Segment {
40 public:
41 // Maximum relative roughness to guarantee non-singularity for Re>=4000 in Haaland friction
42 // factor calculations (see MSWellHelpers in opm-simulators).
43 // cannot be constexpr due to clang not having constexpr std::pow
44 static const double MAX_REL_ROUGHNESS; // = 3.7 * std::pow((1.0 - 1.0e-3) - 6.9/4000.0, 9. / 10.);
45
46
47 enum class SegmentType {
48 REGULAR,
49 SICD,
50 AICD, // Not really supported - just included to complete the enum
51 VALVE,
52 };
53
54 Segment();
55
56 Segment(const Segment& src, double new_depth, double new_length, double new_volume, double new_x, double new_y);
57 Segment(const Segment& src, double new_depth, double new_length, double new_x, double new_y);
58 Segment(const Segment& src, double new_depth, double new_length, double new_volume);
59 Segment(const Segment& src, double new_depth, double new_length);
60 Segment(const Segment& src, double new_volume);
61 Segment(const int segment_number_in,
62 const int branch_in,
63 const int outlet_segment_in,
64 const double length_in,
65 const double depth_in,
66 const double internal_diameter_in,
67 const double roughness_in,
68 const double cross_area_in,
69 const double volume_in,
70 const bool data_ready_in,
71 const double x_in,
72 const double y_in);
73
74 explicit Segment(const RestartIO::RstSegment& rst_segment, const std::string& wname);
75
76 static Segment serializationTestObject();
77
78 int segmentNumber() const;
79 int branchNumber() const;
80 int outletSegment() const;
81 double perfLength() const;
82 double totalLength() const;
83 double node_X() const;
84 double node_Y() const;
85 double depth() const;
86 double internalDiameter() const;
87 double roughness() const;
88 double crossArea() const;
89 double volume() const;
90 bool dataReady() const;
91
92 SegmentType segmentType() const;
93 int ecl_type_id() const;
94
95 const std::vector<int>& inletSegments() const;
96
97 static double invalidValue();
98
99 bool operator==( const Segment& ) const;
100 bool operator!=( const Segment& ) const;
101
102 const SICD& spiralICD() const;
103 const AutoICD& autoICD() const;
104 const Valve& valve() const;
105
106 void updatePerfLength(double perf_length);
107 void updateSpiralICD(const SICD& spiral_icd);
108 void updateAutoICD(const AutoICD& aicd);
109 void updateValve(const Valve& valve, const double segment_length);
110 void updateValve(const Valve& valve);
111 void addInletSegment(const int segment_number);
112
113 bool isRegular() const
114 {
115 return std::holds_alternative<RegularSegment>(this->m_icd);
116 }
117
118 inline bool isSpiralICD() const
119 {
120 return std::holds_alternative<SICD>(this->m_icd);
121 }
122
123 inline bool isAICD() const
124 {
125 return std::holds_alternative<AutoICD>(this->m_icd);
126 }
127
128 inline bool isValve() const
129 {
130 return std::holds_alternative<Valve>(this->m_icd);
131 }
132
133 template<class Serializer>
134 void serializeOp(Serializer& serializer)
135 {
136 serializer(m_segment_number);
137 serializer(m_branch);
138 serializer(m_outlet_segment);
139 serializer(m_inlet_segments);
140 serializer(m_total_length);
141 serializer(m_depth);
142 serializer(m_internal_diameter);
143 serializer(m_roughness);
144 serializer(m_cross_area);
145 serializer(m_volume);
146 serializer(m_data_ready);
147 serializer(m_x);
148 serializer(m_y);
149 serializer(m_perf_length);
150 serializer(m_icd);
151 }
152
153 private:
154 // The current serialization of std::variant<> requires that all the
155 // types in the variant have a serializeOp() method. We introduce
156 // this RegularSegment to meet this requirement. Ideally, the ICD
157 // variant<> would use std::monostate to represent non-ICD segments.
158 struct RegularSegment
159 {
160 template <class Serializer>
161 void serializeOp(Serializer&) {}
162
163 static RegularSegment serializationTestObject() { return {}; }
164
165 bool operator==(const RegularSegment&) const { return true; }
166 };
167
168 // segment number
169 // it should work as a ID.
170 int m_segment_number;
171
172 // branch number
173 // for top segment, it should always be 1
174 int m_branch;
175
176 // the outlet junction segment
177 // for top segment, it should be -1
178 int m_outlet_segment;
179
180 // the segments whose outlet segments are the current segment
181 std::vector<int> m_inlet_segments;
182
183 // length of the segment node to the bhp reference point.
184 // when reading in from deck, with 'INC',
185 // it will be incremental length before processing.
186 // After processing and in the class Well, it always stores the 'ABS' value.
187 // which means the total_length
188 double m_total_length;
189
190 // depth of the nodes to the bhp reference point
191 // when reading in from deck, with 'INC',
192 // it will be the incremental depth before processing.
193 // in the class Well, it always stores the 'ABS' value.
194 // TODO: to check if it is good to use 'ABS' always.
195 double m_depth;
196
197 // tubing internal diameter
198 // or the equivalent diameter for annular cross-sections
199 // for top segment, it is UNDEFINED
200 // we use invalid_value for the top segment
201 double m_internal_diameter;
202
203 // effective roughness of the tubing
204 // used to calculate the Fanning friction factor
205 // for top segment, it is UNDEFINED
206 // we use invalid_value for the top segment
207 double m_roughness;
208
209 // cross-sectional area for fluid flow
210 // not defined for the top segment,
211 // we use invalid_value for the top segment.
212 double m_cross_area;
213
214 // valume of the segment;
215 // it is defined for top segment.
216 // TODO: to check if the definition is the same with other segments.
217 double m_volume;
218
219 // indicate if the data related to 'INC' or 'ABS' is ready
220 // the volume will be updated at a final step.
221 bool m_data_ready;
222
223 // Length of segment projected onto the X axis. Not used in
224 // simulations, but needed for the SEG option in WRFTPLT.
225 double m_x{};
226
227 // Length of segment projected onto the Y axis. Not used in
228 // simulations, but needed for the SEG option in WRFTPLT.
229 double m_y{};
230
231 std::optional<double> m_perf_length;
232 std::variant<RegularSegment, SICD, AutoICD, Valve> m_icd;
233
234 // There are three other properties for the segment pertaining to
235 // thermal conduction. These are not currently supported.
236
237 void updateValve__(Valve& valve, const double segment_length);
238 };
239
240}
241
242#endif
Definition AICD.hpp:44
Definition SICD.hpp:46
Definition Segment.hpp:39
Class for (de-)serializing.
Definition Serializer.hpp:94
Definition Valve.hpp:63
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition segment.hpp:34