24#ifndef OPM_WELLINTERFACE_HEADER_INCLUDED
25#define OPM_WELLINTERFACE_HEADER_INCLUDED
27#include <opm/common/OpmLog/OpmLog.hpp>
28#include <opm/common/ErrorMacros.hpp>
29#include <opm/common/Exceptions.hpp>
31#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
33#include <opm/core/props/BlackoilPhases.hpp>
35#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
36#include <opm/simulators/wells/WellState.hpp>
41 template<
typename TypeTag>
class GasLiftSingleWell;
42 template<
typename TypeTag>
class BlackoilWellModel;
44#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
45#include <opm/simulators/wells/GasLiftSingleWell.hpp>
46#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
47#include <opm/simulators/wells/BlackoilWellModel.hpp>
48#include <opm/simulators/flow/BlackoilModelParametersEbos.hpp>
50#include <opm/simulators/utils/DeferredLogger.hpp>
52#include<dune/common/fmatrix.hh>
53#include<dune/istl/bcrsmatrix.hh>
54#include<dune/istl/matrixmatrix.hh>
56#include <opm/material/densead/Evaluation.hpp>
58#include <opm/simulators/wells/WellInterfaceIndices.hpp>
59#include <opm/simulators/timestepping/ConvergenceReport.hpp>
67class WellInjectionProperties;
68class WellProductionProperties;
70template<
typename TypeTag>
72 GetPropType<TypeTag, Properties::Indices>,
73 GetPropType<TypeTag, Properties::Scalar>>
78 using Grid = GetPropType<TypeTag, Properties::Grid>;
79 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
80 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
81 using Indices = GetPropType<TypeTag, Properties::Indices>;
82 using IntensiveQuantities = GetPropType<TypeTag, Properties::IntensiveQuantities>;
83 using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
84 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
85 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
87 using GLiftOptWells =
typename BlackoilWellModel<TypeTag>::GLiftOptWells;
88 using GLiftProdWells =
typename BlackoilWellModel<TypeTag>::GLiftProdWells;
89 using GLiftWellStateMap =
90 typename BlackoilWellModel<TypeTag>::GLiftWellStateMap;
91 using GLiftSyncGroups =
typename GasLiftSingleWellGeneric::GLiftSyncGroups;
93 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
95 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
96 using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
97 using BVector = Dune::BlockVector<VectorBlockType>;
98 using Eval = DenseAd::Evaluation<Scalar, Indices::numEq>;
99 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
101 using RateConverterType =
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>();
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::EnableEvaporation>();
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>();
123 using FluidState = BlackOilFluidState<Eval,
127 Indices::compositionSwitchIdx >= 0,
132 Indices::numPhases >;
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>& perf_data);
148 virtual void init(
const PhaseUsage* phase_usage_arg,
149 const std::vector<double>& depth_arg,
150 const double gravity_arg,
152 const std::vector< Scalar >& B_avg,
153 const bool changed_to_open_this_step);
155 virtual void initPrimaryVariablesEvaluation() = 0;
159 void assembleWellEq(
const Simulator& ebosSimulator,
165 virtual void computeWellRatesWithBhp(
166 const Simulator& ebosSimulator,
168 std::vector<double>& well_flux,
172 virtual std::optional<double> computeBhpAtThpLimitProdWithAlq(
173 const Simulator& ebos_simulator,
174 const SummaryState& summary_state,
175 const double alq_value,
187 virtual void apply(
const BVector& x, BVector& Ax)
const = 0;
190 virtual void apply(BVector& r)
const = 0;
193 virtual void computeWellPotentials(
const Simulator& ebosSimulator,
195 std::vector<double>& well_potentials,
198 virtual void updateWellStateWithTarget(
const Simulator& ebos_simulator,
203 virtual bool updateWellStateWithTHPTargetProd(
const Simulator& ebos_simulator,
207 enum class IndividualOrGroup { Individual, Group, Both };
208 bool updateWellControl(
const Simulator& ebos_simulator,
209 const IndividualOrGroup iog,
214 virtual void updatePrimaryVariables(
const SummaryState& summary_state,
218 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
222 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
234 virtual void addWellContributions(SparseMatrixAdapter&)
const = 0;
236 virtual void addWellPressureEquations(PressureMatrix& mat,
238 const int pressureVarIndex,
239 const bool use_well_weights,
242 void addCellRates(RateVector& rates,
int cellIdx)
const;
244 Scalar volumetricSurfaceRateForConnection(
int cellIdx,
int phaseIdx)
const;
246 template <
class EvalWell>
247 Eval restrictEval(
const EvalWell& in)
const
250 out.setValue(in.value());
251 for (
int eqIdx = 0; eqIdx < Indices::numEq; ++eqIdx) {
252 out.setDerivative(eqIdx, in.derivative(eqIdx));
259 void wellTesting(
const Simulator& simulator,
260 const double simulation_time,
261 WellState& well_state,
const GroupState& group_state, WellTestState& welltest_state,
262 DeferredLogger& deferred_logger);
264 void checkWellOperability(
const Simulator& ebos_simulator,
const WellState& well_state, DeferredLogger& deferred_logger);
266 void gliftBeginTimeStepWellTestUpdateALQ(
const Simulator& ebos_simulator,
267 WellState& well_state,
268 DeferredLogger& deferred_logger);
272 void updateWellOperability(
const Simulator& ebos_simulator,
273 const WellState& well_state,
274 DeferredLogger& deferred_logger);
278 virtual void updateWaterThroughput(
const double dt, WellState& well_state)
const = 0;
292 void solveWellEquation(
const Simulator& ebosSimulator,
297 const std::vector<RateVector>& connectionRates()
const
299 return connectionRates_;
305 const ModelParameters& param_;
307 std::vector<RateVector> connectionRates_;
309 std::vector< Scalar > B_avg_;
311 bool changed_to_stopped_this_step_ =
false;
313 double wpolymer()
const;
315 double wfoam()
const;
317 double wsalt()
const;
319 double wmicrobes()
const;
321 double woxygen()
const;
323 double wurea()
const;
325 virtual double getRefDensity()
const = 0;
328 const std::vector<double>& compFrac()
const;
330 std::vector<double> initialWellRateFractions(
const Simulator& ebosSimulator,
const WellState& well_state)
const;
333 virtual void checkOperabilityUnderBHPLimit(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger) =0;
336 virtual void checkOperabilityUnderTHPLimit(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger) =0;
338 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const=0;
340 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
342 const WellInjectionControls& inj_controls,
343 const WellProductionControls& prod_controls,
349 virtual bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
351 const WellInjectionControls& inj_controls,
352 const WellProductionControls& prod_controls,
357 bool iterateWellEquations(
const Simulator& ebosSimulator,
363 bool solveWellForTesting(
const Simulator& ebosSimulator,
WellState& well_state,
const GroupState& group_state,
366 Eval getPerfCellPressure(
const FluidState& fs)
const;
372#include "WellInterface_impl.hpp"
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:38
Definition: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:184
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:70
Definition: WellInterfaceFluidSystem.hpp:47
Definition: WellInterfaceIndices.hpp:33
Definition: WellInterface.hpp:74
virtual void apply(const BVector &x, BVector &Ax) const =0
Ax = Ax - C D^-1 B x.
void updateWellStateRates(const Simulator &ebosSimulator, WellState &well_state, DeferredLogger &deferred_logger) const
Modify the well_state's rates if there is only one nonzero rate.
Definition: WellInterface_impl.hpp:1079
virtual std::vector< double > computeCurrentWellRates(const Simulator &ebosSimulator, DeferredLogger &deferred_logger) const =0
Compute well rates based on current reservoir conditions and well variables.
virtual void apply(BVector &r) const =0
r = r - C D^-1 Rw
virtual ~WellInterface()=default
Virtual destructor.
WellInterface(const Well &well, const ParallelWellInfo &pw_info, const int time_step, const ModelParameters ¶m, const RateConverterType &rate_converter, const int pvtRegionIdx, const int num_components, const int num_phases, const int index_of_well, const std::vector< PerforationData > &perf_data)
Constructor.
Definition: WellInterface_impl.hpp:41
virtual bool jacobianContainsWellContributions() const
Wether the Jacobian will also have well contributions in it.
Definition: WellInterface.hpp:228
virtual void recoverWellSolutionAndUpdateWellState(const SummaryState &summary_state, const BVector &x, WellState &well_state, DeferredLogger &deferred_logger)=0
using the solution x to recover the solution xw for wells and applying xw to update Well State
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