23#ifndef OPM_WELLINTERFACE_HEADER_INCLUDED
25#define OPM_WELLINTERFACE_IMPL_HEADER_INCLUDED
29#include <opm/common/Exceptions.hpp>
31#include <opm/input/eclipse/Schedule/ScheduleTypes.hpp>
32#include <opm/input/eclipse/Schedule/Well/WDFAC.hpp>
41#include <dune/common/version.hh>
46#include <fmt/format.h>
52 template<
typename TypeTag>
59 const int pvtRegionIdx,
60 const int num_components,
62 const int index_of_well,
63 const std::vector<PerforationData>& perf_data)
78 if (well.isInjector()) {
79 auto injectorType = this->
well_ecl_.injectorType();
80 if (injectorType == InjectorType::GAS) {
88 template<
typename TypeTag>
92 const std::vector<Scalar>& ,
95 const std::vector<Scalar>& B_avg,
96 const bool changed_to_open_this_step)
98 this->phase_usage_ = phase_usage_arg;
99 this->gravity_ = gravity_arg;
101 this->changed_to_open_this_step_ = changed_to_open_this_step;
107 template<
typename TypeTag>
112 if constexpr (has_polymer) {
113 return this->wpolymer_();
123 template<
typename TypeTag>
128 if constexpr (has_foam) {
129 return this->wfoam_();
137 template<
typename TypeTag>
142 if constexpr (has_brine) {
143 return this->wsalt_();
149 template<
typename TypeTag>
154 if constexpr (has_micp) {
155 return this->wmicrobes_();
161 template<
typename TypeTag>
166 if constexpr (has_micp) {
167 return this->woxygen_();
179 template<
typename TypeTag>
184 if constexpr (has_micp) {
185 return this->wurea_();
191 template<
typename TypeTag>
200 const auto& summary_state = simulator.vanguard().summaryState();
201 if (this->stopppedOrZeroRateTarget(summary_state, well_state)) {
205 const auto& summaryState = simulator.vanguard().summaryState();
206 const auto& schedule = simulator.vanguard().schedule();
207 const auto& well = this->well_ecl_;
208 auto& ws = well_state.
well(this->index_of_well_);
210 if (well.isInjector()) {
211 from = WellInjectorCMode2String(ws.injection_cmode);
213 from = WellProducerCMode2String(ws.production_cmode);
215 bool oscillating = std::count(this->well_control_log_.begin(), this->well_control_log_.end(), from) >= param_.max_number_of_well_switches_;
219 bool output = std::count(this->well_control_log_.begin(), this->well_control_log_.end(), from) == param_.max_number_of_well_switches_;
221 std::ostringstream ss;
222 ss <<
" The control mode for well " << this->name()
223 <<
" is oscillating\n"
224 <<
" We don't allow for more than "
225 << param_.max_number_of_well_switches_
226 <<
" switches. The control is kept at " << from;
227 deferred_logger.
info(ss.str());
229 this->well_control_log_.push_back(from);
233 bool changed =
false;
234 if (iog == IndividualOrGroup::Individual) {
235 changed = this->checkIndividualConstraints(ws, summaryState, deferred_logger);
236 }
else if (iog == IndividualOrGroup::Group) {
237 changed = this->checkGroupConstraints(well_state, group_state, schedule, summaryState, deferred_logger);
239 assert(iog == IndividualOrGroup::Both);
240 changed = this->checkConstraints(well_state, group_state, schedule, summaryState, deferred_logger);
246 if (well.isInjector()) {
247 to = WellInjectorCMode2String(ws.injection_cmode);
249 to = WellProducerCMode2String(ws.production_cmode);
251 std::ostringstream ss;
252 ss <<
" Switching control mode for well " << this->name()
256 ss <<
" on rank " << cc.rank();
258 deferred_logger.
debug(ss.str());
260 this->well_control_log_.push_back(from);
261 updateWellStateWithTarget(simulator, group_state, well_state, deferred_logger);
262 updatePrimaryVariables(summaryState, well_state, deferred_logger);
268 template<
typename TypeTag>
274 const Well::InjectionControls& inj_controls,
275 const Well::ProductionControls& prod_controls,
278 const bool fixed_control,
279 const bool fixed_status)
281 const auto& summary_state = simulator.vanguard().summaryState();
282 const auto& schedule = simulator.vanguard().schedule();
284 if (this->wellUnderZeroRateTarget(summary_state, well_state) || !(this->well_ecl_.getStatus() == WellStatus::OPEN)) {
288 const Scalar sgn = this->isInjector() ? 1.0 : -1.0;
289 if (!this->wellIsStopped()){
290 if (wqTotal*sgn <= 0.0 && !fixed_status){
294 bool changed =
false;
295 if (!fixed_control) {
296 auto& ws = well_state.
well(this->index_of_well_);
297 const bool hasGroupControl = this->isInjector() ? inj_controls.hasControl(Well::InjectorCMode::GRUP) :
298 prod_controls.hasControl(Well::ProducerCMode::GRUP);
300 changed = this->checkIndividualConstraints(ws, summary_state, deferred_logger, inj_controls, prod_controls);
301 if (hasGroupControl) {
302 changed = changed || this->checkGroupConstraints(well_state, group_state, schedule, summary_state,deferred_logger);
306 const bool thp_controlled = this->isInjector() ? ws.injection_cmode == Well::InjectorCMode::THP :
307 ws.production_cmode == Well::ProducerCMode::THP;
308 if (!thp_controlled){
310 updateWellStateWithTarget(simulator, group_state, well_state, deferred_logger);
312 ws.thp = this->getTHPConstraint(summary_state);
314 updatePrimaryVariables(summary_state, well_state, deferred_logger);
319 }
else if (!fixed_status){
321 const Scalar bhp = well_state.
well(this->index_of_well_).bhp;
322 Scalar prod_limit = prod_controls.bhp_limit;
323 Scalar inj_limit = inj_controls.bhp_limit;
324 const bool has_thp = this->wellHasTHPConstraints(summary_state);
326 std::vector<Scalar> rates(this->num_components_);
327 if (this->isInjector()){
329 calculateBhpFromThp(well_state, rates,
332 this->getRefDensity(),
334 inj_limit = std::min(bhp_thp,
static_cast<Scalar>(inj_controls.bhp_limit));
339 calculateMinimumBhpFromThp(well_state,
342 this->getRefDensity());
343 prod_limit = std::max(bhp_min,
static_cast<Scalar>(prod_controls.bhp_limit));
346 const Scalar bhp_diff = (this->isInjector())? inj_limit -
bhp:
bhp - prod_limit;
349 well_state.
well(this->index_of_well_).bhp = (this->isInjector())? inj_limit : prod_limit;
351 well_state.
well(this->index_of_well_).thp = this->getTHPConstraint(summary_state);
362 template<
typename TypeTag>
366 const double simulation_time,
369 WellTestState& well_test_state,
372 deferred_logger.
info(
" well " + this->name() +
" is being tested");
375 auto& ws = well_state_copy.
well(this->indexOfWell());
377 updateWellStateWithTarget(simulator, group_state, well_state_copy, deferred_logger);
378 calculateExplicitQuantities(simulator, well_state_copy, deferred_logger);
379 const auto& summary_state = simulator.vanguard().summaryState();
380 updatePrimaryVariables(summary_state, well_state_copy, deferred_logger);
381 initPrimaryVariablesEvaluation();
383 if (this->isProducer()) {
384 const auto& schedule = simulator.vanguard().schedule();
385 const auto report_step = simulator.episodeIndex();
386 const auto& glo = schedule.glo(report_step);
388 gliftBeginTimeStepWellTestUpdateALQ(simulator, well_state_copy, deferred_logger);
392 WellTestState welltest_state_temp;
394 bool testWell =
true;
399 const std::size_t original_number_closed_completions = welltest_state_temp.num_closed_completions();
400 bool converged = solveWellForTesting(simulator, well_state_copy, group_state, deferred_logger);
402 const auto msg = fmt::format(
"WTEST: Well {} is not solvable (physical)", this->name());
403 deferred_logger.
debug(msg);
408 updateWellOperability(simulator, well_state_copy, deferred_logger);
409 if ( !this->isOperableAndSolvable() ) {
410 const auto msg = fmt::format(
"WTEST: Well {} is not operable (physical)", this->name());
411 deferred_logger.
debug(msg);
414 std::vector<Scalar> potentials;
416 computeWellPotentials(simulator, well_state_copy, potentials, deferred_logger);
417 }
catch (
const std::exception& e) {
418 const std::string msg = fmt::format(
"well {}: computeWellPotentials() "
419 "failed during testing for re-opening: ",
420 this->name(), e.what());
421 deferred_logger.
info(msg);
424 const int np = well_state_copy.
numPhases();
425 for (
int p = 0; p < np; ++p) {
426 ws.well_potentials[p] = std::max(
Scalar{0.0}, potentials[p]);
428 this->updateWellTestState(well_state_copy.
well(this->indexOfWell()),
433 this->closeCompletions(welltest_state_temp);
439 if ( welltest_state_temp.num_closed_wells() > 0 ||
440 (original_number_closed_completions == welltest_state_temp.num_closed_completions()) ) {
446 if (!welltest_state_temp.well_is_closed(this->name())) {
447 well_test_state.open_well(this->name());
449 std::string msg = std::string(
"well ") + this->name() + std::string(
" is re-opened");
450 deferred_logger.
info(msg);
453 for (
auto& completion : this->well_ecl_.getCompletions()) {
454 if (!welltest_state_temp.completion_is_closed(this->name(), completion.first))
455 well_test_state.open_completion(this->name(), completion.first);
459 well_state = well_state_copy;
466 template<
typename TypeTag>
475 const auto& summary_state = simulator.vanguard().summaryState();
476 const auto inj_controls = this->well_ecl_.isInjector() ? this->well_ecl_.injectionControls(summary_state) : Well::InjectionControls(0);
477 const auto prod_controls = this->well_ecl_.isProducer() ? this->well_ecl_.productionControls(summary_state) : Well::ProductionControls(0);
478 bool converged =
false;
481 if (!this->param_.local_well_solver_control_switching_){
482 converged = this->iterateWellEqWithControl(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
484 if (this->param_.use_implicit_ipr_ && this->well_ecl_.isProducer() && this->wellHasTHPConstraints(summary_state) && (this->well_ecl_.getStatus() == WellStatus::OPEN)) {
485 converged = solveWellWithTHPConstraint(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
487 converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
491 }
catch (NumericalProblem& e ) {
492 const std::string msg =
"Inner well iterations failed for well " + this->name() +
" Treat the well as unconverged. ";
493 deferred_logger.
warning(
"INNER_ITERATION_FAILED", msg);
499 template<
typename TypeTag>
504 const Well::InjectionControls& inj_controls,
505 const Well::ProductionControls& prod_controls,
510 const auto& summary_state = simulator.vanguard().summaryState();
511 bool is_operable =
true;
512 bool converged =
true;
513 auto& ws = well_state.
well(this->index_of_well_);
515 if (this->wellIsStopped()) {
517 auto bhp_target = estimateOperableBhp(simulator, dt, well_state, summary_state, deferred_logger);
518 if (!bhp_target.has_value()) {
520 const auto msg = fmt::format(
"estimateOperableBhp: Did not find operable BHP for well {}", this->name());
521 deferred_logger.
debug(msg);
524 solveWellWithZeroRate(simulator, dt, well_state, deferred_logger);
528 ws.thp = this->getTHPConstraint(summary_state);
529 const Scalar bhp = std::max(bhp_target.value(),
530 static_cast<Scalar>(prod_controls.bhp_limit));
531 solveWellWithBhp(simulator, dt,
bhp, well_state, deferred_logger);
536 converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
539 const bool isThp = ws.production_cmode == Well::ProducerCMode::THP;
541 if (converged && !this->stopppedOrZeroRateTarget(summary_state, well_state) && isThp) {
542 auto rates = well_state.
well(this->index_of_well_).surface_rates;
543 this->adaptRatesForVFP(rates);
544 this->updateIPRImplicit(simulator, well_state, deferred_logger);
548 this->operability_status_.use_vfpexplicit =
true;
551 const Scalar reltol = 1e-3;
552 const Scalar cur_bhp = ws.bhp;
553 if (bhp_stable.has_value() && cur_bhp - bhp_stable.value() > cur_bhp*reltol){
554 const auto msg = fmt::format(
"Well {} converged to an unstable solution, re-solving", this->name());
555 deferred_logger.
debug(msg);
556 solveWellWithBhp(simulator, dt, bhp_stable.value(), well_state, deferred_logger);
558 ws.thp = this->getTHPConstraint(summary_state);
559 converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
566 this->operability_status_.use_vfpexplicit =
true;
568 auto bhp_target = estimateOperableBhp(simulator, dt, well_state, summary_state, deferred_logger);
569 if (!bhp_target.has_value()) {
573 converged = solveWellWithZeroRate(simulator, dt, well_state, deferred_logger);
577 const Scalar bhp = std::max(bhp_target.value(),
578 static_cast<Scalar>(prod_controls.bhp_limit));
579 solveWellWithBhp(simulator, dt,
bhp, well_state, deferred_logger);
580 ws.thp = this->getTHPConstraint(summary_state);
581 converged = this->iterateWellEqWithSwitching(simulator, dt,
590 is_operable = is_operable && !this->wellIsStopped();
591 this->operability_status_.can_obtain_bhp_with_thp_limit = is_operable;
592 this->operability_status_.obey_thp_limit_under_bhp_limit = is_operable;
596 template<
typename TypeTag>
597 std::optional<typename WellInterface<TypeTag>::Scalar>
602 const SummaryState& summary_state,
609 const bool converged = solveWellWithBhp(simulator, dt, bhp_min, well_state, deferred_logger);
610 if (!converged || this->wellIsStopped()) {
613 this->updateIPRImplicit(simulator, well_state, deferred_logger);
614 auto rates = well_state.
well(this->index_of_well_).surface_rates;
615 this->adaptRatesForVFP(rates);
619 template<
typename TypeTag>
630 auto inj_controls = Well::InjectionControls(0);
631 auto prod_controls = Well::ProductionControls(0);
632 auto& ws = well_state.
well(this->index_of_well_);
633 auto cmode_inj = ws.injection_cmode;
634 auto cmode_prod = ws.production_cmode;
635 if (this->isInjector()) {
636 inj_controls.addControl(Well::InjectorCMode::BHP);
637 inj_controls.bhp_limit =
bhp;
638 inj_controls.cmode = Well::InjectorCMode::BHP;
639 ws.injection_cmode = Well::InjectorCMode::BHP;
641 prod_controls.addControl(Well::ProducerCMode::BHP);
642 prod_controls.bhp_limit =
bhp;
643 prod_controls.cmode = Well::ProducerCMode::BHP;
644 ws.production_cmode = Well::ProducerCMode::BHP;
649 const bool converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger,
true);
650 ws.injection_cmode = cmode_inj;
651 ws.production_cmode = cmode_prod;
655 template<
typename TypeTag>
664 const auto well_status_orig = this->wellStatus_;
668 auto inj_controls = Well::InjectionControls(0);
669 auto prod_controls = Well::ProductionControls(0);
670 const bool converged = this->iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger,
true,
true);
671 this->wellStatus_ = well_status_orig;
675 template<
typename TypeTag>
685 const double dt = simulator.timeStepSize();
686 const auto& summary_state = simulator.vanguard().summaryState();
687 const bool has_thp_limit = this->wellHasTHPConstraints(summary_state);
690 well_state.
well(this->indexOfWell()).production_cmode = Well::ProducerCMode::THP;
691 converged = gliftBeginTimeStepWellTestIterateWellEquations(
692 simulator, dt, well_state, group_state, deferred_logger);
695 well_state.
well(this->indexOfWell()).production_cmode = Well::ProducerCMode::BHP;
696 converged = iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
699 deferred_logger.
debug(
"WellTest: Well equation for well " + this->name() +
" converged");
702 const int max_iter = param_.max_welleq_iter_;
703 deferred_logger.
debug(
"WellTest: Well equation for well " + this->name() +
" failed converging in "
705 well_state = well_state0;
710 template<
typename TypeTag>
718 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
723 const double dt = simulator.timeStepSize();
724 bool converged = iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
734 auto& ws = well_state.
well(this->indexOfWell());
735 bool thp_control =
false;
736 if (this->well_ecl_.isInjector()) {
737 thp_control = ws.injection_cmode == Well::InjectorCMode::THP;
739 ws.injection_cmode = Well::InjectorCMode::BHP;
740 this->well_control_log_.push_back(WellInjectorCMode2String(Well::InjectorCMode::THP));
743 thp_control = ws.production_cmode == Well::ProducerCMode::THP;
745 ws.production_cmode = Well::ProducerCMode::BHP;
746 this->well_control_log_.push_back(WellProducerCMode2String(Well::ProducerCMode::THP));
750 const std::string msg = std::string(
"The newly opened well ") + this->name()
751 + std::string(
" with THP control did not converge during inner iterations, we try again with bhp control");
752 deferred_logger.
debug(msg);
753 converged = this->iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
758 const int max_iter = param_.max_welleq_iter_;
759 deferred_logger.
debug(
"Compute initial well solution for well " + this->name() +
". Failed to converge in "
761 well_state = well_state0;
767 template <
typename TypeTag>
776 prepareWellBeforeAssembling(simulator, dt, well_state, group_state, deferred_logger);
777 assembleWellEqWithoutIteration(simulator, dt, well_state, group_state, deferred_logger);
782 template <
typename TypeTag>
791 const auto& summary_state = simulator.vanguard().summaryState();
792 const auto inj_controls = this->well_ecl_.isInjector() ? this->well_ecl_.injectionControls(summary_state) : Well::InjectionControls(0);
793 const auto prod_controls = this->well_ecl_.isProducer() ? this->well_ecl_.productionControls(summary_state) : Well::ProductionControls(0);
796 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
801 template<
typename TypeTag>
810 const bool old_well_operable = this->operability_status_.isOperableAndSolvable();
812 if (param_.check_well_operability_iter_)
813 checkWellOperability(simulator, well_state, deferred_logger);
816 const int iteration_idx = simulator.model().newtonMethod().numIterations();
817 if (iteration_idx < param_.max_niter_inner_well_iter_ || this->well_ecl_.isMultiSegment()) {
818 this->operability_status_.solvable =
true;
819 bool converged = this->iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
823 if (param_.shut_unsolvable_wells_)
824 this->operability_status_.solvable =
false;
827 if (this->operability_status_.has_negative_potentials) {
828 auto well_state_copy = well_state;
829 std::vector<Scalar> potentials;
831 computeWellPotentials(simulator, well_state_copy, potentials, deferred_logger);
832 }
catch (
const std::exception& e) {
833 const std::string msg = fmt::format(
"well {}: computeWellPotentials() failed "
834 "during attempt to recompute potentials for well: ",
835 this->name(), e.what());
836 deferred_logger.
info(msg);
837 this->operability_status_.has_negative_potentials =
true;
839 auto& ws = well_state.
well(this->indexOfWell());
841 for (
int p = 0; p < np; ++p) {
842 ws.well_potentials[p] = std::max(
Scalar{0.0}, potentials[p]);
845 this->changed_to_open_this_step_ =
false;
846 const bool well_operable = this->operability_status_.isOperableAndSolvable();
848 if (!well_operable && old_well_operable) {
849 deferred_logger.
info(
" well " + this->name() +
" gets STOPPED during iteration ");
851 changed_to_stopped_this_step_ =
true;
852 }
else if (well_operable && !old_well_operable) {
853 deferred_logger.
info(
" well " + this->name() +
" gets REVIVED during iteration ");
855 changed_to_stopped_this_step_ =
false;
856 this->changed_to_open_this_step_ =
true;
860 template<
typename TypeTag>
864 if(!this->isOperableAndSolvable() && !this->wellIsStopped())
867 for (
int perfIdx = 0; perfIdx < this->number_of_perforations_; ++perfIdx) {
868 if (this->cells()[perfIdx] == cellIdx) {
869 for (
int i = 0; i < RateVector::dimension; ++i) {
870 rates[i] += connectionRates_[perfIdx][i];
876 template<
typename TypeTag>
880 for (
int perfIdx = 0; perfIdx < this->number_of_perforations_; ++perfIdx) {
881 if (this->cells()[perfIdx] == cellIdx) {
882 const unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
883 return connectionRates_[perfIdx][activeCompIdx].value();
887 OPM_THROW(std::invalid_argument,
"The well with name " + this->name()
895 template<
typename TypeTag>
903 if (!param_.check_well_operability_) {
907 if (this->wellIsStopped() && !changed_to_stopped_this_step_) {
911 updateWellOperability(simulator, well_state, deferred_logger);
912 if (!this->operability_status_.isOperableAndSolvable()) {
913 this->operability_status_.use_vfpexplicit =
true;
914 deferred_logger.
debug(
"EXPLICIT_LOOKUP_VFP",
915 "well not operable, trying with explicit vfp lookup: " + this->name());
916 updateWellOperability(simulator, well_state, deferred_logger);
920 template<
typename TypeTag>
929 const auto& well_name = this->name();
930 assert(this->wellHasTHPConstraints(simulator.vanguard().summaryState()));
931 const auto& schedule = simulator.vanguard().schedule();
932 auto report_step_idx = simulator.episodeIndex();
933 const auto& glo = schedule.glo(report_step_idx);
934 if(glo.active() && glo.has_well(well_name)) {
935 const auto increment = glo.gaslift_increment();
936 auto alq = well_state.
getALQ(well_name);
939 well_state.
setALQ(well_name, alq);
941 iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger)))
950 return iterateWellEquations(simulator, dt, well_state, group_state, deferred_logger);
954 template<
typename TypeTag>
961 const auto& summary_state = simulator.vanguard().summaryState();
962 const auto& well_name = this->name();
963 if (!this->wellHasTHPConstraints(summary_state)) {
964 const std::string msg = fmt::format(
"GLIFT WTEST: Well {} does not have THP constraints", well_name);
965 deferred_logger.
info(msg);
968 const auto& schedule = simulator.vanguard().schedule();
969 const auto report_step_idx = simulator.episodeIndex();
970 const auto& glo = schedule.glo(report_step_idx);
971 if (!glo.has_well(well_name)) {
972 const std::string msg = fmt::format(
973 "GLIFT WTEST: Well {} : Gas Lift not activated: "
974 "WLIFTOPT is probably missing. Skipping.", well_name);
975 deferred_logger.
info(msg);
978 const auto& gl_well = glo.well(well_name);
979 auto& max_alq_optional = gl_well.max_rate();
981 if (max_alq_optional) {
982 max_alq = *max_alq_optional;
985 const auto& well_ecl = this->wellEcl();
986 const auto& controls = well_ecl.productionControls(summary_state);
987 const auto& table = this->vfpProperties()->getProd()->getTable(controls.vfp_table_number);
988 const auto& alq_values = table.getALQAxis();
989 max_alq = alq_values.back();
991 well_state.
setALQ(well_name, max_alq);
992 const std::string msg = fmt::format(
993 "GLIFT WTEST: Well {} : Setting ALQ to max value: {}",
995 deferred_logger.
info(msg);
998 template<
typename TypeTag>
1005 if (this->param_.local_well_solver_control_switching_) {
1006 const bool success = updateWellOperabilityFromWellEq(simulator, well_state, deferred_logger);
1010 deferred_logger.
debug(
"Operability check using well equations did not converge for well "
1011 + this->name() +
", reverting to classical approach." );
1014 this->operability_status_.resetOperability();
1016 bool thp_controlled = this->isInjector() ? well_state.
well(this->index_of_well_).injection_cmode == Well::InjectorCMode::THP:
1017 well_state.
well(this->index_of_well_).production_cmode == Well::ProducerCMode::THP;
1018 bool bhp_controlled = this->isInjector() ? well_state.
well(this->index_of_well_).injection_cmode == Well::InjectorCMode::BHP:
1019 well_state.
well(this->index_of_well_).production_cmode == Well::ProducerCMode::BHP;
1023 bool check_thp = thp_controlled || this->operability_status_.thp_limit_violated_but_not_switched;
1024 if (check_thp || bhp_controlled) {
1025 updateIPR(simulator, deferred_logger);
1026 checkOperabilityUnderBHPLimit(well_state, simulator, deferred_logger);
1030 checkOperabilityUnderTHPLimit(simulator, well_state, deferred_logger);
1034 template<
typename TypeTag>
1042 assert(this->param_.local_well_solver_control_switching_);
1043 this->operability_status_.resetOperability();
1045 const auto& group_state = simulator.problem().wellModel().groupState();
1046 const double dt = simulator.timeStepSize();
1048 bool converged = iterateWellEquations(simulator, dt, well_state_copy, group_state, deferred_logger);
1052 template<
typename TypeTag>
1062 const auto& well = this->well_ecl_;
1063 const int well_index = this->index_of_well_;
1064 auto& ws = well_state.
well(well_index);
1067 const auto& summaryState = simulator.vanguard().summaryState();
1068 const auto& schedule = simulator.vanguard().schedule();
1070 if (this->wellIsStopped()) {
1071 for (
int p = 0; p<np; ++p) {
1072 ws.surface_rates[p] = 0;
1078 if (this->isInjector() )
1080 const auto& controls = well.injectionControls(summaryState);
1082 InjectorType injectorType = controls.injector_type;
1084 switch (injectorType) {
1085 case InjectorType::WATER:
1090 case InjectorType::OIL:
1095 case InjectorType::GAS:
1101 OPM_DEFLOG_THROW(std::runtime_error,
"Expected WATER, OIL or GAS as type for injectors " + this->name(), deferred_logger );
1104 const auto current = ws.injection_cmode;
1107 case Well::InjectorCMode::RATE:
1109 ws.surface_rates[phasePos] = (1.0 - this->rsRvInj()) * controls.surface_rate;
1110 if(this->rsRvInj() > 0) {
1111 if (injectorType == InjectorType::OIL && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1113 }
else if (injectorType == InjectorType::GAS && FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
1116 OPM_DEFLOG_THROW(std::runtime_error,
"Expected OIL or GAS as type for injectors when RS/RV (item 10) is non-zero " + this->name(), deferred_logger );
1122 case Well::InjectorCMode::RESV:
1124 std::vector<Scalar> convert_coeff(this->number_of_phases_, 1.0);
1125 this->rateConverter_.calcCoeff( 0, this->pvtRegionIdx_, convert_coeff);
1126 const Scalar coeff = convert_coeff[phasePos];
1127 ws.surface_rates[phasePos] = controls.reservoir_rate/coeff;
1131 case Well::InjectorCMode::THP:
1133 auto rates = ws.surface_rates;
1138 this->getRefDensity(),
1141 ws.thp = this->getTHPConstraint(summaryState);
1146 Scalar total_rate = std::accumulate(rates.begin(), rates.end(), 0.0);
1147 if (total_rate <= 0.0)
1148 ws.surface_rates = ws.well_potentials;
1152 case Well::InjectorCMode::BHP:
1154 ws.bhp = controls.bhp_limit;
1156 for (
int p = 0; p<np; ++p) {
1157 total_rate += ws.surface_rates[p];
1162 if (total_rate <= 0.0)
1163 ws.surface_rates = ws.well_potentials;
1167 case Well::InjectorCMode::GRUP:
1169 assert(well.isAvailableForGroupControl());
1170 const auto& group = schedule.getGroup(well.groupName(), this->currentStep());
1171 const Scalar efficiencyFactor = well.getEfficiencyFactor();
1172 std::optional<Scalar> target =
1173 this->getGroupInjectionTargetRate(group,
1182 ws.surface_rates[phasePos] = *target;
1185 case Well::InjectorCMode::CMODE_UNDEFINED:
1187 OPM_DEFLOG_THROW(std::runtime_error,
"Well control must be specified for well " + this->name(), deferred_logger );
1197 ws.surface_rates[phasePos] = std::max(
Scalar{1.e-7}, ws.surface_rates[phasePos]);
1200 ws.bhp = controls.bhp_limit;
1206 const auto current = ws.production_cmode;
1207 const auto& controls = well.productionControls(summaryState);
1209 case Well::ProducerCMode::ORAT:
1211 Scalar current_rate = -ws.surface_rates[ pu.phase_pos[Oil] ];
1214 if (current_rate > 0.0) {
1215 for (
int p = 0; p<np; ++p) {
1216 ws.surface_rates[p] *= controls.oil_rate/current_rate;
1219 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1220 double control_fraction = fractions[pu.phase_pos[Oil]];
1221 if (control_fraction != 0.0) {
1222 for (
int p = 0; p<np; ++p) {
1223 ws.surface_rates[p] = - fractions[p] * controls.oil_rate/control_fraction;
1229 case Well::ProducerCMode::WRAT:
1231 Scalar current_rate = -ws.surface_rates[ pu.phase_pos[Water] ];
1234 if (current_rate > 0.0) {
1235 for (
int p = 0; p<np; ++p) {
1236 ws.surface_rates[p] *= controls.water_rate/current_rate;
1239 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1240 const Scalar control_fraction = fractions[pu.phase_pos[Water]];
1241 if (control_fraction != 0.0) {
1242 for (
int p = 0; p<np; ++p) {
1243 ws.surface_rates[p] = - fractions[p] * controls.water_rate / control_fraction;
1249 case Well::ProducerCMode::GRAT:
1251 Scalar current_rate = -ws.surface_rates[pu.phase_pos[Gas] ];
1254 if (current_rate > 0.0) {
1255 for (
int p = 0; p<np; ++p) {
1256 ws.surface_rates[p] *= controls.gas_rate/current_rate;
1259 const std::vector<Scalar > fractions = initialWellRateFractions(simulator, well_state);
1260 const Scalar control_fraction = fractions[pu.phase_pos[Gas]];
1261 if (control_fraction != 0.0) {
1262 for (
int p = 0; p<np; ++p) {
1263 ws.surface_rates[p] = - fractions[p] * controls.gas_rate / control_fraction;
1271 case Well::ProducerCMode::LRAT:
1273 Scalar current_rate = - ws.surface_rates[ pu.phase_pos[Water] ]
1274 - ws.surface_rates[ pu.phase_pos[Oil] ];
1277 if (current_rate > 0.0) {
1278 for (
int p = 0; p<np; ++p) {
1279 ws.surface_rates[p] *= controls.liquid_rate/current_rate;
1282 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1283 const Scalar control_fraction = fractions[pu.phase_pos[Water]] + fractions[pu.phase_pos[Oil]];
1284 if (control_fraction != 0.0) {
1285 for (
int p = 0; p<np; ++p) {
1286 ws.surface_rates[p] = - fractions[p] * controls.liquid_rate / control_fraction;
1292 case Well::ProducerCMode::CRAT:
1295 fmt::format(
"CRAT control not supported, well {}", this->name()),
1298 case Well::ProducerCMode::RESV:
1300 std::vector<Scalar> convert_coeff(this->number_of_phases_, 1.0);
1301 this->rateConverter_.calcCoeff( 0, this->pvtRegionIdx_, ws.surface_rates, convert_coeff);
1302 Scalar total_res_rate = 0.0;
1303 for (
int p = 0; p<np; ++p) {
1304 total_res_rate -= ws.surface_rates[p] * convert_coeff[p];
1306 if (controls.prediction_mode) {
1309 if (total_res_rate > 0.0) {
1310 for (
int p = 0; p<np; ++p) {
1311 ws.surface_rates[p] *= controls.resv_rate/total_res_rate;
1314 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1315 for (
int p = 0; p<np; ++p) {
1316 ws.surface_rates[p] = - fractions[p] * controls.resv_rate / convert_coeff[p];
1320 std::vector<Scalar> hrates(this->number_of_phases_,0.);
1321 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
1322 hrates[pu.phase_pos[Water]] = controls.water_rate;
1324 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
1325 hrates[pu.phase_pos[Oil]] = controls.oil_rate;
1327 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1328 hrates[pu.phase_pos[Gas]] = controls.gas_rate;
1330 std::vector<Scalar> hrates_resv(this->number_of_phases_,0.);
1331 this->rateConverter_.calcReservoirVoidageRates( 0, this->pvtRegionIdx_, hrates, hrates_resv);
1332 Scalar target = std::accumulate(hrates_resv.begin(), hrates_resv.end(), 0.0);
1335 if (total_res_rate > 0.0) {
1336 for (
int p = 0; p<np; ++p) {
1337 ws.surface_rates[p] *= target/total_res_rate;
1340 const std::vector<Scalar> fractions = initialWellRateFractions(simulator, well_state);
1341 for (
int p = 0; p<np; ++p) {
1342 ws.surface_rates[p] = - fractions[p] * target / convert_coeff[p];
1348 case Well::ProducerCMode::BHP:
1350 ws.bhp = controls.bhp_limit;
1352 for (
int p = 0; p<np; ++p) {
1353 total_rate -= ws.surface_rates[p];
1358 if (total_rate <= 0.0){
1359 for (
int p = 0; p<np; ++p) {
1360 ws.surface_rates[p] = -ws.well_potentials[p];
1365 case Well::ProducerCMode::THP:
1367 const bool update_success = updateWellStateWithTHPTargetProd(simulator, well_state, deferred_logger);
1369 if (!update_success) {
1373 auto rates = ws.surface_rates;
1374 this->adaptRatesForVFP(rates);
1376 well_state, rates, well, summaryState, this->getRefDensity(), deferred_logger);
1378 ws.thp = this->getTHPConstraint(summaryState);
1382 const Scalar total_rate = -std::accumulate(rates.begin(), rates.end(), 0.0);
1383 if (total_rate <= 0.0) {
1384 for (
int p = 0; p < this->number_of_phases_; ++p) {
1385 ws.surface_rates[p] = -ws.well_potentials[p];
1391 case Well::ProducerCMode::GRUP:
1393 assert(well.isAvailableForGroupControl());
1394 const auto& group = schedule.getGroup(well.groupName(), this->currentStep());
1395 const Scalar efficiencyFactor = well.getEfficiencyFactor();
1396 Scalar scale = this->getGroupProductionTargetRate(group,
1406 for (
int p = 0; p<np; ++p) {
1407 ws.surface_rates[p] *= scale;
1409 ws.trivial_target =
false;
1411 ws.trivial_target =
true;
1415 case Well::ProducerCMode::CMODE_UNDEFINED:
1418 OPM_DEFLOG_THROW(std::runtime_error,
"Well control must be specified for well " + this->name() , deferred_logger);
1424 ws.bhp = controls.bhp_limit;
1429 template<
typename TypeTag>
1430 std::vector<typename WellInterface<TypeTag>::Scalar>
1435 const int np = this->number_of_phases_;
1436 std::vector<Scalar> scaling_factor(np);
1437 const auto& ws = well_state.
well(this->index_of_well_);
1439 Scalar total_potentials = 0.0;
1440 for (
int p = 0; p<np; ++p) {
1441 total_potentials += ws.well_potentials[p];
1443 if (total_potentials > 0) {
1444 for (
int p = 0; p<np; ++p) {
1445 scaling_factor[p] = ws.well_potentials[p] / total_potentials;
1447 return scaling_factor;
1452 const int nperf = this->number_of_perforations_;
1453 for (
int perf = 0; perf < nperf; ++perf) {
1454 total_tw += this->well_index_[perf];
1456 for (
int perf = 0; perf < nperf; ++perf) {
1457 const int cell_idx = this->well_cells_[perf];
1458 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
1459 const auto& fs = intQuants.fluidState();
1460 const Scalar well_tw_fraction = this->well_index_[perf] / total_tw;
1461 Scalar total_mobility = 0.0;
1462 for (
int p = 0; p < np; ++p) {
1463 int modelPhaseIdx = this->flowPhaseToModelPhaseIdx(p);
1464 total_mobility += fs.invB(modelPhaseIdx).value() * intQuants.mobility(modelPhaseIdx).value();
1466 for (
int p = 0; p < np; ++p) {
1467 int modelPhaseIdx = this->flowPhaseToModelPhaseIdx(p);
1468 scaling_factor[p] += well_tw_fraction * fs.invB(modelPhaseIdx).value() * intQuants.mobility(modelPhaseIdx).value() / total_mobility;
1471 return scaling_factor;
1476 template <
typename TypeTag>
1485 auto& ws = well_state.
well(this->index_of_well_);
1486 int nonzero_rate_index = -1;
1487 const Scalar floating_point_error_epsilon = 1e-14;
1488 for (
int p = 0; p < this->number_of_phases_; ++p) {
1489 if (std::abs(ws.surface_rates[p]) > floating_point_error_epsilon) {
1490 if (nonzero_rate_index == -1) {
1491 nonzero_rate_index = p;
1500 std::vector<Scalar> well_q_s = computeCurrentWellRates(simulator, deferred_logger);
1502 if (nonzero_rate_index == -1) {
1505 for (
int p = 0; p < this->number_of_phases_; ++p) {
1506 ws.surface_rates[p] = well_q_s[this->flowPhaseToModelCompIdx(p)];
1512 const Scalar initial_nonzero_rate = ws.surface_rates[nonzero_rate_index];
1513 const int comp_idx_nz = this->flowPhaseToModelCompIdx(nonzero_rate_index);
1514 if (std::abs(well_q_s[comp_idx_nz]) > floating_point_error_epsilon) {
1515 for (
int p = 0; p < this->number_of_phases_; ++p) {
1516 if (p != nonzero_rate_index) {
1517 const int comp_idx = this->flowPhaseToModelCompIdx(p);
1518 Scalar& rate = ws.surface_rates[p];
1519 rate = (initial_nonzero_rate / well_q_s[comp_idx_nz]) * (well_q_s[comp_idx]);
1525 template <
typename TypeTag>
1526 std::vector<typename WellInterface<TypeTag>::Scalar>
1536 auto wi = std::vector<Scalar>
1537 (this->num_components_, this->well_index_[perf] * trans_mult);
1539 if constexpr (! Indices::gasEnabled) {
1543 const auto& wdfac = this->well_ecl_.getWDFAC();
1545 if (! wdfac.useDFactor() || (this->well_index_[perf] == 0.0)) {
1549 const Scalar d = this->computeConnectionDFactor(perf, intQuants, ws);
1556 const auto& connection = this->well_ecl_.getConnections()[ws.
perf_data.ecl_index[perf]];
1557 const Scalar Kh = connection.Kh();
1558 const Scalar scaling = 3.141592653589 * Kh * connection.wpimult();
1559 const unsigned gas_comp_idx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
1562 const Scalar cell_pressure = getValue(intQuants.fluidState().pressure(FluidSystem::gasPhaseIdx));
1563 const Scalar drawdown = cell_pressure - connection_pressure;
1564 const Scalar invB = getValue(intQuants.fluidState().invB(FluidSystem::gasPhaseIdx));
1565 const Scalar mob_g = getValue(intQuants.mobility(FluidSystem::gasPhaseIdx)) * invB;
1567 const Scalar b = 2*scaling/wi[gas_comp_idx];
1568 const Scalar c = -2*scaling*mob_g*drawdown;
1570 Scalar consistent_Q = -1.0e20;
1572 const Scalar r2n = b*b + 4*a*c;
1574 const Scalar rn = std::sqrt(r2n);
1575 const Scalar xn1 = (b-rn)*0.5/a;
1579 const Scalar xn2 = (b+rn)*0.5/a;
1580 if (xn2 <= 0 && xn2 > consistent_Q) {
1586 const Scalar r2p = b*b - 4*a*c;
1588 const Scalar rp = std::sqrt(r2p);
1589 const Scalar xp1 = (rp-b)*0.5/a;
1590 if (xp1 > 0 && xp1 < consistent_Q) {
1593 const Scalar xp2 = -(rp+b)*0.5/a;
1594 if (xp2 > 0 && xp2 < consistent_Q) {
1598 wi[gas_comp_idx] = 1.0/(1.0/(trans_mult * this->well_index_[perf]) + (consistent_Q/2 * d / scaling));
1603 template <
typename TypeTag>
1609 if (! this->well_ecl_.getWDFAC().useDFactor()) {
1613 auto& d_factor = ws.
perf_data.connection_d_factor;
1615 for (
int perf = 0; perf < this->number_of_perforations_; ++perf) {
1616 const int cell_idx = this->well_cells_[perf];
1617 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
1619 d_factor[perf] = this->computeConnectionDFactor(perf, intQuants, ws);
1623 template <
typename TypeTag>
1630 auto rhoGS = [regIdx = this->pvtRegionIdx()]() {
1631 return FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, regIdx);
1635 auto gas_visc = [connection_pressure = ws.
perf_data.pressure[perf],
1637 regIdx = this->pvtRegionIdx(), &intQuants]()
1639 const auto rv = getValue(intQuants.fluidState().Rv());
1641 const auto& gasPvt = FluidSystem::gasPvt();
1646 const Scalar rv_sat = gasPvt.saturatedOilVaporizationFactor
1647 (regIdx, temperature, connection_pressure);
1649 if (! (rv < rv_sat)) {
1650 return gasPvt.saturatedViscosity(regIdx, temperature,
1651 connection_pressure);
1654 return gasPvt.viscosity(regIdx, temperature, connection_pressure,
1655 rv, getValue(intQuants.fluidState().Rvw()));
1658 const auto& connection = this->well_ecl_.getConnections()
1661 return this->well_ecl_.getWDFAC().getDFactor(rhoGS, gas_visc, connection);
1665 template <
typename TypeTag>
1671 auto connCF = [&connIx = std::as_const(ws.
perf_data.ecl_index),
1672 &conns = this->well_ecl_.getConnections()]
1675 return conns[connIx[perf]].CF();
1678 auto& tmult = ws.
perf_data.connection_compaction_tmult;
1679 auto& ctf = ws.
perf_data.connection_transmissibility_factor;
1681 for (
int perf = 0; perf < this->number_of_perforations_; ++perf) {
1682 const int cell_idx = this->well_cells_[perf];
1684 const auto& intQuants = simulator.model()
1685 .intensiveQuantities(cell_idx, 0);
1687 tmult[perf] = simulator.problem()
1688 .template wellTransMultiplier<double>(intQuants, cell_idx);
1690 ctf[perf] = connCF(perf) * tmult[perf];
1695 template<
typename TypeTag>
1699 if constexpr (Indices::oilEnabled) {
1700 return fs.pressure(FluidSystem::oilPhaseIdx);
1701 }
else if constexpr (Indices::gasEnabled) {
1702 return fs.pressure(FluidSystem::gasPhaseIdx);
1704 return fs.pressure(FluidSystem::waterPhaseIdx);
1708 template <
typename TypeTag>
1709 template<
class Value,
class Callback>
1714 std::vector<Value>& mob,
1715 Callback& extendEval,
1718 auto relpermArray = []()
1720 if constexpr (std::is_same_v<Value, Scalar>) {
1721 return std::array<Scalar,3>{};
1723 return std::array<Eval,3>{};
1726 const int cell_idx = this->well_cells_[perf];
1727 assert (
int(mob.size()) == this->num_components_);
1728 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
1729 const auto& materialLawManager = simulator.problem().materialLawManager();
1733 const int satid = this->saturation_table_number_[perf] - 1;
1734 const int satid_elem = materialLawManager->satnumRegionIdx(cell_idx);
1735 if (satid == satid_elem) {
1736 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
1737 if (!FluidSystem::phaseIsActive(phaseIdx)) {
1741 const unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
1742 mob[activeCompIdx] = extendEval(intQuants.mobility(phaseIdx));
1744 if constexpr (has_solvent) {
1745 mob[Indices::contiSolventEqIdx] = extendEval(intQuants.solventMobility());
1748 const auto& paramsCell = materialLawManager->connectionMaterialLawParams(satid, cell_idx);
1749 auto relativePerms = relpermArray();
1750 MaterialLaw::relativePermeabilities(relativePerms, paramsCell, intQuants.fluidState());
1753 materialLawManager->connectionMaterialLawParams(satid_elem, cell_idx);
1756 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
1757 if (!FluidSystem::phaseIsActive(phaseIdx)) {
1761 const unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
1762 mob[activeCompIdx] = extendEval(relativePerms[phaseIdx] / intQuants.fluidState().viscosity(phaseIdx));
1766 if constexpr (has_solvent) {
1767 OPM_DEFLOG_THROW(std::runtime_error,
"individual mobility for wells does not work in combination with solvent", deferred_logger);
1771 if (this->isInjector() && !this->inj_fc_multiplier_.empty()) {
1772 const auto perf_ecl_index = this->perforationData()[perf].ecl_index;
1773 const auto& connections = this->well_ecl_.getConnections();
1774 const auto& connection = connections[perf_ecl_index];
1775 if (connection.filterCakeActive()) {
1776 for (
auto& val : mob) {
1777 val *= this->inj_fc_multiplier_[perf];
1784 template<
typename TypeTag>
1791 const auto& summary_state = simulator.vanguard().summaryState();
1793 auto bhp_at_thp_limit = computeBhpAtThpLimitProdWithAlq(
1794 simulator, summary_state, this->getALQ(well_state), deferred_logger);
1795 if (bhp_at_thp_limit) {
1796 std::vector<Scalar> rates(this->number_of_phases_, 0.0);
1797 if (thp_update_iterations) {
1798 computeWellRatesWithBhpIterations(simulator, *bhp_at_thp_limit,
1799 rates, deferred_logger);
1801 computeWellRatesWithBhp(simulator, *bhp_at_thp_limit,
1802 rates, deferred_logger);
1804 auto& ws = well_state.
well(this->name());
1805 ws.surface_rates = rates;
1806 ws.bhp = *bhp_at_thp_limit;
1807 ws.thp = this->getTHPConstraint(summary_state);
1814 template <
typename TypeTag>
1819 const std::vector<Scalar>& mobility,
1823 const int np = this->number_of_phases_;
1824 for (
int p = 0; p < np; ++p) {
1827 const auto connMob =
1828 mobility[this->flowPhaseToModelCompIdx(p)]
1829 * fs.invB(this->flowPhaseToModelPhaseIdx(p)).value();
1831 connPI[p] = connPICalc(connMob);
1834 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) &&
1835 FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx))
1837 const auto io = pu.phase_pos[Oil];
1838 const auto ig = pu.phase_pos[Gas];
1840 const auto vapoil = connPI[ig] * fs.Rv().value();
1841 const auto disgas = connPI[io] * fs.Rs().value();
1843 connPI[io] += vapoil;
1844 connPI[ig] += disgas;
1849 template <
typename TypeTag>
1853 const Phase preferred_phase,
1855 const std::vector<Scalar>& mobility,
1863 if (preferred_phase == Phase::GAS) {
1864 phase_pos = pu.phase_pos[Gas];
1866 else if (preferred_phase == Phase::OIL) {
1867 phase_pos = pu.phase_pos[Oil];
1869 else if (preferred_phase == Phase::WATER) {
1870 phase_pos = pu.phase_pos[Water];
1874 fmt::format(
"Unsupported Injector Type ({}) "
1875 "for well {} during connection I.I. calculation",
1876 static_cast<int>(preferred_phase), this->name()),
1880 const auto mt = std::accumulate(mobility.begin(), mobility.end(), 0.0);
1881 connII[phase_pos] = connIICalc(mt * fs.invB(this->flowPhaseToModelPhaseIdx(phase_pos)).value());
#define OPM_DEFLOG_THROW(Exception, message, deferred_logger)
Definition: DeferredLoggingErrorHelpers.hpp:45
@ Liquid
Definition: BlackoilPhases.hpp:42
@ Aqua
Definition: BlackoilPhases.hpp:42
@ Vapour
Definition: BlackoilPhases.hpp:42
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: GroupState.hpp:35
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:184
Definition: SingleWellState.hpp:41
PerfData< Scalar > perf_data
Definition: SingleWellState.hpp:109
Scalar temperature
Definition: SingleWellState.hpp:89
Class for computing BHP limits.
Definition: WellBhpThpCalculator.hpp:42
double calculateMinimumBhpFromThp(const WellState< double > &well_state, const Well &well, const SummaryState &summaryState, const double rho) const
bool isStableSolution(const WellState< double > &well_state, const Well &well, const std::vector< double > &rates, const SummaryState &summaryState) const
std::optional< double > estimateStableBhp(const WellState< double > &well_state, const Well &well, const std::vector< double > &rates, const double rho, const SummaryState &summaryState) const
EvalWell calculateBhpFromThp(const WellState< double > &well_state, const std::vector< EvalWell > &rates, const Well &well, const SummaryState &summaryState, const double rho, DeferredLogger &deferred_logger) const
int number_of_perforations_
Definition: WellInterfaceGeneric.hpp:305
Well well_ecl_
Definition: WellInterfaceGeneric.hpp:272
const std::vector< FluidSystem::Scalar > & wellIndex() const
Definition: WellInterfaceGeneric.hpp:140
FluidSystem::Scalar wsolvent_
Definition: WellInterfaceGeneric.hpp:347
Definition: WellInterfaceIndices.hpp:34
bool updateWellControl(const Simulator &simulator, const IndividualOrGroup iog, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:194
Scalar woxygen() const
Definition: WellInterface_impl.hpp:164
IndividualOrGroup
Definition: WellInterface.hpp:233
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:82
void getMobility(const Simulator &simulator, const int perf, std::vector< Value > &mob, Callback &extendEval, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1712
void assembleWellEq(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:770
typename WellInterfaceFluidSystem< FluidSystem >::RateConverterType RateConverterType
Definition: WellInterface.hpp:104
BlackOilFluidState< Eval, FluidSystem, has_temperature, has_energy, Indices::compositionSwitchIdx >=0, has_watVapor, has_brine, has_saltPrecip, has_disgas_in_water, Indices::numPhases > FluidState
Definition: WellInterface.hpp:134
Scalar wfoam() const
Definition: WellInterface_impl.hpp:126
std::vector< RateVector > connectionRates_
Definition: WellInterface.hpp:369
void computeConnLevelProdInd(const FluidState &fs, const std::function< Scalar(const Scalar)> &connPICalc, const std::vector< Scalar > &mobility, Scalar *connPI) const
Definition: WellInterface_impl.hpp:1817
bool updateWellStateWithTHPTargetProd(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1787
void addCellRates(RateVector &rates, int cellIdx) const
Definition: WellInterface_impl.hpp:862
Scalar volumetricSurfaceRateForConnection(int cellIdx, int phaseIdx) const
Definition: WellInterface_impl.hpp:878
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)
Definition: WellInterface_impl.hpp:91
bool iterateWellEquations(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:469
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: WellInterface.hpp:86
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:83
void checkWellOperability(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:898
std::optional< Scalar > estimateOperableBhp(const Simulator &ebos_simulator, const double dt, WellState< Scalar > &well_state, const SummaryState &summary_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:599
void updateWellOperability(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1001
void updateWellStateRates(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1479
virtual void updateWellStateWithTarget(const Simulator &simulator, const GroupState< Scalar > &group_state, WellState< Scalar > &well_state, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1055
Eval getPerfCellPressure(const FluidState &fs) const
Definition: WellInterface_impl.hpp:1697
void updateConnectionTransmissibilityFactor(const Simulator &simulator, SingleWellState< Scalar > &ws) const
Definition: WellInterface_impl.hpp:1668
WellInterface(const Well &well, const ParallelWellInfo &pw_info, const int time_step, const ModelParameters ¶m, const RateConverterType &rate_converter, const int pvtRegionIdx, const int num_components, const int num_phases, const int index_of_well, const std::vector< PerforationData > &perf_data)
Constructor.
Definition: WellInterface_impl.hpp:54
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:84
static constexpr bool has_solvent
Definition: WellInterface.hpp:110
void prepareWellBeforeAssembling(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:804
bool gliftBeginTimeStepWellTestIterateWellEquations(const Simulator &ebos_simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:923
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: WellInterface.hpp:89
void computeConnLevelInjInd(const FluidState &fs, const Phase preferred_phase, const std::function< Scalar(const Scalar)> &connIICalc, const std::vector< Scalar > &mobility, Scalar *connII, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1852
void solveWellEquation(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:713
void updateConnectionDFactor(const Simulator &simulator, SingleWellState< Scalar > &ws) const
Definition: WellInterface_impl.hpp:1606
std::vector< Scalar > initialWellRateFractions(const Simulator &ebosSimulator, const WellState< Scalar > &well_state) const
Definition: WellInterface_impl.hpp:1432
Scalar wsalt() const
Definition: WellInterface_impl.hpp:140
bool updateWellOperabilityFromWellEq(const Simulator &simulator, const WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:1037
void gliftBeginTimeStepWellTestUpdateALQ(const Simulator &simulator, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:957
typename Base::Eval Eval
Definition: WellInterface.hpp:99
void wellTesting(const Simulator &simulator, const double simulation_time, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, WellTestState &welltest_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:365
Scalar wpolymer() const
Definition: WellInterface_impl.hpp:110
Scalar computeConnectionDFactor(const int perf, const IntensiveQuantities &intQuants, const SingleWellState< Scalar > &ws) const
Definition: WellInterface_impl.hpp:1626
bool solveWellWithBhp(const Simulator &simulator, const double dt, const Scalar bhp, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:622
bool solveWellWithZeroRate(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:658
bool solveWellForTesting(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:678
GetPropType< TypeTag, Properties::Indices > Indices
Definition: WellInterface.hpp:85
bool solveWellWithTHPConstraint(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:502
Scalar wurea() const
Definition: WellInterface_impl.hpp:182
Scalar wmicrobes() const
Definition: WellInterface_impl.hpp:152
void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:785
bool updateWellControlAndStatusLocalIteration(const Simulator &simulator, WellState< Scalar > &well_state, const GroupState< Scalar > &group_state, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, const Scalar WQTotal, DeferredLogger &deferred_logger, const bool fixed_control=false, const bool fixed_status=false)
Definition: WellInterface_impl.hpp:271
static constexpr bool has_zFraction
Definition: WellInterface.hpp:111
Definition: WellState.hpp:62
const SingleWellState< Scalar > & well(std::size_t well_index) const
Definition: WellState.hpp:300
void setALQ(const std::string &name, Scalar value)
Definition: WellState.hpp:180
int numPhases() const
The number of phases present.
Definition: WellState.hpp:254
Scalar getALQ(const std::string &name) const
Definition: WellState.hpp:175
@ NONE
Definition: DeferredLogger.hpp:46
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
VFPEvaluation bhp(const VFPProdTable &table, const double aqua, const double liquid, const double vapour, const double thp, const double alq, const double explicit_wfr, const double explicit_gfr, const bool use_vfpexplicit)
Definition: BlackoilPhases.hpp:27
PhaseUsage phaseUsage(const Phases &phases)
Determine the active phases.
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:485
Definition: BlackoilPhases.hpp:46