All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
BlackoilWellModel.hpp
1 /*
2  Copyright 2016 SINTEF ICT, Applied Mathematics.
3  Copyright 2016 - 2017 Statoil ASA.
4  Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5  Copyright 2016 - 2017 IRIS AS
6 
7  This file is part of the Open Porous Media project (OPM).
8 
9  OPM is free software: you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation, either version 3 of the License, or
12  (at your option) any later version.
13 
14  OPM is distributed in the hope that it will be useful,
15  but WITHOUT ANY WARRANTY; without even the implied warranty of
16  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  GNU General Public License for more details.
18 
19  You should have received a copy of the GNU General Public License
20  along with OPM. If not, see <http://www.gnu.org/licenses/>.
21 */
22 
23 
24 #ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
25 #define OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
26 
27 #include <opm/common/OpmLog/OpmLog.hpp>
28 
29 #include <opm/common/utility/platform_dependent/disable_warnings.h>
30 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
31 
32 #include <cassert>
33 #include <tuple>
34 
35 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
36 
37 #include <opm/core/wells.h>
38 #include <opm/core/wells/DynamicListEconLimited.hpp>
39 #include <opm/core/wells/WellCollection.hpp>
40 #include <opm/core/simulator/SimulatorReport.hpp>
41 #include <opm/autodiff/VFPProperties.hpp>
42 #include <opm/autodiff/WellHelpers.hpp>
43 #include <opm/autodiff/BlackoilModelEnums.hpp>
44 #include <opm/autodiff/WellDensitySegmented.hpp>
45 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
46 #include <opm/autodiff/BlackoilDetails.hpp>
47 #include <opm/autodiff/BlackoilModelParameters.hpp>
48 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
50 #include <opm/autodiff/WellInterface.hpp>
51 #include <opm/autodiff/StandardWell.hpp>
52 #include <opm/autodiff/MultisegmentWell.hpp>
53 #include<dune/common/fmatrix.hh>
54 #include<dune/istl/bcrsmatrix.hh>
55 #include<dune/istl/matrixmatrix.hh>
56 
57 #include <opm/material/densead/Math.hpp>
58 
59 #include <opm/simulators/WellSwitchingLogger.hpp>
60 
61 
62 namespace Opm {
63 
65  template<typename TypeTag>
67  public:
68  // --------- Types ---------
71 
72  typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
73  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
74  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
75  typedef typename GET_PROP_TYPE(TypeTag, Indices) BlackoilIndices;
76  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
77  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
78 
79  static const int numEq = BlackoilIndices::numEq;
80  static const int solventSaturationIdx = BlackoilIndices::solventSaturationIdx;
81 
82  // TODO: where we should put these types, WellInterface or Well Model?
83  // or there is some other strategy, like TypeTag
84  typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
85  typedef Dune::BlockVector<VectorBlockType> BVector;
86 
87  typedef Ewoms::BlackOilPolymerModule<TypeTag> PolymerModule;
88 
89  // For the conversion between the surface volume rate and resrevoir voidage rate
92 
93  // --------- Public methods ---------
94  BlackoilWellModel(const Wells* wells_arg,
95  WellCollection* well_collection,
96  const std::vector< const Well* >& wells_ecl,
97  const ModelParameters& param,
98  const RateConverterType& rate_converter,
99  const bool terminal_output,
100  const int current_index,
101  const std::vector<int>& pvt_region_idx);
102 
103  void init(const PhaseUsage phase_usage_arg,
104  const std::vector<bool>& active_arg,
105  const double gravity_arg,
106  const std::vector<double>& depth_arg,
107  long int global_nc,
108  const Grid& grid);
109 
110  void setVFPProperties(const VFPProperties* vfp_properties_arg);
111 
112 
113  SimulatorReport assemble(Simulator& ebosSimulator,
114  const int iterationIdx,
115  const double dt,
116  WellState& well_state);
117 
118  // substract Binv(D)rw from r;
119  void apply( BVector& r) const;
120 
121  // subtract B*inv(D)*C * x from A*x
122  void apply(const BVector& x, BVector& Ax) const;
123 
124  // apply well model with scaling of alpha
125  void applyScaleAdd(const Scalar alpha, const BVector& x, BVector& Ax) const;
126 
127  // using the solution x to recover the solution xw for wells and applying
128  // xw to update Well State
129  void recoverWellSolutionAndUpdateWellState(const BVector& x, WellState& well_state) const;
130 
131  int numWells() const;
132 
134  bool wellsActive() const;
135 
136  void setWellsActive(const bool wells_active);
137 
139  bool localWellsActive() const;
140 
141  bool getWellConvergence(const Simulator& ebosSimulator,
142  const std::vector<Scalar>& B_avg) const;
143 
145  void updateListEconLimited(const Schedule& schedule,
146  const int current_step,
147  const Wells* wells_struct,
148  const WellState& well_state,
149  DynamicListEconLimited& list_econ_limited) const;
150 
151  WellCollection* wellCollection() const;
152 
153 
154  protected:
155  bool wells_active_;
156  const Wells* wells_;
157  const std::vector< const Well* > wells_ecl_;
158 
159  // the number of wells in this process
160  // trying not to use things from Wells struct
161  // TODO: maybe a better name to emphasize it is local?
162  const int number_of_wells_;
163 
164  const int number_of_phases_;
165 
166  const ModelParameters& param_;
167 
168  using WellInterfacePtr = std::unique_ptr<WellInterface<TypeTag> >;
169  // a vector of all the wells.
170  // eventually, the wells_ above should be gone.
171  // the name is just temporary
172  // later, might make share_ptr const later.
173  std::vector<WellInterfacePtr > well_container_;
174 
175  using ConvergenceReport = typename WellInterface<TypeTag>::ConvergenceReport;
176 
177  // create the well container
178  static std::vector<WellInterfacePtr > createWellContainer(const Wells* wells,
179  const std::vector<const Well*>& wells_ecl,
180  const bool use_multisegment_well,
181  const int time_step,
182  const ModelParameters& param);
183 
184  // Well collection is used to enforce the group control
185  WellCollection* well_collection_;
186 
187  bool terminal_output_;
188  bool has_solvent_;
189  bool has_polymer_;
190  int current_timeIdx_;
191 
192  PhaseUsage phase_usage_;
193  std::vector<bool> active_;
194  const RateConverterType& rate_converter_;
195  const std::vector<int>& pvt_region_idx_;
196 
197  // the number of the cells in the local grid
198  int number_of_cells_;
199 
200  long int global_nc_;
201 
202  // used to better efficiency of calcuation
203  mutable BVector scaleAddRes_;
204 
205  void updateWellControls(WellState& xw) const;
206 
207  void updateGroupControls(WellState& well_state) const;
208 
209  // setting the well_solutions_ based on well_state.
210  void updatePrimaryVariables(const WellState& well_state) const;
211 
212  void setupCompressedToCartesian(const int* global_cell, int number_of_cells, std::map<int,int>& cartesian_to_compressed ) const;
213 
214  void computeRepRadiusPerfLength(const Grid& grid);
215 
216 
217  void computeAverageFormationFactor(const Simulator& ebosSimulator,
218  std::vector<double>& B_avg) const;
219 
220  void applyVREPGroupControl(WellState& well_state) const;
221 
222  void computeWellVoidageRates(const WellState& well_state,
223  std::vector<double>& well_voidage_rates,
224  std::vector<double>& voidage_conversion_coeffs) const;
225 
226  // Calculating well potentials for each well
227  // TODO: getBhp() will be refactored to reduce the duplication of the code calculating the bhp from THP.
228  void computeWellPotentials(const Simulator& ebosSimulator,
229  const WellState& well_state,
230  std::vector<double>& well_potentials) const;
231 
232  const std::vector<double>& wellPerfEfficiencyFactors() const;
233 
234  void calculateEfficiencyFactors();
235 
236  // it should be able to go to prepareTimeStep(), however, the updateWellControls() and initPrimaryVariablesEvaluation()
237  // makes it a little more difficult. unless we introduce if (iterationIdx != 0) to avoid doing the above functions
238  // twice at the beginning of the time step
241  void calculateExplicitQuantities(const Simulator& ebosSimulator,
242  const WellState& xw) const;
243 
244  SimulatorReport solveWellEq(Simulator& ebosSimulator,
245  const double dt,
246  WellState& well_state) const;
247 
248  void initPrimaryVariablesEvaluation() const;
249 
250  // The number of components in the model.
251  int numComponents() const
252  {
253  if (numPhases() == 2) {
254  return 2;
255  }
256  int numComp = FluidSystem::numComponents;
257  if (has_solvent_) {
258  numComp ++;
259  }
260 
261  return numComp;
262  }
263 
264  int numPhases() const;
265 
266  int flowPhaseToEbosPhaseIdx( const int phaseIdx ) const;
267 
268  void resetWellControlFromState(const WellState& xw) const;
269 
270  void assembleWellEq(Simulator& ebosSimulator,
271  const double dt,
272  WellState& well_state,
273  bool only_wells) const;
274 
275  // some preparation work, mostly related to group control and RESV,
276  // at the beginning of each time step (Not report step)
277  void prepareTimeStep(const Simulator& ebos_simulator,
278  WellState& well_state) const;
279 
280  void prepareGroupControl(const Simulator& ebos_simulator,
281  WellState& well_state) const;
282 
283  };
284 
285 
286 } // namespace Opm
287 
288 #include "BlackoilWellModel_impl.hpp"
289 #endif
A thin wrapper class that holds one VFPProdProperties and one VFPInjProperties object.
Definition: VFPProperties.hpp:37
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:66
a struct to collect information about the convergence checking
Definition: WellInterface.hpp:117
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilWellModel_impl.hpp:354
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
void updateListEconLimited(const Schedule &schedule, const int current_step, const Wells *wells_struct, const WellState &well_state, DynamicListEconLimited &list_econ_limited) const
upate the dynamic lists related to economic limits
Definition: BlackoilWellModel_impl.hpp:579
bool localWellsActive() const
return true if wells are available on this process
Definition: BlackoilWellModel_impl.hpp:378
void calculateExplicitQuantities(const Simulator &ebosSimulator, const WellState &xw) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:540