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,
101 Properties::Scalar>>
102 {
103 public:
104 // --------- Types ---------
116
118
120
121 static const int numEq = Indices::numEq;
122 static const int solventSaturationIdx = Indices::solventSaturationIdx;
123 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
124 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
125 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
126 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
127
128 // TODO: where we should put these types, WellInterface or Well Model?
129 // or there is some other strategy, like TypeTag
130 using VectorBlockType = Dune::FieldVector<Scalar, numEq>;
131 using BVector = Dune::BlockVector<VectorBlockType>;
132
135
136 // For the conversion between the surface volume rate and reservoir voidage rate
137 using RateConverterType = RateConverter::
138 SurfaceToReservoirVoidage<FluidSystem, std::vector<int> >;
139
140 // For computing average pressured used by gpmaint
141 using AverageRegionalPressureType = RegionAverageCalculator::
142 AverageRegionalPressure<FluidSystem, std::vector<int> >;
143
145
146 void init();
147 void initWellContainer(const int reportStepIdx) override;
148
150 {
151 OPM_TIMEBLOCK(beginEpsiode);
152 beginReportStep(simulator_.episodeIndex());
153 }
154
155 void beginTimeStep();
156
158 {
159 OPM_TIMEBLOCK(beginIteration);
160 assemble(simulator_.model().newtonMethod().numIterations(),
161 simulator_.timeStepSize());
162 }
163
165 { }
166
168 {
169 OPM_TIMEBLOCK(endTimeStep);
170 timeStepSucceeded(simulator_.time(), simulator_.timeStepSize());
171 }
172
174 {
176 }
177
179 unsigned globalIdx) const;
180
181 template <class Context>
183 const Context& context,
184 unsigned spaceIdx,
185 unsigned timeIdx) const;
186
187
188 using WellInterfacePtr = std::shared_ptr<WellInterface<TypeTag> >;
189
191 void initFromRestartFile(const RestartValue& restartValues)
192 {
193 initFromRestartFile(restartValues,
194 this->simulator_.vanguard().transferWTestState(),
195 grid().size(0),
197 this->simulator_.vanguard().enableDistributedWells());
198 }
199
201 void prepareDeserialize(const int report_step)
202 {
203 prepareDeserialize(report_step, grid().size(0),
205 this->simulator_.vanguard().enableDistributedWells());
206 }
207
208 data::Wells wellData() const
209 {
210 auto wsrpt = this->wellState()
211 .report(this->simulator_.vanguard().globalCell().data(),
212 [this](const int well_index)
213 {
214 return this->wasDynamicallyShutThisTimeStep(well_index);
215 });
216
218 .assignWellGuideRates(wsrpt, this->reportStepIndex());
219
220 this->assignWellTracerRates(wsrpt);
221
222 if (const auto& rspec = eclState().runspec();
223 rspec.co2Storage() || rspec.h2Storage())
224 {
225 // The gas reference density (surface condition) is the
226 // same for all PVT regions in CO2STORE/H2STORE runs so,
227 // for simplicity, we use region zero (0) here.
228
229 this->assignMassGasRate(wsrpt, FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, 0));
230 }
231
232 this->assignWellTargets(wsrpt);
233
234 this->assignDynamicWellStatus(wsrpt, this->reportStepIndex());
235
236 // Assigning (a subset of the) property values in shut
237 // connections should be the last step of wellData().
238 this->assignShutConnections(wsrpt, this->reportStepIndex());
239
240 return wsrpt;
241 }
242
243 data::WellBlockAveragePressures wellBlockAveragePressures() const
244 {
246 }
247
248#if COMPILE_GPU_BRIDGE
249 // accumulate the contributions of all Wells in the WellContributions object
250 void getWellContributions(WellContributions<Scalar>& x) const;
251#endif
252
253 // Check if well equations is converged.
254 ConvergenceReport getWellConvergence(const std::vector<Scalar>& B_avg, const bool checkWellGroupControls = false) const;
255
256 const SimulatorReportSingle& lastReport() const;
257
258 void addWellContributions(SparseMatrixAdapter& jacobian) const;
259
260 // add source from wells to the reservoir matrix
262 const std::vector<typename SparseMatrixAdapter::MatrixBlock*>& diagMatAddress) const;
263
264 // called at the beginning of a report step
265 void beginReportStep(const int time_step);
266
267 // it should be able to go to prepareTimeStep(), however, the updateWellControls()
268 // makes it a little more difficult. unless we introduce if (iterationIdx != 0) to avoid doing the above function
269 // twice at the beginning of the time step
272 void calculateExplicitQuantities(DeferredLogger& deferred_logger) const;
273 // some preparation work, mostly related to group control and RESV,
274 // at the beginning of each time step (Not report step)
275 void prepareTimeStep(DeferredLogger& deferred_logger);
276
277 bool
278 updateWellControls(DeferredLogger& deferred_logger);
279
280 std::tuple<bool, Scalar>
281 updateNetworks(const bool mandatory_network_balance, DeferredLogger& deferred_logger, const bool relax_network_tolerance = false);
282
283
284 void updateAndCommunicate(const int reportStepIdx,
285 const int iterationIdx,
286 DeferredLogger& deferred_logger);
287
288 bool updateGroupControls(const Group& group,
289 DeferredLogger& deferred_logger,
290 const int reportStepIdx,
291 const int iterationIdx);
292
293 WellInterfacePtr getWell(const std::string& well_name) const;
294
295 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<Scalar, 1, 1>>;
296
298 const BVector& weights,
299 const bool use_well_weights) const;
302 const BVector& weights,
303 const bool use_well_weights,
304 const int domainIndex) const
305 {
306 if (!nldd_) {
307 OPM_THROW(std::logic_error, "Attempt to access NLDD data without a NLDD solver");
308 }
309 return nldd_->addWellPressureEquations(jacobian,
310 weights,
311 use_well_weights,
312 domainIndex);
313 }
314
316 const std::vector<WellInterfacePtr>& localNonshutWells() const
317 {
318 return well_container_;
319 }
320
322 {
323 if (!nldd_) {
324 OPM_THROW(std::logic_error, "Attempt to access NLDD data without a NLDD solver");
325 }
326 return nldd_->well_local_cells();
327 }
328
329 const std::map<std::string, int>& well_domain() const
330 {
331 if (!nldd_) {
332 OPM_THROW(std::logic_error, "Attempt to access NLDD data without a NLDD solver");
333 }
334
335 return nldd_->well_domain();
336 }
337
338 auto begin() const { return well_container_.begin(); }
339 auto end() const { return well_container_.end(); }
340 bool empty() const { return well_container_.empty(); }
341
344
347
348 int compressedIndexForInterior(int cartesian_cell_idx) const override
349 {
350 return simulator_.vanguard().compressedIndexForInterior(cartesian_cell_idx);
351 }
352
353 int compressedIndexForInteriorLGR(const std::string& lgr_tag, const Connection& conn) const override
354 {
355 return simulator_.vanguard().compressedIndexForInteriorLGR(lgr_tag, conn);
356 }
357
358 // using the solution x to recover the solution xw for wells and applying
359 // xw to update Well State
361
362 // using the solution x to recover the solution xw for wells and applying
363 // xw to update Well State
365 const int domainIdx);
366
367 const Grid& grid() const
368 { return simulator_.vanguard().grid(); }
369
370 const Simulator& simulator() const
371 { return simulator_; }
372
374 { nldd_ = mod; }
375
376#ifdef RESERVOIR_COUPLING_ENABLED
378 return *(this->simulator_.reservoirCouplingMaster());
379 }
381 return *(this->simulator_.reservoirCouplingSlave());
382 }
384 return this->simulator_.reservoirCouplingMaster() != nullptr;
385 }
387 return this->simulator_.reservoirCouplingSlave() != nullptr;
388 }
390 {
391 this->guide_rate_handler_.setReservoirCouplingMaster(master);
392 }
394 {
395 this->guide_rate_handler_.setReservoirCouplingSlave(slave);
396 }
397 #endif
398 protected:
400
401 // a vector of all the wells.
402 std::vector<WellInterfacePtr> well_container_{};
403
404 std::vector<bool> is_cell_perforated_{};
405
406 void initializeWellState(const int timeStepIdx);
407
408 // create the well container
409 void createWellContainer(const int report_step) override;
410
412 createWellPointer(const int wellID,
413 const int report_step) const;
414
415 template <typename WellType>
416 std::unique_ptr<WellType>
417 createTypedWellPointer(const int wellID,
418 const int time_step) const;
419
420 WellInterfacePtr createWellForWellTest(const std::string& well_name,
421 const int report_step,
422 DeferredLogger& deferred_logger) const;
423
425 std::size_t global_num_cells_{};
426 // the number of the cells in the local grid
427 std::size_t local_num_cells_{};
429 std::vector<Scalar> depth_{};
431 std::map<std::string, Scalar> well_group_thp_calc_;
432 std::unique_ptr<RateConverterType> rateConverter_{};
433 std::map<std::string, std::unique_ptr<AverageRegionalPressureType>> regionalAveragePressureCalculator_{};
434
437
438 // Pre-step network solve at static reservoir conditions (group and well states might be updated)
439 void doPreStepNetworkRebalance(DeferredLogger& deferred_logger);
440
441 std::vector<Scalar> B_avg_{};
442
443 const EquilGrid& equilGrid() const
444 { return simulator_.vanguard().equilGrid(); }
445
446 const EclipseState& eclState() const
447 { return simulator_.vanguard().eclState(); }
448
449 // compute the well fluxes and assemble them in to the reservoir equations as source terms
450 // and in the well equations.
451 void assemble(const int iterationIdx,
452 const double dt);
453
454 // well controls and network pressures affect each other and are solved in an iterative manner.
455 // the function handles one iteration of updating well controls and network pressures.
456 // it is possible to decouple the update of well controls and network pressures further.
457 // the returned two booleans are {continue_due_to_network, well_group_control_changed}, respectively
458 std::tuple<bool, bool, Scalar> updateWellControlsAndNetworkIteration(const bool mandatory_network_balance,
459 const bool relax_network_tolerance,
460 const bool optimize_gas_lift,
461 const double dt,
462 DeferredLogger& local_deferredLogger);
463
464 bool updateWellControlsAndNetwork(const bool mandatory_network_balance,
465 const double dt,
466 DeferredLogger& local_deferredLogger);
467
468 bool computeWellGroupThp(const double dt, DeferredLogger& local_deferredLogger);
469
478 void initializeLocalWellStructure(const int reportStepIdx,
479 const bool enableWellPIScaling);
480
484 void initializeGroupStructure(const int reportStepIdx);
485
486 // called at the end of a time step
487 void timeStepSucceeded(const double simulationTime, const double dt);
488
489 // called at the end of a report step
490 void endReportStep();
491
492 // setting the well_solutions_ based on well_state.
493 void updatePrimaryVariables(DeferredLogger& deferred_logger);
494
496
497 void computePotentials(const std::size_t widx,
498 const WellState<Scalar>& well_state_copy,
499 std::string& exc_msg,
500 ExceptionType::ExcEnum& exc_type,
501 DeferredLogger& deferred_logger) override;
502
503 const std::vector<Scalar>& wellPerfEfficiencyFactors() const;
504
505 void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger& deferred_logger) override;
506 void calculateProductivityIndexValues(DeferredLogger& deferred_logger) override;
508 DeferredLogger& deferred_logger);
509
510 // The number of components in the model.
511 int numComponents() const;
512
513 int reportStepIndex() const;
514
515 void assembleWellEq(const double dt, DeferredLogger& deferred_logger);
516
517 void prepareWellsBeforeAssembling(const double dt, DeferredLogger& deferred_logger);
518
519 // TODO: finding a better naming
520 void assembleWellEqWithoutIteration(const double dt, DeferredLogger& deferred_logger);
521
523
524 void extractLegacyDepth_();
525
527 void updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const;
528
529 void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger& deferred_logger);
530
531 void calcResvCoeff(const int fipnum,
532 const int pvtreg,
533 const std::vector<Scalar>& production_rates,
534 std::vector<Scalar>& resv_coeff) override;
535
536 void calcInjResvCoeff(const int fipnum,
537 const int pvtreg,
538 std::vector<Scalar>& resv_coeff) override;
539
541
542 private:
544
546 BlackoilWellModelNldd<TypeTag>* nldd_ = nullptr;
547
548 // These members are used to avoid reallocation in specific functions
549 // instead of using local variables.
550 // Their state is not relevant between function calls, so they can
551 // (and must) be mutable, as the functions using them are const.
552 mutable BVector x_local_;
553
554 void assignWellTracerRates(data::Wells& wsrpt) const;
555 };
556
557} // namespace Opm
558
560
561#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:94
Class for handling the blackoil well model.
Definition: BlackoilWellModelGeneric.hpp:94
BlackoilWellModelWBP< GetPropType< TypeTag, Properties::Scalar > > wbp_
Definition: BlackoilWellModelGeneric.hpp:516
void assignDynamicWellStatus(data::Wells &wsrpt, const int reportStepIdx) const
void assignShutConnections(data::Wells &wsrpt, const int reportStepIndex) const
const WellState< GetPropType< TypeTag, Properties::Scalar > > & wellState() const
Definition: BlackoilWellModelGeneric.hpp:150
WellTestState & wellTestState()
Definition: BlackoilWellModelGeneric.hpp:173
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:79
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:102
void setReservoirCouplingSlave(ReservoirCouplingSlave *slave)
Definition: BlackoilWellModel.hpp:393
void initializeGroupStructure(const int reportStepIdx)
Definition: BlackoilWellModel_impl.hpp:298
void endTimeStep()
Definition: BlackoilWellModel.hpp:167
void prepareTimeStep(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2129
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:1251
bool updateGroupControls(const Group &group, DeferredLogger &deferred_logger, const int reportStepIdx, const int iterationIdx)
Definition: BlackoilWellModel_impl.hpp:1935
WellInterfacePtr createWellPointer(const int wellID, const int report_step) const
Definition: BlackoilWellModel_impl.hpp:1027
const EclipseState & eclState() const
Definition: BlackoilWellModel.hpp:446
WellInterfacePtr getWell(const std::string &well_name) const
Definition: BlackoilWellModel_impl.hpp:2276
void calcInjResvCoeff(const int fipnum, const int pvtreg, std::vector< Scalar > &resv_coeff) override
Definition: BlackoilWellModel_impl.hpp:2316
const std::vector< WellInterfacePtr > & localNonshutWells() const
Get list of local nonshut wells.
Definition: BlackoilWellModel.hpp:316
void prepareDeserialize(const int report_step)
Definition: BlackoilWellModel.hpp:201
void init()
Definition: BlackoilWellModel_impl.hpp:165
void setNlddAdapter(BlackoilWellModelNldd< TypeTag > *mod)
Definition: BlackoilWellModel.hpp:373
void calcResvCoeff(const int fipnum, const int pvtreg, const std::vector< Scalar > &production_rates, std::vector< Scalar > &resv_coeff) override
Definition: BlackoilWellModel_impl.hpp:2305
const Simulator & simulator() const
Definition: BlackoilWellModel.hpp:370
std::vector< Scalar > depth_
Definition: BlackoilWellModel.hpp:429
void updateAndCommunicate(const int reportStepIdx, const int iterationIdx, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1899
std::size_t global_num_cells_
Definition: BlackoilWellModel.hpp:425
std::map< std::string, Scalar > well_group_thp_calc_
Definition: BlackoilWellModel.hpp:431
void doPreStepNetworkRebalance(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1103
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: BlackoilWellModel.hpp:111
SimulatorReportSingle last_report_
Definition: BlackoilWellModel.hpp:435
static const int solventSaturationIdx
Definition: BlackoilWellModel.hpp:122
Scalar gravity_
Definition: BlackoilWellModel.hpp:428
void initWellContainer(const int reportStepIdx) override
Definition: BlackoilWellModel_impl.hpp:184
void beginReportStep(const int time_step)
Definition: BlackoilWellModel_impl.hpp:201
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: BlackoilWellModel.hpp:107
void beginIteration()
Definition: BlackoilWellModel.hpp:157
bool isReservoirCouplingMaster() const
Definition: BlackoilWellModel.hpp:383
Dune::FieldVector< Scalar, numEq > VectorBlockType
Definition: BlackoilWellModel.hpp:130
void initFromRestartFile(const RestartValue &restartValues)
Definition: BlackoilWellModel.hpp:191
const ModelParameters param_
Definition: BlackoilWellModel.hpp:424
void beginEpisode()
Definition: BlackoilWellModel.hpp:149
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
void assembleWellEq(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1491
const std::map< std::string, int > & well_domain() const
Definition: BlackoilWellModel.hpp:329
void setReservoirCouplingMaster(ReservoirCouplingMaster *master)
Definition: BlackoilWellModel.hpp:389
bool updateWellControls(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1761
void endIteration()
Definition: BlackoilWellModel.hpp:164
static constexpr bool has_polymer_
Definition: BlackoilWellModel.hpp:124
int reportStepIndex() const
Definition: BlackoilWellModel_impl.hpp:2293
void calculateProductivityIndexValues(DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2071
void extractLegacyDepth_()
Definition: BlackoilWellModel_impl.hpp:2264
void extractLegacyCellPvtRegionIndex_()
Definition: BlackoilWellModel_impl.hpp:2231
void recoverWellSolutionAndUpdateWellStateDomain(const BVector &x, const int domainIdx)
Definition: BlackoilWellModel_impl.hpp:1682
void updateAverageFormationFactor()
Definition: BlackoilWellModel_impl.hpp:2174
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: BlackoilWellModel.hpp:110
bool isReservoirCouplingSlave() const
Definition: BlackoilWellModel.hpp:386
void initializeWellState(const int timeStepIdx)
Definition: BlackoilWellModel_impl.hpp:761
void assemble(const int iterationIdx, const double dt)
Definition: BlackoilWellModel_impl.hpp:1133
const Grid & grid() const
Definition: BlackoilWellModel.hpp:367
void computeWellTemperature()
Definition: BlackoilWellModel_impl.hpp:2327
static const int numEq
Definition: BlackoilWellModel.hpp:121
void addWellPressureEquations(PressureMatrix &jacobian, const BVector &weights, const bool use_well_weights) const
Definition: BlackoilWellModel_impl.hpp:1582
const SimulatorReportSingle & lastReport() const
Definition: BlackoilWellModel_impl.hpp:635
bool updateWellControlsAndNetwork(const bool mandatory_network_balance, const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1198
void endEpisode()
Definition: BlackoilWellModel.hpp:173
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: BlackoilWellModel.hpp:114
void addWellContributions(SparseMatrixAdapter &jacobian) const
Definition: BlackoilWellModel_impl.hpp:1572
const EquilGrid & equilGrid() const
Definition: BlackoilWellModel.hpp:443
std::tuple< bool, Scalar > updateNetworks(const bool mandatory_network_balance, DeferredLogger &deferred_logger, const bool relax_network_tolerance=false)
Definition: BlackoilWellModel_impl.hpp:1838
WellInterfacePtr createWellForWellTest(const std::string &well_name, const int report_step, DeferredLogger &deferred_logger) const
Definition: BlackoilWellModel_impl.hpp:1078
GetPropType< TypeTag, Properties::Indices > Indices
Definition: BlackoilWellModel.hpp:109
void computePotentials(const std::size_t widx, const WellState< Scalar > &well_state_copy, std::string &exc_msg, ExceptionType::ExcEnum &exc_type, DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2035
int compressedIndexForInteriorLGR(const std::string &lgr_tag, const Connection &conn) const override
Definition: BlackoilWellModel.hpp:353
void assembleWellEqWithoutIteration(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1515
std::vector< WellInterfacePtr > well_container_
Definition: BlackoilWellModel.hpp:402
GuideRateHandler< Scalar > guide_rate_handler_
Definition: BlackoilWellModel.hpp:436
static constexpr std::size_t pressureVarIndex
Definition: BlackoilWellModel.hpp:119
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: BlackoilWellModel.hpp:295
bool empty() const
Definition: BlackoilWellModel.hpp:340
void computeTotalRatesForDof(RateVector &rate, unsigned globalIdx) const
Definition: BlackoilWellModel_impl.hpp:722
bool addMatrixContributions() const
Definition: BlackoilWellModel.hpp:342
void beginTimeStep()
Definition: BlackoilWellModel_impl.hpp:340
std::vector< Scalar > B_avg_
Definition: BlackoilWellModel.hpp:441
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:1746
std::map< std::string, std::unique_ptr< AverageRegionalPressureType > > regionalAveragePressureCalculator_
Definition: BlackoilWellModel.hpp:433
void updatePrimaryVariables(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2222
bool computeWellGroupThp(const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1308
GetPropType< TypeTag, Properties::GlobalEqVector > GlobalEqVector
Definition: BlackoilWellModel.hpp:113
void initializeLocalWellStructure(const int reportStepIdx, const bool enableWellPIScaling)
Definition: BlackoilWellModel_impl.hpp:253
static constexpr bool has_solvent_
Definition: BlackoilWellModel.hpp:123
Dune::BlockVector< VectorBlockType > BVector
Definition: BlackoilWellModel.hpp:131
static constexpr bool has_micp_
Definition: BlackoilWellModel.hpp:126
std::unique_ptr< RateConverterType > rateConverter_
Definition: BlackoilWellModel.hpp:432
BlackoilWellModel(Simulator &simulator)
Definition: BlackoilWellModel_impl.hpp:157
const SparseTable< int > & well_local_cells() const
Definition: BlackoilWellModel.hpp:321
void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:545
std::vector< bool > is_cell_perforated_
Definition: BlackoilWellModel.hpp:404
int numStrictIterations() const
Definition: BlackoilWellModel.hpp:345
ConvergenceReport getWellConvergence(const std::vector< Scalar > &B_avg, const bool checkWellGroupControls=false) const
Definition: BlackoilWellModel_impl.hpp:1695
data::WellBlockAveragePressures wellBlockAveragePressures() const
Definition: BlackoilWellModel.hpp:243
auto begin() const
Definition: BlackoilWellModel.hpp:338
std::size_t local_num_cells_
Definition: BlackoilWellModel.hpp:427
auto end() const
Definition: BlackoilWellModel.hpp:339
bool alternative_well_rate_init_
Definition: BlackoilWellModel.hpp:430
void timeStepSucceeded(const double simulationTime, const double dt)
Definition: BlackoilWellModel_impl.hpp:645
std::unique_ptr< WellType > createTypedWellPointer(const int wellID, const int time_step) const
Definition: BlackoilWellModel_impl.hpp:1047
ReservoirCouplingMaster & reservoirCouplingMaster()
Definition: BlackoilWellModel.hpp:377
Simulator & simulator_
Definition: BlackoilWellModel.hpp:399
std::shared_ptr< WellInterface< TypeTag > > WellInterfacePtr
Definition: BlackoilWellModel.hpp:188
void createWellContainer(const int report_step) override
Definition: BlackoilWellModel_impl.hpp:804
data::Wells wellData() const
Definition: BlackoilWellModel.hpp:208
ReservoirCouplingSlave & reservoirCouplingSlave()
Definition: BlackoilWellModel.hpp:380
void addWellPressureEquationsDomain(PressureMatrix &jacobian, const BVector &weights, const bool use_well_weights, const int domainIndex) const
Definition: BlackoilWellModel.hpp:301
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:1978
void addReservoirSourceTerms(GlobalEqVector &residual, const std::vector< typename SparseMatrixAdapter::MatrixBlock * > &diagMatAddress) const
Definition: BlackoilWellModel_impl.hpp:1604
int compressedIndexForInterior(int cartesian_cell_idx) const override
get compressed index for interior cells (-1, otherwise
Definition: BlackoilWellModel.hpp:348
static constexpr bool has_energy_
Definition: BlackoilWellModel.hpp:125
void recoverWellSolutionAndUpdateWellState(const BVector &x)
Definition: BlackoilWellModel_impl.hpp:1657
void addWellPressureEquationsStruct(PressureMatrix &jacobian) const
Definition: BlackoilWellModel_impl.hpp:1634
void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2085
void endReportStep()
Definition: BlackoilWellModel_impl.hpp:618
void prepareWellsBeforeAssembling(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1503
int numComponents() const
Definition: BlackoilWellModel_impl.hpp:2247
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:54
Definition: RateConverter.hpp:71
Definition: RegionAverageCalculator.hpp:61
Definition: ReservoirCouplingMaster.hpp:35
Definition: ReservoirCouplingSlave.hpp:35
Definition: BlackoilWellModel.hpp:91
Definition: WellConnectionAuxiliaryModule.hpp:39
Definition: WellContributions.hpp:51
Definition: WellInterface.hpp:77
Definition: WellState.hpp:65
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:177
bool matrix_add_well_contributions_
Whether to add influences of wells between cells to the matrix and preconditioner matrix.
Definition: BlackoilModelParameters.hpp:283
bool use_multisegment_well_
Definition: BlackoilModelParameters.hpp:277
int strict_outer_iter_wells_
Newton iteration where wells are stricly convergent.
Definition: BlackoilModelParameters.hpp:233
Definition: BlackoilPhases.hpp:46
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34