3 #ifndef DUNE_INDICESSYNCER_HH
4 #define DUNE_INDICESSYNCER_HH
44 typedef T ParallelIndexSet;
53 typedef typename ParallelIndexSet::LocalIndex::Attribute Attribute;
58 typedef Dune::RemoteIndices<ParallelIndexSet> RemoteIndices;
69 IndicesSyncer(ParallelIndexSet& indexSet,
70 RemoteIndices& remoteIndices);
94 void sync(T1& numberer);
99 ParallelIndexSet& indexSet_;
102 RemoteIndices& remoteIndices_;
108 char* receiveBuffer_;
111 std::size_t* sendBufferSizes_;
114 int receiveBufferSize_;
119 struct MessageInformation
136 class DefaultNumberer
144 std::size_t operator()(
const GlobalIndex& global)
152 MPI_Datatype datatype_;
161 typedef SLList<std::pair<GlobalIndex,Attribute>,
typename RemoteIndices::Allocator> GlobalIndexList;
164 typedef typename GlobalIndexList::ModifyIterator GlobalIndexModifier;
173 typedef std::map<int, GlobalIndexList> GlobalIndicesMap;
183 GlobalIndicesMap globalMap_;
188 typedef SLList<bool, typename RemoteIndices::Allocator> BoolList;
193 typedef typename BoolList::iterator BoolIterator;
196 typedef typename BoolList::ModifyIterator BoolListModifier;
199 typedef std::map<int,BoolList> BoolMap;
208 std::map<int,MessageInformation> infoSend_;
211 typedef typename RemoteIndices::RemoteIndexList RemoteIndexList;
214 typedef typename RemoteIndexList::ModifyIterator RemoteIndexModifier;
217 typedef Dune::RemoteIndex<GlobalIndex,Attribute> RemoteIndex;
220 typedef typename RemoteIndexList::iterator RemoteIndexIterator;
223 typedef typename RemoteIndexList::const_iterator ConstRemoteIndexIterator;
226 typedef std::tuple<RemoteIndexModifier,GlobalIndexModifier,BoolListModifier,
227 const ConstRemoteIndexIterator> IteratorTuple;
238 friend class IndicesSyncer<T>;
249 Iterators(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
260 Iterators& operator++();
267 void insert(
const RemoteIndex& index,
268 const std::pair<GlobalIndex,Attribute>& global);
274 RemoteIndex& remoteIndex()
const;
280 std::pair<GlobalIndex,Attribute>& globalIndexPair()
const;
282 Attribute& attribute()
const;
300 void reset(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
308 bool isNotAtEnd()
const;
315 bool isAtEnd()
const;
327 IteratorTuple iterators_;
331 typedef std::map<int,Iterators> IteratorsMap;
344 IteratorsMap iteratorsMap_;
347 void calculateMessageSizes();
356 void packAndSend(
int destination,
char* buffer, std::size_t bufferSize, MPI_Request& req);
362 template<
typename T1>
363 void recvAndUnpack(T1& numberer);
368 void registerMessageDatatype();
373 void insertIntoRemoteIndexList(
int process,
374 const std::pair<GlobalIndex,Attribute>& global,
380 void resetIteratorsMap();
396 bool checkReset(
const Iterators& iterators, RemoteIndexList& rlist, GlobalIndexList& gList,
400 template<
typename TG,
typename TA>
401 bool operator<(
const IndexPair<TG,ParallelLocalIndex<TA> >& i1,
402 const std::pair<TG,TA>& i2)
404 return i1.global() < i2.first ||
405 (i1.global() == i2.first && i1.local().attribute()<i2.second);
408 template<
typename TG,
typename TA>
409 bool operator<(
const std::pair<TG,TA>& i1,
410 const IndexPair<TG,ParallelLocalIndex<TA> >& i2)
412 return i1.first < i2.global() ||
413 (i1.first == i2.global() && i1.second<i2.local().attribute());
416 template<
typename TG,
typename TA>
417 bool operator==(
const IndexPair<TG,ParallelLocalIndex<TA> >& i1,
418 const std::pair<TG,TA>& i2)
420 return (i1.global() == i2.first && i1.local().attribute()==i2.second);
423 template<
typename TG,
typename TA>
424 bool operator!=(
const IndexPair<TG,ParallelLocalIndex<TA> >& i1,
425 const std::pair<TG,TA>& i2)
427 return (i1.global() != i2.first || i1.local().attribute()!=i2.second);
430 template<
typename TG,
typename TA>
432 const IndexPair<TG,ParallelLocalIndex<TA> >& i1)
434 return (i1.global() == i2.first && i1.local().attribute()==i2.second);
437 template<
typename TG,
typename TA>
439 const IndexPair<TG,ParallelLocalIndex<TA> >& i1)
441 return (i1.global() != i2.first || i1.local().attribute()!=i2.second);
460 template<
typename T,
typename A,
typename A1>
461 void storeGlobalIndicesOfRemoteIndices(std::map<
int,SLList<std::pair<typename T::GlobalIndex, typename T::LocalIndex::Attribute>,A> >& globalMap,
462 const RemoteIndices<T,A1>& remoteIndices)
464 typedef typename RemoteIndices<T,A1>::const_iterator RemoteIterator;
466 for(RemoteIterator remote = remoteIndices.begin(), end =remoteIndices.end(); remote != end; ++remote) {
467 typedef typename RemoteIndices<T,A1>::RemoteIndexList RemoteIndexList;
468 typedef typename RemoteIndexList::const_iterator RemoteIndexIterator;
469 typedef SLList<std::pair<typename T::GlobalIndex, typename T::LocalIndex::Attribute>,A> GlobalIndexList;
470 GlobalIndexList& global = globalMap[remote->first];
471 RemoteIndexList& rList = *(remote->second.first);
473 for(RemoteIndexIterator index = rList.begin(), riEnd = rList.end();
474 index != riEnd; ++index) {
475 global.push_back(std::make_pair(index->localIndexPair().global(),
476 index->localIndexPair().local().attribute()));
489 template<
typename T,
typename A,
typename A1>
490 inline void repairLocalIndexPointers(std::map<
int,
491 SLList<std::pair<
typename T::GlobalIndex,
492 typename T::LocalIndex::Attribute>,A> >& globalMap,
493 RemoteIndices<T,A1>& remoteIndices,
496 typedef typename RemoteIndices<T,A1>::RemoteIndexMap::iterator RemoteIterator;
497 typedef typename RemoteIndices<T,A1>::RemoteIndexList::iterator RemoteIndexIterator;
498 typedef typename T::GlobalIndex GlobalIndex;
499 typedef typename T::LocalIndex::Attribute Attribute;
500 typedef std::pair<GlobalIndex,Attribute> GlobalIndexPair;
501 typedef SLList<GlobalIndexPair,A> GlobalIndexPairList;
502 typedef typename GlobalIndexPairList::iterator GlobalIndexIterator;
504 assert(globalMap.size()==
static_cast<std::size_t
>(remoteIndices.neighbours()));
506 typename std::map<int,GlobalIndexPairList>::iterator global = globalMap.begin();
507 RemoteIterator end = remoteIndices.remoteIndices_.end();
509 for(RemoteIterator remote = remoteIndices.remoteIndices_.begin(); remote != end; ++remote, ++global) {
510 typedef typename T::const_iterator IndexIterator;
512 assert(remote->first==global->first);
513 assert(remote->second.first->size() == global->second.size());
515 RemoteIndexIterator riEnd = remote->second.first->end();
516 RemoteIndexIterator rIndex = remote->second.first->begin();
517 GlobalIndexIterator gIndex = global->second.begin();
518 IndexIterator index = indexSet.begin();
520 assert(rIndex==riEnd || gIndex != global->second.end());
521 while(rIndex != riEnd) {
523 assert(gIndex != global->second.end());
525 while(!(index->global() == gIndex->first
526 && index->local().attribute() == gIndex->second)) {
531 if (index->global() > gIndex->first) {
532 index=indexSet.begin();
536 assert(index != indexSet.end() && *index == *gIndex);
538 rIndex->localIndex_ = &(*index);
544 remoteIndices.sourceSeqNo_ = remoteIndices.source_->seqNo();
545 remoteIndices.destSeqNo_ = remoteIndices.target_->seqNo();
549 IndicesSyncer<T>::IndicesSyncer(ParallelIndexSet& indexSet,
550 RemoteIndices& remoteIndices)
551 : indexSet_(indexSet), remoteIndices_(remoteIndices)
554 assert(remoteIndices.source_ == remoteIndices.target_);
555 assert(remoteIndices.source_ == &indexSet);
556 MPI_Comm_rank(remoteIndices_.communicator(), &rank_);
560 IndicesSyncer<T>::Iterators::Iterators(RemoteIndexList& remoteIndices,
561 GlobalIndexList& globalIndices,
563 : iterators_(remoteIndices.beginModify(), globalIndices.beginModify(),
564 booleans.beginModify(), remoteIndices.end())
568 IndicesSyncer<T>::Iterators::Iterators()
573 inline typename IndicesSyncer<T>::Iterators& IndicesSyncer<T>::Iterators::operator++()
575 ++(std::get<0>(iterators_));
576 ++(std::get<1>(iterators_));
577 ++(std::get<2>(iterators_));
582 inline void IndicesSyncer<T>::Iterators::insert(
const RemoteIndex & index,
583 const std::pair<GlobalIndex,Attribute>& global)
585 std::get<0>(iterators_).insert(index);
586 std::get<1>(iterators_).insert(global);
587 std::get<2>(iterators_).insert(
false);
591 inline typename IndicesSyncer<T>::RemoteIndex&
592 IndicesSyncer<T>::Iterators::remoteIndex()
const
594 return *(std::get<0>(iterators_));
598 inline std::pair<typename IndicesSyncer<T>::GlobalIndex,
typename IndicesSyncer<T>::Attribute>&
599 IndicesSyncer<T>::Iterators::globalIndexPair()
const
601 return *(std::get<1>(iterators_));
605 inline bool IndicesSyncer<T>::Iterators::isOld()
const
607 return *(std::get<2>(iterators_));
611 inline void IndicesSyncer<T>::Iterators::reset(RemoteIndexList& remoteIndices,
612 GlobalIndexList& globalIndices,
615 std::get<0>(iterators_) = remoteIndices.beginModify();
616 std::get<1>(iterators_) = globalIndices.beginModify();
617 std::get<2>(iterators_) = booleans.beginModify();
621 inline bool IndicesSyncer<T>::Iterators::isNotAtEnd()
const
623 return std::get<0>(iterators_) != std::get<3>(iterators_);
627 inline bool IndicesSyncer<T>::Iterators::isAtEnd()
const
629 return std::get<0>(iterators_) == std::get<3>(iterators_);
633 void IndicesSyncer<T>::registerMessageDatatype()
635 MPI_Datatype type[2] = {MPI_INT, MPI_INT};
636 int blocklength[2] = {1,1};
637 MPI_Aint displacement[2];
641 MessageInformation message;
643 MPI_Get_address( &(message.publish), displacement);
644 MPI_Get_address( &(message.pairs), displacement+1);
647 MPI_Get_address(&message, &base);
648 displacement[0] -= base;
649 displacement[1] -= base;
651 MPI_Type_create_struct( 2, blocklength, displacement, type, &datatype_);
652 MPI_Type_commit(&datatype_);
656 void IndicesSyncer<T>::calculateMessageSizes()
658 typedef typename ParallelIndexSet::const_iterator IndexIterator;
659 typedef CollectiveIterator<T,typename RemoteIndices::Allocator> CollectiveIterator;
661 IndexIterator iEnd = indexSet_.end();
662 CollectiveIterator collIter = remoteIndices_.template iterator<true>();
664 for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index) {
665 collIter.advance(index->global(), index->local().attribute());
670 typedef typename CollectiveIterator::iterator ValidIterator;
671 ValidIterator end = collIter.end();
674 for(ValidIterator valid = collIter.begin(); valid != end; ++valid) {
679 Dune::dverb<<rank_<<
": publishing "<<knownRemote<<
" for index "<<index->global()<<
" for processes ";
682 for(ValidIterator valid = collIter.begin(); valid != end; ++valid) {
683 ++(infoSend_[valid.process()].publish);
684 (infoSend_[valid.process()].pairs) += knownRemote;
686 Dune::dverb<<
"(publish="<<infoSend_[valid.process()].publish<<
", pairs="<<infoSend_[valid.process()].pairs
693 typedef typename std::map<int,MessageInformation>::const_iterator
696 const MessageIterator end = infoSend_.end();
701 MessageInformation dummy;
703 MessageIterator messageIter= infoSend_.begin();
705 typedef typename RemoteIndices::RemoteIndexMap::const_iterator RemoteIterator;
706 const RemoteIterator rend = remoteIndices_.end();
709 for(RemoteIterator remote = remoteIndices_.begin(); remote != rend; ++remote, ++neighbour) {
710 MessageInformation* message;
711 MessageInformation recv;
713 if(messageIter != end && messageIter->first==remote->first) {
715 message =
const_cast<MessageInformation*
>(&(messageIter->second));
721 sendBufferSizes_[neighbour]=0;
724 MPI_Pack_size(1, MPI_INT,remoteIndices_.communicator(), &tsize);
725 sendBufferSizes_[neighbour] += tsize;
727 for(
int i=0; i < message->publish; ++i) {
729 MPI_Pack_size(1, MPITraits<GlobalIndex>::getType(), remoteIndices_.communicator(), &tsize);
730 sendBufferSizes_[neighbour] += tsize;
732 MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
733 sendBufferSizes_[neighbour] += tsize;
735 MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
736 sendBufferSizes_[neighbour] += tsize;
738 for(
int i=0; i < message->pairs; ++i) {
740 MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
741 sendBufferSizes_[neighbour] += tsize;
743 MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
744 sendBufferSizes_[neighbour] += tsize;
747 Dune::dverb<<rank_<<
": Buffer (neighbour="<<remote->first<<
") size is "<< sendBufferSizes_[neighbour]<<
" for publish="<<message->publish<<
" pairs="<<message->pairs<<std::endl;
753 inline void IndicesSyncer<T>::sync()
755 DefaultNumberer numberer;
760 template<
typename T1>
761 void IndicesSyncer<T>::sync(T1& numberer)
769 typedef typename RemoteIndices::RemoteIndexMap::const_iterator
772 const RemoteIterator end = remoteIndices_.end();
776 std::size_t noOldNeighbours = remoteIndices_.neighbours();
777 int* oldNeighbours =
new int[noOldNeighbours];
778 sendBufferSizes_ =
new std::size_t[noOldNeighbours];
779 std::size_t neighbourI = 0;
781 for(RemoteIterator remote = remoteIndices_.begin(); remote != end; ++remote, ++neighbourI) {
782 typedef typename RemoteIndices::RemoteIndexList::const_iterator
785 oldNeighbours[neighbourI] = remote->first;
788 assert(remote->second.first==remote->second.second);
790 RemoteIndexList& rList = *(remote->second.first);
793 GlobalIndexList& global = globalMap_[remote->first];
794 BoolList& added = oldMap_[remote->first];
795 RemoteIndexIterator riEnd = rList.end();
797 for(RemoteIndexIterator index = rList.begin();
798 index != riEnd; ++index) {
799 global.push_back(std::make_pair(index->localIndexPair().global(),
800 index->localIndexPair().local().attribute()));
801 added.push_back(
true);
804 Iterators iterators(rList, global, added);
805 iteratorsMap_.insert(std::make_pair(remote->first, iterators));
806 assert(checkReset(iteratorsMap_[remote->first], rList,global,added));
810 calculateMessageSizes();
813 receiveBufferSize_=1;
814 sendBuffers_ =
new char*[noOldNeighbours];
816 for(std::size_t i=0; i<noOldNeighbours; ++i) {
817 sendBuffers_[i] =
new char[sendBufferSizes_[i]];
818 receiveBufferSize_ =
std::max(receiveBufferSize_,
static_cast<int>(sendBufferSizes_[i]));
821 receiveBuffer_=
new char[receiveBufferSize_];
823 indexSet_.beginResize();
827 for(std::size_t i = 0; i<noOldNeighbours; ++i)
832 MPI_Request* requests =
new MPI_Request[noOldNeighbours];
833 MPI_Status* statuses =
new MPI_Status[noOldNeighbours];
836 for(std::size_t i = 0; i<noOldNeighbours; ++i)
837 packAndSend(oldNeighbours[i], sendBuffers_[i], sendBufferSizes_[i], requests[i]);
840 for(std::size_t i = 0; i<noOldNeighbours; ++i)
841 recvAndUnpack(numberer);
848 delete[] receiveBuffer_;
852 if(MPI_SUCCESS!=MPI_Waitall(noOldNeighbours, requests, statuses)) {
853 std::cerr<<
": MPI_Error occurred while sending message"<<std::endl;
854 for(std::size_t i=0; i< noOldNeighbours; i++)
855 if(MPI_SUCCESS!=statuses[i].MPI_ERROR)
856 std::cerr<<
"Destination "<<statuses[i].MPI_SOURCE<<
" error code: "<<statuses[i].MPI_ERROR<<std::endl;
862 for(std::size_t i=0; i<noOldNeighbours; ++i)
863 delete[] sendBuffers_[i];
865 delete[] sendBuffers_;
866 delete[] sendBufferSizes_;
869 iteratorsMap_.clear();
871 indexSet_.endResize();
873 delete[] oldNeighbours;
875 repairLocalIndexPointers(globalMap_, remoteIndices_, indexSet_);
881 remoteIndices_.sourceSeqNo_ = remoteIndices_.destSeqNo_ = indexSet_.seqNo();
885 void IndicesSyncer<T>::packAndSend(
int destination,
char* buffer, std::size_t bufferSize, MPI_Request& request)
887 typedef typename ParallelIndexSet::const_iterator IndexIterator;
889 IndexIterator iEnd = indexSet_.end();
894 assert(checkReset());
897 MPI_Pack(&(infoSend_[destination].publish), 1, MPI_INT, buffer, bufferSize, &bpos,
898 remoteIndices_.communicator());
900 for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index) {
902 typedef typename IteratorsMap::iterator Iterator;
903 Iterator iteratorsEnd = iteratorsMap_.end();
906 for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators) {
907 while(iterators->second.isNotAtEnd() &&
908 iterators->second.globalIndexPair().first < index->global())
909 ++(iterators->second);
910 assert(!iterators->second.isNotAtEnd() || iterators->second.globalIndexPair().first >= index->global());
917 bool knownRemote =
false;
919 for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
921 std::pair<GlobalIndex,Attribute> p;
922 if (iterators->second.isNotAtEnd())
924 p = iterators->second.globalIndexPair();
927 if(iterators->second.isNotAtEnd() && iterators->second.isOld()
928 && iterators->second.globalIndexPair().first == index->global()) {
930 if(destination == iterators->first)
939 Dune::dverb<<rank_<<
": sending "<<indices<<
" for index "<<index->global()<<
" to "<<destination<<std::endl;
943 MPI_Pack(
const_cast<GlobalIndex*
>(&(index->global())), 1, MPITraits<GlobalIndex>::getType(), buffer, bufferSize, &bpos,
944 remoteIndices_.communicator());
946 char attr = index->local().attribute();
947 MPI_Pack(&attr, 1, MPI_CHAR, buffer, bufferSize, &bpos,
948 remoteIndices_.communicator());
951 MPI_Pack(&indices, 1, MPI_INT, buffer, bufferSize, &bpos,
952 remoteIndices_.communicator());
955 for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
956 if(iterators->second.isNotAtEnd() && iterators->second.isOld()
957 && iterators->second.globalIndexPair().first == index->global()) {
958 int process = iterators->first;
961 assert(pairs <= infoSend_[destination].pairs);
962 MPI_Pack(&process, 1, MPI_INT, buffer, bufferSize, &bpos,
963 remoteIndices_.communicator());
964 char attr2 = iterators->second.remoteIndex().attribute();
966 MPI_Pack(&attr2, 1, MPI_CHAR, buffer, bufferSize, &bpos,
967 remoteIndices_.communicator());
972 Dune::dvverb<<
" (publish="<<published<<
", pairs="<<pairs<<
")"<<std::endl;
973 assert(published <= infoSend_[destination].publish);
977 assert(published == infoSend_[destination].publish);
978 assert(pairs == infoSend_[destination].pairs);
981 Dune::dverb << rank_<<
": Sending message of "<<bpos<<
" bytes to "<<destination<<std::endl;
983 MPI_Issend(buffer, bpos, MPI_PACKED, destination, 345, remoteIndices_.communicator(),&request);
987 inline void IndicesSyncer<T>::insertIntoRemoteIndexList(
int process,
988 const std::pair<GlobalIndex,Attribute>& globalPair,
991 Dune::dverb<<
"Inserting from "<<process<<
" "<<globalPair.first<<
", "<<
992 globalPair.second<<
" "<<attribute<<std::endl;
997 typename IteratorsMap::iterator found = iteratorsMap_.find(process);
999 if( found == iteratorsMap_.end() ) {
1000 Dune::dverb<<
"Discovered new neighbour "<<process<<std::endl;
1001 RemoteIndexList* rlist =
new RemoteIndexList();
1002 remoteIndices_.remoteIndices_.insert(std::make_pair(process,std::make_pair(rlist,rlist)));
1003 Iterators iterators = Iterators(*rlist, globalMap_[process], oldMap_[process]);
1004 found = iteratorsMap_.insert(std::make_pair(process, iterators)).first;
1007 Iterators& iterators = found->second;
1010 while(iterators.isNotAtEnd() && iterators.globalIndexPair() < globalPair) {
1016 if(iterators.isAtEnd() || iterators.globalIndexPair() != globalPair) {
1019 iterators.insert(RemoteIndex(Attribute(attribute)),globalPair);
1024 bool indexIsThere=
false;
1025 for(Iterators tmpIterators = iterators;
1026 !tmpIterators.isAtEnd() && tmpIterators.globalIndexPair() == globalPair;
1029 if(tmpIterators.globalIndexPair().second == attribute) {
1037 iterators.insert(RemoteIndex(Attribute(attribute)),globalPair);
1040 template<
typename T>
1041 template<
typename T1>
1042 void IndicesSyncer<T>::recvAndUnpack(T1& numberer)
1044 typedef typename ParallelIndexSet::const_iterator IndexIterator;
1046 IndexIterator iEnd = indexSet_.end();
1047 IndexIterator index = indexSet_.begin();
1051 assert(checkReset());
1056 MPI_Probe(MPI_ANY_SOURCE, 345, remoteIndices_.communicator(), &status);
1058 int source=status.MPI_SOURCE;
1060 MPI_Get_count(&status, MPI_PACKED, &count);
1062 Dune::dvverb<<rank_<<
": Receiving message from "<< source<<
" with "<<count<<
" bytes"<<std::endl;
1064 if(count>receiveBufferSize_) {
1065 receiveBufferSize_=count;
1066 delete[] receiveBuffer_;
1067 receiveBuffer_ =
new char[receiveBufferSize_];
1070 MPI_Recv(receiveBuffer_, count, MPI_PACKED, source, 345, remoteIndices_.communicator(), &status);
1073 MPI_Unpack(receiveBuffer_, count, &bpos, &publish, 1, MPI_INT, remoteIndices_.communicator());
1080 char sourceAttribute;
1083 MPI_Unpack(receiveBuffer_, count, &bpos, &global, 1, MPITraits<GlobalIndex>::getType(),
1084 remoteIndices_.communicator());
1085 MPI_Unpack(receiveBuffer_, count, &bpos, &sourceAttribute, 1, MPI_CHAR,
1086 remoteIndices_.communicator());
1087 MPI_Unpack(receiveBuffer_, count, &bpos, &pairs, 1, MPI_INT,
1088 remoteIndices_.communicator());
1092 SLList<std::pair<int,Attribute> > sourceAttributeList;
1093 sourceAttributeList.push_back(std::make_pair(source,Attribute(sourceAttribute)));
1095 bool foundSelf =
false;
1097 Attribute myAttribute=Attribute();
1100 for(; pairs>0; --pairs) {
1104 MPI_Unpack(receiveBuffer_, count, &bpos, &process, 1, MPI_INT,
1105 remoteIndices_.communicator());
1107 MPI_Unpack(receiveBuffer_, count, &bpos, &attribute, 1, MPI_CHAR,
1108 remoteIndices_.communicator());
1110 if(process==rank_) {
1114 myAttribute=Attribute(attribute);
1118 IndexIterator pos = std::lower_bound(index, iEnd, IndexPair(global));
1120 if(pos == iEnd || pos->global() != global) {
1122 indexSet_.add(global,
1123 ParallelLocalIndex<Attribute>(numberer(global),
1124 myAttribute,
true));
1125 Dune::dvverb <<
"Adding "<<global<<
" "<<myAttribute<<std::endl;
1130 bool indexIsThere =
false;
1133 for(; pos->global()==global; ++pos)
1134 if(pos->local().attribute() == myAttribute) {
1135 Dune::dvverb<<
"found "<<global<<
" "<<myAttribute<<std::endl;
1136 indexIsThere =
true;
1141 indexSet_.add(global,
1142 ParallelLocalIndex<Attribute>(numberer(global),
1143 myAttribute,
true));
1144 Dune::dvverb <<
"Adding "<<global<<
" "<<myAttribute<<std::endl;
1148 sourceAttributeList.push_back(std::make_pair(process,Attribute(attribute)));
1153 typedef typename SLList<std::pair<int,Attribute> >::const_iterator Iter;
1154 for(Iter i=sourceAttributeList.begin(), end=sourceAttributeList.end();
1156 insertIntoRemoteIndexList(i->first, std::make_pair(global, myAttribute),
1161 resetIteratorsMap();
1164 template<
typename T>
1165 void IndicesSyncer<T>::resetIteratorsMap(){
1168 typedef typename IteratorsMap::iterator Iterator;
1169 typedef typename RemoteIndices::RemoteIndexMap::iterator
1171 typedef typename GlobalIndicesMap::iterator GlobalIterator;
1172 typedef typename BoolMap::iterator BoolIterator;
1174 const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
1175 Iterator iterators = iteratorsMap_.begin();
1176 GlobalIterator global = globalMap_.begin();
1177 BoolIterator added = oldMap_.begin();
1179 for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
1180 remote != remoteEnd; ++remote, ++global, ++added, ++iterators) {
1181 iterators->second.reset(*(remote->second.first), global->second, added->second);
1185 template<
typename T>
1186 bool IndicesSyncer<T>::checkReset(
const Iterators& iterators, RemoteIndexList& rList, GlobalIndexList& gList,
1189 if(std::get<0>(iterators.iterators_) != rList.begin())
1191 if(std::get<1>(iterators.iterators_) != gList.begin())
1193 if(std::get<2>(iterators.iterators_) != bList.begin())
1199 template<
typename T>
1200 bool IndicesSyncer<T>::checkReset(){
1203 typedef typename IteratorsMap::iterator Iterator;
1204 typedef typename RemoteIndices::RemoteIndexMap::iterator
1206 typedef typename GlobalIndicesMap::iterator GlobalIterator;
1207 typedef typename BoolMap::iterator BoolIterator;
1209 const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
1210 Iterator iterators = iteratorsMap_.begin();
1211 GlobalIterator global = globalMap_.begin();
1212 BoolIterator added = oldMap_.begin();
1215 for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
1216 remote != remoteEnd; ++remote, ++global, ++added, ++iterators) {
1217 if(!checkReset(iterators->second, *(remote->second.first), global->second,