24#ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
25#define OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
27#include <ebos/eclproblem.hh>
28#include <opm/common/OpmLog/OpmLog.hpp>
38#include <unordered_map>
41#include <opm/input/eclipse/Schedule/Schedule.hpp>
42#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
43#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
44#include <opm/input/eclipse/Schedule/Group/Group.hpp>
46#include <opm/simulators/timestepping/SimulatorReport.hpp>
47#include <opm/simulators/flow/countGlobalCells.hpp>
48#include <opm/simulators/wells/BlackoilWellModelGeneric.hpp>
49#include <opm/simulators/wells/BlackoilWellModelGuideRates.hpp>
50#include <opm/simulators/wells/GasLiftSingleWell.hpp>
51#include <opm/simulators/wells/GasLiftWellState.hpp>
52#include <opm/simulators/wells/GasLiftSingleWellGeneric.hpp>
53#include <opm/simulators/wells/GasLiftStage2.hpp>
54#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
55#include <opm/simulators/wells/PerforationData.hpp>
56#include <opm/simulators/wells/VFPInjProperties.hpp>
57#include <opm/simulators/wells/VFPProdProperties.hpp>
58#include <opm/simulators/wells/WellState.hpp>
59#include <opm/simulators/wells/WGState.hpp>
62#include <opm/simulators/wells/WellInterface.hpp>
63#include <opm/simulators/wells/StandardWell.hpp>
64#include <opm/simulators/wells/MultisegmentWell.hpp>
65#include <opm/simulators/wells/WellGroupHelpers.hpp>
66#include <opm/simulators/wells/WellProdIndexCalculator.hpp>
67#include <opm/simulators/wells/ParallelWellInfo.hpp>
68#include <opm/simulators/timestepping/gatherConvergenceReport.hpp>
69#include <dune/common/fmatrix.hh>
70#include <dune/istl/bcrsmatrix.hh>
71#include <dune/istl/matrixmatrix.hh>
73#include <opm/material/densead/Math.hpp>
75#include <opm/simulators/utils/DeferredLogger.hpp>
77namespace Opm::Properties {
79template<
class TypeTag,
class MyTypeTag>
81 using type = UndefinedProperty;
89 template<
typename TypeTag>
97 using Grid = GetPropType<TypeTag, Properties::Grid>;
98 using EquilGrid = GetPropType<TypeTag, Properties::EquilGrid>;
99 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
100 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
101 using Indices = GetPropType<TypeTag, Properties::Indices>;
102 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
103 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
104 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
105 using GlobalEqVector = GetPropType<TypeTag, Properties::GlobalEqVector>;
106 using SparseMatrixAdapter = GetPropType<TypeTag, Properties::SparseMatrixAdapter>;
108 using GLiftOptWells =
typename BlackoilWellModelGeneric::GLiftOptWells;
109 using GLiftProdWells =
typename BlackoilWellModelGeneric::GLiftProdWells;
110 using GLiftWellStateMap =
111 typename BlackoilWellModelGeneric::GLiftWellStateMap;
112 using GLiftEclWells =
typename GasLiftGroupInfo::GLiftEclWells;
113 using GLiftSyncGroups =
typename GasLiftSingleWellGeneric::GLiftSyncGroups;
114 constexpr static std::size_t pressureVarIndex = GetPropType<TypeTag, Properties::Indices>::pressureSwitchIdx;
115 typedef typename BaseAuxiliaryModule<TypeTag>::NeighborSet NeighborSet;
117 static const int numEq = Indices::numEq;
118 static const int solventSaturationIdx = Indices::solventSaturationIdx;
119 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
120 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
121 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
122 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
126 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
127 typedef Dune::BlockVector<VectorBlockType> BVector;
129 typedef BlackOilPolymerModule<TypeTag> PolymerModule;
130 typedef BlackOilMICPModule<TypeTag> MICPModule;
134 SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;
138 AverageRegionalPressure<FluidSystem, std::vector<int> >;
143 void initWellContainer(
const int reportStepIdx)
override;
148 unsigned numDofs()
const override
152 void addNeighbors(std::vector<NeighborSet>& neighbors)
const override;
154 void applyInitial()
override
157 void linearize(SparseMatrixAdapter& jacobian, GlobalEqVector& res)
override;
159 void postSolve(GlobalEqVector& deltaX)
override
161 recoverWellSolutionAndUpdateWellState(deltaX);
168 template <
class Restarter>
169 void deserialize(Restarter& )
178 template <
class Restarter>
186 OPM_TIMEBLOCK(beginEpsiode);
187 beginReportStep(ebosSimulator_.episodeIndex());
190 void beginTimeStep();
192 void beginIteration()
194 OPM_TIMEBLOCK(beginIteration);
195 assemble(ebosSimulator_.model().newtonMethod().numIterations(),
196 ebosSimulator_.timeStepSize());
204 OPM_TIMEBLOCK(endTimeStep);
205 timeStepSucceeded(ebosSimulator_.time(), ebosSimulator_.timeStepSize());
213 void computeTotalRatesForDof(RateVector& rate,
214 unsigned globalIdx)
const;
216 template <
class Context>
217 void computeTotalRatesForDof(RateVector& rate,
218 const Context& context,
220 unsigned timeIdx)
const;
223 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
225 using BlackoilWellModelGeneric::initFromRestartFile;
226 void initFromRestartFile(
const RestartValue& restartValues)
228 initFromRestartFile(restartValues,
229 this->ebosSimulator_.vanguard().transferWTestState(),
234 using BlackoilWellModelGeneric::prepareDeserialize;
235 void prepareDeserialize(
const int report_step)
237 prepareDeserialize(report_step, grid().size(0),
241 data::Wells wellData()
const
243 auto wsrpt = this->wellState()
244 .report(ebosSimulator_.vanguard().globalCell().data(),
245 [
this](
const int well_index) ->
bool
247 return this->wasDynamicallyShutThisTimeStep(well_index);
250 this->assignWellTracerRates(wsrpt);
252 BlackoilWellModelGuideRates(*this).assignWellGuideRates(wsrpt, this->reportStepIndex());
253 this->assignShutConnections(wsrpt, this->reportStepIndex());
259 void apply( BVector& r)
const;
262 void apply(
const BVector& x, BVector& Ax)
const;
265 void getWellContributions(WellContributions& x)
const;
268 void applyScaleAdd(
const Scalar alpha,
const BVector& x, BVector& Ax)
const;
271 ConvergenceReport getWellConvergence(
const std::vector<Scalar>& B_avg,
const bool checkWellGroupControls =
false)
const;
273 const SimulatorReportSingle& lastReport()
const;
275 void addWellContributions(SparseMatrixAdapter& jacobian)
const;
278 void addReservoirSourceTerms(GlobalEqVector& residual,
279 std::vector<typename SparseMatrixAdapter::MatrixBlock*>& diagMatAddress)
const;
282 void beginReportStep(
const int time_step);
284 void updatePerforationIntensiveQuantities();
293 void prepareTimeStep(DeferredLogger& deferred_logger);
294 void initPrimaryVariablesEvaluation()
const;
295 std::tuple<bool, bool, double> updateWellControls(DeferredLogger& deferred_logger);
297 void updateAndCommunicate(
const int reportStepIdx,
298 const int iterationIdx,
299 DeferredLogger& deferred_logger);
301 bool updateGroupControls(
const Group& group,
302 DeferredLogger& deferred_logger,
303 const int reportStepIdx,
304 const int iterationIdx);
306 WellInterfacePtr getWell(
const std::string& well_name)
const;
307 bool hasWell(
const std::string& well_name)
const;
309 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
311 int numLocalWellsEnd()
const;
313 void addWellPressureEquations(PressureMatrix& jacobian,
const BVector& weights,
const bool use_well_weights)
const;
315 std::vector<std::vector<int>> getMaxWellConnections()
const;
317 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const;
319 void initGliftEclWellMap(GLiftEclWells &ecl_well_map);
324 return well_container_;
327 int numLocalNonshutWells()
const;
330 Simulator& ebosSimulator_;
333 std::vector<WellInterfacePtr> well_container_{};
335 std::vector<bool> is_cell_perforated_{};
337 void initializeWellState(
const int timeStepIdx,
338 const SummaryState& summaryState);
341 void createWellContainer(
const int time_step)
override;
344 createWellPointer(
const int wellID,
345 const int time_step)
const;
347 template <
typename WellType>
348 std::unique_ptr<WellType>
349 createTypedWellPointer(
const int wellID,
350 const int time_step)
const;
352 WellInterfacePtr createWellForWellTest(
const std::string& well_name,
const int report_step, DeferredLogger& deferred_logger)
const;
355 const ModelParameters param_;
356 std::size_t global_num_cells_{};
358 std::size_t local_num_cells_{};
360 std::vector<double> depth_{};
361 bool alternative_well_rate_init_{};
363 std::unique_ptr<RateConverterType> rateConverter_{};
364 std::map<std::string, std::unique_ptr<AverageRegionalPressureType>> regionalAveragePressureCalculator_{};
367 SimulatorReportSingle last_report_{};
370 mutable BVector scaleAddRes_{};
372 std::vector<Scalar> B_avg_{};
374 const Grid& grid()
const
375 {
return ebosSimulator_.vanguard().grid(); }
377 const EquilGrid& equilGrid()
const
378 {
return ebosSimulator_.vanguard().equilGrid(); }
380 const EclipseState& eclState()
const
381 {
return ebosSimulator_.vanguard().eclState(); }
385 void assemble(
const int iterationIdx,
387 bool assembleImpl(
const int iterationIdx,
389 const std::size_t recursion_level,
390 DeferredLogger& local_deferredLogger);
393 void timeStepSucceeded(
const double& simulationTime,
const double dt);
396 void endReportStep();
400 void recoverWellSolutionAndUpdateWellState(
const BVector& x);
403 void updatePrimaryVariables(DeferredLogger& deferred_logger);
405 void updateAverageFormationFactor();
407 void computePotentials(
const std::size_t widx,
408 const WellState& well_state_copy,
409 std::string& exc_msg,
410 ExceptionType::ExcEnum& exc_type,
411 DeferredLogger& deferred_logger)
override;
413 const std::vector<double>& wellPerfEfficiencyFactors()
const;
415 void calculateProductivityIndexValuesShutWells(
const int reportStepIdx, DeferredLogger& deferred_logger)
override;
416 void calculateProductivityIndexValues(DeferredLogger& deferred_logger)
override;
417 void calculateProductivityIndexValues(
const WellInterface<TypeTag>* wellPtr,
418 DeferredLogger& deferred_logger);
421 int numComponents()
const;
423 int reportStepIndex()
const;
425 void assembleWellEq(
const double dt, DeferredLogger& deferred_logger);
427 bool maybeDoGasLiftOptimize(DeferredLogger& deferred_logger);
429 void gasLiftOptimizationStage1(DeferredLogger& deferred_logger,
430 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
431 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map);
434 void gasLiftOptimizationStage1SingleWell(WellInterface<TypeTag> *well,
435 DeferredLogger& deferred_logger,
436 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
437 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map,
438 GLiftSyncGroups& groups_to_sync);
440 void extractLegacyCellPvtRegionIndex_();
442 void extractLegacyDepth_();
445 void updateWellTestState(
const double& simulationTime, WellTestState& wellTestState)
const;
447 void wellTesting(
const int timeStepIdx,
const double simulationTime, DeferredLogger& deferred_logger);
449 void calcRates(
const int fipnum,
451 const std::vector<double>& production_rates,
452 std::vector<double>& resv_coeff)
override;
454 void calcInjRates(
const int fipnum,
456 std::vector<double>& resv_coeff)
override;
458 void computeWellTemperature();
460 void assignWellTracerRates(data::Wells& wsrpt)
const;
463 return ebosSimulator_.vanguard().compressedIndexForInterior(cartesian_cell_idx);
473#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:73
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:92
void serialize(Restarter &)
This method writes the complete state of the well to the harddisk.
Definition: BlackoilWellModel.hpp:179
const std::vector< WellInterfacePtr > & localNonshutWells() const
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:322
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:1438
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1572
int compressedIndexForInterior(int cartesian_cell_idx) const override
get compressed index for interior cells (-1, otherwise
Definition: BlackoilWellModel.hpp:462
Definition: GasLiftSingleWell.hpp:38
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:70
Computes hydrocarbon weighed average pressures over regions.
Definition: RegionAverageCalculator.hpp:65
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:80