From a7664f9ea6345986cf044b343ce66ac329761e59 Mon Sep 17 00:00:00 2001 From: Tor Harald Sandve Date: Fri, 27 Nov 2020 07:57:55 +0100 Subject: [PATCH] Add check for operability of MSW --- .../wells/BlackoilWellModel_impl.hpp | 7 +- opm/simulators/wells/MultisegmentWell.hpp | 29 ++- .../wells/MultisegmentWell_impl.hpp | 242 +++++++++++++++--- opm/simulators/wells/StandardWell.hpp | 42 +-- opm/simulators/wells/StandardWell_impl.hpp | 153 +---------- opm/simulators/wells/WellInterface.hpp | 36 ++- opm/simulators/wells/WellInterface_impl.hpp | 141 +++++++++- 7 files changed, 416 insertions(+), 234 deletions(-) diff --git a/opm/simulators/wells/BlackoilWellModel_impl.hpp b/opm/simulators/wells/BlackoilWellModel_impl.hpp index 45e314d55..ec8a26b1e 100644 --- a/opm/simulators/wells/BlackoilWellModel_impl.hpp +++ b/opm/simulators/wells/BlackoilWellModel_impl.hpp @@ -343,8 +343,11 @@ namespace Opm { // do the initialization for all the wells // TODO: to see whether we can postpone of the intialization of the well containers to // optimize the usage of the following several member variables + std::vector< Scalar > B_avg(numComponents(), Scalar() ); + computeAverageFormationFactor(B_avg); + for (auto& well : well_container_) { - well->init(&phase_usage_, depth_, gravity_, local_num_cells_); + well->init(&phase_usage_, depth_, gravity_, local_num_cells_, B_avg); } // update the updated cell flag @@ -435,7 +438,7 @@ namespace Opm { WellInterfacePtr well = createWellForWellTest(well_name, timeStepIdx, deferred_logger); // some preparation before the well can be used - well->init(&phase_usage_, depth_, gravity_, local_num_cells_); + well->init(&phase_usage_, depth_, gravity_, local_num_cells_, B_avg); const Well& wellEcl = schedule().getWell(well_name, timeStepIdx); double well_efficiency_factor = wellEcl.getEfficiencyFactor(); WellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), timeStepIdx), schedule(), timeStepIdx, well_efficiency_factor); diff --git a/opm/simulators/wells/MultisegmentWell.hpp b/opm/simulators/wells/MultisegmentWell.hpp index b2eaaec08..eed0ddd27 100644 --- a/opm/simulators/wells/MultisegmentWell.hpp +++ b/opm/simulators/wells/MultisegmentWell.hpp @@ -113,7 +113,8 @@ namespace Opm virtual void init(const PhaseUsage* phase_usage_arg, const std::vector& depth_arg, const double gravity_arg, - const int num_cells) override; + const int num_cells, + const std::vector< Scalar >& B_avg) override; virtual void initPrimaryVariablesEvaluation() const override; @@ -230,6 +231,9 @@ namespace Opm using Base::perf_depth_; using Base::num_components_; using Base::connectionRates_; + using Base::ipr_a_; + using Base::ipr_b_; + using Base::changed_to_stopped_this_step_; // protected functions from the Base class using Base::phaseUsage; @@ -239,7 +243,8 @@ namespace Opm using Base::getAllowCrossFlow; using Base::scalingFactor; using Base::wellIsStopped_; - + using Base::updateWellOperability; + using Base::checkWellOperability; // TODO: trying to use the information from the Well opm-parser as much // as possible, it will possibly be re-implemented later for efficiency reason. @@ -431,12 +436,6 @@ namespace Opm // handling the overshooting and undershooting of the fractions void processFractions(const int seg) const; - // checking the operability of the well based on current reservoir condition - // it is not implemented for multisegment well yet - virtual void checkWellOperability(const Simulator& ebos_simulator, - const WellState& well_state, - Opm::DeferredLogger& deferred_logger) override; - void updateWellStateFromPrimaryVariables(WellState& well_state, Opm::DeferredLogger& deferred_logger) const; bool frictionalPressureLossConsidered() const; @@ -458,10 +457,6 @@ namespace Opm WellState& well_state, Opm::DeferredLogger& deferred_logger) override; - virtual void wellTestingPhysical(const Simulator& simulator, const std::vector& B_avg, - const double simulation_time, const int report_step, - WellState& well_state, WellTestState& welltest_state, Opm::DeferredLogger& deferred_logger) override; - virtual void updateWaterThroughput(const double dt, WellState& well_state) const override; EvalWell getSegmentSurfaceVolume(const Simulator& ebos_simulator, const int seg_idx) const; @@ -517,6 +512,16 @@ namespace Opm void assembleValvePressureEq(const int seg, WellState& well_state) const; EvalWell pressureDropValve(const int seg) const; + + // check whether the well is operable under BHP limit with current reservoir condition + virtual void checkOperabilityUnderBHPLimitProducer(const WellState& well_state, const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) override; + + // check whether the well is operable under THP limit with current reservoir condition + virtual void checkOperabilityUnderTHPLimitProducer(const Simulator& ebos_simulator, const WellState& well_state, Opm::DeferredLogger& deferred_logger) override; + + // updating the inflow based on the current reservoir condition + virtual void updateIPR(const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) const override; + }; } diff --git a/opm/simulators/wells/MultisegmentWell_impl.hpp b/opm/simulators/wells/MultisegmentWell_impl.hpp index bfc982191..d972740e6 100644 --- a/opm/simulators/wells/MultisegmentWell_impl.hpp +++ b/opm/simulators/wells/MultisegmentWell_impl.hpp @@ -135,9 +135,10 @@ namespace Opm init(const PhaseUsage* phase_usage_arg, const std::vector& depth_arg, const double gravity_arg, - const int num_cells) + const int num_cells, + const std::vector< Scalar >& B_avg) { - Base::init(phase_usage_arg, depth_arg, gravity_arg, num_cells); + Base::init(phase_usage_arg, depth_arg, gravity_arg, num_cells, B_avg); // TODO: for StandardWell, we need to update the perf depth here using depth_arg. // for MultisegmentWell, it is much more complicated. @@ -264,6 +265,8 @@ namespace Opm Opm::DeferredLogger& deferred_logger) { + checkWellOperability(ebosSimulator, well_state, deferred_logger); + const bool use_inner_iterations = param_.use_inner_iterations_ms_wells_; if (use_inner_iterations) { this->iterateWellEquations(ebosSimulator, B_avg, dt, well_state, deferred_logger); @@ -624,6 +627,8 @@ namespace Opm MultisegmentWell:: apply(const BVector& x, BVector& Ax) const { + if (!this->isOperable() && !this->wellIsStopped()) return; + if ( param_.matrix_add_well_contributions_ ) { // Contributions are already in the matrix itself @@ -649,6 +654,8 @@ namespace Opm MultisegmentWell:: apply(BVector& r) const { + if (!this->isOperable() && !this->wellIsStopped()) return; + // invDrw_ = duneD^-1 * resWell_ const BVectorWell invDrw = mswellhelpers::applyUMFPack(duneD_, duneDSolver_, resWell_); // r = r - duneC_^T * invDrw @@ -726,6 +733,8 @@ namespace Opm WellState& well_state, Opm::DeferredLogger& deferred_logger) const { + if (!this->isOperable() && !this->wellIsStopped()) return; + BVectorWell xw(1); recoverSolutionWell(x, xw); updateWellState(xw, well_state, deferred_logger); @@ -820,7 +829,6 @@ namespace Opm - template void MultisegmentWell:: @@ -856,6 +864,7 @@ namespace Opm well_state_copy.currentProductionControls()[index_of_well_] = Well::ProducerCMode::BHP; } well_state_copy.bhp()[well_copy.index_of_well_] = bhp; + well_copy.calculateExplicitQuantities(ebosSimulator, well_state_copy, deferred_logger); const double dt = ebosSimulator.timeStepSize(); // iterate to get a solution at the given bhp. @@ -932,6 +941,7 @@ namespace Opm { // TODO: to test using rate conversion coefficients to see if it will be better than // this default one + if (!this->isOperable() && !this->wellIsStopped()) return; const Well& well = Base::wellEcl(); @@ -1005,6 +1015,8 @@ namespace Opm MultisegmentWell:: recoverSolutionWell(const BVector& x, BVectorWell& xw) const { + if (!this->isOperable() && !this->wellIsStopped()) return; + BVectorWell resWell = resWell_; // resWell = resWell - B * x duneB_.mmv(x, resWell); @@ -1021,6 +1033,8 @@ namespace Opm MultisegmentWell:: solveEqAndUpdateWellState(WellState& well_state, Opm::DeferredLogger& deferred_logger) { + if (!this->isOperable() && !this->wellIsStopped()) return; + // We assemble the well equations, then we check the convergence, // which is why we do not put the assembleWellEq here. const BVectorWell dx_well = mswellhelpers::applyUMFPack(duneD_, duneDSolver_, resWell_); @@ -1113,6 +1127,8 @@ namespace Opm Opm::DeferredLogger& deferred_logger, const double relaxation_factor) const { + if (!this->isOperable() && !this->wellIsStopped()) return; + const double dFLimit = param_.dwell_fraction_max_; const double max_pressure_change = param_.max_pressure_change_ms_wells_; const std::vector > old_primary_variables = primary_variables_; @@ -2374,33 +2390,199 @@ namespace Opm } - - - - template + template void MultisegmentWell:: - checkWellOperability(const Simulator& /* ebos_simulator */, - const WellState& /* well_state */, - Opm::DeferredLogger& deferred_logger) + checkOperabilityUnderBHPLimitProducer(const WellState& well_state, const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) { - const bool checkOperability = EWOMS_GET_PARAM(TypeTag, bool, EnableWellOperabilityCheck); - if (!checkOperability) { - return; - } + const auto& summaryState = ebos_simulator.vanguard().summaryState(); + const double bhp_limit = Base::mostStrictBhpFromBhpLimits(summaryState); + // Crude but works: default is one atmosphere. + // TODO: a better way to detect whether the BHP is defaulted or not + const bool bhp_limit_not_defaulted = bhp_limit > 1.5 * unit::barsa; + if ( bhp_limit_not_defaulted || !this->wellHasTHPConstraints(summaryState) ) { + // if the BHP limit is not defaulted or the well does not have a THP limit + // we need to check the BHP limit - // focusing on PRODUCER for now + for (int p = 0; p < number_of_phases_; ++p) { + const double temp = ipr_a_[p] - ipr_b_[p] * bhp_limit; + if (temp < 0.) { + this->operability_status_.operable_under_only_bhp_limit = false; + break; + } + } + + // checking whether running under BHP limit will violate THP limit + if (this->operability_status_.operable_under_only_bhp_limit && this->wellHasTHPConstraints(summaryState)) { + // option 1: calculate well rates based on the BHP limit. + // option 2: stick with the above IPR curve + // we use IPR here + std::vector well_rates_bhp_limit; + computeWellRatesWithBhp(ebos_simulator, Base::B_avg_, bhp_limit, well_rates_bhp_limit, deferred_logger); + + const double thp = calculateThpFromBhp(well_rates_bhp_limit, bhp_limit, deferred_logger); + + const double thp_limit = this->getTHPConstraint(summaryState); + + if (thp < thp_limit) { + this->operability_status_.obey_thp_limit_under_bhp_limit = false; + } + } + } else { + // defaulted BHP and there is a THP constraint + // default BHP limit is about 1 atm. + // when applied the hydrostatic pressure correction dp, + // most likely we get a negative value (bhp + dp)to search in the VFP table, + // which is not desirable. + // we assume we can operate under defaulted BHP limit and will violate the THP limit + // when operating under defaulted BHP limit. + this->operability_status_.operable_under_only_bhp_limit = true; + this->operability_status_.obey_thp_limit_under_bhp_limit = false; + } + } + + + + template + void + MultisegmentWell:: + updateIPR(const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) const + { + // TODO: not handling solvent related here for now + + // TODO: it only handles the producers for now + // the formular for the injectors are not formulated yet if (this->isInjector()) { return; } - if (!this->underPredictionMode() ) { - return; - } + // initialize all the values to be zero to begin with + std::fill(ipr_a_.begin(), ipr_a_.end(), 0.); + std::fill(ipr_b_.begin(), ipr_b_.end(), 0.); - const std::string msg = "Support of well operability checking for multisegment wells is not implemented " - "yet, checkWellOperability() for " + name() + " will do nothing"; - deferred_logger.warning("NO_OPERATABILITY_CHECKING_MS_WELLS", msg); + const int nseg = numberOfSegments(); + double seg_bhp_press_diff = 0; + double ref_depth = ref_depth_; + for (int seg = 0; seg < nseg; ++seg) { + // calculating the perforation rate for each perforation that belongs to this segment + const double segment_depth = segmentSet()[seg].depth(); + const double dp = wellhelpers::computeHydrostaticCorrection(ref_depth, segment_depth, segment_densities_[seg].value(), gravity_); + ref_depth = segment_depth; + seg_bhp_press_diff += dp; + for (const int perf : segment_perforations_[seg]) { + //std::vector mob(num_components_, {numWellEq_ + numEq, 0.0}); + std::vector mob(num_components_, 0.0); + + // TODO: mabye we should store the mobility somewhere, so that we only need to calculate it one per iteration + getMobility(ebos_simulator, perf, mob); + + const int cell_idx = well_cells_[perf]; + const auto& int_quantities = *(ebos_simulator.model().cachedIntensiveQuantities(cell_idx, /*timeIdx=*/ 0)); + const auto& fs = int_quantities.fluidState(); + // the pressure of the reservoir grid block the well connection is in + // pressure difference between the segment and the perforation + const double perf_seg_press_diff = gravity_ * segment_densities_[seg].value() * perforation_segment_depth_diffs_[perf]; + // pressure difference between the perforation and the grid cell + const double cell_perf_press_diff = cell_perforation_pressure_diffs_[perf]; + const double pressure_cell = fs.pressure(FluidSystem::oilPhaseIdx).value(); + + // calculating the b for the connection + std::vector b_perf(num_components_); + for (size_t phase = 0; phase < FluidSystem::numPhases; ++phase) { + if (!FluidSystem::phaseIsActive(phase)) { + continue; + } + const unsigned comp_idx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phase)); + b_perf[comp_idx] = fs.invB(phase).value(); + } + + // the pressure difference between the connection and BHP + const double h_perf = cell_perf_press_diff + perf_seg_press_diff + seg_bhp_press_diff; + const double pressure_diff = pressure_cell - h_perf; + + // Let us add a check, since the pressure is calculated based on zero value BHP + // it should not be negative anyway. If it is negative, we might need to re-formulate + // to taking into consideration the crossflow here. + if (pressure_diff <= 0.) { + deferred_logger.warning("NON_POSITIVE_DRAWDOWN_IPR", + "non-positive drawdown found when updateIPR for well " + name()); + } + + // the well index associated with the connection + const double tw_perf = well_index_[perf]*ebos_simulator.problem().template rockCompTransMultiplier(int_quantities, cell_idx); + + // TODO: there might be some indices related problems here + // phases vs components + // ipr values for the perforation + std::vector ipr_a_perf(ipr_a_.size()); + std::vector ipr_b_perf(ipr_b_.size()); + for (int p = 0; p < number_of_phases_; ++p) { + const double tw_mob = tw_perf * mob[p].value() * b_perf[p]; + ipr_a_perf[p] += tw_mob * pressure_diff; + ipr_b_perf[p] += tw_mob; + } + + // we need to handle the rs and rv when both oil and gas are present + if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { + const unsigned oil_comp_idx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); + const unsigned gas_comp_idx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + const double rs = (fs.Rs()).value(); + const double rv = (fs.Rv()).value(); + + const double dis_gas_a = rs * ipr_a_perf[oil_comp_idx]; + const double vap_oil_a = rv * ipr_a_perf[gas_comp_idx]; + + ipr_a_perf[gas_comp_idx] += dis_gas_a; + ipr_a_perf[oil_comp_idx] += vap_oil_a; + + const double dis_gas_b = rs * ipr_b_perf[oil_comp_idx]; + const double vap_oil_b = rv * ipr_b_perf[gas_comp_idx]; + + ipr_b_perf[gas_comp_idx] += dis_gas_b; + ipr_b_perf[oil_comp_idx] += vap_oil_b; + } + + for (int p = 0; p < number_of_phases_; ++p) { + // TODO: double check the indices here + ipr_a_[ebosCompIdxToFlowCompIdx(p)] += ipr_a_perf[p]; + ipr_b_[ebosCompIdxToFlowCompIdx(p)] += ipr_b_perf[p]; + } + } + } + } + + template + void + MultisegmentWell:: + checkOperabilityUnderTHPLimitProducer(const Simulator& ebos_simulator, const WellState& well_state, Opm::DeferredLogger& deferred_logger) + { + const auto& summaryState = ebos_simulator.vanguard().summaryState(); + const auto obtain_bhp = computeBhpAtThpLimitProd(ebos_simulator, Base::B_avg_, summaryState, deferred_logger); + + if (obtain_bhp) { + this->operability_status_.can_obtain_bhp_with_thp_limit = true; + + const double bhp_limit = Base::mostStrictBhpFromBhpLimits(summaryState); + this->operability_status_.obey_bhp_limit_with_thp_limit = (*obtain_bhp >= bhp_limit); + + const double thp_limit = this->getTHPConstraint(summaryState); + if (*obtain_bhp < thp_limit) { + const std::string msg = " obtained bhp " + std::to_string(unit::convert::to(*obtain_bhp, unit::barsa)) + + " bars is SMALLER than thp limit " + + std::to_string(unit::convert::to(thp_limit, unit::barsa)) + + " bars as a producer for well " + name(); + deferred_logger.debug(msg); + } + } else { + this->operability_status_.can_obtain_bhp_with_thp_limit = false; + this->operability_status_.obey_bhp_limit_with_thp_limit = false; + if (!this->wellIsStopped()) { + const double thp_limit = this->getTHPConstraint(summaryState); + deferred_logger.debug(" could not find bhp value at thp limit " + + std::to_string(unit::convert::to(thp_limit, unit::barsa)) + + " bar for well " + name() + ", the well might need to be closed "); + } + } } @@ -2503,6 +2685,8 @@ namespace Opm WellState& well_state, Opm::DeferredLogger& deferred_logger) { + if (!this->isOperable() && !this->wellIsStopped()) return true; + const int max_iter_number = param_.max_inner_iter_ms_wells_; const WellState well_state0 = well_state; const std::vector residuals0 = getWellResiduals(B_avg); @@ -2620,6 +2804,8 @@ namespace Opm Opm::DeferredLogger& deferred_logger) { + if (!this->isOperable() && !this->wellIsStopped()) return; + // update the upwinding segments updateUpwindingSegments(); @@ -2847,20 +3033,6 @@ namespace Opm } - template - void - MultisegmentWell:: - wellTestingPhysical(const Simulator& /* simulator */, const std::vector& /* B_avg */, - const double /* simulation_time */, const int /* report_step */, - WellState& /* well_state */, WellTestState& /* welltest_state */, Opm::DeferredLogger& deferred_logger) - { - const std::string msg = "Support of well testing for physical limits for multisegment wells is not " - "implemented yet, wellTestingPhysical() for " + name() + " will do nothing"; - deferred_logger.warning("NO_WELLTESTPHYSICAL_CHECKING_MS_WELLS", msg); - } - - - template diff --git a/opm/simulators/wells/StandardWell.hpp b/opm/simulators/wells/StandardWell.hpp index de14cf081..fa49e1622 100644 --- a/opm/simulators/wells/StandardWell.hpp +++ b/opm/simulators/wells/StandardWell.hpp @@ -174,7 +174,8 @@ namespace Opm virtual void init(const PhaseUsage* phase_usage_arg, const std::vector& depth_arg, const double gravity_arg, - const int num_cells) override; + const int num_cells, + const std::vector< Scalar >& B_avg) override; virtual void initPrimaryVariablesEvaluation() const override; @@ -337,6 +338,8 @@ namespace Opm using Base::wfoam; using Base::scalingFactor; using Base::mostStrictBhpFromBhpLimits; + using Base::updateWellOperability; + using Base::checkWellOperability; // protected member variables from the Base class using Base::current_step_; @@ -359,6 +362,9 @@ namespace Opm using Base::perf_rep_radius_; using Base::perf_length_; using Base::bore_diameters_; + using Base::ipr_a_; + using Base::ipr_b_; + using Base::changed_to_stopped_this_step_; using Base::wellIsStopped_; @@ -397,14 +403,6 @@ namespace Opm // the saturations in the well bore under surface conditions at the beginning of the time step std::vector F0_; - // the vectors used to describe the inflow performance relationship (IPR) - // Q = IPR_A - BHP * IPR_B - // TODO: it minght need to go to WellInterface, let us implement it in StandardWell first - // it is only updated and used for producers for now - mutable std::vector ipr_a_; - mutable std::vector ipr_b_; - - bool changed_to_stopped_this_step_ = false; // Enable GLIFT debug mode. This will enable output of logging messages. bool glift_debug = false; @@ -523,14 +521,6 @@ namespace Opm // handle the non reasonable fractions due to numerical overshoot void processFractions() const; - // updating the inflow based on the current reservoir condition - void updateIPR(const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) const; - - // update the operability status of the well is operable under the current reservoir condition - // mostly related to BHP limit and THP limit - virtual void checkWellOperability(const Simulator& ebos_simulator, - const WellState& well_state, - Opm::DeferredLogger& deferred_logger) override; virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator, const double dt, @@ -553,17 +543,14 @@ namespace Opm EvalWell& cq_s_zfrac_effective, Opm::DeferredLogger& deferred_logger) const; - // check whether the well is operable under the current reservoir condition - // mostly related to BHP limit and THP limit - void updateWellOperability(const Simulator& ebos_simulator, - const WellState& well_state, - Opm::DeferredLogger& deferred_logger); - // check whether the well is operable under BHP limit with current reservoir condition - void checkOperabilityUnderBHPLimitProducer(const WellState& well_state, const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger); + virtual void checkOperabilityUnderBHPLimitProducer(const WellState& well_state, const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) override; // check whether the well is operable under THP limit with current reservoir condition - void checkOperabilityUnderTHPLimitProducer(const Simulator& ebos_simulator, const WellState& well_state, Opm::DeferredLogger& deferred_logger); + virtual void checkOperabilityUnderTHPLimitProducer(const Simulator& ebos_simulator, const WellState& well_state, Opm::DeferredLogger& deferred_logger) override; + + // updating the inflow based on the current reservoir condition + virtual void updateIPR(const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) const override; // for a well, when all drawdown are in the wrong direction, then this well will not // be able to produce/inject . @@ -594,11 +581,6 @@ namespace Opm static double relaxationFactorRate(const std::vector& primary_variables, const BVectorWell& dwells); - virtual void wellTestingPhysical(const Simulator& simulator, const std::vector& B_avg, - const double simulation_time, const int report_step, - WellState& well_state, WellTestState& welltest_state, - Opm::DeferredLogger& deferred_logger) override; - // calculate the skin pressure based on water velocity, throughput and polymer concentration. // throughput is used to describe the formation damage during water/polymer injection. // calculated skin pressure will be applied to the drawdown during perforation rate calculation diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index 0d81e9546..01592a6c0 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -49,8 +49,6 @@ namespace Opm , perf_pressure_diffs_(number_of_perforations_) , parallelB_(duneB_, pw_info) , F0_(numWellConservationEq) - , ipr_a_(number_of_phases_) - , ipr_b_(number_of_phases_) { assert(num_components_ == numWellConservationEq); @@ -69,9 +67,10 @@ namespace Opm init(const PhaseUsage* phase_usage_arg, const std::vector& depth_arg, const double gravity_arg, - const int num_cells) + const int num_cells, + const std::vector< Scalar >& B_avg) { - Base::init(phase_usage_arg, depth_arg, gravity_arg, num_cells); + Base::init(phase_usage_arg, depth_arg, gravity_arg, num_cells, B_avg); perf_depth_.resize(number_of_perforations_, 0.); for (int perf = 0; perf < number_of_perforations_; ++perf) { @@ -1680,88 +1679,6 @@ namespace Opm } - - - - template - void - StandardWell:: - checkWellOperability(const Simulator& ebos_simulator, - const WellState& well_state, - Opm::DeferredLogger& deferred_logger) - { - - const bool checkOperability = EWOMS_GET_PARAM(TypeTag, bool, EnableWellOperabilityCheck); - if (!checkOperability) { - return; - } - - // focusing on PRODUCER for now - if (this->isInjector()) { - return; - } - - if (!this->underPredictionMode() ) { - return; - } - - if (this->wellIsStopped() && !changed_to_stopped_this_step_) { - return; - } - - const bool old_well_operable = this->operability_status_.isOperable(); - - updateWellOperability(ebos_simulator, well_state, deferred_logger); - - const bool well_operable = this->operability_status_.isOperable(); - - if (!well_operable && old_well_operable) { - if (well_ecl_.getAutomaticShutIn()) { - deferred_logger.info(" well " + name() + " gets SHUT during iteration "); - } else { - if (!this->wellIsStopped()) { - deferred_logger.info(" well " + name() + " gets STOPPED during iteration "); - this->stopWell(); - changed_to_stopped_this_step_ = true; - } - } - } else if (well_operable && !old_well_operable) { - deferred_logger.info(" well " + name() + " gets REVIVED during iteration "); - this->openWell(); - changed_to_stopped_this_step_ = false; - } - } - - - - - - template - void - StandardWell:: - updateWellOperability(const Simulator& ebos_simulator, - const WellState& well_state, - Opm::DeferredLogger& deferred_logger) - { - this->operability_status_.reset(); - - updateIPR(ebos_simulator, deferred_logger); - - // checking the BHP limit related - checkOperabilityUnderBHPLimitProducer(well_state, ebos_simulator, deferred_logger); - - const auto& summaryState = ebos_simulator.vanguard().summaryState(); - - // checking whether the well can operate under the THP constraints. - if (this->wellHasTHPConstraints(summaryState)) { - checkOperabilityUnderTHPLimitProducer(ebos_simulator, well_state, deferred_logger); - } - } - - - - - template void StandardWell:: @@ -3346,70 +3263,6 @@ namespace Opm - - template - void - StandardWell:: - wellTestingPhysical(const Simulator& ebos_simulator, const std::vector& B_avg, - const double /* simulation_time */, const int /* report_step */, - WellState& well_state, WellTestState& welltest_state, - Opm::DeferredLogger& deferred_logger) - { - deferred_logger.info(" well " + name() + " is being tested for physical limits"); - - // some most difficult things are the explicit quantities, since there is no information - // in the WellState to do a decent initialization - - // TODO: Let us assume that the simulator is updated - - // Let us try to do a normal simualtion running, to keep checking the operability status - // If the well is not operable during any of the time. It means it does not pass the physical - // limit test. - - // create a copy of the well_state to use. If the operability checking is sucessful, we use this one - // to replace the original one - WellState well_state_copy = well_state; - - // TODO: well state for this well is kind of all zero status - // we should be able to provide a better initialization - calculateExplicitQuantities(ebos_simulator, well_state_copy, deferred_logger); - - updateWellOperability(ebos_simulator, well_state_copy, deferred_logger); - - if ( !this->isOperable() ) { - const std::string msg = " well " + name() + " is not operable during well testing for physical reason"; - deferred_logger.debug(msg); - return; - } - - updateWellStateWithTarget(ebos_simulator, well_state_copy, deferred_logger); - - calculateExplicitQuantities(ebos_simulator, well_state_copy, deferred_logger); - - const double dt = ebos_simulator.timeStepSize(); - const bool converged = this->iterateWellEquations(ebos_simulator, B_avg, dt, well_state_copy, deferred_logger); - - if (!converged) { - const std::string msg = " well " + name() + " did not get converged during well testing for physical reason"; - deferred_logger.debug(msg); - return; - } - - if (this->isOperable() ) { - welltest_state.openWell(name(), WellTestConfig::PHYSICAL ); - const std::string msg = " well " + name() + " is re-opened through well testing for physical reason"; - deferred_logger.info(msg); - well_state = well_state_copy; - } else { - const std::string msg = " well " + name() + " is not operable during well testing for physical reason"; - deferred_logger.debug(msg); - } - } - - - - - template typename StandardWell::EvalWell StandardWell:: diff --git a/opm/simulators/wells/WellInterface.hpp b/opm/simulators/wells/WellInterface.hpp index a3bff7353..7f87907b1 100644 --- a/opm/simulators/wells/WellInterface.hpp +++ b/opm/simulators/wells/WellInterface.hpp @@ -158,7 +158,8 @@ namespace Opm virtual void init(const PhaseUsage* phase_usage_arg, const std::vector& depth_arg, const double gravity_arg, - const int num_cells); + const int num_cells, + const std::vector< Scalar >& B_avg); virtual void initPrimaryVariablesEvaluation() const = 0; @@ -271,7 +272,14 @@ namespace Opm void updatePerforatedCell(std::vector& is_cell_perforated); - virtual void checkWellOperability(const Simulator& ebos_simulator, const WellState& well_state, Opm::DeferredLogger& deferred_logger) = 0; + void checkWellOperability(const Simulator& ebos_simulator, const WellState& well_state, Opm::DeferredLogger& deferred_logger); + + // check whether the well is operable under the current reservoir condition + // mostly related to BHP limit and THP limit + void updateWellOperability(const Simulator& ebos_simulator, + const WellState& well_state, + Opm::DeferredLogger& deferred_logger); + // whether the well is operable bool isOperable() const; @@ -411,6 +419,17 @@ namespace Opm std::optional dynamic_thp_limit_; + std::vector< Scalar > B_avg_; + + // the vectors used to describe the inflow performance relationship (IPR) + // Q = IPR_A - BHP * IPR_B + // TODO: it minght need to go to WellInterface, let us implement it in StandardWell first + // it is only updated and used for producers for now + mutable std::vector ipr_a_; + mutable std::vector ipr_b_; + + bool changed_to_stopped_this_step_ = false; + const PhaseUsage& phaseUsage() const; int flowPhaseToEbosCompIdx( const int phaseIdx ) const; @@ -476,13 +495,22 @@ namespace Opm OperabilityStatus operability_status_; + // check whether the well is operable under BHP limit with current reservoir condition + virtual void checkOperabilityUnderBHPLimitProducer(const WellState& well_state, const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) =0; + + // check whether the well is operable under THP limit with current reservoir condition + virtual void checkOperabilityUnderTHPLimitProducer(const Simulator& ebos_simulator, const WellState& well_state, Opm::DeferredLogger& deferred_logger) =0; + + virtual void updateIPR(const Simulator& ebos_simulator, Opm::DeferredLogger& deferred_logger) const=0; + + void wellTestingEconomic(const Simulator& simulator, const std::vector& B_avg, const double simulation_time, const WellState& well_state, WellTestState& welltest_state, Opm::DeferredLogger& deferred_logger); - virtual void wellTestingPhysical(const Simulator& simulator, const std::vector& B_avg, + void wellTestingPhysical(const Simulator& simulator, const std::vector& B_avg, const double simulation_time, const int report_step, - WellState& well_state, WellTestState& welltest_state, Opm::DeferredLogger& deferred_logger) = 0; + WellState& well_state, WellTestState& welltest_state, Opm::DeferredLogger& deferred_logger); virtual void assembleWellEqWithoutIteration(const Simulator& ebosSimulator, diff --git a/opm/simulators/wells/WellInterface_impl.hpp b/opm/simulators/wells/WellInterface_impl.hpp index 30a7c5944..3ae225a40 100644 --- a/opm/simulators/wells/WellInterface_impl.hpp +++ b/opm/simulators/wells/WellInterface_impl.hpp @@ -51,6 +51,8 @@ namespace Opm , index_of_well_(index_of_well) , first_perf_(first_perf_index) , perf_data_(&perf_data) + , ipr_a_(number_of_phases_) + , ipr_b_(number_of_phases_) { assert(well.name()==pw_info.name()); assert(std::is_sorted(perf_data.begin(), perf_data.end(), @@ -122,10 +124,12 @@ namespace Opm init(const PhaseUsage* phase_usage_arg, const std::vector& /* depth_arg */, const double gravity_arg, - const int /* num_cells */) + const int /* num_cells */, + const std::vector< Scalar >& B_avg) { phase_usage_ = phase_usage_arg; gravity_ = gravity_arg; + B_avg_ = B_avg; } @@ -1409,6 +1413,141 @@ namespace Opm } + template + void + WellInterface:: + wellTestingPhysical(const Simulator& ebos_simulator, const std::vector& B_avg, + const double /* simulation_time */, const int /* report_step */, + WellState& well_state, WellTestState& welltest_state, + Opm::DeferredLogger& deferred_logger) + { + deferred_logger.info(" well " + name() + " is being tested for physical limits"); + + // some most difficult things are the explicit quantities, since there is no information + // in the WellState to do a decent initialization + + // TODO: Let us assume that the simulator is updated + + // Let us try to do a normal simualtion running, to keep checking the operability status + // If the well is not operable during any of the time. It means it does not pass the physical + // limit test. + + // create a copy of the well_state to use. If the operability checking is sucessful, we use this one + // to replace the original one + WellState well_state_copy = well_state; + + // TODO: well state for this well is kind of all zero status + // we should be able to provide a better initialization + calculateExplicitQuantities(ebos_simulator, well_state_copy, deferred_logger); + + updateWellOperability(ebos_simulator, well_state_copy, deferred_logger); + + if ( !this->isOperable() ) { + const std::string msg = " well " + name() + " is not operable during well testing for physical reason"; + deferred_logger.debug(msg); + return; + } + + updateWellStateWithTarget(ebos_simulator, well_state_copy, deferred_logger); + + calculateExplicitQuantities(ebos_simulator, well_state_copy, deferred_logger); + + const double dt = ebos_simulator.timeStepSize(); + const bool converged = this->iterateWellEquations(ebos_simulator, B_avg, dt, well_state_copy, deferred_logger); + + if (!converged) { + const std::string msg = " well " + name() + " did not get converged during well testing for physical reason"; + deferred_logger.debug(msg); + return; + } + + if (this->isOperable() ) { + welltest_state.openWell(name(), WellTestConfig::PHYSICAL ); + const std::string msg = " well " + name() + " is re-opened through well testing for physical reason"; + deferred_logger.info(msg); + well_state = well_state_copy; + } else { + const std::string msg = " well " + name() + " is not operable during well testing for physical reason"; + deferred_logger.debug(msg); + } + } + + + + template + void + WellInterface:: + checkWellOperability(const Simulator& ebos_simulator, + const WellState& well_state, + Opm::DeferredLogger& deferred_logger) + { + + const bool checkOperability = EWOMS_GET_PARAM(TypeTag, bool, EnableWellOperabilityCheck); + if (!checkOperability) { + return; + } + + // focusing on PRODUCER for now + if (this->isInjector()) { + return; + } + + if (!this->underPredictionMode() ) { + return; + } + + if (this->wellIsStopped() && !changed_to_stopped_this_step_) { + return; + } + + const bool old_well_operable = this->operability_status_.isOperable(); + + updateWellOperability(ebos_simulator, well_state, deferred_logger); + + const bool well_operable = this->operability_status_.isOperable(); + + if (!well_operable && old_well_operable) { + if (well_ecl_.getAutomaticShutIn()) { + deferred_logger.info(" well " + name() + " gets SHUT during iteration "); + } else { + if (!this->wellIsStopped()) { + deferred_logger.info(" well " + name() + " gets STOPPED during iteration "); + this->stopWell(); + changed_to_stopped_this_step_ = true; + } + } + } else if (well_operable && !old_well_operable) { + deferred_logger.info(" well " + name() + " gets REVIVED during iteration "); + this->openWell(); + changed_to_stopped_this_step_ = false; + } + } + + + + + + template + void + WellInterface:: + updateWellOperability(const Simulator& ebos_simulator, + const WellState& well_state, + Opm::DeferredLogger& deferred_logger) + { + this->operability_status_.reset(); + + updateIPR(ebos_simulator, deferred_logger); + + // checking the BHP limit related + checkOperabilityUnderBHPLimitProducer(well_state, ebos_simulator, deferred_logger); + + const auto& summaryState = ebos_simulator.vanguard().summaryState(); + + // checking whether the well can operate under the THP constraints. + if (this->wellHasTHPConstraints(summaryState)) { + checkOperabilityUnderTHPLimitProducer(ebos_simulator, well_state, deferred_logger); + } + }