BlackoilWellModel_impl.hpp
Go to the documentation of this file.
1/*
2 Copyright 2016 - 2019 SINTEF Digital, Mathematics & Cybernetics.
3 Copyright 2016 - 2018 Equinor ASA.
4 Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5 Copyright 2016 - 2018 Norce 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_IMPL_HEADER_INCLUDED
24#define OPM_BLACKOILWELLMODEL_IMPL_HEADER_INCLUDED
25
26// Improve IDE experience
27#ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
28#include <config.h>
30#endif
31
32#include <opm/grid/utility/cartesianToCompressed.hpp>
33#include <opm/common/utility/numeric/RootFinders.hpp>
34
35#include <opm/input/eclipse/Schedule/Network/Balance.hpp>
36#include <opm/input/eclipse/Schedule/Network/ExtNetwork.hpp>
37#include <opm/input/eclipse/Schedule/Well/PAvgDynamicSourceData.hpp>
38#include <opm/input/eclipse/Schedule/Well/WellMatcher.hpp>
39#include <opm/input/eclipse/Schedule/Well/WellTestConfig.hpp>
40
41#include <opm/input/eclipse/Units/UnitSystem.hpp>
42
53
57
58#if COMPILE_GPU_BRIDGE
60#endif
61
62#include <algorithm>
63#include <cassert>
64#include <iomanip>
65#include <utility>
66#include <optional>
67
68#include <fmt/format.h>
69
70namespace Opm {
71 template<typename TypeTag>
73 BlackoilWellModel(Simulator& simulator, const PhaseUsage& phase_usage)
74 : WellConnectionModule(*this, simulator.gridView().comm())
75 , BlackoilWellModelGeneric<Scalar>(simulator.vanguard().schedule(),
76 gaslift_,
77 simulator.vanguard().summaryState(),
78 simulator.vanguard().eclState(),
79 phase_usage,
80 simulator.gridView().comm())
81 , simulator_(simulator)
82 , guide_rate_handler_{
83 *this,
84 simulator.vanguard().schedule(),
85 simulator.vanguard().summaryState(),
86 simulator.vanguard().grid().comm()
87 }
88 , gaslift_(this->terminal_output_, this->phase_usage_)
89 {
90 local_num_cells_ = simulator_.gridView().size(0);
91
92 // Number of cells the global grid view
93 global_num_cells_ = simulator_.vanguard().globalNumCells();
94
95 {
96 auto& parallel_wells = simulator.vanguard().parallelWells();
97
98 this->parallel_well_info_.reserve(parallel_wells.size());
99 for( const auto& name_bool : parallel_wells) {
100 this->parallel_well_info_.emplace_back(name_bool, grid().comm());
101 }
102 }
103
104 this->alternative_well_rate_init_ =
105 Parameters::Get<Parameters::AlternativeWellRateInit>();
106
107 using SourceDataSpan =
108 typename PAvgDynamicSourceData<Scalar>::template SourceDataSpan<Scalar>;
109
110 this->wbp_.initializeSources(
111 [this](const std::size_t globalIndex)
112 { return this->compressedIndexForInterior(globalIndex); },
113 [this](const int localCell, SourceDataSpan sourceTerms)
114 {
115 using Item = typename SourceDataSpan::Item;
116
117 const auto* intQuants = this->simulator_.model()
118 .cachedIntensiveQuantities(localCell, /*timeIndex = */0);
119 const auto& fs = intQuants->fluidState();
120
121 sourceTerms
122 .set(Item::PoreVol, intQuants->porosity().value() *
123 this->simulator_.model().dofTotalVolume(localCell))
124 .set(Item::Depth, this->depth_[localCell]);
125
126 constexpr auto io = FluidSystem::oilPhaseIdx;
127 constexpr auto ig = FluidSystem::gasPhaseIdx;
128 constexpr auto iw = FluidSystem::waterPhaseIdx;
129
130 // Ideally, these would be 'constexpr'.
131 const auto haveOil = FluidSystem::phaseIsActive(io);
132 const auto haveGas = FluidSystem::phaseIsActive(ig);
133 const auto haveWat = FluidSystem::phaseIsActive(iw);
134
135 auto weightedPhaseDensity = [&fs](const auto ip)
136 {
137 return fs.saturation(ip).value() * fs.density(ip).value();
138 };
139
140 if (haveOil) { sourceTerms.set(Item::Pressure, fs.pressure(io).value()); }
141 else if (haveGas) { sourceTerms.set(Item::Pressure, fs.pressure(ig).value()); }
142 else { sourceTerms.set(Item::Pressure, fs.pressure(iw).value()); }
143
144 // Strictly speaking, assumes SUM(s[p]) == 1.
145 auto rho = 0.0;
146 if (haveOil) { rho += weightedPhaseDensity(io); }
147 if (haveGas) { rho += weightedPhaseDensity(ig); }
148 if (haveWat) { rho += weightedPhaseDensity(iw); }
149
150 sourceTerms.set(Item::MixtureDensity, rho);
151 }
152 );
153 }
154
155 template<typename TypeTag>
157 BlackoilWellModel(Simulator& simulator) :
158 BlackoilWellModel(simulator, phaseUsageFromDeck(simulator.vanguard().eclState()))
159 {}
160
161
162 template<typename TypeTag>
163 void
165 init()
166 {
167 extractLegacyCellPvtRegionIndex_();
168 extractLegacyDepth_();
169
170 gravity_ = simulator_.problem().gravity()[2];
171
172 this->initial_step_ = true;
173
174 // add the eWoms auxiliary module for the wells to the list
175 simulator_.model().addAuxiliaryModule(this);
176
177 is_cell_perforated_.resize(local_num_cells_, false);
178 }
179
180
181 template<typename TypeTag>
182 void
184 initWellContainer(const int reportStepIdx)
185 {
186 const uint64_t effective_events_mask = ScheduleEvents::WELL_STATUS_CHANGE
187 + ScheduleEvents::NEW_WELL;
188 const auto& events = this->schedule()[reportStepIdx].wellgroup_events();
189 for (auto& wellPtr : this->well_container_) {
190 const bool well_opened_this_step = this->report_step_starts_ &&
191 events.hasEvent(wellPtr->name(),
192 effective_events_mask);
193 wellPtr->init(&this->phase_usage_, this->depth_, this->gravity_,
194 this->B_avg_, well_opened_this_step);
195 }
196 }
197
198 template<typename TypeTag>
199 void
201 beginReportStep(const int timeStepIdx)
202 {
203 DeferredLogger local_deferredLogger{};
204
205 this->report_step_starts_ = true;
206 this->report_step_start_events_ = this->schedule()[timeStepIdx].wellgroup_events();
207
208 this->rateConverter_ = std::make_unique<RateConverterType>
209 (this->phase_usage_, std::vector<int>(this->local_num_cells_, 0));
210
211 {
212 // WELPI scaling runs at start of report step.
213 const auto enableWellPIScaling = true;
214 this->initializeLocalWellStructure(timeStepIdx, enableWellPIScaling);
215 }
216
217 this->initializeGroupStructure(timeStepIdx);
218
219 const auto& comm = this->simulator_.vanguard().grid().comm();
220
222 {
223 // Create facility for calculating reservoir voidage volumes for
224 // purpose of RESV controls.
225 this->rateConverter_->template defineState<ElementContext>(this->simulator_);
226
227 // Update VFP properties.
228 {
229 const auto& sched_state = this->schedule()[timeStepIdx];
230
231 this->vfp_properties_ = std::make_unique<VFPProperties<Scalar>>
232 (sched_state.vfpinj(), sched_state.vfpprod(), this->wellState());
233 }
234 }
235 OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
236 "beginReportStep() failed: ",
237 this->terminal_output_, comm)
238
239 // Store the current well and group states in order to recover in
240 // the case of failed iterations
241 this->commitWGState();
242
243 this->wellStructureChangedDynamically_ = false;
244 }
245
246
247
248
249
250 template <typename TypeTag>
251 void
253 initializeLocalWellStructure(const int reportStepIdx,
254 const bool enableWellPIScaling)
255 {
256 DeferredLogger local_deferredLogger{};
257
258 const auto& comm = this->simulator_.vanguard().grid().comm();
259
260 // Wells_ecl_ holds this rank's wells, both open and stopped/shut.
261 this->wells_ecl_ = this->getLocalWells(reportStepIdx);
262 this->local_parallel_well_info_ =
263 this->createLocalParallelWellInfo(this->wells_ecl_);
264
265 // At least initializeWellState() might be throw an exception in
266 // UniformTabulated2DFunction. Playing it safe by extending the
267 // scope a bit.
269 {
270 this->initializeWellPerfData();
271 this->initializeWellState(reportStepIdx);
272 this->wbp_.initializeWBPCalculationService();
273
274 if (this->param_.use_multisegment_well_ && this->anyMSWellOpenLocal()) {
275 this->wellState().initWellStateMSWell(this->wells_ecl_, &this->prevWellState());
276 }
277
278 this->initializeWellProdIndCalculators();
279
280 if (enableWellPIScaling && this->schedule()[reportStepIdx].events()
281 .hasEvent(ScheduleEvents::Events::WELL_PRODUCTIVITY_INDEX))
282 {
283 this->runWellPIScaling(reportStepIdx, local_deferredLogger);
284 }
285 }
286 OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
287 "Failed to initialize local well structure: ",
288 this->terminal_output_, comm)
289 }
290
291
292
293
294
295 template <typename TypeTag>
296 void
298 initializeGroupStructure(const int reportStepIdx)
299 {
300 DeferredLogger local_deferredLogger{};
301
302 const auto& comm = this->simulator_.vanguard().grid().comm();
303
305 {
306 const auto& fieldGroup =
307 this->schedule().getGroup("FIELD", reportStepIdx);
308
310 this->schedule(),
311 this->summaryState(),
312 reportStepIdx,
313 this->groupState());
314
315 // Define per region average pressure calculators for use by
316 // pressure maintenance groups (GPMAINT keyword).
317 if (this->schedule()[reportStepIdx].has_gpmaint()) {
319 (fieldGroup,
320 this->schedule(),
321 reportStepIdx,
322 this->eclState_.fieldProps(),
323 this->phase_usage_,
324 this->regionalAveragePressureCalculator_);
325 }
326 }
327 OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
328 "Failed to initialize group structure: ",
329 this->terminal_output_, comm)
330 }
331
332
333
334
335
336 // called at the beginning of a time step
337 template<typename TypeTag>
338 void
341 {
342 OPM_TIMEBLOCK(beginTimeStep);
343
344 this->updateAverageFormationFactor();
345
346 DeferredLogger local_deferredLogger;
347
348 this->switched_prod_groups_.clear();
349 this->switched_inj_groups_.clear();
350
351 if (this->wellStructureChangedDynamically_) {
352 // Something altered the well structure/topology. Possibly
353 // WELSPECS/COMPDAT and/or WELOPEN run from an ACTIONX block.
354 // Reconstruct the local wells to account for the new well
355 // structure.
356 const auto reportStepIdx =
357 this->simulator_.episodeIndex();
358
359 // Disable WELPI scaling when well structure is updated in the
360 // middle of a report step.
361 const auto enableWellPIScaling = false;
362
363 this->initializeLocalWellStructure(reportStepIdx, enableWellPIScaling);
364 this->initializeGroupStructure(reportStepIdx);
365
366 this->commitWGState();
367
368 // Reset topology flag to signal that we've handled this
369 // structure change. That way we don't end up here in
370 // subsequent calls to beginTimeStep() unless there's a new
371 // dynamic change to the well structure during a report step.
372 this->wellStructureChangedDynamically_ = false;
373 }
374
375 this->resetWGState();
376
377 const int reportStepIdx = simulator_.episodeIndex();
378 this->updateAndCommunicateGroupData(reportStepIdx,
379 simulator_.model().newtonMethod().numIterations(),
380 param_.nupcol_group_rate_tolerance_, /*update_wellgrouptarget*/ false,
381 local_deferredLogger);
382
383 this->wellState().updateWellsDefaultALQ(this->schedule(), reportStepIdx, this->summaryState());
384 this->wellState().gliftTimeStepInit();
385
386 const double simulationTime = simulator_.time();
388 {
389 // test wells
390 wellTesting(reportStepIdx, simulationTime, local_deferredLogger);
391
392 // create the well container
393 createWellContainer(reportStepIdx);
394
395 // Wells are active if they are active wells on at least one process.
396 const Grid& grid = simulator_.vanguard().grid();
397 this->wells_active_ = grid.comm().max(!this->well_container_.empty());
398
399 // do the initialization for all the wells
400 // TODO: to see whether we can postpone of the intialization of the well containers to
401 // optimize the usage of the following several member variables
402 this->initWellContainer(reportStepIdx);
403
404 // update the updated cell flag
405 std::fill(is_cell_perforated_.begin(), is_cell_perforated_.end(), false);
406 for (auto& well : well_container_) {
407 well->updatePerforatedCell(is_cell_perforated_);
408 }
409
410 // calculate the efficiency factors for each well
411 this->calculateEfficiencyFactors(reportStepIdx);
412
413 if constexpr (has_polymer_)
414 {
415 if (PolymerModule::hasPlyshlog() || getPropValue<TypeTag, Properties::EnablePolymerMW>() ) {
416 this->setRepRadiusPerfLength();
417 }
418 }
419
420 }
421
422 OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger, "beginTimeStep() failed: ",
423 this->terminal_output_, simulator_.vanguard().grid().comm());
424
425 for (auto& well : well_container_) {
426 well->setVFPProperties(this->vfp_properties_.get());
427 well->setGuideRate(&this->guideRate_);
428 }
429
430 this->updateFiltrationModelsPreStep(local_deferredLogger);
431
432 // Close completions due to economic reasons
433 for (auto& well : well_container_) {
434 well->closeCompletions(this->wellTestState());
435 }
436
437 // we need the inj_multiplier from the previous time step
438 this->initInjMult();
439
440 if (alternative_well_rate_init_) {
441 // Update the well rates of well_state_, if only single-phase rates, to
442 // have proper multi-phase rates proportional to rates at bhp zero.
443 // This is done only for producers, as injectors will only have a single
444 // nonzero phase anyway.
445 for (const auto& well : well_container_) {
446 if (well->isProducer() && !well->wellIsStopped()) {
447 well->initializeProducerWellState(simulator_, this->wellState(), local_deferredLogger);
448 }
449 }
450 }
451
452 for (const auto& well : well_container_) {
453 if (well->isVFPActive(local_deferredLogger)){
454 well->setPrevSurfaceRates(this->wellState(), this->prevWellState());
455 }
456 }
457 try {
458 this->updateWellPotentials(reportStepIdx,
459 /*onlyAfterEvent*/true,
460 simulator_.vanguard().summaryConfig(),
461 local_deferredLogger);
462 } catch ( std::runtime_error& e ) {
463 const std::string msg = "A zero well potential is returned for output purposes. ";
464 local_deferredLogger.warning("WELL_POTENTIAL_CALCULATION_FAILED", msg);
465 }
466 this->guide_rate_handler_.setLogger(&local_deferredLogger);
467#ifdef RESERVOIR_COUPLING_ENABLED
468 if (this->isReservoirCouplingMaster()) {
469 this->guide_rate_handler_.receiveMasterGroupPotentialsFromSlaves();
470 }
471#endif
472 //update guide rates
473 this->guide_rate_handler_.updateGuideRates(
474 reportStepIdx, simulationTime, this->wellState(), this->groupState()
475 );
476#ifdef RESERVOIR_COUPLING_ENABLED
477 if (this->isReservoirCouplingSlave()) {
478 this->guide_rate_handler_.sendSlaveGroupPotentialsToMaster(this->groupState());
479 }
480#endif
481 std::string exc_msg;
482 auto exc_type = ExceptionType::NONE;
483 // update gpmaint targets
484 if (this->schedule_[reportStepIdx].has_gpmaint()) {
485 for (const auto& calculator : regionalAveragePressureCalculator_) {
486 calculator.second->template defineState<ElementContext>(simulator_);
487 }
488 const double dt = simulator_.timeStepSize();
489 const Group& fieldGroup = this->schedule().getGroup("FIELD", reportStepIdx);
491 this->schedule_,
492 regionalAveragePressureCalculator_,
493 reportStepIdx,
494 dt,
495 this->wellState(),
496 this->groupState());
497 }
498
499 this->updateAndCommunicateGroupData(reportStepIdx,
500 simulator_.model().newtonMethod().numIterations(),
501 param_.nupcol_group_rate_tolerance_,
502 /*update_wellgrouptarget*/ true,
503 local_deferredLogger);
504 try {
505 // Compute initial well solution for new wells and injectors that change injection type i.e. WAG.
506 for (auto& well : well_container_) {
507 const uint64_t effective_events_mask = ScheduleEvents::WELL_STATUS_CHANGE
508 + ScheduleEvents::INJECTION_TYPE_CHANGED
509 + ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER
510 + ScheduleEvents::NEW_WELL;
511
512 const auto& events = this->schedule()[reportStepIdx].wellgroup_events();
513 const bool event = this->report_step_starts_ && events.hasEvent(well->name(), effective_events_mask);
514 const bool dyn_status_change = this->wellState().well(well->name()).status
515 != this->prevWellState().well(well->name()).status;
516
517 if (event || dyn_status_change) {
518 try {
519 well->updateWellStateWithTarget(simulator_, this->groupState(), this->wellState(), local_deferredLogger);
520 well->calculateExplicitQuantities(simulator_, this->wellState(), local_deferredLogger);
521 well->solveWellEquation(simulator_, this->wellState(), this->groupState(), local_deferredLogger);
522 } catch (const std::exception& e) {
523 const std::string msg = "Compute initial well solution for new well " + well->name() + " failed. Continue with zero initial rates";
524 local_deferredLogger.warning("WELL_INITIAL_SOLVE_FAILED", msg);
525 }
526 }
527 }
528 }
529 // Catch clauses for all errors setting exc_type and exc_msg
530 OPM_PARALLEL_CATCH_CLAUSE(exc_type, exc_msg);
531
532 if (exc_type != ExceptionType::NONE) {
533 const std::string msg = "Compute initial well solution for new wells failed. Continue with zero initial rates";
534 local_deferredLogger.warning("WELL_INITIAL_SOLVE_FAILED", msg);
535 }
536
537 const auto& comm = simulator_.vanguard().grid().comm();
538 logAndCheckForExceptionsAndThrow(local_deferredLogger,
539 exc_type, "beginTimeStep() failed: " + exc_msg, this->terminal_output_, comm);
540
541 }
542
543 template<typename TypeTag>
544 void
546 const double simulationTime,
547 DeferredLogger& deferred_logger)
548 {
549 for (const std::string& well_name : this->getWellsForTesting(timeStepIdx, simulationTime)) {
550 const Well& wellEcl = this->schedule().getWell(well_name, timeStepIdx);
551 if (wellEcl.getStatus() == Well::Status::SHUT)
552 continue;
553
554 WellInterfacePtr well = createWellForWellTest(well_name, timeStepIdx, deferred_logger);
555 // some preparation before the well can be used
556 well->init(&this->phase_usage_, depth_, gravity_, B_avg_, true);
557
558 Scalar well_efficiency_factor = wellEcl.getEfficiencyFactor() *
559 this->wellState().getGlobalEfficiencyScalingFactor(well_name);
560 WellGroupHelpers<Scalar>::accumulateGroupEfficiencyFactor(this->schedule().getGroup(wellEcl.groupName(),
561 timeStepIdx),
562 this->schedule(),
563 timeStepIdx,
564 well_efficiency_factor);
565
566 well->setWellEfficiencyFactor(well_efficiency_factor);
567 well->setVFPProperties(this->vfp_properties_.get());
568 well->setGuideRate(&this->guideRate_);
569
570 // initialize rates/previous rates to prevent zero fractions in vfp-interpolation
571 if (well->isProducer()) {
572 well->initializeProducerWellState(simulator_, this->wellState(), deferred_logger);
573 }
574 if (well->isVFPActive(deferred_logger)) {
575 well->setPrevSurfaceRates(this->wellState(), this->prevWellState());
576 }
577
578 const auto& network = this->schedule()[timeStepIdx].network();
579 if (network.active() && !this->node_pressures_.empty()) {
580 if (well->isProducer()) {
581 const auto it = this->node_pressures_.find(well->wellEcl().groupName());
582 if (it != this->node_pressures_.end()) {
583 // The well belongs to a group which has a network nodal pressure,
584 // set the dynamic THP constraint based on the network nodal pressure
585 const Scalar nodal_pressure = it->second;
586 well->setDynamicThpLimit(nodal_pressure);
587 }
588 }
589 }
590 try {
591 using GLiftEclWells = typename GasLiftGroupInfo<Scalar>::GLiftEclWells;
592 GLiftEclWells ecl_well_map;
593 gaslift_.initGliftEclWellMap(well_container_, ecl_well_map);
594 well->wellTesting(simulator_,
595 simulationTime,
596 this->wellState(),
597 this->groupState(),
598 this->wellTestState(),
599 this->phase_usage_,
600 ecl_well_map,
601 this->well_open_times_,
602 deferred_logger);
603 } catch (const std::exception& e) {
604 const std::string msg = fmt::format("Exception during testing of well: {}. The well will not open.\n Exception message: {}", wellEcl.name(), e.what());
605 deferred_logger.warning("WELL_TESTING_FAILED", msg);
606 }
607 }
608 }
609
610
611
612
613
614 // called at the end of a report step
615 template<typename TypeTag>
616 void
619 {
620 // Clear the communication data structures for above values.
621 for (auto&& pinfo : this->local_parallel_well_info_)
622 {
623 pinfo.get().clear();
624 }
625 }
626
627
628
629
630
631 // called at the end of a report step
632 template<typename TypeTag>
635 lastReport() const {return last_report_; }
636
637
638
639
640
641 // called at the end of a time step
642 template<typename TypeTag>
643 void
645 timeStepSucceeded(const double simulationTime, const double dt)
646 {
647 this->closed_this_step_.clear();
648
649 // time step is finished and we are not any more at the beginning of an report step
650 this->report_step_starts_ = false;
651 const int reportStepIdx = simulator_.episodeIndex();
652
653 DeferredLogger local_deferredLogger;
654 for (const auto& well : well_container_) {
655 if (getPropValue<TypeTag, Properties::EnablePolymerMW>() && well->isInjector()) {
656 well->updateWaterThroughput(dt, this->wellState());
657 }
658 }
659 // update connection transmissibility factor and d factor (if applicable) in the wellstate
660 for (const auto& well : well_container_) {
661 well->updateConnectionTransmissibilityFactor(simulator_, this->wellState().well(well->indexOfWell()));
662 well->updateConnectionDFactor(simulator_, this->wellState().well(well->indexOfWell()));
663 }
664
665 if (Indices::waterEnabled) {
666 this->updateFiltrationModelsPostStep(dt, FluidSystem::waterPhaseIdx, local_deferredLogger);
667 }
668
669 // WINJMULT: At the end of the time step, update the inj_multiplier saved in WellState for later use
670 this->updateInjMult(local_deferredLogger);
671
672 // report well switching
673 for (const auto& well : well_container_) {
674 well->reportWellSwitching(this->wellState().well(well->indexOfWell()), local_deferredLogger);
675 }
676 // report group switching
677 if (this->terminal_output_) {
678 this->reportGroupSwitching(local_deferredLogger);
679 }
680
681 // update the rate converter with current averages pressures etc in
682 rateConverter_->template defineState<ElementContext>(simulator_);
683
684 // calculate the well potentials
685 try {
686 this->updateWellPotentials(reportStepIdx,
687 /*onlyAfterEvent*/false,
688 simulator_.vanguard().summaryConfig(),
689 local_deferredLogger);
690 } catch ( std::runtime_error& e ) {
691 const std::string msg = "A zero well potential is returned for output purposes. ";
692 local_deferredLogger.warning("WELL_POTENTIAL_CALCULATION_FAILED", msg);
693 }
694
695 updateWellTestState(simulationTime, this->wellTestState());
696
697 // check group sales limits at the end of the timestep
698 const Group& fieldGroup = this->schedule_.getGroup("FIELD", reportStepIdx);
699 this->checkGEconLimits(fieldGroup, simulationTime,
700 simulator_.episodeIndex(), local_deferredLogger);
701 this->checkGconsaleLimits(fieldGroup, this->wellState(),
702 simulator_.episodeIndex(), local_deferredLogger);
703
704 this->calculateProductivityIndexValues(local_deferredLogger);
705
706 this->commitWGState();
707
708 const Opm::Parallel::Communication& comm = grid().comm();
709 DeferredLogger global_deferredLogger = gatherDeferredLogger(local_deferredLogger, comm);
710 if (this->terminal_output_) {
711 global_deferredLogger.logMessages();
712 }
713
714 //reporting output temperatures
715 this->computeWellTemperature();
716 }
717
718
719 template<typename TypeTag>
720 void
723 unsigned elemIdx) const
724 {
725 rate = 0;
726
727 if (!is_cell_perforated_[elemIdx]) {
728 return;
729 }
730
731 for (const auto& well : well_container_)
732 well->addCellRates(rate, elemIdx);
733 }
734
735
736 template<typename TypeTag>
737 template <class Context>
738 void
741 const Context& context,
742 unsigned spaceIdx,
743 unsigned timeIdx) const
744 {
745 rate = 0;
746 int elemIdx = context.globalSpaceIndex(spaceIdx, timeIdx);
747
748 if (!is_cell_perforated_[elemIdx]) {
749 return;
750 }
751
752 for (const auto& well : well_container_)
753 well->addCellRates(rate, elemIdx);
754 }
755
756
757
758 template<typename TypeTag>
759 void
761 initializeWellState(const int timeStepIdx)
762 {
763 const auto pressIx = []()
764 {
765 if (Indices::oilEnabled) { return FluidSystem::oilPhaseIdx; }
766 if (Indices::waterEnabled) { return FluidSystem::waterPhaseIdx; }
767
768 return FluidSystem::gasPhaseIdx;
769 }();
770
771 auto cellPressures = std::vector<Scalar>(this->local_num_cells_, Scalar{0});
772 auto cellTemperatures = std::vector<Scalar>(this->local_num_cells_, Scalar{0});
773
774 auto elemCtx = ElementContext { this->simulator_ };
775 const auto& gridView = this->simulator_.vanguard().gridView();
776
778 for (const auto& elem : elements(gridView, Dune::Partitions::interior)) {
779 elemCtx.updatePrimaryStencil(elem);
780 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
781
782 const auto ix = elemCtx.globalSpaceIndex(/*spaceIdx=*/0, /*timeIdx=*/0);
783 const auto& fs = elemCtx.intensiveQuantities(/*spaceIdx=*/0, /*timeIdx=*/0).fluidState();
784
785 cellPressures[ix] = fs.pressure(pressIx).value();
786 cellTemperatures[ix] = fs.temperature(0).value();
787 }
788 OPM_END_PARALLEL_TRY_CATCH("BlackoilWellModel::initializeWellState() failed: ",
789 this->simulator_.vanguard().grid().comm());
790
791 this->wellState().init(cellPressures, cellTemperatures, this->schedule(), this->wells_ecl_,
792 this->local_parallel_well_info_, timeStepIdx,
793 &this->prevWellState(), this->well_perf_data_,
794 this->summaryState(), simulator_.vanguard().enableDistributedWells());
795 }
796
797
798
799
800
801 template<typename TypeTag>
802 void
804 createWellContainer(const int report_step)
805 {
806 DeferredLogger local_deferredLogger;
807
808 const int nw = this->numLocalWells();
809
810 well_container_.clear();
811
812 if (nw > 0) {
813 well_container_.reserve(nw);
814
815 const auto& wmatcher = this->schedule().wellMatcher(report_step);
816 const auto& wcycle = this->schedule()[report_step].wcycle.get();
817
818 // First loop and check for status changes. This is necessary
819 // as wcycle needs the updated open/close times.
820 std::for_each(this->wells_ecl_.begin(), this->wells_ecl_.end(),
821 [this, &wg_events = this->report_step_start_events_](const auto& well_ecl)
822 {
823 if (!well_ecl.hasConnections()) {
824 // No connections in this well. Nothing to do.
825 return;
826 }
827
828 constexpr auto events_mask = ScheduleEvents::WELL_STATUS_CHANGE |
829 ScheduleEvents::REQUEST_OPEN_WELL;
830 const bool well_status_change =
831 this->report_step_starts_ &&
832 wg_events.hasEvent(well_ecl.name(), events_mask);
833 if (well_status_change) {
834 if (well_ecl.getStatus() == WellStatus::OPEN) {
835 this->well_open_times_.insert_or_assign(well_ecl.name(),
836 this->simulator_.time());
837 this->well_close_times_.erase(well_ecl.name());
838 } else if (well_ecl.getStatus() == WellStatus::SHUT) {
839 this->well_close_times_.insert_or_assign(well_ecl.name(),
840 this->simulator_.time());
841 this->well_open_times_.erase(well_ecl.name());
842 }
843 }
844 });
845
846 // Grab wcycle states. This needs to run before the schedule gets processed
847 const auto cycle_states = wcycle.wellStatus(this->simulator_.time(),
848 wmatcher,
849 this->well_open_times_,
850 this->well_close_times_);
851
852 for (int w = 0; w < nw; ++w) {
853 const Well& well_ecl = this->wells_ecl_[w];
854
855 if (!well_ecl.hasConnections()) {
856 // No connections in this well. Nothing to do.
857 continue;
858 }
859
860 const std::string& well_name = well_ecl.name();
861 const auto well_status = this->schedule()
862 .getWell(well_name, report_step).getStatus();
863
864 if ((well_ecl.getStatus() == Well::Status::SHUT) ||
865 (well_status == Well::Status::SHUT))
866 {
867 // Due to ACTIONX the well might have been closed behind our back.
868 if (well_ecl.getStatus() != Well::Status::SHUT) {
869 this->closed_this_step_.insert(well_name);
870 this->wellState().shutWell(w);
871 }
872
873 this->well_open_times_.erase(well_name);
874 this->well_close_times_.erase(well_name);
875 continue;
876 }
877
878 // A new WCON keywords can re-open a well that was closed/shut due to Physical limit
879 if (this->wellTestState().well_is_closed(well_name)) {
880 // The well was shut this timestep, we are most likely retrying
881 // a timestep without the well in question, after it caused
882 // repeated timestep cuts. It should therefore not be opened,
883 // even if it was new or received new targets this report step.
884 const bool closed_this_step = (this->wellTestState().lastTestTime(well_name) == simulator_.time());
885 // TODO: more checking here, to make sure this standard more specific and complete
886 // maybe there is some WCON keywords will not open the well
887 auto& events = this->wellState().well(w).events;
888 if (events.hasEvent(ScheduleEvents::REQUEST_OPEN_WELL)) {
889 if (!closed_this_step) {
890 this->wellTestState().open_well(well_name);
891 this->wellTestState().open_completions(well_name);
892 this->well_open_times_.insert_or_assign(well_name,
893 this->simulator_.time());
894 this->well_close_times_.erase(well_name);
895 }
896 events.clearEvent(ScheduleEvents::REQUEST_OPEN_WELL);
897 }
898 }
899
900 // TODO: should we do this for all kinds of closing reasons?
901 // something like wellTestState().hasWell(well_name)?
902 bool wellIsStopped = false;
903 if (this->wellTestState().well_is_closed(well_name))
904 {
905 if (well_ecl.getAutomaticShutIn()) {
906 // shut wells are not added to the well container
907 this->wellState().shutWell(w);
908 this->well_close_times_.erase(well_name);
909 this->well_open_times_.erase(well_name);
910 continue;
911 } else {
912 if (!well_ecl.getAllowCrossFlow()) {
913 // stopped wells where cross flow is not allowed
914 // are not added to the well container
915 this->wellState().shutWell(w);
916 this->well_close_times_.erase(well_name);
917 this->well_open_times_.erase(well_name);
918 continue;
919 }
920 // stopped wells are added to the container but marked as stopped
921 this->wellState().stopWell(w);
922 wellIsStopped = true;
923 }
924 }
925
926 // shut wells with zero rante constraints and disallowing
927 if (!well_ecl.getAllowCrossFlow()) {
928 const bool any_zero_rate_constraint = well_ecl.isProducer()
929 ? well_ecl.productionControls(this->summaryState_).anyZeroRateConstraint()
930 : well_ecl.injectionControls(this->summaryState_).anyZeroRateConstraint();
931 if (any_zero_rate_constraint) {
932 // Treat as shut, do not add to container.
933 local_deferredLogger.debug(fmt::format(" Well {} gets shut due to having zero rate constraint and disallowing crossflow ", well_ecl.name()) );
934 this->wellState().shutWell(w);
935 this->well_close_times_.erase(well_name);
936 this->well_open_times_.erase(well_name);
937 continue;
938 }
939 }
940
941 if (well_status == Well::Status::STOP) {
942 this->wellState().stopWell(w);
943 this->well_close_times_.erase(well_name);
944 this->well_open_times_.erase(well_name);
945 wellIsStopped = true;
946 }
947
948 if (!wcycle.empty()) {
949 const auto it = cycle_states.find(well_name);
950 if (it != cycle_states.end()) {
951 if (!it->second) {
952 this->wellState().shutWell(w);
953 continue;
954 } else {
955 this->wellState().openWell(w);
956 }
957 }
958 }
959
960 well_container_.emplace_back(this->createWellPointer(w, report_step));
961
962 if (wellIsStopped) {
963 well_container_.back()->stopWell();
964 this->well_close_times_.erase(well_name);
965 this->well_open_times_.erase(well_name);
966 }
967 }
968
969 if (!wcycle.empty()) {
970 const auto schedule_open =
971 [&wg_events = this->report_step_start_events_](const std::string& name)
972 {
973 return wg_events.hasEvent(name, ScheduleEvents::REQUEST_OPEN_WELL);
974 };
975 for (const auto& [wname, wscale] : wcycle.efficiencyScale(this->simulator_.time(),
976 this->simulator_.timeStepSize(),
977 wmatcher,
978 this->well_open_times_,
979 schedule_open))
980 {
981 this->wellState().updateEfficiencyScalingFactor(wname, wscale);
982 this->schedule_.add_event(ScheduleEvents::WELLGROUP_EFFICIENCY_UPDATE, report_step);
983 }
984 }
985 }
986
987 // Collect log messages and print.
988
989 const Opm::Parallel::Communication& comm = grid().comm();
990 DeferredLogger global_deferredLogger = gatherDeferredLogger(local_deferredLogger, comm);
991 if (this->terminal_output_) {
992 global_deferredLogger.logMessages();
993 }
994
995 this->well_container_generic_.clear();
996 for (auto& w : well_container_)
997 this->well_container_generic_.push_back(w.get());
998
999 const auto& network = this->schedule()[report_step].network();
1000 if (network.active() && !this->node_pressures_.empty()) {
1001 for (auto& well: this->well_container_generic_) {
1002 // Producers only, since we so far only support the
1003 // "extended" network model (properties defined by
1004 // BRANPROP and NODEPROP) which only applies to producers.
1005 if (well->isProducer()) {
1006 const auto it = this->node_pressures_.find(well->wellEcl().groupName());
1007 if (it != this->node_pressures_.end()) {
1008 // The well belongs to a group which has a network nodal pressure,
1009 // set the dynamic THP constraint based on the network nodal pressure
1010 const Scalar nodal_pressure = it->second;
1011 well->setDynamicThpLimit(nodal_pressure);
1012 }
1013 }
1014 }
1015 }
1016
1017 this->wbp_.registerOpenWellsForWBPCalculation();
1018 }
1019
1020
1021
1022
1023
1024 template <typename TypeTag>
1025 typename BlackoilWellModel<TypeTag>::WellInterfacePtr
1027 createWellPointer(const int wellID, const int report_step) const
1028 {
1029 const auto is_multiseg = this->wells_ecl_[wellID].isMultiSegment();
1030
1031 if (! (this->param_.use_multisegment_well_ && is_multiseg)) {
1032 return this->template createTypedWellPointer<StandardWell<TypeTag>>(wellID, report_step);
1033 }
1034 else {
1035 return this->template createTypedWellPointer<MultisegmentWell<TypeTag>>(wellID, report_step);
1036 }
1037 }
1038
1039
1040
1041
1042
1043 template <typename TypeTag>
1044 template <typename WellType>
1045 std::unique_ptr<WellType>
1047 createTypedWellPointer(const int wellID, const int time_step) const
1048 {
1049 // Use the pvtRegionIdx from the top cell
1050 const auto& perf_data = this->well_perf_data_[wellID];
1051
1052 // Cater for case where local part might have no perforations.
1053 const auto pvtreg = perf_data.empty()
1054 ? 0 : this->pvt_region_idx_[perf_data.front().cell_index];
1055
1056 const auto& parallel_well_info = this->local_parallel_well_info_[wellID].get();
1057 const auto global_pvtreg = parallel_well_info.broadcastFirstPerforationValue(pvtreg);
1058
1059 return std::make_unique<WellType>(this->wells_ecl_[wellID],
1060 parallel_well_info,
1061 time_step,
1062 this->param_,
1063 *this->rateConverter_,
1064 global_pvtreg,
1065 this->numComponents(),
1066 this->numPhases(),
1067 wellID,
1068 perf_data);
1069 }
1070
1071
1072
1073
1074
1075 template<typename TypeTag>
1078 createWellForWellTest(const std::string& well_name,
1079 const int report_step,
1080 DeferredLogger& deferred_logger) const
1081 {
1082 // Finding the location of the well in wells_ecl
1083 const auto it = std::find_if(this->wells_ecl_.begin(),
1084 this->wells_ecl_.end(),
1085 [&well_name](const auto& w)
1086 { return well_name == w.name(); });
1087 // It should be able to find in wells_ecl.
1088 if (it == this->wells_ecl_.end()) {
1089 OPM_DEFLOG_THROW(std::logic_error,
1090 fmt::format("Could not find well {} in wells_ecl ", well_name),
1091 deferred_logger);
1092 }
1093
1094 const int pos = static_cast<int>(std::distance(this->wells_ecl_.begin(), it));
1095 return this->createWellPointer(pos, report_step);
1096 }
1097
1098
1099
1100 template<typename TypeTag>
1101 void
1104 {
1105 OPM_TIMEFUNCTION();
1106 const double dt = this->simulator_.timeStepSize();
1107 // TODO: should we also have the group and network backed-up here in case the solution did not get converged?
1108 auto& well_state = this->wellState();
1109
1110 const bool changed_well_group = updateWellControlsAndNetwork(true, dt, deferred_logger);
1111 assembleWellEqWithoutIteration(dt, deferred_logger);
1112 const bool converged = this->getWellConvergence(this->B_avg_, true).converged() && !changed_well_group;
1113
1115 for (auto& well : this->well_container_) {
1116 well->solveEqAndUpdateWellState(simulator_, well_state, deferred_logger);
1117 }
1118 OPM_END_PARALLEL_TRY_CATCH("BlackoilWellModel::doPreStepNetworkRebalance() failed: ",
1119 this->simulator_.vanguard().grid().comm());
1120
1121 if (!converged) {
1122 const std::string msg = fmt::format("Initial (pre-step) network balance did not converge.");
1123 deferred_logger.warning(msg);
1124 }
1125 }
1126
1127
1128
1129
1130 template<typename TypeTag>
1131 void
1133 assemble(const int iterationIdx,
1134 const double dt)
1135 {
1136 OPM_TIMEFUNCTION();
1137 DeferredLogger local_deferredLogger;
1138
1139 this->guide_rate_handler_.setLogger(&local_deferredLogger);
1141 if (gaslift_.terminalOutput()) {
1142 const std::string msg =
1143 fmt::format("assemble() : iteration {}" , iterationIdx);
1144 gaslift_.gliftDebug(msg, local_deferredLogger);
1145 }
1146 }
1147 last_report_ = SimulatorReportSingle();
1148 Dune::Timer perfTimer;
1149 perfTimer.start();
1150 this->closed_offending_wells_.clear();
1151
1152 {
1153 const int episodeIdx = simulator_.episodeIndex();
1154 const auto& network = this->schedule()[episodeIdx].network();
1155 if (!this->wellsActive() && !network.active()) {
1156 return;
1157 }
1158 }
1159
1160 if (iterationIdx == 0 && this->wellsActive()) {
1161 OPM_TIMEBLOCK(firstIterationAssmble);
1162 // try-catch is needed here as updateWellControls
1163 // contains global communication and has either to
1164 // be reached by all processes or all need to abort
1165 // before.
1167 {
1168 calculateExplicitQuantities(local_deferredLogger);
1169 prepareTimeStep(local_deferredLogger);
1170 }
1171 OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
1172 "assemble() failed (It=0): ",
1173 this->terminal_output_, grid().comm());
1174 }
1175
1176 const bool well_group_control_changed = updateWellControlsAndNetwork(false, dt, local_deferredLogger);
1177
1178 // even when there is no wells active, the network nodal pressure still need to be updated through updateWellControlsAndNetwork()
1179 // but there is no need to assemble the well equations
1180 if ( ! this->wellsActive() ) {
1181 return;
1182 }
1183
1184 assembleWellEqWithoutIteration(dt, local_deferredLogger);
1185
1186 // if group or well control changes we don't consider the
1187 // case converged
1188 last_report_.well_group_control_changed = well_group_control_changed;
1189 last_report_.assemble_time_well += perfTimer.stop();
1190 }
1191
1192
1193
1194
1195 template<typename TypeTag>
1196 bool
1198 updateWellControlsAndNetwork(const bool mandatory_network_balance,
1199 const double dt,
1200 DeferredLogger& local_deferredLogger)
1201 {
1202 OPM_TIMEFUNCTION();
1203 // not necessarily that we always need to update once of the network solutions
1204 bool do_network_update = true;
1205 bool well_group_control_changed = false;
1206 Scalar network_imbalance = 0.0;
1207 // after certain number of the iterations, we use relaxed tolerance for the network update
1208 const std::size_t iteration_to_relax = param_.network_max_strict_outer_iterations_;
1209 // after certain number of the iterations, we terminate
1210 const std::size_t max_iteration = param_.network_max_outer_iterations_;
1211 std::size_t network_update_iteration = 0;
1212 while (do_network_update) {
1213 if (network_update_iteration >= max_iteration ) {
1214 // only output to terminal if we at the last newton iterations where we try to balance the network.
1215 const int episodeIdx = simulator_.episodeIndex();
1216 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
1217 if (this->shouldBalanceNetwork(episodeIdx, iterationIdx + 1)) {
1218 const std::string msg = fmt::format("Maximum of {:d} network iterations has been used and we stop the update, \n"
1219 "and try again after the next Newton iteration (imbalance = {:.2e} bar, ctrl_change = {})",
1220 max_iteration, network_imbalance*1.0e-5, well_group_control_changed);
1221 local_deferredLogger.debug(msg);
1222 } else {
1223 if (this->terminal_output_) {
1224 const std::string msg = fmt::format("Maximum of {:d} network iterations has been used and we stop the update. \n"
1225 "The simulator will continue with unconverged network results (imbalance = {:.2e} bar, ctrl_change = {})",
1226 max_iteration, network_imbalance*1.0e-5, well_group_control_changed);
1227 local_deferredLogger.info(msg);
1228 }
1229 }
1230 break;
1231 }
1232 if (this->terminal_output_ && (network_update_iteration == iteration_to_relax) ) {
1233 local_deferredLogger.debug("We begin using relaxed tolerance for network update now after " + std::to_string(iteration_to_relax) + " iterations ");
1234 }
1235 const bool relax_network_balance = network_update_iteration >= iteration_to_relax;
1236 // Never optimize gas lift in last iteration, to allow network convergence (unless max_iter < 2)
1237 const bool optimize_gas_lift = ( (network_update_iteration + 1) < std::max(max_iteration, static_cast<std::size_t>(2)) );
1238 std::tie(well_group_control_changed, do_network_update, network_imbalance) =
1239 updateWellControlsAndNetworkIteration(mandatory_network_balance, relax_network_balance, optimize_gas_lift, dt,local_deferredLogger);
1240 ++network_update_iteration;
1241 }
1242 return well_group_control_changed;
1243 }
1244
1245
1246
1247
1248 template<typename TypeTag>
1249 std::tuple<bool, bool, typename BlackoilWellModel<TypeTag>::Scalar>
1251 updateWellControlsAndNetworkIteration(const bool mandatory_network_balance,
1252 const bool relax_network_tolerance,
1253 const bool optimize_gas_lift,
1254 const double dt,
1255 DeferredLogger& local_deferredLogger)
1256 {
1257 OPM_TIMEFUNCTION();
1258 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
1259 const int reportStepIdx = simulator_.episodeIndex();
1260 this->updateAndCommunicateGroupData(reportStepIdx, iterationIdx,
1261 param_.nupcol_group_rate_tolerance_, /*update_wellgrouptarget*/ true, local_deferredLogger);
1262 const auto [more_inner_network_update, network_imbalance] =
1263 updateNetworks(mandatory_network_balance,
1264 local_deferredLogger,
1265 relax_network_tolerance);
1266
1267 bool well_group_control_changed = updateWellControls(local_deferredLogger);
1268
1269 bool alq_updated = false;
1271 {
1272 if (optimize_gas_lift) alq_updated = gaslift_.maybeDoGasLiftOptimize(simulator_,
1273 well_container_,
1274 this->wellState(),
1275 this->groupState(),
1276 local_deferredLogger);
1277
1278 prepareWellsBeforeAssembling(dt, local_deferredLogger);
1279 }
1280 OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
1281 "updateWellControlsAndNetworkIteration() failed: ",
1282 this->terminal_output_, grid().comm());
1283
1284 // update guide rates
1285 if (alq_updated || BlackoilWellModelGuideRates(*this).
1286 guideRateUpdateIsNeeded(reportStepIdx)) {
1287 const double simulationTime = simulator_.time();
1288 // NOTE: For reservoir coupling: Slave group potentials are only communicated
1289 // at the start of the time step, see beginTimeStep(). Here, we assume those
1290 // potentials remain unchanged during the time step when updating guide rates below.
1291 this->guide_rate_handler_.updateGuideRates(
1292 reportStepIdx, simulationTime, this->wellState(), this->groupState()
1293 );
1294 }
1295 // we need to re-iterate the network when the well group controls changed or gaslift/alq is changed or
1296 // the inner iterations are did not converge
1297 const bool more_network_update = this->shouldBalanceNetwork(reportStepIdx, iterationIdx) &&
1298 (more_inner_network_update || well_group_control_changed || alq_updated);
1299 return {well_group_control_changed, more_network_update, network_imbalance};
1300 }
1301
1302 // This function is to be used for well groups in an extended network that act as a subsea manifold
1303 // The wells of such group should have a common THP and total phase rate(s) obeying (if possible)
1304 // the well group constraint set by GCONPROD
1305 template <typename TypeTag>
1306 bool
1308 computeWellGroupThp(const double dt, DeferredLogger& local_deferredLogger)
1309 {
1310 OPM_TIMEFUNCTION();
1311 const int reportStepIdx = this->simulator_.episodeIndex();
1312 const auto& network = this->schedule()[reportStepIdx].network();
1313 const auto& balance = this->schedule()[reportStepIdx].network_balance();
1314 const Scalar thp_tolerance = balance.thp_tolerance();
1315
1316 if (!network.active()) {
1317 return false;
1318 }
1319
1320 auto& well_state = this->wellState();
1321 auto& group_state = this->groupState();
1322
1323 bool well_group_thp_updated = false;
1324 for (const std::string& nodeName : network.node_names()) {
1325 const bool has_choke = network.node(nodeName).as_choke();
1326 if (has_choke) {
1327 const auto& summary_state = this->simulator_.vanguard().summaryState();
1328 const Group& group = this->schedule().getGroup(nodeName, reportStepIdx);
1329
1330 const auto pu = this->phase_usage_;
1331 //TODO: Auto choke combined with RESV control is not supported
1332 std::vector<Scalar> resv_coeff(pu.num_phases, 1.0);
1333 Scalar gratTargetFromSales = 0.0;
1334 if (group_state.has_grat_sales_target(group.name()))
1335 gratTargetFromSales = group_state.grat_sales_target(group.name());
1336
1337 const auto ctrl = group.productionControls(summary_state);
1338 auto cmode_tmp = ctrl.cmode;
1339 Scalar target_tmp{0.0};
1340 bool fld_none = false;
1341 if (cmode_tmp == Group::ProductionCMode::FLD || cmode_tmp == Group::ProductionCMode::NONE) {
1342 fld_none = true;
1343 // Target is set for an ancestor group. Target for autochoke group to be
1344 // derived via group guide rates
1345 const Scalar efficiencyFactor = 1.0;
1346 const Group& parentGroup = this->schedule().getGroup(group.parent(), reportStepIdx);
1348 group.name(),
1349 parentGroup,
1350 well_state,
1351 group_state,
1352 this->schedule(),
1353 summary_state,
1354 resv_coeff,
1355 efficiencyFactor,
1356 reportStepIdx,
1357 pu,
1358 &this->guideRate_,
1359 local_deferredLogger);
1360 target_tmp = target.first;
1361 cmode_tmp = target.second;
1362 }
1363 const auto cmode = cmode_tmp;
1364 WGHelpers::TargetCalculator tcalc(cmode, pu, resv_coeff,
1365 gratTargetFromSales, nodeName, group_state,
1366 group.has_gpmaint_control(cmode));
1367 if (!fld_none)
1368 {
1369 // Target is set for the autochoke group itself
1370 target_tmp = tcalc.groupTarget(ctrl, local_deferredLogger);
1371 }
1372
1373 const Scalar orig_target = target_tmp;
1374
1375 auto mismatch = [&] (auto group_thp) {
1376 Scalar group_rate(0.0);
1377 Scalar rate(0.0);
1378 for (auto& well : this->well_container_) {
1379 std::string well_name = well->name();
1380 auto& ws = well_state.well(well_name);
1381 if (group.hasWell(well_name)) {
1382 well->setDynamicThpLimit(group_thp);
1383 const Well& well_ecl = this->wells_ecl_[well->indexOfWell()];
1384 const auto inj_controls = Well::InjectionControls(0);
1385 const auto prod_controls = well_ecl.productionControls(summary_state);
1386 well->iterateWellEqWithSwitching(this->simulator_, dt, inj_controls, prod_controls, well_state, group_state, local_deferredLogger, false, false);
1387 rate = -tcalc.calcModeRateFromRates(ws.surface_rates);
1388 group_rate += rate;
1389 }
1390 }
1391 return (group_rate - orig_target)/orig_target;
1392 };
1393
1394 const auto upbranch = network.uptree_branch(nodeName);
1395 const auto it = this->node_pressures_.find((*upbranch).uptree_node());
1396 const Scalar nodal_pressure = it->second;
1397 Scalar well_group_thp = nodal_pressure;
1398
1399 std::optional<Scalar> autochoke_thp;
1400 if (auto iter = this->well_group_thp_calc_.find(nodeName); iter != this->well_group_thp_calc_.end()) {
1401 autochoke_thp = this->well_group_thp_calc_.at(nodeName);
1402 }
1403
1404 //Find an initial bracket
1405 std::array<Scalar, 2> range_initial;
1406 if (!autochoke_thp.has_value()){
1407 Scalar min_thp, max_thp;
1408 // Retrieve the terminal pressure of the associated root of the manifold group
1409 std::string node_name = nodeName;
1410 while (!network.node(node_name).terminal_pressure().has_value()) {
1411 auto branch = network.uptree_branch(node_name).value();
1412 node_name = branch.uptree_node();
1413 }
1414 min_thp = network.node(node_name).terminal_pressure().value();
1416 // Narrow down the bracket
1417 Scalar low1, high1;
1418 std::array<Scalar, 2> range = {Scalar{0.9}*min_thp, Scalar{1.1}*max_thp};
1419 std::optional<Scalar> appr_sol;
1420 WellBhpThpCalculator<Scalar>::bruteForceBracketCommonTHP(mismatch, range, low1, high1, appr_sol, 0.0, local_deferredLogger);
1421 min_thp = low1;
1422 max_thp = high1;
1423 range_initial = {min_thp, max_thp};
1424 }
1425
1426 if (!autochoke_thp.has_value() || autochoke_thp.value() > nodal_pressure) {
1427 // The bracket is based on the initial bracket or on a range based on a previous calculated group thp
1428 std::array<Scalar, 2> range = autochoke_thp.has_value() ?
1429 std::array<Scalar, 2>{Scalar{0.9} * autochoke_thp.value(),
1430 Scalar{1.1} * autochoke_thp.value()} : range_initial;
1431 Scalar low, high;
1432 std::optional<Scalar> approximate_solution;
1433 const Scalar tolerance1 = thp_tolerance;
1434 local_deferredLogger.debug("Using brute force search to bracket the group THP");
1435 const bool finding_bracket = WellBhpThpCalculator<Scalar>::bruteForceBracketCommonTHP(mismatch, range, low, high, approximate_solution, tolerance1, local_deferredLogger);
1436
1437 if (approximate_solution.has_value()) {
1438 autochoke_thp = *approximate_solution;
1439 local_deferredLogger.debug("Approximate group THP value found: " + std::to_string(autochoke_thp.value()));
1440 } else if (finding_bracket) {
1441 const Scalar tolerance2 = thp_tolerance;
1442 const int max_iteration_solve = 100;
1443 int iteration = 0;
1444 autochoke_thp = RegulaFalsiBisection<ThrowOnError>::
1445 solve(mismatch, low, high, max_iteration_solve, tolerance2, iteration);
1446 local_deferredLogger.debug(" bracket = [" + std::to_string(low) + ", " + std::to_string(high) + "], " +
1447 "iteration = " + std::to_string(iteration));
1448 local_deferredLogger.debug("Group THP value = " + std::to_string(autochoke_thp.value()));
1449 } else {
1450 autochoke_thp.reset();
1451 local_deferredLogger.debug("Group THP solve failed due to bracketing failure");
1452 }
1453 }
1454 if (autochoke_thp.has_value()) {
1455 well_group_thp_calc_[nodeName] = autochoke_thp.value();
1456 // Note: The node pressure of the auto-choke node is set to well_group_thp in computeNetworkPressures()
1457 // and must be larger or equal to the pressure of the uptree node of its branch.
1458 well_group_thp = std::max(autochoke_thp.value(), nodal_pressure);
1459 }
1460
1461 for (auto& well : this->well_container_) {
1462 std::string well_name = well->name();
1463
1464 if (well->isInjector() || !well->wellEcl().predictionMode())
1465 continue;
1466
1467 if (group.hasWell(well_name)) {
1468 well->setDynamicThpLimit(well_group_thp);
1469 }
1470 const auto& ws = this->wellState().well(well->indexOfWell());
1471 const bool thp_is_limit = ws.production_cmode == Well::ProducerCMode::THP;
1472 if (thp_is_limit) {
1473 well->prepareWellBeforeAssembling(this->simulator_, dt, this->wellState(), this->groupState(), local_deferredLogger);
1474 }
1475 }
1476
1477 // Use the group THP in computeNetworkPressures().
1478 const auto& current_well_group_thp = group_state.is_autochoke_group(nodeName) ? group_state.well_group_thp(nodeName) : 1e30;
1479 if (std::abs(current_well_group_thp - well_group_thp) > balance.pressure_tolerance()) {
1480 well_group_thp_updated = true;
1481 group_state.update_well_group_thp(nodeName, well_group_thp);
1482 }
1483 }
1484 }
1485 return well_group_thp_updated;
1486 }
1487
1488 template<typename TypeTag>
1489 void
1491 assembleWellEq(const double dt, DeferredLogger& deferred_logger)
1492 {
1493 OPM_TIMEFUNCTION();
1494 for (auto& well : well_container_) {
1495 well->assembleWellEq(simulator_, dt, this->wellState(), this->groupState(), deferred_logger);
1496 }
1497 }
1498
1499
1500 template<typename TypeTag>
1501 void
1503 prepareWellsBeforeAssembling(const double dt, DeferredLogger& deferred_logger)
1504 {
1505 OPM_TIMEFUNCTION();
1506 for (auto& well : well_container_) {
1507 well->prepareWellBeforeAssembling(simulator_, dt, this->wellState(), this->groupState(), deferred_logger);
1508 }
1509 }
1510
1511
1512 template<typename TypeTag>
1513 void
1515 assembleWellEqWithoutIteration(const double dt, DeferredLogger& deferred_logger)
1516 {
1517 OPM_TIMEFUNCTION();
1518 // We make sure that all processes throw in case there is an exception
1519 // on one of them (WetGasPvt::saturationPressure might throw if not converged)
1521
1522 for (auto& well: well_container_) {
1523 well->assembleWellEqWithoutIteration(simulator_, dt, this->wellState(), this->groupState(),
1524 deferred_logger);
1525 }
1526 OPM_END_PARALLEL_TRY_CATCH_LOG(deferred_logger, "BlackoilWellModel::assembleWellEqWithoutIteration failed: ",
1527 this->terminal_output_, grid().comm());
1528
1529 }
1530
1531#if COMPILE_GPU_BRIDGE
1532 template<typename TypeTag>
1533 void
1536 {
1537 // prepare for StandardWells
1539
1540 for(unsigned int i = 0; i < well_container_.size(); i++){
1541 auto& well = well_container_[i];
1542 std::shared_ptr<StandardWell<TypeTag> > derived = std::dynamic_pointer_cast<StandardWell<TypeTag> >(well);
1543 if (derived) {
1544 wellContribs.addNumBlocks(derived->linSys().getNumBlocks());
1545 }
1546 }
1547
1548 // allocate memory for data from StandardWells
1549 wellContribs.alloc();
1550
1551 for(unsigned int i = 0; i < well_container_.size(); i++){
1552 auto& well = well_container_[i];
1553 // maybe WellInterface could implement addWellContribution()
1554 auto derived_std = std::dynamic_pointer_cast<StandardWell<TypeTag>>(well);
1555 if (derived_std) {
1556 derived_std->linSys().extract(derived_std->numStaticWellEq, wellContribs);
1557 } else {
1558 auto derived_ms = std::dynamic_pointer_cast<MultisegmentWell<TypeTag> >(well);
1559 if (derived_ms) {
1560 derived_ms->linSys().extract(wellContribs);
1561 } else {
1562 OpmLog::warning("Warning unknown type of well");
1563 }
1564 }
1565 }
1566 }
1567#endif
1568
1569 template<typename TypeTag>
1570 void
1572 addWellContributions(SparseMatrixAdapter& jacobian) const
1573 {
1574 for ( const auto& well: well_container_ ) {
1575 well->addWellContributions(jacobian);
1576 }
1577 }
1578
1579 template<typename TypeTag>
1580 void
1583 const BVector& weights,
1584 const bool use_well_weights) const
1585 {
1586 int nw = this->numLocalWellsEnd();
1587 int rdofs = local_num_cells_;
1588 for ( int i = 0; i < nw; i++ ) {
1589 int wdof = rdofs + i;
1590 jacobian[wdof][wdof] = 1.0;// better scaling ?
1591 }
1592
1593 for (const auto& well : well_container_) {
1594 well->addWellPressureEquations(jacobian,
1595 weights,
1596 pressureVarIndex,
1597 use_well_weights,
1598 this->wellState());
1599 }
1600 }
1601
1602 template <typename TypeTag>
1604 addReservoirSourceTerms(GlobalEqVector& residual,
1605 const std::vector<typename SparseMatrixAdapter::MatrixBlock*>& diagMatAddress) const
1606 {
1607 // NB this loop may write multiple times to the same element
1608 // if a cell is perforated by more than one well, so it should
1609 // not be OpenMP-parallelized.
1610 for (const auto& well : well_container_) {
1611 if (!well->isOperableAndSolvable() && !well->wellIsStopped()) {
1612 continue;
1613 }
1614 const auto& cells = well->cells();
1615 const auto& rates = well->connectionRates();
1616 for (unsigned perfIdx = 0; perfIdx < rates.size(); ++perfIdx) {
1617 unsigned cellIdx = cells[perfIdx];
1618 auto rate = rates[perfIdx];
1619 rate *= -1.0;
1620 VectorBlockType res(0.0);
1621 using MatrixBlockType = typename SparseMatrixAdapter::MatrixBlock;
1622 MatrixBlockType bMat(0.0);
1623 simulator_.model().linearizer().setResAndJacobi(res, bMat, rate);
1624 residual[cellIdx] += res;
1625 *diagMatAddress[cellIdx] += bMat;
1626 }
1627 }
1628 }
1629
1630
1631 template<typename TypeTag>
1632 void
1635 {
1636 int nw = this->numLocalWellsEnd();
1637 int rdofs = local_num_cells_;
1638 for (int i = 0; i < nw; ++i) {
1639 int wdof = rdofs + i;
1640 jacobian.entry(wdof,wdof) = 1.0;// better scaling ?
1641 }
1642 const auto wellconnections = this->getMaxWellConnections();
1643 for (int i = 0; i < nw; ++i) {
1644 const auto& perfcells = wellconnections[i];
1645 for (int perfcell : perfcells) {
1646 int wdof = rdofs + i;
1647 jacobian.entry(wdof, perfcell) = 0.0;
1648 jacobian.entry(perfcell, wdof) = 0.0;
1649 }
1650 }
1651 }
1652
1653
1654 template<typename TypeTag>
1655 void
1658 {
1659 DeferredLogger local_deferredLogger;
1661 {
1662 for (const auto& well : well_container_) {
1663 const auto& cells = well->cells();
1664 x_local_.resize(cells.size());
1665
1666 for (size_t i = 0; i < cells.size(); ++i) {
1667 x_local_[i] = x[cells[i]];
1668 }
1669 well->recoverWellSolutionAndUpdateWellState(simulator_, x_local_,
1670 this->wellState(), local_deferredLogger);
1671 }
1672 }
1673 OPM_END_PARALLEL_TRY_CATCH_LOG(local_deferredLogger,
1674 "recoverWellSolutionAndUpdateWellState() failed: ",
1675 this->terminal_output_, simulator_.vanguard().grid().comm());
1676 }
1677
1678
1679 template<typename TypeTag>
1680 void
1682 recoverWellSolutionAndUpdateWellStateDomain(const BVector& x, const int domainIdx)
1683 {
1684 if (!nldd_) {
1685 OPM_THROW(std::logic_error, "Attempt to call NLDD method without a NLDD solver");
1686 }
1687
1688 return nldd_->recoverWellSolutionAndUpdateWellState(x, domainIdx);
1689 }
1690
1691
1692 template<typename TypeTag>
1695 getWellConvergence(const std::vector<Scalar>& B_avg, bool checkWellGroupControls) const
1696 {
1697
1698 DeferredLogger local_deferredLogger;
1699 // Get global (from all processes) convergence report.
1700 ConvergenceReport local_report;
1701 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
1702 for (const auto& well : well_container_) {
1703 if (well->isOperableAndSolvable() || well->wellIsStopped()) {
1704 local_report += well->getWellConvergence(
1705 simulator_, this->wellState(), B_avg, local_deferredLogger,
1706 iterationIdx > param_.strict_outer_iter_wells_);
1707 } else {
1708 ConvergenceReport report;
1709 using CR = ConvergenceReport;
1710 report.setWellFailed({CR::WellFailure::Type::Unsolvable, CR::Severity::Normal, -1, well->name()});
1711 local_report += report;
1712 }
1713 }
1714
1715 const Opm::Parallel::Communication comm = grid().comm();
1716 DeferredLogger global_deferredLogger = gatherDeferredLogger(local_deferredLogger, comm);
1717 ConvergenceReport report = gatherConvergenceReport(local_report, comm);
1718
1719 // the well_group_control_changed info is already communicated
1720 if (checkWellGroupControls) {
1721 report.setWellGroupTargetsViolated(this->lastReport().well_group_control_changed);
1722 }
1723
1724 if (this->terminal_output_) {
1725 global_deferredLogger.logMessages();
1726
1727 // Log debug messages for NaN or too large residuals.
1728 for (const auto& f : report.wellFailures()) {
1729 if (f.severity() == ConvergenceReport::Severity::NotANumber) {
1730 OpmLog::debug("NaN residual found with phase " + std::to_string(f.phase()) + " for well " + f.wellName());
1731 } else if (f.severity() == ConvergenceReport::Severity::TooLarge) {
1732 OpmLog::debug("Too large residual found with phase " + std::to_string(f.phase()) + " for well " + f.wellName());
1733 }
1734 }
1735 }
1736 return report;
1737 }
1738
1739
1740
1741
1742
1743 template<typename TypeTag>
1744 void
1746 calculateExplicitQuantities(DeferredLogger& deferred_logger) const
1747 {
1748 // TODO: checking isOperableAndSolvable() ?
1749 for (auto& well : well_container_) {
1750 well->calculateExplicitQuantities(simulator_, this->wellState(), deferred_logger);
1751 }
1752 }
1753
1754
1755
1756
1757
1758 template<typename TypeTag>
1759 bool
1761 updateWellControls(DeferredLogger& deferred_logger)
1762 {
1763 OPM_TIMEFUNCTION();
1764 if (!this->wellsActive()) {
1765 return false;
1766 }
1767 const int episodeIdx = simulator_.episodeIndex();
1768 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
1769 const auto& comm = simulator_.vanguard().grid().comm();
1770 size_t iter = 0;
1771 bool changed_well_group = false;
1772 const Group& fieldGroup = this->schedule().getGroup("FIELD", episodeIdx);
1773 // Check group individual constraints.
1774 // iterate a few times to make sure all constraints are honored
1775 const std::size_t max_iter = param_.well_group_constraints_max_iterations_;
1776 while(!changed_well_group && iter < max_iter) {
1777 changed_well_group = updateGroupControls(fieldGroup, deferred_logger, episodeIdx, iterationIdx);
1778
1779 // Check wells' group constraints and communicate.
1780 bool changed_well_to_group = false;
1781 {
1782 OPM_TIMEBLOCK(UpdateWellControls);
1783 // For MS Wells a linear solve is performed below and the matrix might be singular.
1784 // We need to communicate the exception thrown to the others and rethrow.
1786 for (const auto& well : well_container_) {
1788 const bool changed_well = well->updateWellControl(simulator_, mode, this->wellState(), this->groupState(), deferred_logger);
1789 if (changed_well) {
1790 changed_well_to_group = changed_well || changed_well_to_group;
1791 }
1792 }
1793 OPM_END_PARALLEL_TRY_CATCH("BlackoilWellModel: updating well controls failed: ",
1794 simulator_.gridView().comm());
1795 }
1796
1797 changed_well_to_group = comm.sum(static_cast<int>(changed_well_to_group));
1798 if (changed_well_to_group) {
1799 updateAndCommunicate(episodeIdx, iterationIdx, deferred_logger);
1800 changed_well_group = true;
1801 }
1802
1803 // Check individual well constraints and communicate.
1804 bool changed_well_individual = false;
1805 {
1806 // For MS Wells a linear solve is performed below and the matrix might be singular.
1807 // We need to communicate the exception thrown to the others and rethrow.
1809 for (const auto& well : well_container_) {
1811 const bool changed_well = well->updateWellControl(simulator_, mode, this->wellState(), this->groupState(), deferred_logger);
1812 if (changed_well) {
1813 changed_well_individual = changed_well || changed_well_individual;
1814 }
1815 }
1816 OPM_END_PARALLEL_TRY_CATCH("BlackoilWellModel: updating well controls failed: ",
1817 simulator_.gridView().comm());
1818 }
1819
1820 changed_well_individual = comm.sum(static_cast<int>(changed_well_individual));
1821 if (changed_well_individual) {
1822 updateAndCommunicate(episodeIdx, iterationIdx, deferred_logger);
1823 changed_well_group = true;
1824 }
1825 iter++;
1826 }
1827
1828 // update wsolvent fraction for REIN wells
1829 this->updateWsolvent(fieldGroup, episodeIdx, this->nupcolWellState());
1830
1831 return changed_well_group;
1832 }
1833
1834
1835 template<typename TypeTag>
1836 std::tuple<bool, typename BlackoilWellModel<TypeTag>::Scalar>
1838 updateNetworks(const bool mandatory_network_balance,
1839 DeferredLogger& deferred_logger,
1840 const bool relax_network_tolerance)
1841 {
1842 OPM_TIMEFUNCTION();
1843 const int episodeIdx = simulator_.episodeIndex();
1844 const auto& network = this->schedule()[episodeIdx].network();
1845 if (!this->wellsActive() && !network.active()) {
1846 return {false, 0.0};
1847 }
1848
1849 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
1850 const auto& comm = simulator_.vanguard().grid().comm();
1851
1852 // network related
1853 Scalar network_imbalance = 0.0;
1854 bool more_network_update = false;
1855 if (this->shouldBalanceNetwork(episodeIdx, iterationIdx) || mandatory_network_balance) {
1856 OPM_TIMEBLOCK(BalanceNetwork);
1857 const double dt = this->simulator_.timeStepSize();
1858 // Calculate common THP for subsea manifold well group (item 3 of NODEPROP set to YES)
1859 const bool well_group_thp_updated = computeWellGroupThp(dt, deferred_logger);
1860 const int max_number_of_sub_iterations = param_.network_max_sub_iterations_;
1861 const Scalar network_pressure_update_damping_factor = param_.network_pressure_update_damping_factor_;
1862 const Scalar network_max_pressure_update = param_.network_max_pressure_update_in_bars_ * unit::barsa;
1863 bool more_network_sub_update = false;
1864 for (int i = 0; i < max_number_of_sub_iterations; i++) {
1865 const auto local_network_imbalance = this->updateNetworkPressures(episodeIdx, network_pressure_update_damping_factor, network_max_pressure_update);
1866 network_imbalance = comm.max(local_network_imbalance);
1867 const auto& balance = this->schedule()[episodeIdx].network_balance();
1868 constexpr Scalar relaxation_factor = 10.0;
1869 const Scalar tolerance = relax_network_tolerance ? relaxation_factor * balance.pressure_tolerance() : balance.pressure_tolerance();
1870 more_network_sub_update = this->networkActive() && network_imbalance > tolerance;
1871 if (!more_network_sub_update)
1872 break;
1873
1874 for (const auto& well : well_container_) {
1875 if (well->isInjector() || !well->wellEcl().predictionMode())
1876 continue;
1877
1878 const auto it = this->node_pressures_.find(well->wellEcl().groupName());
1879 if (it != this->node_pressures_.end()) {
1880 const auto& ws = this->wellState().well(well->indexOfWell());
1881 const bool thp_is_limit = ws.production_cmode == Well::ProducerCMode::THP;
1882 if (thp_is_limit) {
1883 well->prepareWellBeforeAssembling(this->simulator_, dt, this->wellState(), this->groupState(), deferred_logger);
1884 }
1885 }
1886 }
1887 this->updateAndCommunicateGroupData(episodeIdx, iterationIdx, param_.nupcol_group_rate_tolerance_,
1888 /*update_wellgrouptarget*/ true, deferred_logger);
1889 }
1890 more_network_update = more_network_sub_update || well_group_thp_updated;
1891 }
1892 return { more_network_update, network_imbalance };
1893 }
1894
1895
1896 template<typename TypeTag>
1897 void
1899 updateAndCommunicate(const int reportStepIdx,
1900 const int iterationIdx,
1901 DeferredLogger& deferred_logger)
1902 {
1903 this->updateAndCommunicateGroupData(reportStepIdx,
1904 iterationIdx,
1905 param_.nupcol_group_rate_tolerance_,
1906 /*update_wellgrouptarget*/ true,
1907 deferred_logger);
1908
1909 // updateWellStateWithTarget might throw for multisegment wells hence we
1910 // have a parallel try catch here to thrown on all processes.
1912 // if a well or group change control it affects all wells that are under the same group
1913 for (const auto& well : well_container_) {
1914 // We only want to update wells under group-control here
1915 const auto& ws = this->wellState().well(well->indexOfWell());
1916 if (ws.production_cmode == Well::ProducerCMode::GRUP ||
1917 ws.injection_cmode == Well::InjectorCMode::GRUP)
1918 {
1919 well->updateWellStateWithTarget(simulator_, this->groupState(),
1920 this->wellState(), deferred_logger);
1921 }
1922 }
1923 OPM_END_PARALLEL_TRY_CATCH("BlackoilWellModel::updateAndCommunicate failed: ",
1924 simulator_.gridView().comm())
1925 this->updateAndCommunicateGroupData(reportStepIdx,
1926 iterationIdx,
1927 param_.nupcol_group_rate_tolerance_,
1928 /*update_wellgrouptarget*/ true,
1929 deferred_logger);
1930 }
1931
1932 template<typename TypeTag>
1933 bool
1935 updateGroupControls(const Group& group,
1936 DeferredLogger& deferred_logger,
1937 const int reportStepIdx,
1938 const int iterationIdx)
1939 {
1940 OPM_TIMEFUNCTION();
1941 bool changed = false;
1942 // restrict the number of group switches but only after nupcol iterations.
1943 const int nupcol = this->schedule()[reportStepIdx].nupcol();
1944 const int max_number_of_group_switches = iterationIdx < nupcol ? 9999 : param_.max_number_of_group_switches_;
1945 bool changed_hc = this->checkGroupHigherConstraints( group, deferred_logger, reportStepIdx, max_number_of_group_switches);
1946 if (changed_hc) {
1947 changed = true;
1948 updateAndCommunicate(reportStepIdx, iterationIdx, deferred_logger);
1949 }
1950
1951 bool changed_individual =
1953 updateGroupIndividualControl(group,
1954 reportStepIdx,
1955 max_number_of_group_switches,
1956 this->switched_inj_groups_,
1957 this->switched_prod_groups_,
1958 this->closed_offending_wells_,
1959 this->groupState(),
1960 this->wellState(),
1961 deferred_logger);
1962
1963 if (changed_individual) {
1964 changed = true;
1965 updateAndCommunicate(reportStepIdx, iterationIdx, deferred_logger);
1966 }
1967 // call recursively down the group hierarchy
1968 for (const std::string& groupName : group.groups()) {
1969 bool changed_this = updateGroupControls( this->schedule().getGroup(groupName, reportStepIdx), deferred_logger, reportStepIdx,iterationIdx);
1970 changed = changed || changed_this;
1971 }
1972 return changed;
1973 }
1974
1975 template<typename TypeTag>
1976 void
1978 updateWellTestState(const double& simulationTime, WellTestState& wellTestState) const
1979 {
1980 OPM_TIMEFUNCTION();
1981 DeferredLogger local_deferredLogger;
1982 for (const auto& well : well_container_) {
1983 const auto& wname = well->name();
1984 const auto wasClosed = wellTestState.well_is_closed(wname);
1985 well->checkWellOperability(simulator_,
1986 this->wellState(),
1987 local_deferredLogger);
1988 const bool under_zero_target =
1989 well->wellUnderZeroGroupRateTarget(this->simulator_,
1990 this->wellState(),
1991 local_deferredLogger);
1992 well->updateWellTestState(this->wellState().well(wname),
1993 simulationTime,
1994 /*writeMessageToOPMLog=*/ true,
1995 under_zero_target,
1996 wellTestState,
1997 local_deferredLogger);
1998
1999 if (!wasClosed && wellTestState.well_is_closed(wname)) {
2000 this->closed_this_step_.insert(wname);
2001 }
2002 }
2003
2004 for (const auto& [group_name, to] : this->closed_offending_wells_) {
2005 if (this->hasOpenLocalWell(to.second) &&
2006 !this->wasDynamicallyShutThisTimeStep(to.second))
2007 {
2008 wellTestState.close_well(to.second,
2009 WellTestConfig::Reason::GROUP,
2010 simulationTime);
2011 this->updateClosedWellsThisStep(to.second);
2012 const std::string msg =
2013 fmt::format("Procedure on exceeding {} limit is WELL for group {}. "
2014 "Well {} is {}.",
2015 to.first,
2016 group_name,
2017 to.second,
2018 "shut");
2019 local_deferredLogger.info(msg);
2020 }
2021 }
2022
2023 const Opm::Parallel::Communication comm = grid().comm();
2024 DeferredLogger global_deferredLogger =
2025 gatherDeferredLogger(local_deferredLogger, comm);
2026
2027 if (this->terminal_output_) {
2028 global_deferredLogger.logMessages();
2029 }
2030 }
2031
2032
2033 template<typename TypeTag>
2034 void
2036 const WellState<Scalar>& well_state_copy,
2037 std::string& exc_msg,
2038 ExceptionType::ExcEnum& exc_type,
2039 DeferredLogger& deferred_logger)
2040 {
2041 OPM_TIMEFUNCTION();
2042 const int np = this->numPhases();
2043 std::vector<Scalar> potentials;
2044 const auto& well = well_container_[widx];
2045 std::string cur_exc_msg;
2046 auto cur_exc_type = ExceptionType::NONE;
2047 try {
2048 well->computeWellPotentials(simulator_, well_state_copy, potentials, deferred_logger);
2049 }
2050 // catch all possible exception and store type and message.
2051 OPM_PARALLEL_CATCH_CLAUSE(cur_exc_type, cur_exc_msg);
2052 if (cur_exc_type != ExceptionType::NONE) {
2053 exc_msg += fmt::format("\nFor well {}: {}", well->name(), cur_exc_msg);
2054 }
2055 exc_type = std::max(exc_type, cur_exc_type);
2056 // Store it in the well state
2057 // potentials is resized and set to zero in the beginning of well->ComputeWellPotentials
2058 // and updated only if sucessfull. i.e. the potentials are zero for exceptions
2059 auto& ws = this->wellState().well(well->indexOfWell());
2060 for (int p = 0; p < np; ++p) {
2061 // make sure the potentials are positive
2062 ws.well_potentials[p] = std::max(Scalar{0.0}, potentials[p]);
2063 }
2064 }
2065
2066
2067
2068 template <typename TypeTag>
2069 void
2072 {
2073 for (const auto& wellPtr : this->well_container_) {
2074 this->calculateProductivityIndexValues(wellPtr.get(), deferred_logger);
2075 }
2076 }
2077
2078
2079
2080
2081
2082 template <typename TypeTag>
2083 void
2085 calculateProductivityIndexValuesShutWells(const int reportStepIdx,
2086 DeferredLogger& deferred_logger)
2087 {
2088 // For the purpose of computing PI/II values, it is sufficient to
2089 // construct StandardWell instances only. We don't need to form
2090 // well objects that honour the 'isMultisegment()' flag of the
2091 // corresponding "this->wells_ecl_[shutWell]".
2092
2093 for (const auto& shutWell : this->local_shut_wells_) {
2094 if (!this->wells_ecl_[shutWell].hasConnections()) {
2095 // No connections in this well. Nothing to do.
2096 continue;
2097 }
2098
2099 auto wellPtr = this->template createTypedWellPointer
2100 <StandardWell<TypeTag>>(shutWell, reportStepIdx);
2101
2102 wellPtr->init(&this->phase_usage_, this->depth_, this->gravity_, this->B_avg_, true);
2103
2104 this->calculateProductivityIndexValues(wellPtr.get(), deferred_logger);
2105 }
2106 }
2107
2108
2109
2110
2111
2112 template <typename TypeTag>
2113 void
2116 DeferredLogger& deferred_logger)
2117 {
2118 wellPtr->updateProductivityIndex(this->simulator_,
2119 this->prod_index_calc_[wellPtr->indexOfWell()],
2120 this->wellState(),
2121 deferred_logger);
2122 }
2123
2124
2125
2126 template<typename TypeTag>
2127 void
2129 prepareTimeStep(DeferredLogger& deferred_logger)
2130 {
2131 // Check if there is a network with active prediction wells at this time step.
2132 const auto episodeIdx = simulator_.episodeIndex();
2133 this->updateNetworkActiveState(episodeIdx);
2134
2135 // Rebalance the network initially if any wells in the network have status changes
2136 // (Need to check this before clearing events)
2137 const bool do_prestep_network_rebalance = param_.pre_solve_network_ && this->needPreStepNetworkRebalance(episodeIdx);
2138
2139 for (const auto& well : well_container_) {
2140 auto& events = this->wellState().well(well->indexOfWell()).events;
2141 if (events.hasEvent(WellState<Scalar>::event_mask)) {
2142 well->updateWellStateWithTarget(simulator_, this->groupState(), this->wellState(), deferred_logger);
2143 well->updatePrimaryVariables(simulator_, this->wellState(), deferred_logger);
2144 // There is no new well control change input within a report step,
2145 // so next time step, the well does not consider to have effective events anymore.
2146 events.clearEvent(WellState<Scalar>::event_mask);
2147 }
2148 // these events only work for the first time step within the report step
2149 if (events.hasEvent(ScheduleEvents::REQUEST_OPEN_WELL)) {
2150 events.clearEvent(ScheduleEvents::REQUEST_OPEN_WELL);
2151 }
2152 // solve the well equation initially to improve the initial solution of the well model
2153 if (param_.solve_welleq_initially_ && well->isOperableAndSolvable()) {
2154 try {
2155 well->solveWellEquation(simulator_, this->wellState(), this->groupState(), deferred_logger);
2156 } catch (const std::exception& e) {
2157 const std::string msg = "Compute initial well solution for " + well->name() + " initially failed. Continue with the previous rates";
2158 deferred_logger.warning("WELL_INITIAL_SOLVE_FAILED", msg);
2159 }
2160 }
2161 // If we're using local well solves that include control switches, they also update
2162 // operability, so reset before main iterations begin
2163 well->resetWellOperability();
2164 }
2165 updatePrimaryVariables(deferred_logger);
2166
2167 // Actually do the pre-step network rebalance, using the updated well states and initial solutions
2168 if (do_prestep_network_rebalance) doPreStepNetworkRebalance(deferred_logger);
2169 }
2170
2171 template<typename TypeTag>
2172 void
2175 {
2176 std::vector< Scalar > B_avg(numComponents(), Scalar() );
2177 const auto& grid = simulator_.vanguard().grid();
2178 const auto& gridView = grid.leafGridView();
2179 ElementContext elemCtx(simulator_);
2180
2182 for (const auto& elem : elements(gridView, Dune::Partitions::interior)) {
2183 elemCtx.updatePrimaryStencil(elem);
2184 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
2185
2186 const auto& intQuants = elemCtx.intensiveQuantities(/*spaceIdx=*/0, /*timeIdx=*/0);
2187 const auto& fs = intQuants.fluidState();
2188
2189 for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx)
2190 {
2191 if (!FluidSystem::phaseIsActive(phaseIdx)) {
2192 continue;
2193 }
2194
2195 const unsigned compIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
2196 auto& B = B_avg[ compIdx ];
2197
2198 B += 1 / fs.invB(phaseIdx).value();
2199 }
2200 if constexpr (has_solvent_) {
2201 auto& B = B_avg[solventSaturationIdx];
2202 B += 1 / intQuants.solventInverseFormationVolumeFactor().value();
2203 }
2204 }
2205 OPM_END_PARALLEL_TRY_CATCH("BlackoilWellModel::updateAverageFormationFactor() failed: ", grid.comm())
2206
2207 // compute global average
2208 grid.comm().sum(B_avg.data(), B_avg.size());
2209 B_avg_.resize(B_avg.size());
2210 std::transform(B_avg.begin(), B_avg.end(), B_avg_.begin(),
2211 [gcells = global_num_cells_](const auto bval)
2212 { return bval / gcells; });
2213 }
2214
2215
2216
2217
2218
2219 template<typename TypeTag>
2220 void
2223 {
2224 for (const auto& well : well_container_) {
2225 well->updatePrimaryVariables(simulator_, this->wellState(), deferred_logger);
2226 }
2227 }
2228
2229 template<typename TypeTag>
2230 void
2232 {
2233 const auto& grid = simulator_.vanguard().grid();
2234 const auto& eclProblem = simulator_.problem();
2235 const unsigned numCells = grid.size(/*codim=*/0);
2236
2237 this->pvt_region_idx_.resize(numCells);
2238 for (unsigned cellIdx = 0; cellIdx < numCells; ++cellIdx) {
2239 this->pvt_region_idx_[cellIdx] =
2240 eclProblem.pvtRegionIndex(cellIdx);
2241 }
2242 }
2243
2244 // The number of components in the model.
2245 template<typename TypeTag>
2246 int
2248 {
2249 // The numComponents here does not reflect the actual number of the components in the system.
2250 // It more or less reflects the number of mass conservation equations for the well equations.
2251 // For example, in the current formulation, we do not have the polymer conservation equation
2252 // in the well equations. As a result, for an oil-water-polymer system, this function will return 2.
2253 // In some way, it makes this function appear to be confusing from its name, and we need
2254 // to revisit/revise this function again when extending the variants of system that flow can simulate.
2255 int numComp = this->numPhases() < 3 ? this->numPhases() : FluidSystem::numComponents;
2256 if constexpr (has_solvent_) {
2257 numComp++;
2258 }
2259 return numComp;
2260 }
2261
2262 template<typename TypeTag>
2263 void
2265 {
2266 const auto& eclProblem = simulator_.problem();
2267 depth_.resize(local_num_cells_);
2268 for (unsigned cellIdx = 0; cellIdx < local_num_cells_; ++cellIdx) {
2269 depth_[cellIdx] = eclProblem.dofCenterDepth(cellIdx);
2270 }
2271 }
2272
2273 template<typename TypeTag>
2276 getWell(const std::string& well_name) const
2277 {
2278 // finding the iterator of the well in wells_ecl
2279 auto well = std::find_if(well_container_.begin(),
2280 well_container_.end(),
2281 [&well_name](const WellInterfacePtr& elem)->bool {
2282 return elem->name() == well_name;
2283 });
2284
2285 assert(well != well_container_.end());
2286
2287 return *well;
2288 }
2289
2290 template <typename TypeTag>
2291 int
2293 reportStepIndex() const
2294 {
2295 return std::max(this->simulator_.episodeIndex(), 0);
2296 }
2297
2298
2299
2300
2301
2302 template<typename TypeTag>
2303 void
2305 calcResvCoeff(const int fipnum,
2306 const int pvtreg,
2307 const std::vector<Scalar>& production_rates,
2308 std::vector<Scalar>& resv_coeff)
2309 {
2310 rateConverter_->calcCoeff(fipnum, pvtreg, production_rates, resv_coeff);
2311 }
2312
2313 template<typename TypeTag>
2314 void
2316 calcInjResvCoeff(const int fipnum,
2317 const int pvtreg,
2318 std::vector<Scalar>& resv_coeff)
2319 {
2320 rateConverter_->calcInjCoeff(fipnum, pvtreg, resv_coeff);
2321 }
2322
2323
2324 template <typename TypeTag>
2325 void
2328 {
2329 if constexpr (has_energy_) {
2330 int np = this->numPhases();
2331 Scalar cellInternalEnergy;
2332 Scalar cellBinv;
2333 Scalar cellDensity;
2334 Scalar perfPhaseRate;
2335 const int nw = this->numLocalWells();
2336 for (auto wellID = 0*nw; wellID < nw; ++wellID) {
2337 const Well& well = this->wells_ecl_[wellID];
2338 auto& ws = this->wellState().well(wellID);
2339 if (well.isInjector()) {
2340 if (ws.status != WellStatus::STOP) {
2341 this->wellState().well(wellID).temperature = well.inj_temperature();
2342 continue;
2343 }
2344 }
2345
2346 std::array<Scalar,2> weighted{0.0,0.0};
2347 auto& [weighted_temperature, total_weight] = weighted;
2348
2349 auto& well_info = this->local_parallel_well_info_[wellID].get();
2350 auto& perf_data = ws.perf_data;
2351 auto& perf_phase_rate = perf_data.phase_rates;
2352
2353 using int_type = decltype(this->well_perf_data_[wellID].size());
2354 for (int_type perf = 0, end_perf = this->well_perf_data_[wellID].size(); perf < end_perf; ++perf) {
2355 const int cell_idx = this->well_perf_data_[wellID][perf].cell_index;
2356 const auto& intQuants = simulator_.model().intensiveQuantities(cell_idx, /*timeIdx=*/0);
2357 const auto& fs = intQuants.fluidState();
2358
2359 // we on only have one temperature pr cell any phaseIdx will do
2360 Scalar cellTemperatures = fs.temperature(/*phaseIdx*/0).value();
2361
2362 Scalar weight_factor = 0.0;
2363 for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
2364 if (!FluidSystem::phaseIsActive(phaseIdx)) {
2365 continue;
2366 }
2367 cellInternalEnergy = fs.enthalpy(phaseIdx).value() -
2368 fs.pressure(phaseIdx).value() / fs.density(phaseIdx).value();
2369 cellBinv = fs.invB(phaseIdx).value();
2370 cellDensity = fs.density(phaseIdx).value();
2371 perfPhaseRate = perf_phase_rate[perf*np + phaseIdx];
2372 weight_factor += cellDensity * perfPhaseRate / cellBinv * cellInternalEnergy / cellTemperatures;
2373 }
2374 weight_factor = std::abs(weight_factor) + 1e-13;
2375 total_weight += weight_factor;
2376 weighted_temperature += weight_factor * cellTemperatures;
2377 }
2378 well_info.communication().sum(weighted.data(), 2);
2379 this->wellState().well(wellID).temperature = weighted_temperature / total_weight;
2380 }
2381 }
2382 }
2383
2384
2385 template <typename TypeTag>
2387 assignWellTracerRates(data::Wells& wsrpt) const
2388 {
2389 const auto reportStepIdx = static_cast<unsigned int>(this->reportStepIndex());
2390 const auto& trMod = this->simulator_.problem().tracerModel();
2391
2392 BlackoilWellModelGeneric<Scalar>::assignWellTracerRates(wsrpt, trMod.getWellTracerRates(), reportStepIdx);
2393 BlackoilWellModelGeneric<Scalar>::assignWellTracerRates(wsrpt, trMod.getWellFreeTracerRates(), reportStepIdx);
2394 BlackoilWellModelGeneric<Scalar>::assignWellTracerRates(wsrpt, trMod.getWellSolTracerRates(), reportStepIdx);
2395
2396 this->assignMswTracerRates(wsrpt, trMod.getMswTracerRates(), reportStepIdx);
2397 }
2398
2399} // namespace Opm
2400
2401#endif // OPM_BLACKOILWELLMODEL_IMPL_HEADER_INCLUDED
#define OPM_END_PARALLEL_TRY_CATCH_LOG(obptc_logger, obptc_prefix, obptc_output, comm)
Catch exception, log, and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:202
#define OPM_DEFLOG_THROW(Exception, message, deferred_logger)
Definition: DeferredLoggingErrorHelpers.hpp:45
#define OPM_END_PARALLEL_TRY_CATCH(prefix, comm)
Catch exception and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:192
#define OPM_PARALLEL_CATCH_CLAUSE(obptc_exc_type, obptc_exc_msg)
Inserts catch classes for the parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:166
#define OPM_BEGIN_PARALLEL_TRY_CATCH()
Macro to setup the try of a parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:158
void logAndCheckForExceptionsAndThrow(Opm::DeferredLogger &deferred_logger, Opm::ExceptionType::ExcEnum exc_type, const std::string &message, const bool terminal_output, Opm::Parallel::Communication comm)
Definition: DeferredLoggingErrorHelpers.hpp:111
Class for handling constraints for the blackoil well model.
Definition: BlackoilWellModelConstraints.hpp:41
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
Class for handling the guide rates in the blackoil well model.
Definition: BlackoilWellModelGuideRates.hpp:46
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:102
void initializeGroupStructure(const int reportStepIdx)
Definition: BlackoilWellModel_impl.hpp:298
void init()
Definition: BlackoilWellModel_impl.hpp:165
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: BlackoilWellModel.hpp:111
void initWellContainer(const int reportStepIdx) override
Definition: BlackoilWellModel_impl.hpp:184
void beginReportStep(const int time_step)
Definition: BlackoilWellModel_impl.hpp:201
Dune::FieldVector< Scalar, numEq > VectorBlockType
Definition: BlackoilWellModel.hpp:130
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: BlackoilWellModel.hpp:108
GetPropType< TypeTag, Properties::Grid > Grid
Definition: BlackoilWellModel.hpp:105
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: BlackoilWellModel.hpp:110
void initializeWellState(const int timeStepIdx)
Definition: BlackoilWellModel_impl.hpp:761
const Grid & grid() const
Definition: BlackoilWellModel.hpp:367
const SimulatorReportSingle & lastReport() const
Definition: BlackoilWellModel_impl.hpp:635
void addWellContributions(SparseMatrixAdapter &jacobian) const
Definition: BlackoilWellModel_impl.hpp:1572
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
void beginTimeStep()
Definition: BlackoilWellModel_impl.hpp:340
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: BlackoilWellModel.hpp:112
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Definition: BlackoilWellModel_impl.hpp:1746
void updatePrimaryVariables(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2222
void initializeLocalWellStructure(const int reportStepIdx, const bool enableWellPIScaling)
Definition: BlackoilWellModel_impl.hpp:253
Dune::BlockVector< VectorBlockType > BVector
Definition: BlackoilWellModel.hpp:131
BlackoilWellModel(Simulator &simulator)
Definition: BlackoilWellModel_impl.hpp:157
void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:545
void timeStepSucceeded(const double simulationTime, const double dt)
Definition: BlackoilWellModel_impl.hpp:645
std::shared_ptr< WellInterface< TypeTag > > WellInterfacePtr
Definition: BlackoilWellModel.hpp:188
void createWellContainer(const int report_step) override
Definition: BlackoilWellModel_impl.hpp:804
void endReportStep()
Definition: BlackoilWellModel_impl.hpp:618
Definition: ConvergenceReport.hpp:38
void setWellFailed(const WellFailure &wf)
Definition: ConvergenceReport.hpp:270
void setWellGroupTargetsViolated(const bool wellGroupTargetsViolated)
Definition: ConvergenceReport.hpp:288
const std::vector< WellFailure > & wellFailures() const
Definition: ConvergenceReport.hpp:369
Definition: DeferredLogger.hpp:57
void info(const std::string &tag, const std::string &message)
void warning(const std::string &tag, const std::string &message)
void debug(const std::string &tag, const std::string &message)
std::map< std::string, std::pair< const Well *, int > > GLiftEclWells
Definition: GasLiftGroupInfo.hpp:70
Definition: StandardWell.hpp:60
virtual void init(const PhaseUsage *phase_usage_arg, const std::vector< Scalar > &depth_arg, const Scalar gravity_arg, const std::vector< Scalar > &B_avg, const bool changed_to_open_this_step) override
Definition: StandardWell_impl.hpp:76
Definition: TargetCalculator.hpp:43
RateType calcModeRateFromRates(const std::vector< RateType > &rates) const
Definition: TargetCalculator.hpp:54
Scalar groupTarget(const std::optional< Group::ProductionControls > &ctrl, DeferredLogger &deferred_logger) const
Class for computing BHP limits.
Definition: WellBhpThpCalculator.hpp:41
Definition: WellContributions.hpp:51
void alloc()
Allocate memory for the StandardWells.
void setBlockSize(unsigned int dim, unsigned int dim_wells)
void addNumBlocks(unsigned int numBlocks)
Class for computing well group controls.
Definition: WellGroupControls.hpp:49
static void setCmodeGroup(const Group &group, const Schedule &schedule, const SummaryState &summaryState, const int reportStepIdx, GroupState< Scalar > &group_state)
static void setRegionAveragePressureCalculator(const Group &group, const Schedule &schedule, const int reportStepIdx, const FieldPropsManager &fp, const PhaseUsage &pu, std::map< std::string, std::unique_ptr< AverageRegionalPressureType > > &regionalAveragePressureCalculator)
static void updateGpMaintTargetForGroups(const Group &group, const Schedule &schedule, const RegionalValues &regional_values, const int reportStepIdx, const double dt, const WellState< Scalar > &well_state, GroupState< Scalar > &group_state)
static void accumulateGroupEfficiencyFactor(const Group &group, const Schedule &schedule, const int reportStepIdx, Scalar &factor)
int indexOfWell() const
Index of well in the wells struct and wellState.
Definition: WellInterface.hpp:77
bool updateWellControl(const Simulator &simulator, const IndividualOrGroup iog, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:191
virtual void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const =0
Definition: WellState.hpp:65
ExcEnum
Definition: DeferredLogger.hpp:45
@ NONE
Definition: DeferredLogger.hpp:46
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
Definition: blackoilboundaryratevector.hh:39
Opm::DeferredLogger gatherDeferredLogger(const Opm::DeferredLogger &local_deferredlogger, Parallel::Communication communicator)
Create a global log combining local logs.
ConvergenceReport gatherConvergenceReport(const ConvergenceReport &local_report, Parallel::Communication communicator)
PhaseUsage phaseUsageFromDeck(const EclipseState &eclipseState)
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34