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
55namespace 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 using typename Base::PressureMatrix;
98
99 // number of the conservation equations
100 static constexpr int numWellConservationEq = Indices::numPhases + Indices::numSolvents;
101 // number of the well control equations
102 static constexpr int numWellControlEq = 1;
103 // number of the well equations that will always be used
104 // based on the solution strategy, there might be other well equations be introduced
105 static constexpr int numStaticWellEq = numWellConservationEq + numWellControlEq;
106
107 // the index for Bhp in primary variables and also the index of well control equation
108 // they both will be the last one in their respective system.
109 // TODO: we should have indices for the well equations and well primary variables separately
110 static constexpr int Bhp = numStaticWellEq - numWellControlEq;
111
112 using typename Base::Scalar;
113
114
115 using Base::name;
116 using Base::Water;
117 using Base::Oil;
118 using Base::Gas;
119
120 using typename Base::BVector;
121
122 using Eval = typename StdWellEval::Eval;
123 using EvalWell = typename StdWellEval::EvalWell;
124 using BVectorWell = typename StdWellEval::BVectorWell;
125 using DiagMatrixBlockWellType = typename StdWellEval::DiagMatrixBlockWellType;
126
127 StandardWell(const Well& well,
128 const ParallelWellInfo& pw_info,
129 const int time_step,
130 const ModelParameters& param,
131 const RateConverterType& rate_converter,
132 const int pvtRegionIdx,
133 const int num_components,
134 const int num_phases,
135 const int index_of_well,
136 const std::vector<PerforationData>& perf_data);
137
138 virtual void init(const PhaseUsage* phase_usage_arg,
139 const std::vector<double>& depth_arg,
140 const double gravity_arg,
141 const int num_cells,
142 const std::vector< Scalar >& B_avg,
143 const bool changed_to_open_this_step) override;
144
145
146 virtual void initPrimaryVariablesEvaluation() const override;
147
149 virtual ConvergenceReport getWellConvergence(const WellState& well_state,
150 const std::vector<double>& B_avg,
151 DeferredLogger& deferred_logger,
152 const bool relax_tolerance = false) const override;
153
155 virtual void apply(const BVector& x, BVector& Ax) const override;
157 virtual void apply(BVector& r) const override;
158
161 virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
162 WellState& well_state,
163 DeferredLogger& deferred_logger) const override;
164
166 virtual void computeWellPotentials(const Simulator& ebosSimulator,
167 const WellState& well_state,
168 std::vector<double>& well_potentials,
169 DeferredLogger& deferred_logger) /* const */ override;
170
171 virtual void updatePrimaryVariables(const WellState& well_state, DeferredLogger& deferred_logger) const override;
172
173 virtual void solveEqAndUpdateWellState(WellState& well_state, DeferredLogger& deferred_logger) override;
174
175 virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
176 const WellState& well_state,
177 DeferredLogger& deferred_logger) override; // should be const?
178
179 virtual void updateProductivityIndex(const Simulator& ebosSimulator,
180 const WellProdIndexCalculator& wellPICalc,
181 WellState& well_state,
182 DeferredLogger& deferred_logger) const override;
183
184 virtual void addWellContributions(SparseMatrixAdapter& mat) const override;
185
186 virtual void addWellPressureEquations(PressureMatrix& mat,
187 const BVector& x,
188 const int pressureVarIndex,
189 const bool use_well_weights,
190 const WellState& well_state) const override;
191
192 // iterate well equations with the specified control until converged
193 bool iterateWellEqWithControl(const Simulator& ebosSimulator,
194 const double dt,
195 const Well::InjectionControls& inj_controls,
196 const Well::ProductionControls& prod_controls,
197 WellState& well_state,
198 const GroupState& group_state,
199 DeferredLogger& deferred_logger) override;
200
202 virtual bool jacobianContainsWellContributions() const override
203 {
204 return this->param_.matrix_add_well_contributions_;
205 }
206
207 /* returns BHP */
208 double computeWellRatesAndBhpWithThpAlqProd(const Simulator &ebos_simulator,
209 const SummaryState &summary_state,
210 DeferredLogger &deferred_logger,
211 std::vector<double> &potentials,
212 double alq) const;
213
214 void computeWellRatesWithThpAlqProd(
215 const Simulator &ebos_simulator,
216 const SummaryState &summary_state,
217 DeferredLogger &deferred_logger,
218 std::vector<double> &potentials,
219 double alq) const;
220
221 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
222 const Simulator& ebos_simulator,
223 const SummaryState& summary_state,
224 DeferredLogger& deferred_logger,
225 double alq_value) const override;
226
227 virtual void computeWellRatesWithBhp(
228 const Simulator& ebosSimulator,
229 const double& bhp,
230 std::vector<double>& well_flux,
231 DeferredLogger& deferred_logger) const override;
232
233 // NOTE: These cannot be protected since they are used by GasLiftRuntime
234 using Base::phaseUsage;
235 using Base::vfp_properties_;
236
237 virtual std::vector<double> computeCurrentWellRates(const Simulator& ebosSimulator,
238 DeferredLogger& deferred_logger) const override;
239
240 void computeConnLevelProdInd(const FluidState& fs,
241 const std::function<double(const double)>& connPICalc,
242 const std::vector<EvalWell>& mobility,
243 double* connPI) const;
244
245 void computeConnLevelInjInd(const typename StandardWell<TypeTag>::FluidState& fs,
246 const Phase preferred_phase,
247 const std::function<double(const double)>& connIICalc,
248 const std::vector<EvalWell>& mobility,
249 double* connII,
250 DeferredLogger& deferred_logger) const;
251
252
253 protected:
254 bool regularize_;
255
256 // xw = inv(D)*(rw - C*x)
257 void recoverSolutionWell(const BVector& x, BVectorWell& xw) const;
258
259 // updating the well_state based on well solution dwells
260 void updateWellState(const BVectorWell& dwells,
261 WellState& well_state,
262 DeferredLogger& deferred_logger) const;
263
264 // calculate the properties for the well connections
265 // to calulate the pressure difference between well connections.
266 void computePropertiesForWellConnectionPressures(const Simulator& ebosSimulator,
267 const WellState& well_state,
268 std::vector<double>& b_perf,
269 std::vector<double>& rsmax_perf,
270 std::vector<double>& rvmax_perf,
271 std::vector<double>& rvwmax_perf,
272 std::vector<double>& surf_dens_perf) const;
273
274 void computeWellConnectionDensitesPressures(const Simulator& ebosSimulator,
275 const WellState& well_state,
276 const std::vector<double>& b_perf,
277 const std::vector<double>& rsmax_perf,
278 const std::vector<double>& rvmax_perf,
279 const std::vector<double>& rvwmax_perf,
280 const std::vector<double>& surf_dens_perf,
281 DeferredLogger& deferred_logger);
282
283 void computeWellConnectionPressures(const Simulator& ebosSimulator,
284 const WellState& well_state,
285 DeferredLogger& deferred_logger);
286
287 void computePerfRateEval(const IntensiveQuantities& intQuants,
288 const std::vector<EvalWell>& mob,
289 const EvalWell& bhp,
290 const double Tw,
291 const int perf,
292 const bool allow_cf,
293 std::vector<EvalWell>& cq_s,
294 double& perf_dis_gas_rate,
295 double& perf_vap_oil_rate,
296 double& perf_vap_wat_rate,
297 DeferredLogger& deferred_logger) const;
298
299 void computePerfRateScalar(const IntensiveQuantities& intQuants,
300 const std::vector<Scalar>& mob,
301 const Scalar& bhp,
302 const double Tw,
303 const int perf,
304 const bool allow_cf,
305 std::vector<Scalar>& cq_s,
306 DeferredLogger& deferred_logger) const;
307
308 template<class Value>
309 void computePerfRate(const std::vector<Value>& mob,
310 const Value& pressure,
311 const Value& bhp,
312 const Value& rs,
313 const Value& rv,
314 const Value& rvw,
315 std::vector<Value>& b_perfcells_dense,
316 const double Tw,
317 const int perf,
318 const bool allow_cf,
319 const Value& skin_pressure,
320 const std::vector<Value>& cmix_s,
321 std::vector<Value>& cq_s,
322 double& perf_dis_gas_rate,
323 double& perf_vap_oil_rate,
324 double& perf_vap_wat_rate,
325 DeferredLogger& deferred_logger) const;
326
327 void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
328 const double& bhp,
329 std::vector<double>& well_flux,
330 DeferredLogger& deferred_logger) const;
331
332 std::vector<double> computeWellPotentialWithTHP(
333 const Simulator& ebosSimulator,
334 DeferredLogger& deferred_logger,
335 const WellState &well_state) const;
336
337
338 virtual double getRefDensity() const override;
339
340 // get the mobility for specific perforation
341 void getMobilityEval(const Simulator& ebosSimulator,
342 const int perf,
343 std::vector<EvalWell>& mob,
344 DeferredLogger& deferred_logger) const;
345
346 // get the mobility for specific perforation
347 void getMobilityScalar(const Simulator& ebosSimulator,
348 const int perf,
349 std::vector<Scalar>& mob,
350 DeferredLogger& deferred_logger) const;
351
352
353 void updateWaterMobilityWithPolymer(const Simulator& ebos_simulator,
354 const int perf,
355 std::vector<EvalWell>& mob_water,
356 DeferredLogger& deferred_logger) const;
357
358 void updatePrimaryVariablesNewton(const BVectorWell& dwells,
359 const WellState& well_state,
360 DeferredLogger& deferred_logger) const;
361
362 // update extra primary vriables if there are any
363 void updateExtraPrimaryVariables(const BVectorWell& dwells) const;
364
365
366 void updateWellStateFromPrimaryVariables(WellState& well_state, DeferredLogger& deferred_logger) const;
367
368 virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator,
369 const double dt,
370 const Well::InjectionControls& inj_controls,
371 const Well::ProductionControls& prod_controls,
372 WellState& well_state,
373 const GroupState& group_state,
374 DeferredLogger& deferred_logger) override;
375
376 void assembleWellEqWithoutIterationImpl(const Simulator& ebosSimulator,
377 const double dt,
378 WellState& well_state,
379 const GroupState& group_state,
380 DeferredLogger& deferred_logger);
381
382 void calculateSinglePerf(const Simulator& ebosSimulator,
383 const int perf,
384 WellState& well_state,
385 std::vector<RateVector>& connectionRates,
386 std::vector<EvalWell>& cq_s,
387 EvalWell& water_flux_s,
388 EvalWell& cq_s_zfrac_effective,
389 DeferredLogger& deferred_logger) const;
390
391 // check whether the well is operable under BHP limit with current reservoir condition
392 virtual void checkOperabilityUnderBHPLimit(const WellState& well_state, const Simulator& ebos_simulator, DeferredLogger& deferred_logger) override;
393
394 // check whether the well is operable under THP limit with current reservoir condition
395 virtual void checkOperabilityUnderTHPLimit(const Simulator& ebos_simulator, const WellState& well_state, DeferredLogger& deferred_logger) override;
396
397 // updating the inflow based on the current reservoir condition
398 virtual void updateIPR(const Simulator& ebos_simulator, DeferredLogger& deferred_logger) const override;
399
400 // for a well, when all drawdown are in the wrong direction, then this well will not
401 // be able to produce/inject .
402 bool allDrawDownWrongDirection(const Simulator& ebos_simulator) const;
403
404 // whether the well can produce / inject based on the current well state (bhp)
405 bool canProduceInjectWithCurrentBhp(const Simulator& ebos_simulator,
406 const WellState& well_state,
407 DeferredLogger& deferred_logger);
408
409 // turn on crossflow to avoid singular well equations
410 // when the well is banned from cross-flow and the BHP is not properly initialized,
411 // we turn on crossflow to avoid singular well equations. It can result in wrong-signed
412 // well rates, it can cause problem for THP calculation
413 // TODO: looking for better alternative to avoid wrong-signed well rates
414 bool openCrossFlowAvoidSingularity(const Simulator& ebos_simulator) const;
415
416 // calculate the skin pressure based on water velocity, throughput and polymer concentration.
417 // throughput is used to describe the formation damage during water/polymer injection.
418 // calculated skin pressure will be applied to the drawdown during perforation rate calculation
419 // to handle the effect from formation damage.
420 EvalWell pskin(const double throuhgput,
421 const EvalWell& water_velocity,
422 const EvalWell& poly_inj_conc,
423 DeferredLogger& deferred_logger) const;
424
425 // calculate the skin pressure based on water velocity, throughput during water injection.
426 EvalWell pskinwater(const double throughput,
427 const EvalWell& water_velocity,
428 DeferredLogger& deferred_logger) const;
429
430 // calculate the injecting polymer molecular weight based on the througput and water velocity
431 EvalWell wpolymermw(const double throughput,
432 const EvalWell& water_velocity,
433 DeferredLogger& deferred_logger) const;
434
435 // modify the water rate for polymer injectivity study
436 void handleInjectivityRate(const Simulator& ebosSimulator,
437 const int perf,
438 std::vector<EvalWell>& cq_s) const;
439
440 // handle the extra equations for polymer injectivity study
441 void handleInjectivityEquations(const Simulator& ebosSimulator,
442 const WellState& well_state,
443 const int perf,
444 const EvalWell& water_flux_s,
445 DeferredLogger& deferred_logger);
446
447 virtual void updateWaterThroughput(const double dt, WellState& well_state) const override;
448
449 // checking convergence of extra equations, if there are any
450 void checkConvergenceExtraEqs(const std::vector<double>& res,
451 ConvergenceReport& report) const;
452
453 // updating the connectionRates_ related polymer molecular weight
454 void updateConnectionRatePolyMW(const EvalWell& cq_s_poly,
455 const IntensiveQuantities& int_quants,
456 const WellState& well_state,
457 const int perf,
458 std::vector<RateVector>& connectionRates,
459 DeferredLogger& deferred_logger) const;
460
461
462 std::optional<double> computeBhpAtThpLimitProd(const WellState& well_state,
463 const Simulator& ebos_simulator,
464 const SummaryState& summary_state,
465 DeferredLogger& deferred_logger) const;
466
467 std::optional<double> computeBhpAtThpLimitInj(const Simulator& ebos_simulator,
468 const SummaryState& summary_state,
469 DeferredLogger& deferred_logger) const;
470
471 };
472
473}
474
475#include "StandardWell_impl.hpp"
476
477#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:243
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:202
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:2709
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:1465
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:1777
virtual void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: StandardWell_impl.hpp:1713
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:1963
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