My Project
Loading...
Searching...
No Matches
WellConnections.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#ifndef CONNECTIONSET_HPP_
21#define CONNECTIONSET_HPP_
22
23#include <opm/input/eclipse/Schedule/Well/Connection.hpp>
24#include <external/resinsight/LibGeometry/cvfBoundingBoxTree.h>
25
26#include <array>
27#include <cstddef>
28#include <optional>
29#include <string>
30#include <utility>
31#include <vector>
32
33#include <stddef.h>
34
35namespace Opm {
36 class ActiveGridCells;
37 class DeckRecord;
38 class EclipseGrid;
39 class FieldPropsManager;
40 class KeywordLocation;
41 class ScheduleGrid;
42 class WDFAC;
43} // namespace Opm
44
45namespace Opm {
46
48 {
49 public:
50 using const_iterator = std::vector<Connection>::const_iterator;
51
52 WellConnections() = default;
53 WellConnections(const Connection::Order ordering, const int headI, const int headJ);
54 WellConnections(const Connection::Order ordering, const int headI, const int headJ,
55 const std::vector<Connection>& connections);
56
57 static WellConnections serializationTestObject();
58
59 // cppcheck-suppress noExplicitConstructor
60 template <class Grid>
61 WellConnections(const WellConnections& src, const Grid& grid)
62 : m_ordering(src.ordering())
63 , headI (src.headI)
64 , headJ (src.headJ)
65 {
66 for (const auto& c : src) {
67 if (grid.isCellActive(c.getI(), c.getJ(), c.getK())) {
68 this->add(c);
69 }
70 }
71 }
72
73 void add(const Connection& conn)
74 {
75 this->m_connections.push_back(conn);
76 }
77
78 void addConnection(const int i, const int j, const int k,
79 const std::size_t global_index,
80 const Connection::State state,
81 const double depth,
82 const Connection::CTFProperties& ctf_props,
83 const int satTableId,
84 const Connection::Direction direction = Connection::Direction::Z,
85 const Connection::CTFKind ctf_kind = Connection::CTFKind::DeckValue,
86 const std::size_t seqIndex = 0,
87 const bool defaultSatTabId = true);
88
89 void loadCOMPDAT(const DeckRecord& record,
90 const ScheduleGrid& grid,
91 const std::string& wname,
92 const WDFAC& wdfac,
93 const KeywordLocation& location);
94
95 void loadCOMPTRAJ(const DeckRecord& record,
96 const ScheduleGrid& grid,
97 const std::string& wname,
98 const KeywordLocation& location,
99 external::cvf::ref<external::cvf::BoundingBoxTree>& cellSearchTree);
100
101 void loadWELTRAJ(const DeckRecord& record,
102 const ScheduleGrid& grid,
103 const std::string& wname,
104 const KeywordLocation& location);
105
106 void applyDFactorCorrelation(const ScheduleGrid& grid,
107 const WDFAC& wdfac);
108
109 int getHeadI() const;
110 int getHeadJ() const;
111 const std::vector<double>& getMD() const;
112 std::size_t size() const;
113 bool empty() const;
114 std::size_t num_open() const;
115 const Connection& operator[](size_t index) const;
116 const Connection& get(size_t index) const;
117 const Connection& getFromIJK(const int i, const int j, const int k) const;
118 const Connection& getFromGlobalIndex(std::size_t global_index) const;
119 const Connection& lowest() const;
120 Connection& getFromIJK(const int i, const int j, const int k);
121 bool hasGlobalIndex(std::size_t global_index) const;
122 double segment_perf_length(int segment) const;
123
124 const_iterator begin() const { return this->m_connections.begin(); }
125 const_iterator end() const { return this->m_connections.end(); }
126 void filter(const ActiveGridCells& grid);
127 bool allConnectionsShut() const;
140 void order();
141
142 bool operator==( const WellConnections& ) const;
143 bool operator!=( const WellConnections& ) const;
144
145 Connection::Order ordering() const { return this->m_ordering; }
146 std::vector<const Connection *> output(const EclipseGrid& grid) const;
147
157
165 void applyWellPIScaling(const double scaleFactor,
166 std::vector<bool>& scalingApplicable);
167
168 template <class Serializer>
169 void serializeOp(Serializer& serializer)
170 {
171 serializer(this->m_ordering);
172 serializer(this->headI);
173 serializer(this->headJ);
174 serializer(this->m_connections);
175 serializer(this->coord);
176 serializer(this->md);
177 }
178
179 private:
180 Connection::Order m_ordering { Connection::Order::TRACK };
181 int headI{0};
182 int headJ{0};
183 std::vector<Connection> m_connections{};
184
185 std::array<std::vector<double>, 3> coord{};
186 std::vector<double> md{};
187
188 void addConnection(const int i, const int j, const int k,
189 const std::size_t global_index,
190 const int complnum,
191 const Connection::State state,
192 const double depth,
193 const Connection::CTFProperties& ctf_props,
194 const int satTableId,
195 const Connection::Direction direction,
196 const Connection::CTFKind ctf_kind,
197 const std::size_t seqIndex,
198 const bool defaultSatTabId);
199
200 size_t findClosestConnection(int oi, int oj, double oz, size_t start_pos);
201 void orderTRACK();
202 void orderMSW();
203 void orderDEPTH();
204 };
205
206 std::optional<int>
207 getCompletionNumberFromGlobalConnectionIndex(const WellConnections& connections,
208 const std::size_t global_index);
209} // namespace Opm
210
211#endif // CONNECTIONSET_HPP_
Simple class capturing active cells of a grid.
Definition ActiveGridCells.hpp:36
Definition Connection.hpp:53
Definition DeckRecord.hpp:32
About cell information and dimension: The actual grid information is held in a pointer to an ERT ecl_...
Definition EclipseGrid.hpp:55
Definition KeywordLocation.hpp:27
Definition ScheduleGrid.hpp:37
Class for (de-)serializing.
Definition Serializer.hpp:84
Definition WDFAC.hpp:40
Definition WellConnections.hpp:48
void order()
Order connections irrespective of input order.
Definition WellConnections.cpp:873
void applyWellPIScaling(const double scaleFactor, std::vector< bool > &scalingApplicable)
Scale pertinent connections' CF value by supplied value.
Definition WellConnections.cpp:301
bool prepareWellPIScaling()
Activate or reactivate WELPI scaling for this connection set.
Definition WellConnections.cpp:291
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Quantities that go into calculating the connection transmissibility factor.
Definition Connection.hpp:90