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/RatioCalculator.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
41
42#include <opm/material/densead/Evaluation.hpp>
43#include <opm/input/eclipse/Schedule/ScheduleTypes.hpp>
44
45#include <opm/simulators/wells/StandardWellEval.hpp>
46
47#include <dune/common/dynvector.hh>
48#include <dune/common/dynmatrix.hh>
49
50#include <memory>
51#include <optional>
52
53namespace Opm
54{
55
56 template<typename TypeTag>
57 class StandardWell : public WellInterface<TypeTag>
58 , public StandardWellEval<GetPropType<TypeTag, Properties::FluidSystem>,
59 GetPropType<TypeTag, Properties::Indices>>
60 {
61
62 public:
66
67 // TODO: some functions working with AD variables handles only with values (double) without
68 // dealing with derivatives. It can be beneficial to make functions can work with either AD or scalar value.
69 // And also, it can also be beneficial to make these functions hanle different types of AD variables.
70 using typename Base::Simulator;
71 using typename Base::IntensiveQuantities;
72 using typename Base::FluidSystem;
73 using typename Base::MaterialLaw;
74 using typename Base::ModelParameters;
75 using typename Base::Indices;
76 using typename Base::RateConverterType;
77 using typename Base::SparseMatrixAdapter;
78 using typename Base::FluidState;
79 using typename Base::RateVector;
80
81 using Base::has_solvent;
82 using Base::has_zFraction;
83 using Base::has_polymer;
84 using Base::has_polymermw;
85 using Base::has_foam;
86 using Base::has_brine;
87 using Base::has_energy;
88 using Base::has_micp;
89
93 using typename Base::PressureMatrix;
94
95 // number of the conservation equations
96 static constexpr int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
97 // number of the well control equations
98 static constexpr int numWellControlEq = 1;
99 // number of the well equations that will always be used
100 // based on the solution strategy, there might be other well equations be introduced
101 static constexpr int numStaticWellEq = numWellConservationEq + numWellControlEq;
102
103 // the index for Bhp in primary variables and also the index of well control equation
104 // they both will be the last one in their respective system.
105 // TODO: we should have indices for the well equations and well primary variables separately
106 static constexpr int Bhp = numStaticWellEq - numWellControlEq;
107
108 using StdWellEval::WQTotal;
109
110 using typename Base::Scalar;
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,
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<Scalar>>& perf_data);
133
134 virtual void init(const PhaseUsage* phase_usage_arg,
135 const std::vector<Scalar>& depth_arg,
136 const Scalar gravity_arg,
137 const std::vector<Scalar>& B_avg,
138 const bool changed_to_open_this_step) override;
139
141 virtual ConvergenceReport getWellConvergence(const Simulator& simulator,
142 const WellState<Scalar>& well_state,
143 const std::vector<Scalar>& B_avg,
145 const bool relax_tolerance) const override;
146
148 virtual void apply(const BVector& x, BVector& Ax) const override;
150 virtual void apply(BVector& r) const override;
151
154 void recoverWellSolutionAndUpdateWellState(const Simulator& simulator,
155 const BVector& x,
156 WellState<Scalar>& well_state,
158
160 void computeWellPotentials(const Simulator& simulator,
161 const WellState<Scalar>& well_state,
162 std::vector<Scalar>& well_potentials,
163 DeferredLogger& deferred_logger) /* const */ override;
164
165 void updatePrimaryVariables(const Simulator& simulator,
166 const WellState<Scalar>& well_state,
168
169 void solveEqAndUpdateWellState(const Simulator& simulator,
170 WellState<Scalar>& well_state,
172
173 void calculateExplicitQuantities(const Simulator& simulator,
174 const WellState<Scalar>& well_state,
175 DeferredLogger& deferred_logger) override; // should be const?
176
177 void updateProductivityIndex(const Simulator& simulator,
179 WellState<Scalar>& well_state,
180 DeferredLogger& deferred_logger) const override;
181
182 Scalar connectionDensity(const int globalConnIdx,
183 const int openConnIdx) const override;
184
185 void addWellContributions(SparseMatrixAdapter& mat) const override;
186
187 void addWellPressureEquations(PressureMatrix& mat,
188 const BVector& x,
189 const int pressureVarIndex,
190 const bool use_well_weights,
191 const WellState<Scalar>& well_state) const override;
192
193 // iterate well equations with the specified control until converged
194 bool iterateWellEqWithControl(const Simulator& simulator,
195 const double dt,
196 const Well::InjectionControls& inj_controls,
197 const Well::ProductionControls& prod_controls,
198 WellState<Scalar>& well_state,
199 const GroupState<Scalar>& group_state,
201
202 // iterate well equations including control switching
203 bool iterateWellEqWithSwitching(const Simulator& simulator,
204 const double dt,
205 const Well::InjectionControls& inj_controls,
206 const Well::ProductionControls& prod_controls,
207 WellState<Scalar>& well_state,
208 const GroupState<Scalar>& group_state,
210 const bool fixed_control = false,
211 const bool fixed_status = false) override;
212
215 {
216 return this->param_.matrix_add_well_contributions_;
217 }
218
219 /* returns BHP */
220 Scalar computeWellRatesAndBhpWithThpAlqProd(const Simulator& ebos_simulator,
223 std::vector<Scalar>& potentials,
224 Scalar alq) const;
225
226 void computeWellRatesWithThpAlqProd(const Simulator& ebos_simulator,
229 std::vector<Scalar>& potentials,
230 Scalar alq) const;
231
232 std::optional<Scalar>
233 computeBhpAtThpLimitProdWithAlq(const Simulator& ebos_simulator,
235 const Scalar alq_value,
237 bool iterate_if_no_solution) const override;
238
239 void updateIPRImplicit(const Simulator& simulator,
240 WellState<Scalar>& well_state,
242
243 void computeWellRatesWithBhp(const Simulator& ebosSimulator,
244 const Scalar& bhp,
245 std::vector<Scalar>& well_flux,
246 DeferredLogger& deferred_logger) const override;
247
248 // NOTE: These cannot be protected since they are used by GasLiftRuntime
249 using Base::phaseUsage;
250 using Base::vfp_properties_;
251
252 std::vector<Scalar>
254 DeferredLogger& deferred_logger) const override;
255
256 std::vector<Scalar> getPrimaryVars() const override;
257
258 int setPrimaryVars(typename std::vector<Scalar>::const_iterator it) override;
259
260 protected:
261 bool regularize_;
262
263 // updating the well_state based on well solution dwells
264 void updateWellState(const Simulator& simulator,
265 const BVectorWell& dwells,
266 WellState<Scalar>& well_state,
268
269 using WellConnectionProps = typename StdWellEval::StdWellConnections::Properties;
270
271 // Compute connection level PVT properties needed to calulate the
272 // pressure difference between well connections.
273 WellConnectionProps
274 computePropertiesForWellConnectionPressures(const Simulator& simulator,
275 const WellState<Scalar>& well_state) const;
276
277 void computeWellConnectionDensitesPressures(const Simulator& simulator,
278 const WellState<Scalar>& well_state,
279 const WellConnectionProps& props,
281
282 void computeWellConnectionPressures(const Simulator& simulator,
283 const WellState<Scalar>& well_state,
285
286 template<class Value>
287 void computePerfRate(const IntensiveQuantities& intQuants,
288 const std::vector<Value>& mob,
289 const Value& bhp,
290 const std::vector<Scalar>& Tw,
291 const int perf,
292 const bool allow_cf,
293 std::vector<Value>& cq_s,
296
297 template<class Value>
298 void computePerfRate(const std::vector<Value>& mob,
299 const Value& pressure,
300 const Value& bhp,
301 const Value& rs,
302 const Value& rv,
303 const Value& rvw,
304 const Value& rsw,
305 std::vector<Value>& b_perfcells_dense,
306 const std::vector<Scalar>& Tw,
307 const int perf,
308 const bool allow_cf,
309 const Value& skin_pressure,
310 const std::vector<Value>& cmix_s,
311 std::vector<Value>& cq_s,
314
315 void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
316 const Scalar& bhp,
317 std::vector<Scalar>& well_flux,
318 DeferredLogger& deferred_logger) const override;
319
320 std::vector<Scalar>
321 computeWellPotentialWithTHP(const Simulator& ebosSimulator,
323 const WellState<Scalar>& well_state) const;
324
325 bool computeWellPotentialsImplicit(const Simulator& ebos_simulator,
326 const WellState<Scalar>& well_state,
327 std::vector<Scalar>& well_potentials,
329
330 Scalar getRefDensity() const override;
331
332 // get the mobility for specific perforation
333 template<class Value>
334 void getMobility(const Simulator& simulator,
335 const int perf,
336 std::vector<Value>& mob,
338
339 void updateWaterMobilityWithPolymer(const Simulator& simulator,
340 const int perf,
341 std::vector<EvalWell>& mob_water,
343
344 void updatePrimaryVariablesNewton(const BVectorWell& dwells,
345 const bool stop_or_zero_rate_target,
347
348 void updateWellStateFromPrimaryVariables(WellState<Scalar>& well_state,
351
352 void assembleWellEqWithoutIteration(const Simulator& simulator,
353 const double dt,
354 const Well::InjectionControls& inj_controls,
355 const Well::ProductionControls& prod_controls,
356 WellState<Scalar>& well_state,
357 const GroupState<Scalar>& group_state,
359
360 void assembleWellEqWithoutIterationImpl(const Simulator& simulator,
361 const double dt,
362 const Well::InjectionControls& inj_controls,
363 const Well::ProductionControls& prod_controls,
364 WellState<Scalar>& well_state,
365 const GroupState<Scalar>& group_state,
367
368 void calculateSinglePerf(const Simulator& simulator,
369 const int perf,
370 WellState<Scalar>& well_state,
371 std::vector<RateVector>& connectionRates,
372 std::vector<EvalWell>& cq_s,
373 EvalWell& water_flux_s,
374 EvalWell& cq_s_zfrac_effective,
376
377 // check whether the well is operable under BHP limit with current reservoir condition
378 void checkOperabilityUnderBHPLimit(const WellState<Scalar>& well_state,
379 const Simulator& simulator,
381
382 // check whether the well is operable under THP limit with current reservoir condition
383 void checkOperabilityUnderTHPLimit(const Simulator& simulator,
384 const WellState<Scalar>& well_state,
386
387 // updating the inflow based on the current reservoir condition
388 void updateIPR(const Simulator& simulator,
389 DeferredLogger& deferred_logger) const override;
390
391 // for a well, when all drawdown are in the wrong direction, then this well will not
392 // be able to produce/inject .
393 bool allDrawDownWrongDirection(const Simulator& simulator) const;
394
395 // whether the well can produce / inject based on the current well state (bhp)
396 bool canProduceInjectWithCurrentBhp(const Simulator& simulator,
397 const WellState<Scalar>& well_state,
399
400 // turn on crossflow to avoid singular well equations
401 // when the well is banned from cross-flow and the BHP is not properly initialized,
402 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
403 // well rates, it can cause problem for THP calculation
404 // TODO: looking for better alternative to avoid wrong-signed well rates
405 bool openCrossFlowAvoidSingularity(const Simulator& simulator) const;
406
407 // calculate the skin pressure based on water velocity, throughput and polymer concentration.
408 // throughput is used to describe the formation damage during water/polymer injection.
409 // calculated skin pressure will be applied to the drawdown during perforation rate calculation
410 // to handle the effect from formation damage.
411 EvalWell pskin(const Scalar throughput,
412 const EvalWell& water_velocity,
413 const EvalWell& poly_inj_conc,
415
416 // calculate the skin pressure based on water velocity, throughput during water injection.
417 EvalWell pskinwater(const Scalar throughput,
418 const EvalWell& water_velocity,
420
421 // calculate the injecting polymer molecular weight based on the througput and water velocity
422 EvalWell wpolymermw(const Scalar throughput,
423 const EvalWell& water_velocity,
425
426 // modify the water rate for polymer injectivity study
427 void handleInjectivityRate(const Simulator& simulator,
428 const int perf,
429 std::vector<EvalWell>& cq_s) const;
430
431 // handle the extra equations for polymer injectivity study
432 void handleInjectivityEquations(const Simulator& simulator,
433 const WellState<Scalar>& well_state,
434 const int perf,
435 const EvalWell& water_flux_s,
437
438 void updateWaterThroughput(const double dt,
439 WellState<Scalar>& well_state) const override;
440
441 // checking convergence of extra equations, if there are any
442 void checkConvergenceExtraEqs(const std::vector<Scalar>& res,
443 ConvergenceReport& report) const;
444
445 // updating the connectionRates_ related polymer molecular weight
446 void updateConnectionRatePolyMW(const EvalWell& cq_s_poly,
447 const IntensiveQuantities& int_quants,
448 const WellState<Scalar>& well_state,
449 const int perf,
450 std::vector<RateVector>& connectionRates,
452
453 std::optional<Scalar>
454 computeBhpAtThpLimitProd(const WellState<Scalar>& well_state,
455 const Simulator& simulator,
458
459 std::optional<Scalar>
460 computeBhpAtThpLimitInj(const Simulator& simulator,
463
464 private:
465 Eval connectionRateEnergy(const Scalar maxOilSaturation,
466 const std::vector<EvalWell>& cq_s,
467 const IntensiveQuantities& intQuants,
469 };
470
471}
472
473#include "StandardWell_impl.hpp"
474
475#endif // OPM_STANDARDWELL_HEADER_INCLUDED
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Contains the classes required to extend the black-oil model by brine.
Contains the classes required to extend the black-oil model by solvent component.
Contains the classes required to extend the black-oil model to include the effects of foam.
Contains the classes required to extend the black-oil model by MICP.
Contains the classes required to extend the black-oil model by polymer.
Contains the classes required to extend the black-oil model by solvents.
Contains the high level supplements required to extend the black oil model by brine.
Definition blackoilbrinemodules.hh:50
Contains the high level supplements required to extend the black oil model to include the effects of ...
Definition blackoilfoammodules.hh:60
Contains the high level supplements required to extend the black oil model by polymer.
Definition blackoilpolymermodules.hh:54
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:43
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:195
Manages the initializing and running of time dependent problems.
Definition simulator.hh:97
Definition StandardWellEval.hpp:46
Definition StandardWell.hpp:60
void computeWellPotentials(const Simulator &simulator, const WellState< Scalar > &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition StandardWell_impl.hpp:1746
std::vector< Scalar > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const override
Compute well rates based on current reservoir conditions and well variables.
Definition StandardWell_impl.hpp:2505
void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellState< Scalar > &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:1456
bool jacobianContainsWellContributions() const override
Wether the Jacobian will also have well contributions in it.
Definition StandardWell.hpp:214
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition StandardWell_impl.hpp:1424
virtual ConvergenceReport getWellConvergence(const Simulator &simulator, const WellState< Scalar > &well_state, const std::vector< Scalar > &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:1200
const std::string & name() const
Well name.
Definition WellInterfaceGeneric.cpp:164
Definition WellInterface.hpp:77
Collect per-connection static information to enable calculating connection-level or well-level produc...
Definition WellProdIndexCalculator.hpp:37
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:66
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
bool matrix_add_well_contributions_
Whether to add influences of wells between cells to the matrix and preconditioner matrix.
Definition BlackoilModelParameters.hpp:280
Static data associated with a well perforation.
Definition PerforationData.hpp:30
Definition PerforationData.hpp:41
Definition BlackoilPhases.hpp:46