dune-common  2.2.1
indicessyncer.hh
Go to the documentation of this file.
1 // $Id$
2 #ifndef DUNE_INDICESSYNCER_HH
3 #define DUNE_INDICESSYNCER_HH
4 
5 #include"indexset.hh"
6 #include"remoteindices.hh"
8 #include<dune/common/tuples.hh>
9 #include<dune/common/sllist.hh>
10 #include<cassert>
11 #include<cmath>
12 #include<limits>
13 #include<algorithm>
14 #include<functional>
15 #include<map>
16 
17 #if HAVE_MPI
18 namespace Dune
19 {
36  template<typename T>
38  {
39  public:
40 
42  typedef T ParallelIndexSet;
43 
46 
49 
51  typedef typename ParallelIndexSet::LocalIndex::Attribute Attribute;
52 
57 
68  RemoteIndices& remoteIndices);
69 
79  void sync();
80 
91  template<typename T1>
92  void sync(T1& numberer);
93 
94  private:
95 
97  ParallelIndexSet& indexSet_;
98 
100  RemoteIndices& remoteIndices_;
101 
103  char** sendBuffers_;
104 
106  char* receiveBuffer_;
107 
109  std::size_t* sendBufferSizes_;
110 
112  int receiveBufferSize_; // int because of MPI
113 
117  struct MessageInformation
118  {
119  MessageInformation()
120  : publish(), pairs()
121  {}
123  int publish;
128  int pairs;
129  };
130 
134  class DefaultNumberer
135  {
136  public:
142  std::size_t operator()(const GlobalIndex& global)
143  {
144  return std::numeric_limits<size_t>::max();
145  }
146  };
147 
149  MPI_Datatype datatype_;
150 
152  int rank_;
153 
158  typedef SLList<std::pair<GlobalIndex,Attribute>, typename RemoteIndices::Allocator> GlobalIndexList;
159 
161  typedef typename GlobalIndexList::ModifyIterator GlobalIndexModifier;
162 
167  GlobalIndexIterator;
168 
170  typedef std::map<int, GlobalIndexList> GlobalIndicesMap;
171 
180  GlobalIndicesMap globalMap_;
181 
186 
190  typedef typename BoolList::iterator BoolIterator;
191 
193  typedef typename BoolList::ModifyIterator BoolListModifier;
194 
196  typedef std::map<int,BoolList> BoolMap;
197 
202  BoolMap oldMap_;
203 
205  std::map<int,MessageInformation> infoSend_;
206 
208  typedef typename RemoteIndices::RemoteIndexList RemoteIndexList;
209 
211  typedef typename RemoteIndexList::ModifyIterator RemoteIndexModifier;
212 
215 
217  typedef typename RemoteIndexList::iterator RemoteIndexIterator;
218 
220  typedef typename RemoteIndexList::const_iterator ConstRemoteIndexIterator;
221 
223  typedef Dune::tuple<RemoteIndexModifier,GlobalIndexModifier,BoolListModifier,
224  const ConstRemoteIndexIterator> IteratorTuple;
225 
233  class Iterators
234  {
235  friend class IndicesSyncer<T>;
236  public:
246  Iterators(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
247  BoolList& booleans);
248 
252  Iterators();
253 
257  Iterators& operator++();
258 
264  void insert(const RemoteIndex& index,
265  const std::pair<GlobalIndex,Attribute>& global);
266 
271  RemoteIndex& remoteIndex() const;
272 
277  std::pair<GlobalIndex,Attribute>& globalIndexPair() const;
278 
279  Attribute& attribute() const;
280 
286  bool isOld() const;
287 
297  void reset(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
298  BoolList& booleans);
299 
305  bool isNotAtEnd() const;
306 
312  bool isAtEnd() const;
313 
314  private:
324  IteratorTuple iterators_;
325  };
326 
328  typedef std::map<int,Iterators> IteratorsMap;
329 
341  IteratorsMap iteratorsMap_;
342 
344  void calculateMessageSizes();
345 
353  void packAndSend(int destination, char* buffer, std::size_t bufferSize, MPI_Request& req);
354 
359  template<typename T1>
360  void recvAndUnpack(T1& numberer);
361 
365  void registerMessageDatatype();
366 
370  void insertIntoRemoteIndexList(int process,
371  const std::pair<GlobalIndex,Attribute>& global,
372  char attribute);
373 
377  void resetIteratorsMap();
378 
383  bool checkReset();
384 
393  bool checkReset(const Iterators& iterators, RemoteIndexList& rlist, GlobalIndexList& gList,
394  BoolList& bList);
395  };
396 
397  template<typename TG, typename TA>
398  bool operator<(const IndexPair<TG,ParallelLocalIndex<TA> >& i1,
399  const std::pair<TG,TA>& i2)
400  {
401  return i1.global() < i2.first ||
402  (i1.global() == i2.first && i1.local().attribute()<i2.second);
403  }
404 
405  template<typename TG, typename TA>
406  bool operator<(const std::pair<TG,TA>& i1,
408  {
409  return i1.first < i2.global() ||
410  (i1.first == i2.global() && i1.second<i2.local().attribute());
411  }
412 
413  template<typename TG, typename TA>
415  const std::pair<TG,TA>& i2)
416  {
417  return (i1.global() == i2.first && i1.local().attribute()==i2.second);
418  }
419 
420  template<typename TG, typename TA>
422  const std::pair<TG,TA>& i2)
423  {
424  return (i1.global() != i2.first || i1.local().attribute()!=i2.second);
425  }
426 
427  template<typename TG, typename TA>
428  bool operator==(const std::pair<TG,TA>& i2,
429  const IndexPair<TG,ParallelLocalIndex<TA> >& i1)
430  {
431  return (i1.global() == i2.first && i1.local().attribute()==i2.second);
432  }
433 
434  template<typename TG, typename TA>
435  bool operator!=(const std::pair<TG,TA>& i2,
436  const IndexPair<TG,ParallelLocalIndex<TA> >& i1)
437  {
438  return (i1.global() != i2.first || i1.local().attribute()!=i2.second);
439  }
440 
457  template<typename T, typename A, typename A1>
458  void storeGlobalIndicesOfRemoteIndices(std::map<int,SLList<std::pair<typename T::GlobalIndex, typename T::LocalIndex::Attribute>,A> >& globalMap,
459  const RemoteIndices<T,A1>& remoteIndices,
460  const T& indexSet)
461  {
462  typedef typename RemoteIndices<T,A1>::const_iterator RemoteIterator;
463 
464  for(RemoteIterator remote = remoteIndices.begin(), end =remoteIndices.end(); remote != end; ++remote){
465  typedef typename RemoteIndices<T,A1>::RemoteIndexList RemoteIndexList;
466  typedef typename RemoteIndexList::const_iterator RemoteIndexIterator;
468  GlobalIndexList& global = globalMap[remote->first];
469  RemoteIndexList& rList = *(remote->second.first);
470 
471  for(RemoteIndexIterator index = rList.begin(), riEnd = rList.end();
472  index != riEnd; ++index){
473  global.push_back(std::make_pair(index->localIndexPair().global(),
474  index->localIndexPair().local().attribute()));
475  }
476  }
477  }
478 
487  template<typename T, typename A, typename A1>
488  inline void repairLocalIndexPointers(std::map<int,
489  SLList<std::pair<typename T::GlobalIndex,
490  typename T::LocalIndex::Attribute>,A> >& globalMap,
491  RemoteIndices<T,A1>& remoteIndices,
492  const T& indexSet)
493  {
494  typedef typename RemoteIndices<T,A1>::RemoteIndexMap::iterator RemoteIterator;
495  typedef typename RemoteIndices<T,A1>::RemoteIndexList::iterator RemoteIndexIterator;
496  typedef typename T::GlobalIndex GlobalIndex;
497  typedef typename T::LocalIndex::Attribute Attribute;
498  typedef std::pair<GlobalIndex,Attribute> GlobalIndexPair;
499  typedef SLList<GlobalIndexPair,A> GlobalIndexPairList;
500  typedef typename GlobalIndexPairList::iterator GlobalIndexIterator;
501 
502  assert(globalMap.size()==static_cast<std::size_t>(remoteIndices.neighbours()));
503  // Repair pointers to index set in remote indices.
504  typename std::map<int,GlobalIndexPairList>::iterator global = globalMap.begin();
505  RemoteIterator end = remoteIndices.remoteIndices_.end();
506 
507  for(RemoteIterator remote = remoteIndices.remoteIndices_.begin(); remote != end; ++remote, ++global){
508  typedef typename T::const_iterator IndexIterator;
509 
510  assert(remote->first==global->first);
511  assert(remote->second.first->size() == global->second.size());
512 
513  RemoteIndexIterator riEnd = remote->second.first->end();
514  RemoteIndexIterator rIndex = remote->second.first->begin();
515  GlobalIndexIterator gIndex = global->second.begin();
516  IndexIterator index = indexSet.begin();
517 
518  assert(rIndex==riEnd || gIndex != global->second.end());
519  while(rIndex != riEnd){
520  // Search for the index in the set.
521  assert(gIndex != global->second.end());
522 
523  while(!(index->global() == gIndex->first
524  && index->local().attribute() == gIndex->second)){
525  ++index;
526  // this is only needed for ALU, where there may exist
527  // more entries with the same global index in the remote index set
528  // than in the index set
529  if (index->global() > gIndex->first){
530  index=indexSet.begin();
531  }
532  }
533 
534  assert(index != indexSet.end() && *index == *gIndex);
535 
536  rIndex->localIndex_ = &(*index);
537  ++index;
538  ++rIndex;
539  ++gIndex;
540  }
541  }
542  remoteIndices.sourceSeqNo_ = remoteIndices.source_->seqNo();
543  remoteIndices.destSeqNo_ = remoteIndices.target_->seqNo();
544  }
545 
546  template<typename T>
548  RemoteIndices& remoteIndices)
549  : indexSet_(indexSet), remoteIndices_(remoteIndices)
550  {
551  // index sets must match.
552  assert(remoteIndices.source_ == remoteIndices.target_);
553  assert(remoteIndices.source_ == &indexSet);
554  MPI_Comm_rank(remoteIndices_.communicator(), &rank_);
555  }
556 
557  template<typename T>
559  GlobalIndexList& globalIndices,
560  BoolList& booleans)
561  : iterators_(remoteIndices.beginModify(), globalIndices.beginModify(),
562  booleans.beginModify(), remoteIndices.end())
563  { }
564 
565  template<typename T>
567  : iterators_()
568  {}
569 
570  template<typename T>
572  {
573  ++(get<0>(iterators_));
574  ++(get<1>(iterators_));
575  ++(get<2>(iterators_));
576  return *this;
577  }
578 
579  template<typename T>
581  const std::pair<GlobalIndex,Attribute>& global)
582  {
583  get<0>(iterators_).insert(index);
584  get<1>(iterators_).insert(global);
585  get<2>(iterators_).insert(false);
586  }
587 
588  template<typename T>
589  inline typename IndicesSyncer<T>::RemoteIndex&
591  {
592  return *(get<0>(iterators_));
593  }
594 
595  template<typename T>
596  inline std::pair<typename IndicesSyncer<T>::GlobalIndex,typename IndicesSyncer<T>::Attribute>&
598  {
599  return *(get<1>(iterators_));
600  }
601 
602  template<typename T>
604  {
605  return *(get<2>(iterators_));
606  }
607 
608  template<typename T>
610  GlobalIndexList& globalIndices,
611  BoolList& booleans)
612  {
613  get<0>(iterators_) = remoteIndices.beginModify();
614  get<1>(iterators_) = globalIndices.beginModify();
615  get<2>(iterators_) = booleans.beginModify();
616  }
617 
618  template<typename T>
620  {
621  return get<0>(iterators_)!=get<3>(iterators_);
622  }
623 
624  template<typename T>
626  {
627  return get<0>(iterators_)==get<3>(iterators_);
628  }
629 
630  template<typename T>
632  {
633  MPI_Datatype type[2] = {MPI_INT, MPI_INT};
634  int blocklength[2] = {1,1};
635  MPI_Aint displacement[2];
636  MPI_Aint base;
637 
638  // Compute displacement
639  MessageInformation message;
640 
641  MPI_Address( &(message.publish), displacement);
642  MPI_Address( &(message.pairs), displacement+1);
643 
644  // Make the displacement relative
645  MPI_Address(&message, &base);
646  displacement[0] -= base;
647  displacement[1] -= base;
648 
649  MPI_Type_struct( 2, blocklength, displacement, type, &datatype_);
650  MPI_Type_commit(&datatype_);
651  }
652 
653  template<typename T>
654  void IndicesSyncer<T>::calculateMessageSizes()
655  {
656  typedef typename ParallelIndexSet::const_iterator IndexIterator;
657  typedef CollectiveIterator<T,typename RemoteIndices::Allocator> CollectiveIterator;
658 
659  IndexIterator iEnd = indexSet_.end();
660  CollectiveIterator collIter = remoteIndices_.template iterator<true>();
661 
662  for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index){
663  collIter.advance(index->global(), index->local().attribute());
664  if(collIter.empty())
665  break;
666  int knownRemote=0;
667 
668  typedef typename CollectiveIterator::iterator ValidIterator;
669  ValidIterator end = collIter.end();
670 
671  // Count the remote indices we know.
672  for(ValidIterator valid = collIter.begin(); valid != end; ++valid){
673  ++knownRemote;
674  }
675 
676  if(knownRemote>0){
677  Dune::dverb<<rank_<<": publishing "<<knownRemote<<" for index "<<index->global()<< " for processes ";
678 
679  // Update MessageInformation
680  for(ValidIterator valid = collIter.begin(); valid != end; ++valid){
681  ++(infoSend_[valid.process()].publish);
682  (infoSend_[valid.process()].pairs) += knownRemote;
683  Dune::dverb<<valid.process()<<" ";
684  Dune::dverb<<"(publish="<<infoSend_[valid.process()].publish<<", pairs="<<infoSend_[valid.process()].pairs
685  <<") ";
686  }
687  Dune::dverb<<std::endl;
688  }
689  }
690 
691  typedef typename std::map<int,MessageInformation>::const_iterator
692  MessageIterator;
693 
694  const MessageIterator end = infoSend_.end();
695 
696  // registerMessageDatatype();
697 
698  // Now determine the buffersizes needed for each neighbour using MPI_Pack_size
699  MessageInformation dummy;
700 
701  MessageIterator messageIter= infoSend_.begin();
702 
703  typedef typename RemoteIndices::RemoteIndexMap::const_iterator RemoteIterator;
704  const RemoteIterator rend = remoteIndices_.end();
705  int neighbour=0;
706 
707  for(RemoteIterator remote = remoteIndices_.begin(); remote != rend; ++remote, ++neighbour){
708  MessageInformation* message;
709  MessageInformation recv;
710 
711  if(messageIter != end && messageIter->first==remote->first){
712  // We want to send message information to that process
713  message = const_cast<MessageInformation*>(&(messageIter->second));
714  ++messageIter;
715  }else
716  // We do not want to send information but the other process might.
717  message = &dummy;
718 
719  sendBufferSizes_[neighbour]=0;
720  int tsize;
721  // The number of indices published
722  MPI_Pack_size(1, MPI_INT,remoteIndices_.communicator(), &tsize);
723  sendBufferSizes_[neighbour] += tsize;
724 
725  for(int i=0; i < message->publish; ++i){
726  // The global index
727  MPI_Pack_size(1, MPITraits<GlobalIndex>::getType(), remoteIndices_.communicator(), &tsize);
728  sendBufferSizes_[neighbour] += tsize;
729  // The attribute in the local index
730  MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
731  sendBufferSizes_[neighbour] += tsize;
732  // The number of corresponding remote indices
733  MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
734  sendBufferSizes_[neighbour] += tsize;
735  }
736  for(int i=0; i < message->pairs; ++i){
737  // The process of the remote index
738  MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
739  sendBufferSizes_[neighbour] += tsize;
740  // The attribute of the remote index
741  MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
742  sendBufferSizes_[neighbour] += tsize;
743  }
744 
745  Dune::dverb<<rank_<<": Buffer (neighbour="<<remote->first<<") size is "<< sendBufferSizes_[neighbour]<<" for publish="<<message->publish<<" pairs="<<message->pairs<<std::endl;
746  }
747 
748  }
749 
750  template<typename T>
752  {
753  DefaultNumberer numberer;
754  sync(numberer);
755  }
756 
757  template<typename T>
758  template<typename T1>
759  void IndicesSyncer<T>::sync(T1& numberer)
760  {
761 
762  // The pointers to the local indices in the remote indices
763  // will become invalid due to the resorting of the index set.
764  // Therefore store the corresponding global indices.
765  // Mark all indices as not added
766 
768  RemoteIterator;
769 
770  const RemoteIterator end = remoteIndices_.end();
771 
772  // Number of neighbours might change during the syncing.
773  // save the old neighbours
774  std::size_t noOldNeighbours = remoteIndices_.neighbours();
775  int* oldNeighbours = new int[noOldNeighbours];
776  sendBufferSizes_ = new std::size_t[noOldNeighbours];
777  std::size_t i=0;
778 
779  for(RemoteIterator remote = remoteIndices_.begin(); remote != end; ++remote, ++i){
782 
783  oldNeighbours[i]=remote->first;
784 
785  // Make sure we only have one remote index list.
786  assert(remote->second.first==remote->second.second);
787 
788  RemoteIndexList& rList = *(remote->second.first);
789 
790  // Store the corresponding global indices.
791  GlobalIndexList& global = globalMap_[remote->first];
792  BoolList& added = oldMap_[remote->first];
793  RemoteIndexIterator riEnd = rList.end();
794 
795  for(RemoteIndexIterator index = rList.begin();
796  index != riEnd; ++index){
797  global.push_back(std::make_pair(index->localIndexPair().global(),
798  index->localIndexPair().local().attribute()));
799  added.push_back(true);
800  }
801 
802  Iterators iterators(rList, global, added);
803  iteratorsMap_.insert(std::make_pair(remote->first, iterators));
804  assert(checkReset(iteratorsMap_[remote->first], rList,global,added));
805  }
806 
807  // Exchange indices with each neighbour
808  calculateMessageSizes();
809 
810  // Allocate the buffers
811  receiveBufferSize_=1;
812  sendBuffers_ = new char*[noOldNeighbours];
813 
814  for(std::size_t i=0; i<noOldNeighbours; ++i){
815  sendBuffers_[i] = new char[sendBufferSizes_[i]];
816  receiveBufferSize_ = std::max(receiveBufferSize_, static_cast<int>(sendBufferSizes_[i]));
817  }
818 
819  receiveBuffer_=new char[receiveBufferSize_];
820 
821  indexSet_.beginResize();
822 
823  Dune::dverb<<rank_<<": Neighbours: ";
824 
825  for(i = 0; i<noOldNeighbours; ++i)
826  Dune::dverb<<oldNeighbours[i]<<" ";
827 
828  Dune::dverb<<std::endl;
829 
830  MPI_Request* requests = new MPI_Request[noOldNeighbours];
831  MPI_Status* statuses = new MPI_Status[noOldNeighbours];
832 
833  // Pack Message data and start the sends
834  for(i = 0; i<noOldNeighbours; ++i)
835  packAndSend(oldNeighbours[i], sendBuffers_[i], sendBufferSizes_[i], requests[i]);
836 
837  // Probe for incoming messages, receive and unpack them
838  for(i = 0; i<noOldNeighbours; ++i)
839  recvAndUnpack(numberer);
840 // }else{
841 // recvAndUnpack(oldNeighbours[i], numberer);
842 // packAndSend(oldNeighbours[i]);
843 // }
844 // }
845 
846  delete[] receiveBuffer_;
847 
848  // Wait for the completion of the sends
849  // Wait for completion of sends
850  if(MPI_SUCCESS!=MPI_Waitall(noOldNeighbours, requests, statuses)){
851  std::cerr<<": MPI_Error occurred while sending message"<<std::endl;
852  for(i=0;i< noOldNeighbours;i++)
853  if(MPI_SUCCESS!=statuses[i].MPI_ERROR)
854  std::cerr<<"Destination "<<statuses[i].MPI_SOURCE<<" error code: "<<statuses[i].MPI_ERROR<<std::endl;
855  }
856 
857  delete[] statuses;
858  delete[] requests;
859 
860  for(std::size_t i=0; i<noOldNeighbours; ++i)
861  delete[] sendBuffers_[i];
862 
863  delete[] sendBuffers_;
864  delete[] sendBufferSizes_;
865 
866  // No need for the iterator tuples any more
867  iteratorsMap_.clear();
868 
869  indexSet_.endResize();
870 
871  delete[] oldNeighbours;
872 
873  repairLocalIndexPointers(globalMap_, remoteIndices_, indexSet_);
874 
875  oldMap_.clear();
876  globalMap_.clear();
877 
878  // update the sequence number
879  remoteIndices_.sourceSeqNo_ = remoteIndices_.destSeqNo_ = indexSet_.seqNo();
880  }
881 
882  template<typename T>
883  void IndicesSyncer<T>::packAndSend(int destination, char* buffer, std::size_t bufferSize, MPI_Request& request)
884  {
885  typedef typename ParallelIndexSet::const_iterator IndexIterator;
886 
887  IndexIterator iEnd = indexSet_.end();
888  int bpos = 0;
889  int published = 0;
890  int pairs = 0;
891 
892  assert(checkReset());
893 
894  // Pack the number of indices we publish
895  MPI_Pack(&(infoSend_[destination].publish), 1, MPI_INT, buffer, bufferSize, &bpos,
896  remoteIndices_.communicator());
897 
898  for(IndexIterator index = indexSet_.begin(); index != iEnd; ++index){
899  // Search for corresponding remote indices in all iterator tuples
900  typedef typename IteratorsMap::iterator Iterator;
901  Iterator iteratorsEnd = iteratorsMap_.end();
902 
903  // advance all iterators to a position with global index >= index->global()
904  for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators){
905  while(iterators->second.isNotAtEnd() &&
906  iterators->second.globalIndexPair().first < index->global())
907  ++(iterators->second);
908  assert(!iterators->second.isNotAtEnd() || iterators->second.globalIndexPair().first >= index->global());
909  }
910 
911  // Add all remote indices positioned at global which were already present before calling sync
912  // to the message.
913  // Count how many remote indices we will send
914  int indices = 0;
915  bool knownRemote = false; // Is the remote process supposed to know this index?
916 
917  for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
918  {
919  std::pair<GlobalIndex,Attribute> p;
920  if (iterators->second.isNotAtEnd())
921  {
922  p = iterators->second.globalIndexPair();
923  }
924 
925  if(iterators->second.isNotAtEnd() && iterators->second.isOld()
926  && iterators->second.globalIndexPair().first == index->global()){
927  indices++;
928  if(destination == iterators->first)
929  knownRemote = true;
930  }
931  }
932 
933  if(!knownRemote)
934  // We do not need to send any indices
935  continue;
936 
937  Dune::dverb<<rank_<<": sending "<<indices<<" for index "<<index->global()<<" to "<<destination<<std::endl;
938 
939 
940  // Pack the global index, the attribute and the number
941  MPI_Pack(const_cast<GlobalIndex*>(&(index->global())), 1, MPITraits<GlobalIndex>::getType(), buffer, bufferSize, &bpos,
942  remoteIndices_.communicator());
943 
944  char attr = index->local().attribute();
945  MPI_Pack(&attr, 1, MPI_CHAR, buffer, bufferSize, &bpos,
946  remoteIndices_.communicator());
947 
948  // Pack the number of remote indices we send.
949  MPI_Pack(&indices, 1, MPI_INT, buffer, bufferSize, &bpos,
950  remoteIndices_.communicator());
951 
952  // Pack the information about the remote indices
953  for(Iterator iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
954  if(iterators->second.isNotAtEnd() && iterators->second.isOld()
955  && iterators->second.globalIndexPair().first == index->global()){
956  int process = iterators->first;
957 
958  ++pairs;
959  assert(pairs <= infoSend_[destination].pairs);
960  MPI_Pack(&process, 1, MPI_INT, buffer, bufferSize, &bpos,
961  remoteIndices_.communicator());
962  char attr = iterators->second.remoteIndex().attribute();
963 
964  MPI_Pack(&attr, 1, MPI_CHAR, buffer, bufferSize, &bpos,
965  remoteIndices_.communicator());
966  --indices;
967  }
968  assert(indices==0);
969  ++published;
970  Dune::dvverb<<" (publish="<<published<<", pairs="<<pairs<<")"<<std::endl;
971  assert(published <= infoSend_[destination].publish);
972  }
973 
974  // Make sure we send all expected entries
975  assert(published == infoSend_[destination].publish);
976  assert(pairs == infoSend_[destination].pairs);
977  resetIteratorsMap();
978 
979  Dune::dverb << rank_<<": Sending message of "<<bpos<<" bytes to "<<destination<<std::endl;
980 
981  MPI_Issend(buffer, bpos, MPI_PACKED, destination, 345, remoteIndices_.communicator(),&request);
982  }
983 
984  template<typename T>
985  inline void IndicesSyncer<T>::insertIntoRemoteIndexList(int process,
986  const std::pair<GlobalIndex,Attribute>& globalPair,
987  char attribute)
988  {
989  Dune::dverb<<"Inserting from "<<process<<" "<<globalPair.first<<", "<<
990  globalPair.second<<" "<<attribute<<std::endl;
991 
992  resetIteratorsMap();
993 
994  // There might be cases where there no remote indices for that process yet
995  typename IteratorsMap::iterator found = iteratorsMap_.find(process);
996 
997  if( found == iteratorsMap_.end() ){
998  Dune::dverb<<"Discovered new neighbour "<<process<<std::endl;
999  RemoteIndexList* rlist = new RemoteIndexList();
1000  remoteIndices_.remoteIndices_.insert(std::make_pair(process,std::make_pair(rlist,rlist)));
1001  Iterators iterators = Iterators(*rlist, globalMap_[process], oldMap_[process]);
1002  found = iteratorsMap_.insert(std::make_pair(process, iterators)).first;
1003  }
1004 
1005  Iterators& iterators = found->second;
1006 
1007  // Search for the remote index
1008  while(iterators.isNotAtEnd() && iterators.globalIndexPair() < globalPair){
1009  // Increment all iterators
1010  ++iterators;
1011 
1012  }
1013 
1014  if(iterators.isAtEnd() || iterators.globalIndexPair() != globalPair){
1015  // The entry is not yet known
1016  // Insert in the the list and do not change the first iterator.
1017  iterators.insert(RemoteIndex(Attribute(attribute)),globalPair);
1018  return;
1019  }
1020 
1021  // Global indices match
1022  bool indexIsThere=false;
1023  for(Iterators tmpIterators = iterators;
1024  !tmpIterators.isAtEnd() && tmpIterators.globalIndexPair() == globalPair;
1025  ++tmpIterators)
1026  //entry already exists with the same attribute
1027  if(tmpIterators.globalIndexPair().second == attribute){
1028  indexIsThere=true;
1029  break;
1030  }
1031 
1032  if(!indexIsThere)
1033  // The entry is not yet known
1034  // Insert in the the list and do not change the first iterator.
1035  iterators.insert(RemoteIndex(Attribute(attribute)),globalPair);
1036  }
1037 
1038  template<typename T>
1039  template<typename T1>
1040  void IndicesSyncer<T>::recvAndUnpack(T1& numberer)
1041  {
1042  typedef typename ParallelIndexSet::const_iterator IndexIterator;
1043 
1044  IndexIterator iEnd = indexSet_.end();
1045  IndexIterator index = indexSet_.begin();
1046  int bpos = 0;
1047  int publish;
1048 
1049  assert(checkReset());
1050 
1051  MPI_Status status;
1052 
1053  // We have to determine the message size and source before the receive
1054  MPI_Probe(MPI_ANY_SOURCE, 345, remoteIndices_.communicator(), &status);
1055 
1056  int source=status.MPI_SOURCE;
1057  int count;
1058  MPI_Get_count(&status, MPI_PACKED, &count);
1059 
1060  Dune::dvverb<<rank_<<": Receiving message from "<< source<<" with "<<count<<" bytes"<<std::endl;
1061 
1062  if(count>receiveBufferSize_){
1063  receiveBufferSize_=count;
1064  delete[] receiveBuffer_;
1065  receiveBuffer_ = new char[receiveBufferSize_];
1066  }
1067 
1068  MPI_Recv(receiveBuffer_, count, MPI_PACKED, source, 345, remoteIndices_.communicator(), &status);
1069 
1070  // How many global entries were published?
1071  MPI_Unpack(receiveBuffer_, count, &bpos, &publish, 1, MPI_INT, remoteIndices_.communicator());
1072 
1073  // Now unpack the remote indices and add them.
1074  while(publish>0){
1075 
1076  // Unpack information about the local index on the source process
1077  GlobalIndex global; // global index of the current entry
1078  char sourceAttribute; // Attribute on the source process
1079  int pairs;
1080 
1081  MPI_Unpack(receiveBuffer_, count, &bpos, &global, 1, MPITraits<GlobalIndex>::getType(),
1082  remoteIndices_.communicator());
1083  MPI_Unpack(receiveBuffer_, count, &bpos, &sourceAttribute, 1, MPI_CHAR,
1084  remoteIndices_.communicator());
1085  MPI_Unpack(receiveBuffer_, count, &bpos, &pairs, 1, MPI_INT,
1086  remoteIndices_.communicator());
1087 
1088  // Insert the entry on the remote process to our
1089  // remote index list
1090  SLList<std::pair<int,Attribute> > sourceAttributeList;
1091  sourceAttributeList.push_back(std::make_pair(source,Attribute(sourceAttribute)));
1092 #ifndef NDEBUG
1093  bool foundSelf = false;
1094 #endif
1095  Attribute myAttribute=Attribute();
1096 
1097  // Unpack the remote indices
1098  for(;pairs>0;--pairs){
1099  // Unpack the process id that knows the index
1100  int process;
1101  char attribute;
1102  MPI_Unpack(receiveBuffer_, count, &bpos, &process, 1, MPI_INT,
1103  remoteIndices_.communicator());
1104  // Unpack the attribute
1105  MPI_Unpack(receiveBuffer_, count, &bpos, &attribute, 1, MPI_CHAR,
1106  remoteIndices_.communicator());
1107 
1108  if(process==rank_){
1109 #ifndef NDEBUG
1110  foundSelf=true;
1111 #endif
1112  myAttribute=Attribute(attribute);
1113  // Now we know the local attribute of the global index
1114  //Only add the index if it is unknown.
1115  // Do we know that global index already?
1116  IndexIterator pos = std::lower_bound(index, iEnd, IndexPair(global));
1117 
1118  if(pos == iEnd || pos->global() != global){
1119  // no entry with this global index
1120  indexSet_.add(global,
1121  ParallelLocalIndex<Attribute>(numberer(global),
1122  myAttribute, true));
1123  Dune::dvverb << "Adding "<<global<<" "<<myAttribute<<std::endl;
1124  continue;
1125  }
1126 
1127  // because of above the global indices match. Add only if the attribute is different
1128  bool indexIsThere = false;
1129  index=pos;
1130 
1131  for(;pos->global()==global;++pos)
1132  if(pos->local().attribute() == myAttribute){
1133  Dune::dvverb<<"found "<<global<<" "<<myAttribute<<std::endl;
1134  indexIsThere = true;
1135  break;
1136  }
1137 
1138  if(!indexIsThere){
1139  indexSet_.add(global,
1140  ParallelLocalIndex<Attribute>(numberer(global),
1141  myAttribute, true));
1142  Dune::dvverb << "Adding "<<global<<" "<<myAttribute<<std::endl;
1143  }
1144 
1145  }else{
1146  sourceAttributeList.push_back(std::make_pair(process,Attribute(attribute)));
1147  }
1148  }
1149  assert(foundSelf);
1150  // Insert remote indices
1151  typedef typename SLList<std::pair<int,Attribute> >::const_iterator Iter;
1152  for(Iter i=sourceAttributeList.begin(), end=sourceAttributeList.end();
1153  i!=end; ++i)
1154  insertIntoRemoteIndexList(i->first, std::make_pair(global, myAttribute),
1155  i->second);
1156  --publish;
1157  }
1158 
1159  resetIteratorsMap();
1160  }
1161 
1162  template<typename T>
1163  void IndicesSyncer<T>::resetIteratorsMap(){
1164 
1165  // Reset iterators in all tuples.
1166  typedef typename IteratorsMap::iterator Iterator;
1167  typedef typename RemoteIndices::RemoteIndexMap::iterator
1168  RemoteIterator;
1169  typedef typename GlobalIndicesMap::iterator GlobalIterator;
1170  typedef typename BoolMap::iterator BoolIterator;
1171 
1172  const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
1173  Iterator iterators = iteratorsMap_.begin();
1174  GlobalIterator global = globalMap_.begin();
1175  BoolIterator added = oldMap_.begin();
1176 
1177  for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
1178  remote != remoteEnd; ++remote, ++global, ++added, ++iterators){
1179  iterators->second.reset(*(remote->second.first), global->second, added->second);
1180  }
1181  }
1182 
1183  template<typename T>
1184  bool IndicesSyncer<T>::checkReset(const Iterators& iterators, RemoteIndexList& rList, GlobalIndexList& gList,
1185  BoolList& bList){
1186 
1187  if(get<0>(iterators.iterators_)!=rList.begin())
1188  return false;
1189  if(get<1>(iterators.iterators_)!=gList.begin())
1190  return false;
1191  if(get<2>(iterators.iterators_)!=bList.begin())
1192  return false;
1193  return true;
1194  }
1195 
1196 
1197  template<typename T>
1198  bool IndicesSyncer<T>::checkReset(){
1199 
1200  // Reset iterators in all tuples.
1201  typedef typename IteratorsMap::iterator Iterator;
1202  typedef typename RemoteIndices::RemoteIndexMap::iterator
1203  RemoteIterator;
1204  typedef typename GlobalIndicesMap::iterator GlobalIterator;
1205  typedef typename BoolMap::iterator BoolIterator;
1206 
1207  const RemoteIterator remoteEnd = remoteIndices_.remoteIndices_.end();
1208  Iterator iterators = iteratorsMap_.begin();
1209  GlobalIterator global = globalMap_.begin();
1210  BoolIterator added = oldMap_.begin();
1211  bool ret = true;
1212 
1213  for(RemoteIterator remote = remoteIndices_.remoteIndices_.begin();
1214  remote != remoteEnd; ++remote, ++global, ++added, ++iterators){
1215  if(!checkReset(iterators->second, *(remote->second.first), global->second,
1216  added->second))
1217  ret=false;
1218  }
1219  return ret;
1220  }
1221 }
1222 
1223 #endif
1224 #endif