24#ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED
25#define OPM_BLACKOILWELLMODEL_IMPL_HEADER_INCLUDED
30#include <opm/grid/utility/cartesianToCompressed.hpp>
31#include <opm/common/utility/numeric/RootFinders.hpp>
33#include <opm/input/eclipse/Schedule/Network/Balance.hpp>
34#include <opm/input/eclipse/Schedule/Network/ExtNetwork.hpp>
35#include <opm/input/eclipse/Schedule/Well/PAvgDynamicSourceData.hpp>
36#include <opm/input/eclipse/Schedule/Well/WellTestConfig.hpp>
38#include <opm/input/eclipse/Units/UnitSystem.hpp>
66#include <fmt/format.h>
69 template<
typename TypeTag>
72 : BlackoilWellModelGeneric<Scalar>(simulator.vanguard().schedule(),
73 simulator.vanguard().summaryState(),
74 simulator.vanguard().eclState(),
76 simulator.gridView().comm())
77 , simulator_(simulator)
79 this->terminal_output_ = (simulator.gridView().comm().rank() == 0)
80 && Parameters::Get<Parameters::EnableTerminalOutput>();
82 local_num_cells_ = simulator_.gridView().size(0);
85 global_num_cells_ = simulator_.vanguard().globalNumCells();
88 auto& parallel_wells = simulator.vanguard().parallelWells();
90 this->parallel_well_info_.reserve(parallel_wells.size());
91 for(
const auto& name_bool : parallel_wells) {
92 this->parallel_well_info_.emplace_back(name_bool, grid().comm());
96 this->alternative_well_rate_init_ =
97 Parameters::Get<Parameters::AlternativeWellRateInit>();
99 using SourceDataSpan =
100 typename PAvgDynamicSourceData<Scalar>::template SourceDataSpan<Scalar>;
102 this->wbpCalculationService_
103 .localCellIndex([
this](
const std::size_t globalIndex)
104 {
return this->compressedIndexForInterior(globalIndex); })
105 .evalCellSource([
this](
const int localCell, SourceDataSpan sourceTerms)
107 using Item =
typename SourceDataSpan::Item;
109 const auto* intQuants = this->simulator_.model()
110 .cachedIntensiveQuantities(localCell, 0);
111 const auto& fs = intQuants->fluidState();
113 sourceTerms.set(Item::PoreVol, intQuants->porosity().value() *
114 this->simulator_.model().dofTotalVolume(localCell));
116 constexpr auto io = FluidSystem::oilPhaseIdx;
117 constexpr auto ig = FluidSystem::gasPhaseIdx;
118 constexpr auto iw = FluidSystem::waterPhaseIdx;
121 const auto haveOil = FluidSystem::phaseIsActive(io);
122 const auto haveGas = FluidSystem::phaseIsActive(ig);
123 const auto haveWat = FluidSystem::phaseIsActive(iw);
125 auto weightedPhaseDensity = [&fs](
const auto ip)
127 return fs.saturation(ip).value() * fs.density(ip).value();
130 if (haveOil) { sourceTerms.set(Item::Pressure, fs.pressure(io).value()); }
131 else if (haveGas) { sourceTerms.set(Item::Pressure, fs.pressure(ig).value()); }
132 else { sourceTerms.set(Item::Pressure, fs.pressure(iw).value()); }
136 if (haveOil) { rho += weightedPhaseDensity(io); }
137 if (haveGas) { rho += weightedPhaseDensity(ig); }
138 if (haveWat) { rho += weightedPhaseDensity(iw); }
140 sourceTerms.set(Item::MixtureDensity, rho);
144 template<
typename TypeTag>
151 template<
typename TypeTag>
156 extractLegacyCellPvtRegionIndex_();
157 extractLegacyDepth_();
159 gravity_ = simulator_.problem().gravity()[2];
161 this->initial_step_ =
true;
164 simulator_.model().addAuxiliaryModule(
this);
166 is_cell_perforated_.resize(local_num_cells_,
false);
170 template<
typename TypeTag>
175 const uint64_t effective_events_mask = ScheduleEvents::WELL_STATUS_CHANGE
176 + ScheduleEvents::NEW_WELL;
177 const auto& events = this->schedule()[reportStepIdx].wellgroup_events();
178 for (
auto& wellPtr : this->well_container_) {
179 const bool well_opened_this_step = this->report_step_starts_ && events.hasEvent(wellPtr->name(), effective_events_mask);
180 wellPtr->init(&this->phase_usage_, this->depth_, this->gravity_,
181 this->local_num_cells_, this->B_avg_, well_opened_this_step);
186 template<
typename TypeTag>
191 if (!param_.matrix_add_well_contributions_) {
196 const auto& schedule_wells = this->schedule().getWellsatEnd();
197 auto possibleFutureConnections = this->schedule().getPossibleFutureConnections();
201 const auto& comm = this->simulator_.vanguard().grid().comm();
203 ser.
broadcast(possibleFutureConnections);
206 for (
const auto& well : schedule_wells)
208 std::vector<int> wellCells = this->getCellsForConnections(well);
210 const auto possibleFutureConnectionSetIt = possibleFutureConnections.find(well.name());
211 if (possibleFutureConnectionSetIt != possibleFutureConnections.end()) {
212 for (
auto& global_index : possibleFutureConnectionSetIt->second) {
213 int compressed_idx = compressedIndexForInterior(global_index);
214 if (compressed_idx >= 0) {
215 wellCells.push_back(compressed_idx);
219 for (
int cellIdx : wellCells) {
220 neighbors[cellIdx].insert(wellCells.begin(),
227 template<
typename TypeTag>
230 linearize(SparseMatrixAdapter& jacobian, GlobalEqVector& res)
233 for (
const auto& well: well_container_) {
236 if (param_.matrix_add_well_contributions_) {
237 well->addWellContributions(jacobian);
244 simulator_.gridView().comm());
248 template<
typename TypeTag>
256 for (
const auto& well: well_container_) {
257 if (well_domain_.at(well->name()) == domain.
index) {
260 if (param_.matrix_add_well_contributions_) {
261 well->addWellContributions(jacobian);
271 template<
typename TypeTag>
278 this->report_step_starts_ =
true;
281 this->rateConverter_ = std::make_unique<RateConverterType>
282 (this->phase_usage_, std::vector<int>(this->local_num_cells_, 0));
286 const auto enableWellPIScaling =
true;
287 this->initializeLocalWellStructure(timeStepIdx, enableWellPIScaling);
290 this->initializeGroupStructure(timeStepIdx);
292 const auto& comm = this->simulator_.vanguard().grid().comm();
298 this->rateConverter_->template defineState<ElementContext>(this->simulator_);
302 const auto& sched_state = this->schedule()[timeStepIdx];
304 this->vfp_properties_ = std::make_unique<VFPProperties<Scalar>>
305 (sched_state.vfpinj(), sched_state.vfpprod(), this->wellState());
309 "beginReportStep() failed: ",
310 this->terminal_output_, comm)
314 this->commitWGState();
316 this->wellStructureChangedDynamically_ =
false;
323 template <
typename TypeTag>
327 const bool enableWellPIScaling)
331 const auto& comm = this->simulator_.vanguard().grid().comm();
334 this->wells_ecl_ = this->getLocalWells(reportStepIdx);
335 this->local_parallel_well_info_ =
336 this->createLocalParallelWellInfo(this->wells_ecl_);
343 this->initializeWellPerfData();
344 this->initializeWellState(reportStepIdx);
345 this->initializeWBPCalculationService();
347 if (this->param_.use_multisegment_well_ && this->anyMSWellOpenLocal()) {
348 this->wellState().initWellStateMSWell(this->wells_ecl_, &this->prevWellState());
351 this->initializeWellProdIndCalculators();
353 if (enableWellPIScaling && this->schedule()[reportStepIdx].events()
354 .hasEvent(ScheduleEvents::Events::WELL_PRODUCTIVITY_INDEX))
356 this->runWellPIScaling(reportStepIdx, local_deferredLogger);
360 "Failed to initialize local well structure: ",
361 this->terminal_output_, comm)
368 template <
typename TypeTag>
375 const auto& comm = this->simulator_.vanguard().grid().comm();
379 const auto& fieldGroup =
380 this->schedule().getGroup(
"FIELD", reportStepIdx);
384 this->summaryState(),
390 if (this->schedule()[reportStepIdx].has_gpmaint()) {
395 this->eclState_.fieldProps(),
397 this->regionalAveragePressureCalculator_);
401 "Failed to initialize group structure: ",
402 this->terminal_output_, comm)
410 template<
typename TypeTag>
415 OPM_TIMEBLOCK(beginTimeStep);
417 this->updateAverageFormationFactor();
421 this->switched_prod_groups_.clear();
422 this->switched_inj_groups_.clear();
424 if (this->wellStructureChangedDynamically_) {
429 const auto reportStepIdx =
430 this->simulator_.episodeIndex();
434 const auto enableWellPIScaling =
false;
436 this->initializeLocalWellStructure(reportStepIdx, enableWellPIScaling);
437 this->initializeGroupStructure(reportStepIdx);
439 this->commitWGState();
445 this->wellStructureChangedDynamically_ =
false;
448 this->resetWGState();
450 const int reportStepIdx = simulator_.episodeIndex();
451 this->updateAndCommunicateGroupData(reportStepIdx,
452 simulator_.model().newtonMethod().numIterations());
454 this->wellState().updateWellsDefaultALQ(this->schedule(), reportStepIdx, this->summaryState());
455 this->wellState().gliftTimeStepInit();
457 const double simulationTime = simulator_.time();
461 wellTesting(reportStepIdx, simulationTime, local_deferredLogger);
464 createWellContainer(reportStepIdx);
467 const Grid& grid = simulator_.vanguard().grid();
468 this->wells_active_ = grid.comm().max(!this->well_container_.empty());
473 this->initWellContainer(reportStepIdx);
476 std::fill(is_cell_perforated_.begin(), is_cell_perforated_.end(),
false);
477 for (
auto& well : well_container_) {
478 well->updatePerforatedCell(is_cell_perforated_);
482 this->calculateEfficiencyFactors(reportStepIdx);
484 if constexpr (has_polymer_)
486 if (PolymerModule::hasPlyshlog() || getPropValue<TypeTag, Properties::EnablePolymerMW>() ) {
487 this->setRepRadiusPerfLength();
493 this->terminal_output_, simulator_.vanguard().grid().comm());
495 for (
auto& well : well_container_) {
496 well->setVFPProperties(this->vfp_properties_.get());
497 well->setGuideRate(&this->guideRate_);
500 this->updateInjFCMult(local_deferredLogger);
503 for (
auto& well : well_container_) {
504 well->closeCompletions(this->wellTestState());
510 const auto& summaryState = simulator_.vanguard().summaryState();
511 if (alternative_well_rate_init_) {
516 for (
auto& well : well_container_) {
517 const bool zero_target = well->stoppedOrZeroRateTarget(simulator_, this->wellState(), local_deferredLogger);
518 if (well->isProducer() && !zero_target) {
519 well->updateWellStateRates(simulator_, this->wellState(), local_deferredLogger);
524 for (
auto& well : well_container_) {
525 if (well->isVFPActive(local_deferredLogger)){
526 well->setPrevSurfaceRates(this->wellState(), this->prevWellState());
532 this->updateWellPotentials(reportStepIdx,
534 simulator_.vanguard().summaryConfig(),
535 local_deferredLogger);
536 }
catch ( std::runtime_error& e ) {
537 const std::string msg =
"A zero well potential is returned for output purposes. ";
538 local_deferredLogger.
warning(
"WELL_POTENTIAL_CALCULATION_FAILED", msg);
542 const auto& comm = simulator_.vanguard().grid().comm();
543 std::vector<Scalar> pot(this->numPhases(), 0.0);
544 const Group& fieldGroup = this->schedule().getGroup(
"FIELD", reportStepIdx);
556 local_deferredLogger);
560 if (this->schedule_[reportStepIdx].has_gpmaint()) {
561 for (
auto& calculator : regionalAveragePressureCalculator_) {
562 calculator.second->template defineState<ElementContext>(simulator_);
564 const double dt = simulator_.timeStepSize();
567 regionalAveragePressureCalculator_,
575 for (
auto& well : well_container_) {
576 const uint64_t effective_events_mask = ScheduleEvents::WELL_STATUS_CHANGE
577 + ScheduleEvents::INJECTION_TYPE_CHANGED
578 + ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER
579 + ScheduleEvents::NEW_WELL;
581 const auto& events = this->schedule()[reportStepIdx].wellgroup_events();
582 const bool event = this->report_step_starts_ && events.hasEvent(well->name(), effective_events_mask);
583 const bool dyn_status_change = this->wellState().well(well->name()).status
584 != this->prevWellState().well(well->name()).status;
586 if (event || dyn_status_change) {
588 well->updateWellStateWithTarget(simulator_, this->groupState(), this->wellState(), local_deferredLogger);
589 well->calculateExplicitQuantities(simulator_, this->wellState(), local_deferredLogger);
590 well->solveWellEquation(simulator_, this->wellState(), this->groupState(), local_deferredLogger);
591 }
catch (
const std::exception& e) {
592 const std::string msg =
"Compute initial well solution for new well " + well->name() +
" failed. Continue with zero initial rates";
593 local_deferredLogger.
warning(
"WELL_INITIAL_SOLVE_FAILED", msg);
602 const std::string msg =
"Compute initial well solution for new wells failed. Continue with zero initial rates";
603 local_deferredLogger.
warning(
"WELL_INITIAL_SOLVE_FAILED", msg);
607 exc_type,
"beginTimeStep() failed: " + exc_msg, this->terminal_output_, comm);
611 template<
typename TypeTag>
614 const double simulationTime,
617 for (
const std::string& well_name : this->getWellsForTesting(timeStepIdx, simulationTime)) {
618 const Well& wellEcl = this->schedule().getWell(well_name, timeStepIdx);
619 if (wellEcl.getStatus() == Well::Status::SHUT)
622 WellInterfacePtr well = createWellForWellTest(well_name, timeStepIdx, deferred_logger);
624 well->init(&this->phase_usage_, depth_, gravity_, local_num_cells_, B_avg_,
true);
626 Scalar well_efficiency_factor = wellEcl.getEfficiencyFactor();
631 well_efficiency_factor);
633 well->setWellEfficiencyFactor(well_efficiency_factor);
634 well->setVFPProperties(this->vfp_properties_.get());
635 well->setGuideRate(&this->guideRate_);
638 if (well->isProducer()) {
639 well->updateWellStateRates(simulator_, this->wellState(), deferred_logger);
641 if (well->isVFPActive(deferred_logger)) {
642 well->setPrevSurfaceRates(this->wellState(), this->prevWellState());
646 well->wellTesting(simulator_, simulationTime, this->wellState(),
647 this->groupState(), this->wellTestState(), deferred_logger);
648 }
catch (
const std::exception& e) {
649 const std::string msg = fmt::format(
"Exception during testing of well: {}. The well will not open.\n Exception message: {}", wellEcl.name(), e.what());
650 deferred_logger.
warning(
"WELL_TESTING_FAILED", msg);
660 template<
typename TypeTag>
666 for (
auto&& pinfo : this->local_parallel_well_info_)
677 template<
typename TypeTag>
687 template<
typename TypeTag>
692 this->closed_this_step_.clear();
695 this->report_step_starts_ =
false;
696 const int reportStepIdx = simulator_.episodeIndex();
699 for (
const auto& well : well_container_) {
700 if (getPropValue<TypeTag, Properties::EnablePolymerMW>() && well->isInjector()) {
701 well->updateWaterThroughput(dt, this->wellState());
705 for (
const auto& well : well_container_) {
706 well->updateConnectionTransmissibilityFactor(simulator_, this->wellState().well(well->indexOfWell()));
707 well->updateConnectionDFactor(simulator_, this->wellState().well(well->indexOfWell()));
710 if (Indices::waterEnabled) {
711 this->updateFiltrationParticleVolume(dt, FluidSystem::waterPhaseIdx);
715 this->updateInjMult(local_deferredLogger);
718 for (
const auto& well : well_container_) {
719 well->reportWellSwitching(this->wellState().well(well->indexOfWell()), local_deferredLogger);
722 if (this->terminal_output_) {
724 for (
const auto& [name, to] : this->switched_prod_groups_) {
725 const Group::ProductionCMode& oldControl = this->prevWGState().group_state.production_control(name);
726 std::string from = Group::ProductionCMode2String(oldControl);
728 std::string msg =
" Production Group " + name
729 +
" control mode changed from ";
732 local_deferredLogger.
info(msg);
735 for (
const auto& [key, to] : this->switched_inj_groups_) {
736 const std::string& name = key.first;
737 const Opm::Phase& phase = key.second;
739 const Group::InjectionCMode& oldControl = this->prevWGState().group_state.injection_control(name, phase);
740 std::string from = Group::InjectionCMode2String(oldControl);
742 std::string msg =
" Injection Group " + name
743 +
" control mode changed from ";
746 local_deferredLogger.
info(msg);
752 rateConverter_->template defineState<ElementContext>(simulator_);
756 this->updateWellPotentials(reportStepIdx,
758 simulator_.vanguard().summaryConfig(),
759 local_deferredLogger);
760 }
catch ( std::runtime_error& e ) {
761 const std::string msg =
"A zero well potential is returned for output purposes. ";
762 local_deferredLogger.
warning(
"WELL_POTENTIAL_CALCULATION_FAILED", msg);
765 updateWellTestState(simulationTime, this->wellTestState());
768 const Group& fieldGroup = this->schedule_.getGroup(
"FIELD", reportStepIdx);
769 this->checkGEconLimits(fieldGroup, simulationTime,
770 simulator_.episodeIndex(), local_deferredLogger);
771 this->checkGconsaleLimits(fieldGroup, this->wellState(),
772 simulator_.episodeIndex(), local_deferredLogger);
774 this->calculateProductivityIndexValues(local_deferredLogger);
776 this->commitWGState();
780 if (this->terminal_output_) {
785 this->computeWellTemperature();
789 template<
typename TypeTag>
793 unsigned elemIdx)
const
797 if (!is_cell_perforated_[elemIdx])
800 for (
const auto& well : well_container_)
801 well->addCellRates(rate, elemIdx);
805 template<
typename TypeTag>
806 template <
class Context>
810 const Context& context,
812 unsigned timeIdx)
const
815 int elemIdx = context.globalSpaceIndex(spaceIdx, timeIdx);
817 if (!is_cell_perforated_[elemIdx])
820 for (
const auto& well : well_container_)
821 well->addCellRates(rate, elemIdx);
826 template<
typename TypeTag>
831 const auto pressIx = []()
833 if (Indices::oilEnabled) {
return FluidSystem::oilPhaseIdx; }
834 if (Indices::waterEnabled) {
return FluidSystem::waterPhaseIdx; }
836 return FluidSystem::gasPhaseIdx;
839 auto cellPressures = std::vector<Scalar>(this->local_num_cells_,
Scalar{0});
840 auto cellTemperatures = std::vector<Scalar>(this->local_num_cells_,
Scalar{0});
843 const auto& gridView = this->simulator_.vanguard().gridView();
846 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
847 elemCtx.updatePrimaryStencil(elem);
848 elemCtx.updatePrimaryIntensiveQuantities(0);
850 const auto ix = elemCtx.globalSpaceIndex(0, 0);
851 const auto& fs = elemCtx.intensiveQuantities(0, 0).fluidState();
853 cellPressures[ix] = fs.pressure(pressIx).value();
854 cellTemperatures[ix] = fs.temperature(0).value();
857 this->simulator_.vanguard().grid().comm());
859 this->wellState().init(cellPressures, cellTemperatures, this->schedule(), this->wells_ecl_,
860 this->local_parallel_well_info_, timeStepIdx,
861 &this->prevWellState(), this->well_perf_data_,
862 this->summaryState());
869 template<
typename TypeTag>
876 const int nw = this->numLocalWells();
878 well_container_.clear();
881 well_container_.reserve(nw);
883 for (
int w = 0; w < nw; ++w) {
884 const Well& well_ecl = this->wells_ecl_[w];
886 if (!well_ecl.hasConnections()) {
891 const std::string& well_name = well_ecl.name();
892 const auto well_status = this->schedule()
893 .getWell(well_name, report_step).getStatus();
895 if ((well_ecl.getStatus() == Well::Status::SHUT) ||
896 (well_status == Well::Status::SHUT))
899 if (well_ecl.getStatus() != Well::Status::SHUT) {
900 this->closed_this_step_.insert(well_name);
901 this->wellState().shutWell(w);
908 if (this->wellTestState().well_is_closed(well_name)) {
913 const bool closed_this_step = (this->wellTestState().lastTestTime(well_name) == simulator_.time());
916 auto& events = this->wellState().well(w).events;
917 if (events.hasEvent(ScheduleEvents::REQUEST_OPEN_WELL)) {
918 if (!closed_this_step) {
919 this->wellTestState().open_well(well_name);
920 this->wellTestState().open_completions(well_name);
922 events.clearEvent(ScheduleEvents::REQUEST_OPEN_WELL);
928 bool wellIsStopped =
false;
929 if (this->wellTestState().well_is_closed(well_name))
931 if (well_ecl.getAutomaticShutIn()) {
933 this->wellState().shutWell(w);
936 if (!well_ecl.getAllowCrossFlow()) {
939 this->wellState().shutWell(w);
943 this->wellState().stopWell(w);
944 wellIsStopped =
true;
949 if (!well_ecl.getAllowCrossFlow()) {
950 const bool any_zero_rate_constraint = well_ecl.isProducer()
951 ? well_ecl.productionControls(this->summaryState_).anyZeroRateConstraint()
952 : well_ecl.injectionControls(this->summaryState_).anyZeroRateConstraint();
953 if (any_zero_rate_constraint) {
955 local_deferredLogger.
debug(fmt::format(
" Well {} gets shut due to having zero rate constraint and disallowing crossflow ", well_ecl.name()) );
956 this->wellState().shutWell(w);
961 if (well_status == Well::Status::STOP) {
962 this->wellState().stopWell(w);
963 wellIsStopped =
true;
966 well_container_.emplace_back(this->createWellPointer(w, report_step));
969 well_container_.back()->stopWell();
977 if (this->terminal_output_) {
981 this->well_container_generic_.clear();
982 for (
auto& w : well_container_)
983 this->well_container_generic_.push_back(w.get());
985 const auto& network = this->schedule()[report_step].network();
986 if (network.active() && !this->node_pressures_.empty()) {
987 for (
auto& well: this->well_container_generic_) {
991 if (well->isProducer()) {
992 const auto it = this->node_pressures_.find(well->wellEcl().groupName());
993 if (it != this->node_pressures_.end()) {
996 const Scalar nodal_pressure = it->second;
997 well->setDynamicThpLimit(nodal_pressure);
1003 this->registerOpenWellsForWBPCalculation();
1010 template <
typename TypeTag>
1015 const auto is_multiseg = this->wells_ecl_[wellID].isMultiSegment();
1017 if (! (this->param_.use_multisegment_well_ && is_multiseg)) {
1018 return this->
template createTypedWellPointer<StandardWell<TypeTag>>(wellID, report_step);
1021 return this->
template createTypedWellPointer<MultisegmentWell<TypeTag>>(wellID, report_step);
1029 template <
typename TypeTag>
1030 template <
typename WellType>
1031 std::unique_ptr<WellType>
1036 const auto& perf_data = this->well_perf_data_[wellID];
1039 const auto pvtreg = perf_data.empty()
1040 ? 0 : this->pvt_region_idx_[perf_data.front().cell_index];
1042 const auto& parallel_well_info = this->local_parallel_well_info_[wellID].get();
1043 const auto global_pvtreg = parallel_well_info.broadcastFirstPerforationValue(pvtreg);
1045 return std::make_unique<WellType>(this->wells_ecl_[wellID],
1049 *this->rateConverter_,
1051 this->numComponents(),
1061 template<
typename TypeTag>
1065 const int report_step,
1069 const int nw_wells_ecl = this->wells_ecl_.size();
1070 int index_well_ecl = 0;
1071 for (; index_well_ecl < nw_wells_ecl; ++index_well_ecl) {
1072 if (well_name == this->wells_ecl_[index_well_ecl].name()) {
1077 if (index_well_ecl == nw_wells_ecl) {
1079 fmt::format(
"Could not find well {} in wells_ecl ", well_name),
1083 return this->createWellPointer(index_well_ecl, report_step);
1088 template<
typename TypeTag>
1092 const double dt = this->simulator_.timeStepSize();
1094 auto& well_state = this->wellState();
1095 const std::size_t max_iter = param_.network_max_iterations_;
1096 bool converged =
false;
1097 std::size_t iter = 0;
1098 bool changed_well_group =
false;
1100 changed_well_group = updateWellControlsAndNetwork(
true, dt, deferred_logger);
1101 assembleWellEqWithoutIteration(dt, deferred_logger);
1102 converged = this->getWellConvergence(this->B_avg_,
true).converged() && !changed_well_group;
1107 for (
auto& well : this->well_container_) {
1108 well->solveEqAndUpdateWellState(simulator_, well_state, deferred_logger);
1110 this->initPrimaryVariablesEvaluation();
1111 }
while (iter < max_iter);
1114 const std::string msg = fmt::format(
"Initial (pre-step) network balance did not get converged with {} iterations, "
1115 "unconverged network balance result will be used", max_iter);
1118 const std::string msg = fmt::format(
"Initial (pre-step) network balance converged with {} iterations", iter);
1119 deferred_logger.
debug(msg);
1126 template<
typename TypeTag>
1134 if (this->glift_debug) {
1135 const std::string msg = fmt::format(
1136 "assemble() : iteration {}" , iterationIdx);
1137 this->gliftDebug(msg, local_deferredLogger);
1140 Dune::Timer perfTimer;
1142 this->closed_offending_wells_.clear();
1145 const int episodeIdx = simulator_.episodeIndex();
1146 const auto& network = this->schedule()[episodeIdx].network();
1147 if (!this->wellsActive() && !network.active()) {
1152 if (iterationIdx == 0 && this->wellsActive()) {
1159 calculateExplicitQuantities(local_deferredLogger);
1160 prepareTimeStep(local_deferredLogger);
1163 "assemble() failed (It=0): ",
1164 this->terminal_output_, grid().comm());
1167 const bool well_group_control_changed = updateWellControlsAndNetwork(
false, dt, local_deferredLogger);
1171 if ( ! this->wellsActive() ) {
1175 assembleWellEqWithoutIteration(dt, local_deferredLogger);
1179 last_report_.well_group_control_changed = well_group_control_changed;
1180 last_report_.assemble_time_well += perfTimer.stop();
1186 template<
typename TypeTag>
1192 bool do_network_update =
true;
1193 bool well_group_control_changed =
false;
1195 const std::size_t iteration_to_relax = param_.network_max_strict_iterations_;
1197 const std::size_t max_iteration = param_.network_max_iterations_;
1198 std::size_t network_update_iteration = 0;
1199 while (do_network_update) {
1200 if (this->terminal_output_ && (network_update_iteration == iteration_to_relax) ) {
1201 local_deferredLogger.
info(
" we begin using relaxed tolerance for network update now after " +
std::to_string(iteration_to_relax) +
" iterations ");
1203 const bool relax_network_balance = network_update_iteration >= iteration_to_relax;
1204 std::tie(do_network_update, well_group_control_changed) =
1205 updateWellControlsAndNetworkIteration(mandatory_network_balance, relax_network_balance, dt,local_deferredLogger);
1206 ++network_update_iteration;
1208 if (network_update_iteration >= max_iteration ) {
1209 if (this->terminal_output_) {
1210 local_deferredLogger.
info(
"maximum of " +
std::to_string(max_iteration) +
" iterations has been used, we stop the network update now. "
1211 "The simulation will continue with unconverged network results");
1216 return well_group_control_changed;
1222 template<
typename TypeTag>
1223 std::pair<bool, bool>
1226 const bool relax_network_tolerance,
1230 auto [well_group_control_changed, more_network_update] =
1231 updateWellControls(mandatory_network_balance, local_deferredLogger, relax_network_tolerance);
1233 bool alq_updated =
false;
1237 initPrimaryVariablesEvaluation();
1239 alq_updated = maybeDoGasLiftOptimize(local_deferredLogger);
1241 prepareWellsBeforeAssembling(dt, local_deferredLogger);
1244 this->terminal_output_, grid().comm());
1247 const int reportStepIdx = simulator_.episodeIndex();
1249 guideRateUpdateIsNeeded(reportStepIdx)) {
1250 const double simulationTime = simulator_.time();
1251 const auto& comm = simulator_.vanguard().grid().comm();
1252 const auto& summaryState = simulator_.vanguard().summaryState();
1253 std::vector<Scalar> pot(this->numPhases(), 0.0);
1254 const Group& fieldGroup = this->schedule().getGroup(
"FIELD", reportStepIdx);
1266 local_deferredLogger);
1269 return {more_network_update, well_group_control_changed};
1275 template <
typename TypeTag>
1280 const int reportStepIdx = this->simulator_.episodeIndex();
1281 const auto& network = this->schedule()[reportStepIdx].network();
1282 const auto& balance = this->schedule()[reportStepIdx].network_balance();
1283 const Scalar thp_tolerance = balance.thp_tolerance();
1285 if (!network.active()) {
1289 auto& well_state = this->wellState();
1290 auto& group_state = this->groupState();
1292 for (
const std::string& nodeName : network.node_names()) {
1293 const bool has_choke = network.node(nodeName).as_choke();
1295 const auto& summary_state = this->simulator_.vanguard().summaryState();
1296 const Group& group = this->schedule().getGroup(nodeName, reportStepIdx);
1297 const auto ctrl = group.productionControls(summary_state);
1298 const auto cmode = ctrl.cmode;
1299 const auto pu = this->phase_usage_;
1301 std::vector<Scalar> resv_coeff(pu.num_phases, 1.0);
1302 Scalar gratTargetFromSales = 0.0;
1303 if (group_state.has_grat_sales_target(group.name()))
1304 gratTargetFromSales = group_state.grat_sales_target(group.name());
1307 gratTargetFromSales, nodeName, group_state,
1308 group.has_gpmaint_control(cmode));
1311 auto mismatch = [&] (
auto group_thp) {
1314 for (
auto& well : this->well_container_) {
1315 std::string well_name = well->name();
1316 auto& ws = well_state.well(well_name);
1317 if (group.hasWell(well_name)) {
1318 well->setDynamicThpLimit(group_thp);
1319 const Well& well_ecl = this->wells_ecl_[well->indexOfWell()];
1320 const auto inj_controls = Well::InjectionControls(0);
1321 const auto prod_controls = well_ecl.productionControls(summary_state);
1322 well->iterateWellEqWithSwitching(this->simulator_, dt, inj_controls, prod_controls, well_state, group_state, local_deferredLogger,
false,
false);
1327 return (group_rate - orig_target)/orig_target;
1330 const auto upbranch = network.uptree_branch(nodeName);
1331 const auto it = this->node_pressures_.find((*upbranch).uptree_node());
1332 const Scalar nodal_pressure = it->second;
1333 Scalar well_group_thp = nodal_pressure;
1335 std::optional<Scalar> autochoke_thp;
1336 if (
auto iter = this->well_group_thp_calc_.find(nodeName); iter != this->well_group_thp_calc_.end()) {
1337 autochoke_thp = this->well_group_thp_calc_.at(nodeName);
1341 std::array<Scalar, 2> range_initial;
1342 if (!autochoke_thp.has_value()){
1345 std::string node_name = nodeName;
1346 while (!network.node(node_name).terminal_pressure().has_value()) {
1347 auto branch = network.uptree_branch(node_name).value();
1348 node_name = branch.uptree_node();
1350 min_thp = network.node(node_name).terminal_pressure().value();
1354 std::array<Scalar, 2> range = {
Scalar{0.9}*min_thp,
Scalar{1.1}*max_thp};
1355 std::optional<Scalar> appr_sol;
1359 range_initial = {min_thp, max_thp};
1362 if (!autochoke_thp.has_value() || autochoke_thp.value() > nodal_pressure) {
1364 std::array<Scalar, 2> range = autochoke_thp.has_value() ?
1365 std::array<Scalar, 2>{
Scalar{0.9} * autochoke_thp.value(),
1366 Scalar{1.1} * autochoke_thp.value()} : range_initial;
1368 std::optional<Scalar> approximate_solution;
1369 const Scalar tolerance1 = thp_tolerance;
1370 local_deferredLogger.
debug(
"Using brute force search to bracket the group THP");
1373 if (approximate_solution.has_value()) {
1374 autochoke_thp = *approximate_solution;
1375 local_deferredLogger.
debug(
"Approximate group THP value found: " +
std::to_string(autochoke_thp.value()));
1376 }
else if (finding_bracket) {
1377 const Scalar tolerance2 = thp_tolerance;
1378 const int max_iteration_solve = 100;
1380 autochoke_thp = RegulaFalsiBisection<ThrowOnError>::
1381 solve(mismatch, low, high, max_iteration_solve, tolerance2, iteration);
1386 autochoke_thp.reset();
1387 local_deferredLogger.
debug(
"Group THP solve failed due to bracketing failure");
1390 if (autochoke_thp.has_value()) {
1391 well_group_thp_calc_[nodeName] = autochoke_thp.value();
1394 well_group_thp = std::max(autochoke_thp.value(), nodal_pressure);
1397 for (
auto& well : this->well_container_) {
1398 std::string well_name = well->name();
1399 if (group.hasWell(well_name)) {
1400 well->setDynamicThpLimit(well_group_thp);
1405 group_state.update_well_group_thp(nodeName, well_group_thp);
1410 template<
typename TypeTag>
1418 Dune::Timer perfTimer;
1422 const int episodeIdx = simulator_.episodeIndex();
1423 const auto& network = this->schedule()[episodeIdx].network();
1424 if (!this->wellsActive() && !network.active()) {
1435 updateWellControlsDomain(local_deferredLogger, domain);
1436 initPrimaryVariablesEvaluationDomain(domain);
1437 assembleWellEqDomain(dt, domain, local_deferredLogger);
1441 if (this->terminal_output_) {
1445 last_report_.converged =
true;
1446 last_report_.assemble_time_well += perfTimer.stop();
1450 template<
typename TypeTag>
1455 bool do_glift_optimization =
false;
1456 int num_wells_changed = 0;
1457 const double simulation_time = simulator_.time();
1458 const Scalar min_wait = simulator_.vanguard().schedule().glo(simulator_.episodeIndex()).min_wait();
1463 if ( simulation_time == this->last_glift_opt_time_ || simulation_time >= (this->last_glift_opt_time_ + min_wait)) {
1464 do_glift_optimization =
true;
1465 this->last_glift_opt_time_ = simulation_time;
1468 if (do_glift_optimization) {
1480 initGliftEclWellMap(ecl_well_map);
1483 simulator_.vanguard().schedule(),
1484 simulator_.vanguard().summaryState(),
1485 simulator_.episodeIndex(),
1486 simulator_.model().newtonMethod().numIterations(),
1491 simulator_.vanguard().grid().comm(),
1494 group_info.initialize();
1495 gasLiftOptimizationStage1(deferred_logger, prod_wells, glift_wells,
1496 group_info, state_map);
1497 this->gasLiftOptimizationStage2(deferred_logger, prod_wells, glift_wells,
1498 group_info, state_map, simulator_.episodeIndex());
1499 if (this->glift_debug) {
1500 this->gliftDebugShowALQ(deferred_logger);
1502 num_wells_changed = glift_wells.size();
1504 num_wells_changed = this->comm_.sum(num_wells_changed);
1505 return num_wells_changed > 0;
1508 template<
typename TypeTag>
1517 auto comm = simulator_.vanguard().grid().comm();
1518 int num_procs = comm.size();
1544 for (
int i = 0; i< num_procs; i++) {
1545 int num_rates_to_sync = 0;
1547 if (comm.rank() == i) {
1549 for (
const auto& well : well_container_) {
1551 if (group_info.
hasWell(well->name())) {
1552 gasLiftOptimizationStage1SingleWell(
1553 well.get(), deferred_logger, prod_wells, glift_wells,
1554 group_info, state_map, groups_to_sync
1558 num_rates_to_sync = groups_to_sync.size();
1560 num_rates_to_sync = comm.sum(num_rates_to_sync);
1561 if (num_rates_to_sync > 0) {
1562 std::vector<int> group_indexes;
1563 group_indexes.reserve(num_rates_to_sync);
1564 std::vector<Scalar> group_alq_rates;
1565 group_alq_rates.reserve(num_rates_to_sync);
1566 std::vector<Scalar> group_oil_rates;
1567 group_oil_rates.reserve(num_rates_to_sync);
1568 std::vector<Scalar> group_gas_rates;
1569 group_gas_rates.reserve(num_rates_to_sync);
1570 std::vector<Scalar> group_water_rates;
1571 group_water_rates.reserve(num_rates_to_sync);
1572 if (comm.rank() == i) {
1573 for (
auto idx : groups_to_sync) {
1574 auto [oil_rate, gas_rate, water_rate, alq] = group_info.
getRates(idx);
1575 group_indexes.push_back(idx);
1576 group_oil_rates.push_back(oil_rate);
1577 group_gas_rates.push_back(gas_rate);
1578 group_water_rates.push_back(water_rate);
1579 group_alq_rates.push_back(alq);
1582 group_indexes.resize(num_rates_to_sync);
1583 group_oil_rates.resize(num_rates_to_sync);
1584 group_gas_rates.resize(num_rates_to_sync);
1585 group_water_rates.resize(num_rates_to_sync);
1586 group_alq_rates.resize(num_rates_to_sync);
1590 ser.
broadcast(i, group_indexes, group_oil_rates,
1591 group_gas_rates, group_water_rates, group_alq_rates);
1593 if (comm.rank() != i) {
1594 for (
int j=0; j<num_rates_to_sync; j++) {
1596 group_oil_rates[j], group_gas_rates[j], group_water_rates[j], group_alq_rates[j]);
1599 if (this->glift_debug) {
1601 if (comm.rank() == i) {
1602 counter = this->wellState().gliftGetDebugCounter();
1604 counter = comm.sum(counter);
1605 if (comm.rank() != i) {
1606 this->wellState().gliftSetDebugCounter(counter);
1616 template<
typename TypeTag>
1627 const auto& summary_state = simulator_.vanguard().summaryState();
1628 std::unique_ptr<GasLiftSingleWell> glift
1629 = std::make_unique<GasLiftSingleWell>(
1630 *well, simulator_, summary_state,
1631 deferred_logger, this->wellState(), this->groupState(),
1632 group_info, sync_groups, this->comm_, this->glift_debug);
1633 auto state = glift->runOptimize(
1634 simulator_.model().newtonMethod().numIterations());
1636 state_map.insert({well->
name(), std::move(state)});
1637 glift_wells.insert({well->
name(), std::move(glift)});
1640 prod_wells.insert({well->
name(), well});
1644 template<
typename TypeTag>
1649 for (
const auto& well: well_container_ ) {
1650 ecl_well_map.try_emplace(
1651 well->name(), &(well->wellEcl()), well->indexOfWell());
1656 template<
typename TypeTag>
1661 for (
auto& well : well_container_) {
1662 well->assembleWellEq(simulator_, dt, this->wellState(), this->groupState(), deferred_logger);
1667 template<
typename TypeTag>
1672 for (
auto& well : well_container_) {
1673 if (well_domain_.at(well->name()) == domain.
index) {
1674 well->assembleWellEq(simulator_, dt, this->wellState(), this->groupState(), deferred_logger);
1680 template<
typename TypeTag>
1685 for (
auto& well : well_container_) {
1686 well->prepareWellBeforeAssembling(simulator_, dt, this->wellState(), this->groupState(), deferred_logger);
1691 template<
typename TypeTag>
1700 for (
auto& well: well_container_) {
1701 well->assembleWellEqWithoutIteration(simulator_, dt, this->wellState(), this->groupState(),
1705 this->terminal_output_, grid().comm());
1710 template<
typename TypeTag>
1715 for (
auto& well : well_container_) {
1722 template<
typename TypeTag>
1727 for (
auto& well : well_container_) {
1732#if COMPILE_BDA_BRIDGE
1733 template<
typename TypeTag>
1741 for(
unsigned int i = 0; i < well_container_.size(); i++){
1742 auto& well = well_container_[i];
1743 std::shared_ptr<StandardWell<TypeTag> > derived = std::dynamic_pointer_cast<StandardWell<TypeTag> >(well);
1745 wellContribs.
addNumBlocks(derived->linSys().getNumBlocks());
1750 wellContribs.
alloc();
1752 for(
unsigned int i = 0; i < well_container_.size(); i++){
1753 auto& well = well_container_[i];
1755 auto derived_std = std::dynamic_pointer_cast<StandardWell<TypeTag>>(well);
1757 derived_std->linSys().extract(derived_std->numStaticWellEq, wellContribs);
1759 auto derived_ms = std::dynamic_pointer_cast<MultisegmentWell<TypeTag> >(well);
1761 derived_ms->linSys().extract(wellContribs);
1763 OpmLog::warning(
"Warning unknown type of well");
1771 template<
typename TypeTag>
1776 if (this->well_container_.empty()) {
1780 if( scaleAddRes_.size() != Ax.size() ) {
1781 scaleAddRes_.resize( Ax.size() );
1786 apply( x, scaleAddRes_ );
1788 Ax.axpy( alpha, scaleAddRes_ );
1791 template<
typename TypeTag>
1796 for (
const auto& well: well_container_ ) {
1797 well->addWellContributions(jacobian);
1801 template<
typename TypeTag>
1806 int nw = this->numLocalWellsEnd();
1807 int rdofs = local_num_cells_;
1808 for (
int i = 0; i < nw; i++ ){
1809 int wdof = rdofs + i;
1810 jacobian[wdof][wdof] = 1.0;
1813 for (
const auto& well : well_container_ ) {
1814 well->addWellPressureEquations(jacobian, weights, pressureVarIndex, use_well_weights, this->wellState());
1818 template <
typename TypeTag>
1821 std::vector<typename SparseMatrixAdapter::MatrixBlock*>& diagMatAddress)
const
1826 for (
const auto& well : well_container_) {
1827 if (!well->isOperableAndSolvable() && !well->wellIsStopped()) {
1830 const auto& cells = well->cells();
1831 const auto& rates = well->connectionRates();
1832 for (
unsigned perfIdx = 0; perfIdx < rates.size(); ++perfIdx) {
1833 unsigned cellIdx = cells[perfIdx];
1834 auto rate = rates[perfIdx];
1837 using MatrixBlockType =
typename SparseMatrixAdapter::MatrixBlock;
1838 MatrixBlockType bMat(0.0);
1839 simulator_.model().linearizer().setResAndJacobi(res, bMat, rate);
1840 residual[cellIdx] += res;
1841 *diagMatAddress[cellIdx] += bMat;
1847 template<
typename TypeTag>
1852 int nw = this->numLocalWellsEnd();
1853 int rdofs = local_num_cells_;
1854 for(
int i=0; i < nw; i++){
1855 int wdof = rdofs + i;
1856 jacobian.entry(wdof,wdof) = 1.0;
1858 std::vector<std::vector<int>> wellconnections = this->getMaxWellConnections();
1859 for(
int i=0; i < nw; i++){
1860 const auto& perfcells = wellconnections[i];
1861 for(
int perfcell : perfcells){
1862 int wdof = rdofs + i;
1863 jacobian.entry(wdof,perfcell) = 0.0;
1864 jacobian.entry(perfcell, wdof) = 0.0;
1870 template<
typename TypeTag>
1878 for (
auto& well : well_container_) {
1879 well->recoverWellSolutionAndUpdateWellState(simulator_, x, this->wellState(), local_deferredLogger);
1883 "recoverWellSolutionAndUpdateWellState() failed: ",
1884 this->terminal_output_, simulator_.vanguard().grid().comm());
1888 template<
typename TypeTag>
1897 for (
auto& well : well_container_) {
1898 if (well_domain_.at(well->name()) == domain.
index) {
1899 well->recoverWellSolutionAndUpdateWellState(simulator_, x,
1901 local_deferredLogger);
1906 if (this->terminal_output_) {
1912 template<
typename TypeTag>
1917 for (
auto& well : well_container_) {
1918 well->initPrimaryVariablesEvaluation();
1923 template<
typename TypeTag>
1928 for (
auto& well : well_container_) {
1929 if (well_domain_.at(well->name()) == domain.
index) {
1930 well->initPrimaryVariablesEvaluation();
1940 template<
typename TypeTag>
1944 const std::vector<Scalar>& B_avg,
1947 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
1948 const bool relax_tolerance = iterationIdx > param_.strict_outer_iter_wells_;
1951 for (
const auto& well : well_container_) {
1952 if ((well_domain_.at(well->name()) == domain.
index)) {
1953 if (well->isOperableAndSolvable() || well->wellIsStopped()) {
1954 report += well->getWellConvergence(simulator_,
1957 local_deferredLogger,
1962 xreport.
setWellFailed({CR::WellFailure::Type::Unsolvable, CR::Severity::Normal, -1, well->name()});
1969 if (this->terminal_output_) {
1972 local_deferredLogger.
debug(
"NaN residual found with phase " +
std::to_string(f.phase()) +
" for well " + f.wellName());
1974 local_deferredLogger.
debug(
"Too large residual found with phase " +
std::to_string(f.phase()) +
" for well " + f.wellName());
1985 template<
typename TypeTag>
1988 getWellConvergence(
const std::vector<Scalar>& B_avg,
bool checkWellGroupControls)
const
1994 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
1995 for (
const auto& well : well_container_) {
1996 if (well->isOperableAndSolvable() || well->wellIsStopped()) {
1997 local_report += well->getWellConvergence(
1998 simulator_, this->wellState(), B_avg, local_deferredLogger,
1999 iterationIdx > param_.strict_outer_iter_wells_);
2003 report.
setWellFailed({CR::WellFailure::Type::Unsolvable, CR::Severity::Normal, -1, well->name()});
2004 local_report += report;
2013 if (checkWellGroupControls) {
2017 if (this->terminal_output_) {
2021 if (this->terminal_output_) {
2024 OpmLog::debug(
"NaN residual found with phase " +
std::to_string(f.phase()) +
" for well " + f.wellName());
2026 OpmLog::debug(
"Too large residual found with phase " +
std::to_string(f.phase()) +
" for well " + f.wellName());
2037 template<
typename TypeTag>
2043 for (
auto& well : well_container_) {
2044 well->calculateExplicitQuantities(simulator_, this->wellState(), deferred_logger);
2052 template<
typename TypeTag>
2053 std::pair<bool, bool>
2057 const int episodeIdx = simulator_.episodeIndex();
2058 const auto& network = this->schedule()[episodeIdx].network();
2059 if (!this->wellsActive() && !network.active()) {
2060 return {
false,
false};
2063 const int iterationIdx = simulator_.model().newtonMethod().numIterations();
2064 const auto& comm = simulator_.vanguard().grid().comm();
2065 this->updateAndCommunicateGroupData(episodeIdx, iterationIdx);
2068 bool more_network_update =
false;
2069 if (this->shouldBalanceNetwork(episodeIdx, iterationIdx) || mandatory_network_balance) {
2070 const double dt = this->simulator_.timeStepSize();
2072 computeWellGroupThp(dt, deferred_logger);
2073 const auto local_network_imbalance = this->updateNetworkPressures(episodeIdx);
2074 const Scalar network_imbalance = comm.max(local_network_imbalance);
2075 const auto& balance = this->schedule()[episodeIdx].network_balance();
2076 constexpr Scalar relaxation_factor = 10.0;
2077 const Scalar tolerance = relax_network_tolerance ? relaxation_factor * balance.pressure_tolerance() : balance.pressure_tolerance();
2078 more_network_update = this->networkActive() && network_imbalance > tolerance;
2081 bool changed_well_group =
false;
2083 const int nupcol = this->schedule()[episodeIdx].nupcol();
2086 if (iterationIdx <= nupcol) {
2087 const Group& fieldGroup = this->schedule().getGroup(
"FIELD", episodeIdx);
2088 changed_well_group = updateGroupControls(fieldGroup, deferred_logger, episodeIdx, iterationIdx);
2091 bool changed_well_to_group =
false;
2096 for (
const auto& well : well_container_) {
2098 const bool changed_well = well->
updateWellControl(simulator_, mode, this->wellState(), this->groupState(), deferred_logger);
2100 changed_well_to_group = changed_well || changed_well_to_group;
2104 simulator_.gridView().comm());
2107 changed_well_to_group = comm.sum(
static_cast<int>(changed_well_to_group));
2108 if (changed_well_to_group) {
2109 updateAndCommunicate(episodeIdx, iterationIdx, deferred_logger);
2110 changed_well_group =
true;
2114 bool changed_well_individual =
false;
2119 for (
const auto& well : well_container_) {
2121 const bool changed_well = well->
updateWellControl(simulator_, mode, this->wellState(), this->groupState(), deferred_logger);
2123 changed_well_individual = changed_well || changed_well_individual;
2127 simulator_.gridView().comm());
2130 changed_well_individual = comm.sum(
static_cast<int>(changed_well_individual));
2131 if (changed_well_individual) {
2132 updateAndCommunicate(episodeIdx, iterationIdx, deferred_logger);
2133 changed_well_group =
true;
2137 const Group& fieldGroup = this->schedule().getGroup(
"FIELD", episodeIdx);
2138 this->updateWsolvent(fieldGroup, episodeIdx, this->nupcolWellState());
2140 return { changed_well_group, more_network_update };
2143 template<
typename TypeTag>
2148 if ( !this->wellsActive() ) return ;
2154 for (
const auto& well : well_container_) {
2155 if (well_domain_.at(well->name()) == domain.
index) {
2157 well->
updateWellControl(simulator_, mode, this->wellState(), this->groupState(), deferred_logger);
2166 template <
typename TypeTag>
2171 this->wbpCalcMap_.clear();
2172 this->wbpCalcMap_.resize(this->wells_ecl_.size());
2174 this->registerOpenWellsForWBPCalculation();
2176 auto wellID = std::size_t{0};
2177 for (
const auto& well : this->wells_ecl_) {
2178 this->wbpCalcMap_[wellID].wbpCalcIdx_ = this->wbpCalculationService_
2179 .createCalculator(well,
2180 this->local_parallel_well_info_[wellID],
2181 this->conn_idx_map_[wellID].local(),
2182 this->makeWellSourceEvaluatorFactory(wellID));
2187 this->wbpCalculationService_.defineCommunication();
2194 template <
typename TypeTag>
2195 data::WellBlockAveragePressures
2199 auto wbpResult = data::WellBlockAveragePressures{};
2201 using Calculated =
typename PAvgCalculator<Scalar>::Result::WBPMode;
2202 using Output = data::WellBlockAvgPress::Quantity;
2204 this->wbpCalculationService_.collectDynamicValues();
2206 const auto numWells = this->wells_ecl_.size();
2207 for (
auto wellID = 0*numWells; wellID < numWells; ++wellID) {
2208 const auto calcIdx = this->wbpCalcMap_[wellID].wbpCalcIdx_;
2209 const auto& well = this->wells_ecl_[wellID];
2211 if (! well.hasRefDepth()) {
2217 this->wbpCalculationService_
2218 .inferBlockAveragePressures(calcIdx, well.pavg(),
2220 well.getWPaveRefDepth());
2222 const auto& result = this->wbpCalculationService_
2223 .averagePressures(calcIdx);
2225 auto& reported = wbpResult.values[well.name()];
2227 reported[Output::WBP] = result.value(Calculated::WBP);
2228 reported[Output::WBP4] = result.value(Calculated::WBP4);
2229 reported[Output::WBP5] = result.value(Calculated::WBP5);
2230 reported[Output::WBP9] = result.value(Calculated::WBP9);
2240 template <
typename TypeTag>
2245 using Span =
typename PAvgDynamicSourceData<Scalar>::template SourceDataSpan<Scalar>;
2246 using Item =
typename Span::Item;
2250 if (! this->wbpCalcMap_[wellIdx].openWellIdx_.has_value()) {
2252 return []([[maybe_unused]]
const int connIdx, Span sourceTerm)
2258 .set(Item::Pressure , 0.0)
2259 .set(Item::PoreVol , 0.0)
2260 .set(Item::MixtureDensity, 0.0);
2265 return [
this, wellPtr = this->well_container_[*this->wbpCalcMap_[wellIdx].openWellIdx_].get()]
2266 (
const int connIdx, Span sourceTerm)
2272 const auto& connIdxMap =
2273 this->conn_idx_map_[wellPtr->indexOfWell()];
2275 const auto rho = wellPtr->
2276 connectionDensity(connIdxMap.global(connIdx),
2277 connIdxMap.open(connIdx));
2280 .set(Item::Pressure , 0.0)
2281 .set(Item::PoreVol , 0.0)
2282 .set(Item::MixtureDensity, rho);
2291 template <
typename TypeTag>
2296 assert (this->wbpCalcMap_.size() == this->wells_ecl_.size());
2298 for (
auto& wbpCalc : this->wbpCalcMap_) {
2299 wbpCalc.openWellIdx_.reset();
2302 auto openWellIdx =
typename std::vector<WellInterfacePtr>::size_type{0};
2303 for (
const auto* openWell : this->well_container_generic_) {
2304 this->wbpCalcMap_[openWell->indexOfWell()].openWellIdx_ = openWellIdx++;
2312 template<
typename TypeTag>
2316 const int iterationIdx,
2319 this->updateAndCommunicateGroupData(reportStepIdx, iterationIdx);
2325 for (
const auto& well : well_container_) {
2326 well->updateWellStateWithTarget(simulator_, this->groupState(),
2327 this->wellState(), deferred_logger);
2330 simulator_.gridView().comm())
2331 this->updateAndCommunicateGroupData(reportStepIdx, iterationIdx);
2334 template<
typename TypeTag>
2339 const int reportStepIdx,
2340 const int iterationIdx)
2342 bool changed =
false;
2343 bool changed_hc = this->checkGroupHigherConstraints( group, deferred_logger, reportStepIdx);
2346 updateAndCommunicate(reportStepIdx, iterationIdx, deferred_logger);
2349 bool changed_individual =
2351 updateGroupIndividualControl(group,
2353 this->switched_inj_groups_,
2354 this->switched_prod_groups_,
2355 this->closed_offending_wells_,
2360 if (changed_individual) {
2362 updateAndCommunicate(reportStepIdx, iterationIdx, deferred_logger);
2365 for (
const std::string& groupName : group.groups()) {
2366 bool changed_this = updateGroupControls( this->schedule().getGroup(groupName, reportStepIdx), deferred_logger, reportStepIdx,iterationIdx);
2367 changed = changed || changed_this;
2372 template<
typename TypeTag>
2378 for (
const auto& well : well_container_) {
2379 const auto& wname = well->name();
2380 const auto wasClosed = wellTestState.well_is_closed(wname);
2381 well->checkWellOperability(simulator_, this->wellState(), local_deferredLogger);
2382 const bool under_zero_target = well->wellUnderZeroGroupRateTarget(this->simulator_, this->wellState(), local_deferredLogger);
2383 well->updateWellTestState(this->wellState().well(wname), simulationTime,
true, under_zero_target, wellTestState, local_deferredLogger);
2385 if (!wasClosed && wellTestState.well_is_closed(wname)) {
2386 this->closed_this_step_.insert(wname);
2393 for (
const auto& [group_name, to] : this->closed_offending_wells_) {
2394 if (this->hasWell(to.second) && !this->wasDynamicallyShutThisTimeStep(to.second)) {
2395 wellTestState.close_well(to.second, WellTestConfig::Reason::GROUP, simulationTime);
2396 this->updateClosedWellsThisStep(to.second);
2397 const std::string msg =
2398 fmt::format(
"Procedure on exceeding {} limit is WELL for group {}. Well {} is {}.",
2403 global_deferredLogger.
info(msg);
2407 if (this->terminal_output_) {
2413 template<
typename TypeTag>
2417 std::string& exc_msg,
2421 const int np = this->numPhases();
2422 std::vector<Scalar> potentials;
2423 const auto& well = well_container_[widx];
2424 std::string cur_exc_msg;
2427 well->computeWellPotentials(simulator_, well_state_copy, potentials, deferred_logger);
2432 exc_msg += fmt::format(
"\nFor well {}: {}", well->name(), cur_exc_msg);
2434 exc_type = std::max(exc_type, cur_exc_type);
2438 auto& ws = this->wellState().well(well->indexOfWell());
2439 for (
int p = 0; p < np; ++p) {
2441 ws.well_potentials[p] = std::max(
Scalar{0.0}, potentials[p]);
2447 template <
typename TypeTag>
2452 for (
const auto& wellPtr : this->well_container_) {
2453 this->calculateProductivityIndexValues(wellPtr.get(), deferred_logger);
2461 template <
typename TypeTag>
2472 for (
const auto& shutWell : this->local_shut_wells_) {
2473 if (!this->wells_ecl_[shutWell].hasConnections()) {
2478 auto wellPtr = this->
template createTypedWellPointer
2481 wellPtr->
init(&this->phase_usage_, this->depth_, this->gravity_,
2482 this->local_num_cells_, this->B_avg_,
true);
2484 this->calculateProductivityIndexValues(wellPtr.get(), deferred_logger);
2492 template <
typename TypeTag>
2506 template<
typename TypeTag>
2512 const auto episodeIdx = simulator_.episodeIndex();
2513 this->updateNetworkActiveState(episodeIdx);
2517 const bool do_prestep_network_rebalance = this->needPreStepNetworkRebalance(episodeIdx);
2519 for (
const auto& well : well_container_) {
2520 auto& events = this->wellState().well(well->indexOfWell()).events;
2522 well->updateWellStateWithTarget(simulator_, this->groupState(), this->wellState(), deferred_logger);
2523 well->updatePrimaryVariables(simulator_, this->wellState(), deferred_logger);
2524 well->initPrimaryVariablesEvaluation();
2530 if (events.hasEvent(ScheduleEvents::REQUEST_OPEN_WELL)) {
2531 events.clearEvent(ScheduleEvents::REQUEST_OPEN_WELL);
2534 if (param_.solve_welleq_initially_ && well->isOperableAndSolvable()) {
2536 well->solveWellEquation(simulator_, this->wellState(), this->groupState(), deferred_logger);
2537 }
catch (
const std::exception& e) {
2538 const std::string msg =
"Compute initial well solution for " + well->name() +
" initially failed. Continue with the previous rates";
2539 deferred_logger.
warning(
"WELL_INITIAL_SOLVE_FAILED", msg);
2544 well->resetWellOperability();
2546 updatePrimaryVariables(deferred_logger);
2549 if (do_prestep_network_rebalance) doPreStepNetworkRebalance(deferred_logger);
2552 template<
typename TypeTag>
2557 std::vector< Scalar > B_avg(numComponents(),
Scalar() );
2558 const auto& grid = simulator_.vanguard().grid();
2559 const auto& gridView = grid.leafGridView();
2563 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
2564 elemCtx.updatePrimaryStencil(elem);
2565 elemCtx.updatePrimaryIntensiveQuantities(0);
2567 const auto& intQuants = elemCtx.intensiveQuantities(0, 0);
2568 const auto& fs = intQuants.fluidState();
2570 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx)
2572 if (!FluidSystem::phaseIsActive(phaseIdx)) {
2576 const unsigned compIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
2577 auto& B = B_avg[ compIdx ];
2579 B += 1 / fs.invB(phaseIdx).value();
2581 if constexpr (has_solvent_) {
2582 auto& B = B_avg[solventSaturationIdx];
2583 B += 1 / intQuants.solventInverseFormationVolumeFactor().value();
2589 grid.comm().sum(B_avg.data(), B_avg.size());
2590 for (
auto& bval : B_avg)
2592 bval /= global_num_cells_;
2601 template<
typename TypeTag>
2606 for (
const auto& well : well_container_) {
2607 well->updatePrimaryVariables(simulator_, this->wellState(), deferred_logger);
2611 template<
typename TypeTag>
2615 const auto& grid = simulator_.vanguard().grid();
2616 const auto& eclProblem = simulator_.problem();
2617 const unsigned numCells = grid.size(0);
2619 this->pvt_region_idx_.resize(numCells);
2620 for (
unsigned cellIdx = 0; cellIdx < numCells; ++cellIdx) {
2621 this->pvt_region_idx_[cellIdx] =
2622 eclProblem.pvtRegionIndex(cellIdx);
2627 template<
typename TypeTag>
2637 int numComp = this->numPhases() < 3 ? this->numPhases() : FluidSystem::numComponents;
2638 if constexpr (has_solvent_) {
2644 template<
typename TypeTag>
2648 const auto& eclProblem = simulator_.problem();
2649 depth_.resize(local_num_cells_);
2650 for (
unsigned cellIdx = 0; cellIdx < local_num_cells_; ++cellIdx) {
2651 depth_[cellIdx] = eclProblem.dofCenterDepth(cellIdx);
2655 template<
typename TypeTag>
2658 getWell(
const std::string& well_name)
const
2661 auto well = std::find_if(well_container_.begin(),
2662 well_container_.end(),
2664 return elem->name() == well_name;
2667 assert(well != well_container_.end());
2672 template<
typename TypeTag>
2675 hasWell(
const std::string& well_name)
const
2677 return std::any_of(well_container_.begin(), well_container_.end(),
2680 return elem->name() == well_name;
2687 template <
typename TypeTag>
2692 return std::max(this->simulator_.episodeIndex(), 0);
2699 template<
typename TypeTag>
2704 const std::vector<Scalar>& production_rates,
2705 std::vector<Scalar>& resv_coeff)
2707 rateConverter_->calcCoeff(fipnum, pvtreg, production_rates, resv_coeff);
2710 template<
typename TypeTag>
2715 std::vector<Scalar>& resv_coeff)
2717 rateConverter_->calcInjCoeff(fipnum, pvtreg, resv_coeff);
2721 template <
typename TypeTag>
2729 int np = this->numPhases();
2730 Scalar cellInternalEnergy;
2734 const int nw = this->numLocalWells();
2735 for (
auto wellID = 0*nw; wellID < nw; ++wellID) {
2736 const Well& well = this->wells_ecl_[wellID];
2737 auto& ws = this->wellState().well(wellID);
2738 if (well.isInjector()){
2739 if( !(ws.status == WellStatus::STOP)){
2740 this->wellState().well(wellID).temperature = well.inj_temperature();
2745 std::array<Scalar,2> weighted{0.0,0.0};
2746 auto& [weighted_temperature, total_weight] = weighted;
2748 auto& well_info = this->local_parallel_well_info_[wellID].get();
2749 auto& perf_data = ws.perf_data;
2750 auto& perf_phase_rate = perf_data.phase_rates;
2752 using int_type =
decltype(this->well_perf_data_[wellID].size());
2753 for (int_type perf = 0, end_perf = this->well_perf_data_[wellID].size(); perf < end_perf; ++perf) {
2754 const int cell_idx = this->well_perf_data_[wellID][perf].cell_index;
2755 const auto& intQuants = simulator_.model().intensiveQuantities(cell_idx, 0);
2756 const auto& fs = intQuants.fluidState();
2759 Scalar cellTemperatures = fs.temperature(0).value();
2761 Scalar weight_factor = 0.0;
2762 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx)
2764 if (!FluidSystem::phaseIsActive(phaseIdx)) {
2767 cellInternalEnergy = fs.enthalpy(phaseIdx).value() - fs.pressure(phaseIdx).value() / fs.density(phaseIdx).value();
2768 cellBinv = fs.invB(phaseIdx).value();
2769 cellDensity = fs.density(phaseIdx).value();
2770 perfPhaseRate = perf_phase_rate[ perf*np + phaseIdx ];
2771 weight_factor += cellDensity * perfPhaseRate/cellBinv * cellInternalEnergy/cellTemperatures;
2773 weight_factor = std::abs(weight_factor)+1e-13;
2774 total_weight += weight_factor;
2775 weighted_temperature += weight_factor * cellTemperatures;
2777 well_info.communication().sum(weighted.data(), 2);
2778 this->wellState().well(wellID).temperature = weighted_temperature/total_weight;
2783 template <
typename TypeTag>
2788 std::ostringstream os;
2789 for (
const auto& w : well_container_) {
2790 os << w->name() <<
":";
2791 auto pv = w->getPrimaryVars();
2792 for (
const Scalar v : pv) {
2797 OpmLog::debug(os.str());
2802 template <
typename TypeTag>
2803 std::vector<typename BlackoilWellModel<TypeTag>::Scalar>
2807 std::vector<Scalar> ret;
2808 for (
const auto& well : well_container_) {
2809 if (well_domain_.at(well->name()) == domain.
index) {
2810 const auto& pv = well->getPrimaryVars();
2811 ret.insert(ret.end(), pv.begin(), pv.end());
2819 template <
typename TypeTag>
2824 std::size_t offset = 0;
2825 for (
auto& well : well_container_) {
2826 if (well_domain_.at(well->name()) == domain.
index) {
2827 int num_pri_vars = well->setPrimaryVars(vars.begin() + offset);
2828 offset += num_pri_vars;
2831 assert(offset == vars.size());
2836 template <
typename TypeTag>
2846 for (
const auto& wellPtr : this->well_container_) {
2847 const int first_well_cell = wellPtr->cells().front();
2848 for (
const auto& domain : domains) {
2849 auto cell_present = [&domain](
const auto cell)
2851 return std::binary_search(domain.cells.begin(),
2852 domain.cells.end(), cell);
2855 if (cell_present(first_well_cell)) {
2858 well_domain_[wellPtr->name()] = domain.index;
2861 for (
int well_cell : wellPtr->cells()) {
2862 if (! cell_present(well_cell)) {
2863 OPM_THROW(std::runtime_error,
2864 fmt::format(
"Well '{}' found on multiple domains.",
2872 simulator_.gridView().comm());
2876 const int rank = comm.rank();
2878 if (!well_domain_.empty()) {
2879 std::ostringstream os;
2880 os <<
"Well name Rank Domain\n";
2881 for (
const auto& [wname, domain] : well_domain_) {
2882 os << wname << std::setw(19 - wname.size()) << rank << std::setw(12) << domain <<
'\n';
2884 local_log.
debug(os.str());
2887 if (this->terminal_output_) {
2888 global_log.logMessages();
#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 guide rates in the blackoil well model.
Definition: BlackoilWellModelGuideRates.hpp:46
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:98
void initializeGroupStructure(const int reportStepIdx)
Definition: BlackoilWellModel_impl.hpp:371
void apply(BVector &r) const
Definition: BlackoilWellModel_impl.hpp:1713
void prepareTimeStep(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2509
bool updateGroupControls(const Group &group, DeferredLogger &deferred_logger, const int reportStepIdx, const int iterationIdx)
Definition: BlackoilWellModel_impl.hpp:2337
WellInterfacePtr createWellPointer(const int wellID, const int report_step) const
Definition: BlackoilWellModel_impl.hpp:1013
WellInterfacePtr getWell(const std::string &well_name) const
Definition: BlackoilWellModel_impl.hpp:2658
void calcInjResvCoeff(const int fipnum, const int pvtreg, std::vector< Scalar > &resv_coeff) override
Definition: BlackoilWellModel_impl.hpp:2713
void init()
Definition: BlackoilWellModel_impl.hpp:154
typename BlackoilWellModelGeneric< Scalar >::GLiftProdWells GLiftProdWells
Definition: BlackoilWellModel.hpp:113
void calcResvCoeff(const int fipnum, const int pvtreg, const std::vector< Scalar > &production_rates, std::vector< Scalar > &resv_coeff) override
Definition: BlackoilWellModel_impl.hpp:2702
void updateAndCommunicate(const int reportStepIdx, const int iterationIdx, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2315
void doPreStepNetworkRebalance(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1091
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: BlackoilWellModel.hpp:107
std::vector< Scalar > getPrimaryVarsDomain(const Domain &domain) const
Definition: BlackoilWellModel_impl.hpp:2805
typename GasLiftGroupInfo< Scalar >::GLiftEclWells GLiftEclWells
Definition: BlackoilWellModel.hpp:116
void initWellContainer(const int reportStepIdx) override
Definition: BlackoilWellModel_impl.hpp:173
void addReservoirSourceTerms(GlobalEqVector &residual, std::vector< typename SparseMatrixAdapter::MatrixBlock * > &diagMatAddress) const
Definition: BlackoilWellModel_impl.hpp:1820
void beginReportStep(const int time_step)
Definition: BlackoilWellModel_impl.hpp:274
void assembleWellEqDomain(const double dt, const Domain &domain, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1670
void applyScaleAdd(const Scalar alpha, const BVector &x, BVector &Ax) const
Definition: BlackoilWellModel_impl.hpp:1774
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: BlackoilWellModel.hpp:104
GetPropType< TypeTag, Properties::Grid > Grid
Definition: BlackoilWellModel.hpp:101
void updateWellControlsDomain(DeferredLogger &deferred_logger, const Domain &domain)
Definition: BlackoilWellModel_impl.hpp:2146
void assembleWellEq(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1659
void gasLiftOptimizationStage1SingleWell(WellInterface< TypeTag > *well, DeferredLogger &deferred_logger, GLiftProdWells &prod_wells, GLiftOptWells &glift_wells, GasLiftGroupInfo< Scalar > &group_info, GLiftWellStateMap &state_map, GLiftSyncGroups &groups_to_sync)
Definition: BlackoilWellModel_impl.hpp:1619
void linearizeDomain(const Domain &domain, SparseMatrixAdapter &jacobian, GlobalEqVector &res)
Definition: BlackoilWellModel_impl.hpp:251
int reportStepIndex() const
Definition: BlackoilWellModel_impl.hpp:2690
void calculateProductivityIndexValues(DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2450
void extractLegacyDepth_()
Definition: BlackoilWellModel_impl.hpp:2646
void extractLegacyCellPvtRegionIndex_()
Definition: BlackoilWellModel_impl.hpp:2613
void updateAverageFormationFactor()
Definition: BlackoilWellModel_impl.hpp:2555
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: BlackoilWellModel.hpp:106
void initializeWellState(const int timeStepIdx)
Definition: BlackoilWellModel_impl.hpp:829
void assemble(const int iterationIdx, const double dt)
Definition: BlackoilWellModel_impl.hpp:1129
void computeWellTemperature()
Definition: BlackoilWellModel_impl.hpp:2724
void addWellPressureEquations(PressureMatrix &jacobian, const BVector &weights, const bool use_well_weights) const
Definition: BlackoilWellModel_impl.hpp:1804
void logPrimaryVars() const
Definition: BlackoilWellModel_impl.hpp:2786
const SimulatorReportSingle & lastReport() const
Definition: BlackoilWellModel_impl.hpp:680
bool updateWellControlsAndNetwork(const bool mandatory_network_balance, const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1189
void addWellContributions(SparseMatrixAdapter &jacobian) const
Definition: BlackoilWellModel_impl.hpp:1794
WellInterfacePtr createWellForWellTest(const std::string &well_name, const int report_step, DeferredLogger &deferred_logger) const
Definition: BlackoilWellModel_impl.hpp:1064
typename BlackoilWellModelGeneric< Scalar >::GLiftOptWells GLiftOptWells
Definition: BlackoilWellModel.hpp:112
void computePotentials(const std::size_t widx, const WellState< Scalar > &well_state_copy, std::string &exc_msg, ExceptionType::ExcEnum &exc_type, DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2415
std::pair< bool, bool > updateWellControls(const bool mandatory_network_balance, DeferredLogger &deferred_logger, const bool relax_network_tolerance=false)
Definition: BlackoilWellModel_impl.hpp:2055
void recoverWellSolutionAndUpdateWellStateDomain(const BVector &x, const Domain &domain)
Definition: BlackoilWellModel_impl.hpp:1891
typename GasLiftSingleWellGeneric< Scalar >::GLiftSyncGroups GLiftSyncGroups
Definition: BlackoilWellModel.hpp:117
void assembleDomain(const int iterationIdx, const double dt, const Domain &domain)
Definition: BlackoilWellModel_impl.hpp:1413
bool hasWell(const std::string &well_name) const
Definition: BlackoilWellModel_impl.hpp:2675
ParallelWBPCalculation< Scalar >::EvaluatorFactory makeWellSourceEvaluatorFactory(const std::vector< Well >::size_type wellIdx) const
Definition: BlackoilWellModel_impl.hpp:2243
void assembleWellEqWithoutIteration(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1694
data::WellBlockAveragePressures computeWellBlockAveragePressures() const
Definition: BlackoilWellModel_impl.hpp:2197
void initPrimaryVariablesEvaluation() const
Definition: BlackoilWellModel_impl.hpp:1915
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: BlackoilWellModel.hpp:360
void computeTotalRatesForDof(RateVector &rate, unsigned globalIdx) const
Definition: BlackoilWellModel_impl.hpp:792
void beginTimeStep()
Definition: BlackoilWellModel_impl.hpp:413
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: BlackoilWellModel.hpp:108
void calculateExplicitQuantities(DeferredLogger &deferred_logger) const
Definition: BlackoilWellModel_impl.hpp:2040
Dune::FieldVector< Scalar, numEq > VectorBlockType
Definition: BlackoilWellModel.hpp:133
void updatePrimaryVariables(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:2604
void initializeLocalWellStructure(const int reportStepIdx, const bool enableWellPIScaling)
Definition: BlackoilWellModel_impl.hpp:326
std::pair< bool, bool > updateWellControlsAndNetworkIteration(const bool mandatory_network_balance, const bool relax_network_tolerance, const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1225
void addNeighbors(std::vector< NeighborSet > &neighbors) const override
Specify the additional neighboring correlations caused by the auxiliary module.
Definition: BlackoilWellModel_impl.hpp:189
BlackoilWellModel(Simulator &simulator)
Definition: BlackoilWellModel_impl.hpp:146
void wellTesting(const int timeStepIdx, const double simulationTime, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:613
ConvergenceReport getWellConvergence(const std::vector< Scalar > &B_avg, const bool checkWellGroupControls=false) const
Definition: BlackoilWellModel_impl.hpp:1988
void registerOpenWellsForWBPCalculation()
Definition: BlackoilWellModel_impl.hpp:2294
void linearize(SparseMatrixAdapter &jacobian, GlobalEqVector &res) override
Linearize the auxiliary equation.
Definition: BlackoilWellModel_impl.hpp:230
void initGliftEclWellMap(GLiftEclWells &ecl_well_map)
Definition: BlackoilWellModel_impl.hpp:1647
void gasLiftOptimizationStage1(DeferredLogger &deferred_logger, GLiftProdWells &prod_wells, GLiftOptWells &glift_wells, GasLiftGroupInfo< Scalar > &group_info, GLiftWellStateMap &state_map)
Definition: BlackoilWellModel_impl.hpp:1511
void timeStepSucceeded(const double simulationTime, const double dt)
Definition: BlackoilWellModel_impl.hpp:690
std::unique_ptr< WellType > createTypedWellPointer(const int wellID, const int time_step) const
Definition: BlackoilWellModel_impl.hpp:1033
void setPrimaryVarsDomain(const Domain &domain, const std::vector< Scalar > &vars)
Definition: BlackoilWellModel_impl.hpp:2822
std::shared_ptr< WellInterface< TypeTag > > WellInterfacePtr
Definition: BlackoilWellModel.hpp:238
void createWellContainer(const int report_step) override
Definition: BlackoilWellModel_impl.hpp:872
Dune::BlockVector< VectorBlockType > BVector
Definition: BlackoilWellModel.hpp:134
bool maybeDoGasLiftOptimize(DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1453
void computeWellGroupThp(const double dt, DeferredLogger &local_deferredLogger)
Definition: BlackoilWellModel_impl.hpp:1278
void updateWellTestState(const double &simulationTime, WellTestState &wellTestState) const
upate the wellTestState related to economic limits
Definition: BlackoilWellModel_impl.hpp:2375
ConvergenceReport getDomainWellConvergence(const Domain &domain, const std::vector< Scalar > &B_avg, DeferredLogger &local_deferredLogger) const
Definition: BlackoilWellModel_impl.hpp:1943
void initializeWBPCalculationService()
Definition: BlackoilWellModel_impl.hpp:2169
void initPrimaryVariablesEvaluationDomain(const Domain &domain) const
Definition: BlackoilWellModel_impl.hpp:1926
void recoverWellSolutionAndUpdateWellState(const BVector &x)
Definition: BlackoilWellModel_impl.hpp:1873
void setupDomains(const std::vector< Domain > &domains)
Definition: BlackoilWellModel_impl.hpp:2839
void addWellPressureEquationsStruct(PressureMatrix &jacobian) const
Definition: BlackoilWellModel_impl.hpp:1850
void calculateProductivityIndexValuesShutWells(const int reportStepIdx, DeferredLogger &deferred_logger) override
Definition: BlackoilWellModel_impl.hpp:2464
void endReportStep()
Definition: BlackoilWellModel_impl.hpp:663
void prepareWellsBeforeAssembling(const double dt, DeferredLogger &deferred_logger)
Definition: BlackoilWellModel_impl.hpp:1683
typename BlackoilWellModelGeneric< Scalar >::GLiftWellStateMap GLiftWellStateMap
Definition: BlackoilWellModel.hpp:115
int numComponents() const
Definition: BlackoilWellModel_impl.hpp:2629
Definition: ConvergenceReport.hpp:38
void setWellFailed(const WellFailure &wf)
Definition: ConvergenceReport.hpp:197
void setWellGroupTargetsViolated(const bool wellGroupTargetsViolated)
Definition: ConvergenceReport.hpp:209
const std::vector< WellFailure > & wellFailures() const
Definition: ConvergenceReport.hpp:289
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)
Definition: GasLiftGroupInfo.hpp:46
std::tuple< Scalar, Scalar, Scalar, Scalar > getRates(const int group_idx) const
bool hasWell(const std::string &well_name)
void updateRate(int idx, Scalar oil_rate, Scalar gas_rate, Scalar water_rate, Scalar alq)
Class for serializing and broadcasting data using MPI.
Definition: MPISerializer.hpp:31
void broadcast(T &data, int root=0)
Serialize and broadcast on root process, de-serialize on others.
Definition: MPISerializer.hpp:46
Definition: ParallelWBPCalculation.hpp:51
typename ParallelPAvgDynamicSourceData< Scalar >::Evaluator Evaluator
Callback for evaluating WBPn source terms on the current MPI rank.
Definition: ParallelWBPCalculation.hpp:58
Definition: StandardWell.hpp:59
virtual void init(const PhaseUsage *phase_usage_arg, const std::vector< Scalar > &depth_arg, const Scalar gravity_arg, const int num_cells, const std::vector< Scalar > &B_avg, const bool changed_to_open_this_step) override
Definition: StandardWell_impl.hpp:96
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
static bool bruteForceBracketCommonTHP(const std::function< Scalar(const Scalar)> &eq, const std::array< Scalar, 2 > &range, Scalar &low, Scalar &high, std::optional< Scalar > &approximate_solution, const Scalar &limit, DeferredLogger &deferred_logger)
Find limits using brute-force solver.
Definition: WellContributions.hpp:53
void alloc()
Allocate memory for the StandardWells.
void setBlockSize(unsigned int dim, unsigned int dim_wells)
void addNumBlocks(unsigned int numBlocks)
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 > > ®ionalAveragePressureCalculator)
static void updateGuideRates(const Group &group, const Schedule &schedule, const SummaryState &summary_state, const PhaseUsage &pu, int report_step, double sim_time, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, const Parallel::Communication &comm, GuideRate *guide_rate, std::vector< Scalar > &pot, DeferredLogger &deferred_logger)
static void updateGpMaintTargetForGroups(const Group &group, const Schedule &schedule, const RegionalValues ®ional_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)
const std::string & name() const
Well name.
int indexOfWell() const
Index of well in the wells struct and wellState.
Definition: WellInterface.hpp:73
bool updateWellControl(const Simulator &simulator, const IndividualOrGroup iog, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:195
virtual void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const =0
Definition: WellState.hpp:62
ExcEnum
Definition: DeferredLogger.hpp:45
@ NONE
Definition: DeferredLogger.hpp:46
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
Definition: blackoilboundaryratevector.hh:37
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
Definition: SubDomain.hpp:62
int index
Definition: SubDomain.hpp:65