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