My Project
Loading...
Searching...
No Matches
WellInterface.hpp
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4 Copyright 2017 IRIS
5 Copyright 2019 Norce
6
7 This file is part of the Open Porous Media project (OPM).
8
9 OPM is free software: you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation, either version 3 of the License, or
12 (at your option) any later version.
13
14 OPM is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18
19 You should have received a copy of the GNU General Public License
20 along with OPM. If not, see <http://www.gnu.org/licenses/>.
21*/
22
23#ifndef OPM_WELLINTERFACE_HEADER_INCLUDED
24#define OPM_WELLINTERFACE_HEADER_INCLUDED
25
26// NOTE: GasLiftSingleWell.hpp includes StandardWell.hpp which includes ourself
27// (WellInterface.hpp), so we need to forward declare GasLiftSingleWell
28// for it to be defined in this file. Similar for BlackoilWellModel
29namespace Opm {
30 template<typename TypeTag> class GasLiftSingleWell;
31 template<typename TypeTag> class BlackoilWellModel;
32}
33
34#include <opm/common/OpmLog/OpmLog.hpp>
35#include <opm/common/ErrorMacros.hpp>
36#include <opm/common/Exceptions.hpp>
37
38#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
39
40#include <opm/material/fluidstates/BlackOilFluidState.hpp>
41
43
45
46#include <opm/simulators/wells/BlackoilWellModel.hpp>
47#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
48#include <opm/simulators/wells/GasLiftSingleWell.hpp>
49#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
50#include <opm/simulators/wells/PerforationData.hpp>
51#include <opm/simulators/wells/WellInterfaceIndices.hpp>
52#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
53#include <opm/simulators/wells/WellState.hpp>
54
55#include <opm/simulators/timestepping/ConvergenceReport.hpp>
56
57#include <opm/simulators/utils/BlackoilPhases.hpp>
58#include <opm/simulators/utils/DeferredLogger.hpp>
59
60#include <dune/common/fmatrix.hh>
61#include <dune/istl/bcrsmatrix.hh>
62#include <dune/istl/matrixmatrix.hh>
63
64#include <opm/material/densead/Evaluation.hpp>
65
66#include <vector>
67
68namespace Opm
69{
70
73
74template<typename TypeTag>
75class WellInterface : public WellInterfaceIndices<GetPropType<TypeTag, Properties::FluidSystem>,
76 GetPropType<TypeTag, Properties::Indices>>
77{
80public:
91 using GLiftEclWells = typename GasLiftGroupInfo<Scalar>::GLiftEclWells;
92
93 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
94 using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
95 using Eval = typename Base::Eval;
96 using BVector = Dune::BlockVector<VectorBlockType>;
97 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<Scalar, 1, 1>>;
98
99 using RateConverterType =
101
102 using WellInterfaceFluidSystem<FluidSystem>::Gas;
103 using WellInterfaceFluidSystem<FluidSystem>::Oil;
104 using WellInterfaceFluidSystem<FluidSystem>::Water;
105
106 using ModelParameters = typename Base::ModelParameters;
107
108 static constexpr bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
109 static constexpr bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
110 static constexpr bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
111 static constexpr bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
112 static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
113 // flag for polymer molecular weight related
114 static constexpr bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
115 static constexpr bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
116 static constexpr bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
117 static constexpr bool has_watVapor = getPropValue<TypeTag, Properties::EnableVapwat>();
118 static constexpr bool has_disgas_in_water = getPropValue<TypeTag, Properties::EnableDisgasInWater>();
119 static constexpr bool has_saltPrecip = getPropValue<TypeTag, Properties::EnableSaltPrecipitation>();
120 static constexpr bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
121
122 // For the conversion between the surface volume rate and reservoir voidage rate
123 using FluidState = BlackOilFluidState<Eval,
124 FluidSystem,
125 has_temperature,
126 has_energy,
127 Indices::compositionSwitchIdx >= 0,
128 has_watVapor,
129 has_brine,
130 has_saltPrecip,
131 has_disgas_in_water,
132 Indices::numPhases >;
134 WellInterface(const Well& well,
136 const int time_step,
137 const ModelParameters& param,
138 const RateConverterType& rate_converter,
139 const int pvtRegionIdx,
140 const int num_components,
141 const int num_phases,
142 const int index_of_well,
143 const std::vector<PerforationData<Scalar>>& perf_data);
144
146 virtual ~WellInterface() = default;
147
148 virtual void init(const PhaseUsage* phase_usage_arg,
149 const std::vector<Scalar>& depth_arg,
150 const Scalar gravity_arg,
151 const std::vector<Scalar>& B_avg,
152 const bool changed_to_open_this_step);
153
154 virtual ConvergenceReport getWellConvergence(const Simulator& simulator,
155 const WellState<Scalar>& well_state,
156 const std::vector<Scalar>& B_avg,
158 const bool relax_tolerance) const = 0;
159
160 virtual void solveEqAndUpdateWellState(const Simulator& simulator,
161 WellState<Scalar>& well_state,
163
164 void assembleWellEq(const Simulator& simulator,
165 const double dt,
166 WellState<Scalar>& well_state,
167 const GroupState<Scalar>& group_state,
169
170 void assembleWellEqWithoutIteration(const Simulator& simulator,
171 const double dt,
172 WellState<Scalar>& well_state,
173 const GroupState<Scalar>& group_state,
175
176 // TODO: better name or further refactoring the function to make it more clear
177 void prepareWellBeforeAssembling(const Simulator& simulator,
178 const double dt,
179 WellState<Scalar>& well_state,
180 const GroupState<Scalar>& group_state,
182
183
184 virtual void computeWellRatesWithBhp(const Simulator& ebosSimulator,
185 const Scalar& bhp,
186 std::vector<Scalar>& well_flux,
188
189 virtual std::optional<Scalar>
190 computeBhpAtThpLimitProdWithAlq(const Simulator& ebos_simulator,
192 const Scalar alq_value,
194 bool iterate_if_no_solution) const = 0;
195
198 virtual void recoverWellSolutionAndUpdateWellState(const Simulator& simulator,
199 const BVector& x,
200 WellState<Scalar>& well_state,
202
204 virtual void apply(const BVector& x, BVector& Ax) const = 0;
205
207 virtual void apply(BVector& r) const = 0;
208
209 // TODO: before we decide to put more information under mutable, this function is not const
210 virtual void computeWellPotentials(const Simulator& simulator,
211 const WellState<Scalar>& well_state,
212 std::vector<Scalar>& well_potentials,
214
215 virtual void updateWellStateWithTarget(const Simulator& simulator,
216 const GroupState<Scalar>& group_state,
217 WellState<Scalar>& well_state,
219
220 virtual void computeWellRatesWithBhpIterations(const Simulator& simulator,
221 const Scalar& bhp,
222 std::vector<Scalar>& well_flux,
224
225 bool wellUnderZeroRateTarget(const Simulator& simulator,
226 const WellState<Scalar>& well_state,
228
229 bool wellUnderZeroGroupRateTarget(const Simulator& simulator,
230 const WellState<Scalar>& well_state,
232 std::optional<bool> group_control = std::nullopt) const;
233
234 bool stoppedOrZeroRateTarget(const Simulator& simulator,
235 const WellState<Scalar>& well_state,
237
238 bool updateWellStateWithTHPTargetProd(const Simulator& simulator,
239 WellState<Scalar>& well_state,
241
242 enum class IndividualOrGroup { Individual, Group, Both };
243 bool updateWellControl(const Simulator& simulator,
244 const IndividualOrGroup iog,
245 WellState<Scalar>& well_state,
246 const GroupState<Scalar>& group_state,
247 DeferredLogger& deferred_logger) /* const */;
248
249 bool updateWellControlAndStatusLocalIteration(const Simulator& simulator,
250 WellState<Scalar>& well_state,
251 const GroupState<Scalar>& group_state,
252 const Well::InjectionControls& inj_controls,
253 const Well::ProductionControls& prod_controls,
254 const Scalar WQTotal,
256 const bool fixed_control = false,
257 const bool fixed_status = false);
258
259 virtual void updatePrimaryVariables(const Simulator& simulator,
260 const WellState<Scalar>& well_state,
262
263 virtual void calculateExplicitQuantities(const Simulator& simulator,
264 const WellState<Scalar>& well_state,
265 DeferredLogger& deferred_logger) = 0; // should be const?
266
267 virtual void updateProductivityIndex(const Simulator& simulator,
269 WellState<Scalar>& well_state,
271
274 {
275 return false;
276 }
277
278 // Add well contributions to matrix
279 virtual void addWellContributions(SparseMatrixAdapter&) const = 0;
280
281 virtual void addWellPressureEquations(PressureMatrix& mat,
282 const BVector& x,
283 const int pressureVarIndex,
284 const bool use_well_weights,
285 const WellState<Scalar>& well_state) const = 0;
286
287 void addCellRates(RateVector& rates, int cellIdx) const;
288
289 Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const;
290
291 // TODO: theoretically, it should be a const function
292 // Simulator is not const is because that assembleWellEq is non-const Simulator
293 void wellTesting(const Simulator& simulator,
294 const double simulation_time,
295 /* const */ WellState<Scalar>& well_state,
296 const GroupState<Scalar>& group_state,
298 const PhaseUsage& phase_usage,
299 GLiftEclWells& ecl_well_map,
300 std::map<std::string, double>& open_times,
302
303 void checkWellOperability(const Simulator& simulator,
304 const WellState<Scalar>& well_state,
306
307 bool gliftBeginTimeStepWellTestIterateWellEquations(const Simulator& ebos_simulator,
308 const double dt,
309 WellState<Scalar>& well_state,
310 const GroupState<Scalar>& group_state,
312
313 void gliftBeginTimeStepWellTestUpdateALQ(const Simulator& simulator,
314 WellState<Scalar>& well_state,
315 const GroupState<Scalar>& group_state,
316 const PhaseUsage& phase_usage,
317 GLiftEclWells& ecl_well_map,
319
320 // check whether the well is operable under the current reservoir condition
321 // mostly related to BHP limit and THP limit
322 void updateWellOperability(const Simulator& simulator,
323 const WellState<Scalar>& well_state,
325
326 bool updateWellOperabilityFromWellEq(const Simulator& simulator,
327 const WellState<Scalar>& well_state,
329
330 // update perforation water throughput based on solved water rate
331 virtual void updateWaterThroughput(const double dt,
332 WellState<Scalar>& well_state) const = 0;
333
336 virtual std::vector<Scalar>
337 computeCurrentWellRates(const Simulator& simulator,
339
343 void updateWellStateRates(const Simulator& simulator,
344 WellState<Scalar>& well_state,
346
347 void solveWellEquation(const Simulator& simulator,
348 WellState<Scalar>& well_state,
349 const GroupState<Scalar>& group_state,
351
352 const std::vector<RateVector>& connectionRates() const
353 {
354 return connectionRates_;
355 }
356
357 std::vector<Scalar> wellIndex(const int perf,
358 const IntensiveQuantities& intQuants,
359 const Scalar trans_mult,
360 const SingleWellState<Scalar>& ws) const;
361
362 void updateConnectionDFactor(const Simulator& simulator,
364
365 void updateConnectionTransmissibilityFactor(const Simulator& simulator,
367
368 virtual bool iterateWellEqWithSwitching(const Simulator& simulator,
369 const double dt,
372 WellState<Scalar>& well_state,
373 const GroupState<Scalar>& group_state,
375 const bool fixed_control = false,
376 const bool fixed_status = false) = 0;
377protected:
378 // simulation parameters
379 std::vector<RateVector> connectionRates_;
380 std::vector<Scalar> B_avg_;
381 bool changed_to_stopped_this_step_ = false;
382 bool thp_update_iterations = false;
383
384 Scalar wpolymer() const;
385 Scalar wfoam() const;
386 Scalar wsalt() const;
387 Scalar wmicrobes() const;
388 Scalar woxygen() const;
389 Scalar wurea() const;
390
391 virtual Scalar getRefDensity() const = 0;
392
393 // Component fractions for each phase for the well
394 const std::vector<Scalar>& compFrac() const;
395
396 std::vector<Scalar>
397 initialWellRateFractions(const Simulator& ebosSimulator,
398 const WellState<Scalar>& well_state) const;
399
400 // check whether the well is operable under BHP limit with current reservoir condition
401 virtual void checkOperabilityUnderBHPLimit(const WellState<Scalar>& well_state,
402 const Simulator& simulator,
404
405 // check whether the well is operable under THP limit with current reservoir condition
406 virtual void checkOperabilityUnderTHPLimit(const Simulator& simulator,
407 const WellState<Scalar>& well_state,
409
410 virtual void updateIPR(const Simulator& simulator,
412
413 virtual void assembleWellEqWithoutIteration(const Simulator& simulator,
414 const double dt,
417 WellState<Scalar>& well_state,
418 const GroupState<Scalar>& group_state,
420
421 // iterate well equations with the specified control until converged
422 virtual bool iterateWellEqWithControl(const Simulator& simulator,
423 const double dt,
426 WellState<Scalar>& well_state,
427 const GroupState<Scalar>& group_state,
429
430 virtual void updateIPRImplicit(const Simulator& simulator,
431 WellState<Scalar>& well_state,
433
434 bool iterateWellEquations(const Simulator& simulator,
435 const double dt,
436 WellState<Scalar>& well_state,
437 const GroupState<Scalar>& group_state,
439
440 bool solveWellWithTHPConstraint(const Simulator& simulator,
441 const double dt,
442 const Well::InjectionControls& inj_controls,
443 const Well::ProductionControls& prod_controls,
444 WellState<Scalar>& well_state,
445 const GroupState<Scalar>& group_state,
447
448 std::optional<Scalar>
449 estimateOperableBhp(const Simulator& ebos_simulator,
450 const double dt,
451 WellState<Scalar>& well_state,
454
455 bool solveWellWithBhp(const Simulator& simulator,
456 const double dt,
457 const Scalar bhp,
458 WellState<Scalar>& well_state,
460
461 bool solveWellWithZeroRate(const Simulator& simulator,
462 const double dt,
463 WellState<Scalar>& well_state,
465
466 bool solveWellForTesting(const Simulator& simulator,
467 WellState<Scalar>& well_state,
468 const GroupState<Scalar>& group_state,
470
471
472 template<class GasLiftSingleWell>
473 std::unique_ptr<GasLiftSingleWell> initializeGliftWellTest_(const Simulator& simulator,
474 WellState<Scalar>& well_state,
475 const GroupState<Scalar>& group_state,
476 const PhaseUsage& phase_usage,
477 GLiftEclWells& ecl_well_map,
479
480 Eval getPerfCellPressure(const FluidState& fs) const;
481
482 // get the mobility for specific perforation
483 template<class Value, class Callback>
484 void getMobility(const Simulator& simulator,
485 const int perf,
486 std::vector<Value>& mob,
487 Callback& extendEval,
489
490 void computeConnLevelProdInd(const FluidState& fs,
491 const std::function<Scalar(const Scalar)>& connPICalc,
492 const std::vector<Scalar>& mobility,
493 Scalar* connPI) const;
494
495 void computeConnLevelInjInd(const FluidState& fs,
496 const Phase preferred_phase,
497 const std::function<Scalar(const Scalar)>& connIICalc,
498 const std::vector<Scalar>& mobility,
499 Scalar* connII,
501
502 Scalar computeConnectionDFactor(const int perf,
503 const IntensiveQuantities& intQuants,
504 const SingleWellState<Scalar>& ws) const;
505};
506
507} // namespace Opm
508
509#include "WellInterface_impl.hpp"
510
511#endif // OPM_WELLINTERFACE_HEADER_INCLUDED
Declares the properties required by the black oil model.
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 GasLiftSingleWell.hpp:37
Definition GroupState.hpp:43
Class encapsulating some information about parallel wells.
Definition ParallelWellInfo.hpp:195
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition RateConverter.hpp:71
Manages the initializing and running of time dependent problems.
Definition simulator.hh:97
Definition SingleWellState.hpp:42
Definition WellInterfaceFluidSystem.hpp:51
Definition WellInterfaceIndices.hpp:34
Definition WellInterface.hpp:77
virtual std::vector< Scalar > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const =0
Compute well rates based on current reservoir conditions and well variables.
virtual void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)=0
using the solution x to recover the solution xw for wells and applying xw to update Well State
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
void updateWellStateRates(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Modify the well_state's rates if there is only one nonzero rate.
Definition WellInterface_impl.hpp:1619
virtual void apply(BVector &r) const =0
r = r - C D^-1 Rw
virtual ~WellInterface()=default
Virtual destructor.
virtual bool jacobianContainsWellContributions() const
Wether the Jacobian will also have well contributions in it.
Definition WellInterface.hpp:273
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
Declares the properties required by the black oil model.
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
Static data associated with a well perforation.
Definition PerforationData.hpp:30
Definition BlackoilPhases.hpp:46