5#include <dune/common/timer.hh>
8#include <amdis/common/parallel/Collective.hpp>
9#include <amdis/common/parallel/RequestOperations.hpp>
16template <
class GlobalId,
class LI>
17 template <
class Basis>
18void PetscParallelIndexDistribution<GlobalId,LI>::
19update(Basis
const& basis)
25 buildParallelIndexSet(basis, *indexSet_);
27 remoteIndices_->setIndexSets(*indexSet_, *indexSet_, world_);
28 remoteIndices_->template rebuild<true>();
31 globalIndices_.resize(indexSet_->size(), 0);
32 owner_.resize(indexSet_->size(),
false);
35 for (
auto const& ip : *indexSet_) {
36 if (ip.local().attribute() == Attribute::owner) {
37 size_type idx = ip.local();
38 globalIndices_[idx] = localSize_++;
46 Mpi::all_gather(world_, localSize_, sizes_);
49 starts_.resize(world_.size() + 1);
51 std::partial_sum(sizes_.begin(), sizes_.end(), starts_.begin()+1);
52 globalSize_ = starts_.back();
55 for (
auto& i : globalIndices_)
56 i += starts_[world_.rank()];
62 using GlobalAssoc = std::pair<typename IndexSet::GlobalIndex, size_type>;
63 std::vector<std::vector<GlobalAssoc>> sendList(world_.size());
64 std::vector<std::size_t> receiveList(world_.size(), 0);
67 for (
auto const& rim : *remoteIndices_) {
70 auto* sourceRemoteIndexList = rim.second.first;
71 auto* targetRemoteIndexList = rim.second.second;
74 for (
auto const& ri : *sourceRemoteIndexList) {
75 auto const& lip = ri.localIndexPair();
76 Attribute remoteAttr = ri.attribute();
77 Attribute myAttr = lip.local().attribute();
78 if (myAttr == Attribute::owner && remoteAttr != Attribute::owner) {
79 size_type globalIndex = globalIndices_[size_type(lip.local())];
80 sendList[p].push_back({lip.global(), globalIndex});
85 for (
auto const& ri : *targetRemoteIndexList) {
86 auto const& lip = ri.localIndexPair();
87 Attribute remoteAttr = ri.attribute();
88 Attribute myAttr = lip.local().attribute();
89 if (myAttr != Attribute::owner && remoteAttr == Attribute::owner) {
95 assert(ghostSize_ == std::accumulate(receiveList.begin(), receiveList.end(), 0u));
98 std::vector<Mpi::Request> sendRequests;
99 for (
int p = 0; p < world_.size(); ++p) {
100 if (!sendList[p].empty()) {
101 sendRequests.emplace_back( world_.isend(sendList[p], p, tag_) );
106 std::vector<Mpi::Request> recvRequests;
107 std::vector<std::vector<GlobalAssoc>> recvData(world_.size());
108 for (
int p = 0; p < world_.size(); ++p) {
109 if (receiveList[p] > 0)
110 recvRequests.emplace_back( world_.irecv(recvData[p], p, tag_) );
113 Mpi::wait_all(recvRequests.begin(), recvRequests.end());
115 ghostIndices_.reserve(ghostSize_);
116 ghostLocalIndices_.resize(indexSet_->size(), LocalIndex(-1));
119 std::size_t counter = 0;
120 for (
int p = 0; p < world_.size(); ++p) {
121 auto const& data = recvData[p];
122 assert(data.size() == receiveList[p]);
123 for (
auto const& d : data) {
124 typename IndexSet::IndexPair
const& l = indexSet_->at(d.first);
125 assert(!owner_[size_type(l.local())]);
127 globalIndices_[size_type(l.local())] = d.second;
128 ghostIndices_.push_back(d.second);
129 ghostLocalIndices_[size_type(l.local())] = counter++;
132 assert(counter == ghostSize_);
133 assert(ghostSize_ + localSize_ == indexSet_->size());
135 Mpi::wait_all(sendRequests.begin(), sendRequests.end());
136 msg(
"update IndexDistribution need {} sec", t.elapsed());
140template <
class GlobalId,
class LI>
141void PetscParallelIndexDistribution<GlobalId,LI>::
144 int p = world_.rank();
145 std::cout <<
"[" << p <<
"] sizes_.size()=" << sizes_.size() <<
", my_size=" << sizes_[p] << std::endl;
146 std::cout <<
"[" << p <<
"] starts_.size()=" << starts_.size() <<
", my_start=" << starts_[p] << std::endl;
147 std::cout <<
"[" << p <<
"] localSize_=" << localSize_ <<
", globalSize_=" << globalSize_ <<
", ghostSize_=" << ghostSize_ << std::endl;
148 std::cout <<
"[" << p <<
"] globalIndices_.size()=" << globalIndices_.size() << std::endl;
149 std::cout <<
"[" << p <<
"] ghostIndices_.size()=" << ghostIndices_.size() << std::endl;
150 std::cout <<
"[" << p <<
"] ghostLocalIndices_.size()=" << ghostLocalIndices_.size() << std::endl;
151 std::cout <<
"[" << p <<
"] owner_.size()=" << owner_.size() << std::endl;