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>
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>;
109 using GLiftOptWells =
typename BlackoilWellModelGeneric::GLiftOptWells;
110 using GLiftProdWells =
typename BlackoilWellModelGeneric::GLiftProdWells;
111 using GLiftWellStateMap =
112 typename BlackoilWellModelGeneric::GLiftWellStateMap;
113 using GLiftEclWells =
typename GasLiftGroupInfo::GLiftEclWells;
114 using GLiftSyncGroups =
typename GasLiftSingleWellGeneric::GLiftSyncGroups;
116 typedef typename BaseAuxiliaryModule<TypeTag>::NeighborSet NeighborSet;
118 static const int numEq = Indices::numEq;
119 static const int solventSaturationIdx = Indices::solventSaturationIdx;
120 static constexpr
bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
121 static constexpr
bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
122 static constexpr
bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
123 static constexpr
bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
127 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
128 typedef Dune::BlockVector<VectorBlockType> BVector;
130 typedef Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
132 typedef BlackOilPolymerModule<TypeTag> PolymerModule;
133 typedef BlackOilMICPModule<TypeTag> MICPModule;
137 SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;
141 AverageRegionalPressure<FluidSystem, std::vector<int> >;
146 void initWellContainer()
override;
151 unsigned numDofs()
const override
155 void addNeighbors(std::vector<NeighborSet>& neighbors)
const override;
157 void applyInitial()
override
160 void linearize(SparseMatrixAdapter& jacobian, GlobalEqVector& res)
override;
162 void postSolve(GlobalEqVector& deltaX)
override
164 recoverWellSolutionAndUpdateWellState(deltaX);
171 template <
class Restarter>
172 void deserialize(Restarter& )
181 template <
class Restarter>
189 beginReportStep(ebosSimulator_.episodeIndex());
192 void beginTimeStep();
194 void beginIteration()
196 assemble(ebosSimulator_.model().newtonMethod().numIterations(),
197 ebosSimulator_.timeStepSize());
205 timeStepSucceeded(ebosSimulator_.time(), ebosSimulator_.timeStepSize());
213 template <
class Context>
214 void computeTotalRatesForDof(RateVector& rate,
215 const Context& context,
217 unsigned timeIdx)
const;
220 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
222 using BlackoilWellModelGeneric::initFromRestartFile;
223 void initFromRestartFile(
const RestartValue& restartValues)
225 initFromRestartFile(restartValues,
226 this->ebosSimulator_.vanguard().transferWTestState(),
227 UgGridHelpers::numCells(grid()),
231 data::Wells wellData()
const
233 auto wsrpt = this->wellState()
234 .report(UgGridHelpers::globalCell(this->grid()),
235 [
this](
const int well_index) ->
bool
237 return this->wasDynamicallyShutThisTimeStep(well_index);
240 this->assignWellTracerRates(wsrpt);
242 this->assignWellGuideRates(wsrpt, this->reportStepIndex());
243 this->assignShutConnections(wsrpt, this->reportStepIndex());
249 void apply( BVector& r)
const;
252 void apply(
const BVector& x, BVector& Ax)
const;
254 #if HAVE_CUDA || HAVE_OPENCL
256 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 checkGroupConvergence =
false)
const;
265 const SimulatorReportSingle& lastReport()
const;
267 void addWellContributions(SparseMatrixAdapter& jacobian)
const
269 for (
const auto& well: well_container_ ) {
270 well->addWellContributions(jacobian);
275 void beginReportStep(
const int time_step);
277 void updatePerforationIntensiveQuantities();
286 void prepareTimeStep(DeferredLogger& deferred_logger);
287 void initPrimaryVariablesEvaluation()
const;
288 void updateWellControls(DeferredLogger& deferred_logger,
const bool checkGroupControls);
290 void updateAndCommunicate(
const int reportStepIdx,
291 const int iterationIdx,
292 DeferredLogger& deferred_logger);
294 WellInterfacePtr getWell(
const std::string& well_name)
const;
295 bool hasWell(
const std::string& well_name)
const;
297 void initGliftEclWellMap(GLiftEclWells &ecl_well_map);
302 return well_container_;
306 Simulator& ebosSimulator_;
309 std::vector<WellInterfacePtr > well_container_{};
311 std::vector<bool> is_cell_perforated_{};
313 void initializeWellState(
const int timeStepIdx,
314 const SummaryState& summaryState);
317 void createWellContainer(
const int time_step)
override;
320 createWellPointer(
const int wellID,
321 const int time_step)
const;
323 template <
typename WellType>
324 std::unique_ptr<WellType>
325 createTypedWellPointer(
const int wellID,
326 const int time_step)
const;
328 WellInterfacePtr createWellForWellTest(
const std::string& well_name,
const int report_step, DeferredLogger& deferred_logger)
const;
331 const ModelParameters param_;
332 size_t global_num_cells_{};
334 size_t local_num_cells_{};
336 std::vector<double> depth_{};
337 bool alternative_well_rate_init_{};
339 std::unique_ptr<RateConverterType> rateConverter_{};
340 std::unique_ptr<AverageRegionalPressureType> regionalAveragePressureCalculator_{};
343 SimulatorReportSingle last_report_{};
346 mutable BVector scaleAddRes_{};
348 std::vector<Scalar> B_avg_{};
350 const Grid& grid()
const
351 {
return ebosSimulator_.vanguard().grid(); }
353 const EclipseState& eclState()
const
354 {
return ebosSimulator_.vanguard().eclState(); }
358 void assemble(
const int iterationIdx,
362 void timeStepSucceeded(
const double& simulationTime,
const double dt);
365 void endReportStep();
369 void recoverWellSolutionAndUpdateWellState(
const BVector& x);
372 void updatePrimaryVariables(DeferredLogger& deferred_logger);
374 void updateAverageFormationFactor();
376 void computePotentials(
const std::size_t widx,
377 const WellState& well_state_copy,
378 std::string& exc_msg,
379 ExceptionType::ExcEnum& exc_type,
380 DeferredLogger& deferred_logger)
override;
382 const std::vector<double>& wellPerfEfficiencyFactors()
const;
384 void calculateProductivityIndexValuesShutWells(
const int reportStepIdx, DeferredLogger& deferred_logger)
override;
385 void calculateProductivityIndexValues(DeferredLogger& deferred_logger)
override;
386 void calculateProductivityIndexValues(
const WellInterface<TypeTag>* wellPtr,
387 DeferredLogger& deferred_logger);
390 int numComponents()
const;
392 int reportStepIndex()
const;
394 void assembleWellEq(
const double dt, DeferredLogger& deferred_logger);
396 bool maybeDoGasLiftOptimize(DeferredLogger& deferred_logger);
398 void gasLiftOptimizationStage1(DeferredLogger& deferred_logger,
399 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
400 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map);
403 void gasLiftOptimizationStage1SingleWell(WellInterface<TypeTag> *well,
404 DeferredLogger& deferred_logger,
405 GLiftProdWells &prod_wells, GLiftOptWells &glift_wells,
406 GasLiftGroupInfo &group_info, GLiftWellStateMap &state_map,
407 GLiftSyncGroups& groups_to_sync);
409 void extractLegacyCellPvtRegionIndex_();
411 void extractLegacyDepth_();
414 void updateWellTestState(
const double& simulationTime, WellTestState& wellTestState)
const;
416 void wellTesting(
const int timeStepIdx,
const double simulationTime, DeferredLogger& deferred_logger);
418 void calcRates(
const int fipnum,
420 std::vector<double>& resv_coeff)
override;
422 void calcInjRates(
const int fipnum,
424 std::vector<double>& resv_coeff)
override;
426 void computeWellTemperature();
428 void assignWellTracerRates(data::Wells& wsrpt)
const;
431 return ebosSimulator_.vanguard().compressedIndexForInterior(cartesian_cell_idx);
441 #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:182
const std::vector< WellInterfacePtr > & localNonshutWells()
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:300
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:1273
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1376
int compressedIndexForInterior(int cartesian_cell_idx) const override
get compressed index for interior cells (-1, otherwise
Definition: BlackoilWellModel.hpp:430
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