24#ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
25#define OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
27#include <ebos/eclproblem.hh>
28#include <opm/common/OpmLog/OpmLog.hpp>
37#include <unordered_map>
42#include <opm/input/eclipse/EclipseState/Runspec.hpp>
44#include <opm/input/eclipse/Schedule/Schedule.hpp>
45#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
46#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
47#include <opm/input/eclipse/Schedule/Group/Group.hpp>
49#include <opm/simulators/timestepping/SimulatorReport.hpp>
50#include <opm/simulators/flow/countGlobalCells.hpp>
51#include <opm/simulators/wells/BlackoilWellModelGeneric.hpp>
52#include <opm/simulators/wells/GasLiftSingleWell.hpp>
53#include <opm/simulators/wells/GasLiftWellState.hpp>
54#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
55#include <opm/simulators/wells/GasLiftStage2.hpp>
56#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
57#include <opm/simulators/wells/PerforationData.hpp>
58#include <opm/simulators/wells/VFPInjProperties.hpp>
59#include <opm/simulators/wells/VFPProdProperties.hpp>
60#include <opm/simulators/wells/WellState.hpp>
61#include <opm/simulators/wells/WGState.hpp>
64#include <opm/simulators/wells/WellInterface.hpp>
65#include <opm/simulators/wells/StandardWell.hpp>
66#include <opm/simulators/wells/MultisegmentWell.hpp>
67#include <opm/simulators/wells/WellGroupHelpers.hpp>
68#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
69#include <opm/simulators/wells/ParallelWellInfo.hpp>
70#include <opm/simulators/timestepping/gatherConvergenceReport.hpp>
71#include <dune/common/fmatrix.hh>
72#include <dune/istl/bcrsmatrix.hh>
73#include <dune/istl/matrixmatrix.hh>
75#include <opm/material/densead/Math.hpp>
77#include <opm/simulators/utils/DeferredLogger.hpp>
79namespace Opm::Properties {
81template<
class TypeTag,
class MyTypeTag>
83 using type = UndefinedProperty;
91 template<
typename TypeTag>
99 using Grid = GetPropType<TypeTag, Properties::Grid>;
100 using EquilGrid = GetPropType<TypeTag, Properties::EquilGrid>;
101 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
102 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
103 using Indices = GetPropType<TypeTag, Properties::Indices>;
104 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
105 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
106 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
107 using GlobalEqVector = GetPropType<TypeTag, Properties::GlobalEqVector>;
108 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
110 using GLiftOptWells =
typename BlackoilWellModelGeneric::GLiftOptWells;
111 using GLiftProdWells =
typename BlackoilWellModelGeneric::GLiftProdWells;
112 using GLiftWellStateMap =
113 typename BlackoilWellModelGeneric::GLiftWellStateMap;
114 using GLiftEclWells =
typename GasLiftGroupInfo::GLiftEclWells;
115 using GLiftSyncGroups =
typename GasLiftSingleWellGeneric::GLiftSyncGroups;
116 constexpr static std::size_t pressureVarIndex = GetPropType<TypeTag, Properties::Indices>::pressureSwitchIdx;
117 typedef typename BaseAuxiliaryModule<TypeTag>::NeighborSet NeighborSet;
119 static const int numEq = Indices::numEq;
120 static const int solventSaturationIdx = Indices::solventSaturationIdx;
121 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
122 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
123 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
124 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
128 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
129 typedef Dune::BlockVector<VectorBlockType> BVector;
131 typedef BlackOilPolymerModule<TypeTag> PolymerModule;
132 typedef BlackOilMICPModule<TypeTag> MICPModule;
136 SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;
140 AverageRegionalPressure<FluidSystem, std::vector<int> >;
145 void initWellContainer(
const int reportStepIdx)
override;
150 unsigned numDofs()
const override
154 void addNeighbors(std::vector<NeighborSet>& neighbors)
const override;
156 void applyInitial()
override
159 void linearize(SparseMatrixAdapter& jacobian, GlobalEqVector& res)
override;
161 void postSolve(GlobalEqVector& deltaX)
override
163 recoverWellSolutionAndUpdateWellState(deltaX);
170 template <
class Restarter>
171 void deserialize(Restarter& )
180 template <
class Restarter>
188 beginReportStep(ebosSimulator_.episodeIndex());
191 void beginTimeStep();
193 void beginIteration()
195 assemble(ebosSimulator_.model().newtonMethod().numIterations(),
196 ebosSimulator_.timeStepSize());
204 timeStepSucceeded(ebosSimulator_.time(), ebosSimulator_.timeStepSize());
212 void computeTotalRatesForDof(RateVector& rate,
213 unsigned globalIdx)
const;
215 template <
class Context>
216 void computeTotalRatesForDof(RateVector& rate,
217 const Context& context,
219 unsigned timeIdx)
const;
222 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
224 using BlackoilWellModelGeneric::initFromRestartFile;
225 void initFromRestartFile(
const RestartValue& restartValues)
227 initFromRestartFile(restartValues,
228 this->ebosSimulator_.vanguard().transferWTestState(),
233 data::Wells wellData()
const
235 auto wsrpt = this->wellState()
236 .report(ebosSimulator_.vanguard().globalCell().data(),
237 [
this](
const int well_index) ->
bool
239 return this->wasDynamicallyShutThisTimeStep(well_index);
242 this->assignWellTracerRates(wsrpt);
244 this->assignWellGuideRates(wsrpt, this->reportStepIndex());
245 this->assignShutConnections(wsrpt, this->reportStepIndex());
251 void apply( BVector& r)
const;
254 void apply(
const BVector& x, BVector& Ax)
const;
257 void getWellContributions(WellContributions& x)
const;
260 void applyScaleAdd(
const Scalar alpha,
const BVector& x, BVector& Ax)
const;
263 ConvergenceReport getWellConvergence(
const std::vector<Scalar>& B_avg,
const bool checkWellGroupControls =
false)
const;
265 const SimulatorReportSingle& lastReport()
const;
267 void addWellContributions(SparseMatrixAdapter& jacobian)
const;
270 void beginReportStep(
const int time_step);
272 void updatePerforationIntensiveQuantities();
281 void prepareTimeStep(DeferredLogger& deferred_logger);
282 void initPrimaryVariablesEvaluation()
const;
283 bool shouldBalanceNetwork(
const int reportStepIndex,
const int iterationIdx)
const;
284 std::tuple<bool, bool, double> updateWellControls(DeferredLogger& deferred_logger);
286 void updateAndCommunicate(
const int reportStepIdx,
287 const int iterationIdx,
288 DeferredLogger& deferred_logger);
290 bool updateGroupControls(
const Group& group,
291 DeferredLogger& deferred_logger,
292 const int reportStepIdx,
293 const int iterationIdx);
295 WellInterfacePtr getWell(
const std::string& well_name)
const;
296 bool hasWell(
const std::string& well_name)
const;
298 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
300 int numLocalWellsEnd()
const;
302 void addWellPressureEquations(PressureMatrix& jacobian,
const BVector& weights,
const bool use_well_weights)
const;
304 std::vector<std::vector<int>> getMaxWellConnections()
const;
306 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const;
308 void initGliftEclWellMap(GLiftEclWells &ecl_well_map);
313 return well_container_;
316 int numLocalNonshutWells()
const;
319 Simulator& ebosSimulator_;
322 std::vector<WellInterfacePtr > well_container_{};
324 std::vector<bool> is_cell_perforated_{};
326 void initializeWellState(
const int timeStepIdx,
327 const SummaryState& summaryState);
330 void createWellContainer(
const int time_step)
override;
333 createWellPointer(
const int wellID,
334 const int time_step)
const;
336 template <
typename WellType>
337 std::unique_ptr<WellType>
338 createTypedWellPointer(
const int wellID,
339 const int time_step)
const;
341 WellInterfacePtr createWellForWellTest(
const std::string& well_name,
const int report_step, DeferredLogger& deferred_logger)
const;
344 const ModelParameters param_;
345 size_t global_num_cells_{};
347 size_t local_num_cells_{};
349 std::vector<double> depth_{};
350 bool alternative_well_rate_init_{};
352 std::unique_ptr<RateConverterType> rateConverter_{};
353 std::unique_ptr<AverageRegionalPressureType> regionalAveragePressureCalculator_{};
356 SimulatorReportSingle last_report_{};
359 mutable BVector scaleAddRes_{};
361 std::vector<Scalar> B_avg_{};
363 const Grid& grid()
const
364 {
return ebosSimulator_.vanguard().grid(); }
366 const EquilGrid& equilGrid()
const
367 {
return ebosSimulator_.vanguard().equilGrid(); }
369 const EclipseState& eclState()
const
370 {
return ebosSimulator_.vanguard().eclState(); }
374 void assemble(
const int iterationIdx,
376 bool assembleImpl(
const int iterationIdx,
378 const std::size_t recursion_level,
379 DeferredLogger& local_deferredLogger);
382 void timeStepSucceeded(
const double& simulationTime,
const double dt);
385 void endReportStep();
389 void recoverWellSolutionAndUpdateWellState(
const BVector& x);
392 void updatePrimaryVariables(DeferredLogger& deferred_logger);
394 void updateAverageFormationFactor();
396 void computePotentials(
const std::size_t widx,
397 const WellState& well_state_copy,
398 std::string& exc_msg,
399 ExceptionType::ExcEnum& exc_type,
400 DeferredLogger& deferred_logger)
override;
402 const std::vector<double>& wellPerfEfficiencyFactors()
const;
404 void calculateProductivityIndexValuesShutWells(
const int reportStepIdx, DeferredLogger& deferred_logger)
override;
405 void calculateProductivityIndexValues(DeferredLogger& deferred_logger)
override;
406 void calculateProductivityIndexValues(
const WellInterface<TypeTag>* wellPtr,
407 DeferredLogger& deferred_logger);
410 int numComponents()
const;
412 int reportStepIndex()
const;
414 void assembleWellEq(
const double dt, DeferredLogger& deferred_logger);
416 bool maybeDoGasLiftOptimize(DeferredLogger& deferred_logger);
418 void gasLiftOptimizationStage1(DeferredLogger& deferred_logger,
419 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
420 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map);
423 void gasLiftOptimizationStage1SingleWell(WellInterface<TypeTag> *well,
424 DeferredLogger& deferred_logger,
425 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
426 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map,
427 GLiftSyncGroups& groups_to_sync);
429 void extractLegacyCellPvtRegionIndex_();
431 void extractLegacyDepth_();
434 void updateWellTestState(
const double& simulationTime, WellTestState& wellTestState)
const;
436 void wellTesting(
const int timeStepIdx,
const double simulationTime, DeferredLogger& deferred_logger);
438 void calcRates(
const int fipnum,
440 std::vector<double>& resv_coeff)
override;
442 void calcInjRates(
const int fipnum,
444 std::vector<double>& resv_coeff)
override;
446 void computeWellTemperature();
448 void assignWellTracerRates(data::Wells& wsrpt)
const;
451 return ebosSimulator_.vanguard().compressedIndexForInterior(cartesian_cell_idx);
461#include "BlackoilWellModel_impl.hpp"
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Class for handling the blackoil well model.
Definition: BlackoilWellModelGeneric.hpp:72
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:94
void serialize(Restarter &)
This method writes the complete state of the well to the harddisk.
Definition: BlackoilWellModel.hpp:181
const std::vector< WellInterfacePtr > & localNonshutWells() const
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:311
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:1445
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1595
int compressedIndexForInterior(int cartesian_cell_idx) const override
get compressed index for interior cells (-1, otherwise
Definition: BlackoilWellModel.hpp:450
Definition: GasLiftSingleWell.hpp:39
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:68
Computes hydrocarbon weighed average pressures over regions.
Definition: RegionAverageCalculator.hpp:62
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 use_multisegment_well_
Whether to use MultisegmentWell to handle multisegment wells it is something temporary before the mul...
Definition: BlackoilModelParametersEbos.hpp:408
Definition: BlackoilPhases.hpp:46
Definition: BlackoilWellModel.hpp:82