00001
00002 #ifndef DUNE_INDICESSYNCER_HH
00003 #define DUNE_INDICESSYNCER_HH
00004
00005 #include"indexset.hh"
00006 #include"remoteindices.hh"
00007 #include<dune/common/stdstreams.hh>
00008 #include<dune/common/tuples.hh>
00009 #include<dune/common/sllist.hh>
00010 #include<cassert>
00011 #include<cmath>
00012 #include<limits>
00013 #include<algorithm>
00014 #include<functional>
00015
00016 #if HAVE_MPI
00017 namespace Dune
00018 {
00035 template<typename T>
00036 class IndicesSyncer
00037 {
00038 public:
00039
00041 typedef T ParallelIndexSet;
00042
00044 typedef typename ParallelIndexSet::IndexPair IndexPair;
00045
00047 typedef typename ParallelIndexSet::GlobalIndex GlobalIndex;
00048
00050 typedef typename ParallelIndexSet::LocalIndex::Attribute Attribute;
00051
00055 typedef RemoteIndices<ParallelIndexSet> RemoteIndices;
00056
00066 IndicesSyncer(ParallelIndexSet& indexSet,
00067 RemoteIndices& remoteIndices);
00068
00077 void sync();
00078
00088 template<typename T1>
00089 void sync(T1& numberer);
00090
00091 private:
00092
00094 ParallelIndexSet& indexSet_;
00095
00097 RemoteIndices& remoteIndices_;
00098
00100 char* buffer_;
00101
00103 int bufferSize_;
00104
00108 struct MessageInformation
00109 {
00110 MessageInformation()
00111 {
00112 publish=0;
00113 pairs=0;
00114 }
00116 int publish;
00121 int pairs;
00122 };
00123
00127 class DefaultNumberer
00128 {
00129 public:
00135 size_t operator()(const GlobalIndex& global)
00136 {
00137 return std::numeric_limits<size_t>::max();
00138 }
00139 };
00140
00142 MPI_Datatype datatype_;
00143
00145 int rank_;
00146
00151 typedef SLList<GlobalIndex, typename RemoteIndices::Allocator> GlobalIndexList;
00152
00154 typedef typename GlobalIndexList::ModifyIterator GlobalIndexModifier;
00155
00159 typedef typename SLList<GlobalIndex, typename RemoteIndices::Allocator>::iterator
00160 GlobalIndexIterator;
00161
00163 typedef std::map<int, GlobalIndexList> GlobalIndicesMap;
00164
00173 GlobalIndicesMap globalMap_;
00174
00178 typedef SLList<bool, typename RemoteIndices::Allocator> BoolList;
00179
00183 typedef typename BoolList::iterator BoolIterator;
00184
00186 typedef typename BoolList::ModifyIterator BoolListModifier;
00187
00189 typedef std::map<int,BoolList> BoolMap;
00190
00195 BoolMap oldMap_;
00196
00198 std::map<int,MessageInformation> infoSend_;
00199
00201 typedef typename RemoteIndices::RemoteIndexList RemoteIndexList;
00202
00204 typedef typename RemoteIndexList::ModifyIterator RemoteIndexModifier;
00205
00207 typedef RemoteIndex<GlobalIndex,Attribute> RemoteIndex;
00208
00210 typedef typename RemoteIndexList::iterator RemoteIndexIterator;
00211
00213 typedef typename RemoteIndexList::const_iterator ConstRemoteIndexIterator;
00214
00216 typedef Tuple<RemoteIndexModifier,GlobalIndexModifier,BoolListModifier,
00217 const ConstRemoteIndexIterator> IteratorTuple;
00218
00226 class Iterators
00227 {
00228 friend class IndicesSyncer<T>;
00229 public:
00239 Iterators(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
00240 BoolList& booleans);
00241
00245 Iterators();
00246
00250 Iterators& operator++();
00251
00257 void insert(const RemoteIndex& index, const GlobalIndex& global);
00258
00263 RemoteIndex& remoteIndex() const;
00264
00269 GlobalIndex& globalIndex() const;
00270
00276 bool isOld() const;
00277
00287 void reset(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
00288 BoolList& booleans);
00289
00295 bool isNotAtEnd() const;
00296
00302 bool isAtEnd() const;
00303
00304 private:
00314 IteratorTuple iterators_;
00315 };
00316
00318 typedef std::map<int,Iterators> IteratorsMap;
00319
00331 IteratorsMap iteratorsMap_;
00332
00334 void calculateMessageSizes();
00335
00340 void packAndSend(int destination);
00341
00347 template<typename T1>
00348 void recvAndUnpack(int source, T1& numberer);
00349
00353 void registerMessageDatatype();
00354
00358 void insertIntoRemoteIndexList(int process, const GlobalIndex& global, char attribute);
00359
00363 void resetIteratorsMap();
00364
00369 bool checkReset();
00370
00379 bool checkReset(const Iterators& iterators, RemoteIndexList& rlist, GlobalIndexList& gList,
00380 BoolList& bList);
00381 };
00382
00391 template<typename T, typename A>
00392 inline void repairLocalIndexPointers(std::map<int,SLList<typename T::GlobalIndex,A> >& globalMap,
00393 RemoteIndices<T>& remoteIndices,
00394 const T& indexSet)
00395 {
00396 typedef typename RemoteIndices<T>::RemoteIndexMap::iterator RemoteIterator;
00397 typedef typename RemoteIndices<T>::RemoteIndexList::iterator RemoteIndexIterator;
00398 typedef typename T::GlobalIndex GlobalIndex;
00399 typedef typename SLList<GlobalIndex,A>::iterator GlobalIndexIterator;
00400
00401 assert(globalMap.size()==remoteIndices.remoteIndices_.size());
00402
00403 typename std::map<int,SLList<GlobalIndex,A> >::iterator global = globalMap.begin();
00404 RemoteIterator end = remoteIndices.remoteIndices_.end();
00405
00406 for(RemoteIterator remote = remoteIndices.remoteIndices_.begin(); remote != end; ++remote, ++global){
00407 typedef typename T::const_iterator IndexIterator;
00408
00409 assert(remote->first==global->first);
00410 assert(remote->second.first->size() == global->second.size());
00411
00412 RemoteIndexIterator riEnd = remote->second.first->end();
00413 RemoteIndexIterator rIndex = remote->second.first->begin();
00414 GlobalIndexIterator gIndex = global->second.begin();
00415 GlobalIndexIterator gEnd = global->second.end();
00416 IndexIterator index = indexSet.begin();
00417
00418 assert(rIndex==riEnd || gIndex != global->second.end());
00419 while(rIndex != riEnd){
00420
00421 assert(gIndex != global->second.end());
00422
00423 while(index->global() < *gIndex)
00424 ++index;
00425
00426 assert(index != indexSet.end() && index->global() == *gIndex);
00427
00428 rIndex->localIndex_ = &(*index);
00429
00430 ++rIndex;
00431 ++gIndex;
00432 }
00433 }
00434
00435 }
00436
00437 template<typename T>
00438 IndicesSyncer<T>::IndicesSyncer(ParallelIndexSet& indexSet,
00439 RemoteIndices& remoteIndices)
00440 : indexSet_(indexSet), remoteIndices_(remoteIndices)
00441 {
00442
00443 assert(remoteIndices.source_ == remoteIndices.target_);
00444 assert(remoteIndices.source_ == &indexSet);
00445 MPI_Comm_rank(remoteIndices_.communicator(), &rank_);
00446 }
00447
00448 template<typename T>
00449 IndicesSyncer<T>::Iterators::Iterators(RemoteIndexList& remoteIndices,
00450 GlobalIndexList& globalIndices,
00451 BoolList& booleans)
00452 : iterators_(remoteIndices.beginModify(), globalIndices.beginModify(),
00453 booleans.beginModify(), remoteIndices.end())
00454 { }
00455
00456 template<typename T>
00457 IndicesSyncer<T>::Iterators::Iterators()
00458 : iterators_()
00459 {}
00460
00461 template<typename T>
00462 inline typename IndicesSyncer<T>::Iterators& IndicesSyncer<T>::Iterators::operator++()
00463 {
00464 ++(Element<0>::get(iterators_));
00465 ++(Element<1>::get(iterators_));
00466 ++(Element<2>::get(iterators_));
00467 return *this;
00468 }
00469
00470 template<typename T>
00471 inline void IndicesSyncer<T>::Iterators::insert(const RemoteIndex& index,
00472 const GlobalIndex& global)
00473 {
00474 Element<0>::get(iterators_).insert(index);
00475 Element<1>::get(iterators_).insert(global);
00476 Element<2>::get(iterators_).insert(false);
00477 }
00478
00479 template<typename T>
00480 inline typename IndicesSyncer<T>::RemoteIndex&
00481 IndicesSyncer<T>::Iterators::remoteIndex() const
00482 {
00483 return *(Element<0>::get(iterators_));
00484 }
00485
00486 template<typename T>
00487 inline typename IndicesSyncer<T>::GlobalIndex& IndicesSyncer<T>::Iterators::globalIndex() const
00488 {
00489 return *(Element<1>::get(iterators_));
00490 }
00491
00492 template<typename T>
00493 inline bool IndicesSyncer<T>::Iterators::isOld() const
00494 {
00495 return *(Element<2>::get(iterators_));
00496 }
00497
00498 template<typename T>
00499 inline void IndicesSyncer<T>::Iterators::reset(RemoteIndexList& remoteIndices,
00500 GlobalIndexList& globalIndices,
00501 BoolList& booleans)
00502 {
00503 Element<0>::get(iterators_) = remoteIndices.beginModify();
00504 Element<1>::get(iterators_) = globalIndices.beginModify();
00505 Element<2>::get(iterators_) = booleans.beginModify();
00506 }
00507
00508 template<typename T>
00509 inline bool IndicesSyncer<T>::Iterators::isNotAtEnd() const
00510 {
00511 return Element<0>::get(iterators_)!=Element<3>::get(iterators_);
00512 }
00513
00514 template<typename T>
00515 inline bool IndicesSyncer<T>::Iterators::isAtEnd() const
00516 {
00517 return Element<0>::get(iterators_)==Element<3>::get(iterators_);
00518 }
00519
00520 template<typename T>
00521 void IndicesSyncer<T>::registerMessageDatatype()
00522 {
00523 MPI_Datatype type[2] = {MPI_INT, MPI_INT};
00524 int blocklength[2] = {1,1};
00525 MPI_Aint displacement[2];
00526 MPI_Aint base;
00527
00528
00529 MessageInformation message;
00530
00531 MPI_Address( &(message.publish), displacement);
00532 MPI_Address( &(message.pairs), displacement+1);
00533
00534
00535 MPI_Address(&message, &base);
00536 displacement[0] -= base;
00537 displacement[1] -= base;
00538
00539 MPI_Type_struct( 2, blocklength, displacement, type, &datatype_);
00540 MPI_Type_commit(&datatype_);
00541 }
00542
00543 template<typename T>
00544 void IndicesSyncer<T>::calculateMessageSizes()
00545 {
00546 typedef typename ParallelIndexSet::const_iterator IndexIterator;
00547 typedef CollectiveIterator<T> CollectiveIterator;
00548
00549 IndexIterator iEnd = indexSet_.end();
00550 CollectiveIterator collIter = remoteIndices_.template iterator<true>();
00551
00552 for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index){
00553 collIter.advance(index->global());
00554 if(collIter.empty())
00555 break;
00556 int knownRemote=0;
00557
00558 typedef typename CollectiveIterator::iterator ValidIterator;
00559 ValidIterator end = collIter.end();
00560
00561
00562 for(ValidIterator valid = collIter.begin(); valid != end; ++valid)
00563 ++knownRemote;
00564
00565 if(knownRemote>0){
00566 Dune::dverb<<rank_<<": publishing "<<knownRemote<<" for index "<<index->global()<< " for processes ";
00567
00568
00569 for(ValidIterator valid = collIter.begin(); valid != end; ++valid){
00570 ++(infoSend_[valid.process()].publish);
00571 (infoSend_[valid.process()].pairs) += knownRemote;
00572 Dune::dverb<<valid.process()<<" ";
00573 Dune::dverb<<"(publish="<<infoSend_[valid.process()].publish<<", pairs="<<infoSend_[valid.process()].pairs
00574 <<") ";
00575 }
00576 Dune::dverb<<std::endl;
00577 }
00578 }
00579
00580 typedef typename std::map<int,MessageInformation>::const_iterator
00581 MessageIterator;
00582
00583 const MessageIterator end = infoSend_.end();
00584
00585 registerMessageDatatype();
00586
00587 MessageInformation maxSize, dummy;
00588
00589 MessageIterator messageIter= infoSend_.begin();
00590
00591 typedef typename RemoteIndices::RemoteIndexMap::const_iterator RemoteIterator;
00592 const RemoteIterator rend = remoteIndices_.end();
00593
00594 for(RemoteIterator remote = remoteIndices_.begin(); remote != rend; ++remote){
00595 MessageInformation* message;
00596 MessageInformation recv;
00597 MPI_Status status;
00598
00599 if(messageIter != end && messageIter->first==remote->first){
00600
00601 message = const_cast<MessageInformation*>(&(messageIter->second));
00602 ++messageIter;
00603 }else
00604
00605 message = &dummy;
00606
00607 if(remote->first < rank_){
00608 MPI_Send(message, 1, datatype_, remote->first, 122, remoteIndices_.communicator());
00609 MPI_Recv(&recv, 1, datatype_, remote->first, 122, remoteIndices_.communicator(), &status);
00610 }else{
00611 MPI_Recv(&recv, 1, datatype_, remote->first, 122, remoteIndices_.communicator(), &status);
00612 MPI_Send(message, 1, datatype_, remote->first, 122, remoteIndices_.communicator());
00613 }
00614
00615
00616
00617
00618
00619 maxSize.publish = std::max(maxSize.publish, message->publish);
00620 maxSize.pairs = std::max(maxSize.pairs, message->pairs);
00621 maxSize.publish = std::max(maxSize.publish, recv.publish);
00622 maxSize.pairs = std::max(maxSize.pairs, recv.pairs);
00623 }
00624
00625 MPI_Type_free(&datatype_);
00626
00627
00628 bufferSize_=0;
00629 int tsize;
00630
00631 MPI_Pack_size(1, MPI_INT,remoteIndices_.communicator(), &tsize);
00632 bufferSize_ += tsize;
00633
00634 for(int i=0; i < maxSize.publish; ++i){
00635
00636 MPI_Pack_size(1, MPITraits<GlobalIndex>::getType(), remoteIndices_.communicator(), &tsize);
00637 bufferSize_ += tsize;
00638
00639 MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
00640 bufferSize_ += tsize;
00641
00642 MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
00643 bufferSize_ += tsize;
00644 }
00645 for(int i=0; i < maxSize.pairs; ++i){
00646
00647 MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
00648 bufferSize_ += tsize;
00649
00650 MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
00651 bufferSize_ += tsize;
00652 }
00653
00654
00655 buffer_ = new char[bufferSize_];
00656
00657 Dune::dverb<<rank_<<": Buffer size is "<< bufferSize_<<" for publish="<<maxSize.publish<<" pairs="<<maxSize.pairs<<std::endl;
00658 }
00659
00660 template<typename T>
00661 inline void IndicesSyncer<T>::sync()
00662 {
00663 DefaultNumberer numberer;
00664 sync(numberer);
00665 }
00666
00667 template<typename T>
00668 template<typename T1>
00669 void IndicesSyncer<T>::sync(T1& numberer)
00670 {
00671
00672
00673
00674
00675
00676
00677 typedef typename RemoteIndices::RemoteIndexMap::const_iterator
00678 RemoteIterator;
00679
00680 const RemoteIterator end = remoteIndices_.end();
00681
00682
00683
00684 int noOldNeighbours = remoteIndices_.neighbours();
00685 int* oldNeighbours = new int[noOldNeighbours];
00686 int i=0;
00687
00688 for(RemoteIterator remote = remoteIndices_.begin(); remote != end; ++remote, ++i){
00689 typedef typename RemoteIndices::RemoteIndexList::const_iterator
00690 RemoteIndexIterator;
00691
00692 oldNeighbours[i]=remote->first;
00693
00694
00695 assert(remote->second.first==remote->second.second);
00696
00697 RemoteIndexList& rList = *(remote->second.first);
00698
00699
00700 GlobalIndexList& global = globalMap_[remote->first];
00701 BoolList& added = oldMap_[remote->first];
00702 RemoteIndexIterator riEnd = rList.end();
00703
00704 for(RemoteIndexIterator index = rList.begin();
00705 index != riEnd; ++index){
00706 global.push_back(index->localIndexPair().global());
00707 added.push_back(true);
00708 }
00709
00710 Iterators iterators(rList, global, added);
00711 iteratorsMap_.insert(std::make_pair(remote->first, iterators));
00712 assert(checkReset(iteratorsMap_[remote->first], rList,global,added));
00713 }
00714
00715
00716 const RemoteIterator rend = remoteIndices_.end();
00717
00718 calculateMessageSizes();
00719
00720 indexSet_.beginResize();
00721
00722 Dune::dverb<<rank_<<": Neighbours: ";
00723
00724 for(i = 0; i<noOldNeighbours; ++i)
00725 Dune::dverb<<oldNeighbours[i]<<" ";
00726
00727 Dune::dverb<<std::endl;
00728
00729 for(i = 0; i<noOldNeighbours; ++i){
00730 if(oldNeighbours[i] < rank_){
00731 packAndSend(oldNeighbours[i]);
00732 recvAndUnpack(oldNeighbours[i], numberer);
00733 }else{
00734 recvAndUnpack(oldNeighbours[i], numberer);
00735 packAndSend(oldNeighbours[i]);
00736 }
00737 }
00738
00739 delete[] buffer_;
00740
00741
00742 iteratorsMap_.clear();
00743
00744 indexSet_.endResize();
00745
00746 delete[] oldNeighbours;
00747
00748 repairLocalIndexPointers<T,typename RemoteIndices::Allocator>(globalMap_, remoteIndices_, indexSet_);
00749
00750 oldMap_.clear();
00751 globalMap_.clear();
00752
00753
00754 remoteIndices_.sourceSeqNo_ = remoteIndices_.destSeqNo_ = indexSet_.seqNo();
00755 }
00756
00757
00758 template<typename T>
00759 void IndicesSyncer<T>::packAndSend(int destination)
00760 {
00761 typedef typename ParallelIndexSet::const_iterator IndexIterator;
00762 typedef typename RemoteIndexList::const_iterator RemoteIndexIterator;
00763 typedef typename GlobalIndexList::const_iterator GlobalIterator;
00764 typedef typename BoolList::const_iterator BoolIterator;
00765
00766 IndexIterator iEnd = indexSet_.end();
00767 int bpos = 0;
00768 int published = 0;
00769 int pairs = 0;
00770
00771 assert(checkReset());
00772
00773
00774 MPI_Pack(&(infoSend_[destination].publish), 1, MPI_INT, buffer_, bufferSize_, &bpos,
00775 remoteIndices_.communicator());
00776
00777 for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index){
00778
00779 typedef typename IteratorsMap::iterator Iterator;
00780 Iterator iteratorsEnd = iteratorsMap_.end();
00781
00782
00783 for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
00784 while(iterators->second.isNotAtEnd() && iterators->second.globalIndex() < index->global())
00785 ++(iterators->second);
00786
00787
00788
00789
00790 int indices = 0;
00791 bool knownRemote = false;
00792
00793 for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
00794 if(iterators->second.isNotAtEnd() && iterators->second.isOld()
00795 && iterators->second.globalIndex() == index->global()){
00796 indices++;
00797 if(destination == iterators->first)
00798 knownRemote = true;
00799 }
00800
00801 if(!knownRemote || indices==0)
00802
00803 continue;
00804
00805 Dune::dverb<<rank_<<": sending "<<indices<<" for index "<<index->global()<<" to "<<destination<<std::endl;
00806
00807 pairs+=indices;
00808 assert(pairs <= infoSend_[destination].pairs);
00809
00810
00811 MPI_Pack(const_cast<GlobalIndex*>(&(index->global())), 1, MPITraits<GlobalIndex>::getType(), buffer_, bufferSize_, &bpos,
00812 remoteIndices_.communicator());
00813
00814 char attr = index->local().attribute();
00815 MPI_Pack(&attr, 1, MPI_CHAR, buffer_, bufferSize_, &bpos,
00816 remoteIndices_.communicator());
00817
00818
00819 MPI_Pack(&indices, 1, MPI_INT, buffer_, bufferSize_, &bpos,
00820 remoteIndices_.communicator());
00821
00822
00823 for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
00824 if(iterators->second.isNotAtEnd() && iterators->second.isOld()
00825 && iterators->second.globalIndex() == index->global()){
00826 int process = iterators->first;
00827 MPI_Pack(&process, 1, MPI_INT, buffer_, bufferSize_, &bpos,
00828 remoteIndices_.communicator());
00829 char attr = iterators->second.remoteIndex().attribute();
00830
00831 MPI_Pack(&attr, 1, MPI_CHAR, buffer_, bufferSize_, &bpos,
00832 remoteIndices_.communicator());
00833 }
00834 ++published;
00835 assert(published <= infoSend_[destination].publish);
00836 }
00837
00838
00839 assert(published == infoSend_[destination].publish);
00840
00841 resetIteratorsMap();
00842
00843 Dune::dverb << rank_<<": Sending message of "<<bpos<<" bytes to "<<destination<<std::endl;
00844
00845 MPI_Send(buffer_, bpos, MPI_PACKED, destination, 111, remoteIndices_.communicator());
00846 }
00847
00848 template<typename T>
00849 inline void IndicesSyncer<T>::insertIntoRemoteIndexList(int process, const GlobalIndex& global,
00850 char attribute)
00851 {
00852 Dune::dvverb<<"Inserting from "<<process<<" "<<global<<" "<<attribute<<std::endl;
00853
00854
00855 typename IteratorsMap::iterator found = iteratorsMap_.find(process);
00856
00857 if( found == iteratorsMap_.end() ){
00858 Dune::dverb<<"Discovered new neighbour "<<process<<std::endl;
00859 RemoteIndexList* rlist = new RemoteIndexList();
00860 remoteIndices_.remoteIndices_.insert(std::make_pair(process,std::make_pair(rlist,rlist)));
00861 Iterators iterators = Iterators(*rlist, globalMap_[process], oldMap_[process]);
00862 found = iteratorsMap_.insert(std::make_pair(process, iterators)).first;
00863 }
00864
00865 Iterators& iterators = found->second;
00866
00867
00868 while(iterators.isNotAtEnd() && iterators.globalIndex() < global){
00869
00870 ++iterators;
00871
00872 }
00873
00874 if(iterators.isAtEnd() || iterators.globalIndex() != global){
00875
00876
00877 iterators.insert(RemoteIndex(Attribute(attribute)),global);
00878 }else if(iterators.isNotAtEnd()){
00879
00880 assert(iterators.remoteIndex().attribute() == attribute);
00881 }
00882 }
00883
00884 template<typename T>
00885 template<typename T1>
00886 void IndicesSyncer<T>::recvAndUnpack(int source, T1& numberer)
00887 {
00888 typedef typename ParallelIndexSet::const_iterator IndexIterator;
00889 typedef typename RemoteIndexList::iterator RemoteIndexIterator;
00890 typedef typename GlobalIndexList::iterator GlobalIndexIterator;
00891
00892 IndexIterator iEnd = indexSet_.end();
00893 IndexIterator index = indexSet_.begin();
00894 int bpos = 0;
00895 int publish;
00896
00897 assert(checkReset());
00898
00899 Dune::dvverb<<rank_<<": Waiting for message from "<< source<<std::endl;
00900
00901
00902 MPI_Status status;
00903
00904
00905 MPI_Probe(source, 111, remoteIndices_.communicator(), &status);
00906
00907 int count;
00908 MPI_Get_count(&status, MPI_PACKED, &count);
00909
00910 if(count>bufferSize_){
00911 delete[] buffer_;
00912 buffer_ = new char[count];
00913 bufferSize_ = count;
00914 }
00915
00916 MPI_Recv(buffer_, count, MPI_PACKED, source, 111, remoteIndices_.communicator(), &status);
00917
00918
00919 MPI_Unpack(buffer_, count, &bpos, &publish, 1, MPI_INT, remoteIndices_.communicator());
00920
00921
00922 while(publish>0){
00923
00924
00925 GlobalIndex global;
00926 char sourceAttribute;
00927 int pairs;
00928
00929 MPI_Unpack(buffer_, count, &bpos, &global, 1, MPITraits<GlobalIndex>::getType(),
00930 remoteIndices_.communicator());
00931 MPI_Unpack(buffer_, count, &bpos, &sourceAttribute, 1, MPI_CHAR,
00932 remoteIndices_.communicator());
00933 MPI_Unpack(buffer_, count, &bpos, &pairs, 1, MPI_INT,
00934 remoteIndices_.communicator());
00935
00936
00937
00938 insertIntoRemoteIndexList(source, global, sourceAttribute);
00939
00940
00941 while(pairs>0){
00942
00943 int process;
00944 char attribute;
00945 MPI_Unpack(buffer_, count, &bpos, &process, 1, MPI_INT,
00946 remoteIndices_.communicator());
00947
00948 MPI_Unpack(buffer_, count, &bpos, &attribute, 1, MPI_CHAR,
00949 remoteIndices_.communicator());
00950
00951 if(process==rank_){
00952
00953
00954 IndexIterator pos = std::lower_bound(index, iEnd, IndexPair(global));
00955
00956 if(pos == iEnd || pos->global() != global){
00957
00958 indexSet_.add(global,ParallelLocalIndex<Attribute>(numberer(global),
00959 Attribute(attribute), true));
00960 }else{
00961
00962 assert(attribute==pos->local().attribute());
00963 }
00964 index=pos;
00965 }else{
00966 insertIntoRemoteIndexList(process, global, attribute);
00967 }
00968
00969 --pairs;
00970 }
00971 --publish;
00972 }
00973
00974 resetIteratorsMap();
00975 }
00976
00977 template<typename T>
00978 void IndicesSyncer<T>::resetIteratorsMap(){
00979
00980
00981 typedef typename IteratorsMap::iterator Iterator;
00982 typedef typename RemoteIndices::RemoteIndexMap::iterator
00983 RemoteIterator;
00984 typedef typename GlobalIndicesMap::iterator GlobalIterator;
00985 typedef typename BoolMap::iterator BoolIterator;
00986
00987 const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
00988 Iterator iterators = iteratorsMap_.begin();
00989 GlobalIterator global = globalMap_.begin();
00990 BoolIterator added = oldMap_.begin();
00991
00992 for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
00993 remote != remoteEnd; ++remote, ++global, ++added, ++iterators){
00994 iterators->second.reset(*(remote->second.first), global->second, added->second);
00995 }
00996 }
00997
00998 template<typename T>
00999 bool IndicesSyncer<T>::checkReset(const Iterators& iterators, RemoteIndexList& rList, GlobalIndexList& gList,
01000 BoolList& bList){
01001
01002 if(Element<0>::get(iterators.iterators_)!=rList.begin())
01003 return false;
01004 if(Element<1>::get(iterators.iterators_)!=gList.begin())
01005 return false;
01006 if(Element<2>::get(iterators.iterators_)!=bList.begin())
01007 return false;
01008 return true;
01009 }
01010
01011
01012 template<typename T>
01013 bool IndicesSyncer<T>::checkReset(){
01014
01015
01016 typedef typename IteratorsMap::iterator Iterator;
01017 typedef typename RemoteIndices::RemoteIndexMap::iterator
01018 RemoteIterator;
01019 typedef typename GlobalIndicesMap::iterator GlobalIterator;
01020 typedef typename BoolMap::iterator BoolIterator;
01021
01022 const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
01023 Iterator iterators = iteratorsMap_.begin();
01024 GlobalIterator global = globalMap_.begin();
01025 BoolIterator added = oldMap_.begin();
01026 bool ret = true;
01027
01028 for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
01029 remote != remoteEnd; ++remote, ++global, ++added, ++iterators){
01030 if(!checkReset(iterators->second, *(remote->second.first), global->second,
01031 added->second))
01032 ret=false;
01033 }
01034 return ret;
01035 }
01036 }
01037
01038 #endif
01039 #endif