diff --git a/opm/simulators/wells/GasLiftStage2.hpp b/opm/simulators/wells/GasLiftStage2.hpp index aded0b167..7f33d5b11 100644 --- a/opm/simulators/wells/GasLiftStage2.hpp +++ b/opm/simulators/wells/GasLiftStage2.hpp @@ -99,7 +99,7 @@ namespace Opm void displayWarning_(const std::string &msg); std::tuple getCurrentGroupRates_( const Opm::Group &group); - std::tuple getCurrentGroupRatesRecursive_( + std::array getCurrentGroupRatesRecursive_( const Opm::Group &group); std::tuple getCurrentWellRates_( const std::string &well_name, const std::string &group_name); @@ -127,6 +127,9 @@ namespace Opm const std::string &name, GradInfo &grad, bool increase); void updateGradVector_( const std::string &name, std::vector &grads, double grad); + std::vector updateGlobalGradVector_(const std::vector &grads_global) const; + std::vector localToGlobalGradVector_(const std::vector &grads_local) const; + DeferredLogger &deferred_logger_; const Simulator &ebos_simulator_; @@ -213,7 +216,7 @@ namespace Opm bool checkEcoGradient(const std::string &well_name, double eco_grad); bool checkGasTarget(); bool checkOilTarget(); - void updateRates(const std::string &name, const GradInfo &gi); + void updateRates(const std::string &name); private: }; }; diff --git a/opm/simulators/wells/GasLiftStage2_impl.hpp b/opm/simulators/wells/GasLiftStage2_impl.hpp index 2213d00fa..1ceedadde 100644 --- a/opm/simulators/wells/GasLiftStage2_impl.hpp +++ b/opm/simulators/wells/GasLiftStage2_impl.hpp @@ -100,8 +100,13 @@ void GasLiftStage2:: addOrRemoveALQincrement_(GradMap &grad_map, const std::string well_name, bool add) { + // only applies to wells in the well_state_map (i.e. wells on this rank) + auto it = this->well_state_map_.find(well_name); + if (it == this->well_state_map_.end()) + return; + + GLiftWellState &state = *(it->second.get()); const GradInfo &gi = grad_map.at(well_name); - GLiftWellState &state = *(this->well_state_map_.at(well_name).get()); if (this->debug_) { auto new_alq = gi.alq; auto old_alq = state.alq(); @@ -122,6 +127,11 @@ GasLiftStage2:: calcIncOrDecGrad_( const std::string well_name, const GasLiftSingleWell &gs_well, bool increase) { + + // only applies to wells in the well_state_map (i.e. wells on this rank) + if(this->well_state_map_.count(well_name) == 0) + return std::nullopt; + GLiftWellState &state = *(this->well_state_map_.at(well_name).get()); if (checkRateAlreadyLimited_(state, increase)) { /* @@ -273,23 +283,23 @@ std::tuple GasLiftStage2:: getCurrentGroupRates_(const Opm::Group &group) { + auto rates = getCurrentGroupRatesRecursive_(group); + const auto& comm = ebos_simulator_.vanguard().grid().comm(); + comm.sum(rates.data(), rates.size()); + auto [oil_rate, gas_rate, alq] = rates; if (this->debug_) { - auto [oil_rate, gas_rate, alq] = getCurrentGroupRatesRecursive_(group); const std::string msg = fmt::format( "Current group rates for {} : oil: {}, gas: {}, alq: {}", group.name(), oil_rate, gas_rate, alq); displayDebugMessage2B_(msg); - return {oil_rate, gas_rate, alq}; - } - else { - return getCurrentGroupRatesRecursive_(group); } + + return {oil_rate, gas_rate, alq}; } template -std::tuple -GasLiftStage2:: +std::array GasLiftStage2:: getCurrentGroupRatesRecursive_(const Opm::Group &group) { double oil_rate = 0.0; @@ -304,6 +314,7 @@ getCurrentGroupRatesRecursive_(const Opm::Group &group) gas_rate += sw_gas_rate; alq += sw_alq; } + } else { for (const std::string& group_name : group.groups()) { @@ -317,8 +328,8 @@ getCurrentGroupRatesRecursive_(const Opm::Group &group) // they are added to the lift gas supply rate of the // parent group. const auto gefac = sub_group.getGroupEfficiencyFactor(); - auto [sg_oil_rate, sg_gas_rate, sg_alq] = - getCurrentGroupRatesRecursive_(sub_group); + auto rates = getCurrentGroupRatesRecursive_(sub_group); + auto [sg_oil_rate, sg_gas_rate, sg_alq] = rates; oil_rate += (gefac * sg_oil_rate); gas_rate += (gefac * sg_gas_rate); alq += (gefac * sg_alq); @@ -504,17 +515,23 @@ recalculateGradientAndUpdateData_( // NOTE: We make a copy of the name string instead of taking a reference // since we may have to erase grad_itr (in the "else" condition below) const std::string name = grad_itr->first; - GasLiftSingleWell &gs_well = *(this->stage1_wells_.at(name).get()); - auto grad = calcIncOrDecGrad_(name, gs_well, increase); std::optional old_grad = std::nullopt; - if (grad) { - grad_itr->second = grad->grad; - old_grad = updateGrad_(name, *grad, increase); - } - else { - grads.erase(grad_itr); // NOTE: this invalidates grad_itr - old_grad = deleteGrad_(name, increase); + + // only applies to wells in the well_state_map (i.e. wells on this rank) + // the grads and other grads are synchronized later + if(this->stage1_wells_.count(name) > 0) { + GasLiftSingleWell &gs_well = *(this->stage1_wells_.at(name).get()); + auto grad = calcIncOrDecGrad_(name, gs_well, increase); + if (grad) { + grad_itr->second = grad->grad; + old_grad = updateGrad_(name, *grad, increase); + } + else { + grads.erase(grad_itr); // NOTE: this invalidates grad_itr + old_grad = deleteGrad_(name, increase); + } } + if (old_grad) { // NOTE: Either creates a new item or reassigns // The old incremental gradient becomes the new decremental gradient @@ -574,9 +591,16 @@ redistributeALQ_(std::vector &wells, const Opm::Group &gro std::vector &inc_grads, std::vector &dec_grads) { OptimizeState state {*this, group}; - inc_grads.reserve(wells.size()); - dec_grads.reserve(wells.size()); - state.calculateEcoGradients(wells, inc_grads, dec_grads); + std::vector inc_grads_local; + std::vector dec_grads_local; + inc_grads_local.reserve(wells.size()); + dec_grads_local.reserve(wells.size()); + state.calculateEcoGradients(wells, inc_grads_local, dec_grads_local); + + // the gradients needs to be communicated to all ranks + dec_grads = localToGlobalGradVector_(dec_grads_local); + inc_grads = localToGlobalGradVector_(inc_grads_local); + if (!state.checkAtLeastTwoWells(wells)) { // NOTE: Even though we here in redistributeALQ_() do not use the // economic gradient if there is only a single well, we still @@ -584,7 +608,6 @@ redistributeALQ_(std::vector &wells, const Opm::Group &gro // and will be used by removeSurplusALQ_() later. return; } - assert(wells.size() >= 2); bool stop_iteration = false; while (!stop_iteration && (state.it++ <= this->max_iterations_)) { state.debugShowIterationInfo(); @@ -645,12 +668,13 @@ removeSurplusALQ_(const Opm::Group &group, } SurplusState state {*this, group, oil_rate, gas_rate, alq, min_eco_grad, controls.oil_target, controls.gas_target, max_glift }; + while (!stop_iteration) { if (dec_grads.size() >= 2) { sortGradients_(dec_grads); } auto dec_grad_itr = dec_grads.begin(); - const auto &well_name = dec_grad_itr->first; + const auto well_name = dec_grad_itr->first; auto eco_grad = dec_grad_itr->second; bool remove = false; if (state.checkOilTarget() || state.checkGasTarget() || state.checkALQlimit()) { @@ -665,11 +689,14 @@ removeSurplusALQ_(const Opm::Group &group, if (state.checkEcoGradient(well_name, eco_grad)) remove = true; } if (remove) { - const GradInfo &gi = this->dec_grads_.at(well_name); - state.updateRates(well_name, gi); + state.updateRates(well_name); state.addOrRemoveALQincrement( this->dec_grads_, well_name, /*add=*/false); recalculateGradientAndUpdateData_( - dec_grad_itr, /*increase=*/false, dec_grads, inc_grads); + dec_grad_itr, /*increase=*/false, dec_grads, inc_grads); + + // The dec_grads and inc_grads needs to be syncronized across ranks + dec_grads = updateGlobalGradVector_(dec_grads); + inc_grads = updateGlobalGradVector_(inc_grads); // NOTE: recalculateGradientAndUpdateData_() will remove the current gradient // from dec_grads if it cannot calculate a new decremental gradient. // This will invalidate dec_grad_itr and well_name @@ -766,6 +793,51 @@ updateGradVector_(const std::string &name, std::vector &grads, double // later in getEcoGradients() } +template +std::vector::GradPair> +GasLiftStage2:: +localToGlobalGradVector_(const std::vector &grads_local) const +{ + const auto& comm = ebos_simulator_.vanguard().grid().comm(); + if (comm.size() == 1) + return grads_local; + + std::vector sizes_(comm.size()); + std::vector displ_(comm.size() + 1, 0); + int mySize = grads_local.size(); + comm.allgather(&mySize, 1, sizes_.data()); + std::partial_sum(sizes_.begin(), sizes_.end(), displ_.begin()+1); + using Pair = std::pair; + std::vector grads_global(displ_.back()); + std::vector grads_local_tmp(grads_local.size()); + for (size_t i = 0; i < grads_local.size(); ++i) { + grads_local_tmp[i] = std::make_pair(well_state_.wellNameToGlobalIdx(grads_local[i].first), grads_local[i].second); + } + comm.allgatherv(grads_local_tmp.data(), grads_local_tmp.size(), grads_global.data(), sizes_.data(), displ_.data()); + std::vector grads(grads_global.size()); + for (size_t i = 0; i < grads_global.size(); ++i) { + grads[i] = std::make_pair(well_state_.globalIdxToWellName(grads_global[i].first), grads_global[i].second); + } + return grads; +} + +template +std::vector::GradPair> +GasLiftStage2:: +updateGlobalGradVector_(const std::vector &grads_global) const +{ + const auto& comm = ebos_simulator_.vanguard().grid().comm(); + if (comm.size() == 1) + return grads_global; + + std::vector grads_local; + for (auto itr = grads_global.begin(); itr != grads_global.end(); itr++) { + if (well_state_map_.count(itr->first) > 0) { + grads_local.push_back(*itr); + } + } + return localToGlobalGradVector_(grads_local); +} /*********************************************** * Public methods declared in OptimizeState @@ -799,9 +871,12 @@ bool GasLiftStage2::OptimizeState:: checkAtLeastTwoWells(std::vector &wells) { - if (wells.size() < 2) { + int numberOfwells = wells.size(); + const auto& comm = this->parent.ebos_simulator_.vanguard().grid().comm(); + numberOfwells = comm.sum(numberOfwells); + if (numberOfwells < 2) { const std::string msg = fmt::format( - "skipping: too few wells ({}) to optimize.", wells.size()); + "skipping: too few wells ({}) to optimize.", numberOfwells); displayDebugMessage_(msg); return false; } @@ -871,6 +946,10 @@ recalculateGradients( max_inc_grad_itr, /*increase=*/true, inc_grads, dec_grads); this->parent.recalculateGradientAndUpdateData_( min_dec_grad_itr, /*increase=*/false, dec_grads, inc_grads); + + // The dec_grads and inc_grads needs to be syncronized across ranks + dec_grads = this->parent.updateGlobalGradVector_(dec_grads); + inc_grads = this->parent.updateGlobalGradVector_(inc_grads); } // Take one ALQ increment from well1, and give it to well2 @@ -1010,17 +1089,33 @@ checkOilTarget() template void GasLiftStage2::SurplusState:: -updateRates(const std::string &well_name, const GradInfo &gi) +updateRates(const std::string &well_name) { - GLiftWellState &state = *(this->parent.well_state_map_.at(well_name).get()); - GasLiftSingleWell &gs_well = *(this->parent.stage1_wells_.at(well_name).get()); - const WellInterface &well = gs_well.getStdWell(); - const auto &well_ecl = well.wellEcl(); - double factor = well_ecl.getEfficiencyFactor(); - auto delta_oil = gi.new_oil_rate - state.oilRate(); - this->oil_rate += factor * delta_oil; - auto delta_gas = gi.new_gas_rate - state.gasRate(); - this->gas_rate += factor * delta_gas; - auto delta_alq = gi.alq - state.alq(); - this->alq += factor * delta_alq; + std::array delta = {0.0,0.0,0.0}; + + // compute the delta on wells on own rank + if(this->parent.well_state_map_.count(well_name) > 0 ) { + const GradInfo &gi = this->parent.dec_grads_.at(well_name); + GLiftWellState &state = *(this->parent.well_state_map_.at(well_name).get()); + GasLiftSingleWell &gs_well = *(this->parent.stage1_wells_.at(well_name).get()); + const WellInterface &well = gs_well.getStdWell(); + const auto &well_ecl = well.wellEcl(); + double factor = well_ecl.getEfficiencyFactor(); + auto& [delta_oil, delta_gas, delta_alq] = delta; + delta_oil = factor * (gi.new_oil_rate - state.oilRate()); + delta_gas = factor * (gi.new_gas_rate - state.gasRate()); + delta_alq = factor * (gi.alq - state.alq()); + } + + // and communicate the results + const auto& comm = this->parent.ebos_simulator_.vanguard().grid().comm(); + comm.sum(delta.data(), delta.size()); + + // and update + const auto& [delta_oil, delta_gas, delta_alq] = delta; + this->oil_rate += delta_oil; + this->gas_rate += delta_gas; + this->alq += delta_alq; } + + diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index b96956b09..76d0d5043 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -2692,13 +2692,6 @@ namespace Opm { gliftDebug("checking if GLIFT should be done..", deferred_logger); - std::size_t num_procs = ebos_simulator.gridView().comm().size(); - if (num_procs > 1u) { - const std::string msg = fmt::format(" GLIFT: skipping optimization. " - "Parallel run not supported yet: num procs = {}", num_procs); - deferred_logger.warning(msg); - return false; - } if (!well_state.gliftOptimizationEnabled()) { gliftDebug("Optimization disabled in WellState", deferred_logger); return false; diff --git a/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp b/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp index 16e8d71b2..db876de2d 100644 --- a/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp +++ b/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp @@ -1351,6 +1351,26 @@ namespace Opm do_glift_optimization_ = true; } + int wellNameToGlobalIdx(const std::string &name) { + auto it = wellNameToGlobalIdx_.find(name); + if (it == wellNameToGlobalIdx_.end()) + OPM_THROW(std::logic_error, "Could not find global well idx for well " << name); + + return it->second; + } + + std::string globalIdxToWellName(const int index) { + using Pair = std::pair; + auto it = find_if(wellNameToGlobalIdx_.begin(), wellNameToGlobalIdx_.end(), [index](const Pair & p) { + return p.second == index; + }); + if (it == wellNameToGlobalIdx_.end() ) + { + OPM_THROW(std::logic_error, "Could not find well name for global idx " << index); + } + return it->first; + } + private: std::vector perfphaserates_; std::vector is_producer_; // Size equal to number of local wells.