My Project
StandardWell.hpp
1 /*
2  Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3  Copyright 2017 Statoil ASA.
4  Copyright 2016 - 2017 IRIS AS.
5 
6  This file is part of the Open Porous Media project (OPM).
7 
8  OPM is free software: you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation, either version 3 of the License, or
11  (at your option) any later version.
12 
13  OPM is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with OPM. If not, see <http://www.gnu.org/licenses/>.
20 */
21 
22 
23 #ifndef OPM_STANDARDWELL_HEADER_INCLUDED
24 #define OPM_STANDARDWELL_HEADER_INCLUDED
25 
26 #include <opm/simulators/timestepping/ConvergenceReport.hpp>
28 #include <opm/simulators/wells/StandardWellGeneric.hpp>
29 #include <opm/simulators/wells/VFPInjProperties.hpp>
30 #include <opm/simulators/wells/VFPProdProperties.hpp>
31 #include <opm/simulators/wells/WellInterface.hpp>
32 #include <opm/simulators/wells/WellProdIndexCalculator.hpp>
33 #include <opm/simulators/wells/ParallelWellInfo.hpp>
34 #include <opm/simulators/wells/GasLiftSingleWell.hpp>
35 #include <opm/simulators/wells/GasLiftGroupInfo.hpp>
36 
37 #include <opm/models/blackoil/blackoilpolymermodules.hh>
38 #include <opm/models/blackoil/blackoilsolventmodules.hh>
39 #include <opm/models/blackoil/blackoilextbomodules.hh>
40 #include <opm/models/blackoil/blackoilfoammodules.hh>
41 #include <opm/models/blackoil/blackoilbrinemodules.hh>
42 #include <opm/models/blackoil/blackoilmicpmodules.hh>
43 
44 #include <opm/material/densead/DynamicEvaluation.hpp>
45 #include <opm/parser/eclipse/EclipseState/Runspec.hpp>
46 #include <opm/parser/eclipse/EclipseState/Schedule/ScheduleTypes.hpp>
47 
48 #include <opm/simulators/wells/StandardWellEval.hpp>
49 
50 #include <dune/common/dynvector.hh>
51 #include <dune/common/dynmatrix.hh>
52 
53 #include <memory>
54 #include <optional>
55 #include <fmt/format.h>
56 
57 namespace Opm
58 {
59 
60  template<typename TypeTag>
61  class StandardWell : public WellInterface<TypeTag>
62  , public StandardWellEval<GetPropType<TypeTag, Properties::FluidSystem>,
63  GetPropType<TypeTag, Properties::Indices>,
64  GetPropType<TypeTag, Properties::Scalar>>
65  {
66 
67  public:
70  GetPropType<TypeTag, Properties::Indices>,
71  GetPropType<TypeTag, Properties::Scalar>>;
72 
73  // TODO: some functions working with AD variables handles only with values (double) without
74  // dealing with derivatives. It can be beneficial to make functions can work with either AD or scalar value.
75  // And also, it can also be beneficial to make these functions hanle different types of AD variables.
76  using typename Base::Simulator;
77  using typename Base::IntensiveQuantities;
78  using typename Base::FluidSystem;
79  using typename Base::MaterialLaw;
80  using typename Base::ModelParameters;
81  using typename Base::Indices;
82  using typename Base::RateConverterType;
83  using typename Base::SparseMatrixAdapter;
84  using typename Base::FluidState;
85  using typename Base::RateVector;
86  using typename Base::GasLiftSingleWell;
87  using typename Base::GLiftOptWells;
88  using typename Base::GLiftProdWells;
89  using typename Base::GLiftWellStateMap;
90  using typename Base::GLiftSyncGroups;
91 
92  using Base::has_solvent;
93  using Base::has_zFraction;
94  using Base::has_polymer;
95  using Base::has_polymermw;
96  using Base::has_foam;
97  using Base::has_brine;
98  using Base::has_energy;
99  using Base::has_micp;
100 
101  using PolymerModule = BlackOilPolymerModule<TypeTag>;
102  using FoamModule = BlackOilFoamModule<TypeTag>;
103  using BrineModule = BlackOilBrineModule<TypeTag>;
104 
105  // number of the conservation equations
106  static constexpr int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
107  // number of the well control equations
108  static constexpr int numWellControlEq = 1;
109  // number of the well equations that will always be used
110  // based on the solution strategy, there might be other well equations be introduced
111  static constexpr int numStaticWellEq = numWellConservationEq + numWellControlEq;
112 
113  // the index for Bhp in primary variables and also the index of well control equation
114  // they both will be the last one in their respective system.
115  // TODO: we should have indices for the well equations and well primary variables separately
116  static constexpr int Bhp = numStaticWellEq - numWellControlEq;
117 
118  using typename Base::Scalar;
119 
120 
121  using Base::name;
122  using Base::Water;
123  using Base::Oil;
124  using Base::Gas;
125 
126  using typename Base::BVector;
127 
128  using Eval = typename StdWellEval::Eval;
129  using EvalWell = typename StdWellEval::EvalWell;
130  using BVectorWell = typename StdWellEval::BVectorWell;
131 
132  StandardWell(const Well& well,
133  const ParallelWellInfo& pw_info,
134  const int time_step,
135  const ModelParameters& param,
136  const RateConverterType& rate_converter,
137  const int pvtRegionIdx,
138  const int num_components,
139  const int num_phases,
140  const int index_of_well,
141  const std::vector<PerforationData>& perf_data);
142 
143  virtual void init(const PhaseUsage* phase_usage_arg,
144  const std::vector<double>& depth_arg,
145  const double gravity_arg,
146  const int num_cells,
147  const std::vector< Scalar >& B_avg) override;
148 
149 
150  virtual void initPrimaryVariablesEvaluation() const override;
151 
153  virtual ConvergenceReport getWellConvergence(const WellState& well_state,
154  const std::vector<double>& B_avg,
155  DeferredLogger& deferred_logger,
156  const bool relax_tolerance = false) const override;
157 
159  virtual void apply(const BVector& x, BVector& Ax) const override;
161  virtual void apply(BVector& r) const override;
162 
165  virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
166  WellState& well_state,
167  DeferredLogger& deferred_logger) const override;
168 
170  virtual void computeWellPotentials(const Simulator& ebosSimulator,
171  const WellState& well_state,
172  std::vector<double>& well_potentials,
173  DeferredLogger& deferred_logger) /* const */ override;
174 
175  virtual void updatePrimaryVariables(const WellState& well_state, DeferredLogger& deferred_logger) const override;
176 
177  virtual void solveEqAndUpdateWellState(WellState& well_state, DeferredLogger& deferred_logger) override;
178 
179  virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
180  const WellState& well_state,
181  DeferredLogger& deferred_logger) override; // should be const?
182 
183  virtual void updateProductivityIndex(const Simulator& ebosSimulator,
184  const WellProdIndexCalculator& wellPICalc,
185  WellState& well_state,
186  DeferredLogger& deferred_logger) const override;
187 
188  virtual void addWellContributions(SparseMatrixAdapter& mat) const override;
189 
190  // iterate well equations with the specified control until converged
191  bool iterateWellEqWithControl(const Simulator& ebosSimulator,
192  const double dt,
193  const Well::InjectionControls& inj_controls,
194  const Well::ProductionControls& prod_controls,
195  WellState& well_state,
196  const GroupState& group_state,
197  DeferredLogger& deferred_logger) override;
198 
200  virtual bool jacobianContainsWellContributions() const override
201  {
202  return this->param_.matrix_add_well_contributions_;
203  }
204 
205  virtual void gasLiftOptimizationStage1 (
206  WellState& well_state,
207  const GroupState& group_state,
208  const Simulator& ebosSimulator,
209  DeferredLogger& deferred_logger,
210  GLiftProdWells &prod_wells,
211  GLiftOptWells &glift_wells,
212  GLiftWellStateMap &state_map,
213  GasLiftGroupInfo &group_info,
214  GLiftSyncGroups &sync_groups
215  ) const override;
216 
217  /* returns BHP */
218  double computeWellRatesAndBhpWithThpAlqProd(const Simulator &ebos_simulator,
219  const SummaryState &summary_state,
220  DeferredLogger &deferred_logger,
221  std::vector<double> &potentials,
222  double alq) const;
223 
224  void computeWellRatesWithThpAlqProd(
225  const Simulator &ebos_simulator,
226  const SummaryState &summary_state,
227  DeferredLogger &deferred_logger,
228  std::vector<double> &potentials,
229  double alq) const;
230 
231  // NOTE: Cannot be protected since it is used by GasLiftRuntime
232  std::optional<double> computeBhpAtThpLimitProdWithAlq(
233  const Simulator& ebos_simulator,
234  const SummaryState& summary_state,
235  DeferredLogger& deferred_logger,
236  double alq_value) const;
237 
238  // NOTE: Cannot be protected since it is used by GasLiftRuntime
239  void computeWellRatesWithBhp(
240  const Simulator& ebosSimulator,
241  const double& bhp,
242  std::vector<double>& well_flux,
243  DeferredLogger& deferred_logger) const;
244 
245  // NOTE: These cannot be protected since they are used by GasLiftRuntime
246  using Base::phaseUsage;
247  using Base::vfp_properties_;
248 
249  virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
250  DeferredLogger& deferred_logger) const override;
251 
252  void computeConnLevelProdInd(const FluidState& fs,
253  const std::function<double(const double)>& connPICalc,
254  const std::vector<EvalWell>& mobility,
255  double* connPI) const;
256 
257  void computeConnLevelInjInd(const typename StandardWell<TypeTag>::FluidState& fs,
258  const Phase preferred_phase,
259  const std::function<double(const double)>& connIICalc,
260  const std::vector<EvalWell>& mobility,
261  double* connII,
262  DeferredLogger& deferred_logger) const;
263 
264 
265  protected:
266  Eval getPerfCellPressure(const FluidState& fs) const;
267 
268  // xw = inv(D)*(rw - C*x)
269  void recoverSolutionWell(const BVector& x, BVectorWell& xw) const;
270 
271  // updating the well_state based on well solution dwells
272  void updateWellState(const BVectorWell& dwells,
273  WellState& well_state,
274  DeferredLogger& deferred_logger) const;
275 
276  // calculate the properties for the well connections
277  // to calulate the pressure difference between well connections.
278  void computePropertiesForWellConnectionPressures(const Simulator& ebosSimulator,
279  const WellState& well_state,
280  std::vector<double>& b_perf,
281  std::vector<double>& rsmax_perf,
282  std::vector<double>& rvmax_perf,
283  std::vector<double>& surf_dens_perf) const;
284 
285  void computeWellConnectionDensitesPressures(const Simulator& ebosSimulator,
286  const WellState& well_state,
287  const std::vector<double>& b_perf,
288  const std::vector<double>& rsmax_perf,
289  const std::vector<double>& rvmax_perf,
290  const std::vector<double>& surf_dens_perf);
291 
292  void computeWellConnectionPressures(const Simulator& ebosSimulator,
293  const WellState& well_state);
294 
295  void computePerfRateEval(const IntensiveQuantities& intQuants,
296  const std::vector<EvalWell>& mob,
297  const EvalWell& bhp,
298  const double Tw,
299  const int perf,
300  const bool allow_cf,
301  std::vector<EvalWell>& cq_s,
302  double& perf_dis_gas_rate,
303  double& perf_vap_oil_rate,
304  DeferredLogger& deferred_logger) const;
305 
306  void computePerfRateScalar(const IntensiveQuantities& intQuants,
307  const std::vector<Scalar>& mob,
308  const Scalar& bhp,
309  const double Tw,
310  const int perf,
311  const bool allow_cf,
312  std::vector<Scalar>& cq_s,
313  DeferredLogger& deferred_logger) const;
314 
315  template<class Value>
316  void computePerfRate(const std::vector<Value>& mob,
317  const Value& pressure,
318  const Value& bhp,
319  const Value& rs,
320  const Value& rv,
321  std::vector<Value>& b_perfcells_dense,
322  const double Tw,
323  const int perf,
324  const bool allow_cf,
325  const Value& skin_pressure,
326  const std::vector<Value>& cmix_s,
327  std::vector<Value>& cq_s,
328  double& perf_dis_gas_rate,
329  double& perf_vap_oil_rate,
330  DeferredLogger& deferred_logger) const;
331 
332  void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
333  const double& bhp,
334  std::vector<double>& well_flux,
335  DeferredLogger& deferred_logger) const;
336 
337  std::vector<double> computeWellPotentialWithTHP(
338  const Simulator& ebosSimulator,
339  DeferredLogger& deferred_logger,
340  const WellState &well_state) const;
341 
342 
343  virtual double getRefDensity() const override;
344 
345  // get the mobility for specific perforation
346  void getMobilityEval(const Simulator& ebosSimulator,
347  const int perf,
348  std::vector<EvalWell>& mob,
349  DeferredLogger& deferred_logger) const;
350 
351  // get the mobility for specific perforation
352  void getMobilityScalar(const Simulator& ebosSimulator,
353  const int perf,
354  std::vector<Scalar>& mob,
355  DeferredLogger& deferred_logger) const;
356 
357 
358  void updateWaterMobilityWithPolymer(const Simulator& ebos_simulator,
359  const int perf,
360  std::vector<EvalWell>& mob_water,
361  DeferredLogger& deferred_logger) const;
362 
363  void updatePrimaryVariablesNewton(const BVectorWell& dwells,
364  const WellState& well_state,
365  DeferredLogger& deferred_logger) const;
366 
367  // update extra primary vriables if there are any
368  void updateExtraPrimaryVariables(const BVectorWell& dwells) const;
369 
370 
371  void updateWellStateFromPrimaryVariables(WellState& well_state, DeferredLogger& deferred_logger) const;
372 
373  virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
374  const double dt,
375  const Well::InjectionControls& inj_controls,
376  const Well::ProductionControls& prod_controls,
377  WellState& well_state,
378  const GroupState& group_state,
379  DeferredLogger& deferred_logger) override;
380 
381  void assembleWellEqWithoutIterationImpl(const Simulator& ebosSimulator,
382  const double dt,
383  WellState& well_state,
384  const GroupState& group_state,
385  DeferredLogger& deferred_logger);
386 
387  void calculateSinglePerf(const Simulator& ebosSimulator,
388  const int perf,
389  WellState& well_state,
390  std::vector<RateVector>& connectionRates,
391  std::vector<EvalWell>& cq_s,
392  EvalWell& water_flux_s,
393  EvalWell& cq_s_zfrac_effective,
394  DeferredLogger& deferred_logger) const;
395 
396  // check whether the well is operable under BHP limit with current reservoir condition
397  virtual void checkOperabilityUnderBHPLimitProducer(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) override;
398 
399  // check whether the well is operable under THP limit with current reservoir condition
400  virtual void checkOperabilityUnderTHPLimitProducer(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) override;
401 
402  // updating the inflow based on the current reservoir condition
403  virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const override;
404 
405  // for a well, when all drawdown are in the wrong direction, then this well will not
406  // be able to produce/inject .
407  bool allDrawDownWrongDirection(const Simulator& ebos_simulator) const;
408 
409  // whether the well can produce / inject based on the current well state (bhp)
410  bool canProduceInjectWithCurrentBhp(const Simulator& ebos_simulator,
411  const WellState& well_state,
412  DeferredLogger& deferred_logger);
413 
414  // turn on crossflow to avoid singular well equations
415  // when the well is banned from cross-flow and the BHP is not properly initialized,
416  // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
417  // well rates, it can cause problem for THP calculation
418  // TODO: looking for better alternative to avoid wrong-signed well rates
419  bool openCrossFlowAvoidSingularity(const Simulator& ebos_simulator) const;
420 
421  // calculate the skin pressure based on water velocity, throughput and polymer concentration.
422  // throughput is used to describe the formation damage during water/polymer injection.
423  // calculated skin pressure will be applied to the drawdown during perforation rate calculation
424  // to handle the effect from formation damage.
425  EvalWell pskin(const double throuhgput,
426  const EvalWell& water_velocity,
427  const EvalWell& poly_inj_conc,
428  DeferredLogger& deferred_logger) const;
429 
430  // calculate the skin pressure based on water velocity, throughput during water injection.
431  EvalWell pskinwater(const double throughput,
432  const EvalWell& water_velocity,
433  DeferredLogger& deferred_logger) const;
434 
435  // calculate the injecting polymer molecular weight based on the througput and water velocity
436  EvalWell wpolymermw(const double throughput,
437  const EvalWell& water_velocity,
438  DeferredLogger& deferred_logger) const;
439 
440  // modify the water rate for polymer injectivity study
441  void handleInjectivityRate(const Simulator& ebosSimulator,
442  const int perf,
443  std::vector<EvalWell>& cq_s) const;
444 
445  // handle the extra equations for polymer injectivity study
446  void handleInjectivityEquations(const Simulator& ebosSimulator,
447  const WellState& well_state,
448  const int perf,
449  const EvalWell& water_flux_s,
450  DeferredLogger& deferred_logger);
451 
452  virtual void updateWaterThroughput(const double dt, WellState& well_state) const override;
453 
454  // checking convergence of extra equations, if there are any
455  void checkConvergenceExtraEqs(const std::vector<double>& res,
456  ConvergenceReport& report) const;
457 
458  // updating the connectionRates_ related polymer molecular weight
459  void updateConnectionRatePolyMW(const EvalWell& cq_s_poly,
460  const IntensiveQuantities& int_quants,
461  const WellState& well_state,
462  const int perf,
463  std::vector<RateVector>& connectionRates,
464  DeferredLogger& deferred_logger) const;
465 
466 
467  std::optional<double> computeBhpAtThpLimitProd(const WellState& well_state,
468  const Simulator& ebos_simulator,
469  const SummaryState& summary_state,
470  DeferredLogger& deferred_logger) const;
471 
472  std::optional<double> computeBhpAtThpLimitInj(const Simulator& ebos_simulator,
473  const SummaryState& summary_state,
474  DeferredLogger& deferred_logger) const;
475 
476  };
477 
478 }
479 
480 #include "StandardWell_impl.hpp"
481 
482 #endif // OPM_STANDARDWELL_HEADER_INCLUDED
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
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
Definition: StandardWellEval.hpp:47
Definition: StandardWell.hpp:65
virtual bool jacobianContainsWellContributions() const override
Wether the Jacobian will also have well contributions in it.
Definition: StandardWell.hpp:200
virtual std::vector< double > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const override
Compute well rates based on current reservoir conditions and well variables.
Definition: StandardWell_impl.hpp:2426
virtual ConvergenceReport getWellConvergence(const WellState &well_state, const std::vector< double > &B_avg, DeferredLogger &deferred_logger, const bool relax_tolerance=false) const override
check whether the well equations get converged for this well
Definition: StandardWell_impl.hpp:1344
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state, DeferredLogger &deferred_logger) const override
using the solution x to recover the solution xw for wells and applying xw to update Well State
Definition: StandardWell_impl.hpp:1642
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1578
virtual void computeWellPotentials(const Simulator &ebosSimulator, const WellState &well_state, std::vector< double > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition: StandardWell_impl.hpp:1858
const std::string & name() const
Well name.
Definition: WellInterfaceGeneric.cpp:108
Definition: WellInterface.hpp:71
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
bool matrix_add_well_contributions_
Whether to add influences of wells between cells to the matrix and preconditioner matrix.
Definition: BlackoilModelParametersEbos.hpp:414
Definition: BlackoilPhases.hpp:45