All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
MultisegmentWells.hpp
1 /*
2  Copyright 2016 SINTEF ICT, Applied Mathematics.
3  Copyright 2016 Statoil ASA.
4 
5  This file is part of the Open Porous Media project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #ifndef OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
23 #define OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
24 
25 #include <dune/common/parallel/mpihelper.hh>
26 
27 #include <opm/common/utility/platform_dependent/disable_warnings.h>
28 #include <Eigen/Eigen>
29 #include <Eigen/Sparse>
30 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
31 
32 #include <cassert>
33 
34 #include <opm/core/props/BlackoilPhases.hpp>
35 #include <opm/core/wells/WellCollection.hpp>
36 
37 #include <opm/autodiff/AutoDiffBlock.hpp>
38 #include <opm/autodiff/AutoDiffHelpers.hpp>
39 #include <opm/autodiff/BlackoilModelEnums.hpp>
40 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
41 #include <opm/autodiff/LinearisedBlackoilResidual.hpp>
42 #include <opm/autodiff/WellHelpers.hpp>
43 #include <opm/autodiff/VFPProperties.hpp>
44 
45 #include <opm/autodiff/WellMultiSegment.hpp>
46 #include <opm/autodiff/WellDensitySegmented.hpp>
47 #include <opm/simulators/WellSwitchingLogger.hpp>
48 
49 
50 
51 namespace Opm {
52 
53 
56  public:
57  // --------- Types ---------
58  using ADB = AutoDiffBlock<double>;
59  using Vector = ADB::V;
60 
61  // Well operations and data needed.
63  explicit MultisegmentWellOps(const std::vector<WellMultiSegmentConstPtr>& wells_ms);
64  Eigen::SparseMatrix<double> w2p; // well -> perf (scatter)
65  Eigen::SparseMatrix<double> p2w; // perf -> well (gather)
66  Eigen::SparseMatrix<double> w2s; // well -> segment (scatter)
67  Eigen::SparseMatrix<double> s2w; // segment -> well (gather)
68  Eigen::SparseMatrix<double> s2p; // segment -> perf (scatter)
69  Eigen::SparseMatrix<double> p2s; // perf -> segment (gather)
70  Eigen::SparseMatrix<double> s2s_inlets; // segment -> its inlet segments
71  Eigen::SparseMatrix<double> s2s_outlet; // segment -> its outlet segment
72  Eigen::SparseMatrix<double> topseg2w; // top segment -> well
73  AutoDiffMatrix eliminate_topseg; // change the top segment related to be zero
74  std::vector<int> well_cells; // the set of perforated cells
75  Vector conn_trans_factors; // connection transmissibility factors
76  bool has_multisegment_wells; // flag indicating whether there is any muli-segment well
77  };
78 
79  // copied from BlackoilModelBase
80  // should put to somewhere better
81  using DataBlock = Eigen::Array<double,
82  Eigen::Dynamic,
83  Eigen::Dynamic,
84  Eigen::RowMajor>;
85  using Communication =
86  Dune::CollectiveCommunication<typename Dune::MPIHelper
87  ::MPICommunicator>;
88 
89  // --------- Public methods ---------
90  // TODO: using a vector of WellMultiSegmentConstPtr for now
91  // TODO: it should use const Wells or something else later.
92  MultisegmentWells(const Wells* wells_arg,
93  WellCollection* well_collection,
94  const std::vector< const Well* >& wells_ecl,
95  const int time_step);
96 
97  std::vector<WellMultiSegmentConstPtr> createMSWellVector(const Wells* wells_arg,
98  const std::vector< const Well* >& wells_ecl,
99  const int time_step);
100 
101  void init(const BlackoilPropsAdFromDeck* fluid_arg,
102  const std::vector<bool>* active_arg,
103  const std::vector<PhasePresence>* pc_arg,
104  const VFPProperties* vfp_properties_arg,
105  const double gravity_arg,
106  const Vector& depth_arg);
107 
108  const std::vector<WellMultiSegmentConstPtr>& msWells() const;
109  const MultisegmentWellOps& wellOps() const;
110 
111  const Wells& wells() const;
112 
113  const Wells* wellsPointer() const;
114 
115  int numPhases() const { return num_phases_; };
116 
117  int numWells() const { return wells_multisegment_.size(); }
118 
119  int numSegment() const { return nseg_total_; };
120 
121  int numPerf() const { return nperf_total_; };
122 
123  bool wellsActive() const { return wells_active_; };
124 
125  void setWellsActive(const bool wells_active) { wells_active_ = wells_active; };
126 
127  bool localWellsActive() const { return ! wells_multisegment_.empty(); };
128 
129  int numWellVars() const { return (num_phases_ + 1) * nseg_total_; };
130 
131  template <class ReservoirResidualQuant, class SolutionState>
132  void extractWellPerfProperties(const SolutionState& state,
133  const std::vector<ReservoirResidualQuant>& rq,
134  std::vector<ADB>& mob_perfcells,
135  std::vector<ADB>& b_perfcells) const;
136 
137  Vector& wellPerforationCellPressureDiffs() { return well_perforation_cell_pressure_diffs_; };
138 
139  Vector& wellSegmentPerforationDepthDiffs() { return well_segment_perforation_depth_diffs_; };
140 
141  const Vector& wellPerforationCellDensities() const { return well_perforation_cell_densities_; };
142  Vector& wellPerforationCellDensities() { return well_perforation_cell_densities_; };
143 
144  const std::vector<Vector>& segmentCompSurfVolumeInitial() const { return segment_comp_surf_volume_initial_; };
145  std::vector<Vector>& segmentCompSurfVolumeInitial() { return segment_comp_surf_volume_initial_; };
146 
147  const std::vector<ADB>& segmentCompSurfVolumeCurrent() const { return segment_comp_surf_volume_current_; };
148 
149  const std::vector<int>& topWellSegments() const { return top_well_segments_; };
150  std::vector<int>& topWellSegments() { return top_well_segments_; };
151 
152  Vector& segVDt() { return segvdt_; };
153 
154  const Vector& wellPerforationDensities() const { return well_perforation_densities_; };
155  Vector& wellPerforationDensities() { return well_perforation_densities_; };
156 
157  const Vector& wellPerforationPressureDiffs() const { return well_perforation_pressure_diffs_; };
158  Vector& wellPerforationPressureDiffs() { return well_perforation_pressure_diffs_; };
159 
160  template <class WellState>
161  void
162  updateWellState(const Vector& dwells,
163  const double dpmaxrel,
164  WellState& well_state) const;
165 
166  // TODO: some arguments can be removed later
167  // TODO: compi will be required in the multisegment wells
168  template <class SolutionState>
169  void
170  computeWellFlux(const SolutionState& state,
171  const std::vector<ADB>& mob_perfcells,
172  const std::vector<ADB>& b_perfcells,
173  Vector& aliveWells,
174  std::vector<ADB>& cq_s) const;
175 
176  template <class SolutionState, class WellState>
177  void updatePerfPhaseRatesAndPressures(const std::vector<ADB>& cq_s,
178  const SolutionState& state,
179  WellState& xw) const;
180 
181 
182  // Calculate the density of the mixture in the segments
183  // And the surface volume of the components in the segments by dt
184  template <class SolutionState>
185  void
186  computeSegmentFluidProperties(const SolutionState& state);
187 
188  void
189  computeSegmentPressuresDelta(const double grav);
190 
191  template <class SolutionState>
192  void
193  addWellFluxEq(const std::vector<ADB>& cq_s,
194  const SolutionState& state,
195  LinearisedBlackoilResidual& residual);
196 
197  template <class SolutionState, class WellState>
198  void
199  addWellControlEq(const SolutionState& state,
200  const WellState& xw,
201  const Vector& aliveWells,
202  LinearisedBlackoilResidual& residual);
203 
204  template <class WellState>
205  void
206  updateWellControls(WellState& xw) const;
207 
208  // TODO: these code are same with the StandardWells
209  // to find a better solution later.
210  void
211  variableStateWellIndices(std::vector<int>& indices,
212  int& next) const;
213 
214  template <class SolutionState>
215  void
216  variableStateExtractWellsVars(const std::vector<int>& indices,
217  std::vector<ADB>& vars,
218  SolutionState& state) const;
219 
220  std::vector<int>
221  variableWellStateIndices() const;
222 
223  template <class WellState>
224  void
225  variableWellStateInitials(const WellState& xw,
226  std::vector<Vector>& vars0) const;
227 
228  template <class SolutionState, class WellState>
229  void computeWellConnectionPressures(const SolutionState& state,
230  const WellState& xw,
231  const std::vector<ADB>& kr_adb,
232  const std::vector<ADB>& fluid_density);
233 
234  WellCollection* wellCollection() const;
235 
236  void calculateEfficiencyFactors();
237 
238  const Vector& wellPerfEfficiencyFactors() const;
239 
240 
241  protected:
242  // TODO: probably a wells_active_ will be required here.
243  bool wells_active_;
244  std::vector<WellMultiSegmentConstPtr> wells_multisegment_;
245  MultisegmentWellOps wops_ms_;
246  // It will probably need to be updated during running time.
247  WellCollection* well_collection_;
248 
249  // The efficiency factor for each connection
250  // It is specified based on wells and groups
251  // By default, they should all be one.
252  Vector well_perforation_efficiency_factors_;
253 
254  const int num_phases_;
255  int nseg_total_;
256  int nperf_total_;
257 
258  // TODO: put the Wells here to simplify the interface
259  // TODO: at the moment, both wells_ and wells_multisegment_
260  // TODO: acutally contain all the wells
261  // TODO: they should be split eventually.
262  const Wells* wells_;
263 
264  const BlackoilPropsAdFromDeck* fluid_;
265  const std::vector<bool>* active_;
266  const std::vector<PhasePresence>* phase_condition_;
267  const VFPProperties* vfp_properties_;
268  double gravity_;
269  // The depth of the all the cell centers
270  // It is different from the perforation depth in MultisegmentWells
271  Vector perf_cell_depth_;
272 
273  // Pressure correction due to the different depth of the perforation
274  // and the cell center of the grid block
275  // For the non-segmented wells, since the perforation are forced to be
276  // at the center of the grid cell, it should be ZERO.
277  // It only applies to the mutli-segmented wells.
278  Vector well_perforation_cell_pressure_diffs_;
279 
280  // Pressure correction due to the depth differennce between segment depth and perforation depth.
281  ADB well_segment_perforation_pressure_diffs_;
282 
283  // The depth difference between segment nodes and perforations
284  Vector well_segment_perforation_depth_diffs_;
285 
286  // The depth difference between the perforations and the perforation cells.
287  Vector perf_cell_depth_diffs_;
288 
289  // the average of the fluid densities in the grid block
290  // which is used to calculate the hydrostatic head correction due to the depth difference of the perforation
291  // and the cell center of the grid block
292  Vector well_perforation_cell_densities_;
293 
294  // the density of the fluid mixture in the segments
295  // which is calculated in an implicit way
296  ADB well_segment_densities_;
297 
298  // the hydrostatic pressure drop between segment nodes
299  // calculated with the above density of fluid mixtures
300  // for the top segment, they should always be zero for the moment.
301  ADB well_segment_pressures_delta_;
302 
303  // the surface volume of components in the segments
304  // the initial value at the beginning of the time step
305  std::vector<Vector> segment_comp_surf_volume_initial_;
306 
307  // the value within the current iteration.
308  std::vector<ADB> segment_comp_surf_volume_current_;
309 
310  // the mass flow rate in the segments
311  ADB segment_mass_flow_rates_;
312 
313  // the viscosity of the fluid mixture in the segments
314  // TODO: it is only used to calculate the Reynolds number as we know
315  // maybe it is not better just to store the Reynolds number?
316  ADB segment_viscosities_;
317 
318  std::vector<int> top_well_segments_;
319 
320  // segment volume by dt (time step)
321  // to handle the volume effects of the segment
322  Vector segvdt_;
323 
324  // technically, they are only useful for standard wells
325  // since at the moment, we are handling both the standard
326  // wells and the multi-segment wells under the MultisegmentWells
327  // we need them to remove the dependency on StandardWells
328  Vector well_perforation_densities_;
329  Vector well_perforation_pressure_diffs_;
330 
331  };
332 
333 } // namespace Opm
334 
335 #include "MultisegmentWells_impl.hpp"
336 
337 #endif // OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
A thin wrapper class that holds one VFPProdProperties and one VFPInjProperties object.
Definition: VFPProperties.hpp:37
AutoDiffMatrix is a wrapper class that optimizes matrix operations.
Definition: AutoDiffMatrix.hpp:43
Definition: MultisegmentWells.hpp:62
A class for forward-mode automatic differentiation with vector values and sparse jacobian matrices...
Definition: AutoDiffBlock.hpp:95
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
Eigen::Array< double, Eigen::Dynamic, 1 > V
Underlying type for values.
Definition: AutoDiffBlock.hpp:99
Class for handling the multi-segment well model.
Definition: MultisegmentWells.hpp:55