All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
WellMultiSegment.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 OPM_WELLMULTISEGMENT_HEADER_INCLUDED
21 #define OPM_WELLMULTISEGMENT_HEADER_INCLUDED
22 
23 
24 #include <opm/common/utility/platform_dependent/disable_warnings.h>
25 #include <Eigen/Eigen>
26 #include <Eigen/Sparse>
27 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
28 
29 #include <opm/core/wells.h>
30 #include <opm/core/well_controls.h>
31 #include <opm/core/simulator/WellState.hpp>
32 #include <opm/common/ErrorMacros.hpp>
33 #include <opm/parser/eclipse/EclipseState/Schedule/ScheduleEnums.hpp>
34 #include <opm/parser/eclipse/EclipseState/Schedule/MSW/SegmentSet.hpp>
35 #include <opm/parser/eclipse/EclipseState/Schedule/Well.hpp>
36 #include <vector>
37 #include <cassert>
38 #include <string>
39 #include <utility>
40 #include <map>
41 #include <algorithm>
42 #include <array>
43 
44 namespace Opm
45 {
46 
48  {
49  public:
50 
51  typedef Eigen::SparseMatrix<double> Matrix;
52 
57  WellMultiSegment(const Well* well, size_t time_step, const Wells* wells);
58 
60  const std::string& name() const;
61 
63  bool isMultiSegmented() const;
64 
66  int numberOfPerforations() const;
67 
69  int numberOfSegments() const;
70 
75  std::string compPressureDrop() const;
76 
78  const WellControls* wellControls() const;
79 
81  const std::vector<double>& compFrac() const;
82 
84  int numberOfPhases() const;
85 
87  WellType wellType() const;
88 
90  const std::vector<double>& wellIndex() const;
91 
93  const std::vector<double>& perfDepth() const;
94 
96  const std::vector<int>& wellCells() const;
97 
99  const std::vector<int>& segmentCells() const;
100 
103  const std::vector<int>& outletSegment() const;
104 
106  const std::vector<std::vector<int>>& inletSegments() const;
107 
109  const std::vector<double>& segmentLength() const;
110 
112  const std::vector<double>& segmentDepth() const;
113 
115  const std::vector<double>& segmentDiameter() const;
116 
118  const std::vector<double>& segmentCrossArea() const;
119 
121  const std::vector<double>& segmentRoughness() const;
122 
124  const std::vector<double>& segmentVolume() const;
125 
127  const std::vector<std::vector<int>>& segmentPerforations() const;
128 
131  struct WellOps {
132  Matrix s2p; // segment -> perf (scatter)
133  Matrix p2s; // perf -> segment (gather)
134  Matrix p2s_average; // perf -> segment (avarage)
135  Matrix s2s_gather; // segment -> segment (in an accumlative way)
136  // means the outlet segments will gather all the contribution
137  // from all the inlet segments in a recurisive way
138  Matrix p2s_gather; // perforation -> segment (in an accumative way)
139  Matrix s2s_inlets; // segment -> its inlet segments
140  Matrix s2s_outlet; // segment -> its outlet segment
141  };
142 
144  const WellOps& wellOps() const;
145 
146  private:
147  // for the moment, we use the information from wells.
148  // TODO: remove the dependency on wells from opm-core.
149  void initMultiSegmentWell(const Well* well, size_t time_step, const Wells* wells);
150  void initNonMultiSegmentWell(const Well* well, size_t time_step, const Wells* wells);
151  void updateWellOps();
152 
153  private:
154  // name of the well
155  std::string m_well_name_;
156  // flag to indicate if this well is a
157  // multi-segmented well
158  bool m_is_multi_segment_;
159  // well type
160  // INJECTOR or PRODUCER
161  enum WellType m_well_type_;
162  // number of phases
163  int m_number_of_phases_;
164  // component fractions for each well
165  std::vector<double> m_comp_frac_;
166  // controls for this well
167  // using pointer for temporary
168  // changing it when figuring out how to using it
169  struct WellControls *m_well_controls_;
170  // components of the pressure drop to be included
171  WellSegment::CompPressureDropEnum m_comp_pressure_drop_;
172  // multi-phase flow model
173  WellSegment::MultiPhaseModelEnum m_multiphase_model_;
174  // number of perforation for this well
175  int m_number_of_perforations_;
176  // number of segments for this well
177  int m_number_of_segments_;
178  // well index for each completion
179  std::vector<double> m_well_index_;
180  // depth for each completion // form the keyword COMPSEGS?
181  std::vector<double> m_perf_depth_;
182  // well cell for each completion
183  std::vector<int> m_well_cell_;
184  // cell for each segment
185  std::vector<int> m_segment_cell_;
186  // how to organize the segment structure here?
187  // indicate the outlet segment for each segment
188  // maybe here we can use the location in the vector
189  // at the moment, we still use the ID number
190  // then a mapping from the ID number to the actual location will be required
191  // The ID is already changed to the location now.
192  std::vector<int> m_outlet_segment_;
193  // for convinience, we store the inlet segments for each segment
194  std::vector<std::vector<int>> m_inlet_segments_;
195  // this one is not necessary any more, since the segment number is its location.
196  // std::map<int, int> m_number_to_location_;
197  // has not decided to use the absolute length from the well head
198  // or the length of this single segment
199  // using the absolute length first
200  std::vector<double> m_segment_length_;
201  // the depth of the segmnet node
202  std::vector<double> m_segment_depth_;
203  // the internal diameter of the segment
204  std::vector<double> m_segment_internal_diameter_;
205  // the roughness of the segment
206  std::vector<double> m_segment_roughness_;
207  // the cross area of the segment
208  std::vector<double> m_segment_cross_area_;
209  // the volume of the segment
210  std::vector<double> m_segment_volume_;
211  // the completions that is related to each segment
212  // the completions's ids are their location in the vector m_well_index_
213  // m_well_cell_
214  // This is also assuming the order of the completions in Well is the same with
215  // the order of the completions in wells.
216  std::vector<std::vector<int>> m_segment_perforations_;
217 
218  WellOps m_wops_;
219  };
220 
221  typedef std::shared_ptr<WellMultiSegment> WellMultiSegmentPtr;
222  typedef std::shared_ptr<const WellMultiSegment> WellMultiSegmentConstPtr;
223 
224 } // namespace Opm
225 
226 
227 #endif // OPM_WELLMULTISEGMENT_HEADER_INCLUDE
Struct for the well operator matrices.
Definition: WellMultiSegment.hpp:131
const WellOps & wellOps() const
Well operator matrics.
Definition: WellMultiSegment.cpp:399
const std::vector< double > & perfDepth() const
Depth of the perforations.
Definition: WellMultiSegment.cpp:351
bool isMultiSegmented() const
Flag indicating if the well is a multi-segment well.
Definition: WellMultiSegment.cpp:315
WellMultiSegment(const Well *well, size_t time_step, const Wells *wells)
Constructor of WellMultiSegment.
Definition: WellMultiSegment.cpp:29
const std::vector< double > & segmentDepth() const
The depth of the segment nodes.
Definition: WellMultiSegment.cpp:375
const std::vector< double > & segmentDiameter() const
Tubing internal diameter.
Definition: WellMultiSegment.cpp:379
const std::vector< double > & segmentRoughness() const
Effective absolute roughness of the tube.
Definition: WellMultiSegment.cpp:387
const std::vector< double > & wellIndex() const
Well productivity index.
Definition: WellMultiSegment.cpp:347
std::string compPressureDrop() const
Components of the pressure drop invloved.
Definition: WellMultiSegment.cpp:335
int numberOfSegments() const
Number of the segments.
Definition: WellMultiSegment.cpp:331
WellType wellType() const
Well type.
Definition: WellMultiSegment.cpp:319
const WellControls * wellControls() const
Well control.
Definition: WellMultiSegment.cpp:323
const std::vector< std::vector< int > > & inletSegments() const
Inlet segments. a segment can have more than one inlet segments.
Definition: WellMultiSegment.cpp:367
int numberOfPerforations() const
Number of the perforations.
Definition: WellMultiSegment.cpp:327
const std::vector< double > & segmentCrossArea() const
Cross-sectional area.
Definition: WellMultiSegment.cpp:383
const std::vector< int > & segmentCells() const
Indices of the gird blocks that segments locate at.
Definition: WellMultiSegment.cpp:359
const std::vector< double > & segmentLength() const
The length of the segment nodes down the wellbore from the reference point.
Definition: WellMultiSegment.cpp:371
const std::vector< double > & compFrac() const
Component fractions for each well.
Definition: WellMultiSegment.cpp:339
Definition: WellMultiSegment.hpp:47
int numberOfPhases() const
Number of phases.
Definition: WellMultiSegment.cpp:343
const std::vector< int > & wellCells() const
Indices of the grid blocks that perforations are completed in.
Definition: WellMultiSegment.cpp:355
const std::vector< std::vector< int > > & segmentPerforations() const
Perforations related to each segment.
Definition: WellMultiSegment.cpp:395
const std::string & name() const
Well name.
Definition: WellMultiSegment.cpp:311
const std::vector< double > & segmentVolume() const
Volume of segment.
Definition: WellMultiSegment.cpp:391
const std::vector< int > & outletSegment() const
Outlet segments, a segment (except top segment) can only have one outlet segment. ...
Definition: WellMultiSegment.cpp:363