BlackoilWellModel.hpp
Go to the documentation of this file.
1/*
2 Copyright 2016 SINTEF ICT, Applied Mathematics.
3 Copyright 2016 - 2017 Statoil ASA.
4 Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5 Copyright 2016 - 2018 IRIS AS
6
7 This file is part of the Open Porous Media project (OPM).
8
9 OPM is free software: you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation, either version 3 of the License, or
12 (at your option) any later version.
13
14 OPM is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18
19 You should have received a copy of the GNU General Public License
20 along with OPM. If not, see <http://www.gnu.org/licenses/>.
21*/
22
23#ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
24#define OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
25
26#if HAVE_MPI
27#define RESERVOIR_COUPLING_ENABLED
28#endif
29#ifdef RESERVOIR_COUPLING_ENABLED
32#endif
33
34#include <dune/common/fmatrix.hh>
35#include <dune/istl/bcrsmatrix.hh>
36#include <dune/istl/matrixmatrix.hh>
37
38#include <opm/common/OpmLog/OpmLog.hpp>
39
40#include <opm/input/eclipse/Schedule/Group/Group.hpp>
41#include <opm/input/eclipse/Schedule/Group/GuideRate.hpp>
42#include <opm/input/eclipse/Schedule/Schedule.hpp>
43#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
44
45#include <opm/material/densead/Math.hpp>
46
49
51
54
56
80
81#include <cstddef>
82#include <map>
83#include <memory>
84#include <string>
85#include <tuple>
86#include <vector>
87
88namespace Opm {
89
90template<class Scalar> class BlackoilWellModelNldd;
91template<class T> class SparseTable;
92
93#if COMPILE_GPU_BRIDGE
94template<class Scalar> class WellContributions;
95#endif
96
98 template<typename TypeTag>
99 class BlackoilWellModel : public WellConnectionAuxiliaryModule<TypeTag, BlackoilWellModel<TypeTag>>
100 , public BlackoilWellModelGeneric<GetPropType<TypeTag, Properties::Scalar>,
101 typename GetPropType<TypeTag, Properties::FluidSystem>::IndexTraitsType>
102 {
103 public:
104 // --------- Types ---------
116
118
119 using IndexTraits = typename FluidSystem::IndexTraitsType;
120
122
123 static const int numEq = Indices::numEq;
124 static const int solventSaturationIdx = Indices::solventSaturationIdx;
125 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
126 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
127 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
128 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
129
130 // TODO: where we should put these types, WellInterface or Well Model?
131 // or there is some other strategy, like TypeTag
132 using VectorBlockType = Dune::FieldVector<Scalar, numEq>;
133 using BVector = Dune::BlockVector<VectorBlockType>;
134
137
138 // For the conversion between the surface volume rate and reservoir voidage rate
139 using RateConverterType = RateConverter::
140 SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;
141
142 // For computing average pressured used by gpmaint
143 using AverageRegionalPressureType = RegionAverageCalculator::
144 AverageRegionalPressure<FluidSystem, std::vector<int> >;
145
147
149
150 void init();
151 void initWellContainer(const int reportStepIdx) override;
152
154 {
155 OPM_TIMEBLOCK(beginEpsiode);
156 beginReportStep(simulator_.episodeIndex());
157 }
158
159 void beginTimeStep();
160
162 {
163 OPM_TIMEBLOCK(beginIteration);
164 assemble(simulator_.model().newtonMethod().numIterations(),
165 simulator_.timeStepSize());
166 }
167
169 { }
170
172 {
173 OPM_TIMEBLOCK(endTimeStep);
174 timeStepSucceeded(simulator_.time(), simulator_.timeStepSize());
175 }
176
178 {
180 }
181
183 unsigned globalIdx) const;
184
185 template <class Context>
187 const Context& context,
188 unsigned spaceIdx,
189 unsigned timeIdx) const;
190
191
192 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
193
195 void initFromRestartFile(const RestartValue& restartValues)
196 {
197 initFromRestartFile(restartValues,
198 this->simulator_.vanguard().transferWTestState(),
199 grid().size(0),
201 this->simulator_.vanguard().enableDistributedWells());
202 }
203
205 void prepareDeserialize(const int report_step)
206 {
207 prepareDeserialize(report_step, grid().size(0),
209 this->simulator_.vanguard().enableDistributedWells());
210 }
211
212 data::Wells wellData() const
213 {
214 auto wsrpt = this->wellState()
215 .report(this->simulator_.vanguard().globalCell().data(),
216 [this](const int well_index)
217 {
218 return this->wasDynamicallyShutThisTimeStep(well_index);
219 });
220
222 .assignWellGuideRates(wsrpt, this->reportStepIndex());
223
224 this->assignWellTracerRates(wsrpt);
225
226 if (const auto& rspec = eclState().runspec();
227 rspec.co2Storage() || rspec.h2Storage())
228 {
229 // The gas reference density (surface condition) is the
230 // same for all PVT regions in CO2STORE/H2STORE runs so,
231 // for simplicity, we use region zero (0) here.
232
233 this->assignMassGasRate(wsrpt, FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, 0));
234 }
235
236 this->assignWellTargets(wsrpt);
237
238 this->assignDynamicWellStatus(wsrpt);
239
240 // Assigning (a subset of the) property values in shut
241 // connections should be the last step of wellData().
242 this->assignShutConnections(wsrpt, this->reportStepIndex());
243
244 return wsrpt;
245 }
246
247 data::WellBlockAveragePressures wellBlockAveragePressures() const
248 {
250 }
251
252#if COMPILE_GPU_BRIDGE
253 // accumulate the contributions of all Wells in the WellContributions object
254 void getWellContributions(WellContributions<Scalar>& x) const;
255#endif
256
257 // Check if well equations is converged.
258 ConvergenceReport getWellConvergence(const std::vector<Scalar>& B_avg, const bool checkWellGroupControls = false) const;
259
260 const SimulatorReportSingle& lastReport() const;
261
262 void addWellContributions(SparseMatrixAdapter& jacobian) const;
263
264 // add source from wells to the reservoir matrix
266 const std::vector<typename SparseMatrixAdapter::MatrixBlock*>& diagMatAddress) const;
267
268 // called at the beginning of a report step
269 void beginReportStep(const int time_step);
270
271 // it should be able to go to prepareTimeStep(), however, the updateWellControls()
272 // makes it a little more difficult. unless we introduce if (iterationIdx != 0) to avoid doing the above function
273 // twice at the beginning of the time step
276 void calculateExplicitQuantities(DeferredLogger& deferred_logger) const;
277 // some preparation work, mostly related to group control and RESV,
278 // at the beginning of each time step (Not report step)
279 void prepareTimeStep(DeferredLogger& deferred_logger);
280
281 bool
282 updateWellControls(DeferredLogger& deferred_logger);
283
284 std::tuple<bool, Scalar>
285 updateNetworks(const bool mandatory_network_balance, DeferredLogger& deferred_logger, const bool relax_network_tolerance = false);
286
287
288 void updateAndCommunicate(const int reportStepIdx,
289 const int iterationIdx,
290 DeferredLogger& deferred_logger);
291
292 bool updateGroupControls(const Group& group,
293 DeferredLogger& deferred_logger,
294 const int reportStepIdx,
295 const int iterationIdx);
296
297 WellInterfacePtr getWell(const std::string& well_name) const;
298
299 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<Scalar, 1, 1>>;
300
302 const BVector& weights,
303 const bool use_well_weights) const;
306 const BVector& weights,
307 const bool use_well_weights,
308 const int domainIndex) const
309 {
310 if (!nldd_) {
311 OPM_THROW(std::logic_error, "Attempt to access NLDD data without a NLDD solver");
312 }
313 return nldd_->addWellPressureEquations(jacobian,
314 weights,
315 use_well_weights,
316 domainIndex);
317 }
318
320 const std::vector<WellInterfacePtr>& localNonshutWells() const
321 {
322 return well_container_;
323 }
324
326 {
327 if (!nldd_) {
328 OPM_THROW(std::logic_error, "Attempt to access NLDD data without a NLDD solver");
329 }
330 return nldd_->well_local_cells();
331 }
332
333 const std::map<std::string, int>& well_domain() const
334 {
335 if (!nldd_) {
336 OPM_THROW(std::logic_error, "Attempt to access NLDD data without a NLDD solver");
337 }
338
339 return nldd_->well_domain();
340 }
341
342 auto begin() const { return well_container_.begin(); }
343 auto end() const { return well_container_.end(); }
344 bool empty() const { return well_container_.empty(); }
345
348
351
352 int compressedIndexForInterior(int cartesian_cell_idx) const override
353 {
354 return simulator_.vanguard().compressedIndexForInterior(cartesian_cell_idx);
355 }
356
357 int compressedIndexForInteriorLGR(const std::string& lgr_tag, const Connection& conn) const override
358 {
359 return simulator_.vanguard().compressedIndexForInteriorLGR(lgr_tag, conn);
360 }
361
362 // using the solution x to recover the solution xw for wells and applying
363 // xw to update Well State
365
366 // using the solution x to recover the solution xw for wells and applying
367 // xw to update Well State
369 const int domainIdx);
370
371 const Grid& grid() const
372 { return simulator_.vanguard().grid(); }
373
374 const Simulator& simulator() const
375 { return simulator_; }
376
378 { nldd_ = mod; }
379
380#ifdef RESERVOIR_COUPLING_ENABLED
382 return *(this->simulator_.reservoirCouplingMaster());
383 }
385 return *(this->simulator_.reservoirCouplingSlave());
386 }
388 return this->simulator_.reservoirCouplingMaster() != nullptr;
389 }
391 return this->simulator_.reservoirCouplingSlave() != nullptr;
392 }
394 {
395 this->guide_rate_handler_.setReservoirCouplingMaster(master);
396 }
398 {
399 this->guide_rate_handler_.setReservoirCouplingSlave(slave);
400 }
401 #endif
402 protected:
404
405 // a vector of all the wells.
406 std::vector<WellInterfacePtr> well_container_{};
407
408 std::vector<bool> is_cell_perforated_{};
409
410 void initializeWellState(const int timeStepIdx);
411
412 // create the well container
413 void createWellContainer(const int report_step) override;
414
416 createWellPointer(const int wellID,
417 const int report_step) const;
418
419 template <typename WellType>
420 std::unique_ptr<WellType>
421 createTypedWellPointer(const int wellID,
422 const int time_step) const;
423
424 WellInterfacePtr createWellForWellTest(const std::string& well_name,
425 const int report_step,
426 DeferredLogger& deferred_logger) const;
427
429 std::size_t global_num_cells_{};
430 // the number of the cells in the local grid
431 std::size_t local_num_cells_{};
433 std::vector<Scalar> depth_{};
435 std::map<std::string, Scalar> well_group_thp_calc_;
436 std::unique_ptr<RateConverterType> rateConverter_{};
437 std::map<std::string, std::unique_ptr<AverageRegionalPressureType>> regionalAveragePressureCalculator_{};
438
441
442 // Pre-step network solve at static reservoir conditions (group and well states might be updated)
443 void doPreStepNetworkRebalance(DeferredLogger& deferred_logger);
444
445 std::vector<Scalar> B_avg_{};
446
447 const EquilGrid& equilGrid() const
448 { return simulator_.vanguard().equilGrid(); }
449
450 const EclipseState& eclState() const
451 { return simulator_.vanguard().eclState(); }
452
453 // compute the well fluxes and assemble them in to the reservoir equations as source terms
454 // and in the well equations.
455 void assemble(const int iterationIdx,
456 const double dt);
457
458 // well controls and network pressures affect each other and are solved in an iterative manner.
459 // the function handles one iteration of updating well controls and network pressures.
460 // it is possible to decouple the update of well controls and network pressures further.
461 // the returned two booleans are {continue_due_to_network, well_group_control_changed}, respectively
462 std::tuple<bool, bool, Scalar> updateWellControlsAndNetworkIteration(const bool mandatory_network_balance,
463 const bool relax_network_tolerance,
464 const bool optimize_gas_lift,
465 const double dt,
466 DeferredLogger& local_deferredLogger);
467
468 bool updateWellControlsAndNetwork(const bool mandatory_network_balance,
469 const double dt,
470 DeferredLogger& local_deferredLogger);
471
472 bool computeWellGroupThp(const double dt, DeferredLogger& local_deferredLogger);
473
482 void initializeLocalWellStructure(const int reportStepIdx,
483 const bool enableWellPIScaling);
484
488 void initializeGroupStructure(const int reportStepIdx);
489
490 // called at the end of a time step
491 void timeStepSucceeded(const double simulationTime, const double dt);
492
493 // called at the end of a report step
494 void endReportStep();
495
496 // setting the well_solutions_ based on well_state.
497 void updatePrimaryVariables(DeferredLogger& deferred_logger);
498
500
501 void computePotentials(const std::size_t widx,
502 const WellState<Scalar, IndexTraits>& well_state_copy,
503 std::string& exc_msg,
504 ExceptionType::ExcEnum& exc_type,
505 DeferredLogger& deferred_logger) override;
506
507 const std::vector<Scalar>& wellPerfEfficiencyFactors() const;
508
509 void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger& deferred_logger) override;
510 void calculateProductivityIndexValues(DeferredLogger& deferred_logger) override;
512 DeferredLogger& deferred_logger);
513
514 // The number of conservation quantities.
515 int numConservationQuantities() const;
516
517 int reportStepIndex() const;
518
519 void assembleWellEq(const double dt, DeferredLogger& deferred_logger);
520
521 void prepareWellsBeforeAssembling(const double dt, DeferredLogger& deferred_logger);
522
523 // TODO: finding a better naming
524 void assembleWellEqWithoutIteration(const double dt, DeferredLogger& deferred_logger);
525
527
528 void extractLegacyDepth_();
529
531 void updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const;
532
533 void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger& deferred_logger);
534
535 void calcResvCoeff(const int fipnum,
536 const int pvtreg,
537 const std::vector<Scalar>& production_rates,
538 std::vector<Scalar>& resv_coeff) override;
539
540 void calcInjResvCoeff(const int fipnum,
541 const int pvtreg,
542 std::vector<Scalar>& resv_coeff) override;
543
545
546 private:
548 BlackoilWellModelNldd<TypeTag>* nldd_ = nullptr;
549
550 // These members are used to avoid reallocation in specific functions
551 // instead of using local variables.
552 // Their state is not relevant between function calls, so they can
553 // (and must) be mutable, as the functions using them are const.
554 mutable BVector x_local_;
555
556 void assignWellTracerRates(data::Wells& wsrpt) const;
557 };
558
559} // namespace Opm
560
562
563#endif // OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
Contains the high level supplements required to extend the black oil model by MICP.
Definition: blackoilmicpmodules.hh:54
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
Class for handling the gaslift in the blackoil well model.
Definition: BlackoilWellModelGasLift.hpp:96
Class for handling the blackoil well model.
Definition: BlackoilWellModelGeneric.hpp:95
BlackoilWellModelWBP< GetPropType< TypeTag, Properties::Scalar >, GetPropType< TypeTag, Properties::FluidSystem >::IndexTraitsType > wbp_
Definition: BlackoilWellModelGeneric.hpp:517
const WellState< GetPropType< TypeTag, Properties::Scalar >, GetPropType< TypeTag, Properties::FluidSystem >::IndexTraitsType > & wellState() const
Definition: BlackoilWellModelGeneric.hpp:152
void assignMassGasRate(data::Wells &wsrpt, const GetPropType< TypeTag, Properties::Scalar > gasDensity) const
Class for handling the guide rates in the blackoil well model.
Definition: BlackoilWellModelGuideRates.hpp:46
void assignWellGuideRates(data::Wells &wsrpt, const int reportStepIdx) const
Assign well guide rates.
Class for handling the blackoil well model in a NLDD solver.
Definition: BlackoilWellModelNldd.hpp:80
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:102
void setReservoirCouplingSlave(ReservoirCouplingSlave *slave)
Definition: BlackoilWellModel.hpp:397
void initializeGroupStructure(const int reportStepIdx)
Definition: BlackoilWellModel_impl.hpp:290
void endTimeStep()
Definition: BlackoilWellModel.hpp:171
void prepareTimeStep(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2126
std::tuple< bool, bool, Scalar > updateWellControlsAndNetworkIteration(const bool mandatory_network_balance, const bool relax_network_tolerance, const bool optimize_gas_lift, const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1240
bool updateGroupControls(const Group &group, DeferredLogger &deferred_logger, const int reportStepIdx, const int iterationIdx)
Definition: BlackoilWellModel_impl.hpp:1930
WellInterfacePtr createWellPointer(const int wellID, const int report_step) const
Definition: BlackoilWellModel_impl.hpp:1016
GuideRateHandler< Scalar, IndexTraits > guide_rate_handler_
Definition: BlackoilWellModel.hpp:440
const EclipseState & eclState() const
Definition: BlackoilWellModel.hpp:450
WellInterfacePtr getWell(const std::string &well_name) const
Definition: BlackoilWellModel_impl.hpp:2272
void calcInjResvCoeff(const int fipnum, const int pvtreg, std::vector< Scalar > &resv_coeff) override
Definition: BlackoilWellModel_impl.hpp:2312
const std::vector< WellInterfacePtr > & localNonshutWells() const
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:320
void prepareDeserialize(const int report_step)
Definition: BlackoilWellModel.hpp:205
void init()
Definition: BlackoilWellModel_impl.hpp:157
void setNlddAdapter(BlackoilWellModelNldd< TypeTag > *mod)
Definition: BlackoilWellModel.hpp:377
void calcResvCoeff(const int fipnum, const int pvtreg, const std::vector< Scalar > &production_rates, std::vector< Scalar > &resv_coeff) override
Definition: BlackoilWellModel_impl.hpp:2301
const Simulator & simulator() const
Definition: BlackoilWellModel.hpp:374
std::vector< Scalar > depth_
Definition: BlackoilWellModel.hpp:433
void updateAndCommunicate(const int reportStepIdx, const int iterationIdx, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1894
std::size_t global_num_cells_
Definition: BlackoilWellModel.hpp:429
std::map< std::string, Scalar > well_group_thp_calc_
Definition: BlackoilWellModel.hpp:435
void doPreStepNetworkRebalance(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1092
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: BlackoilWellModel.hpp:111
SimulatorReportSingle last_report_
Definition: BlackoilWellModel.hpp:439
static const int solventSaturationIdx
Definition: BlackoilWellModel.hpp:124
Scalar gravity_
Definition: BlackoilWellModel.hpp:432
void initWellContainer(const int reportStepIdx) override
Definition: BlackoilWellModel_impl.hpp:176
void beginReportStep(const int time_step)
Definition: BlackoilWellModel_impl.hpp:193
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: BlackoilWellModel.hpp:107
void beginIteration()
Definition: BlackoilWellModel.hpp:161
bool isReservoirCouplingMaster() const
Definition: BlackoilWellModel.hpp:387
Dune::FieldVector< Scalar, numEq > VectorBlockType
Definition: BlackoilWellModel.hpp:132
void initFromRestartFile(const RestartValue &restartValues)
Definition: BlackoilWellModel.hpp:195
const ModelParameters param_
Definition: BlackoilWellModel.hpp:428
void beginEpisode()
Definition: BlackoilWellModel.hpp:153
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: BlackoilWellModel.hpp:108
GetPropType< TypeTag, Properties::Grid > Grid
Definition: BlackoilWellModel.hpp:105
GetPropType< TypeTag, Properties::EquilGrid > EquilGrid
Definition: BlackoilWellModel.hpp:106
int numConservationQuantities() const
Definition: BlackoilWellModel_impl.hpp:2244
void assembleWellEq(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1486
const std::map< std::string, int > & well_domain() const
Definition: BlackoilWellModel.hpp:333
void setReservoirCouplingMaster(ReservoirCouplingMaster *master)
Definition: BlackoilWellModel.hpp:393
bool updateWellControls(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1756
void endIteration()
Definition: BlackoilWellModel.hpp:168
static constexpr bool has_polymer_
Definition: BlackoilWellModel.hpp:126
int reportStepIndex() const
Definition: BlackoilWellModel_impl.hpp:2289
void calculateProductivityIndexValues(DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2068
void extractLegacyDepth_()
Definition: BlackoilWellModel_impl.hpp:2260
void extractLegacyCellPvtRegionIndex_()
Definition: BlackoilWellModel_impl.hpp:2228
void recoverWellSolutionAndUpdateWellStateDomain(const BVector &x, const int domainIdx)
Definition: BlackoilWellModel_impl.hpp:1677
void updateAverageFormationFactor()
Definition: BlackoilWellModel_impl.hpp:2171
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: BlackoilWellModel.hpp:110
bool isReservoirCouplingSlave() const
Definition: BlackoilWellModel.hpp:390
void initializeWellState(const int timeStepIdx)
Definition: BlackoilWellModel_impl.hpp:750
void assemble(const int iterationIdx, const double dt)
Definition: BlackoilWellModel_impl.hpp:1122
const Grid & grid() const
Definition: BlackoilWellModel.hpp:371
void computeWellTemperature()
Definition: BlackoilWellModel_impl.hpp:2323
static const int numEq
Definition: BlackoilWellModel.hpp:123
void addWellPressureEquations(PressureMatrix &jacobian, const BVector &weights, const bool use_well_weights) const
Definition: BlackoilWellModel_impl.hpp:1577
const SimulatorReportSingle & lastReport() const
Definition: BlackoilWellModel_impl.hpp:624
bool updateWellControlsAndNetwork(const bool mandatory_network_balance, const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1187
void endEpisode()
Definition: BlackoilWellModel.hpp:177
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: BlackoilWellModel.hpp:114
void addWellContributions(SparseMatrixAdapter &jacobian) const
Definition: BlackoilWellModel_impl.hpp:1567
const EquilGrid & equilGrid() const
Definition: BlackoilWellModel.hpp:447
std::tuple< bool, Scalar > updateNetworks(const bool mandatory_network_balance, DeferredLogger &deferred_logger, const bool relax_network_tolerance=false)
Definition: BlackoilWellModel_impl.hpp:1833
WellInterfacePtr createWellForWellTest(const std::string &well_name, const int report_step, DeferredLogger &deferred_logger) const
Definition: BlackoilWellModel_impl.hpp:1067
GetPropType< TypeTag, Properties::Indices > Indices
Definition: BlackoilWellModel.hpp:109
int compressedIndexForInteriorLGR(const std::string &lgr_tag, const Connection &conn) const override
Definition: BlackoilWellModel.hpp:357
void assembleWellEqWithoutIteration(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1510
std::vector< WellInterfacePtr > well_container_
Definition: BlackoilWellModel.hpp:406
static constexpr std::size_t pressureVarIndex
Definition: BlackoilWellModel.hpp:121
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: BlackoilWellModel.hpp:299
bool empty() const
Definition: BlackoilWellModel.hpp:344
void computeTotalRatesForDof(RateVector &rate, unsigned globalIdx) const
Definition: BlackoilWellModel_impl.hpp:711
bool addMatrixContributions() const
Definition: BlackoilWellModel.hpp:346
void beginTimeStep()
Definition: BlackoilWellModel_impl.hpp:330
std::vector< Scalar > B_avg_
Definition: BlackoilWellModel.hpp:445
const std::vector< Scalar > & wellPerfEfficiencyFactors() const
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: BlackoilWellModel.hpp:112
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Definition: BlackoilWellModel_impl.hpp:1741
std::map< std::string, std::unique_ptr< AverageRegionalPressureType > > regionalAveragePressureCalculator_
Definition: BlackoilWellModel.hpp:437
void updatePrimaryVariables(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2219
bool computeWellGroupThp(const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1303
GetPropType< TypeTag, Properties::GlobalEqVector > GlobalEqVector
Definition: BlackoilWellModel.hpp:113
void initializeLocalWellStructure(const int reportStepIdx, const bool enableWellPIScaling)
Definition: BlackoilWellModel_impl.hpp:245
static constexpr bool has_solvent_
Definition: BlackoilWellModel.hpp:125
Dune::BlockVector< VectorBlockType > BVector
Definition: BlackoilWellModel.hpp:133
static constexpr bool has_micp_
Definition: BlackoilWellModel.hpp:128
std::unique_ptr< RateConverterType > rateConverter_
Definition: BlackoilWellModel.hpp:436
BlackoilWellModel(Simulator &simulator)
Definition: BlackoilWellModel_impl.hpp:72
const SparseTable< int > & well_local_cells() const
Definition: BlackoilWellModel.hpp:325
void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:535
std::vector< bool > is_cell_perforated_
Definition: BlackoilWellModel.hpp:408
int numStrictIterations() const
Definition: BlackoilWellModel.hpp:349
ConvergenceReport getWellConvergence(const std::vector< Scalar > &B_avg, const bool checkWellGroupControls=false) const
Definition: BlackoilWellModel_impl.hpp:1690
typename FluidSystem::IndexTraitsType IndexTraits
Definition: BlackoilWellModel.hpp:119
data::WellBlockAveragePressures wellBlockAveragePressures() const
Definition: BlackoilWellModel.hpp:247
auto begin() const
Definition: BlackoilWellModel.hpp:342
std::size_t local_num_cells_
Definition: BlackoilWellModel.hpp:431
auto end() const
Definition: BlackoilWellModel.hpp:343
bool alternative_well_rate_init_
Definition: BlackoilWellModel.hpp:434
void timeStepSucceeded(const double simulationTime, const double dt)
Definition: BlackoilWellModel_impl.hpp:634
std::unique_ptr< WellType > createTypedWellPointer(const int wellID, const int time_step) const
Definition: BlackoilWellModel_impl.hpp:1036
ReservoirCouplingMaster & reservoirCouplingMaster()
Definition: BlackoilWellModel.hpp:381
Simulator & simulator_
Definition: BlackoilWellModel.hpp:403
std::shared_ptr< WellInterface< TypeTag > > WellInterfacePtr
Definition: BlackoilWellModel.hpp:192
void createWellContainer(const int report_step) override
Definition: BlackoilWellModel_impl.hpp:793
data::Wells wellData() const
Definition: BlackoilWellModel.hpp:212
ReservoirCouplingSlave & reservoirCouplingSlave()
Definition: BlackoilWellModel.hpp:384
void addWellPressureEquationsDomain(PressureMatrix &jacobian, const BVector &weights, const bool use_well_weights, const int domainIndex) const
Definition: BlackoilWellModel.hpp:305
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1975
void addReservoirSourceTerms(GlobalEqVector &residual, const std::vector< typename SparseMatrixAdapter::MatrixBlock * > &diagMatAddress) const
Definition: BlackoilWellModel_impl.hpp:1599
int compressedIndexForInterior(int cartesian_cell_idx) const override
get compressed index for interior cells (-1, otherwise
Definition: BlackoilWellModel.hpp:352
static constexpr bool has_energy_
Definition: BlackoilWellModel.hpp:127
void recoverWellSolutionAndUpdateWellState(const BVector &x)
Definition: BlackoilWellModel_impl.hpp:1652
void addWellPressureEquationsStruct(PressureMatrix &jacobian) const
Definition: BlackoilWellModel_impl.hpp:1629
void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2082
void computePotentials(const std::size_t widx, const WellState< Scalar, IndexTraits > &well_state_copy, std::string &exc_msg, ExceptionType::ExcEnum &exc_type, DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2032
void endReportStep()
Definition: BlackoilWellModel_impl.hpp:607
void prepareWellsBeforeAssembling(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1498
data::WellBlockAveragePressures computeWellBlockAveragePressures(const Scalar gravity) const
Definition: ConvergenceReport.hpp:38
Definition: DeferredLogger.hpp:57
Handles computation and reporting of guide rates for wells and groups.
Definition: GuideRateHandler.hpp:53
Definition: RateConverter.hpp:70
Definition: RegionAverageCalculator.hpp:60
Definition: ReservoirCouplingMaster.hpp:35
Definition: ReservoirCouplingSlave.hpp:35
Definition: BlackoilWellModel.hpp:91
Definition: WellConnectionAuxiliaryModule.hpp:39
Definition: WellContributions.hpp:51
Definition: WellGroupHelpers.hpp:51
Definition: WellInterface.hpp:76
Definition: WellState.hpp:66
data::Wells report(const int *globalCellIdxMap, const std::function< bool(const int)> &wasDynamicallyClosed) const
ExcEnum
Definition: DeferredLogger.hpp:45
Definition: blackoilboundaryratevector.hh:39
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:179
bool matrix_add_well_contributions_
Whether to add influences of wells between cells to the matrix and preconditioner matrix.
Definition: BlackoilModelParameters.hpp:285
bool use_multisegment_well_
Definition: BlackoilModelParameters.hpp:279
int strict_outer_iter_wells_
Newton iteration where wells are stricly convergent.
Definition: BlackoilModelParameters.hpp:235
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34