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