My Project
WellInterface.hpp
1 /*
2  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3  Copyright 2017 Statoil ASA.
4  Copyright 2017 IRIS
5  Copyright 2019 Norce
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_WELLINTERFACE_HEADER_INCLUDED
25 #define OPM_WELLINTERFACE_HEADER_INCLUDED
26 
27 #include <opm/common/OpmLog/OpmLog.hpp>
28 #include <opm/common/ErrorMacros.hpp>
29 #include <opm/common/Exceptions.hpp>
30 
31 #include <opm/parser/eclipse/EclipseState/Schedule/Well/Well.hpp>
32 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WellTestState.hpp>
33 
34 #include <opm/core/props/BlackoilPhases.hpp>
35 
36 #include <opm/simulators/wells/WellProdIndexCalculator.hpp>
37 #include <opm/simulators/wells/WellState.hpp>
38 // NOTE: GasLiftSingleWell.hpp includes StandardWell.hpp which includes ourself
39 // (WellInterface.hpp), so we need to forward declare GasLiftSingleWell
40 // for it to be defined in this file. Similar for BlackoilWellModel
41 namespace Opm {
42  template<typename TypeTag> class GasLiftSingleWell;
43  template<typename TypeTag> class BlackoilWellModel;
44 }
45 #include <opm/simulators/wells/GasLiftGroupInfo.hpp>
46 #include <opm/simulators/wells/GasLiftSingleWell.hpp>
47 #include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
48 #include <opm/simulators/wells/BlackoilWellModel.hpp>
49 #include <opm/simulators/flow/BlackoilModelParametersEbos.hpp>
50 
51 #include <opm/simulators/utils/DeferredLogger.hpp>
52 
53 #include<dune/common/fmatrix.hh>
54 #include<dune/istl/bcrsmatrix.hh>
55 #include<dune/istl/matrixmatrix.hh>
56 
57 #include <opm/material/densead/Evaluation.hpp>
58 
59 #include <opm/simulators/wells/WellInterfaceIndices.hpp>
60 
61 #include <cassert>
62 #include <vector>
63 
64 namespace Opm
65 {
66 
67 template<typename TypeTag>
68 class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
69  GetPropType<TypeTag, Properties::Indices>,
70  GetPropType<TypeTag, Properties::Scalar>>
71 {
72 public:
73 
75 
76  using Grid = GetPropType<TypeTag, Properties::Grid>;
77  using Simulator = GetPropType<TypeTag, Properties::Simulator>;
78  using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
79  using Indices = GetPropType<TypeTag, Properties::Indices>;
80  using IntensiveQuantities = GetPropType<TypeTag, Properties::IntensiveQuantities>;
81  using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
82  using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
83  using RateVector = GetPropType<TypeTag, Properties::RateVector>;
85  using GLiftOptWells = typename BlackoilWellModel<TypeTag>::GLiftOptWells;
86  using GLiftProdWells = typename BlackoilWellModel<TypeTag>::GLiftProdWells;
87  using GLiftWellStateMap =
88  typename BlackoilWellModel<TypeTag>::GLiftWellStateMap;
89  using GLiftSyncGroups = typename GasLiftSingleWellGeneric::GLiftSyncGroups;
90 
91  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
92 
93  using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
94  using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
95  using BVector = Dune::BlockVector<VectorBlockType>;
96  using Eval = DenseAd::Evaluation<Scalar, /*size=*/Indices::numEq>;
97 
98  using RateConverterType =
100 
104  using RatioLimitCheckReport = typename WellInterfaceFluidSystem<FluidSystem>::RatioLimitCheckReport;
105 
106  static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
107  static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
108  static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
109  static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
110  static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
111  // flag for polymer molecular weight related
112  static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
113  static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
114  static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
115  static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
116 
117  // For the conversion between the surface volume rate and reservoir voidage rate
118  using FluidState = BlackOilFluidState<Eval,
119  FluidSystem,
120  has_temperature,
121  has_energy,
122  Indices::compositionSwitchIdx >= 0,
123  has_brine,
124  Indices::numPhases >;
126  WellInterface(const Well& well,
127  const ParallelWellInfo& pw_info,
128  const int time_step,
129  const ModelParameters& param,
130  const RateConverterType& rate_converter,
131  const int pvtRegionIdx,
132  const int num_components,
133  const int num_phases,
134  const int index_of_well,
135  const std::vector<PerforationData>& perf_data);
136 
138  virtual ~WellInterface() = default;
139 
140  virtual void init(const PhaseUsage* phase_usage_arg,
141  const std::vector<double>& depth_arg,
142  const double gravity_arg,
143  const int num_cells,
144  const std::vector< Scalar >& B_avg);
145 
146  virtual void initPrimaryVariablesEvaluation() const = 0;
147 
148  virtual ConvergenceReport getWellConvergence(const WellState& well_state, const std::vector<double>& B_avg, DeferredLogger& deferred_logger, const bool relax_tolerance) const = 0;
149 
150  virtual void solveEqAndUpdateWellState(WellState& well_state, DeferredLogger& deferred_logger) = 0;
151 
152  void assembleWellEq(const Simulator& ebosSimulator,
153  const double dt,
154  WellState& well_state,
155  const GroupState& group_state,
156  DeferredLogger& deferred_logger);
157 
158  virtual void gasLiftOptimizationStage1 (
159  WellState& well_state,
160  const GroupState& group_state,
161  const Simulator& ebosSimulator,
162  DeferredLogger& deferred_logger,
163  GLiftProdWells& prod_wells,
164  GLiftOptWells& glift_wells,
165  GLiftWellStateMap& state_map,
166  GasLiftGroupInfo &group_info,
167  GLiftSyncGroups &sync_groups
168  ) const = 0;
169 
172  virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
173  WellState& well_state,
174  DeferredLogger& deferred_logger) const = 0;
175 
177  virtual void apply(const BVector& x, BVector& Ax) const = 0;
178 
180  virtual void apply(BVector& r) const = 0;
181 
182  // TODO: before we decide to put more information under mutable, this function is not const
183  virtual void computeWellPotentials(const Simulator& ebosSimulator,
184  const WellState& well_state,
185  std::vector<double>& well_potentials,
186  DeferredLogger& deferred_logger) = 0;
187 
188  virtual void updateWellStateWithTarget(const Simulator& ebos_simulator,
189  const GroupState& group_state,
190  WellState& well_state,
191  DeferredLogger& deferred_logger) const;
192 
193  enum class IndividualOrGroup { Individual, Group, Both };
194  bool updateWellControl(const Simulator& ebos_simulator,
195  const IndividualOrGroup iog,
196  WellState& well_state,
197  const GroupState& group_state,
198  DeferredLogger& deferred_logger) /* const */;
199 
200  virtual void updatePrimaryVariables(const WellState& well_state, DeferredLogger& deferred_logger) const = 0;
201 
202  virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
203  const WellState& well_state,
204  DeferredLogger& deferred_logger) = 0; // should be const?
205 
206  virtual void updateProductivityIndex(const Simulator& ebosSimulator,
207  const WellProdIndexCalculator& wellPICalc,
208  WellState& well_state,
209  DeferredLogger& deferred_logger) const = 0;
210 
213  {
214  return false;
215  }
216 
217  // Add well contributions to matrix
218  virtual void addWellContributions(SparseMatrixAdapter&) const = 0;
219 
220  void addCellRates(RateVector& rates, int cellIdx) const;
221 
222  Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
223 
224  template <class EvalWell>
225  Eval restrictEval(const EvalWell& in) const
226  {
227  Eval out = 0.0;
228  out.setValue(in.value());
229  for (int eqIdx = 0; eqIdx < Indices::numEq; ++eqIdx) {
230  out.setDerivative(eqIdx, in.derivative(eqIdx));
231  }
232  return out;
233  }
234 
235  // TODO: theoretically, it should be a const function
236  // Simulator is not const is because that assembleWellEq is non-const Simulator
237  void wellTesting(const Simulator& simulator,
238  const double simulation_time,
239  /* const */ WellState& well_state, const GroupState& group_state, WellTestState& welltest_state,
240  DeferredLogger& deferred_logger);
241 
242  void checkWellOperability(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger);
243 
244  // check whether the well is operable under the current reservoir condition
245  // mostly related to BHP limit and THP limit
246  void updateWellOperability(const Simulator& ebos_simulator,
247  const WellState& well_state,
248  DeferredLogger& deferred_logger);
249 
250 
251  // update perforation water throughput based on solved water rate
252  virtual void updateWaterThroughput(const double dt, WellState& well_state) const = 0;
253 
256  virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
257  DeferredLogger& deferred_logger) const = 0;
258 
262  void updateWellStateRates(const Simulator& ebosSimulator,
263  WellState& well_state,
264  DeferredLogger& deferred_logger) const;
265 
266  void solveWellEquation(const Simulator& ebosSimulator,
267  WellState& well_state,
268  const GroupState& group_state,
269  DeferredLogger& deferred_logger);
270 
271 protected:
272 
273  // simulation parameters
274  const ModelParameters& param_;
275 
276  std::vector<RateVector> connectionRates_;
277 
278  std::vector< Scalar > B_avg_;
279 
280  bool changed_to_stopped_this_step_ = false;
281 
282  double wpolymer() const;
283 
284  double wfoam() const;
285 
286  double wsalt() const;
287 
288  double wmicrobes() const;
289 
290  double woxygen() const;
291 
292  double wurea() const;
293 
294  virtual double getRefDensity() const = 0;
295 
296  // Component fractions for each phase for the well
297  const std::vector<double>& compFrac() const;
298 
299  std::vector<double> initialWellRateFractions(const Simulator& ebosSimulator, const WellState& well_state) const;
300 
301  // check whether the well is operable under BHP limit with current reservoir condition
302  virtual void checkOperabilityUnderBHPLimitProducer(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) =0;
303 
304  // check whether the well is operable under THP limit with current reservoir condition
305  virtual void checkOperabilityUnderTHPLimitProducer(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) =0;
306 
307  virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const=0;
308 
309  virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
310  const double dt,
311  const Well::InjectionControls& inj_controls,
312  const Well::ProductionControls& prod_controls,
313  WellState& well_state,
314  const GroupState& group_state,
315  DeferredLogger& deferred_logger) = 0;
316 
317  // iterate well equations with the specified control until converged
318  virtual bool iterateWellEqWithControl(const Simulator& ebosSimulator,
319  const double dt,
320  const Well::InjectionControls& inj_controls,
321  const Well::ProductionControls& prod_controls,
322  WellState& well_state,
323  const GroupState& group_state,
324  DeferredLogger& deferred_logger) = 0;
325 
326  bool iterateWellEquations(const Simulator& ebosSimulator,
327  const double dt,
328  WellState& well_state,
329  const GroupState& group_state,
330  DeferredLogger& deferred_logger);
331 
332  bool solveWellForTesting(const Simulator& ebosSimulator, WellState& well_state, const GroupState& group_state,
333  DeferredLogger& deferred_logger);
334 
335  bool shutUnsolvableWells() const;
336 };
337 
338 }
339 
340 #include "WellInterface_impl.hpp"
341 
342 #endif // OPM_WELLINTERFACE_HEADER_INCLUDED
Represents the convergence status of the whole simulator, to make it possible to query and store the ...
Definition: ConvergenceReport.hpp:36
Definition: DeferredLogger.hpp:57
Definition: GasLiftGroupInfo.hpp:46
Definition: GasLiftSingleWell.hpp:46
Definition: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:252
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:68
Definition: WellInterfaceFluidSystem.hpp:46
Definition: WellInterfaceIndices.hpp:35
Definition: WellInterface.hpp:71
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
void updateWellStateRates(const Simulator &ebosSimulator, WellState &well_state, DeferredLogger &deferred_logger) const
Modify the well_state's rates if there is only one nonzero rate.
Definition: WellInterface_impl.hpp:991
virtual std::vector< double > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const =0
Compute well rates based on current reservoir conditions and well variables.
virtual void apply(BVector &r) const =0
r = r - C D^-1 Rw
virtual ~WellInterface()=default
Virtual destructor.
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state, DeferredLogger &deferred_logger) const =0
using the solution x to recover the solution xw for wells and applying xw to update Well State
WellInterface(const Well &well, const ParallelWellInfo &pw_info, const int time_step, const ModelParameters &param, const RateConverterType &rate_converter, const int pvtRegionIdx, const int num_components, const int num_phases, const int index_of_well, const std::vector< PerforationData > &perf_data)
Constructor.
Definition: WellInterface_impl.hpp:35
virtual bool jacobianContainsWellContributions() const
Wether the Jacobian will also have well contributions in it.
Definition: WellInterface.hpp:212
Collect per-connection static information to enable calculating connection-level or well-level produc...
Definition: WellProdIndexCalculator.hpp:36
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:26
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParametersEbos.hpp:327
Definition: BlackoilPhases.hpp:45
Definition: WellInterfaceFluidSystem.hpp:126