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/parser/eclipse/EclipseState/Runspec.hpp>
44 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
45 #include <opm/parser/eclipse/EclipseState/Schedule/Well/WellTestState.hpp>
46 #include <opm/parser/eclipse/EclipseState/Schedule/Group/GuideRate.hpp>
47 #include <opm/parser/eclipse/EclipseState/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>
79 namespace Opm::Properties {
81 template<
class TypeTag,
class MyTypeTag>
83 using type = UndefinedProperty;
91 template<
typename TypeTag>
99 using Grid = GetPropType<TypeTag, Properties::Grid>;
100 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
101 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
102 using Indices = GetPropType<TypeTag, Properties::Indices>;
103 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
104 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
105 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
106 using GlobalEqVector = GetPropType<TypeTag, Properties::GlobalEqVector>;
107 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;
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 Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
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()
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 template <
class Context>
213 void computeTotalRatesForDof(RateVector& rate,
214 const Context& context,
216 unsigned timeIdx)
const;
219 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
221 using BlackoilWellModelGeneric::initFromRestartFile;
222 void initFromRestartFile(
const RestartValue& restartValues)
224 initFromRestartFile(restartValues,
225 this->ebosSimulator_.vanguard().transferWTestState(),
226 UgGridHelpers::numCells(grid()),
230 data::Wells wellData()
const
232 auto wsrpt = this->wellState()
233 .report(UgGridHelpers::globalCell(grid()),
234 [
this](
const int well_ndex) ->
bool
236 return this->wasDynamicallyShutThisTimeStep(well_ndex);
239 this->assignWellTracerRates(wsrpt);
241 this->assignWellGuideRates(wsrpt);
242 this->assignShutConnections(wsrpt, this->reportStepIndex());
248 void apply( BVector& r)
const;
251 void apply(
const BVector& x, BVector& Ax)
const;
253 #if HAVE_CUDA || HAVE_OPENCL
255 void getWellContributions(WellContributions& x)
const;
259 void applyScaleAdd(
const Scalar alpha,
const BVector& x, BVector& Ax)
const;
262 ConvergenceReport getWellConvergence(
const std::vector<Scalar>& B_avg,
const bool checkGroupConvergence =
false)
const;
264 const SimulatorReportSingle& lastReport()
const;
266 void addWellContributions(SparseMatrixAdapter& jacobian)
const
268 for (
const auto& well: well_container_ ) {
269 well->addWellContributions(jacobian);
274 void beginReportStep(
const int time_step);
276 void updatePerforationIntensiveQuantities();
285 void prepareTimeStep(DeferredLogger& deferred_logger);
286 void initPrimaryVariablesEvaluation()
const;
287 void updateWellControls(DeferredLogger& deferred_logger,
const bool checkGroupControls);
289 void updateAndCommunicate(
const int reportStepIdx,
290 const int iterationIdx,
291 DeferredLogger& deferred_logger);
293 WellInterfacePtr getWell(
const std::string& well_name)
const;
294 bool hasWell(
const std::string& well_name)
const;
296 void initGliftEclWellMap(GLiftEclWells &ecl_well_map);
301 return well_container_;
305 Simulator& ebosSimulator_;
308 std::vector<WellInterfacePtr > well_container_{};
310 std::vector<bool> is_cell_perforated_{};
312 void initializeWellState(
const int timeStepIdx,
313 const SummaryState& summaryState);
316 void createWellContainer(
const int time_step)
override;
319 createWellPointer(
const int wellID,
320 const int time_step)
const;
322 template <
typename WellType>
323 std::unique_ptr<WellType>
324 createTypedWellPointer(
const int wellID,
325 const int time_step)
const;
327 WellInterfacePtr createWellForWellTest(
const std::string& well_name,
const int report_step, DeferredLogger& deferred_logger)
const;
330 const ModelParameters param_;
331 size_t global_num_cells_{};
333 size_t local_num_cells_{};
335 std::vector<double> depth_{};
336 bool alternative_well_rate_init_{};
338 std::unique_ptr<RateConverterType> rateConverter_{};
339 std::unique_ptr<AverageRegionalPressureType> regionalAveragePressureCalculator_{};
342 SimulatorReportSingle last_report_{};
345 mutable BVector scaleAddRes_{};
347 std::vector<Scalar> B_avg_{};
349 const Grid& grid()
const
350 {
return ebosSimulator_.vanguard().grid(); }
352 const EclipseState& eclState()
const
353 {
return ebosSimulator_.vanguard().eclState(); }
357 void assemble(
const int iterationIdx,
361 void timeStepSucceeded(
const double& simulationTime,
const double dt);
364 void endReportStep();
368 void recoverWellSolutionAndUpdateWellState(
const BVector& x);
371 void updatePrimaryVariables(DeferredLogger& deferred_logger);
373 void setupCartesianToCompressed_(
const int* global_cell,
int local_num__cells);
375 void updateAverageFormationFactor();
377 void computePotentials(
const std::size_t widx,
378 const WellState& well_state_copy,
379 std::string& exc_msg,
380 ExceptionType::ExcEnum& exc_type,
381 DeferredLogger& deferred_logger)
override;
383 const std::vector<double>& wellPerfEfficiencyFactors()
const;
385 void calculateProductivityIndexValuesShutWells(
const int reportStepIdx, DeferredLogger& deferred_logger)
override;
386 void calculateProductivityIndexValues(DeferredLogger& deferred_logger)
override;
387 void calculateProductivityIndexValues(
const WellInterface<TypeTag>* wellPtr,
388 DeferredLogger& deferred_logger);
391 int numComponents()
const;
393 int reportStepIndex()
const;
395 void assembleWellEq(
const double dt, DeferredLogger& deferred_logger);
397 void maybeDoGasLiftOptimize(DeferredLogger& deferred_logger);
399 void gasLiftOptimizationStage1(DeferredLogger& deferred_logger,
400 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
401 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map);
403 void extractLegacyCellPvtRegionIndex_();
405 void extractLegacyDepth_();
408 void updateWellTestState(
const double& simulationTime, WellTestState& wellTestState)
const;
410 void wellTesting(
const int timeStepIdx,
const double simulationTime, DeferredLogger& deferred_logger);
412 void calcRates(
const int fipnum,
414 std::vector<double>& resv_coeff)
override;
416 void calcInjRates(
const int fipnum,
418 std::vector<double>& resv_coeff)
override;
420 void computeWellTemperature();
422 void assignWellTracerRates(data::Wells& wsrpt)
const;
425 BlackoilWellModel(Simulator& ebosSimulator,
const PhaseUsage& pu);
431 #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:69
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()
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:299
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:1198
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1301
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:26
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: BlackoilWellModel.hpp:82