From 35218bf0427010796799eca4b9f20af30eee913a Mon Sep 17 00:00:00 2001 From: Markus Blatt Date: Tue, 8 Dec 2020 19:02:33 +0100 Subject: [PATCH] Added possibility to communicate values from perforations below. --- .../wells/BlackoilWellModel_impl.hpp | 2 +- opm/simulators/wells/ParallelWellInfo.cpp | 88 ++++++++++++++----- opm/simulators/wells/ParallelWellInfo.hpp | 44 +++++++--- tests/test_parallelwellinfo.cpp | 68 +++++++++----- 4 files changed, 147 insertions(+), 55 deletions(-) diff --git a/opm/simulators/wells/BlackoilWellModel_impl.hpp b/opm/simulators/wells/BlackoilWellModel_impl.hpp index 02fac730e..9b395f5c4 100644 --- a/opm/simulators/wells/BlackoilWellModel_impl.hpp +++ b/opm/simulators/wells/BlackoilWellModel_impl.hpp @@ -467,7 +467,7 @@ namespace Opm { // Clear the communication data structures for above values. for(auto&& pinfo : local_parallel_well_info_) { - pinfo->clearCommunicateAbove(); + pinfo->clearCommunicateAboveBelow(); } } diff --git a/opm/simulators/wells/ParallelWellInfo.cpp b/opm/simulators/wells/ParallelWellInfo.cpp index 5e8964406..03b879623 100644 --- a/opm/simulators/wells/ParallelWellInfo.cpp +++ b/opm/simulators/wells/ParallelWellInfo.cpp @@ -43,13 +43,13 @@ struct CommPolicy namespace Opm { -CommunicateAbove::CommunicateAbove([[maybe_unused]] const Communication& comm) +CommunicateAboveBelow::CommunicateAboveBelow([[maybe_unused]] const Communication& comm) #if HAVE_MPI : comm_(comm), interface_(comm_) #endif {} -void CommunicateAbove::clear() +void CommunicateAboveBelow::clear() { #if HAVE_MPI above_indices_ = {}; @@ -60,7 +60,7 @@ void CommunicateAbove::clear() count_ = 0; } -void CommunicateAbove::beginReset() +void CommunicateAboveBelow::beginReset() { clear(); #if HAVE_MPI @@ -72,7 +72,7 @@ void CommunicateAbove::beginReset() #endif } -void CommunicateAbove::endReset() +void CommunicateAboveBelow::endReset() { #if HAVE_MPI if (comm_.size() > 1) @@ -104,9 +104,9 @@ struct CopyGatherScatter }; -std::vector CommunicateAbove::communicate(double first_above, - const double* current, - std::size_t size) +std::vector CommunicateAboveBelow::communicateAbove(double first_above, + const double* current, + std::size_t size) { std::vector above(size, first_above); @@ -131,10 +131,37 @@ std::vector CommunicateAbove::communicate(double first_above, } return above; } +std::vector CommunicateAboveBelow::communicateBelow(double last_below, + const double* current, + std::size_t size) +{ + std::vector below(size, last_below); -void CommunicateAbove::pushBackEclIndex([[maybe_unused]] int above, - [[maybe_unused]] int current, - [[maybe_unused]] bool isOwner) +#if HAVE_MPI + if (comm_.size() > 1) + { + auto belowData = below.data(); + // Ugly const_cast needed as my compiler says, that + // passing const double*& and double* as parameter is + // incompatible with function decl template backward(Data&, const Data&) + // That would need the first argument to be double* const& + communicator_.backward(belowData, const_cast(current)); + } + else +#endif + { + if (below.size() > 1) + { + // No comunication needed, just copy. + std::copy(current+1, current + below.size(), below.begin()); + } + } + return below; +} + +void CommunicateAboveBelow::pushBackEclIndex([[maybe_unused]] int above, + [[maybe_unused]] int current, + [[maybe_unused]] bool isOwner) { #if HAVE_MPI if (comm_.size() > 1) @@ -175,8 +202,8 @@ ParallelWellInfo::ParallelWellInfo(const std::string& name, : name_(name), hasLocalCells_ (hasLocalCells), isOwner_(true), rankWithFirstPerf_(-1), comm_(new Communication(Dune::MPIHelper::getLocalCommunicator())), - commAbove_(new CommunicateAbove(*comm_)) - {} + commAboveBelow_(new CommunicateAboveBelow(*comm_)) +{} ParallelWellInfo::ParallelWellInfo(const std::pair& well_info, [[maybe_unused]] Communication allComm) @@ -191,7 +218,7 @@ ParallelWellInfo::ParallelWellInfo(const std::pair& well_info, #else comm_.reset(new Communication(Dune::MPIHelper::getLocalCommunicator())); #endif - commAbove_.reset(new CommunicateAbove(*comm_)); + commAboveBelow_.reset(new CommunicateAboveBelow(*comm_)); isOwner_ = (comm_->rank() == 0); } @@ -207,38 +234,53 @@ void ParallelWellInfo::communicateFirstPerforation(bool hasFirst) void ParallelWellInfo::pushBackEclIndex(int above, int current) { - commAbove_->pushBackEclIndex(above, current); + commAboveBelow_->pushBackEclIndex(above, current); } void ParallelWellInfo::beginReset() { - commAbove_->beginReset(); + commAboveBelow_->beginReset(); } void ParallelWellInfo::endReset() { - commAbove_->beginReset(); + commAboveBelow_->beginReset(); } -void ParallelWellInfo::clearCommunicateAbove() +void ParallelWellInfo::clearCommunicateAboveBelow() { - commAbove_->clear(); + commAboveBelow_->clear(); } std::vector ParallelWellInfo::communicateAboveValues(double zero_value, const double* current_values, std::size_t size) const { - return commAbove_->communicate(zero_value, current_values, - size); + return commAboveBelow_->communicateAbove(zero_value, current_values, + size); } std::vector ParallelWellInfo::communicateAboveValues(double zero_value, const std::vector& current_values) const { - return commAbove_->communicate(zero_value, current_values.data(), - current_values.size()); + return commAboveBelow_->communicateAbove(zero_value, current_values.data(), + current_values.size()); +} + +std::vector ParallelWellInfo::communicateBelowValues(double last_value, + const double* current_values, + std::size_t size) const +{ + return commAboveBelow_->communicateBelow(last_value, current_values, + size); +} + +std::vector ParallelWellInfo::communicateBelowValues(double last_value, + const std::vector& current_values) const +{ + return commAboveBelow_->communicateBelow(last_value, current_values.data(), + current_values.size()); } bool operator<(const ParallelWellInfo& well1, const ParallelWellInfo& well2) @@ -294,7 +336,7 @@ bool operator!=(const ParallelWellInfo& well, const std::pair } CheckDistributedWellConnections::CheckDistributedWellConnections(const Well& well, - const ParallelWellInfo& info) + const ParallelWellInfo& info) : well_(well), pwinfo_(info) { foundConnections_.resize(well.getConnections().size(), 0); diff --git a/opm/simulators/wells/ParallelWellInfo.hpp b/opm/simulators/wells/ParallelWellInfo.hpp index 62a29b6dd..1fbbcdaa0 100644 --- a/opm/simulators/wells/ParallelWellInfo.hpp +++ b/opm/simulators/wells/ParallelWellInfo.hpp @@ -33,9 +33,9 @@ namespace Opm { -/// \brief Class to facilitate getting values associated with the above perforation +/// \brief Class to facilitate getting values associated with the above/below perforation /// -class CommunicateAbove +class CommunicateAboveBelow { enum Attribute { owner=1, overlap=2 @@ -53,7 +53,7 @@ class CommunicateAbove #endif public: - explicit CommunicateAbove(const Communication& comm); + explicit CommunicateAboveBelow(const Communication& comm); /// \brief Adds information about original index of the perforations in ECL Schedule. /// /// \warning Theses indices need to be push in the same order as they @@ -81,10 +81,18 @@ public: /// \param current C-array of the values at the perforations /// \param size The size of the C-array and the returned vector /// \return a vector containing the values for the perforation above. - std::vector communicate(double first_value, - const double* current, - std::size_t size); + std::vector communicateAbove(double first_value, + const double* current, + std::size_t size); + /// \brief Creates an array of values for the perforation below. + /// \param first_value Value to use for above of the first perforation + /// \param current C-array of the values at the perforations + /// \param size The size of the C-array and the returned vector + /// \return a vector containing the values for the perforation above. + std::vector communicateBelow(double first_value, + const double* current, + std::size_t size); private: #if HAVE_MPI Communication comm_; @@ -157,15 +165,31 @@ public: /// \param current C-array of the values at the perforations /// \param size The size of the C-array and the returned vector /// \return a vector containing the values for the perforation above. - std::vector communicateAboveValues(double zero_value, + std::vector communicateAboveValues(double first_value, const double* current, std::size_t size) const; /// \brief Creates an array of values for the perforation above. /// \param first_value Value to use for above of the first perforation /// \param current vector of current values - std::vector communicateAboveValues(double zero_value, + std::vector communicateAboveValues(double first_value, const std::vector& current) const; + + /// \brief Creates an array of values for the perforation below. + /// \param last_value Value to use for below of the last perforation + /// \param current C-array of the values at the perforations + /// \param size The size of the C-array and the returned vector + /// \return a vector containing the values for the perforation above. + std::vector communicateBelowValues(double last_value, + const double* current, + std::size_t size) const; + + /// \brief Creates an array of values for the perforation above. + /// \param last_value Value to use for below of the last perforation + /// \param current vector of current values + std::vector communicateBelowValues(double last_value, + const std::vector& current) const; + /// \brief Adds information about the ecl indices of the perforations. /// /// \warning Theses indices need to be push in the same order as they @@ -210,7 +234,7 @@ public: } /// \brief Free data of communication data structures. - void clearCommunicateAbove(); + void clearCommunicateAboveBelow(); private: /// \brief Deleter that also frees custom MPI communicators @@ -234,7 +258,7 @@ private: std::unique_ptr comm_; /// \brief used to communicate the values for the perforation above. - std::unique_ptr commAbove_; + std::unique_ptr commAboveBelow_; }; /// \brief Class checking that all connections are on active cells diff --git a/tests/test_parallelwellinfo.cpp b/tests/test_parallelwellinfo.cpp index 0cfe06469..8e6fa2766 100644 --- a/tests/test_parallelwellinfo.cpp +++ b/tests/test_parallelwellinfo.cpp @@ -153,64 +153,76 @@ BOOST_AUTO_TEST_CASE(ParallelWellComparison) } -BOOST_AUTO_TEST_CASE(CommunicateAboveSelf) +BOOST_AUTO_TEST_CASE(CommunicateAboveBelowSelf) { auto comm = Dune::MPIHelper::getLocalCommunicator(); - Opm::CommunicateAbove commAbove{ comm }; + Opm::CommunicateAboveBelow commAboveBelow{ comm }; for(std::size_t count=0; count < 2; ++count) { std::vector eclIndex = {0, 1, 2, 3, 7 , 8, 10, 11}; std::vector current(eclIndex.size()); std::transform(eclIndex.begin(), eclIndex.end(), current.begin(), [](double v){ return 1+10.0*v;}); - commAbove.beginReset(); + commAboveBelow.beginReset(); for (std::size_t i = 0; i < current.size(); ++i) { if (i==0) - commAbove.pushBackEclIndex(-1, eclIndex[i]); + commAboveBelow.pushBackEclIndex(-1, eclIndex[i]); else - commAbove.pushBackEclIndex(eclIndex[i-1], eclIndex[i]); + commAboveBelow.pushBackEclIndex(eclIndex[i-1], eclIndex[i]); } - commAbove.endReset(); - auto above = commAbove.communicate(-10.0, current.data(), current.size()); + commAboveBelow.endReset(); + auto above = commAboveBelow.communicateAbove(-10.0, current.data(), current.size()); BOOST_CHECK(above[0]==-10.0); BOOST_CHECK(above.size() == current.size()); auto a = above.begin()+1; std::for_each(current.begin(), current.begin() + (current.size()-1), [&a](double v){ BOOST_CHECK(*(a++) == v);}); + auto below = commAboveBelow.communicateBelow(-10.0, current.data(), current.size()); + BOOST_CHECK(below.back() == -10.0); + BOOST_CHECK(below.size() == current.size()); + auto b = below.begin(); + std::for_each(current.begin()+1, current.end(), + [&b](double v){ BOOST_CHECK(*(b++) == v);}); } } -BOOST_AUTO_TEST_CASE(CommunicateAboveSelf1) +BOOST_AUTO_TEST_CASE(CommunicateAboveBelowSelf1) { auto comm = Dune::MPIHelper::getLocalCommunicator(); - Opm::CommunicateAbove commAbove{ comm }; + Opm::CommunicateAboveBelow commAboveBelow{ comm }; for(std::size_t count=0; count < 2; ++count) { std::vector eclIndex = {0}; std::vector current(eclIndex.size()); std::transform(eclIndex.begin(), eclIndex.end(), current.begin(), [](double v){ return 1+10.0*v;}); - commAbove.beginReset(); + commAboveBelow.beginReset(); for (std::size_t i = 0; i < current.size(); ++i) { if (i==0) - commAbove.pushBackEclIndex(-1, eclIndex[i]); + commAboveBelow.pushBackEclIndex(-1, eclIndex[i]); else - commAbove.pushBackEclIndex(eclIndex[i-1], eclIndex[i]); + commAboveBelow.pushBackEclIndex(eclIndex[i-1], eclIndex[i]); } - commAbove.endReset(); - auto above = commAbove.communicate(-10.0, current.data(), current.size()); + commAboveBelow.endReset(); + auto above = commAboveBelow.communicateAbove(-10.0, current.data(), current.size()); BOOST_CHECK(above[0]==-10.0); BOOST_CHECK(above.size() == current.size()); auto a = above.begin()+1; std::for_each(current.begin(), current.begin() + (current.size()-1), [&a](double v){ BOOST_CHECK(*(a++) == v);}); + auto below = commAboveBelow.communicateBelow(-10.0, current.data(), current.size()); + BOOST_CHECK(below.back() == -10.0); + BOOST_CHECK(below.size() == current.size()); + auto b = below.begin(); + std::for_each(current.begin()+1, current.end(), + [&b](double v){ BOOST_CHECK(*(b++) == v);}); } } -BOOST_AUTO_TEST_CASE(CommunicateAboveParalle) +BOOST_AUTO_TEST_CASE(CommunicateAboveBelowParallel) { using MPIComm = typename Dune::MPIHelper::MPICommunicator; #if DUNE_VERSION_NEWER(DUNE_COMMON, 2, 7) @@ -220,7 +232,7 @@ BOOST_AUTO_TEST_CASE(CommunicateAboveParalle) #endif auto comm = Communication(Dune::MPIHelper::getCommunicator()); - Opm::CommunicateAbove commAbove{ comm }; + Opm::CommunicateAboveBelow commAboveBelow{ comm }; for(std::size_t count=0; count < 2; ++count) { std::vector globalEclIndex = {0, 1, 2, 3, 7 , 8, 10, 11}; @@ -244,23 +256,23 @@ BOOST_AUTO_TEST_CASE(CommunicateAboveParalle) std::vector current(3); - commAbove.beginReset(); + commAboveBelow.beginReset(); for (std::size_t i = 0; i < current.size(); ++i) { auto gi = comm.rank() + comm.size() * i; if (gi==0) { - commAbove.pushBackEclIndex(-1, globalEclIndex[gi]); + commAboveBelow.pushBackEclIndex(-1, globalEclIndex[gi]); } else { - commAbove.pushBackEclIndex(globalEclIndex[gi-1], globalEclIndex[gi]); + commAboveBelow.pushBackEclIndex(globalEclIndex[gi-1], globalEclIndex[gi]); } current[i] = globalCurrent[gi]; } - commAbove.endReset(); - auto above = commAbove.communicate(-10.0, current.data(), current.size()); + commAboveBelow.endReset(); + auto above = commAboveBelow.communicateAbove(-10.0, current.data(), current.size()); if (comm.rank() == 0) BOOST_CHECK(above[0]==-10.0); @@ -274,5 +286,19 @@ BOOST_AUTO_TEST_CASE(CommunicateAboveParalle) BOOST_CHECK(above[i]==globalCurrent[gi-1]); } } + auto below = commAboveBelow.communicateBelow(-10.0, current.data(), current.size()); + if (comm.rank() == comm.size() - 1) + BOOST_CHECK(below.back() == -10.0); + + BOOST_CHECK(below.size() == current.size()); + + for (std::size_t i = 0; i < current.size(); ++i) + { + auto gi = comm.rank() + comm.size() * i; + if (gi < globalCurrent.size() - 1) + { + BOOST_CHECK(below[i]==globalCurrent[gi+1]); + } + } } }