WellInterface.hpp
1 /*
2  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3  Copyright 2017 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_WELLINTERFACE_HEADER_INCLUDED
23 #define OPM_WELLINTERFACE_HEADER_INCLUDED
24 
25 #include <opm/common/OpmLog/OpmLog.hpp>
26 
27 
28 #include <opm/parser/eclipse/EclipseState/Schedule/Well.hpp>
29 #include <opm/core/wells.h>
30 #include <opm/core/well_controls.h>
31 #include <opm/core/props/BlackoilPhases.hpp>
32 #include <opm/core/wells/WellsManager.hpp>
33 
34 #include <opm/autodiff/VFPProperties.hpp>
35 #include <opm/autodiff/VFPInjProperties.hpp>
36 #include <opm/autodiff/VFPProdProperties.hpp>
37 #include <opm/autodiff/WellHelpers.hpp>
38 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
39 #include <opm/autodiff/BlackoilModelParameters.hpp>
40 
41 #include <opm/simulators/WellSwitchingLogger.hpp>
42 
43 #include<dune/common/fmatrix.hh>
44 #include<dune/istl/bcrsmatrix.hh>
45 #include<dune/istl/matrixmatrix.hh>
46 
47 #include <opm/material/densead/Math.hpp>
48 #include <opm/material/densead/Evaluation.hpp>
49 
50 #include <string>
51 #include <memory>
52 #include <vector>
53 #include <cassert>
54 
55 namespace Opm
56 {
57 
58 
59  template<typename TypeTag>
61  {
62  public:
63 
65 
67  typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
68  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
69  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
70  typedef typename GET_PROP_TYPE(TypeTag, Indices) BlackoilIndices;
71  typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
72  typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
73 
74  static const int numEq = BlackoilIndices::numEq;
75  typedef double Scalar;
76 
77  typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
78  typedef Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
79  typedef Dune::BCRSMatrix <MatrixBlockType> Mat;
80  typedef Dune::BlockVector<VectorBlockType> BVector;
81  typedef DenseAd::Evaluation<double, /*size=*/numEq> Eval;
82 
83  typedef Ewoms::BlackOilPolymerModule<TypeTag> PolymerModule;
84 
85  static const bool has_solvent = GET_PROP_VALUE(TypeTag, EnableSolvent);
86  static const bool has_polymer = GET_PROP_VALUE(TypeTag, EnablePolymer);
87 
89  WellInterface(const Well* well, const int time_step, const Wells* wells, const ModelParameters& param);
90 
92  virtual ~WellInterface() {}
93 
95  const std::string& name() const;
96 
98  const std::vector<int>& cells() {return well_cells_; }
99 
101  WellType wellType() const;
102 
104  WellControls* wellControls() const;
105 
106  void setVFPProperties(const VFPProperties* vfp_properties_arg);
107 
108  virtual void init(const PhaseUsage* phase_usage_arg,
109  const std::vector<bool>* active_arg,
110  const std::vector<double>& depth_arg,
111  const double gravity_arg,
112  const int num_cells);
113 
114  virtual void initPrimaryVariablesEvaluation() const = 0;
115 
118  struct ProblemWell {
119  std::string well_name;
120  std::string phase_name;
121  };
122  bool converged = true;
123  bool nan_residual_found = false;
124  std::vector<ProblemWell> nan_residual_wells;
125  // We consider Inf is large residual here
126  bool too_large_residual_found = false;
127  std::vector<ProblemWell> too_large_residual_wells;
128 
129  ConvergenceReport& operator+=(const ConvergenceReport& rhs) {
130  converged = converged && rhs.converged;
131  nan_residual_found = nan_residual_found || rhs.nan_residual_found;
132  if (rhs.nan_residual_found) {
133  for (const ProblemWell& well : rhs.nan_residual_wells) {
134  nan_residual_wells.push_back(well);
135  }
136  }
137  too_large_residual_found = too_large_residual_found || rhs.too_large_residual_found;
138  if (rhs.too_large_residual_found) {
139  for (const ProblemWell& well : rhs.too_large_residual_wells) {
140  too_large_residual_wells.push_back(well);
141  }
142  }
143  return *this;
144  }
145  };
146 
147  virtual ConvergenceReport getWellConvergence(const std::vector<double>& B_avg) const = 0;
148 
149  virtual void solveEqAndUpdateWellState(WellState& well_state) = 0;
150 
151  virtual void assembleWellEq(Simulator& ebosSimulator,
152  const double dt,
153  WellState& well_state,
154  bool only_wells) = 0;
155 
156  void updateListEconLimited(const WellState& well_state,
157  DynamicListEconLimited& list_econ_limited) const;
158 
159  void setWellEfficiencyFactor(const double efficiency_factor);
160 
161  void computeRepRadiusPerfLength(const Grid& grid, const std::map<int, int>& cartesian_to_compressed);
162 
165  virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
166  WellState& well_state) const = 0;
167 
169  virtual void apply(const BVector& x, BVector& Ax) const = 0;
170 
172  virtual void apply(BVector& r) const = 0;
173 
174  // TODO: before we decide to put more information under mutable, this function is not const
175  virtual void computeWellPotentials(const Simulator& ebosSimulator,
176  const WellState& well_state,
177  std::vector<double>& well_potentials) = 0;
178 
179  virtual void updateWellStateWithTarget(const int current,
180  WellState& xw) const = 0;
181 
182  void updateWellControl(WellState& xw,
183  wellhelpers::WellSwitchingLogger& logger) const;
184 
185  virtual void updatePrimaryVariables(const WellState& well_state) const = 0;
186 
187  virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
188  const WellState& xw) = 0; // should be const?
189 
190  protected:
191 
192  // to indicate a invalid connection
193  static const int INVALIDCONNECTION = -100000;
194 
195  const Well* well_ecl_;
196 
197  const int current_step_;
198 
199  // the index of well in Wells struct
200  int index_of_well_;
201 
202  // simulation parameters
203  const ModelParameters& param_;
204 
205  // well type
206  // INJECTOR or PRODUCER
207  enum WellType well_type_;
208 
209  // number of phases
210  int number_of_phases_;
211 
212  // component fractions for each well
213  // typically, it should apply to injection wells
214  std::vector<double> comp_frac_;
215 
216  // controls for this well
217  struct WellControls* well_controls_;
218 
219  // number of the perforations for this well
220  int number_of_perforations_;
221 
222  // record the index of the first perforation
223  // of states of individual well.
224  int first_perf_;
225 
226  // well index for each perforation
227  std::vector<double> well_index_;
228 
229  // depth for each perforation
230  std::vector<double> perf_depth_;
231 
232  // reference depth for the BHP
233  double ref_depth_;
234 
235  double well_efficiency_factor_;
236 
237  // cell index for each well perforation
238  std::vector<int> well_cells_;
239 
240  // saturation table nubmer for each well perforation
241  std::vector<int> saturation_table_number_;
242 
243  // representative radius of the perforations, used in shear calculation
244  std::vector<double> perf_rep_radius_;
245 
246  // length of the perforations, use in shear calculation
247  std::vector<double> perf_length_;
248 
249  // well bore diameter
250  std::vector<double> bore_diameters_;
251 
252  const PhaseUsage* phase_usage_;
253 
254  bool getAllowCrossFlow() const;
255 
256  const std::vector<bool>* active_;
257 
258  const VFPProperties* vfp_properties_;
259 
260  double gravity_;
261 
262  const std::vector<bool>& active() const;
263 
264  const PhaseUsage& phaseUsage() const;
265 
266  int flowPhaseToEbosCompIdx( const int phaseIdx ) const;
267 
268  int flowPhaseToEbosPhaseIdx( const int phaseIdx ) const;
269 
270  // TODO: it is dumplicated with BlackoilWellModel
271  int numComponents() const;
272 
273  double wsolvent() const;
274 
275  double wpolymer() const;
276 
277  bool checkRateEconLimits(const WellEconProductionLimits& econ_production_limits,
278  const WellState& well_state) const;
279 
280  bool wellHasTHPConstraints() const;
281 
282  // Component fractions for each phase for the well
283  const std::vector<double>& compFrac() const;
284 
285  double mostStrictBhpFromBhpLimits() const;
286 
287  // a tuple type for ratio limit check.
288  // first value indicates whether ratio limit is violated, when the ratio limit is not violated, the following three
289  // values should not be used.
290  // second value indicates whehter there is only one connection left.
291  // third value indicates the indx of the worst-offending connection.
292  // the last value indicates the extent of the violation for the worst-offending connection, which is defined by
293  // the ratio of the actual value to the value of the violated limit.
294  using RatioCheckTuple = std::tuple<bool, bool, int, double>;
295 
296  RatioCheckTuple checkMaxWaterCutLimit(const WellEconProductionLimits& econ_production_limits,
297  const WellState& well_state) const;
298 
299  RatioCheckTuple checkRatioEconLimits(const WellEconProductionLimits& econ_production_limits,
300  const WellState& well_state) const;
301 
302  };
303 
304 }
305 
306 #include "WellInterface_impl.hpp"
307 
308 #endif // OPM_WELLINTERFACE_HEADER_INCLUDED
A thin wrapper class that holds one VFPProdProperties and one VFPInjProperties object.
Definition: VFPProperties.hpp:37
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state) const =0
using the solution x to recover the solution xw for wells and applying xw to update Well State ...
a struct to collect information about the convergence checking
Definition: WellInterface.hpp:117
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
Definition: WellInterface.hpp:60
const std::string & name() const
Well name.
Definition: WellInterface_impl.hpp:138
virtual ~WellInterface()
Virutal destructor.
Definition: WellInterface.hpp:92
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
WellType wellType() const
Well type, INJECTOR or PRODUCER.
Definition: WellInterface_impl.hpp:150
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
WellInterface(const Well *well, const int time_step, const Wells *wells, const ModelParameters &param)
Constructor.
Definition: WellInterface_impl.hpp:28
const std::vector< int > & cells()
Well cells.
Definition: WellInterface.hpp:98
WellControls * wellControls() const
Well controls.
Definition: WellInterface_impl.hpp:162