22#ifndef OPM_MULTISEGMENTWELL_IMPL_HEADER_INCLUDED
23#define OPM_MULTISEGMENTWELL_IMPL_HEADER_INCLUDED
25#ifndef OPM_MULTISEGMENTWELL_HEADER_INCLUDED
30#include <opm/common/Exceptions.hpp>
31#include <opm/common/OpmLog/OpmLog.hpp>
33#include <opm/input/eclipse/Schedule/MSW/Segment.hpp>
34#include <opm/input/eclipse/Schedule/MSW/Valve.hpp>
35#include <opm/input/eclipse/Schedule/MSW/WellSegments.hpp>
36#include <opm/input/eclipse/Schedule/Well/Connection.hpp>
37#include <opm/input/eclipse/Schedule/Well/WellConnections.hpp>
39#include <opm/input/eclipse/Units/Units.hpp>
41#include <opm/material/densead/EvaluationFormat.hpp>
52#if COMPILE_GPU_BRIDGE && (HAVE_CUDA || HAVE_OPENCL)
60 template <
typename TypeTag>
67 const int pvtRegionIdx,
68 const int num_conservation_quantities,
70 const int index_of_well,
72 :
Base(well, pw_info, time_step, param, rate_converter, pvtRegionIdx, num_conservation_quantities, num_phases, index_of_well, perf_data)
75 , segment_fluid_initial_(this->numberOfSegments(), std::vector<
Scalar>(this->num_conservation_quantities_, 0.0))
79 OPM_THROW(std::runtime_error,
"solvent is not supported by multisegment well yet");
83 OPM_THROW(std::runtime_error,
"polymer is not supported by multisegment well yet");
87 OPM_THROW(std::runtime_error,
"energy is not supported by multisegment well yet");
91 OPM_THROW(std::runtime_error,
"foam is not supported by multisegment well yet");
95 OPM_THROW(std::runtime_error,
"brine is not supported by multisegment well yet");
99 OPM_THROW(std::runtime_error,
"water evaporation is not supported by multisegment well yet");
103 OPM_THROW(std::runtime_error,
"MICP is not supported by multisegment well yet");
107 OPM_THROW(std::runtime_error,
108 "dissolved gas/ vapporized oil in injected oil/gas not supported by multisegment well yet."
109 " \n See (WCONINJE item 10 / WCONHIST item 8)");
119 template <
typename TypeTag>
122 init(
const std::vector<Scalar>& depth_arg,
124 const std::vector< Scalar >& B_avg,
125 const bool changed_to_open_this_step)
127 Base::init(depth_arg, gravity_arg, B_avg, changed_to_open_this_step);
139 this->initMatrixAndVectors();
142 for (
int local_perf_index = 0; local_perf_index < this->number_of_local_perforations_; ++local_perf_index) {
144 const int cell_idx = this->well_cells_[local_perf_index];
146 this->cell_perforation_depth_diffs_[local_perf_index] = depth_arg[cell_idx] - this->perf_depth_[this->pw_info_.localPerfToActivePerf(local_perf_index)];
154 template <
typename TypeTag>
161 const bool stop_or_zero_rate_target = this->stoppedOrZeroRateTarget(simulator, well_state, deferred_logger);
162 this->primary_variables_.update(well_state, stop_or_zero_rate_target);
170 template <
typename TypeTag>
175 this->scaleSegmentRatesWithWellRates(this->segments_.inlets(),
176 this->segments_.perforations(),
178 this->scaleSegmentPressuresWithBhp(well_state);
181 template <
typename TypeTag>
189 Base::updateWellStateWithTarget(simulator, group_state, well_state, deferred_logger);
192 this->scaleSegmentRatesWithWellRates(this->segments_.inlets(),
193 this->segments_.perforations(),
195 this->scaleSegmentPressuresWithBhp(well_state);
201 template <
typename TypeTag>
206 const std::vector<Scalar>& B_avg,
208 const bool relax_tolerance)
const
210 return this->MSWEval::getWellConvergence(well_state,
213 this->param_.max_residual_allowed_,
214 this->param_.tolerance_wells_,
215 this->param_.relaxed_tolerance_flow_well_,
216 this->param_.tolerance_pressure_ms_wells_,
217 this->param_.relaxed_tolerance_pressure_ms_well_,
219 this->wellIsStopped());
227 template <
typename TypeTag>
232 if (!this->isOperableAndSolvable() && !this->wellIsStopped()) {
236 if (this->param_.matrix_add_well_contributions_) {
241 this->linSys_.apply(x, Ax);
248 template <
typename TypeTag>
253 if (!this->isOperableAndSolvable() && !this->wellIsStopped()) {
257 this->linSys_.apply(r);
262 template <
typename TypeTag>
270 if (!this->isOperableAndSolvable() && !this->wellIsStopped()) {
276 this->linSys_.recoverSolutionWell(x, xw);
278 updateWellState(simulator, xw, well_state, deferred_logger);
280 catch (
const NumericalProblem& exp) {
284 deferred_logger.
problem(
"In MultisegmentWell::recoverWellSolutionAndUpdateWellState for well "
285 + this->name() +
": "+exp.what());
294 template <
typename TypeTag>
299 std::vector<Scalar>& well_potentials,
302 const auto [compute_potential, bhp_controlled_well] =
305 if (!compute_potential) {
309 debug_cost_counter_ = 0;
310 bool converged_implicit =
false;
311 if (this->param_.local_well_solver_control_switching_) {
312 converged_implicit = computeWellPotentialsImplicit(simulator, well_state, well_potentials, deferred_logger);
313 if (!converged_implicit) {
314 deferred_logger.
debug(
"Implicit potential calculations failed for well "
315 + this->name() +
", reverting to original aproach.");
318 if (!converged_implicit) {
320 const auto& summaryState = simulator.vanguard().summaryState();
321 if (!Base::wellHasTHPConstraints(summaryState) || bhp_controlled_well) {
322 computeWellRatesAtBhpLimit(simulator, well_potentials, deferred_logger);
324 well_potentials = computeWellPotentialWithTHP(
325 well_state, simulator, deferred_logger);
328 deferred_logger.
debug(
"Cost in iterations of finding well potential for well "
331 this->checkNegativeWellPotentials(well_potentials,
332 this->param_.check_well_operability_,
339 template<
typename TypeTag>
343 std::vector<Scalar>& well_flux,
346 if (this->well_ecl_.isInjector()) {
347 const auto controls = this->well_ecl_.injectionControls(simulator.vanguard().summaryState());
348 computeWellRatesWithBhpIterations(simulator, controls.bhp_limit, well_flux, deferred_logger);
350 const auto controls = this->well_ecl_.productionControls(simulator.vanguard().summaryState());
351 computeWellRatesWithBhpIterations(simulator, controls.bhp_limit, well_flux, deferred_logger);
355 template<
typename TypeTag>
360 std::vector<Scalar>& well_flux,
363 const int np = this->number_of_phases_;
365 well_flux.resize(np, 0.0);
366 const bool allow_cf = this->getAllowCrossFlow();
367 const int nseg = this->numberOfSegments();
368 const WellStateType& well_state = simulator.problem().wellModel().wellState();
369 const auto& ws = well_state.
well(this->indexOfWell());
370 auto segments_copy = ws.segments;
371 segments_copy.scale_pressure(bhp);
372 const auto& segment_pressure = segments_copy.pressure;
373 for (
int seg = 0; seg < nseg; ++seg) {
374 for (
const int perf : this->segments_.perforations()[seg]) {
375 const int local_perf_index = this->pw_info_.activePerfToLocalPerf(perf);
376 if (local_perf_index < 0)
378 const int cell_idx = this->well_cells_[local_perf_index];
379 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
381 std::vector<Scalar> mob(this->num_conservation_quantities_, 0.);
382 getMobility(simulator, local_perf_index, mob, deferred_logger);
383 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(intQuants, cell_idx);
384 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
385 const std::vector<Scalar> Tw = this->wellIndex(local_perf_index, intQuants, trans_mult, wellstate_nupcol);
386 const Scalar seg_pressure = segment_pressure[seg];
387 std::vector<Scalar> cq_s(this->num_conservation_quantities_, 0.);
390 computePerfRate(intQuants, mob, Tw, seg, perf, seg_pressure,
391 allow_cf, cq_s, perf_press, perf_rates, deferred_logger);
393 for(
int p = 0; p < np; ++p) {
394 well_flux[FluidSystem::activeCompToActivePhaseIdx(p)] += cq_s[p];
398 this->parallel_well_info_.communication().sum(well_flux.data(), well_flux.size());
402 template<
typename TypeTag>
407 std::vector<Scalar>& well_flux,
419 WellStateType well_state_copy = simulator.problem().wellModel().wellState();
420 const auto& group_state = simulator.problem().wellModel().groupState();
421 auto& ws = well_state_copy.
well(this->index_of_well_);
424 const auto& summary_state = simulator.vanguard().summaryState();
425 auto inj_controls = well_copy.
well_ecl_.isInjector()
426 ? well_copy.
well_ecl_.injectionControls(summary_state)
427 : Well::InjectionControls(0);
428 auto prod_controls = well_copy.
well_ecl_.isProducer()
429 ? well_copy.
well_ecl_.productionControls(summary_state) :
430 Well::ProductionControls(0);
434 inj_controls.bhp_limit = bhp;
435 ws.injection_cmode = Well::InjectorCMode::BHP;
437 prod_controls.bhp_limit = bhp;
438 ws.production_cmode = Well::ProducerCMode::BHP;
444 const int np = this->number_of_phases_;
446 for (
int phase = 0; phase < np; ++phase){
447 trivial = trivial && (ws.well_potentials[phase] == 0.0) ;
451 for (
int phase = 0; phase < np; ++phase) {
452 ws.surface_rates[phase] = sign * ws.well_potentials[phase];
456 this->segments_.perforations(),
460 const double dt = simulator.timeStepSize();
467 well_flux.resize(np, 0.0);
468 for (
int compIdx = 0; compIdx < this->num_conservation_quantities_; ++compIdx) {
470 well_flux[FluidSystem::activeCompToActivePhaseIdx(compIdx)] = rate.value();
477 template<
typename TypeTag>
478 std::vector<typename MultisegmentWell<TypeTag>::Scalar>
484 std::vector<Scalar> potentials(this->number_of_phases_, 0.0);
485 const auto& summary_state = simulator.vanguard().summaryState();
487 const auto& well = this->well_ecl_;
488 if (well.isInjector()) {
489 auto bhp_at_thp_limit = computeBhpAtThpLimitInj(simulator, summary_state, deferred_logger);
490 if (bhp_at_thp_limit) {
491 const auto& controls = well.injectionControls(summary_state);
492 const Scalar bhp = std::min(*bhp_at_thp_limit,
493 static_cast<Scalar>(controls.bhp_limit));
494 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
495 deferred_logger.
debug(
"Converged thp based potential calculation for well "
498 deferred_logger.
warning(
"FAILURE_GETTING_CONVERGED_POTENTIAL",
499 "Failed in getting converged thp based potential calculation for well "
500 + this->name() +
". Instead the bhp based value is used");
501 const auto& controls = well.injectionControls(summary_state);
502 const Scalar bhp = controls.bhp_limit;
503 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
506 auto bhp_at_thp_limit = computeBhpAtThpLimitProd(
507 well_state, simulator, summary_state, deferred_logger);
508 if (bhp_at_thp_limit) {
509 const auto& controls = well.productionControls(summary_state);
510 const Scalar bhp = std::max(*bhp_at_thp_limit,
511 static_cast<Scalar>(controls.bhp_limit));
512 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
513 deferred_logger.
debug(
"Converged thp based potential calculation for well "
516 deferred_logger.
warning(
"FAILURE_GETTING_CONVERGED_POTENTIAL",
517 "Failed in getting converged thp based potential calculation for well "
518 + this->name() +
". Instead the bhp based value is used");
519 const auto& controls = well.productionControls(summary_state);
520 const Scalar bhp = controls.bhp_limit;
521 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
528 template<
typename TypeTag>
533 std::vector<Scalar>& well_potentials,
544 const auto& group_state = simulator.problem().wellModel().groupState();
545 auto& ws = well_state_copy.
well(this->index_of_well_);
548 const auto& summary_state = simulator.vanguard().summaryState();
549 auto inj_controls = well_copy.
well_ecl_.isInjector()
550 ? well_copy.
well_ecl_.injectionControls(summary_state)
551 : Well::InjectionControls(0);
552 auto prod_controls = well_copy.
well_ecl_.isProducer()
553 ? well_copy.
well_ecl_.productionControls(summary_state)
554 : Well::ProductionControls(0);
562 const int np = this->number_of_phases_;
564 for (
int phase = 0; phase < np; ++phase){
565 trivial = trivial && (ws.well_potentials[phase] == 0.0) ;
569 for (
int phase = 0; phase < np; ++phase) {
570 ws.surface_rates[phase] = sign * ws.well_potentials[phase];
574 this->segments_.perforations(),
578 const double dt = simulator.timeStepSize();
580 bool converged =
false;
581 if (this->well_ecl_.isProducer()) {
582 converged = well_copy.
solveWellWithOperabilityCheck(simulator, dt, inj_controls, prod_controls, well_state_copy, group_state, deferred_logger);
584 converged = well_copy.
iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state_copy, group_state, deferred_logger);
588 well_potentials.clear();
589 well_potentials.resize(np, 0.0);
590 for (
int compIdx = 0; compIdx < this->num_conservation_quantities_; ++compIdx) {
592 well_potentials[FluidSystem::activeCompToActivePhaseIdx(compIdx)] = rate.value();
598 template <
typename TypeTag>
605 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return;
610 const BVectorWell dx_well = this->linSys_.solve();
611 updateWellState(simulator, dx_well, well_state, deferred_logger);
613 catch(
const NumericalProblem& exp) {
617 deferred_logger.
problem(
"In MultisegmentWell::solveEqAndUpdateWellState for well "
618 + this->name() +
": "+exp.what());
627 template <
typename TypeTag>
634 for (
int local_perf_index = 0; local_perf_index < this->number_of_local_perforations_; ++local_perf_index) {
637 std::vector<Scalar> kr(this->number_of_phases_, 0.0);
638 std::vector<Scalar> density(this->number_of_phases_, 0.0);
640 const int cell_idx = this->well_cells_[local_perf_index];
641 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
642 const auto& fs = intQuants.fluidState();
646 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
647 const int water_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::waterPhaseIdx);
648 kr[water_pos] = intQuants.relativePermeability(FluidSystem::waterPhaseIdx).value();
649 sum_kr += kr[water_pos];
650 density[water_pos] = fs.density(FluidSystem::waterPhaseIdx).value();
653 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
654 const int oil_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::oilPhaseIdx);
655 kr[oil_pos] = intQuants.relativePermeability(FluidSystem::oilPhaseIdx).value();
656 sum_kr += kr[oil_pos];
657 density[oil_pos] = fs.density(FluidSystem::oilPhaseIdx).value();
660 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
661 const int gas_pos = FluidSystem::canonicalToActivePhaseIdx(FluidSystem::gasPhaseIdx);
662 kr[gas_pos] = intQuants.relativePermeability(FluidSystem::gasPhaseIdx).value();
663 sum_kr += kr[gas_pos];
664 density[gas_pos] = fs.density(FluidSystem::gasPhaseIdx).value();
667 assert(sum_kr != 0.);
670 Scalar average_density = 0.;
671 for (
int p = 0; p < this->number_of_phases_; ++p) {
672 average_density += kr[p] * density[p];
674 average_density /= sum_kr;
676 this->cell_perforation_pressure_diffs_[local_perf_index] = this->gravity_ * average_density * this->cell_perforation_depth_diffs_[local_perf_index];
684 template <
typename TypeTag>
689 for (
int seg = 0; seg < this->numberOfSegments(); ++seg) {
691 const Scalar surface_volume = getSegmentSurfaceVolume(simulator, seg).value();
692 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) {
693 segment_fluid_initial_[seg][comp_idx] = surface_volume * this->primary_variables_.surfaceVolumeFraction(seg, comp_idx).value();
702 template <
typename TypeTag>
706 const BVectorWell& dwells,
709 const Scalar relaxation_factor)
711 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return;
713 const Scalar dFLimit = this->param_.dwell_fraction_max_;
714 const Scalar max_pressure_change = this->param_.max_pressure_change_ms_wells_;
715 const bool stop_or_zero_rate_target =
716 this->stoppedOrZeroRateTarget(simulator, well_state, deferred_logger);
717 this->primary_variables_.updateNewton(dwells,
720 stop_or_zero_rate_target,
721 max_pressure_change);
723 const auto& summary_state = simulator.vanguard().summaryState();
724 this->primary_variables_.copyToWellState(*
this, getRefDensity(),
730 auto& ws = well_state.
well(this->index_of_well_);
731 this->segments_.copyPhaseDensities(ws.segments);
734 Base::calculateReservoirRates(simulator.vanguard().eclState().runspec().co2Storage(), well_state.
well(this->index_of_well_));
741 template <
typename TypeTag>
748 updatePrimaryVariables(simulator, well_state, deferred_logger);
749 computePerfCellPressDiffs(simulator);
750 computeInitialSegmentFluids(simulator);
757 template<
typename TypeTag>
765 auto fluidState = [&simulator,
this](
const int local_perf_index)
767 const auto cell_idx = this->well_cells_[local_perf_index];
768 return simulator.model()
769 .intensiveQuantities(cell_idx, 0).fluidState();
772 const int np = this->number_of_phases_;
773 auto setToZero = [np](
Scalar* x) ->
void
775 std::fill_n(x, np, 0.0);
778 auto addVector = [np](
const Scalar* src,
Scalar* dest) ->
void
780 std::transform(src, src + np, dest, dest, std::plus<>{});
783 auto& ws = well_state.
well(this->index_of_well_);
784 auto& perf_data = ws.perf_data;
785 auto* connPI = perf_data.prod_index.data();
786 auto* wellPI = ws.productivity_index.data();
790 const auto preferred_phase = this->well_ecl_.getPreferredPhase();
791 auto subsetPerfID = 0;
793 for (
const auto& perf : *this->perf_data_){
794 auto allPerfID = perf.ecl_index;
796 auto connPICalc = [&wellPICalc, allPerfID](
const Scalar mobility) ->
Scalar
801 std::vector<Scalar> mob(this->num_conservation_quantities_, 0.0);
806 getMobility(simulator,
static_cast<int>(subsetPerfID), mob, deferred_logger);
808 const auto& fs = fluidState(subsetPerfID);
811 if (this->isInjector()) {
812 this->computeConnLevelInjInd(fs, preferred_phase, connPICalc,
813 mob, connPI, deferred_logger);
816 this->computeConnLevelProdInd(fs, connPICalc, mob, connPI);
819 addVector(connPI, wellPI);
826 const auto& comm = this->parallel_well_info_.communication();
827 if (comm.size() > 1) {
828 comm.sum(wellPI, np);
831 assert (
static_cast<int>(subsetPerfID) == this->number_of_local_perforations_ &&
832 "Internal logic error in processing connections for PI/II");
839 template<
typename TypeTag>
843 [[maybe_unused]]
const int openConnIdx)
const
848 const auto segNum = this->wellEcl()
849 .getConnections()[globalConnIdx].segment();
851 const auto segIdx = this->wellEcl()
852 .getSegments().segmentNumberToIndex(segNum);
854 return this->segments_.density(segIdx).value();
861 template<
typename TypeTag>
866 if (this->number_of_local_perforations_ == 0) {
870 this->linSys_.extract(jacobian);
874 template<
typename TypeTag>
879 const int pressureVarIndex,
880 const bool use_well_weights,
883 if (this->number_of_local_perforations_ == 0) {
888 this->linSys_.extractCPRPressureMatrix(jacobian,
898 template<
typename TypeTag>
899 template<
class Value>
905 const std::vector<Value>& b_perfcells,
906 const std::vector<Value>& mob_perfcells,
907 const std::vector<Scalar>& Tw,
909 const Value& segment_pressure,
910 const Value& segment_density,
911 const bool& allow_cf,
912 const std::vector<Value>& cmix_s,
913 std::vector<Value>& cq_s,
918 const int local_perf_index = this->pw_info_.activePerfToLocalPerf(perf);
919 if (local_perf_index < 0)
923 const Value perf_seg_press_diff = this->gravity() * segment_density *
924 this->segments_.local_perforation_depth_diff(local_perf_index);
926 const Scalar cell_perf_press_diff = this->cell_perforation_pressure_diffs_[local_perf_index];
930 perf_press = segment_pressure + perf_seg_press_diff;
933 const Value cell_press_at_perf = pressure_cell - cell_perf_press_diff;
936 const Value drawdown = cell_press_at_perf - perf_press;
939 if (drawdown > 0.0) {
941 if (!allow_cf && this->isInjector()) {
946 for (
int comp_idx = 0; comp_idx < this->numConservationQuantities(); ++comp_idx) {
947 const Value cq_p = - Tw[comp_idx] * (mob_perfcells[comp_idx] * drawdown);
948 cq_s[comp_idx] = b_perfcells[comp_idx] * cq_p;
951 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
952 const unsigned oilCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx);
953 const unsigned gasCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx);
954 const Value cq_s_oil = cq_s[oilCompIdx];
955 const Value cq_s_gas = cq_s[gasCompIdx];
956 cq_s[gasCompIdx] += rs * cq_s_oil;
957 cq_s[oilCompIdx] += rv * cq_s_gas;
961 if (!allow_cf && this->isProducer()) {
966 Value total_mob = mob_perfcells[0];
967 for (
int comp_idx = 1; comp_idx < this->numConservationQuantities(); ++comp_idx) {
968 total_mob += mob_perfcells[comp_idx];
972 Value volume_ratio = 0.0;
973 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
974 const unsigned waterCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::waterCompIdx);
975 volume_ratio += cmix_s[waterCompIdx] / b_perfcells[waterCompIdx];
978 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
979 const unsigned oilCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx);
980 const unsigned gasCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx);
985 const Value d = 1.0 - rv * rs;
987 if (getValue(d) == 0.0) {
989 fmt::format(
"Zero d value obtained for well {} "
990 "during flux calculation with rs {} and rv {}",
991 this->name(), rs, rv),
995 const Value tmp_oil = (cmix_s[oilCompIdx] - rv * cmix_s[gasCompIdx]) / d;
996 volume_ratio += tmp_oil / b_perfcells[oilCompIdx];
998 const Value tmp_gas = (cmix_s[gasCompIdx] - rs * cmix_s[oilCompIdx]) / d;
999 volume_ratio += tmp_gas / b_perfcells[gasCompIdx];
1001 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
1002 const unsigned oilCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx);
1003 volume_ratio += cmix_s[oilCompIdx] / b_perfcells[oilCompIdx];
1005 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1006 const unsigned gasCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx);
1007 volume_ratio += cmix_s[gasCompIdx] / b_perfcells[gasCompIdx];
1011 for (
int componentIdx = 0; componentIdx < this->numConservationQuantities(); ++componentIdx) {
1012 const Value cqt_i = - Tw[componentIdx] * (total_mob * drawdown);
1013 Value cqt_is = cqt_i / volume_ratio;
1014 cq_s[componentIdx] = cmix_s[componentIdx] * cqt_is;
1019 if (this->isProducer()) {
1020 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1021 const unsigned oilCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx);
1022 const unsigned gasCompIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx);
1031 const Scalar d = 1.0 - getValue(rv) * getValue(rs);
1034 perf_rates.
vap_oil = getValue(rv) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d;
1037 perf_rates.
dis_gas = getValue(rs) * (getValue(cq_s[oilCompIdx]) - getValue(rv) * getValue(cq_s[gasCompIdx])) / d;
1042 template <
typename TypeTag>
1043 template<
class Value>
1047 const std::vector<Value>& mob_perfcells,
1048 const std::vector<Scalar>& Tw,
1051 const Value& segment_pressure,
1052 const bool& allow_cf,
1053 std::vector<Value>& cq_s,
1059 auto obtain = [
this](
const Eval& value)
1061 if constexpr (std::is_same_v<Value, Scalar>) {
1062 static_cast<void>(
this);
1063 return getValue(value);
1065 return this->extendEval(value);
1068 auto obtainN = [](
const auto& value)
1070 if constexpr (std::is_same_v<Value, Scalar>) {
1071 return getValue(value);
1076 const auto& fs = int_quants.fluidState();
1078 const Value pressure_cell = obtain(this->getPerfCellPressure(fs));
1079 const Value rs = obtain(fs.Rs());
1080 const Value rv = obtain(fs.Rv());
1083 std::vector<Value> b_perfcells(this->num_conservation_quantities_, 0.0);
1085 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
1086 if (!FluidSystem::phaseIsActive(phaseIdx)) {
1090 const unsigned compIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::solventComponentIndex(phaseIdx));
1091 b_perfcells[compIdx] = obtain(fs.invB(phaseIdx));
1094 std::vector<Value> cmix_s(this->numConservationQuantities(), 0.0);
1095 for (
int comp_idx = 0; comp_idx < this->numConservationQuantities(); ++comp_idx) {
1096 cmix_s[comp_idx] = obtainN(this->primary_variables_.surfaceVolumeFraction(seg, comp_idx));
1099 this->computePerfRate(pressure_cell,
1107 obtainN(this->segments_.density(seg)),
1116 template <
typename TypeTag>
1136 auto info = this->getFirstPerforationFluidStateInfo(simulator);
1137 temperature.setValue(std::get<0>(info));
1138 saltConcentration = this->extendEval(std::get<1>(info));
1140 this->segments_.computeFluidProperties(temperature,
1142 this->primary_variables_,
1147 template <
typename TypeTag>
1148 template<
class Value>
1152 const int local_perf_index,
1153 std::vector<Value>& mob,
1156 auto obtain = [
this](
const Eval& value)
1158 if constexpr (std::is_same_v<Value, Scalar>) {
1159 static_cast<void>(
this);
1160 return getValue(value);
1162 return this->extendEval(value);
1169 const auto perf_ecl_index = this->perforationData()[local_perf_index].ecl_index;
1170 const Connection& con = this->well_ecl_.getConnections()[perf_ecl_index];
1171 const int seg = this->segmentNumberToIndex(con.segment());
1175 const Scalar segment_pres = this->primary_variables_.getSegmentPressure(seg).value();
1176 const Scalar perf_seg_press_diff = this->gravity() * this->segments_.density(seg).value()
1177 * this->segments_.local_perforation_depth_diff(local_perf_index);
1178 const Scalar perf_press = segment_pres + perf_seg_press_diff;
1179 const Scalar multiplier = this->getInjMult(local_perf_index, segment_pres, perf_press, deferred_logger);
1180 for (std::size_t i = 0; i < mob.size(); ++i) {
1181 mob[i] *= multiplier;
1188 template<
typename TypeTag>
1193 return this->segments_.getRefDensity();
1196 template<
typename TypeTag>
1203 const auto& summaryState = simulator.vanguard().summaryState();
1207 const bool bhp_limit_not_defaulted = bhp_limit > 1.5 * unit::barsa;
1208 if ( bhp_limit_not_defaulted || !this->wellHasTHPConstraints(summaryState) ) {
1211 Scalar total_ipr_mass_rate = 0.0;
1212 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx)
1214 if (!FluidSystem::phaseIsActive(phaseIdx)) {
1218 const unsigned compIdx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::solventComponentIndex(phaseIdx));
1219 const Scalar ipr_rate = this->ipr_a_[compIdx] - this->ipr_b_[compIdx] * bhp_limit;
1221 const Scalar rho = FluidSystem::referenceDensity( phaseIdx, Base::pvtRegionIdx() );
1222 total_ipr_mass_rate += ipr_rate * rho;
1224 if ( (this->isProducer() && total_ipr_mass_rate < 0.) || (this->isInjector() && total_ipr_mass_rate > 0.) ) {
1225 this->operability_status_.operable_under_only_bhp_limit =
false;
1229 if (this->operability_status_.operable_under_only_bhp_limit && this->wellHasTHPConstraints(summaryState)) {
1233 std::vector<Scalar> well_rates_bhp_limit;
1234 computeWellRatesWithBhp(simulator, bhp_limit, well_rates_bhp_limit, deferred_logger);
1236 const Scalar thp_limit = this->getTHPConstraint(summaryState);
1239 this->getRefDensity(),
1240 this->wellEcl().alq_value(summaryState),
1243 if ( (this->isProducer() && thp < thp_limit) || (this->isInjector() && thp > thp_limit) ) {
1244 this->operability_status_.obey_thp_limit_under_bhp_limit =
false;
1255 this->operability_status_.operable_under_only_bhp_limit =
true;
1256 this->operability_status_.obey_thp_limit_under_bhp_limit =
false;
1262 template<
typename TypeTag>
1270 std::fill(this->ipr_a_.begin(), this->ipr_a_.end(), 0.);
1271 std::fill(this->ipr_b_.begin(), this->ipr_b_.end(), 0.);
1273 const int nseg = this->numberOfSegments();
1274 std::vector<Scalar> seg_dp(nseg, 0.0);
1275 for (
int seg = 0; seg < nseg; ++seg) {
1277 const Scalar dp = this->getSegmentDp(seg,
1278 this->segments_.density(seg).value(),
1281 for (
const int perf : this->segments_.perforations()[seg]) {
1282 const int local_perf_index = this->pw_info_.activePerfToLocalPerf(perf);
1283 if (local_perf_index < 0)
1285 std::vector<Scalar> mob(this->num_conservation_quantities_, 0.0);
1288 getMobility(simulator, local_perf_index, mob, deferred_logger);
1290 const int cell_idx = this->well_cells_[local_perf_index];
1291 const auto& int_quantities = simulator.model().intensiveQuantities(cell_idx, 0);
1292 const auto& fs = int_quantities.fluidState();
1294 const Scalar perf_seg_press_diff = this->segments_.getPressureDiffSegLocalPerf(seg, local_perf_index);
1296 const Scalar cell_perf_press_diff = this->cell_perforation_pressure_diffs_[local_perf_index];
1297 const Scalar pressure_cell = this->getPerfCellPressure(fs).value();
1300 std::vector<Scalar> b_perf(this->num_conservation_quantities_);
1301 for (std::size_t phase = 0; phase < FluidSystem::numPhases; ++phase) {
1302 if (!FluidSystem::phaseIsActive(phase)) {
1305 const unsigned comp_idx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::solventComponentIndex(phase));
1306 b_perf[comp_idx] = fs.invB(phase).value();
1310 const Scalar h_perf = cell_perf_press_diff + perf_seg_press_diff + dp;
1311 const Scalar pressure_diff = pressure_cell - h_perf;
1314 if ( (this->isProducer() && pressure_diff < 0.) || (this->isInjector() && pressure_diff > 0.) ) {
1315 deferred_logger.
debug(
"CROSSFLOW_IPR",
1316 "cross flow found when updateIPR for well " + this->name());
1320 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(int_quantities, cell_idx);
1321 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
1322 const std::vector<Scalar> tw_perf = this->wellIndex(local_perf_index, int_quantities, trans_mult, wellstate_nupcol);
1323 std::vector<Scalar> ipr_a_perf(this->ipr_a_.size());
1324 std::vector<Scalar> ipr_b_perf(this->ipr_b_.size());
1325 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) {
1326 const Scalar tw_mob = tw_perf[comp_idx] * mob[comp_idx] * b_perf[comp_idx];
1327 ipr_a_perf[comp_idx] += tw_mob * pressure_diff;
1328 ipr_b_perf[comp_idx] += tw_mob;
1332 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1333 const unsigned oil_comp_idx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx);
1334 const unsigned gas_comp_idx = FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx);
1335 const Scalar rs = (fs.Rs()).value();
1336 const Scalar rv = (fs.Rv()).value();
1338 const Scalar dis_gas_a = rs * ipr_a_perf[oil_comp_idx];
1339 const Scalar vap_oil_a = rv * ipr_a_perf[gas_comp_idx];
1341 ipr_a_perf[gas_comp_idx] += dis_gas_a;
1342 ipr_a_perf[oil_comp_idx] += vap_oil_a;
1344 const Scalar dis_gas_b = rs * ipr_b_perf[oil_comp_idx];
1345 const Scalar vap_oil_b = rv * ipr_b_perf[gas_comp_idx];
1347 ipr_b_perf[gas_comp_idx] += dis_gas_b;
1348 ipr_b_perf[oil_comp_idx] += vap_oil_b;
1351 for (std::size_t comp_idx = 0; comp_idx < ipr_a_perf.size(); ++comp_idx) {
1352 this->ipr_a_[comp_idx] += ipr_a_perf[comp_idx];
1353 this->ipr_b_[comp_idx] += ipr_b_perf[comp_idx];
1357 this->parallel_well_info_.communication().sum(this->ipr_a_.data(), this->ipr_a_.size());
1358 this->parallel_well_info_.communication().sum(this->ipr_b_.data(), this->ipr_b_.size());
1361 template<
typename TypeTag>
1375 auto rates = well_state.
well(this->index_of_well_).surface_rates;
1377 for (std::size_t p = 0; p < rates.size(); ++p) {
1378 zero_rates &= rates[p] == 0.0;
1380 auto& ws = well_state.
well(this->index_of_well_);
1382 const auto msg = fmt::format(
"updateIPRImplicit: Well {} has zero rate, IPRs might be problematic", this->name());
1383 deferred_logger.
debug(msg);
1395 const auto& group_state = simulator.problem().wellModel().groupState();
1397 std::fill(ws.implicit_ipr_a.begin(), ws.implicit_ipr_a.end(), 0.);
1398 std::fill(ws.implicit_ipr_b.begin(), ws.implicit_ipr_b.end(), 0.);
1400 auto inj_controls = Well::InjectionControls(0);
1401 auto prod_controls = Well::ProductionControls(0);
1402 prod_controls.addControl(Well::ProducerCMode::BHP);
1403 prod_controls.bhp_limit = well_state.
well(this->index_of_well_).bhp;
1406 const auto cmode = ws.production_cmode;
1407 ws.production_cmode = Well::ProducerCMode::BHP;
1408 const double dt = simulator.timeStepSize();
1409 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
1411 BVectorWell rhs(this->numberOfSegments());
1413 rhs[0][SPres] = -1.0;
1415 const BVectorWell x_well = this->linSys_.solve(rhs);
1416 constexpr int num_eq = MSWEval::numWellEq;
1417 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx){
1418 const EvalWell comp_rate = this->primary_variables_.getQs(comp_idx);
1419 const int idx = FluidSystem::activeCompToActivePhaseIdx(comp_idx);
1420 for (
size_t pvIdx = 0; pvIdx < num_eq; ++pvIdx) {
1422 ws.implicit_ipr_b[idx] -= x_well[0][pvIdx]*comp_rate.derivative(pvIdx+Indices::numEq);
1424 ws.implicit_ipr_a[idx] = ws.implicit_ipr_b[idx]*ws.bhp - comp_rate.value();
1427 ws.production_cmode = cmode;
1430 template<
typename TypeTag>
1437 const auto& summaryState = simulator.vanguard().summaryState();
1438 const auto obtain_bhp = this->isProducer()
1439 ? computeBhpAtThpLimitProd(
1440 well_state, simulator, summaryState, deferred_logger)
1441 : computeBhpAtThpLimitInj(simulator, summaryState, deferred_logger);
1444 this->operability_status_.can_obtain_bhp_with_thp_limit =
true;
1447 this->operability_status_.obey_bhp_limit_with_thp_limit = (*obtain_bhp >= bhp_limit);
1449 const Scalar thp_limit = this->getTHPConstraint(summaryState);
1450 if (this->isProducer() && *obtain_bhp < thp_limit) {
1451 const std::string msg =
" obtained bhp " +
std::to_string(unit::convert::to(*obtain_bhp, unit::barsa))
1452 +
" bars is SMALLER than thp limit "
1454 +
" bars as a producer for well " + this->name();
1455 deferred_logger.
debug(msg);
1457 else if (this->isInjector() && *obtain_bhp > thp_limit) {
1458 const std::string msg =
" obtained bhp " +
std::to_string(unit::convert::to(*obtain_bhp, unit::barsa))
1459 +
" bars is LARGER than thp limit "
1461 +
" bars as a injector for well " + this->name();
1462 deferred_logger.
debug(msg);
1467 this->operability_status_.can_obtain_bhp_with_thp_limit =
false;
1468 this->operability_status_.obey_bhp_limit_with_thp_limit =
false;
1469 if (!this->wellIsStopped()) {
1470 const Scalar thp_limit = this->getTHPConstraint(summaryState);
1471 deferred_logger.
debug(
" could not find bhp value at thp limit "
1473 +
" bar for well " + this->name() +
", the well might need to be closed ");
1482 template<
typename TypeTag>
1487 const Well::InjectionControls& inj_controls,
1488 const Well::ProductionControls& prod_controls,
1493 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return true;
1495 const int max_iter_number = this->param_.max_inner_iter_ms_wells_;
1499 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1504 updatePrimaryVariables(simulator, well_state, deferred_logger);
1506 std::vector<std::vector<Scalar> > residual_history;
1507 std::vector<Scalar> measure_history;
1510 Scalar relaxation_factor = 1.;
1511 bool converged =
false;
1512 bool relax_convergence =
false;
1513 this->regularize_ =
false;
1514 for (; it < max_iter_number; ++it, ++debug_cost_counter_) {
1516 if (it > this->param_.strict_inner_iter_wells_) {
1517 relax_convergence =
true;
1518 this->regularize_ =
true;
1521 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls,
1522 well_state, group_state, deferred_logger);
1524 const auto report = getWellConvergence(simulator, well_state, Base::B_avg_, deferred_logger, relax_convergence);
1525 if (report.converged()) {
1532 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1536 residual_history.push_back(residuals);
1537 measure_history.push_back(this->getResidualMeasureValue(well_state,
1538 residual_history[it],
1539 this->param_.tolerance_wells_,
1540 this->param_.tolerance_pressure_ms_wells_,
1543 bool min_relaxation_reached = this->update_relaxation_factor(measure_history, relaxation_factor, this->regularize_, deferred_logger);
1544 if (min_relaxation_reached || this->repeatedStagnation(measure_history, this->regularize_, deferred_logger)) {
1546 const auto reportStag = getWellConvergence(simulator, well_state, Base::B_avg_, deferred_logger,
true);
1547 if (reportStag.converged()) {
1549 std::string message = fmt::format(
"Well stagnates/oscillates but {} manages to get converged with relaxed tolerances in {} inner iterations."
1551 deferred_logger.
debug(message);
1558 BVectorWell dx_well;
1560 dx_well = this->linSys_.solve();
1561 updateWellState(simulator, dx_well, well_state, deferred_logger, relaxation_factor);
1563 catch(
const NumericalProblem& exp) {
1567 deferred_logger.
problem(
"In MultisegmentWell::iterateWellEqWithControl for well "
1568 + this->name() +
": "+exp.what());
1575 std::ostringstream sstr;
1576 sstr <<
" Well " << this->name() <<
" converged in " << it <<
" inner iterations.";
1577 if (relax_convergence)
1578 sstr <<
" (A relaxed tolerance was used after "<< this->param_.strict_inner_iter_wells_ <<
" iterations)";
1582 deferred_logger.
debug(sstr.str(), OpmLog::defaultDebugVerbosityLevel + (it == 0));
1584 std::ostringstream sstr;
1585 sstr <<
" Well " << this->name() <<
" did not converge in " << it <<
" inner iterations.";
1586#define EXTRA_DEBUG_MSW 0
1588 sstr <<
"***** Outputting the residual history for well " << this->name() <<
" during inner iterations:";
1589 for (
int i = 0; i < it; ++i) {
1590 const auto& residual = residual_history[i];
1591 sstr <<
" residual at " << i <<
"th iteration ";
1592 for (
const auto& res : residual) {
1595 sstr <<
" " << measure_history[i] <<
" \n";
1598#undef EXTRA_DEBUG_MSW
1599 deferred_logger.
debug(sstr.str());
1606 template<
typename TypeTag>
1611 const Well::InjectionControls& inj_controls,
1612 const Well::ProductionControls& prod_controls,
1616 const bool fixed_control ,
1617 const bool fixed_status )
1619 const int max_iter_number = this->param_.max_inner_iter_ms_wells_;
1623 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1628 updatePrimaryVariables(simulator, well_state, deferred_logger);
1630 std::vector<std::vector<Scalar> > residual_history;
1631 std::vector<Scalar> measure_history;
1634 Scalar relaxation_factor = 1.;
1635 bool converged =
false;
1636 bool relax_convergence =
false;
1637 this->regularize_ =
false;
1638 const auto& summary_state = simulator.vanguard().summaryState();
1643 const int min_its_after_switch = 3;
1645 const int max_status_switch = this->param_.max_well_status_switch_inner_iter_;
1646 int its_since_last_switch = min_its_after_switch;
1647 int switch_count= 0;
1648 int status_switch_count = 0;
1650 const auto well_status_orig = this->wellStatus_;
1651 const auto operability_orig = this->operability_status_;
1652 auto well_status_cur = well_status_orig;
1654 const bool allow_open = well_state.
well(this->index_of_well_).status == WellStatus::OPEN;
1656 const bool allow_switching = !this->wellUnderZeroRateTarget(simulator, well_state, deferred_logger) &&
1657 (!fixed_control || !fixed_status) && allow_open;
1658 bool final_check =
false;
1660 this->operability_status_.resetOperability();
1661 this->operability_status_.solvable =
true;
1663 for (; it < max_iter_number; ++it, ++debug_cost_counter_) {
1664 ++its_since_last_switch;
1665 if (allow_switching && its_since_last_switch >= min_its_after_switch && status_switch_count < max_status_switch){
1666 const Scalar wqTotal = this->primary_variables_.getWQTotal().value();
1667 bool changed = this->updateWellControlAndStatusLocalIteration(simulator, well_state, group_state,
1668 inj_controls, prod_controls, wqTotal,
1669 deferred_logger, fixed_control,
1672 its_since_last_switch = 0;
1674 if (well_status_cur != this->wellStatus_) {
1675 well_status_cur = this->wellStatus_;
1676 status_switch_count++;
1679 if (!changed && final_check) {
1682 final_check =
false;
1684 if (status_switch_count == max_status_switch) {
1685 this->wellStatus_ = well_status_orig;
1689 if (it > this->param_.strict_inner_iter_wells_) {
1690 relax_convergence =
true;
1691 this->regularize_ =
true;
1694 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls,
1695 well_state, group_state, deferred_logger);
1698 const auto report = getWellConvergence(simulator, well_state, Base::B_avg_, deferred_logger, relax_convergence);
1699 converged = report.converged();
1700 if (this->parallel_well_info_.communication().size() > 1 &&
1701 this->parallel_well_info_.communication().max(converged) != this->parallel_well_info_.communication().min(converged)) {
1702 OPM_THROW(std::runtime_error, fmt::format(
"Misalignment of the parallel simulation run in iterateWellEqWithSwitching - the well calculation for well {} succeeded some ranks but failed on other ranks.", this->name()));
1707 if (switch_count > 0 && its_since_last_switch < min_its_after_switch) {
1709 its_since_last_switch = min_its_after_switch;
1717 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1721 residual_history.push_back(residuals);
1725 measure_history.push_back(this->getResidualMeasureValue(well_state,
1726 residual_history[it],
1727 this->param_.tolerance_wells_,
1728 this->param_.tolerance_pressure_ms_wells_,
1730 bool min_relaxation_reached = this->update_relaxation_factor(measure_history, relaxation_factor, this->regularize_, deferred_logger);
1731 if (min_relaxation_reached || this->repeatedStagnation(measure_history, this->regularize_, deferred_logger)) {
1737 const BVectorWell dx_well = this->linSys_.solve();
1738 updateWellState(simulator, dx_well, well_state, deferred_logger, relaxation_factor);
1740 catch(
const NumericalProblem& exp) {
1744 deferred_logger.
problem(
"In MultisegmentWell::iterateWellEqWithSwitching for well "
1745 + this->name() +
": "+exp.what());
1751 if (allow_switching){
1753 const bool is_stopped = this->wellIsStopped();
1754 if (this->wellHasTHPConstraints(summary_state)){
1755 this->operability_status_.can_obtain_bhp_with_thp_limit = !is_stopped;
1756 this->operability_status_.obey_thp_limit_under_bhp_limit = !is_stopped;
1758 this->operability_status_.operable_under_only_bhp_limit = !is_stopped;
1761 std::string message = fmt::format(
" Well {} converged in {} inner iterations ("
1762 "{} control/status switches).", this->name(), it, switch_count);
1763 if (relax_convergence) {
1764 message.append(fmt::format(
" (A relaxed tolerance was used after {} iterations)",
1765 this->param_.strict_inner_iter_wells_));
1767 deferred_logger.
debug(message, OpmLog::defaultDebugVerbosityLevel + ((it == 0) && (switch_count == 0)));
1769 this->wellStatus_ = well_status_orig;
1770 this->operability_status_ = operability_orig;
1771 const std::string message = fmt::format(
" Well {} did not converge in {} inner iterations ("
1772 "{} switches, {} status changes).", this->name(), it, switch_count, status_switch_count);
1773 deferred_logger.
debug(message);
1774 this->primary_variables_.outputLowLimitPressureSegments(deferred_logger);
1781 template<
typename TypeTag>
1786 const Well::InjectionControls& inj_controls,
1787 const Well::ProductionControls& prod_controls,
1792 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return;
1795 this->segments_.updateUpwindingSegments(this->primary_variables_);
1798 computeSegmentFluidProperties(simulator, deferred_logger);
1801 this->linSys_.clear();
1803 auto& ws = well_state.
well(this->index_of_well_);
1804 ws.phase_mixing_rates.fill(0.0);
1811 const bool allow_cf = this->getAllowCrossFlow() || openCrossFlowAvoidSingularity(simulator);
1813 const int nseg = this->numberOfSegments();
1815 const Scalar rhow = FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) ?
1816 FluidSystem::referenceDensity( FluidSystem::waterPhaseIdx, Base::pvtRegionIdx() ) : 0.0;
1817 const unsigned watCompIdx = FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) ?
1818 FluidSystem::canonicalToActiveCompIdx(FluidSystem::waterCompIdx) : 0;
1820 for (
int seg = 0; seg < nseg; ++seg) {
1822 const EvalWell seg_pressure = this->primary_variables_.getSegmentPressure(seg);
1823 auto& perf_data = ws.perf_data;
1824 auto& perf_rates = perf_data.phase_rates;
1825 auto& perf_press_state = perf_data.pressure;
1826 for (
const int perf : this->segments_.perforations()[seg]) {
1827 const int local_perf_index = this->pw_info_.activePerfToLocalPerf(perf);
1828 if (local_perf_index < 0)
1830 const int cell_idx = this->well_cells_[local_perf_index];
1831 const auto& int_quants = simulator.model().intensiveQuantities(cell_idx, 0);
1832 std::vector<EvalWell> mob(this->num_conservation_quantities_, 0.0);
1833 getMobility(simulator, local_perf_index, mob, deferred_logger);
1834 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(int_quants, cell_idx);
1835 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
1836 const std::vector<Scalar> Tw = this->wellIndex(local_perf_index, int_quants, trans_mult, wellstate_nupcol);
1837 std::vector<EvalWell> cq_s(this->num_conservation_quantities_, 0.0);
1840 computePerfRate(int_quants, mob, Tw, seg, perf, seg_pressure,
1841 allow_cf, cq_s, perf_press, perfRates, deferred_logger);
1844 if (this->isProducer()) {
1845 ws.phase_mixing_rates[ws.dissolved_gas] += perfRates.
dis_gas;
1846 ws.phase_mixing_rates[ws.vaporized_oil] += perfRates.
vap_oil;
1847 perf_data.phase_mixing_rates[local_perf_index][ws.dissolved_gas] = perfRates.
dis_gas;
1848 perf_data.phase_mixing_rates[local_perf_index][ws.vaporized_oil] = perfRates.
vap_oil;
1852 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) {
1853 perf_rates[local_perf_index*this->number_of_phases_ + FluidSystem::activeCompToActivePhaseIdx(comp_idx)] = cq_s[comp_idx].value();
1855 perf_press_state[local_perf_index] = perf_press.value();
1858 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
1859 perf_data.wat_mass_rates[local_perf_index] = cq_s[watCompIdx].value() * rhow;
1862 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) {
1864 const EvalWell cq_s_effective = cq_s[comp_idx] * this->well_efficiency_factor_;
1866 this->connectionRates_[local_perf_index][comp_idx] = Base::restrictEval(cq_s_effective);
1869 assemblePerforationEq(seg, local_perf_index, comp_idx, cq_s_effective, this->linSys_);
1875 const auto& comm = this->parallel_well_info_.communication();
1876 comm.sum(ws.phase_mixing_rates.data(), ws.phase_mixing_rates.size());
1879 if (this->parallel_well_info_.communication().size() > 1) {
1881 this->linSys_.sumDistributed(this->parallel_well_info_.communication());
1883 for (
int seg = 0; seg < nseg; ++seg) {
1887 const EvalWell segment_surface_volume = getSegmentSurfaceVolume(simulator, seg);
1892 const Scalar regularization_factor = this->regularize_? this->param_.regularization_factor_wells_ : 1.0;
1894 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) {
1895 const EvalWell accumulation_term = regularization_factor * (segment_surface_volume * this->primary_variables_.surfaceVolumeFraction(seg, comp_idx)
1896 - segment_fluid_initial_[seg][comp_idx]) / dt;
1898 assembleAccumulationTerm(seg, comp_idx, accumulation_term, this->linSys_);
1903 const int seg_upwind = this->segments_.upwinding_segment(seg);
1904 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) {
1906 this->primary_variables_.getSegmentRateUpwinding(seg,
1909 this->well_efficiency_factor_;
1911 assembleOutflowTerm(seg, seg_upwind, comp_idx, segment_rate, this->linSys_);
1917 for (
const int inlet : this->segments_.inlets()[seg]) {
1918 const int inlet_upwind = this->segments_.upwinding_segment(inlet);
1919 for (
int comp_idx = 0; comp_idx < this->num_conservation_quantities_; ++comp_idx) {
1921 this->primary_variables_.getSegmentRateUpwinding(inlet,
1924 this->well_efficiency_factor_;
1926 assembleInflowTerm(seg, inlet, inlet_upwind, comp_idx, inlet_rate, this->linSys_);
1933 const auto& summaryState = simulator.vanguard().summaryState();
1934 const Schedule& schedule = simulator.vanguard().schedule();
1935 const bool stopped_or_zero_target = this->stoppedOrZeroRateTarget(simulator, well_state, deferred_logger);
1937 assembleControlEq(well_state,
1944 this->primary_variables_,
1946 stopped_or_zero_target,
1949 const UnitSystem& unit_system = simulator.vanguard().eclState().getDeckUnitSystem();
1950 const auto& summary_state = simulator.vanguard().summaryState();
1951 this->assemblePressureEq(seg, unit_system, well_state, summary_state, this->param_.use_average_density_ms_wells_, deferred_logger);
1955 this->parallel_well_info_.communication().sum(this->ipr_a_.data(), this->ipr_a_.size());
1956 this->linSys_.createSolver();
1962 template<
typename TypeTag>
1967 return !this->getAllowCrossFlow() && allDrawDownWrongDirection(simulator);
1971 template<
typename TypeTag>
1976 bool all_drawdown_wrong_direction =
true;
1977 const int nseg = this->numberOfSegments();
1979 for (
int seg = 0; seg < nseg; ++seg) {
1980 const EvalWell segment_pressure = this->primary_variables_.getSegmentPressure(seg);
1981 for (
const int perf : this->segments_.perforations()[seg]) {
1982 const int local_perf_index = this->pw_info_.activePerfToLocalPerf(perf);
1983 if (local_perf_index < 0)
1986 const int cell_idx = this->well_cells_[local_perf_index];
1987 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
1988 const auto& fs = intQuants.fluidState();
1991 const EvalWell perf_seg_press_diff = this->segments_.getPressureDiffSegLocalPerf(seg, local_perf_index);
1993 const Scalar cell_perf_press_diff = this->cell_perforation_pressure_diffs_[local_perf_index];
1995 const Scalar pressure_cell = this->getPerfCellPressure(fs).value();
1996 const Scalar perf_press = pressure_cell - cell_perf_press_diff;
1999 const EvalWell drawdown = perf_press - (segment_pressure + perf_seg_press_diff);
2004 if ( (drawdown < 0. && this->isInjector()) ||
2005 (drawdown > 0. && this->isProducer()) ) {
2006 all_drawdown_wrong_direction =
false;
2011 const auto& comm = this->parallel_well_info_.communication();
2012 if (comm.size() > 1)
2014 all_drawdown_wrong_direction =
2015 (comm.min(all_drawdown_wrong_direction ? 1 : 0) == 1);
2018 return all_drawdown_wrong_direction;
2024 template<
typename TypeTag>
2035 template<
typename TypeTag>
2043 auto info = this->getFirstPerforationFluidStateInfo(simulator);
2044 temperature.setValue(std::get<0>(info));
2045 saltConcentration = this->extendEval(std::get<1>(info));
2047 return this->segments_.getSurfaceVolume(temperature,
2049 this->primary_variables_,
2055 template<
typename TypeTag>
2056 std::optional<typename MultisegmentWell<TypeTag>::Scalar>
2060 const SummaryState& summary_state,
2066 this->getALQ(well_state),
2073 template<
typename TypeTag>
2074 std::optional<typename MultisegmentWell<TypeTag>::Scalar>
2077 const SummaryState& summary_state,
2080 bool iterate_if_no_solution)
const
2084 auto frates = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2090 std::vector<Scalar> rates(3);
2091 computeWellRatesWithBhp(simulator, bhp, rates, deferred_logger);
2096 computeBhpAtThpLimitProd(frates,
2098 this->maxPerfPress(simulator),
2099 this->getRefDensity(),
2101 this->getTHPConstraint(summary_state),
2107 if (!iterate_if_no_solution)
2108 return std::nullopt;
2110 auto fratesIter = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2114 std::vector<Scalar> rates(3);
2115 computeWellRatesWithBhpIterations(simulator, bhp, rates, deferred_logger);
2120 computeBhpAtThpLimitProd(fratesIter,
2122 this->maxPerfPress(simulator),
2123 this->getRefDensity(),
2125 this->getTHPConstraint(summary_state),
2129 template<
typename TypeTag>
2130 std::optional<typename MultisegmentWell<TypeTag>::Scalar>
2133 const SummaryState& summary_state,
2137 auto frates = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2143 std::vector<Scalar> rates(3);
2144 computeWellRatesWithBhp(simulator, bhp, rates, deferred_logger);
2149 computeBhpAtThpLimitInj(frates,
2151 this->getRefDensity(),
2160 auto fratesIter = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2164 std::vector<Scalar> rates(3);
2165 computeWellRatesWithBhpIterations(simulator, bhp, rates, deferred_logger);
2170 computeBhpAtThpLimitInj(fratesIter,
2172 this->getRefDensity(),
2183 template<
typename TypeTag>
2188 Scalar max_pressure = 0.0;
2189 const int nseg = this->numberOfSegments();
2190 for (
int seg = 0; seg < nseg; ++seg) {
2191 for (
const int perf : this->segments_.perforations()[seg]) {
2192 const int local_perf_index = this->pw_info_.activePerfToLocalPerf(perf);
2193 if (local_perf_index < 0)
2196 const int cell_idx = this->well_cells_[local_perf_index];
2197 const auto& int_quants = simulator.model().intensiveQuantities(cell_idx, 0);
2198 const auto& fs = int_quants.fluidState();
2199 Scalar pressure_cell = this->getPerfCellPressure(fs).value();
2200 max_pressure = std::max(max_pressure, pressure_cell);
2203 max_pressure = this->parallel_well_info_.communication().max(max_pressure);
2204 return max_pressure;
2211 template<
typename TypeTag>
2212 std::vector<typename MultisegmentWell<TypeTag>::Scalar>
2218 std::vector<Scalar> well_q_s(this->num_conservation_quantities_, 0.0);
2219 const bool allow_cf = this->getAllowCrossFlow() || openCrossFlowAvoidSingularity(simulator);
2220 const int nseg = this->numberOfSegments();
2221 for (
int seg = 0; seg < nseg; ++seg) {
2223 const Scalar seg_pressure = getValue(this->primary_variables_.getSegmentPressure(seg));
2224 for (
const int perf : this->segments_.perforations()[seg]) {
2225 const int local_perf_index = this->pw_info_.activePerfToLocalPerf(perf);
2226 if (local_perf_index < 0)
2229 const int cell_idx = this->well_cells_[local_perf_index];
2230 const auto& int_quants = simulator.model().intensiveQuantities(cell_idx, 0);
2231 std::vector<Scalar> mob(this->num_conservation_quantities_, 0.0);
2232 getMobility(simulator, local_perf_index, mob, deferred_logger);
2233 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(int_quants, cell_idx);
2234 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
2235 const std::vector<Scalar> Tw = this->wellIndex(local_perf_index, int_quants, trans_mult, wellstate_nupcol);
2236 std::vector<Scalar> cq_s(this->num_conservation_quantities_, 0.0);
2239 computePerfRate(int_quants, mob, Tw, seg, perf, seg_pressure,
2240 allow_cf, cq_s, perf_press, perf_rates, deferred_logger);
2241 for (
int comp = 0; comp < this->num_conservation_quantities_; ++comp) {
2242 well_q_s[comp] += cq_s[comp];
2246 const auto& comm = this->parallel_well_info_.communication();
2247 if (comm.size() > 1)
2249 comm.sum(well_q_s.data(), well_q_s.size());
2255 template <
typename TypeTag>
2256 std::vector<typename MultisegmentWell<TypeTag>::Scalar>
2260 const int num_seg = this->numberOfSegments();
2261 constexpr int num_eq = MSWEval::numWellEq;
2262 std::vector<Scalar> retval(num_seg * num_eq);
2263 for (
int ii = 0; ii < num_seg; ++ii) {
2264 const auto& pv = this->primary_variables_.value(ii);
2265 std::copy(pv.begin(), pv.end(), retval.begin() + ii * num_eq);
2273 template <
typename TypeTag>
2278 const int num_seg = this->numberOfSegments();
2279 constexpr int num_eq = MSWEval::numWellEq;
2280 std::array<Scalar, num_eq> tmp;
2281 for (
int ii = 0; ii < num_seg; ++ii) {
2282 const auto start = it + ii * num_eq;
2283 std::copy(start, start + num_eq, tmp.begin());
2284 this->primary_variables_.setValue(ii, tmp);
2286 return num_seg * num_eq;
2289 template <
typename TypeTag>
2293 Scalar fsTemperature = 0.0;
2294 using SaltConcType =
typename std::decay<
decltype(std::declval<
decltype(simulator.model().intensiveQuantities(0, 0).fluidState())>().saltConcentration())>::type;
2295 SaltConcType fsSaltConcentration{};
2296 int pvt_region_index = 0;
2299 if (this->well_cells_.size() > 0) {
2302 const int cell_idx = this->well_cells_[0];
2303 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
2304 const auto& fs = intQuants.fluidState();
2306 fsTemperature = fs.temperature(FluidSystem::oilPhaseIdx).value();
2307 fsSaltConcentration = fs.saltConcentration();
2308 pvt_region_index = fs.pvtRegionIndex();
2311 auto info = std::make_tuple(fsTemperature, fsSaltConcentration, pvt_region_index);
2315 return this->parallel_well_info_.communication().size() == 1 ? info : this->pw_info_.broadcastFirstPerforationValue(info);
#define OPM_DEFLOG_PROBLEM(Exception, message, deferred_logger)
Definition: DeferredLoggingErrorHelpers.hpp:61
Definition: ConvergenceReport.hpp:38
Definition: DeferredLogger.hpp:57
void warning(const std::string &tag, const std::string &message)
void problem(const std::string &tag, const std::string &message)
void debug(const std::string &tag, const std::string &message)
Definition: GroupState.hpp:41
Class handling assemble of the equation system for MultisegmentWell.
Definition: MultisegmentWellAssemble.hpp:44
typename PrimaryVariables::EvalWell EvalWell
Definition: MultisegmentWellEval.hpp:66
PrimaryVariables primary_variables_
The primary variables.
Definition: MultisegmentWellEval.hpp:143
void scaleSegmentRatesWithWellRates(const std::vector< std::vector< int > > &segment_inlets, const std::vector< std::vector< int > > &segment_perforations, WellState< Scalar, IndexTraits > &well_state) const
void scaleSegmentPressuresWithBhp(WellState< Scalar, IndexTraits > &well_state) const
Definition: MultisegmentWell.hpp:38
bool computeWellPotentialsImplicit(const Simulator &simulator, const WellStateType &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:531
void updateWellState(const Simulator &simulator, const BVectorWell &dwells, WellStateType &well_state, DeferredLogger &deferred_logger, const Scalar relaxation_factor=1.0)
Definition: MultisegmentWell_impl.hpp:705
void updateWaterThroughput(const double dt, WellStateType &well_state) const override
Definition: MultisegmentWell_impl.hpp:2027
void addWellPressureEquations(PressureMatrix &mat, const BVector &x, const int pressureVarIndex, const bool use_well_weights, const WellStateType &well_state) const override
Definition: MultisegmentWell_impl.hpp:877
std::optional< Scalar > computeBhpAtThpLimitProdWithAlq(const Simulator &simulator, const SummaryState &summary_state, const Scalar alq_value, DeferredLogger &deferred_logger, bool iterate_if_no_solution) const override
Definition: MultisegmentWell_impl.hpp:2076
Scalar connectionDensity(const int globalConnIdx, const int openConnIdx) const override
Definition: MultisegmentWell_impl.hpp:842
void addWellContributions(SparseMatrixAdapter &jacobian) const override
Definition: MultisegmentWell_impl.hpp:864
bool iterateWellEqWithSwitching(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger, const bool fixed_control=false, const bool fixed_status=false) override
Definition: MultisegmentWell_impl.hpp:1609
void solveEqAndUpdateWellState(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:601
Scalar getRefDensity() const override
Definition: MultisegmentWell_impl.hpp:1191
std::vector< Scalar > getPrimaryVars() const override
Definition: MultisegmentWell_impl.hpp:2258
void computeWellRatesWithBhpIterations(const Simulator &simulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:405
std::vector< Scalar > computeCurrentWellRates(const Simulator &simulator, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:2214
void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition: MultisegmentWell_impl.hpp:230
MultisegmentWell(const Well &well, const ParallelWellInfo< Scalar > &pw_info, const int time_step, const ModelParameters ¶m, const RateConverterType &rate_converter, const int pvtRegionIdx, const int num_conservation_quantities, const int num_phases, const int index_of_well, const std::vector< PerforationData< Scalar > > &perf_data)
Definition: MultisegmentWell_impl.hpp:62
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:82
void checkOperabilityUnderBHPLimit(const WellStateType &well_state, const Simulator &ebos_simulator, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1199
void getMobility(const Simulator &simulator, const int local_perf_index, std::vector< Value > &mob, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:1151
void computePerfRate(const IntensiveQuantities &int_quants, const std::vector< Value > &mob_perfcells, const std::vector< Scalar > &Tw, const int seg, const int perf, const Value &segment_pressure, const bool &allow_cf, std::vector< Value > &cq_s, Value &perf_press, PerforationRates< Scalar > &perf_rates, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:1046
void updateIPRImplicit(const Simulator &simulator, WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1364
std::optional< Scalar > computeBhpAtThpLimitProd(const WellStateType &well_state, const Simulator &ebos_simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:2058
bool openCrossFlowAvoidSingularity(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:1965
void computeSegmentFluidProperties(const Simulator &simulator, DeferredLogger &deferred_logger)
Definition: MultisegmentWell_impl.hpp:1119
void checkOperabilityUnderTHPLimit(const Simulator &ebos_simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1433
bool iterateWellEqWithControl(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1485
int setPrimaryVars(typename std::vector< Scalar >::const_iterator it) override
Definition: MultisegmentWell_impl.hpp:2276
void computeWellRatesWithBhp(const Simulator &simulator, const Scalar &bhp, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:358
EvalWell getSegmentSurfaceVolume(const Simulator &simulator, const int seg_idx) const
Definition: MultisegmentWell_impl.hpp:2038
void scaleSegmentRatesAndPressure(WellStateType &well_state) const override
updating the segment pressure and rates based the current bhp and well rates
Definition: MultisegmentWell_impl.hpp:173
bool allDrawDownWrongDirection(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:1974
int debug_cost_counter_
Definition: MultisegmentWell.hpp:179
ConvergenceReport getWellConvergence(const Simulator &simulator, const WellStateType &well_state, const std::vector< Scalar > &B_avg, DeferredLogger &deferred_logger, const bool relax_tolerance) const override
check whether the well equations get converged for this well
Definition: MultisegmentWell_impl.hpp:204
void updateProductivityIndex(const Simulator &simulator, const WellProdIndexCalculator< Scalar > &wellPICalc, WellStateType &well_state, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:760
static constexpr bool has_polymer
Definition: WellInterface.hpp:113
static constexpr bool has_solvent
Definition: WellInterface.hpp:111
void computePerfCellPressDiffs(const Simulator &simulator)
Definition: MultisegmentWell_impl.hpp:630
void computeWellPotentials(const Simulator &simulator, const WellStateType &well_state, std::vector< Scalar > &well_potentials, DeferredLogger &deferred_logger) override
computing the well potentials for group control
Definition: MultisegmentWell_impl.hpp:297
void recoverWellSolutionAndUpdateWellState(const Simulator &simulator, const BVector &x, WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:265
void updateWellStateWithTarget(const Simulator &simulator, const GroupState< Scalar > &group_state, WellStateType &well_state, DeferredLogger &deferred_logger) const override
updating the well state based the current control mode
Definition: MultisegmentWell_impl.hpp:184
void updateIPR(const Simulator &ebos_simulator, DeferredLogger &deferred_logger) const override
Definition: MultisegmentWell_impl.hpp:1265
std::tuple< Scalar, typename std::decay< decltype(std::declval< decltype(std::declval< const Simulator & >().model().intensiveQuantities(0, 0).fluidState())>().saltConcentration())>::type, int > FSInfo
Definition: MultisegmentWell.hpp:73
FSInfo getFirstPerforationFluidStateInfo(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:2291
void computeWellRatesAtBhpLimit(const Simulator &simulator, std::vector< Scalar > &well_flux, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:342
void init(const std::vector< Scalar > &depth_arg, const Scalar gravity_arg, const std::vector< Scalar > &B_avg, const bool changed_to_open_this_step) override
Definition: MultisegmentWell_impl.hpp:122
void assembleWellEqWithoutIteration(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:1784
Scalar maxPerfPress(const Simulator &simulator) const
Definition: MultisegmentWell_impl.hpp:2186
void updatePrimaryVariables(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:157
std::optional< Scalar > computeBhpAtThpLimitInj(const Simulator &ebos_simulator, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:2132
void calculateExplicitQuantities(const Simulator &simulator, const WellStateType &well_state, DeferredLogger &deferred_logger) override
Definition: MultisegmentWell_impl.hpp:744
std::vector< Scalar > computeWellPotentialWithTHP(const WellStateType &well_state, const Simulator &simulator, DeferredLogger &deferred_logger) const
Definition: MultisegmentWell_impl.hpp:480
void computeInitialSegmentFluids(const Simulator &simulator)
Definition: MultisegmentWell_impl.hpp:687
EvalWell getQs(const int comp_idx) const
Returns scaled rate for a component.
Class encapsulating some information about parallel wells.
Definition: ParallelWellInfo.hpp:198
Class for computing BHP limits.
Definition: WellBhpThpCalculator.hpp:41
Scalar calculateThpFromBhp(const std::vector< Scalar > &rates, const Scalar bhp, const Scalar rho, const std::optional< Scalar > &alq, const Scalar thp_limit, DeferredLogger &deferred_logger) const
Calculates THP from BHP.
Scalar mostStrictBhpFromBhpLimits(const SummaryState &summaryState) const
Obtain the most strict BHP from BHP limits.
Well well_ecl_
Definition: WellInterfaceGeneric.hpp:304
void onlyKeepBHPandTHPcontrols(const SummaryState &summary_state, WellStateType &well_state, Well::InjectionControls &inj_controls, Well::ProductionControls &prod_controls) const
FluidSystem::Scalar rsRvInj() const
void resetDampening()
Definition: WellInterfaceGeneric.hpp:247
std::pair< bool, bool > computeWellPotentials(std::vector< Scalar > &well_potentials, const WellStateType &well_state)
Definition: WellInterfaceIndices.hpp:34
Definition: WellInterface.hpp:76
static constexpr bool has_brine
Definition: WellInterface.hpp:119
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: WellInterface.hpp:81
typename WellInterfaceFluidSystem< FluidSystem >::RateConverterType RateConverterType
Definition: WellInterface.hpp:103
Dune::BCRSMatrix< Opm::MatrixBlock< Scalar, 1, 1 > > PressureMatrix
Definition: WellInterface.hpp:97
void getMobility(const Simulator &simulator, const int local_perf_index, std::vector< Value > &mob, Callback &extendEval, DeferredLogger &deferred_logger) const
Definition: WellInterface_impl.hpp:1960
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: WellInterface.hpp:86
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: WellInterface.hpp:82
static constexpr bool has_watVapor
Definition: WellInterface.hpp:120
Dune::BlockVector< VectorBlockType > BVector
Definition: WellInterface.hpp:96
typename Base::ModelParameters ModelParameters
Definition: WellInterface.hpp:109
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: WellInterface.hpp:83
static constexpr bool has_foam
Definition: WellInterface.hpp:118
static constexpr bool has_micp
Definition: WellInterface.hpp:124
typename Base::Eval Eval
Definition: WellInterface.hpp:95
static constexpr bool has_energy
Definition: WellInterface.hpp:114
bool solveWellWithOperabilityCheck(const Simulator &simulator, const double dt, const Well::InjectionControls &inj_controls, const Well::ProductionControls &prod_controls, WellStateType &well_state, const GroupState< Scalar > &group_state, DeferredLogger &deferred_logger)
Definition: WellInterface_impl.hpp:566
GetPropType< TypeTag, Properties::Indices > Indices
Definition: WellInterface.hpp:85
bool thp_update_iterations
Definition: WellInterface.hpp:372
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: WellInterface.hpp:88
Definition: WellProdIndexCalculator.hpp:37
Scalar connectionProdIndStandard(const std::size_t connIdx, const Scalar connMobility) const
Definition: WellState.hpp:66
const SingleWellState< Scalar, IndexTraits > & well(std::size_t well_index) const
Definition: WellState.hpp:290
@ NONE
Definition: DeferredLogger.hpp:46
Definition: blackoilbioeffectsmodules.hh:43
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
Static data associated with a well perforation.
Definition: PerforationData.hpp:30
Definition: PerforationData.hpp:41
Scalar dis_gas
Definition: PerforationData.hpp:42
Scalar vap_oil
Definition: PerforationData.hpp:44