diff --git a/opm/simulators/wells/VFPProdProperties.cpp b/opm/simulators/wells/VFPProdProperties.cpp index 971a0e366a..dd99cceae8 100644 --- a/opm/simulators/wells/VFPProdProperties.cpp +++ b/opm/simulators/wells/VFPProdProperties.cpp @@ -134,116 +134,6 @@ bhpwithflo(const std::vector& flos, } - - - -double -VFPProdProperties:: -calculateBhpWithTHPTarget(const std::vector& ipr_a, - const std::vector& ipr_b, - const double bhp_limit, - const double thp_table_id, - const double thp_limit, - const double alq, - const double dp) const -{ - // For producers, bhp_safe_limit is the highest BHP value that can still produce based on IPR - double bhp_safe_limit = 1.e100; - for (size_t i = 0; i < ipr_a.size(); ++i) { - if (ipr_b[i] == 0.) continue; - - const double bhp = ipr_a[i] / ipr_b[i]; - if (bhp < bhp_safe_limit) { - bhp_safe_limit = bhp; - } - } - - // Here, we use the middle point between the bhp_limit and bhp_safe_limit to calculate the ratio of the flow - // and the middle point serves one of the two points to describe inflow performance relationship line - const double bhp_middle = (bhp_limit + bhp_safe_limit) / 2.0; - - // FLO is the rate based on the type specified with the VFP table - // The two points correspond to the bhp values of bhp_limit, and the middle of bhp_limit and bhp_safe_limit - // for producers, the rates are negative - std::vector rates_bhp_limit(ipr_a.size()); - std::vector rates_bhp_middle(ipr_a.size()); - for (size_t i = 0; i < rates_bhp_limit.size(); ++i) { - rates_bhp_limit[i] = bhp_limit * ipr_b[i] - ipr_a[i]; - rates_bhp_middle[i] = bhp_middle * ipr_b[i] - ipr_a[i]; - } - - // TODO: we need to be careful that there is nothings wrong related to the indices here - const int Water = BlackoilPhases::Aqua; - const int Oil = BlackoilPhases::Liquid; - const int Gas = BlackoilPhases::Vapour; - - const VFPProdTable& table = detail::getTable(m_tables, thp_table_id); - const double aqua_bhp_limit = rates_bhp_limit[Water]; - const double liquid_bhp_limit = rates_bhp_limit[Oil]; - const double vapour_bhp_limit = rates_bhp_limit[Gas]; - const double flo_bhp_limit = detail::getFlo(table, aqua_bhp_limit, liquid_bhp_limit, vapour_bhp_limit ); - - const double aqua_bhp_middle = rates_bhp_middle[Water]; - const double liquid_bhp_middle = rates_bhp_middle[Oil]; - const double vapour_bhp_middle = rates_bhp_middle[Gas]; - const double flo_bhp_middle = detail::getFlo(table, aqua_bhp_middle, liquid_bhp_middle, vapour_bhp_middle ); - - // we use the ratios based on the middle value of bhp_limit and bhp_safe_limit - const double wfr = detail::getWFR(table, aqua_bhp_middle, liquid_bhp_middle, vapour_bhp_middle); - const double gfr = detail::getGFR(table, aqua_bhp_middle, liquid_bhp_middle, vapour_bhp_middle); - - // we get the flo sampling points from the table, - // then extend it with zero and rate under bhp_limit for extrapolation - std::vector flo_samples = table.getFloAxis(); - - if (flo_samples[0] > 0.) { - flo_samples.insert(flo_samples.begin(), 0.); - } - - if (flo_samples.back() < std::abs(flo_bhp_limit)) { - flo_samples.push_back(std::abs(flo_bhp_limit)); - } - - // kind of unncessarily following the tradation that producers should have negative rates - // the key is here that it should be consistent with the function bhpwithflo - for (double& value : flo_samples) { - value = -value; - } - - // get the bhp sampling values based on the flo sample values - const std::vector bhp_flo_samples = bhpwithflo(flo_samples, thp_table_id, wfr, gfr, thp_limit, alq, dp); - - std::vector ratebhp_samples; - for (size_t i = 0; i < flo_samples.size(); ++i) { - ratebhp_samples.push_back( detail::RateBhpPair{flo_samples[i], bhp_flo_samples[i]} ); - } - - const std::array ratebhp_twopoints_ipr {detail::RateBhpPair{flo_bhp_middle, bhp_middle}, - detail::RateBhpPair{flo_bhp_limit, bhp_limit} }; - - double obtain_bhp = 0.; - const bool can_obtain_bhp_with_thp_limit = detail::findIntersectionForBhp(ratebhp_samples, ratebhp_twopoints_ipr, obtain_bhp); - - // \Note: assuming that negative BHP does not make sense - if (can_obtain_bhp_with_thp_limit && obtain_bhp > 0.) { - // getting too high bhp that might cause negative rates (rates in the undesired direction) - if (obtain_bhp >= bhp_safe_limit) { - const std::string msg (" We are getting a too high BHP value from the THP constraint, which may " - " cause problems later "); - OpmLog::info("TOO_HIGH_BHP_FOUND_THP_TARGET", msg); - - const std::string debug_msg = " obtain_bhp " + std::to_string(obtain_bhp) - + " bhp_safe_limit " + std::to_string(bhp_safe_limit) - + " thp limit " + std::to_string(thp_limit); - OpmLog::debug(debug_msg); - } - return obtain_bhp; - } else { - OpmLog::warning("NO_BHP_FOUND_THP_TARGET", " we could not find a bhp value with thp target."); - return -100.; - } -} - void VFPProdProperties::addTable(const VFPProdTable& new_table) { this->m_tables.emplace( new_table.getTableNum(), new_table ); } diff --git a/opm/simulators/wells/VFPProdProperties.hpp b/opm/simulators/wells/VFPProdProperties.hpp index e624dc2cd1..17dfc36189 100644 --- a/opm/simulators/wells/VFPProdProperties.hpp +++ b/opm/simulators/wells/VFPProdProperties.hpp @@ -150,19 +150,6 @@ public: } - /** - * Calculate the Bhp value from the THP target/constraint value - * based on inflow performance relationship and VFP curves - */ - double - calculateBhpWithTHPTarget(const std::vector& ipr_a, - const std::vector& ipr_b, - const double bhp_limit, - const double thp_table_id, - const double thp_limit, - const double alq, - const double dp) const; - protected: // calculate a group bhp values with a group of flo rate values std::vector bhpwithflo(const std::vector& flos,