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 GetPropType<TypeTag, Properties::Scalar>>
60 {
61
62 public:
67
68 // TODO: some functions working with AD variables handles only with values (double) without
69 // dealing with derivatives. It can be beneficial to make functions can work with either AD or scalar value.
70 // And also, it can also be beneficial to make these functions hanle different types of AD variables.
77 using typename Base::RateConverterType;
79 using typename Base::FluidState;
81
82 using Base::has_solvent;
83 using Base::has_zFraction;
84 using Base::has_polymer;
85 using Base::has_polymermw;
86 using Base::has_foam;
87 using Base::has_brine;
88 using Base::has_energy;
89 using Base::has_micp;
90
94 using typename Base::PressureMatrix;
95
96 // number of the conservation equations
97 static constexpr int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
98 // number of the well control equations
99 static constexpr int numWellControlEq = 1;
100 // number of the well equations that will always be used
101 // based on the solution strategy, there might be other well equations be introduced
102 static constexpr int numStaticWellEq = numWellConservationEq + numWellControlEq;
103
104 // the index for Bhp in primary variables and also the index of well control equation
105 // they both will be the last one in their respective system.
106 // TODO: we should have indices for the well equations and well primary variables separately
107 static constexpr int Bhp = numStaticWellEq - numWellControlEq;
108
109 using StdWellEval::WQTotal;
110
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,
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,
141 const bool changed_to_open_this_step) override;
142
143
144 void initPrimaryVariablesEvaluation() override;
145
147 virtual ConvergenceReport getWellConvergence(const WellState& well_state,
148 const std::vector<double>& B_avg,
150 const bool relax_tolerance = false) const override;
151
153 virtual void apply(const BVector& x, BVector& Ax) const override;
155 virtual void apply(BVector& r) const override;
156
160 const BVector& x,
161 WellState& well_state,
163
165 virtual void computeWellPotentials(const Simulator& ebosSimulator,
166 const WellState& well_state,
167 std::vector<double>& well_potentials,
168 DeferredLogger& deferred_logger) /* const */ override;
169
170 void updatePrimaryVariables(const SummaryState& summary_state,
171 const WellState& well_state,
173
174 void solveEqAndUpdateWellState(const SummaryState& summary_state,
175 WellState& well_state,
177
178 virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
179 const WellState& well_state,
180 DeferredLogger& deferred_logger) override; // should be const?
181
182 virtual void updateProductivityIndex(const Simulator& ebosSimulator,
184 WellState& well_state,
185 DeferredLogger& deferred_logger) const override;
186
187 virtual void addWellContributions(SparseMatrixAdapter& mat) const override;
188
189 virtual void addWellPressureEquations(PressureMatrix& mat,
190 const BVector& x,
191 const int pressureVarIndex,
192 const bool use_well_weights,
193 const WellState& well_state) const override;
194
195 // iterate well equations with the specified control until converged
196 bool iterateWellEqWithControl(const Simulator& ebosSimulator,
197 const double dt,
198 const Well::InjectionControls& inj_controls,
199 const Well::ProductionControls& prod_controls,
200 WellState& well_state,
201 const GroupState& group_state,
203
205 virtual bool jacobianContainsWellContributions() const override
206 {
207 return this->param_.matrix_add_well_contributions_;
208 }
209
210 /* returns BHP */
211 double computeWellRatesAndBhpWithThpAlqProd(const Simulator &ebos_simulator,
214 std::vector<double> &potentials,
215 double alq) const;
216
217 void computeWellRatesWithThpAlqProd(
218 const Simulator &ebos_simulator,
221 std::vector<double> &potentials,
222 double alq) const;
223
224 std::optional<double> computeBhpAtThpLimitProdWithAlq(
225 const Simulator& ebos_simulator,
227 const double alq_value,
228 DeferredLogger& deferred_logger) const override;
229
230 virtual void computeWellRatesWithBhp(
231 const Simulator& ebosSimulator,
232 const double& bhp,
233 std::vector<double>& well_flux,
234 DeferredLogger& deferred_logger) const override;
235
236 // NOTE: These cannot be protected since they are used by GasLiftRuntime
237 using Base::phaseUsage;
238 using Base::vfp_properties_;
239
240 virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
241 DeferredLogger& deferred_logger) const override;
242
243 void computeConnLevelProdInd(const FluidState& fs,
244 const std::function<double(const double)>& connPICalc,
245 const std::vector<EvalWell>& mobility,
246 double* connPI) const;
247
248 void computeConnLevelInjInd(const typename StandardWell<TypeTag>::FluidState& fs,
250 const std::function<double(const double)>& connIICalc,
251 const std::vector<EvalWell>& mobility,
252 double* connII,
254
255
256 protected:
257 bool regularize_;
258
259 // updating the well_state based on well solution dwells
260 void updateWellState(const SummaryState& summary_state,
261 const BVectorWell& dwells,
262 WellState& well_state,
264
265 // calculate the properties for the well connections
266 // to calulate the pressure difference between well connections.
267 void computePropertiesForWellConnectionPressures(const Simulator& ebosSimulator,
268 const WellState& well_state,
269 std::vector<double>& b_perf,
270 std::vector<double>& rsmax_perf,
271 std::vector<double>& rvmax_perf,
272 std::vector<double>& rvwmax_perf,
273 std::vector<double>& rswmax_perf,
274 std::vector<double>& surf_dens_perf) const;
275
276 void computeWellConnectionDensitesPressures(const Simulator& ebosSimulator,
277 const WellState& well_state,
278 const std::vector<double>& b_perf,
279 const std::vector<double>& rsmax_perf,
280 const std::vector<double>& rvmax_perf,
281 const std::vector<double>& rvwmax_perf,
282 const std::vector<double>& rswmax_perf,
283 const std::vector<double>& surf_dens_perf,
285
286 void computeWellConnectionPressures(const Simulator& ebosSimulator,
287 const WellState& well_state,
289
290 void computePerfRateEval(const IntensiveQuantities& intQuants,
291 const std::vector<EvalWell>& mob,
292 const EvalWell& bhp,
293 const double Tw,
294 const int perf,
295 const bool allow_cf,
296 std::vector<EvalWell>& cq_s,
297 double& perf_dis_gas_rate,
299 double& perf_vap_oil_rate,
300 double& perf_vap_wat_rate,
302
303 void computePerfRateScalar(const IntensiveQuantities& intQuants,
304 const std::vector<Scalar>& mob,
305 const Scalar& bhp,
306 const double Tw,
307 const int perf,
308 const bool allow_cf,
309 std::vector<Scalar>& cq_s,
311
312 template<class Value>
313 void computePerfRate(const std::vector<Value>& mob,
314 const Value& pressure,
315 const Value& bhp,
316 const Value& rs,
317 const Value& rv,
318 const Value& rvw,
319 const Value& rsw,
320 std::vector<Value>& b_perfcells_dense,
321 const double Tw,
322 const int perf,
323 const bool allow_cf,
324 const Value& skin_pressure,
325 const std::vector<Value>& cmix_s,
326 std::vector<Value>& cq_s,
327 double& perf_dis_gas_rate,
329 double& perf_vap_oil_rate,
330 double& perf_vap_wat_rate,
332
333 void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
334 const double& bhp,
335 std::vector<double>& well_flux,
337
338 std::vector<double> computeWellPotentialWithTHP(
339 const Simulator& ebosSimulator,
341 const WellState &well_state) const;
342
343 bool updateWellStateWithTHPTargetProd(const Simulator& ebos_simulator,
344 WellState& well_state,
345 DeferredLogger& deferred_logger) const override;
346
347 virtual double getRefDensity() const override;
348
349 // get the mobility for specific perforation
350 void getMobilityEval(const Simulator& ebosSimulator,
351 const int perf,
352 std::vector<EvalWell>& mob,
354
355 // get the mobility for specific perforation
356 void getMobilityScalar(const Simulator& ebosSimulator,
357 const int perf,
358 std::vector<Scalar>& mob,
360
361
362 void updateWaterMobilityWithPolymer(const Simulator& ebos_simulator,
363 const int perf,
364 std::vector<EvalWell>& mob_water,
366
367 void updatePrimaryVariablesNewton(const BVectorWell& dwells,
368 const bool stop_or_zero_rate_target,
370
371 void updateWellStateFromPrimaryVariables(const bool stop_or_zero_rate_target,
372 WellState& well_state,
374
375 virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
376 const double dt,
377 const Well::InjectionControls& inj_controls,
378 const Well::ProductionControls& prod_controls,
379 WellState& well_state,
380 const GroupState& group_state,
382
383 void assembleWellEqWithoutIterationImpl(const Simulator& ebosSimulator,
384 const double dt,
385 const Well::InjectionControls& inj_controls,
386 const Well::ProductionControls& prod_controls,
387 WellState& well_state,
388 const GroupState& group_state,
390
391 void calculateSinglePerf(const Simulator& ebosSimulator,
392 const int perf,
393 WellState& well_state,
394 std::vector<RateVector>& connectionRates,
395 std::vector<EvalWell>& cq_s,
396 EvalWell& water_flux_s,
397 EvalWell& cq_s_zfrac_effective,
399
400 // check whether the well is operable under BHP limit with current reservoir condition
401 virtual void checkOperabilityUnderBHPLimit(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) override;
402
403 // check whether the well is operable under THP limit with current reservoir condition
404 virtual void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) override;
405
406 // updating the inflow based on the current reservoir condition
407 virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const override;
408
409 // for a well, when all drawdown are in the wrong direction, then this well will not
410 // be able to produce/inject .
411 bool allDrawDownWrongDirection(const Simulator& ebos_simulator) const;
412
413 // whether the well can produce / inject based on the current well state (bhp)
414 bool canProduceInjectWithCurrentBhp(const Simulator& ebos_simulator,
415 const WellState& well_state,
417
418 // turn on crossflow to avoid singular well equations
419 // when the well is banned from cross-flow and the BHP is not properly initialized,
420 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
421 // well rates, it can cause problem for THP calculation
422 // TODO: looking for better alternative to avoid wrong-signed well rates
423 bool openCrossFlowAvoidSingularity(const Simulator& ebos_simulator) const;
424
425 // calculate the skin pressure based on water velocity, throughput and polymer concentration.
426 // throughput is used to describe the formation damage during water/polymer injection.
427 // calculated skin pressure will be applied to the drawdown during perforation rate calculation
428 // to handle the effect from formation damage.
429 EvalWell pskin(const double throuhgput,
430 const EvalWell& water_velocity,
431 const EvalWell& poly_inj_conc,
433
434 // calculate the skin pressure based on water velocity, throughput during water injection.
435 EvalWell pskinwater(const double throughput,
436 const EvalWell& water_velocity,
438
439 // calculate the injecting polymer molecular weight based on the througput and water velocity
440 EvalWell wpolymermw(const double throughput,
441 const EvalWell& water_velocity,
443
444 // modify the water rate for polymer injectivity study
445 void handleInjectivityRate(const Simulator& ebosSimulator,
446 const int perf,
447 std::vector<EvalWell>& cq_s) const;
448
449 // handle the extra equations for polymer injectivity study
450 void handleInjectivityEquations(const Simulator& ebosSimulator,
451 const WellState& well_state,
452 const int perf,
453 const EvalWell& water_flux_s,
455
456 virtual void updateWaterThroughput(const double dt, WellState& well_state) const override;
457
458 // checking convergence of extra equations, if there are any
459 void checkConvergenceExtraEqs(const std::vector<double>& res,
460 ConvergenceReport& report) const;
461
462 // updating the connectionRates_ related polymer molecular weight
463 void updateConnectionRatePolyMW(const EvalWell& cq_s_poly,
464 const IntensiveQuantities& int_quants,
465 const WellState& well_state,
466 const int perf,
467 std::vector<RateVector>& connectionRates,
469
470
471 std::optional<double> computeBhpAtThpLimitProd(const WellState& well_state,
472 const Simulator& ebos_simulator,
475
476 std::optional<double> computeBhpAtThpLimitInj(const Simulator& ebos_simulator,
479
480 };
481
482}
483
484#include "StandardWell_impl.hpp"
485
486#endif // OPM_STANDARDWELL_HEADER_INCLUDED
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Definition AquiferInterface.hpp:35
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:49
Definition StandardWell.hpp:60
virtual bool jacobianContainsWellContributions() const override
Wether the Jacobian will also have well contributions in it.
Definition StandardWell.hpp:205
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:1659
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:2544
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:1407
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition StandardWell_impl.hpp:1627
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:1885
const std::string & name() const
Well name.
Definition WellInterfaceGeneric.cpp:146
Definition WellInterface.hpp:74
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:60
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