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/parser/eclipse/EclipseState/Schedule/Well/Well.hpp>
32 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WellTestState.hpp>
34 #include <opm/core/props/BlackoilPhases.hpp>
36 #include <opm/simulators/wells/WellProdIndexCalculator.hpp>
37 #include <opm/simulators/wells/WellState.hpp>
42 template<
typename TypeTag>
class GasLiftSingleWell;
43 template<
typename TypeTag>
class BlackoilWellModel;
45 #include <opm/simulators/wells/GasLiftGroupInfo.hpp>
46 #include <opm/simulators/wells/GasLiftSingleWell.hpp>
47 #include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
48 #include <opm/simulators/wells/BlackoilWellModel.hpp>
49 #include <opm/simulators/flow/BlackoilModelParametersEbos.hpp>
51 #include <opm/simulators/utils/DeferredLogger.hpp>
53 #include<dune/common/fmatrix.hh>
54 #include<dune/istl/bcrsmatrix.hh>
55 #include<dune/istl/matrixmatrix.hh>
57 #include <opm/material/densead/Evaluation.hpp>
59 #include <opm/simulators/wells/WellInterfaceIndices.hpp>
67 template<
typename TypeTag>
69 GetPropType<TypeTag, Properties::Indices>,
70 GetPropType<TypeTag, Properties::Scalar>>
76 using Grid = GetPropType<TypeTag, Properties::Grid>;
77 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
78 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
79 using Indices = GetPropType<TypeTag, Properties::Indices>;
80 using IntensiveQuantities = GetPropType<TypeTag, Properties::IntensiveQuantities>;
81 using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
82 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
83 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
85 using GLiftOptWells =
typename BlackoilWellModel<TypeTag>::GLiftOptWells;
86 using GLiftProdWells =
typename BlackoilWellModel<TypeTag>::GLiftProdWells;
87 using GLiftWellStateMap =
88 typename BlackoilWellModel<TypeTag>::GLiftWellStateMap;
89 using GLiftSyncGroups =
typename GasLiftSingleWellGeneric::GLiftSyncGroups;
91 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
93 using VectorBlockType = Dune::FieldVector<Scalar, Indices::numEq>;
94 using MatrixBlockType = Dune::FieldMatrix<Scalar, Indices::numEq, Indices::numEq>;
95 using BVector = Dune::BlockVector<VectorBlockType>;
96 using Eval = DenseAd::Evaluation<Scalar, Indices::numEq>;
98 using RateConverterType =
106 static constexpr
bool has_solvent = getPropValue<TypeTag, Properties::EnableSolvent>();
107 static constexpr
bool has_zFraction = getPropValue<TypeTag, Properties::EnableExtbo>();
108 static constexpr
bool has_polymer = getPropValue<TypeTag, Properties::EnablePolymer>();
109 static constexpr
bool has_energy = getPropValue<TypeTag, Properties::EnableEnergy>();
110 static const bool has_temperature = getPropValue<TypeTag, Properties::EnableTemperature>();
112 static constexpr
bool has_polymermw = getPropValue<TypeTag, Properties::EnablePolymerMW>();
113 static constexpr
bool has_foam = getPropValue<TypeTag, Properties::EnableFoam>();
114 static constexpr
bool has_brine = getPropValue<TypeTag, Properties::EnableBrine>();
115 static constexpr
bool has_micp = getPropValue<TypeTag, Properties::EnableMICP>();
118 using FluidState = BlackOilFluidState<Eval,
122 Indices::compositionSwitchIdx >= 0,
124 Indices::numPhases >;
130 const RateConverterType& rate_converter,
131 const int pvtRegionIdx,
132 const int num_components,
133 const int num_phases,
134 const int index_of_well,
135 const std::vector<PerforationData>& perf_data);
140 virtual void init(
const PhaseUsage* phase_usage_arg,
141 const std::vector<double>& depth_arg,
142 const double gravity_arg,
144 const std::vector< Scalar >& B_avg);
146 virtual void initPrimaryVariablesEvaluation()
const = 0;
152 void assembleWellEq(
const Simulator& ebosSimulator,
158 virtual void gasLiftOptimizationStage1 (
161 const Simulator& ebosSimulator,
163 GLiftProdWells& prod_wells,
164 GLiftOptWells& glift_wells,
165 GLiftWellStateMap& state_map,
167 GLiftSyncGroups &sync_groups
177 virtual void apply(
const BVector& x, BVector& Ax)
const = 0;
180 virtual void apply(BVector& r)
const = 0;
183 virtual void computeWellPotentials(
const Simulator& ebosSimulator,
185 std::vector<double>& well_potentials,
188 virtual void updateWellStateWithTarget(
const Simulator& ebos_simulator,
193 enum class IndividualOrGroup { Individual, Group, Both };
194 bool updateWellControl(
const Simulator& ebos_simulator,
195 const IndividualOrGroup iog,
202 virtual void calculateExplicitQuantities(
const Simulator& ebosSimulator,
206 virtual void updateProductivityIndex(
const Simulator& ebosSimulator,
218 virtual void addWellContributions(SparseMatrixAdapter&)
const = 0;
220 void addCellRates(RateVector& rates,
int cellIdx)
const;
222 Scalar volumetricSurfaceRateForConnection(
int cellIdx,
int phaseIdx)
const;
224 template <
class EvalWell>
225 Eval restrictEval(
const EvalWell& in)
const
228 out.setValue(in.value());
229 for (
int eqIdx = 0; eqIdx < Indices::numEq; ++eqIdx) {
230 out.setDerivative(eqIdx, in.derivative(eqIdx));
237 void wellTesting(
const Simulator& simulator,
238 const double simulation_time,
239 WellState& well_state,
const GroupState& group_state, WellTestState& welltest_state,
240 DeferredLogger& deferred_logger);
242 void checkWellOperability(
const Simulator& ebos_simulator,
const WellState& well_state, DeferredLogger& deferred_logger);
246 void updateWellOperability(
const Simulator& ebos_simulator,
247 const WellState& well_state,
248 DeferredLogger& deferred_logger);
252 virtual void updateWaterThroughput(
const double dt, WellState& well_state)
const = 0;
266 void solveWellEquation(
const Simulator& ebosSimulator,
276 std::vector<RateVector> connectionRates_;
278 std::vector< Scalar > B_avg_;
280 bool changed_to_stopped_this_step_ =
false;
282 double wpolymer()
const;
284 double wfoam()
const;
286 double wsalt()
const;
288 double wmicrobes()
const;
290 double woxygen()
const;
292 double wurea()
const;
294 virtual double getRefDensity()
const = 0;
297 const std::vector<double>& compFrac()
const;
299 std::vector<double> initialWellRateFractions(
const Simulator& ebosSimulator,
const WellState& well_state)
const;
302 virtual void checkOperabilityUnderBHPLimitProducer(
const WellState& well_state,
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger) =0;
305 virtual void checkOperabilityUnderTHPLimitProducer(
const Simulator& ebos_simulator,
const WellState& well_state,
DeferredLogger& deferred_logger) =0;
307 virtual void updateIPR(
const Simulator& ebos_simulator,
DeferredLogger& deferred_logger)
const=0;
309 virtual void assembleWellEqWithoutIteration(
const Simulator& ebosSimulator,
311 const Well::InjectionControls& inj_controls,
312 const Well::ProductionControls& prod_controls,
318 virtual bool iterateWellEqWithControl(
const Simulator& ebosSimulator,
320 const Well::InjectionControls& inj_controls,
321 const Well::ProductionControls& prod_controls,
326 bool iterateWellEquations(
const Simulator& ebosSimulator,
332 bool solveWellForTesting(
const Simulator& ebosSimulator,
WellState& well_state,
const GroupState& group_state,
335 bool shutUnsolvableWells()
const;
340 #include "WellInterface_impl.hpp"
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: GasLiftGroupInfo.hpp:46
Definition: GasLiftSingleWell.hpp:46
Definition: GroupState.hpp:34
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:252
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:68
Definition: WellInterfaceFluidSystem.hpp:46
Definition: WellInterfaceIndices.hpp:35
Definition: WellInterface.hpp:71
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:991
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.
virtual void recoverWellSolutionAndUpdateWellState(const BVector &x, WellState &well_state, DeferredLogger &deferred_logger) const =0
using the solution x to recover the solution xw for wells and applying xw to update Well State
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:35
virtual bool jacobianContainsWellContributions() const
Wether the Jacobian will also have well contributions in it.
Definition: WellInterface.hpp:212
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:26
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParametersEbos.hpp:327
Definition: BlackoilPhases.hpp:45
Definition: WellInterfaceFluidSystem.hpp:126