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/input/eclipse/Schedule/Well/Well.hpp>
32 #include <opm/input/eclipse/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 #include <opm/simulators/timestepping/ConvergenceReport.hpp>
61 
62 #include <cassert>
63 #include <vector>
64 
65 namespace Opm
66 {
67 
68 template<typename TypeTag>
69 class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
70  GetPropType<TypeTag, Properties::Indices>,
71  GetPropType<TypeTag, Properties::Scalar>>
72 {
73 public:
74 
76 
77  using Grid = GetPropType<TypeTag, Properties::Grid>;
78  using Simulator = GetPropType<TypeTag, Properties::Simulator>;
79  using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
80  using Indices = GetPropType<TypeTag, Properties::Indices>;
81  using IntensiveQuantities = GetPropType<TypeTag, Properties::IntensiveQuantities>;
82  using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
83  using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
84  using RateVector = GetPropType<TypeTag, Properties::RateVector>;
86  using GLiftOptWells = typename BlackoilWellModel<TypeTag>::GLiftOptWells;
87  using GLiftProdWells = typename BlackoilWellModel<TypeTag>::GLiftProdWells;
88  using GLiftWellStateMap =
89  typename BlackoilWellModel<TypeTag>::GLiftWellStateMap;
90  using GLiftSyncGroups = typename GasLiftSingleWellGeneric::GLiftSyncGroups;
91 
92  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
93 
94  using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
95  using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
96  using BVector = Dune::BlockVector<VectorBlockType>;
97  using Eval = DenseAd::Evaluation<Scalar, /*size=*/Indices::numEq>;
98 
99  using RateConverterType =
101 
105  using RatioLimitCheckReport = typename WellInterfaceFluidSystem<FluidSystem>::RatioLimitCheckReport;
106 
107  static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
108  static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
109  static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
110  static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
111  static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
112  // flag for polymer molecular weight related
113  static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
114  static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
115  static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
116  static constexpr bool has_watVapor = getPropValue<TypeTag, Properties::EnableEvaporation>();
117  static constexpr bool has_saltPrecip = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
118  static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
119 
120  // For the conversion between the surface volume rate and reservoir voidage rate
121  using FluidState = BlackOilFluidState<Eval,
122  FluidSystem,
123  has_temperature,
124  has_energy,
125  Indices::compositionSwitchIdx >= 0,
126  has_watVapor,
127  has_brine,
128  has_saltPrecip,
129  Indices::numPhases >;
131  WellInterface(const Well& well,
132  const ParallelWellInfo& pw_info,
133  const int time_step,
134  const ModelParameters& param,
135  const RateConverterType& rate_converter,
136  const int pvtRegionIdx,
137  const int num_components,
138  const int num_phases,
139  const int index_of_well,
140  const std::vector<PerforationData>& perf_data);
141 
143  virtual ~WellInterface() = default;
144 
145  virtual void init(const PhaseUsage* phase_usage_arg,
146  const std::vector<double>& depth_arg,
147  const double gravity_arg,
148  const int num_cells,
149  const std::vector< Scalar >& B_avg);
150 
151  virtual void initPrimaryVariablesEvaluation() const = 0;
152 
153  virtual ConvergenceReport getWellConvergence(const WellState& well_state, const std::vector<double>& B_avg, DeferredLogger& deferred_logger, const bool relax_tolerance) const = 0;
154 
155  virtual void solveEqAndUpdateWellState(WellState& well_state, DeferredLogger& deferred_logger) = 0;
156 
157  void assembleWellEq(const Simulator& ebosSimulator,
158  const double dt,
159  WellState& well_state,
160  const GroupState& group_state,
161  DeferredLogger& deferred_logger);
162 
163  virtual void computeWellRatesWithBhp(
164  const Simulator& ebosSimulator,
165  const double& bhp,
166  std::vector<double>& well_flux,
167  DeferredLogger& deferred_logger
168  ) const = 0;
169 
170  virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
171  const Simulator& ebos_simulator,
172  const SummaryState& summary_state,
173  DeferredLogger& deferred_logger,
174  double alq_value
175  ) const = 0;
176 
179  virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
180  WellState& well_state,
181  DeferredLogger& deferred_logger) const = 0;
182 
184  virtual void apply(const BVector& x, BVector& Ax) const = 0;
185 
187  virtual void apply(BVector& r) const = 0;
188 
189  // TODO: before we decide to put more information under mutable, this function is not const
190  virtual void computeWellPotentials(const Simulator& ebosSimulator,
191  const WellState& well_state,
192  std::vector<double>& well_potentials,
193  DeferredLogger& deferred_logger) = 0;
194 
195  virtual void updateWellStateWithTarget(const Simulator& ebos_simulator,
196  const GroupState& group_state,
197  WellState& well_state,
198  DeferredLogger& deferred_logger) const;
199 
200  enum class IndividualOrGroup { Individual, Group, Both };
201  bool updateWellControl(const Simulator& ebos_simulator,
202  const IndividualOrGroup iog,
203  WellState& well_state,
204  const GroupState& group_state,
205  DeferredLogger& deferred_logger) /* const */;
206 
207  virtual void updatePrimaryVariables(const WellState& well_state, DeferredLogger& deferred_logger) const = 0;
208 
209  virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
210  const WellState& well_state,
211  DeferredLogger& deferred_logger) = 0; // should be const?
212 
213  virtual void updateProductivityIndex(const Simulator& ebosSimulator,
214  const WellProdIndexCalculator& wellPICalc,
215  WellState& well_state,
216  DeferredLogger& deferred_logger) const = 0;
217 
220  {
221  return false;
222  }
223 
224  // Add well contributions to matrix
225  virtual void addWellContributions(SparseMatrixAdapter&) const = 0;
226 
227  void addCellRates(RateVector& rates, int cellIdx) const;
228 
229  Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
230 
231  template <class EvalWell>
232  Eval restrictEval(const EvalWell& in) const
233  {
234  Eval out = 0.0;
235  out.setValue(in.value());
236  for (int eqIdx = 0; eqIdx < Indices::numEq; ++eqIdx) {
237  out.setDerivative(eqIdx, in.derivative(eqIdx));
238  }
239  return out;
240  }
241 
242  // TODO: theoretically, it should be a const function
243  // Simulator is not const is because that assembleWellEq is non-const Simulator
244  void wellTesting(const Simulator& simulator,
245  const double simulation_time,
246  /* const */ WellState& well_state, const GroupState& group_state, WellTestState& welltest_state,
247  DeferredLogger& deferred_logger);
248 
249  void checkWellOperability(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger);
250 
251  // check whether the well is operable under the current reservoir condition
252  // mostly related to BHP limit and THP limit
253  void updateWellOperability(const Simulator& ebos_simulator,
254  const WellState& well_state,
255  DeferredLogger& deferred_logger);
256 
257 
258  // update perforation water throughput based on solved water rate
259  virtual void updateWaterThroughput(const double dt, WellState& well_state) const = 0;
260 
263  virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
264  DeferredLogger& deferred_logger) const = 0;
265 
269  void updateWellStateRates(const Simulator& ebosSimulator,
270  WellState& well_state,
271  DeferredLogger& deferred_logger) const;
272 
273  void solveWellEquation(const Simulator& ebosSimulator,
274  WellState& well_state,
275  const GroupState& group_state,
276  DeferredLogger& deferred_logger);
277 
278 protected:
279 
280  // simulation parameters
281  const ModelParameters& param_;
282 
283  std::vector<RateVector> connectionRates_;
284 
285  std::vector< Scalar > B_avg_;
286 
287  bool changed_to_stopped_this_step_ = false;
288 
289  double wpolymer() const;
290 
291  double wfoam() const;
292 
293  double wsalt() const;
294 
295  double wmicrobes() const;
296 
297  double woxygen() const;
298 
299  double wurea() const;
300 
301  virtual double getRefDensity() const = 0;
302 
303  // Component fractions for each phase for the well
304  const std::vector<double>& compFrac() const;
305 
306  std::vector<double> initialWellRateFractions(const Simulator& ebosSimulator, const WellState& well_state) const;
307 
308  // check whether the well is operable under BHP limit with current reservoir condition
309  virtual void checkOperabilityUnderBHPLimit(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) =0;
310 
311  // check whether the well is operable under THP limit with current reservoir condition
312  virtual void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) =0;
313 
314  virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const=0;
315 
316  virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
317  const double dt,
318  const Well::InjectionControls& inj_controls,
319  const Well::ProductionControls& prod_controls,
320  WellState& well_state,
321  const GroupState& group_state,
322  DeferredLogger& deferred_logger) = 0;
323 
324  // iterate well equations with the specified control until converged
325  virtual bool iterateWellEqWithControl(const Simulator& ebosSimulator,
326  const double dt,
327  const Well::InjectionControls& inj_controls,
328  const Well::ProductionControls& prod_controls,
329  WellState& well_state,
330  const GroupState& group_state,
331  DeferredLogger& deferred_logger) = 0;
332 
333  bool iterateWellEquations(const Simulator& ebosSimulator,
334  const double dt,
335  WellState& well_state,
336  const GroupState& group_state,
337  DeferredLogger& deferred_logger);
338 
339  bool solveWellForTesting(const Simulator& ebosSimulator, WellState& well_state, const GroupState& group_state,
340  DeferredLogger& deferred_logger);
341 
342  Eval getPerfCellPressure(const FluidState& fs) const;
343 
344 
345 };
346 
347 }
348 
349 #include "WellInterface_impl.hpp"
350 
351 #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: GasLiftSingleWell.hpp:39
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:72
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:1003
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:219
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:27
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParametersEbos.hpp:327
Definition: BlackoilPhases.hpp:46
Definition: WellInterfaceFluidSystem.hpp:129