22#ifndef OPM_MULTISEGMENTWELL_IMPL_HEADER_INCLUDED
23#define OPM_MULTISEGMENTWELL_IMPL_HEADER_INCLUDED
25#ifndef OPM_MULTISEGMENTWELL_HEADER_INCLUDED
27#include <opm/simulators/wells/MultisegmentWell.hpp>
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>
43#include <opm/simulators/wells/MultisegmentWellAssemble.hpp>
44#include <opm/simulators/wells/WellBhpThpCalculator.hpp>
45#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
51#if COMPILE_GPU_BRIDGE && (HAVE_CUDA || HAVE_OPENCL)
52#include <opm/simulators/linalg/gpubridge/WellContributions.hpp>
59 template <
typename TypeTag>
60 MultisegmentWell<TypeTag>::
61 MultisegmentWell(
const Well& well,
64 const ModelParameters& param,
66 const int pvtRegionIdx,
74 , segment_fluid_initial_(
this->numberOfSegments(), std::
vector<Scalar>(
this->num_components_, 0.0))
77 if constexpr (has_solvent) {
78 OPM_THROW(std::runtime_error,
"solvent is not supported by multisegment well yet");
81 if constexpr (has_polymer) {
82 OPM_THROW(std::runtime_error,
"polymer is not supported by multisegment well yet");
85 if constexpr (Base::has_energy) {
86 OPM_THROW(std::runtime_error,
"energy is not supported by multisegment well yet");
89 if constexpr (Base::has_foam) {
90 OPM_THROW(std::runtime_error,
"foam is not supported by multisegment well yet");
93 if constexpr (Base::has_brine) {
94 OPM_THROW(std::runtime_error,
"brine is not supported by multisegment well yet");
97 if constexpr (Base::has_watVapor) {
98 OPM_THROW(std::runtime_error,
"water evaporation is not supported by multisegment well yet");
101 if(this->rsRvInj() > 0) {
102 OPM_THROW(std::runtime_error,
103 "dissolved gas/ vapporized oil in injected oil/gas not supported by multisegment well yet."
104 " \n See (WCONINJE item 10 / WCONHIST item 8)");
107 this->thp_update_iterations =
true;
114 template <
typename TypeTag>
116 MultisegmentWell<TypeTag>::
117 init(
const PhaseUsage* phase_usage_arg,
118 const std::vector<Scalar>& depth_arg,
119 const Scalar gravity_arg,
120 const std::vector< Scalar >& B_avg,
121 const bool changed_to_open_this_step)
123 Base::init(phase_usage_arg, depth_arg, gravity_arg, B_avg, changed_to_open_this_step);
135 this->initMatrixAndVectors();
138 for (
int local_perf_index = 0; local_perf_index < this->number_of_local_perforations_; ++local_perf_index) {
140 const int cell_idx = this->well_cells_[local_perf_index];
142 this->cell_perforation_depth_diffs_[local_perf_index] = depth_arg[cell_idx] - this->perf_depth_[this->pw_info_.localToGlobal(local_perf_index)];
150 template <
typename TypeTag>
152 MultisegmentWell<TypeTag>::
153 updatePrimaryVariables(
const Simulator& simulator,
154 const WellState<Scalar>& well_state,
155 DeferredLogger& deferred_logger)
157 const bool stop_or_zero_rate_target = this->stoppedOrZeroRateTarget(simulator, well_state, deferred_logger);
158 this->primary_variables_.update(well_state, stop_or_zero_rate_target);
166 template <
typename TypeTag>
174 Base::updateWellStateWithTarget(simulator, group_state, well_state, deferred_logger);
177 this->scaleSegmentRatesWithWellRates(this->segments_.inlets(),
178 this->segments_.perforations(),
180 this->scaleSegmentPressuresWithBhp(well_state);
187 template <
typename TypeTag>
192 const std::vector<Scalar>& B_avg,
194 const bool relax_tolerance)
const
196 return this->MSWEval::getWellConvergence(well_state,
199 this->param_.max_residual_allowed_,
200 this->param_.tolerance_wells_,
201 this->param_.relaxed_tolerance_flow_well_,
202 this->param_.tolerance_pressure_ms_wells_,
203 this->param_.relaxed_tolerance_pressure_ms_well_,
205 this->wellIsStopped());
213 template <
typename TypeTag>
216 apply(
const BVector& x, BVector& Ax)
const
218 if (!this->isOperableAndSolvable() && !this->wellIsStopped()) {
222 if (this->param_.matrix_add_well_contributions_) {
227 this->linSys_.
apply(x, Ax);
234 template <
typename TypeTag>
237 apply(BVector& r)
const
239 if (!this->isOperableAndSolvable() && !this->wellIsStopped()) {
243 this->linSys_.
apply(r);
248 template <
typename TypeTag>
256 if (!this->isOperableAndSolvable() && !this->wellIsStopped()) {
261 this->linSys_.recoverSolutionWell(x, xw);
262 updateWellState(simulator, xw, well_state, deferred_logger);
269 template <
typename TypeTag>
274 std::vector<Scalar>& well_potentials,
277 const auto [compute_potential, bhp_controlled_well] =
280 if (!compute_potential) {
284 debug_cost_counter_ = 0;
285 bool converged_implicit =
false;
286 if (this->param_.local_well_solver_control_switching_) {
287 converged_implicit = computeWellPotentialsImplicit(simulator, well_state, well_potentials, deferred_logger);
288 if (!converged_implicit) {
289 deferred_logger.debug(
"Implicit potential calculations failed for well "
290 + this->name() +
", reverting to original aproach.");
293 if (!converged_implicit) {
295 const auto& summaryState = simulator.vanguard().summaryState();
296 if (!Base::wellHasTHPConstraints(summaryState) || bhp_controlled_well) {
297 computeWellRatesAtBhpLimit(simulator, well_potentials, deferred_logger);
299 well_potentials = computeWellPotentialWithTHP(
300 well_state, simulator, deferred_logger);
303 deferred_logger.debug(
"Cost in iterations of finding well potential for well "
304 + this->name() +
": " + std::to_string(debug_cost_counter_));
306 this->checkNegativeWellPotentials(well_potentials,
307 this->param_.check_well_operability_,
314 template<
typename TypeTag>
318 std::vector<Scalar>& well_flux,
321 if (this->well_ecl_.isInjector()) {
322 const auto controls = this->well_ecl_.injectionControls(simulator.
vanguard().summaryState());
323 computeWellRatesWithBhpIterations(simulator, controls.bhp_limit, well_flux, deferred_logger);
325 const auto controls = this->well_ecl_.productionControls(simulator.
vanguard().summaryState());
326 computeWellRatesWithBhpIterations(simulator, controls.bhp_limit, well_flux, deferred_logger);
330 template<
typename TypeTag>
332 MultisegmentWell<TypeTag>::
333 computeWellRatesWithBhp(
const Simulator& simulator,
335 std::vector<Scalar>& well_flux,
336 DeferredLogger& deferred_logger)
const
338 const int np = this->number_of_phases_;
340 well_flux.resize(np, 0.0);
341 const bool allow_cf = this->getAllowCrossFlow();
342 const int nseg = this->numberOfSegments();
343 const WellState<Scalar>& well_state = simulator.problem().wellModel().wellState();
344 const auto& ws = well_state.well(this->indexOfWell());
345 auto segments_copy = ws.segments;
346 segments_copy.scale_pressure(bhp);
347 const auto& segment_pressure = segments_copy.pressure;
348 for (
int seg = 0; seg < nseg; ++seg) {
349 for (
const int perf : this->segments_.perforations()[seg]) {
350 const int local_perf_index = this->pw_info_.globalToLocal(perf);
351 if (local_perf_index < 0)
353 const int cell_idx = this->well_cells_[local_perf_index];
354 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
356 std::vector<Scalar> mob(this->num_components_, 0.);
357 getMobility(simulator, local_perf_index, mob, deferred_logger);
358 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(intQuants, cell_idx);
359 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
360 const std::vector<Scalar> Tw = this->wellIndex(local_perf_index, intQuants, trans_mult, wellstate_nupcol);
361 const Scalar seg_pressure = segment_pressure[seg];
362 std::vector<Scalar> cq_s(this->num_components_, 0.);
363 Scalar perf_press = 0.0;
364 PerforationRates<Scalar> perf_rates;
365 computePerfRate(intQuants, mob, Tw, seg, perf, seg_pressure,
366 allow_cf, cq_s, perf_press, perf_rates, deferred_logger);
368 for(
int p = 0; p < np; ++p) {
369 well_flux[this->modelCompIdxToFlowCompIdx(p)] += cq_s[p];
373 this->parallel_well_info_.communication().sum(well_flux.data(), well_flux.size());
377 template<
typename TypeTag>
379 MultisegmentWell<TypeTag>::
380 computeWellRatesWithBhpIterations(
const Simulator& simulator,
382 std::vector<Scalar>& well_flux,
383 DeferredLogger& deferred_logger)
const
388 MultisegmentWell<TypeTag> well_copy(*
this);
389 well_copy.resetDampening();
391 well_copy.debug_cost_counter_ = 0;
394 WellState<Scalar> well_state_copy = simulator.problem().wellModel().wellState();
395 const auto& group_state = simulator.problem().wellModel().groupState();
396 auto& ws = well_state_copy.well(this->index_of_well_);
399 const auto& summary_state = simulator.vanguard().summaryState();
400 auto inj_controls = well_copy.well_ecl_.isInjector()
401 ? well_copy.well_ecl_.injectionControls(summary_state)
402 : Well::InjectionControls(0);
403 auto prod_controls = well_copy.well_ecl_.isProducer()
404 ? well_copy.well_ecl_.productionControls(summary_state) :
405 Well::ProductionControls(0);
408 if (well_copy.well_ecl_.isInjector()) {
409 inj_controls.bhp_limit = bhp;
410 ws.injection_cmode = Well::InjectorCMode::BHP;
412 prod_controls.bhp_limit = bhp;
413 ws.production_cmode = Well::ProducerCMode::BHP;
416 well_copy.scaleSegmentPressuresWithBhp(well_state_copy);
419 const int np = this->number_of_phases_;
421 for (
int phase = 0; phase < np; ++phase){
422 trivial = trivial && (ws.well_potentials[phase] == 0.0) ;
425 const Scalar sign = well_copy.well_ecl_.isInjector() ? 1.0 : -1.0;
426 for (
int phase = 0; phase < np; ++phase) {
427 ws.surface_rates[phase] = sign * ws.well_potentials[phase];
430 well_copy.scaleSegmentRatesWithWellRates(this->segments_.inlets(),
431 this->segments_.perforations(),
434 well_copy.calculateExplicitQuantities(simulator, well_state_copy, deferred_logger);
435 const double dt = simulator.timeStepSize();
437 well_copy.iterateWellEqWithControl(simulator, dt, inj_controls, prod_controls, well_state_copy, group_state,
442 well_flux.resize(np, 0.0);
443 for (
int compIdx = 0; compIdx < this->num_components_; ++compIdx) {
444 const EvalWell rate = well_copy.primary_variables_.getQs(compIdx);
445 well_flux[this->modelCompIdxToFlowCompIdx(compIdx)] = rate.value();
447 debug_cost_counter_ += well_copy.debug_cost_counter_;
452 template<
typename TypeTag>
453 std::vector<typename MultisegmentWell<TypeTag>::Scalar>
454 MultisegmentWell<TypeTag>::
455 computeWellPotentialWithTHP(
const WellState<Scalar>& well_state,
456 const Simulator& simulator,
457 DeferredLogger& deferred_logger)
const
459 std::vector<Scalar> potentials(this->number_of_phases_, 0.0);
460 const auto& summary_state = simulator.vanguard().summaryState();
462 const auto& well = this->well_ecl_;
463 if (well.isInjector()) {
464 auto bhp_at_thp_limit = computeBhpAtThpLimitInj(simulator, summary_state, deferred_logger);
465 if (bhp_at_thp_limit) {
466 const auto& controls = well.injectionControls(summary_state);
467 const Scalar bhp = std::min(*bhp_at_thp_limit,
468 static_cast<Scalar
>(controls.bhp_limit));
469 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
470 deferred_logger.debug(
"Converged thp based potential calculation for well "
471 + this->name() +
", at bhp = " + std::to_string(bhp));
473 deferred_logger.warning(
"FAILURE_GETTING_CONVERGED_POTENTIAL",
474 "Failed in getting converged thp based potential calculation for well "
475 + this->name() +
". Instead the bhp based value is used");
476 const auto& controls = well.injectionControls(summary_state);
477 const Scalar bhp = controls.bhp_limit;
478 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
481 auto bhp_at_thp_limit = computeBhpAtThpLimitProd(
482 well_state, simulator, summary_state, deferred_logger);
483 if (bhp_at_thp_limit) {
484 const auto& controls = well.productionControls(summary_state);
485 const Scalar bhp = std::max(*bhp_at_thp_limit,
486 static_cast<Scalar
>(controls.bhp_limit));
487 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
488 deferred_logger.debug(
"Converged thp based potential calculation for well "
489 + this->name() +
", at bhp = " + std::to_string(bhp));
491 deferred_logger.warning(
"FAILURE_GETTING_CONVERGED_POTENTIAL",
492 "Failed in getting converged thp based potential calculation for well "
493 + this->name() +
". Instead the bhp based value is used");
494 const auto& controls = well.productionControls(summary_state);
495 const Scalar bhp = controls.bhp_limit;
496 computeWellRatesWithBhpIterations(simulator, bhp, potentials, deferred_logger);
503 template<
typename TypeTag>
505 MultisegmentWell<TypeTag>::
506 computeWellPotentialsImplicit(
const Simulator& simulator,
507 const WellState<Scalar>& well_state,
508 std::vector<Scalar>& well_potentials,
509 DeferredLogger& deferred_logger)
const
514 MultisegmentWell<TypeTag> well_copy(*
this);
515 well_copy.debug_cost_counter_ = 0;
518 WellState<Scalar> well_state_copy = well_state;
519 const auto& group_state = simulator.problem().wellModel().groupState();
520 auto& ws = well_state_copy.well(this->index_of_well_);
523 const auto& summary_state = simulator.vanguard().summaryState();
524 auto inj_controls = well_copy.well_ecl_.isInjector()
525 ? well_copy.well_ecl_.injectionControls(summary_state)
526 : Well::InjectionControls(0);
527 auto prod_controls = well_copy.well_ecl_.isProducer()
528 ? well_copy.well_ecl_.productionControls(summary_state)
529 : Well::ProductionControls(0);
532 well_copy.prepareForPotentialCalculations(summary_state, well_state_copy, inj_controls, prod_controls);
534 well_copy.scaleSegmentPressuresWithBhp(well_state_copy);
537 const int np = this->number_of_phases_;
539 for (
int phase = 0; phase < np; ++phase){
540 trivial = trivial && (ws.well_potentials[phase] == 0.0) ;
543 const Scalar sign = well_copy.well_ecl_.isInjector() ? 1.0 : -1.0;
544 for (
int phase = 0; phase < np; ++phase) {
545 ws.surface_rates[phase] = sign * ws.well_potentials[phase];
548 well_copy.scaleSegmentRatesWithWellRates(this->segments_.inlets(),
549 this->segments_.perforations(),
552 well_copy.calculateExplicitQuantities(simulator, well_state_copy, deferred_logger);
553 const double dt = simulator.timeStepSize();
555 bool converged =
false;
556 if (this->well_ecl_.isProducer() && this->wellHasTHPConstraints(summary_state)) {
557 converged = well_copy.solveWellWithTHPConstraint(simulator, dt, inj_controls, prod_controls, well_state_copy, group_state, deferred_logger);
559 converged = well_copy.iterateWellEqWithSwitching(simulator, dt, inj_controls, prod_controls, well_state_copy, group_state, deferred_logger);
563 well_potentials.clear();
564 well_potentials.resize(np, 0.0);
565 for (
int compIdx = 0; compIdx < this->num_components_; ++compIdx) {
566 const EvalWell rate = well_copy.primary_variables_.getQs(compIdx);
567 well_potentials[this->modelCompIdxToFlowCompIdx(compIdx)] = rate.value();
569 debug_cost_counter_ += well_copy.debug_cost_counter_;
573 template <
typename TypeTag>
575 MultisegmentWell<TypeTag>::
576 solveEqAndUpdateWellState(
const Simulator& simulator,
577 WellState<Scalar>& well_state,
578 DeferredLogger& deferred_logger)
580 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return;
585 const BVectorWell dx_well = this->linSys_.solve();
587 updateWellState(simulator, dx_well, well_state, deferred_logger);
589 catch(
const NumericalProblem& exp) {
593 deferred_logger.problem(
"In MultisegmentWell::solveEqAndUpdateWellState for well "
594 + this->name() +
": "+exp.what());
603 template <
typename TypeTag>
605 MultisegmentWell<TypeTag>::
606 computePerfCellPressDiffs(
const Simulator& simulator)
610 for (
int local_perf_index = 0; local_perf_index < this->number_of_local_perforations_; ++local_perf_index) {
613 std::vector<Scalar> kr(this->number_of_phases_, 0.0);
614 std::vector<Scalar> density(this->number_of_phases_, 0.0);
616 const int cell_idx = this->well_cells_[local_perf_index];
617 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
618 const auto& fs = intQuants.fluidState();
623 if (pu.phase_used[Water]) {
624 const int water_pos = pu.phase_pos[Water];
625 kr[water_pos] = intQuants.relativePermeability(FluidSystem::waterPhaseIdx).value();
626 sum_kr += kr[water_pos];
627 density[water_pos] = fs.density(FluidSystem::waterPhaseIdx).value();
630 if (pu.phase_used[Oil]) {
631 const int oil_pos = pu.phase_pos[Oil];
632 kr[oil_pos] = intQuants.relativePermeability(FluidSystem::oilPhaseIdx).value();
633 sum_kr += kr[oil_pos];
634 density[oil_pos] = fs.density(FluidSystem::oilPhaseIdx).value();
637 if (pu.phase_used[Gas]) {
638 const int gas_pos = pu.phase_pos[Gas];
639 kr[gas_pos] = intQuants.relativePermeability(FluidSystem::gasPhaseIdx).value();
640 sum_kr += kr[gas_pos];
641 density[gas_pos] = fs.density(FluidSystem::gasPhaseIdx).value();
644 assert(sum_kr != 0.);
647 Scalar average_density = 0.;
648 for (
int p = 0; p < this->number_of_phases_; ++p) {
649 average_density += kr[p] * density[p];
651 average_density /= sum_kr;
653 this->cell_perforation_pressure_diffs_[local_perf_index] = this->gravity_ * average_density * this->cell_perforation_depth_diffs_[local_perf_index];
661 template <
typename TypeTag>
663 MultisegmentWell<TypeTag>::
664 computeInitialSegmentFluids(
const Simulator& simulator)
666 for (
int seg = 0; seg < this->numberOfSegments(); ++seg) {
668 const Scalar surface_volume = getSegmentSurfaceVolume(simulator, seg).value();
669 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
670 segment_fluid_initial_[seg][comp_idx] = surface_volume * this->primary_variables_.surfaceVolumeFraction(seg, comp_idx).value();
679 template <
typename TypeTag>
681 MultisegmentWell<TypeTag>::
682 updateWellState(
const Simulator& simulator,
683 const BVectorWell& dwells,
684 WellState<Scalar>& well_state,
685 DeferredLogger& deferred_logger,
686 const Scalar relaxation_factor)
688 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return;
690 const Scalar dFLimit = this->param_.dwell_fraction_max_;
691 const Scalar max_pressure_change = this->param_.max_pressure_change_ms_wells_;
692 const bool stop_or_zero_rate_target =
693 this->stoppedOrZeroRateTarget(simulator, well_state, deferred_logger);
694 this->primary_variables_.updateNewton(dwells,
697 stop_or_zero_rate_target,
698 max_pressure_change);
700 const auto& summary_state = simulator.vanguard().summaryState();
701 this->primary_variables_.copyToWellState(*
this, getRefDensity(),
707 auto& ws = well_state.well(this->index_of_well_);
708 this->segments_.copyPhaseDensities(ws.pu, ws.segments);
711 Base::calculateReservoirRates(simulator.vanguard().eclState().runspec().co2Storage(), well_state.well(this->index_of_well_));
718 template <
typename TypeTag>
720 MultisegmentWell<TypeTag>::
721 calculateExplicitQuantities(
const Simulator& simulator,
722 const WellState<Scalar>& well_state,
723 DeferredLogger& deferred_logger)
725 updatePrimaryVariables(simulator, well_state, deferred_logger);
726 computePerfCellPressDiffs(simulator);
727 computeInitialSegmentFluids(simulator);
734 template<
typename TypeTag>
736 MultisegmentWell<TypeTag>::
737 updateProductivityIndex(
const Simulator& simulator,
738 const WellProdIndexCalculator<Scalar>& wellPICalc,
739 WellState<Scalar>& well_state,
740 DeferredLogger& deferred_logger)
const
742 auto fluidState = [&simulator,
this](
const int local_perf_index)
744 const auto cell_idx = this->well_cells_[local_perf_index];
745 return simulator.model()
746 .intensiveQuantities(cell_idx, 0).fluidState();
749 const int np = this->number_of_phases_;
750 auto setToZero = [np](Scalar* x) ->
void
752 std::fill_n(x, np, 0.0);
755 auto addVector = [np](
const Scalar* src, Scalar* dest) ->
void
757 std::transform(src, src + np, dest, dest, std::plus<>{});
760 auto& ws = well_state.well(this->index_of_well_);
761 auto& perf_data = ws.perf_data;
762 auto* connPI = perf_data.prod_index.data();
763 auto* wellPI = ws.productivity_index.data();
767 const auto preferred_phase = this->well_ecl_.getPreferredPhase();
768 auto subsetPerfID = 0;
770 for (
const auto& perf : *this->perf_data_){
771 auto allPerfID = perf.ecl_index;
773 auto connPICalc = [&wellPICalc, allPerfID](
const Scalar mobility) -> Scalar
775 return wellPICalc.connectionProdIndStandard(allPerfID, mobility);
778 std::vector<Scalar> mob(this->num_components_, 0.0);
783 getMobility(simulator,
static_cast<int>(subsetPerfID), mob, deferred_logger);
785 const auto& fs = fluidState(subsetPerfID);
788 if (this->isInjector()) {
789 this->computeConnLevelInjInd(fs, preferred_phase, connPICalc,
790 mob, connPI, deferred_logger);
793 this->computeConnLevelProdInd(fs, connPICalc, mob, connPI);
796 addVector(connPI, wellPI);
803 const auto& comm = this->parallel_well_info_.communication();
804 if (comm.size() > 1) {
805 comm.sum(wellPI, np);
808 assert (
static_cast<int>(subsetPerfID) == this->number_of_local_perforations_ &&
809 "Internal logic error in processing connections for PI/II");
816 template<
typename TypeTag>
817 typename MultisegmentWell<TypeTag>::Scalar
818 MultisegmentWell<TypeTag>::
819 connectionDensity(
const int globalConnIdx,
820 [[maybe_unused]]
const int openConnIdx)
const
825 const auto segNum = this->wellEcl()
826 .getConnections()[globalConnIdx].segment();
828 const auto segIdx = this->wellEcl()
829 .getSegments().segmentNumberToIndex(segNum);
831 return this->segments_.density(segIdx).value();
838 template<
typename TypeTag>
840 MultisegmentWell<TypeTag>::
841 addWellContributions(SparseMatrixAdapter& jacobian)
const
843 this->linSys_.extract(jacobian);
847 template<
typename TypeTag>
849 MultisegmentWell<TypeTag>::
850 addWellPressureEquations(PressureMatrix& jacobian,
851 const BVector& weights,
852 const int pressureVarIndex,
853 const bool use_well_weights,
854 const WellState<Scalar>& well_state)
const
857 this->linSys_.extractCPRPressureMatrix(jacobian,
867 template<
typename TypeTag>
868 template<
class Value>
870 MultisegmentWell<TypeTag>::
871 computePerfRate(
const Value& pressure_cell,
874 const std::vector<Value>& b_perfcells,
875 const std::vector<Value>& mob_perfcells,
876 const std::vector<Scalar>& Tw,
878 const Value& segment_pressure,
879 const Value& segment_density,
880 const bool& allow_cf,
881 const std::vector<Value>& cmix_s,
882 std::vector<Value>& cq_s,
884 PerforationRates<Scalar>& perf_rates,
885 DeferredLogger& deferred_logger)
const
887 const int local_perf_index = this->pw_info_.globalToLocal(perf);
888 if (local_perf_index < 0)
892 const Value perf_seg_press_diff = this->gravity() * segment_density *
893 this->segments_.local_perforation_depth_diff(local_perf_index);
895 const Scalar cell_perf_press_diff = this->cell_perforation_pressure_diffs_[local_perf_index];
899 perf_press = segment_pressure + perf_seg_press_diff;
902 const Value cell_press_at_perf = pressure_cell - cell_perf_press_diff;
905 const Value drawdown = cell_press_at_perf - perf_press;
908 if (drawdown > 0.0) {
910 if (!allow_cf && this->isInjector()) {
915 for (
int comp_idx = 0; comp_idx < this->numComponents(); ++comp_idx) {
916 const Value cq_p = - Tw[comp_idx] * (mob_perfcells[comp_idx] * drawdown);
917 cq_s[comp_idx] = b_perfcells[comp_idx] * cq_p;
920 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
921 const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
922 const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
923 const Value cq_s_oil = cq_s[oilCompIdx];
924 const Value cq_s_gas = cq_s[gasCompIdx];
925 cq_s[gasCompIdx] += rs * cq_s_oil;
926 cq_s[oilCompIdx] += rv * cq_s_gas;
930 if (!allow_cf && this->isProducer()) {
935 Value total_mob = mob_perfcells[0];
936 for (
int comp_idx = 1; comp_idx < this->numComponents(); ++comp_idx) {
937 total_mob += mob_perfcells[comp_idx];
941 Value volume_ratio = 0.0;
942 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
943 const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
944 volume_ratio += cmix_s[waterCompIdx] / b_perfcells[waterCompIdx];
947 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
948 const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
949 const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
954 const Value d = 1.0 - rv * rs;
956 if (getValue(d) == 0.0) {
957 OPM_DEFLOG_PROBLEM(NumericalProblem,
958 fmt::format(
"Zero d value obtained for well {} "
959 "during flux calculation with rs {} and rv {}",
960 this->name(), rs, rv),
964 const Value tmp_oil = (cmix_s[oilCompIdx] - rv * cmix_s[gasCompIdx]) / d;
965 volume_ratio += tmp_oil / b_perfcells[oilCompIdx];
967 const Value tmp_gas = (cmix_s[gasCompIdx] - rs * cmix_s[oilCompIdx]) / d;
968 volume_ratio += tmp_gas / b_perfcells[gasCompIdx];
970 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
971 const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
972 volume_ratio += cmix_s[oilCompIdx] / b_perfcells[oilCompIdx];
974 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
975 const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
976 volume_ratio += cmix_s[gasCompIdx] / b_perfcells[gasCompIdx];
980 for (
int componentIdx = 0; componentIdx < this->numComponents(); ++componentIdx) {
981 const Value cqt_i = - Tw[componentIdx] * (total_mob * drawdown);
982 Value cqt_is = cqt_i / volume_ratio;
983 cq_s[componentIdx] = cmix_s[componentIdx] * cqt_is;
988 if (this->isProducer()) {
989 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
990 const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
991 const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
1000 const Scalar d = 1.0 - getValue(rv) * getValue(rs);
1003 perf_rates.vap_oil = getValue(rv) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d;
1006 perf_rates.dis_gas = getValue(rs) * (getValue(cq_s[oilCompIdx]) - getValue(rv) * getValue(cq_s[gasCompIdx])) / d;
1011 template <
typename TypeTag>
1012 template<
class Value>
1014 MultisegmentWell<TypeTag>::
1015 computePerfRate(
const IntensiveQuantities& int_quants,
1016 const std::vector<Value>& mob_perfcells,
1017 const std::vector<Scalar>& Tw,
1020 const Value& segment_pressure,
1021 const bool& allow_cf,
1022 std::vector<Value>& cq_s,
1024 PerforationRates<Scalar>& perf_rates,
1025 DeferredLogger& deferred_logger)
const
1028 auto obtain = [
this](
const Eval& value)
1030 if constexpr (std::is_same_v<Value, Scalar>) {
1031 static_cast<void>(
this);
1032 return getValue(value);
1034 return this->extendEval(value);
1037 auto obtainN = [](
const auto& value)
1039 if constexpr (std::is_same_v<Value, Scalar>) {
1040 return getValue(value);
1045 const auto& fs = int_quants.fluidState();
1047 const Value pressure_cell = obtain(this->getPerfCellPressure(fs));
1048 const Value rs = obtain(fs.Rs());
1049 const Value rv = obtain(fs.Rv());
1052 std::vector<Value> b_perfcells(this->num_components_, 0.0);
1054 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
1055 if (!FluidSystem::phaseIsActive(phaseIdx)) {
1059 const unsigned compIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
1060 b_perfcells[compIdx] = obtain(fs.invB(phaseIdx));
1063 std::vector<Value> cmix_s(this->numComponents(), 0.0);
1064 for (
int comp_idx = 0; comp_idx < this->numComponents(); ++comp_idx) {
1065 cmix_s[comp_idx] = obtainN(this->primary_variables_.surfaceVolumeFraction(seg, comp_idx));
1068 this->computePerfRate(pressure_cell,
1076 obtainN(this->segments_.density(seg)),
1085 template <
typename TypeTag>
1087 MultisegmentWell<TypeTag>::
1088 computeSegmentFluidProperties(
const Simulator& simulator, DeferredLogger& deferred_logger)
1097 EvalWell temperature;
1098 EvalWell saltConcentration;
1104 int pvt_region_index;
1107 const int cell_idx = this->well_cells_[0];
1108 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
1109 const auto& fs = intQuants.fluidState();
1113 auto fsTemperature = fs.temperature(FluidSystem::oilPhaseIdx).value();
1114 fsTemperature = this->pw_info_.broadcastFirstPerforationValue(fsTemperature);
1115 temperature.setValue(fsTemperature);
1117 auto fsSaltConcentration = fs.saltConcentration();
1118 fsSaltConcentration = this->pw_info_.broadcastFirstPerforationValue(fsSaltConcentration);
1119 saltConcentration = this->extendEval(fsSaltConcentration);
1121 pvt_region_index = fs.pvtRegionIndex();
1122 pvt_region_index = this->pw_info_.broadcastFirstPerforationValue(pvt_region_index);
1125 this->segments_.computeFluidProperties(temperature,
1127 this->primary_variables_,
1132 template <
typename TypeTag>
1133 template<
class Value>
1135 MultisegmentWell<TypeTag>::
1136 getMobility(
const Simulator& simulator,
1137 const int local_perf_index,
1138 std::vector<Value>& mob,
1139 DeferredLogger& deferred_logger)
const
1141 auto obtain = [
this](
const Eval& value)
1143 if constexpr (std::is_same_v<Value, Scalar>) {
1144 static_cast<void>(
this);
1145 return getValue(value);
1147 return this->extendEval(value);
1151 WellInterface<TypeTag>::getMobility(simulator, local_perf_index, mob, obtain, deferred_logger);
1153 if (this->isInjector() && this->well_ecl_.getInjMultMode() != Well::InjMultMode::NONE) {
1154 const auto perf_ecl_index = this->perforationData()[local_perf_index].ecl_index;
1155 const Connection& con = this->well_ecl_.getConnections()[perf_ecl_index];
1156 const int seg = this->segmentNumberToIndex(con.segment());
1160 const Scalar segment_pres = this->primary_variables_.getSegmentPressure(seg).value();
1161 const Scalar perf_seg_press_diff = this->gravity() * this->segments_.density(seg).value()
1162 * this->segments_.local_perforation_depth_diff(local_perf_index);
1163 const Scalar perf_press = segment_pres + perf_seg_press_diff;
1164 const Scalar multiplier = this->getInjMult(local_perf_index, segment_pres, perf_press, deferred_logger);
1165 for (std::size_t i = 0; i < mob.size(); ++i) {
1166 mob[i] *= multiplier;
1173 template<
typename TypeTag>
1174 typename MultisegmentWell<TypeTag>::Scalar
1175 MultisegmentWell<TypeTag>::
1176 getRefDensity()
const
1178 return this->segments_.getRefDensity();
1181 template<
typename TypeTag>
1183 MultisegmentWell<TypeTag>::
1184 checkOperabilityUnderBHPLimit(
const WellState<Scalar>& ,
1185 const Simulator& simulator,
1186 DeferredLogger& deferred_logger)
1188 const auto& summaryState = simulator.vanguard().summaryState();
1189 const Scalar bhp_limit = WellBhpThpCalculator(*this).mostStrictBhpFromBhpLimits(summaryState);
1192 const bool bhp_limit_not_defaulted = bhp_limit > 1.5 * unit::barsa;
1193 if ( bhp_limit_not_defaulted || !this->wellHasTHPConstraints(summaryState) ) {
1196 Scalar total_ipr_mass_rate = 0.0;
1197 for (
unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx)
1199 if (!FluidSystem::phaseIsActive(phaseIdx)) {
1203 const unsigned compIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
1204 const Scalar ipr_rate = this->ipr_a_[compIdx] - this->ipr_b_[compIdx] * bhp_limit;
1206 const Scalar rho = FluidSystem::referenceDensity( phaseIdx, Base::pvtRegionIdx() );
1207 total_ipr_mass_rate += ipr_rate * rho;
1209 if ( (this->isProducer() && total_ipr_mass_rate < 0.) || (this->isInjector() && total_ipr_mass_rate > 0.) ) {
1210 this->operability_status_.operable_under_only_bhp_limit =
false;
1214 if (this->operability_status_.operable_under_only_bhp_limit && this->wellHasTHPConstraints(summaryState)) {
1218 std::vector<Scalar> well_rates_bhp_limit;
1219 computeWellRatesWithBhp(simulator, bhp_limit, well_rates_bhp_limit, deferred_logger);
1221 const Scalar thp_limit = this->getTHPConstraint(summaryState);
1222 const Scalar thp = WellBhpThpCalculator(*this).calculateThpFromBhp(well_rates_bhp_limit,
1224 this->getRefDensity(),
1225 this->wellEcl().alq_value(summaryState),
1228 if ( (this->isProducer() && thp < thp_limit) || (this->isInjector() && thp > thp_limit) ) {
1229 this->operability_status_.obey_thp_limit_under_bhp_limit =
false;
1240 this->operability_status_.operable_under_only_bhp_limit =
true;
1241 this->operability_status_.obey_thp_limit_under_bhp_limit =
false;
1247 template<
typename TypeTag>
1249 MultisegmentWell<TypeTag>::
1250 updateIPR(
const Simulator& simulator, DeferredLogger& deferred_logger)
const
1255 std::fill(this->ipr_a_.begin(), this->ipr_a_.end(), 0.);
1256 std::fill(this->ipr_b_.begin(), this->ipr_b_.end(), 0.);
1258 const int nseg = this->numberOfSegments();
1259 std::vector<Scalar> seg_dp(nseg, 0.0);
1260 for (
int seg = 0; seg < nseg; ++seg) {
1262 const Scalar dp = this->getSegmentDp(seg,
1263 this->segments_.density(seg).value(),
1266 for (
const int perf : this->segments_.perforations()[seg]) {
1267 const int local_perf_index = this->pw_info_.globalToLocal(perf);
1268 if (local_perf_index < 0)
1270 std::vector<Scalar> mob(this->num_components_, 0.0);
1273 getMobility(simulator, local_perf_index, mob, deferred_logger);
1275 const int cell_idx = this->well_cells_[local_perf_index];
1276 const auto& int_quantities = simulator.model().intensiveQuantities(cell_idx, 0);
1277 const auto& fs = int_quantities.fluidState();
1279 const Scalar perf_seg_press_diff = this->segments_.getPressureDiffSegLocalPerf(seg, local_perf_index);
1281 const Scalar cell_perf_press_diff = this->cell_perforation_pressure_diffs_[local_perf_index];
1282 const Scalar pressure_cell = this->getPerfCellPressure(fs).value();
1285 std::vector<Scalar> b_perf(this->num_components_);
1286 for (std::size_t phase = 0; phase < FluidSystem::numPhases; ++phase) {
1287 if (!FluidSystem::phaseIsActive(phase)) {
1290 const unsigned comp_idx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phase));
1291 b_perf[comp_idx] = fs.invB(phase).value();
1295 const Scalar h_perf = cell_perf_press_diff + perf_seg_press_diff + dp;
1296 const Scalar pressure_diff = pressure_cell - h_perf;
1299 if ( (this->isProducer() && pressure_diff < 0.) || (this->isInjector() && pressure_diff > 0.) ) {
1300 deferred_logger.debug(
"CROSSFLOW_IPR",
1301 "cross flow found when updateIPR for well " + this->name());
1305 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(int_quantities, cell_idx);
1306 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
1307 const std::vector<Scalar> tw_perf = this->wellIndex(local_perf_index, int_quantities, trans_mult, wellstate_nupcol);
1308 std::vector<Scalar> ipr_a_perf(this->ipr_a_.size());
1309 std::vector<Scalar> ipr_b_perf(this->ipr_b_.size());
1310 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
1311 const Scalar tw_mob = tw_perf[comp_idx] * mob[comp_idx] * b_perf[comp_idx];
1312 ipr_a_perf[comp_idx] += tw_mob * pressure_diff;
1313 ipr_b_perf[comp_idx] += tw_mob;
1317 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
1318 const unsigned oil_comp_idx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx);
1319 const unsigned gas_comp_idx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx);
1320 const Scalar rs = (fs.Rs()).value();
1321 const Scalar rv = (fs.Rv()).value();
1323 const Scalar dis_gas_a = rs * ipr_a_perf[oil_comp_idx];
1324 const Scalar vap_oil_a = rv * ipr_a_perf[gas_comp_idx];
1326 ipr_a_perf[gas_comp_idx] += dis_gas_a;
1327 ipr_a_perf[oil_comp_idx] += vap_oil_a;
1329 const Scalar dis_gas_b = rs * ipr_b_perf[oil_comp_idx];
1330 const Scalar vap_oil_b = rv * ipr_b_perf[gas_comp_idx];
1332 ipr_b_perf[gas_comp_idx] += dis_gas_b;
1333 ipr_b_perf[oil_comp_idx] += vap_oil_b;
1336 for (std::size_t comp_idx = 0; comp_idx < ipr_a_perf.size(); ++comp_idx) {
1337 this->ipr_a_[comp_idx] += ipr_a_perf[comp_idx];
1338 this->ipr_b_[comp_idx] += ipr_b_perf[comp_idx];
1342 this->parallel_well_info_.communication().sum(this->ipr_a_.data(), this->ipr_a_.size());
1343 this->parallel_well_info_.communication().sum(this->ipr_b_.data(), this->ipr_b_.size());
1346 template<
typename TypeTag>
1348 MultisegmentWell<TypeTag>::
1349 updateIPRImplicit(
const Simulator& simulator,
1350 WellState<Scalar>& well_state,
1351 DeferredLogger& deferred_logger)
1360 auto rates = well_state.well(this->index_of_well_).surface_rates;
1362 for (std::size_t p = 0; p < rates.size(); ++p) {
1363 zero_rates &= rates[p] == 0.0;
1365 auto& ws = well_state.well(this->index_of_well_);
1367 const auto msg = fmt::format(
"updateIPRImplicit: Well {} has zero rate, IPRs might be problematic", this->name());
1368 deferred_logger.debug(msg);
1380 const auto& group_state = simulator.problem().wellModel().groupState();
1382 std::fill(ws.implicit_ipr_a.begin(), ws.implicit_ipr_a.end(), 0.);
1383 std::fill(ws.implicit_ipr_b.begin(), ws.implicit_ipr_b.end(), 0.);
1385 auto inj_controls = Well::InjectionControls(0);
1386 auto prod_controls = Well::ProductionControls(0);
1387 prod_controls.addControl(Well::ProducerCMode::BHP);
1388 prod_controls.bhp_limit = well_state.well(this->index_of_well_).bhp;
1391 const auto cmode = ws.production_cmode;
1392 ws.production_cmode = Well::ProducerCMode::BHP;
1393 const double dt = simulator.timeStepSize();
1394 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
1396 BVectorWell rhs(this->numberOfSegments());
1398 rhs[0][SPres] = -1.0;
1400 const BVectorWell x_well = this->linSys_.solve(rhs);
1401 constexpr int num_eq = MSWEval::numWellEq;
1402 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx){
1403 const EvalWell comp_rate = this->primary_variables_.getQs(comp_idx);
1404 const int idx = this->modelCompIdxToFlowCompIdx(comp_idx);
1405 for (
size_t pvIdx = 0; pvIdx < num_eq; ++pvIdx) {
1407 ws.implicit_ipr_b[idx] -= x_well[0][pvIdx]*comp_rate.derivative(pvIdx+Indices::numEq);
1409 ws.implicit_ipr_a[idx] = ws.implicit_ipr_b[idx]*ws.bhp - comp_rate.value();
1412 ws.production_cmode = cmode;
1415 template<
typename TypeTag>
1417 MultisegmentWell<TypeTag>::
1418 checkOperabilityUnderTHPLimit(
const Simulator& simulator,
1419 const WellState<Scalar>& well_state,
1420 DeferredLogger& deferred_logger)
1422 const auto& summaryState = simulator.vanguard().summaryState();
1423 const auto obtain_bhp = this->isProducer()
1424 ? computeBhpAtThpLimitProd(
1425 well_state, simulator, summaryState, deferred_logger)
1426 : computeBhpAtThpLimitInj(simulator, summaryState, deferred_logger);
1429 this->operability_status_.can_obtain_bhp_with_thp_limit =
true;
1431 const Scalar bhp_limit = WellBhpThpCalculator(*this).mostStrictBhpFromBhpLimits(summaryState);
1432 this->operability_status_.obey_bhp_limit_with_thp_limit = (*obtain_bhp >= bhp_limit);
1434 const Scalar thp_limit = this->getTHPConstraint(summaryState);
1435 if (this->isProducer() && *obtain_bhp < thp_limit) {
1436 const std::string msg =
" obtained bhp " + std::to_string(unit::convert::to(*obtain_bhp, unit::barsa))
1437 +
" bars is SMALLER than thp limit "
1438 + std::to_string(unit::convert::to(thp_limit, unit::barsa))
1439 +
" bars as a producer for well " + this->name();
1440 deferred_logger.debug(msg);
1442 else if (this->isInjector() && *obtain_bhp > thp_limit) {
1443 const std::string msg =
" obtained bhp " + std::to_string(unit::convert::to(*obtain_bhp, unit::barsa))
1444 +
" bars is LARGER than thp limit "
1445 + std::to_string(unit::convert::to(thp_limit, unit::barsa))
1446 +
" bars as a injector for well " + this->name();
1447 deferred_logger.debug(msg);
1452 this->operability_status_.can_obtain_bhp_with_thp_limit =
false;
1453 this->operability_status_.obey_bhp_limit_with_thp_limit =
false;
1454 if (!this->wellIsStopped()) {
1455 const Scalar thp_limit = this->getTHPConstraint(summaryState);
1456 deferred_logger.debug(
" could not find bhp value at thp limit "
1457 + std::to_string(unit::convert::to(thp_limit, unit::barsa))
1458 +
" bar for well " + this->name() +
", the well might need to be closed ");
1467 template<
typename TypeTag>
1469 MultisegmentWell<TypeTag>::
1470 iterateWellEqWithControl(
const Simulator& simulator,
1472 const Well::InjectionControls& inj_controls,
1473 const Well::ProductionControls& prod_controls,
1474 WellState<Scalar>& well_state,
1475 const GroupState<Scalar>& group_state,
1476 DeferredLogger& deferred_logger)
1478 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return true;
1480 const int max_iter_number = this->param_.max_inner_iter_ms_wells_;
1484 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1489 std::vector<std::vector<Scalar> > residual_history;
1490 std::vector<Scalar> measure_history;
1493 Scalar relaxation_factor = 1.;
1494 bool converged =
false;
1495 bool relax_convergence =
false;
1496 this->regularize_ =
false;
1497 for (; it < max_iter_number; ++it, ++debug_cost_counter_) {
1499 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls, well_state, group_state, deferred_logger);
1501 BVectorWell dx_well;
1503 dx_well = this->linSys_.solve();
1505 catch(
const NumericalProblem& exp) {
1509 deferred_logger.problem(
"In MultisegmentWell::iterateWellEqWithControl for well "
1510 + this->name() +
": "+exp.what());
1514 if (it > this->param_.strict_inner_iter_wells_) {
1515 relax_convergence =
true;
1516 this->regularize_ =
true;
1519 const auto report = getWellConvergence(simulator, well_state, Base::B_avg_, deferred_logger, relax_convergence);
1520 if (report.converged()) {
1527 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1531 residual_history.push_back(residuals);
1532 measure_history.push_back(this->getResidualMeasureValue(well_state,
1533 residual_history[it],
1534 this->param_.tolerance_wells_,
1535 this->param_.tolerance_pressure_ms_wells_,
1538 bool min_relaxation_reached = this->update_relaxation_factor(measure_history, relaxation_factor, this->regularize_, deferred_logger);
1539 if (min_relaxation_reached || this->repeatedStagnation(measure_history, this->regularize_, deferred_logger)) {
1541 const auto reportStag = getWellConvergence(simulator, well_state, Base::B_avg_, deferred_logger,
true);
1542 if (reportStag.converged()) {
1544 std::string message = fmt::format(
"Well stagnates/oscillates but {} manages to get converged with relaxed tolerances in {} inner iterations."
1546 deferred_logger.debug(message);
1552 updateWellState(simulator, dx_well, well_state, deferred_logger, relaxation_factor);
1557 std::ostringstream sstr;
1558 sstr <<
" Well " << this->name() <<
" converged in " << it <<
" inner iterations.";
1559 if (relax_convergence)
1560 sstr <<
" (A relaxed tolerance was used after "<< this->param_.strict_inner_iter_wells_ <<
" iterations)";
1561 deferred_logger.debug(sstr.str());
1563 std::ostringstream sstr;
1564 sstr <<
" Well " << this->name() <<
" did not converge in " << it <<
" inner iterations.";
1565#define EXTRA_DEBUG_MSW 0
1567 sstr <<
"***** Outputting the residual history for well " << this->name() <<
" during inner iterations:";
1568 for (
int i = 0; i < it; ++i) {
1569 const auto& residual = residual_history[i];
1570 sstr <<
" residual at " << i <<
"th iteration ";
1571 for (
const auto& res : residual) {
1574 sstr <<
" " << measure_history[i] <<
" \n";
1577#undef EXTRA_DEBUG_MSW
1578 deferred_logger.debug(sstr.str());
1585 template<
typename TypeTag>
1587 MultisegmentWell<TypeTag>::
1588 iterateWellEqWithSwitching(
const Simulator& simulator,
1590 const Well::InjectionControls& inj_controls,
1591 const Well::ProductionControls& prod_controls,
1592 WellState<Scalar>& well_state,
1593 const GroupState<Scalar>& group_state,
1594 DeferredLogger& deferred_logger,
1595 const bool fixed_control ,
1596 const bool fixed_status )
1598 const int max_iter_number = this->param_.max_inner_iter_ms_wells_;
1602 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1607 std::vector<std::vector<Scalar> > residual_history;
1608 std::vector<Scalar> measure_history;
1611 Scalar relaxation_factor = 1.;
1612 bool converged =
false;
1613 bool relax_convergence =
false;
1614 this->regularize_ =
false;
1615 const auto& summary_state = simulator.vanguard().summaryState();
1620 const int min_its_after_switch = 3;
1621 int its_since_last_switch = min_its_after_switch;
1622 int switch_count= 0;
1624 const auto well_status_orig = this->wellStatus_;
1625 const auto operability_orig = this->operability_status_;
1627 const bool allow_open = this->well_ecl_.getStatus() == WellStatus::OPEN &&
1628 well_state.well(this->index_of_well_).status == WellStatus::OPEN;
1630 const bool allow_switching = !this->wellUnderZeroRateTarget(simulator, well_state, deferred_logger) &&
1631 (!fixed_control || !fixed_status) && allow_open;
1632 bool final_check =
false;
1634 this->operability_status_.resetOperability();
1635 this->operability_status_.solvable =
true;
1637 for (; it < max_iter_number; ++it, ++debug_cost_counter_) {
1638 ++its_since_last_switch;
1639 if (allow_switching && its_since_last_switch >= min_its_after_switch){
1640 const Scalar wqTotal = this->primary_variables_.getWQTotal().value();
1641 bool changed = this->updateWellControlAndStatusLocalIteration(simulator, well_state, group_state,
1642 inj_controls, prod_controls, wqTotal,
1643 deferred_logger, fixed_control,
1646 its_since_last_switch = 0;
1649 if (!changed && final_check) {
1652 final_check =
false;
1656 assembleWellEqWithoutIteration(simulator, dt, inj_controls, prod_controls,
1657 well_state, group_state, deferred_logger);
1659 const BVectorWell dx_well = this->linSys_.solve();
1661 if (it > this->param_.strict_inner_iter_wells_) {
1662 relax_convergence =
true;
1663 this->regularize_ =
true;
1666 const auto report = getWellConvergence(simulator, well_state, Base::B_avg_, deferred_logger, relax_convergence);
1667 converged = report.converged();
1668 if (this->parallel_well_info_.communication().size() > 1 && converged != this->parallel_well_info_.communication().min(converged)) {
1669 OPM_THROW(std::runtime_error, fmt::format(
"Misalignment of the parallel simulation run in iterateWellEqWithSwitching - the well calculation succeeded on rank {} but failed on other ranks.", this->parallel_well_info_.communication().rank()));
1674 if (switch_count > 0 && its_since_last_switch < min_its_after_switch) {
1676 its_since_last_switch = min_its_after_switch;
1684 const auto& [isFinite, residuals] = this->getFiniteWellResiduals(Base::B_avg_, deferred_logger);
1688 residual_history.push_back(residuals);
1692 measure_history.push_back(this->getResidualMeasureValue(well_state,
1693 residual_history[it],
1694 this->param_.tolerance_wells_,
1695 this->param_.tolerance_pressure_ms_wells_,
1697 bool min_relaxation_reached = this->update_relaxation_factor(measure_history, relaxation_factor, this->regularize_, deferred_logger);
1698 if (min_relaxation_reached || this->repeatedStagnation(measure_history, this->regularize_, deferred_logger)) {
1703 updateWellState(simulator, dx_well, well_state, deferred_logger, relaxation_factor);
1707 if (allow_switching){
1709 const bool is_stopped = this->wellIsStopped();
1710 if (this->wellHasTHPConstraints(summary_state)){
1711 this->operability_status_.can_obtain_bhp_with_thp_limit = !is_stopped;
1712 this->operability_status_.obey_thp_limit_under_bhp_limit = !is_stopped;
1714 this->operability_status_.operable_under_only_bhp_limit = !is_stopped;
1717 std::string message = fmt::format(
" Well {} converged in {} inner iterations ("
1718 "{} control/status switches).", this->name(), it, switch_count);
1719 if (relax_convergence) {
1720 message.append(fmt::format(
" (A relaxed tolerance was used after {} iterations)",
1721 this->param_.strict_inner_iter_wells_));
1723 deferred_logger.debug(message);
1725 this->wellStatus_ = well_status_orig;
1726 this->operability_status_ = operability_orig;
1727 const std::string message = fmt::format(
" Well {} did not converge in {} inner iterations ("
1728 "{} control/status switches).", this->name(), it, switch_count);
1729 deferred_logger.debug(message);
1730 this->primary_variables_.outputLowLimitPressureSegments(deferred_logger);
1737 template<
typename TypeTag>
1739 MultisegmentWell<TypeTag>::
1740 assembleWellEqWithoutIteration(
const Simulator& simulator,
1742 const Well::InjectionControls& inj_controls,
1743 const Well::ProductionControls& prod_controls,
1744 WellState<Scalar>& well_state,
1745 const GroupState<Scalar>& group_state,
1746 DeferredLogger& deferred_logger)
1748 if (!this->isOperableAndSolvable() && !this->wellIsStopped())
return;
1751 this->segments_.updateUpwindingSegments(this->primary_variables_);
1754 computeSegmentFluidProperties(simulator, deferred_logger);
1757 this->linSys_.clear();
1759 auto& ws = well_state.well(this->index_of_well_);
1760 ws.phase_mixing_rates.fill(0.0);
1767 const bool allow_cf = this->getAllowCrossFlow() || openCrossFlowAvoidSingularity(simulator);
1769 const int nseg = this->numberOfSegments();
1771 for (
int seg = 0; seg < nseg; ++seg) {
1773 const EvalWell seg_pressure = this->primary_variables_.getSegmentPressure(seg);
1774 auto& perf_data = ws.perf_data;
1775 auto& perf_rates = perf_data.phase_rates;
1776 auto& perf_press_state = perf_data.pressure;
1777 for (
const int perf : this->segments_.perforations()[seg]) {
1778 const int local_perf_index = this->pw_info_.globalToLocal(perf);
1779 if (local_perf_index < 0)
1781 const int cell_idx = this->well_cells_[local_perf_index];
1782 const auto& int_quants = simulator.model().intensiveQuantities(cell_idx, 0);
1783 std::vector<EvalWell> mob(this->num_components_, 0.0);
1784 getMobility(simulator, local_perf_index, mob, deferred_logger);
1785 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(int_quants, cell_idx);
1786 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
1787 const std::vector<Scalar> Tw = this->wellIndex(local_perf_index, int_quants, trans_mult, wellstate_nupcol);
1788 std::vector<EvalWell> cq_s(this->num_components_, 0.0);
1789 EvalWell perf_press;
1790 PerforationRates<Scalar> perfRates;
1791 computePerfRate(int_quants, mob, Tw, seg, perf, seg_pressure,
1792 allow_cf, cq_s, perf_press, perfRates, deferred_logger);
1795 if (this->isProducer()) {
1796 ws.phase_mixing_rates[ws.dissolved_gas] += perfRates.dis_gas;
1797 ws.phase_mixing_rates[ws.vaporized_oil] += perfRates.vap_oil;
1798 perf_data.phase_mixing_rates[local_perf_index][ws.dissolved_gas] = perfRates.dis_gas;
1799 perf_data.phase_mixing_rates[local_perf_index][ws.vaporized_oil] = perfRates.vap_oil;
1803 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
1804 perf_rates[local_perf_index*this->number_of_phases_ + this->modelCompIdxToFlowCompIdx(comp_idx)] = cq_s[comp_idx].value();
1806 perf_press_state[local_perf_index] = perf_press.value();
1808 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
1810 const EvalWell cq_s_effective = cq_s[comp_idx] * this->well_efficiency_factor_;
1812 this->connectionRates_[local_perf_index][comp_idx] = Base::restrictEval(cq_s_effective);
1814 MultisegmentWellAssemble(*this).
1815 assemblePerforationEq(seg, local_perf_index, comp_idx, cq_s_effective, this->linSys_);
1821 const auto& comm = this->parallel_well_info_.communication();
1822 comm.sum(ws.phase_mixing_rates.data(), ws.phase_mixing_rates.size());
1825 if (this->parallel_well_info_.communication().size() > 1) {
1827 this->linSys_.sumDistributed(this->parallel_well_info_.communication());
1829 for (
int seg = 0; seg < nseg; ++seg) {
1833 const EvalWell segment_surface_volume = getSegmentSurfaceVolume(simulator, seg);
1838 const Scalar regularization_factor = this->regularize_? this->param_.regularization_factor_wells_ : 1.0;
1840 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
1841 const EvalWell accumulation_term = regularization_factor * (segment_surface_volume * this->primary_variables_.surfaceVolumeFraction(seg, comp_idx)
1842 - segment_fluid_initial_[seg][comp_idx]) / dt;
1843 MultisegmentWellAssemble(*this).
1844 assembleAccumulationTerm(seg, comp_idx, accumulation_term, this->linSys_);
1849 const int seg_upwind = this->segments_.upwinding_segment(seg);
1850 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
1851 const EvalWell segment_rate =
1852 this->primary_variables_.getSegmentRateUpwinding(seg,
1855 this->well_efficiency_factor_;
1856 MultisegmentWellAssemble(*this).
1857 assembleOutflowTerm(seg, seg_upwind, comp_idx, segment_rate, this->linSys_);
1863 for (
const int inlet : this->segments_.inlets()[seg]) {
1864 const int inlet_upwind = this->segments_.upwinding_segment(inlet);
1865 for (
int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
1866 const EvalWell inlet_rate =
1867 this->primary_variables_.getSegmentRateUpwinding(inlet,
1870 this->well_efficiency_factor_;
1871 MultisegmentWellAssemble(*this).
1872 assembleInflowTerm(seg, inlet, inlet_upwind, comp_idx, inlet_rate, this->linSys_);
1879 const auto& summaryState = simulator.vanguard().summaryState();
1880 const Schedule& schedule = simulator.vanguard().schedule();
1881 const bool stopped_or_zero_target = this->stoppedOrZeroRateTarget(simulator, well_state, deferred_logger);
1882 MultisegmentWellAssemble(*this).
1883 assembleControlEq(well_state,
1890 this->primary_variables_,
1892 stopped_or_zero_target,
1895 const UnitSystem& unit_system = simulator.vanguard().eclState().getDeckUnitSystem();
1896 const auto& summary_state = simulator.vanguard().summaryState();
1897 this->assemblePressureEq(seg, unit_system, well_state, summary_state, this->param_.use_average_density_ms_wells_, deferred_logger);
1901 this->parallel_well_info_.communication().sum(this->ipr_a_.data(), this->ipr_a_.size());
1902 this->linSys_.createSolver();
1908 template<
typename TypeTag>
1910 MultisegmentWell<TypeTag>::
1911 openCrossFlowAvoidSingularity(
const Simulator& simulator)
const
1913 return !this->getAllowCrossFlow() && allDrawDownWrongDirection(simulator);
1917 template<
typename TypeTag>
1919 MultisegmentWell<TypeTag>::
1920 allDrawDownWrongDirection(
const Simulator& simulator)
const
1922 bool all_drawdown_wrong_direction =
true;
1923 const int nseg = this->numberOfSegments();
1925 for (
int seg = 0; seg < nseg; ++seg) {
1926 const EvalWell segment_pressure = this->primary_variables_.getSegmentPressure(seg);
1927 for (
const int perf : this->segments_.perforations()[seg]) {
1928 const int local_perf_index = this->pw_info_.globalToLocal(perf);
1929 if (local_perf_index < 0)
1932 const int cell_idx = this->well_cells_[local_perf_index];
1933 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
1934 const auto& fs = intQuants.fluidState();
1937 const EvalWell perf_seg_press_diff = this->segments_.getPressureDiffSegLocalPerf(seg, local_perf_index);
1939 const Scalar cell_perf_press_diff = this->cell_perforation_pressure_diffs_[local_perf_index];
1941 const Scalar pressure_cell = this->getPerfCellPressure(fs).value();
1942 const Scalar perf_press = pressure_cell - cell_perf_press_diff;
1945 const EvalWell drawdown = perf_press - (segment_pressure + perf_seg_press_diff);
1950 if ( (drawdown < 0. && this->isInjector()) ||
1951 (drawdown > 0. && this->isProducer()) ) {
1952 all_drawdown_wrong_direction =
false;
1957 const auto& comm = this->parallel_well_info_.communication();
1958 if (comm.size() > 1)
1960 all_drawdown_wrong_direction =
1961 (comm.min(all_drawdown_wrong_direction ? 1 : 0) == 1);
1964 return all_drawdown_wrong_direction;
1970 template<
typename TypeTag>
1972 MultisegmentWell<TypeTag>::
1973 updateWaterThroughput(
const double , WellState<Scalar>& )
const
1981 template<
typename TypeTag>
1982 typename MultisegmentWell<TypeTag>::EvalWell
1983 MultisegmentWell<TypeTag>::
1984 getSegmentSurfaceVolume(
const Simulator& simulator,
const int seg_idx)
const
1986 EvalWell temperature;
1987 EvalWell saltConcentration;
1988 int pvt_region_index;
1992 const int cell_idx = this->well_cells_[0];
1993 const auto& intQuants = simulator.model().intensiveQuantities(cell_idx, 0);
1994 const auto& fs = intQuants.fluidState();
1998 auto fsTemperature = fs.temperature(FluidSystem::oilPhaseIdx).value();
1999 fsTemperature = this->pw_info_.broadcastFirstPerforationValue(fsTemperature);
2000 temperature.setValue(fsTemperature);
2002 auto fsSaltConcentration = fs.saltConcentration();
2003 fsSaltConcentration = this->pw_info_.broadcastFirstPerforationValue(fsSaltConcentration);
2004 saltConcentration = this->extendEval(fsSaltConcentration);
2006 pvt_region_index = fs.pvtRegionIndex();
2007 pvt_region_index = this->pw_info_.broadcastFirstPerforationValue(pvt_region_index);
2010 return this->segments_.getSurfaceVolume(temperature,
2012 this->primary_variables_,
2018 template<
typename TypeTag>
2019 std::optional<typename MultisegmentWell<TypeTag>::Scalar>
2020 MultisegmentWell<TypeTag>::
2021 computeBhpAtThpLimitProd(
const WellState<Scalar>& well_state,
2022 const Simulator& simulator,
2023 const SummaryState& summary_state,
2024 DeferredLogger& deferred_logger)
const
2026 return this->MultisegmentWell<TypeTag>::computeBhpAtThpLimitProdWithAlq(
2029 this->getALQ(well_state),
2036 template<
typename TypeTag>
2037 std::optional<typename MultisegmentWell<TypeTag>::Scalar>
2038 MultisegmentWell<TypeTag>::
2039 computeBhpAtThpLimitProdWithAlq(
const Simulator& simulator,
2040 const SummaryState& summary_state,
2041 const Scalar alq_value,
2042 DeferredLogger& deferred_logger,
2043 bool iterate_if_no_solution)
const
2047 auto frates = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2053 std::vector<Scalar> rates(3);
2054 computeWellRatesWithBhp(simulator, bhp, rates, deferred_logger);
2058 auto bhpAtLimit = WellBhpThpCalculator(*this).
2059 computeBhpAtThpLimitProd(frates,
2061 this->maxPerfPress(simulator),
2062 this->getRefDensity(),
2064 this->getTHPConstraint(summary_state),
2070 if (!iterate_if_no_solution)
2071 return std::nullopt;
2073 auto fratesIter = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2077 std::vector<Scalar> rates(3);
2078 computeWellRatesWithBhpIterations(simulator, bhp, rates, deferred_logger);
2082 return WellBhpThpCalculator(*this).
2083 computeBhpAtThpLimitProd(fratesIter,
2085 this->maxPerfPress(simulator),
2086 this->getRefDensity(),
2088 this->getTHPConstraint(summary_state),
2092 template<
typename TypeTag>
2093 std::optional<typename MultisegmentWell<TypeTag>::Scalar>
2094 MultisegmentWell<TypeTag>::
2095 computeBhpAtThpLimitInj(
const Simulator& simulator,
2096 const SummaryState& summary_state,
2097 DeferredLogger& deferred_logger)
const
2100 auto frates = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2106 std::vector<Scalar> rates(3);
2107 computeWellRatesWithBhp(simulator, bhp, rates, deferred_logger);
2111 auto bhpAtLimit = WellBhpThpCalculator(*this).
2112 computeBhpAtThpLimitInj(frates,
2114 this->getRefDensity(),
2123 auto fratesIter = [
this, &simulator, &deferred_logger](
const Scalar bhp) {
2127 std::vector<Scalar> rates(3);
2128 computeWellRatesWithBhpIterations(simulator, bhp, rates, deferred_logger);
2132 return WellBhpThpCalculator(*this).
2133 computeBhpAtThpLimitInj(fratesIter,
2135 this->getRefDensity(),
2146 template<
typename TypeTag>
2147 typename MultisegmentWell<TypeTag>::Scalar
2148 MultisegmentWell<TypeTag>::
2149 maxPerfPress(
const Simulator& simulator)
const
2151 Scalar max_pressure = 0.0;
2152 const int nseg = this->numberOfSegments();
2153 for (
int seg = 0; seg < nseg; ++seg) {
2154 for (
const int perf : this->segments_.perforations()[seg]) {
2155 const int local_perf_index = this->pw_info_.globalToLocal(perf);
2156 if (local_perf_index < 0)
2159 const int cell_idx = this->well_cells_[local_perf_index];
2160 const auto& int_quants = simulator.model().intensiveQuantities(cell_idx, 0);
2161 const auto& fs = int_quants.fluidState();
2162 Scalar pressure_cell = this->getPerfCellPressure(fs).value();
2163 max_pressure = std::max(max_pressure, pressure_cell);
2166 return max_pressure;
2173 template<
typename TypeTag>
2174 std::vector<typename MultisegmentWell<TypeTag>::Scalar>
2180 std::vector<Scalar> well_q_s(this->num_components_, 0.0);
2181 const bool allow_cf = this->getAllowCrossFlow() || openCrossFlowAvoidSingularity(simulator);
2182 const int nseg = this->numberOfSegments();
2183 for (
int seg = 0; seg < nseg; ++seg) {
2185 const Scalar seg_pressure = getValue(this->primary_variables_.getSegmentPressure(seg));
2186 for (
const int perf : this->segments_.perforations()[seg]) {
2187 const int local_perf_index = this->pw_info_.globalToLocal(perf);
2188 if (local_perf_index < 0)
2191 const int cell_idx = this->well_cells_[local_perf_index];
2192 const auto& int_quants = simulator.model().intensiveQuantities(cell_idx, 0);
2193 std::vector<Scalar> mob(this->num_components_, 0.0);
2194 getMobility(simulator, local_perf_index, mob, deferred_logger);
2195 const Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(int_quants, cell_idx);
2196 const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
2197 const std::vector<Scalar> Tw = this->wellIndex(local_perf_index, int_quants, trans_mult, wellstate_nupcol);
2198 std::vector<Scalar> cq_s(this->num_components_, 0.0);
2199 Scalar perf_press = 0.0;
2201 computePerfRate(int_quants, mob, Tw, seg, perf, seg_pressure,
2202 allow_cf, cq_s, perf_press, perf_rates, deferred_logger);
2203 for (
int comp = 0; comp < this->num_components_; ++comp) {
2204 well_q_s[comp] += cq_s[comp];
2208 const auto& comm = this->parallel_well_info_.communication();
2209 if (comm.size() > 1)
2211 comm.sum(well_q_s.data(), well_q_s.size());
2217 template <
typename TypeTag>
2218 std::vector<typename MultisegmentWell<TypeTag>::Scalar>
2222 const int num_seg = this->numberOfSegments();
2223 constexpr int num_eq = MSWEval::numWellEq;
2224 std::vector<Scalar> retval(num_seg * num_eq);
2225 for (
int ii = 0; ii < num_seg; ++ii) {
2226 const auto& pv = this->primary_variables_.value(ii);
2227 std::copy(pv.begin(), pv.end(), retval.begin() + ii * num_eq);
2235 template <
typename TypeTag>
2237 MultisegmentWell<TypeTag>::
2238 setPrimaryVars(
typename std::vector<Scalar>::const_iterator it)
2240 const int num_seg = this->numberOfSegments();
2241 constexpr int num_eq = MSWEval::numWellEq;
2242 std::array<Scalar, num_eq> tmp;
2243 for (
int ii = 0; ii < num_seg; ++ii) {
2244 const auto start = it + ii * num_eq;
2245 std::copy(start, start + num_eq, tmp.begin());
2246 this->primary_variables_.setValue(ii, tmp);
2248 return num_seg * num_eq;
Represents the convergence status of the whole simulator, to make it possible to query and store the ...
Definition ConvergenceReport.hpp:38
Definition DeferredLogger.hpp:57
Definition GroupState.hpp:43
Definition MultisegmentWell.hpp:38
void apply(const BVector &x, BVector &Ax) const override
Ax = Ax - C D^-1 B x.
Definition MultisegmentWell_impl.hpp:216
Manages the initializing and running of time dependent problems.
Definition simulator.hh:92
Vanguard & vanguard()
Return a reference to the grid manager of simulation.
Definition simulator.hh:273
Definition WellInterfaceGeneric.hpp:53
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:66
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
PhaseUsage phaseUsage(const Phases &phases)
Determine the active phases.
Definition phaseUsageFromDeck.cpp:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Definition PerforationData.hpp:41