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