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 <cstddef>
27#include <optional>
28#include <string>
29#include <vector>
30
31#include <stddef.h>
32
33namespace Opm {
34 class ActiveGridCells;
35 class DeckRecord;
36 class FieldPropsManager;
37 class KeywordLocation;
38 class ScheduleGrid;
39 class EclipseGrid;
40} // namespace Opm
41
42namespace Opm {
43
45 {
46 public:
47 using const_iterator = std::vector<Connection>::const_iterator;
48
49 WellConnections() = default;
50 WellConnections(const Connection::Order ordering, const int headI, const int headJ);
51 WellConnections(const Connection::Order ordering, const int headI, const int headJ,
52 const std::vector<Connection>& connections);
53
54 static WellConnections serializationTestObject();
55
56 // cppcheck-suppress noExplicitConstructor
57 template <class Grid>
58 WellConnections(const WellConnections& src, const Grid& grid)
59 : m_ordering(src.ordering())
60 , headI (src.headI)
61 , headJ (src.headJ)
62 {
63 for (const auto& c : src) {
64 if (grid.isCellActive(c.getI(), c.getJ(), c.getK())) {
65 this->add(c);
66 }
67 }
68 }
69
70 void addConnection(const int i, const int j, const int k,
71 const std::size_t global_index,
72 const double depth,
73 const Connection::State state,
74 const double CF,
75 const double Kh,
76 const double rw,
77 const double r0,
78 const double re,
79 const double connection_length,
80 const double skin_factor,
81 const int satTableId,
82 const Connection::Direction direction = Connection::Direction::Z,
83 const Connection::CTFKind ctf_kind = Connection::CTFKind::DeckValue,
84 const std::size_t seqIndex = 0,
85 const bool defaultSatTabId = true);
86
87 void loadCOMPDAT(const DeckRecord& record,
88 const ScheduleGrid& grid,
89 const std::string& wname,
90 const KeywordLocation& location);
91
92 void loadCOMPTRAJ(const DeckRecord& record, const ScheduleGrid& grid, const std::string& wname, const KeywordLocation& location, external::cvf::ref<external::cvf::BoundingBoxTree>& cellSearchTree);
93
94 void loadWELTRAJ(const DeckRecord& record, const ScheduleGrid& grid, const std::string& wname, const KeywordLocation& location);
95
96 int getHeadI() const;
97 int getHeadJ() const;
98 const std::vector<double>& getMD() const;
99 void add(Connection);
100 std::size_t size() const;
101 bool empty() const;
102 std::size_t num_open() const;
103 const Connection& operator[](size_t index) const;
104 const Connection& get(size_t index) const;
105 const Connection& getFromIJK(const int i, const int j, const int k) const;
106 const Connection& getFromGlobalIndex(std::size_t global_index) const;
107 const Connection& lowest() const;
108 Connection& getFromIJK(const int i, const int j, const int k);
109 bool hasGlobalIndex(std::size_t global_index) const;
110 double segment_perf_length(int segment) const;
111
112 const_iterator begin() const { return this->m_connections.begin(); }
113 const_iterator end() const { return this->m_connections.end(); }
114 void filter(const ActiveGridCells& grid);
115 bool allConnectionsShut() const;
128 void order();
129
130 bool operator==( const WellConnections& ) const;
131 bool operator!=( const WellConnections& ) const;
132
133 Connection::Order ordering() const { return this->m_ordering; }
134 std::vector<const Connection *> output(const EclipseGrid& grid) const;
135
145
153 void applyWellPIScaling(const double scaleFactor,
154 std::vector<bool>& scalingApplicable);
155
156 template <class Serializer>
157 void serializeOp(Serializer& serializer)
158 {
159 serializer(this->m_ordering);
160 serializer(this->headI);
161 serializer(this->headJ);
162 serializer(this->m_connections);
163 serializer(this->coord);
164 serializer(this->md);
165 }
166 private:
167 Connection::Order m_ordering { Connection::Order::TRACK };
168 int headI{0};
169 int headJ{0};
170 std::vector<Connection> m_connections{};
171 std::vector<std::vector<double>> coord{3, std::vector<double>(0, 0.0) };
172 std::vector<double> md{};
173
174 void addConnection(const int i, const int j, const int k,
175 const std::size_t global_index,
176 const int complnum,
177 const double depth,
178 const Connection::State state,
179 const double CF,
180 const double Kh,
181 const double rw,
182 const double r0,
183 const double re,
184 const double connection_length,
185 const double skin_factor,
186 const int satTableId,
187 const Connection::Direction direction = Connection::Direction::Z,
188 const Connection::CTFKind ctf_kind = Connection::CTFKind::DeckValue,
189 const std::size_t seqIndex = 0,
190 const bool defaultSatTabId = true);
191
192 size_t findClosestConnection(int oi, int oj, double oz, size_t start_pos);
193 void orderTRACK();
194 void orderMSW();
195 void orderDEPTH();
196
197 };
198
199 std::optional<int>
200 getCompletionNumberFromGlobalConnectionIndex(const WellConnections& connections,
201 const std::size_t global_index);
202} // namespace Opm
203
204#endif // CONNECTIONSET_HPP_
Simple class capturing active cells of a grid.
Definition ActiveGridCells.hpp:36
Definition Connection.hpp:47
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:29
Class for (de-)serializing.
Definition Serializer.hpp:84
Definition WellConnections.hpp:45
void order()
Order connections irrespective of input order.
void applyWellPIScaling(const double scaleFactor, std::vector< bool > &scalingApplicable)
Scale pertinent connections' CF value by supplied value.
bool prepareWellPIScaling()
Activate or reactivate WELPI scaling for this connection set.
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30