indicessyncer.hh

Go to the documentation of this file.
00001 // $Id: indicessyncer.hh 960 2008-11-12 12:58:35Z mblatt $
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 Dune::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** sendBuffers_;
00101 
00103     char* receiveBuffer_;
00104     
00106     std::size_t* sendBufferSizes_;
00107     
00109     int receiveBufferSize_; // int because of MPI
00110     
00114     struct MessageInformation
00115     {
00116       MessageInformation()
00117         : publish(), pairs()
00118       {}
00120       int publish;
00125       int pairs;
00126     };
00127 
00131     class DefaultNumberer
00132     {
00133     public:
00139       std::size_t operator()(const GlobalIndex& global)
00140       {
00141         return std::numeric_limits<size_t>::max();
00142       }
00143     };
00144     
00146     MPI_Datatype datatype_;
00147 
00149     int rank_;
00150     
00155     typedef SLList<GlobalIndex, typename RemoteIndices::Allocator> GlobalIndexList;
00156     
00158     typedef typename GlobalIndexList::ModifyIterator GlobalIndexModifier;
00159     
00163     typedef typename SLList<GlobalIndex, typename RemoteIndices::Allocator>::iterator 
00164     GlobalIndexIterator;
00165 
00167     typedef std::map<int, GlobalIndexList> GlobalIndicesMap;
00168     
00177     GlobalIndicesMap globalMap_;
00178 
00182     typedef SLList<bool, typename RemoteIndices::Allocator> BoolList;
00183 
00187     typedef typename BoolList::iterator BoolIterator;
00188 
00190     typedef typename BoolList::ModifyIterator BoolListModifier;
00191     
00193     typedef std::map<int,BoolList> BoolMap;
00194 
00199     BoolMap oldMap_;
00200     
00202     std::map<int,MessageInformation> infoSend_;
00203     
00205     typedef typename RemoteIndices::RemoteIndexList RemoteIndexList;
00206 
00208     typedef typename RemoteIndexList::ModifyIterator RemoteIndexModifier;
00209     
00211     typedef Dune::RemoteIndex<GlobalIndex,Attribute> RemoteIndex;
00212     
00214     typedef typename RemoteIndexList::iterator RemoteIndexIterator;
00215 
00217     typedef typename RemoteIndexList::const_iterator ConstRemoteIndexIterator;
00218 
00220     typedef Dune::tuple<RemoteIndexModifier,GlobalIndexModifier,BoolListModifier,
00221                   const ConstRemoteIndexIterator> IteratorTuple;
00222 
00230     class Iterators
00231     {
00232       friend class IndicesSyncer<T>;
00233     public:
00243       Iterators(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
00244                 BoolList& booleans);
00245       
00249       Iterators();
00250       
00254       Iterators& operator++();
00255       
00261       void insert(const RemoteIndex& index, const GlobalIndex& global);
00262       
00267       RemoteIndex& remoteIndex() const;
00268       
00273       GlobalIndex& globalIndex() const;
00274       
00280       bool isOld() const;
00281       
00291       void reset(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
00292                  BoolList& booleans);
00293 
00299       bool isNotAtEnd() const;
00300       
00306       bool isAtEnd() const;
00307       
00308     private:
00318       IteratorTuple iterators_;
00319     };
00320     
00322     typedef std::map<int,Iterators> IteratorsMap;
00323     
00335     IteratorsMap iteratorsMap_;
00336         
00338     void calculateMessageSizes();
00339 
00347     void packAndSend(int destination, char* buffer, std::size_t bufferSize, MPI_Request& req);
00348     
00353     template<typename T1>
00354     void recvAndUnpack(T1& numberer);
00355 
00359     void registerMessageDatatype();
00360     
00364     void insertIntoRemoteIndexList(int process, const GlobalIndex& global, char attribute);
00365 
00369     void resetIteratorsMap();
00370     
00375     bool checkReset();
00376     
00385     bool checkReset(const Iterators& iterators, RemoteIndexList& rlist, GlobalIndexList& gList,
00386                     BoolList& bList);
00387   };
00388   
00397   template<typename T, typename A>
00398   inline void repairLocalIndexPointers(std::map<int,SLList<typename T::GlobalIndex,A> >& globalMap,
00399                                        RemoteIndices<T>& remoteIndices,
00400                                        const T& indexSet)
00401   {
00402     typedef typename RemoteIndices<T>::RemoteIndexMap::iterator RemoteIterator;
00403     typedef typename RemoteIndices<T>::RemoteIndexList::iterator RemoteIndexIterator;
00404     typedef typename T::GlobalIndex GlobalIndex;
00405     typedef typename SLList<GlobalIndex,A>::iterator GlobalIndexIterator;
00406 
00407     assert(globalMap.size()==remoteIndices.remoteIndices_.size());
00408     // Repair pointers to index set in remote indices.
00409     typename std::map<int,SLList<GlobalIndex,A> >::iterator global = globalMap.begin();
00410     RemoteIterator end = remoteIndices.remoteIndices_.end();
00411     
00412     for(RemoteIterator remote = remoteIndices.remoteIndices_.begin(); remote != end; ++remote, ++global){
00413       typedef typename T::const_iterator IndexIterator;
00414 
00415       assert(remote->first==global->first);
00416       assert(remote->second.first->size() == global->second.size());
00417       
00418       RemoteIndexIterator riEnd  = remote->second.first->end();
00419       RemoteIndexIterator rIndex = remote->second.first->begin();
00420       GlobalIndexIterator gIndex = global->second.begin();
00421       GlobalIndexIterator gEnd   = global->second.end();
00422       IndexIterator       index  = indexSet.begin();
00423 
00424       assert(rIndex==riEnd || gIndex != global->second.end());
00425       while(rIndex != riEnd){
00426         // Search for the index in the set.
00427         assert(gIndex != global->second.end());
00428         
00429         while(index->global() < *gIndex)
00430           ++index;
00431         
00432         assert(index != indexSet.end() && index->global() == *gIndex);
00433         
00434         rIndex->localIndex_ = &(*index);
00435         
00436         ++rIndex;
00437         ++gIndex;
00438       }
00439     }
00440     
00441   }
00442   
00443   template<typename T>
00444   IndicesSyncer<T>::IndicesSyncer(ParallelIndexSet& indexSet, 
00445                                       RemoteIndices& remoteIndices)
00446     : indexSet_(indexSet), remoteIndices_(remoteIndices)
00447   {
00448     // index sets must match.
00449     assert(remoteIndices.source_ == remoteIndices.target_);
00450     assert(remoteIndices.source_ == &indexSet);
00451     MPI_Comm_rank(remoteIndices_.communicator(), &rank_);
00452   }
00453     
00454   template<typename T>
00455   IndicesSyncer<T>::Iterators::Iterators(RemoteIndexList& remoteIndices, 
00456                                                GlobalIndexList& globalIndices,
00457                                                BoolList& booleans)
00458     : iterators_(remoteIndices.beginModify(), globalIndices.beginModify(),
00459                     booleans.beginModify(), remoteIndices.end())
00460   { }
00461 
00462   template<typename T>
00463   IndicesSyncer<T>::Iterators::Iterators()
00464     : iterators_()
00465   {}
00466   
00467   template<typename T>
00468   inline typename IndicesSyncer<T>::Iterators& IndicesSyncer<T>::Iterators::operator++()
00469   {
00470       ++(get<0>(iterators_));
00471       ++(get<1>(iterators_));
00472       ++(get<2>(iterators_));
00473       return *this;
00474   }
00475 
00476   template<typename T>
00477   inline void IndicesSyncer<T>::Iterators::insert(const RemoteIndex& index, 
00478                                                  const GlobalIndex& global)
00479   {
00480     get<0>(iterators_).insert(index);
00481     get<1>(iterators_).insert(global);
00482     get<2>(iterators_).insert(false);
00483   }
00484   
00485   template<typename T>
00486   inline typename IndicesSyncer<T>::RemoteIndex& 
00487   IndicesSyncer<T>::Iterators::remoteIndex() const
00488   {
00489     return *(get<0>(iterators_));
00490   }
00491   
00492   template<typename T>
00493   inline typename IndicesSyncer<T>::GlobalIndex&  IndicesSyncer<T>::Iterators::globalIndex() const
00494   {
00495     return *(get<1>(iterators_));
00496   }
00497   
00498   template<typename T>
00499   inline bool IndicesSyncer<T>::Iterators::isOld() const
00500   {
00501     return *(get<2>(iterators_));
00502   }
00503   
00504   template<typename T>
00505   inline void IndicesSyncer<T>::Iterators::reset(RemoteIndexList& remoteIndices, 
00506                                                 GlobalIndexList& globalIndices,
00507                                                 BoolList& booleans)
00508   {
00509     get<0>(iterators_) = remoteIndices.beginModify();
00510     get<1>(iterators_) = globalIndices.beginModify();
00511     get<2>(iterators_) = booleans.beginModify();
00512   }
00513   
00514   template<typename T>
00515   inline bool IndicesSyncer<T>::Iterators::isNotAtEnd() const
00516   {
00517     return get<0>(iterators_)!=get<3>(iterators_);
00518   }
00519   
00520   template<typename T>
00521   inline bool IndicesSyncer<T>::Iterators::isAtEnd() const
00522   {
00523     return get<0>(iterators_)==get<3>(iterators_);
00524   }
00525 
00526   template<typename T>
00527   void IndicesSyncer<T>::registerMessageDatatype()
00528   {
00529     MPI_Datatype type[2] = {MPI_INT, MPI_INT};
00530     int blocklength[2] = {1,1};
00531     MPI_Aint displacement[2];
00532     MPI_Aint base;
00533     
00534     // Compute displacement
00535     MessageInformation message;
00536     
00537     MPI_Address( &(message.publish), displacement);
00538     MPI_Address( &(message.pairs), displacement+1);
00539     
00540     // Make the displacement relative
00541     MPI_Address(&message, &base);
00542     displacement[0] -= base;
00543     displacement[1] -= base;
00544     
00545     MPI_Type_struct( 2, blocklength, displacement, type, &datatype_); 
00546     MPI_Type_commit(&datatype_);
00547   }
00548   
00549   template<typename T>
00550   void IndicesSyncer<T>::calculateMessageSizes()
00551   {    
00552     typedef typename ParallelIndexSet::const_iterator IndexIterator;
00553     typedef CollectiveIterator<T> CollectiveIterator;
00554     
00555     IndexIterator iEnd = indexSet_.end();
00556     CollectiveIterator collIter = remoteIndices_.template iterator<true>();
00557     
00558     for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index){
00559       collIter.advance(index->global());
00560       if(collIter.empty())
00561         break;
00562       int knownRemote=0;
00563 
00564       typedef typename CollectiveIterator::iterator ValidIterator;
00565       ValidIterator end = collIter.end();
00566 
00567       // Count the remote indices we know.
00568       for(ValidIterator valid = collIter.begin(); valid != end; ++valid)
00569         ++knownRemote;
00570       
00571       if(knownRemote>0){
00572         Dune::dverb<<rank_<<": publishing "<<knownRemote<<" for index "<<index->global()<< " for processes ";
00573         
00574         // Update MessageInformation
00575         for(ValidIterator valid = collIter.begin(); valid != end; ++valid){
00576           ++(infoSend_[valid.process()].publish);
00577           (infoSend_[valid.process()].pairs) += knownRemote;
00578           Dune::dverb<<valid.process()<<" ";
00579           Dune::dverb<<"(publish="<<infoSend_[valid.process()].publish<<", pairs="<<infoSend_[valid.process()].pairs
00580                <<") ";
00581         }
00582         Dune::dverb<<std::endl;
00583       }
00584     }
00585     
00586     typedef typename std::map<int,MessageInformation>::const_iterator 
00587       MessageIterator;
00588     
00589     const MessageIterator end = infoSend_.end();
00590 
00591     // registerMessageDatatype();
00592     
00593     // Now determine the buffersizes needed for each neighbour using MPI_Pack_size
00594     MessageInformation dummy;
00595     
00596     MessageIterator messageIter= infoSend_.begin();
00597 
00598     typedef typename RemoteIndices::RemoteIndexMap::const_iterator RemoteIterator;
00599     const RemoteIterator rend = remoteIndices_.end();
00600     int neighbour=0;
00601     
00602     for(RemoteIterator remote = remoteIndices_.begin(); remote != rend; ++remote, ++neighbour){
00603       MessageInformation* message;
00604       MessageInformation recv;
00605       
00606       if(messageIter != end && messageIter->first==remote->first){
00607         // We want to send message information to that process
00608         message = const_cast<MessageInformation*>(&(messageIter->second));
00609         ++messageIter;
00610       }else
00611         // We do not want to send information but the other process might.
00612         message = &dummy;
00613 
00614       sendBufferSizes_[neighbour]=0;
00615       int tsize;
00616       // The number of indices published
00617       MPI_Pack_size(1, MPI_INT,remoteIndices_.communicator(), &tsize);
00618       sendBufferSizes_[neighbour] += tsize;
00619 
00620       for(int i=0; i < message->publish; ++i){
00621         // The global index
00622         MPI_Pack_size(1, MPITraits<GlobalIndex>::getType(), remoteIndices_.communicator(), &tsize);
00623         sendBufferSizes_[neighbour] += tsize;
00624         // The attribute in the local index
00625         MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
00626         sendBufferSizes_[neighbour] += tsize;
00627         // The number of corresponding remote indices
00628         MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
00629         sendBufferSizes_[neighbour] += tsize;
00630       }
00631       for(int i=0; i < message->pairs; ++i){
00632         // The process of the remote index
00633         MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
00634         sendBufferSizes_[neighbour] += tsize;
00635         // The attribute of the remote index
00636         MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
00637         sendBufferSizes_[neighbour] += tsize;
00638       }
00639           
00640       Dune::dverb<<rank_<<": Buffer (neighbour="<<remote->first<<") size is "<< sendBufferSizes_[neighbour]<<" for publish="<<message->publish<<" pairs="<<message->pairs<<std::endl;
00641     }
00642 
00643   }
00644   
00645   template<typename T>
00646   inline void IndicesSyncer<T>::sync()
00647   {
00648     DefaultNumberer numberer;
00649     sync(numberer);
00650   }
00651   
00652   template<typename T>
00653   template<typename T1>
00654   void IndicesSyncer<T>::sync(T1& numberer)
00655   {
00656     
00657     // The pointers to the local indices in the remote indices
00658     // will become invalid due to the resorting of the index set.
00659     // Therefore store the corresponding global indices.
00660     // Mark all indices as not added
00661     
00662     typedef typename RemoteIndices::RemoteIndexMap::const_iterator 
00663       RemoteIterator;
00664 
00665     const RemoteIterator end = remoteIndices_.end();
00666     
00667     // Number of neighbours might change during the syncing.
00668     // save the old neighbours
00669     std::size_t noOldNeighbours = remoteIndices_.neighbours();
00670     int* oldNeighbours = new int[noOldNeighbours];
00671     sendBufferSizes_ = new std::size_t[noOldNeighbours];
00672     std::size_t i=0;
00673     
00674     for(RemoteIterator remote = remoteIndices_.begin(); remote != end; ++remote, ++i){
00675       typedef typename RemoteIndices::RemoteIndexList::const_iterator
00676         RemoteIndexIterator;
00677 
00678       oldNeighbours[i]=remote->first;
00679 
00680       // Make sure we only have one remote index list.
00681       assert(remote->second.first==remote->second.second);
00682 
00683       RemoteIndexList& rList = *(remote->second.first);
00684             
00685       // Store the corresponding global indices.
00686       GlobalIndexList& global = globalMap_[remote->first];
00687       BoolList& added = oldMap_[remote->first];
00688       RemoteIndexIterator riEnd = rList.end();
00689       
00690       for(RemoteIndexIterator index = rList.begin();
00691           index != riEnd; ++index){
00692         global.push_back(index->localIndexPair().global());
00693         added.push_back(true);
00694       }
00695 
00696       Iterators iterators(rList, global, added);
00697       iteratorsMap_.insert(std::make_pair(remote->first, iterators));
00698       assert(checkReset(iteratorsMap_[remote->first], rList,global,added));
00699     }
00700     
00701     // Exchange indices with each neighbour
00702     const RemoteIterator rend = remoteIndices_.end();
00703 
00704     calculateMessageSizes();
00705     
00706     // Allocate the buffers
00707     receiveBufferSize_=1;
00708     sendBuffers_ = new char*[noOldNeighbours];
00709     
00710     for(std::size_t i=0; i<noOldNeighbours; ++i){
00711       sendBuffers_[i] = new char[sendBufferSizes_[i]];
00712       receiveBufferSize_ = std::max(receiveBufferSize_, static_cast<int>(sendBufferSizes_[i]));
00713     }
00714     
00715     receiveBuffer_=new char[receiveBufferSize_];
00716     
00717     indexSet_.beginResize();
00718     
00719     Dune::dverb<<rank_<<": Neighbours: ";
00720     
00721     for(i = 0; i<noOldNeighbours; ++i)
00722       Dune::dverb<<oldNeighbours[i]<<" ";
00723     
00724     Dune::dverb<<std::endl;
00725 
00726     MPI_Request* requests = new MPI_Request[noOldNeighbours];
00727     MPI_Status* statuses = new MPI_Status[noOldNeighbours];
00728 
00729     // Pack Message data and start the sends
00730     for(i = 0; i<noOldNeighbours; ++i)
00731         packAndSend(oldNeighbours[i], sendBuffers_[i], sendBufferSizes_[i], requests[i]);
00732 
00733     // Probe for incoming messages, receive and unpack them
00734     for(i = 0; i<noOldNeighbours; ++i)
00735       recvAndUnpack(numberer);
00736 //       }else{
00737 //      recvAndUnpack(oldNeighbours[i], numberer);
00738 //      packAndSend(oldNeighbours[i]);
00739 //       }
00740 //     }
00741 
00742     delete[] receiveBuffer_;
00743     
00744     // Wait for the completion of the sends
00745     // Wait for completion of sends
00746     if(MPI_SUCCESS!=MPI_Waitall(noOldNeighbours, requests, statuses)){
00747       std::cerr<<": MPI_Error occurred while sending message"<<std::endl;
00748       for(i=0;i< noOldNeighbours;i++)
00749         if(MPI_SUCCESS!=statuses[i].MPI_ERROR)
00750           std::cerr<<"Destination "<<statuses[i].MPI_SOURCE<<" error code: "<<statuses[i].MPI_ERROR<<std::endl;
00751     }
00752   
00753     delete[] statuses;
00754     delete[] requests;
00755     
00756     for(std::size_t i=0; i<noOldNeighbours; ++i)
00757       delete[] sendBuffers_[i];
00758       
00759     delete[] sendBuffers_;
00760     
00761     // No need for the iterator tuples any more
00762     iteratorsMap_.clear();
00763     
00764     indexSet_.endResize();
00765     
00766     delete[] oldNeighbours;
00767     
00768     repairLocalIndexPointers<T,typename RemoteIndices::Allocator>(globalMap_, remoteIndices_, indexSet_);
00769     
00770     oldMap_.clear();    
00771     globalMap_.clear();
00772         
00773     // update the sequence number
00774     remoteIndices_.sourceSeqNo_ = remoteIndices_.destSeqNo_ = indexSet_.seqNo();    
00775   }
00776   
00777   template<typename T>
00778   void IndicesSyncer<T>::packAndSend(int destination, char* buffer, std::size_t bufferSize, MPI_Request& request)
00779   {
00780     typedef typename ParallelIndexSet::const_iterator IndexIterator;
00781     typedef typename RemoteIndexList::const_iterator RemoteIndexIterator;
00782     typedef typename GlobalIndexList::const_iterator GlobalIterator;
00783     typedef typename BoolList::const_iterator BoolIterator;
00784         
00785     IndexIterator      iEnd       = indexSet_.end();
00786     int                bpos       = 0;
00787     int                published  = 0;
00788     int                pairs      = 0;
00789     
00790     assert(checkReset());
00791     
00792     // Pack the number of indices we publish
00793     MPI_Pack(&(infoSend_[destination].publish), 1, MPI_INT, buffer, bufferSize, &bpos, 
00794              remoteIndices_.communicator());
00795 
00796     for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index){
00797       // Search for corresponding remote indices in all iterator tuples
00798       typedef typename IteratorsMap::iterator Iterator;
00799       Iterator iteratorsEnd = iteratorsMap_.end();
00800       
00801       // advance all iterators to a position with global index >= index->global()
00802       for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
00803         while(iterators->second.isNotAtEnd() && iterators->second.globalIndex() < index->global())
00804           ++(iterators->second);
00805           
00806       // Add all remote indices positioned at global which were already present before calling sync
00807       // to the message.
00808       // Count how many remote indices we will send
00809       int indices = 0;
00810       bool knownRemote = false; // Is the remote process supposed to know this index?
00811       
00812       for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
00813         if(iterators->second.isNotAtEnd() && iterators->second.isOld() 
00814            && iterators->second.globalIndex() == index->global()){
00815           indices++;
00816           if(destination == iterators->first)
00817             knownRemote = true;
00818         }
00819       
00820       if(!knownRemote || indices==0)
00821         // We do not need to send any indices
00822         continue;
00823       
00824       Dune::dverb<<rank_<<": sending "<<indices<<" for index "<<index->global()<<" to "<<destination<<std::endl;
00825 
00826       pairs+=indices;
00827       assert(pairs <= infoSend_[destination].pairs);
00828 
00829       // Pack the global index, the attribute and the number
00830       MPI_Pack(const_cast<GlobalIndex*>(&(index->global())), 1, MPITraits<GlobalIndex>::getType(), buffer, bufferSize, &bpos, 
00831                remoteIndices_.communicator());
00832       
00833       char attr = index->local().attribute();
00834       MPI_Pack(&attr, 1, MPI_CHAR, buffer, bufferSize, &bpos, 
00835                remoteIndices_.communicator());
00836 
00837       // Pack the number of remote indices we send.
00838       MPI_Pack(&indices, 1, MPI_INT, buffer, bufferSize, &bpos, 
00839                remoteIndices_.communicator());
00840         
00841       // Pack the information about the remote indices
00842       for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
00843         if(iterators->second.isNotAtEnd() && iterators->second.isOld() 
00844            && iterators->second.globalIndex() == index->global()){
00845           int process = iterators->first;
00846           MPI_Pack(&process, 1, MPI_INT, buffer, bufferSize, &bpos, 
00847                    remoteIndices_.communicator());
00848           char attr = iterators->second.remoteIndex().attribute();
00849           
00850           MPI_Pack(&attr, 1, MPI_CHAR, buffer, bufferSize, &bpos, 
00851                    remoteIndices_.communicator());
00852         }
00853       ++published;
00854       assert(published <= infoSend_[destination].publish);
00855     }
00856 
00857     // Make sure we send all expected entries
00858     assert(published == infoSend_[destination].publish);
00859 
00860     resetIteratorsMap();
00861 
00862     Dune::dverb << rank_<<": Sending message of "<<bpos<<" bytes to "<<destination<<std::endl;
00863     
00864     MPI_Issend(buffer, bpos, MPI_PACKED, destination, 345, remoteIndices_.communicator(),&request);
00865   }
00866 
00867   template<typename T>
00868   inline void IndicesSyncer<T>::insertIntoRemoteIndexList(int process, const GlobalIndex& global, 
00869                                                                 char attribute)
00870   {
00871     Dune::dvverb<<"Inserting from "<<process<<" "<<global<<" "<<attribute<<std::endl;
00872     
00873     // There might be cases where there no remote indices for that process yet
00874     typename IteratorsMap::iterator found = iteratorsMap_.find(process);
00875     
00876     if( found == iteratorsMap_.end() ){
00877       Dune::dverb<<"Discovered new neighbour "<<process<<std::endl;
00878       RemoteIndexList* rlist = new RemoteIndexList();
00879       remoteIndices_.remoteIndices_.insert(std::make_pair(process,std::make_pair(rlist,rlist)));
00880       Iterators iterators = Iterators(*rlist, globalMap_[process], oldMap_[process]);
00881       found = iteratorsMap_.insert(std::make_pair(process, iterators)).first;
00882     }
00883     
00884     Iterators& iterators = found->second;
00885     
00886     // Search for the remote index
00887     while(iterators.isNotAtEnd() && iterators.globalIndex() < global){
00888       // Increment all iterators
00889       ++iterators;
00890       
00891     }
00892     
00893     if(iterators.isAtEnd() || iterators.globalIndex() != global){
00894       // The entry is not yet known
00895       // Insert in the the list and do not change the first iterator.
00896       iterators.insert(RemoteIndex(Attribute(attribute)),global);
00897     }else if(iterators.isNotAtEnd()){
00898       // Assert that the attributes match 
00899       assert(iterators.remoteIndex().attribute() == attribute);
00900     }
00901   }
00902   
00903   template<typename T>
00904   template<typename T1>
00905   void IndicesSyncer<T>::recvAndUnpack(T1& numberer)
00906   {
00907     typedef typename ParallelIndexSet::const_iterator IndexIterator;
00908     typedef typename RemoteIndexList::iterator RemoteIndexIterator;
00909     typedef typename GlobalIndexList::iterator GlobalIndexIterator;
00910 
00911     IndexIterator    iEnd   = indexSet_.end();
00912     IndexIterator    index  = indexSet_.begin();
00913     int              bpos   = 0;
00914     int              publish;
00915     
00916     assert(checkReset());
00917 
00918     MPI_Status status;
00919 
00920     // We have to determine the message size and source before the receive
00921     MPI_Probe(MPI_ANY_SOURCE, 345, remoteIndices_.communicator(), &status);
00922     
00923     int source=status.MPI_SOURCE;
00924     int count;
00925     MPI_Get_count(&status, MPI_PACKED, &count);
00926 
00927     Dune::dvverb<<rank_<<": Receiving message from "<< source<<" with "<<count<<" bytes"<<std::endl;
00928 
00929     if(count>receiveBufferSize_){
00930       receiveBufferSize_=count;
00931       delete[] receiveBuffer_;
00932       receiveBuffer_ = new char[receiveBufferSize_];
00933     }
00934 
00935     MPI_Recv(receiveBuffer_, count, MPI_PACKED, source, 345, remoteIndices_.communicator(), &status);
00936     
00937     // How many global entries were published?
00938     MPI_Unpack(receiveBuffer_,  count, &bpos, &publish, 1, MPI_INT, remoteIndices_.communicator());
00939 
00940     // Now unpack the remote indices and add them.
00941     while(publish>0){
00942       
00943       // Unpack information about the local index on the source process
00944       GlobalIndex  global;          // global index of the current entry
00945       char sourceAttribute; // Attribute on the source process
00946       int pairs;
00947       
00948       MPI_Unpack(receiveBuffer_,  count, &bpos, &global, 1, MPITraits<GlobalIndex>::getType(), 
00949                  remoteIndices_.communicator());
00950       MPI_Unpack(receiveBuffer_,  count, &bpos, &sourceAttribute, 1, MPI_CHAR, 
00951                  remoteIndices_.communicator());
00952       MPI_Unpack(receiveBuffer_,  count, &bpos, &pairs, 1, MPI_INT, 
00953                  remoteIndices_.communicator());
00954     
00955       // Insert the entry on the remote process to our
00956       // remote index list
00957       insertIntoRemoteIndexList(source, global, sourceAttribute);
00958 
00959       // Unpack the remote indices
00960       while(pairs>0){
00961         // Unpack the process id that knows the index
00962         int process;
00963         char attribute;
00964         MPI_Unpack(receiveBuffer_,  count, &bpos, &process, 1, MPI_INT, 
00965                  remoteIndices_.communicator());
00966         // Unpack the attribute
00967         MPI_Unpack(receiveBuffer_,  count, &bpos, &attribute, 1, MPI_CHAR, 
00968                    remoteIndices_.communicator());
00969         
00970         if(process==rank_){
00971           // Now we know the local attribute of the global index
00972           // Do we know that global index already?
00973           IndexIterator pos = std::lower_bound(index, iEnd, IndexPair(global));
00974                   
00975           if(pos == iEnd || pos->global() != global){
00976             // No, we do not. Add it!
00977             indexSet_.add(global,ParallelLocalIndex<Attribute>(numberer(global),
00978                                                         Attribute(attribute), true));
00979           }else{
00980             // Attributes have to match!
00981             assert(attribute==pos->local().attribute());
00982           }
00983           index=pos;
00984         }else{
00985           insertIntoRemoteIndexList(process, global, attribute);
00986         }
00987         
00988         --pairs;
00989       }
00990       --publish;
00991     }
00992     
00993     resetIteratorsMap();
00994   }
00995   
00996   template<typename T>
00997   void IndicesSyncer<T>::resetIteratorsMap(){
00998     
00999     // Reset iterators in all tuples.
01000     typedef typename IteratorsMap::iterator Iterator;
01001     typedef typename RemoteIndices::RemoteIndexMap::iterator 
01002       RemoteIterator;
01003     typedef typename GlobalIndicesMap::iterator GlobalIterator;
01004     typedef typename BoolMap::iterator BoolIterator;
01005     
01006     const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
01007     Iterator iterators = iteratorsMap_.begin();
01008     GlobalIterator global = globalMap_.begin();
01009     BoolIterator added = oldMap_.begin();
01010     
01011     for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
01012         remote != remoteEnd; ++remote, ++global, ++added, ++iterators){
01013       iterators->second.reset(*(remote->second.first), global->second, added->second);
01014     }
01015   }
01016  
01017   template<typename T>
01018   bool IndicesSyncer<T>::checkReset(const Iterators& iterators, RemoteIndexList& rList, GlobalIndexList& gList,
01019                                           BoolList& bList){
01020     
01021     if(get<0>(iterators.iterators_)!=rList.begin())
01022       return false;
01023     if(get<1>(iterators.iterators_)!=gList.begin())
01024       return false;
01025     if(get<2>(iterators.iterators_)!=bList.begin())
01026       return false;
01027     return true;
01028   }
01029   
01030 
01031   template<typename T>
01032   bool IndicesSyncer<T>::checkReset(){
01033     
01034     // Reset iterators in all tuples.
01035     typedef typename IteratorsMap::iterator Iterator;
01036     typedef typename RemoteIndices::RemoteIndexMap::iterator 
01037       RemoteIterator;
01038     typedef typename GlobalIndicesMap::iterator GlobalIterator;
01039     typedef typename BoolMap::iterator BoolIterator;
01040     
01041     const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
01042     Iterator iterators = iteratorsMap_.begin();
01043     GlobalIterator global = globalMap_.begin();
01044     BoolIterator added = oldMap_.begin();
01045     bool ret = true;
01046     
01047     for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
01048         remote != remoteEnd; ++remote, ++global, ++added, ++iterators){
01049       if(!checkReset(iterators->second, *(remote->second.first), global->second,
01050                      added->second))
01051         ret=false;
01052     }
01053     return ret;
01054   } 
01055 }
01056 
01057 #endif
01058 #endif

Generated on Sun Nov 15 22:29:36 2009 for dune-istl by  doxygen 1.5.6