dune-common  2.8.0
indicessyncer.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_INDICESSYNCER_HH
4 #define DUNE_INDICESSYNCER_HH
5 
6 #include "indexset.hh"
7 #include "remoteindices.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 #include <tuple>
17 
18 #if HAVE_MPI
19 namespace Dune
20 {
37  template<typename T>
39  {
40  public:
41 
43  typedef T ParallelIndexSet;
44 
47 
50 
52  typedef typename ParallelIndexSet::LocalIndex::Attribute Attribute;
53 
58 
69  RemoteIndices& remoteIndices);
70 
80  void sync();
81 
92  template<typename T1>
93  void sync(T1& numberer);
94 
95  private:
96 
98  ParallelIndexSet& indexSet_;
99 
101  RemoteIndices& remoteIndices_;
102 
104  char** sendBuffers_;
105 
107  char* receiveBuffer_;
108 
110  std::size_t* sendBufferSizes_;
111 
113  int receiveBufferSize_; // int because of MPI
114 
118  struct MessageInformation
119  {
121  : publish(), pairs()
122  {}
124  int publish;
129  int pairs;
130  };
131 
135  class DefaultNumberer
136  {
137  public:
143  std::size_t operator()([[maybe_unused]] const GlobalIndex& global)
144  {
146  }
147  };
148 
150  MPI_Datatype datatype_;
151 
153  int rank_;
154 
159  typedef SLList<std::pair<GlobalIndex,Attribute>, typename RemoteIndices::Allocator> GlobalIndexList;
160 
162  typedef typename GlobalIndexList::ModifyIterator GlobalIndexModifier;
163 
168  GlobalIndexIterator;
169 
171  typedef std::map<int, GlobalIndexList> GlobalIndicesMap;
172 
181  GlobalIndicesMap globalMap_;
182 
187 
191  typedef typename BoolList::iterator BoolIterator;
192 
194  typedef typename BoolList::ModifyIterator BoolListModifier;
195 
197  typedef std::map<int,BoolList> BoolMap;
198 
203  BoolMap oldMap_;
204 
206  std::map<int,MessageInformation> infoSend_;
207 
209  typedef typename RemoteIndices::RemoteIndexList RemoteIndexList;
210 
212  typedef typename RemoteIndexList::ModifyIterator RemoteIndexModifier;
213 
216 
218  typedef typename RemoteIndexList::iterator RemoteIndexIterator;
219 
221  typedef typename RemoteIndexList::const_iterator ConstRemoteIndexIterator;
222 
224  typedef std::tuple<RemoteIndexModifier,GlobalIndexModifier,BoolListModifier,
225  const ConstRemoteIndexIterator> IteratorTuple;
226 
234  class Iterators
235  {
236  friend class IndicesSyncer<T>;
237  public:
247  Iterators(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
248  BoolList& booleans);
249 
253  Iterators();
254 
258  Iterators& operator++();
259 
265  void insert(const RemoteIndex& index,
266  const std::pair<GlobalIndex,Attribute>& global);
267 
272  RemoteIndex& remoteIndex() const;
273 
278  std::pair<GlobalIndex,Attribute>& globalIndexPair() const;
279 
281 
287  bool isOld() const;
288 
298  void reset(RemoteIndexList& remoteIndices, GlobalIndexList& globalIndices,
299  BoolList& booleans);
300 
306  bool isNotAtEnd() const;
307 
313  bool isAtEnd() const;
314 
315  private:
325  IteratorTuple iterators_;
326  };
327 
329  typedef std::map<int,Iterators> IteratorsMap;
330 
342  IteratorsMap iteratorsMap_;
343 
345  void calculateMessageSizes();
346 
354  void packAndSend(int destination, char* buffer, std::size_t bufferSize, MPI_Request& req);
355 
360  template<typename T1>
361  void recvAndUnpack(T1& numberer);
362 
366  void registerMessageDatatype();
367 
371  void insertIntoRemoteIndexList(int process,
372  const std::pair<GlobalIndex,Attribute>& global,
373  char attribute);
374 
378  void resetIteratorsMap();
379 
384  bool checkReset();
385 
394  bool checkReset(const Iterators& iterators, RemoteIndexList& rlist, GlobalIndexList& gList,
395  BoolList& bList);
396  };
397 
398  template<typename TG, typename TA>
400  const std::pair<TG,TA>& i2)
401  {
402  return i1.global() < i2.first ||
403  (i1.global() == i2.first && i1.local().attribute()<i2.second);
404  }
405 
406  template<typename TG, typename TA>
407  bool operator<(const std::pair<TG,TA>& i1,
408  const IndexPair<TG,ParallelLocalIndex<TA> >& i2)
409  {
410  return i1.first < i2.global() ||
411  (i1.first == i2.global() && i1.second<i2.local().attribute());
412  }
413 
414  template<typename TG, typename TA>
416  const std::pair<TG,TA>& i2)
417  {
418  return (i1.global() == i2.first && i1.local().attribute()==i2.second);
419  }
420 
421  template<typename TG, typename TA>
423  const std::pair<TG,TA>& i2)
424  {
425  return (i1.global() != i2.first || i1.local().attribute()!=i2.second);
426  }
427 
428  template<typename TG, typename TA>
429  bool operator==(const std::pair<TG,TA>& i2,
430  const IndexPair<TG,ParallelLocalIndex<TA> >& i1)
431  {
432  return (i1.global() == i2.first && i1.local().attribute()==i2.second);
433  }
434 
435  template<typename TG, typename TA>
436  bool operator!=(const std::pair<TG,TA>& i2,
437  const IndexPair<TG,ParallelLocalIndex<TA> >& i1)
438  {
439  return (i1.global() != i2.first || i1.local().attribute()!=i2.second);
440  }
441 
458  template<typename T, typename A, typename A1>
459  void storeGlobalIndicesOfRemoteIndices(std::map<int,SLList<std::pair<typename T::GlobalIndex, typename T::LocalIndex::Attribute>,A> >& globalMap,
460  const RemoteIndices<T,A1>& remoteIndices)
461  {
462  for(auto remote = remoteIndices.begin(), end =remoteIndices.end(); remote != end; ++remote) {
463  typedef typename RemoteIndices<T,A1>::RemoteIndexList RemoteIndexList;
465  GlobalIndexList& global = globalMap[remote->first];
466  RemoteIndexList& rList = *(remote->second.first);
467 
468  for(auto index = rList.begin(), riEnd = rList.end();
469  index != riEnd; ++index) {
470  global.push_back(std::make_pair(index->localIndexPair().global(),
471  index->localIndexPair().local().attribute()));
472  }
473  }
474  }
475 
484  template<typename T, typename A, typename A1>
485  inline void repairLocalIndexPointers(std::map<int,
486  SLList<std::pair<typename T::GlobalIndex,
487  typename T::LocalIndex::Attribute>,A> >& globalMap,
488  RemoteIndices<T,A1>& remoteIndices,
489  const T& indexSet)
490  {
491  assert(globalMap.size()==static_cast<std::size_t>(remoteIndices.neighbours()));
492  // Repair pointers to index set in remote indices.
493  auto global = globalMap.begin();
494  auto end = remoteIndices.remoteIndices_.end();
495 
496  for(auto remote = remoteIndices.remoteIndices_.begin(); remote != end; ++remote, ++global) {
497  assert(remote->first==global->first);
498  assert(remote->second.first->size() == global->second.size());
499 
500  auto riEnd = remote->second.first->end();
501  auto rIndex = remote->second.first->begin();
502  auto gIndex = global->second.begin();
503  auto index = indexSet.begin();
504 
505  assert(rIndex==riEnd || gIndex != global->second.end());
506  while(rIndex != riEnd) {
507  // Search for the index in the set.
508  assert(gIndex != global->second.end());
509 
510  while(!(index->global() == gIndex->first
511  && index->local().attribute() == gIndex->second)) {
512  ++index;
513  // this is only needed for ALU, where there may exist
514  // more entries with the same global index in the remote index set
515  // than in the index set
516  if (index->global() > gIndex->first) {
517  index=indexSet.begin();
518  }
519  }
520 
521  assert(index != indexSet.end() && *index == *gIndex);
522 
523  rIndex->localIndex_ = &(*index);
524  ++index;
525  ++rIndex;
526  ++gIndex;
527  }
528  }
529  remoteIndices.sourceSeqNo_ = remoteIndices.source_->seqNo();
530  remoteIndices.destSeqNo_ = remoteIndices.target_->seqNo();
531  }
532 
533  template<typename T>
535  RemoteIndices& remoteIndices)
536  : indexSet_(indexSet), remoteIndices_(remoteIndices)
537  {
538  // index sets must match.
539  assert(remoteIndices.source_ == remoteIndices.target_);
540  assert(remoteIndices.source_ == &indexSet);
541  MPI_Comm_rank(remoteIndices_.communicator(), &rank_);
542  }
543 
544  template<typename T>
546  GlobalIndexList& globalIndices,
547  BoolList& booleans)
548  : iterators_(remoteIndices.beginModify(), globalIndices.beginModify(),
549  booleans.beginModify(), remoteIndices.end())
550  { }
551 
552  template<typename T>
554  : iterators_()
555  {}
556 
557  template<typename T>
559  {
560  ++(std::get<0>(iterators_));
561  ++(std::get<1>(iterators_));
562  ++(std::get<2>(iterators_));
563  return *this;
564  }
565 
566  template<typename T>
568  const std::pair<GlobalIndex,Attribute>& global)
569  {
570  std::get<0>(iterators_).insert(index);
571  std::get<1>(iterators_).insert(global);
572  std::get<2>(iterators_).insert(false);
573  }
574 
575  template<typename T>
576  inline typename IndicesSyncer<T>::RemoteIndex&
578  {
579  return *(std::get<0>(iterators_));
580  }
581 
582  template<typename T>
583  inline std::pair<typename IndicesSyncer<T>::GlobalIndex,typename IndicesSyncer<T>::Attribute>&
585  {
586  return *(std::get<1>(iterators_));
587  }
588 
589  template<typename T>
591  {
592  return *(std::get<2>(iterators_));
593  }
594 
595  template<typename T>
597  GlobalIndexList& globalIndices,
598  BoolList& booleans)
599  {
600  std::get<0>(iterators_) = remoteIndices.beginModify();
601  std::get<1>(iterators_) = globalIndices.beginModify();
602  std::get<2>(iterators_) = booleans.beginModify();
603  }
604 
605  template<typename T>
607  {
608  return std::get<0>(iterators_) != std::get<3>(iterators_);
609  }
610 
611  template<typename T>
613  {
614  return std::get<0>(iterators_) == std::get<3>(iterators_);
615  }
616 
617  template<typename T>
619  {
620  MPI_Datatype type[2] = {MPI_INT, MPI_INT};
621  int blocklength[2] = {1,1};
622  MPI_Aint displacement[2];
623  MPI_Aint base;
624 
625  // Compute displacement
626  MessageInformation message;
627 
628  MPI_Get_address( &(message.publish), displacement);
629  MPI_Get_address( &(message.pairs), displacement+1);
630 
631  // Make the displacement relative
632  MPI_Get_address(&message, &base);
633  displacement[0] -= base;
634  displacement[1] -= base;
635 
636  MPI_Type_create_struct( 2, blocklength, displacement, type, &datatype_);
637  MPI_Type_commit(&datatype_);
638  }
639 
640  template<typename T>
641  void IndicesSyncer<T>::calculateMessageSizes()
642  {
643  auto iEnd = indexSet_.end();
644  auto collIter = remoteIndices_.template iterator<true>();
645 
646  for(auto index = indexSet_.begin(); index != iEnd; ++index) {
647  collIter.advance(index->global(), index->local().attribute());
648  if(collIter.empty())
649  break;
650  int knownRemote=0;
651  auto end = collIter.end();
652 
653  // Count the remote indices we know.
654  for(auto valid = collIter.begin(); valid != end; ++valid) {
655  ++knownRemote;
656  }
657 
658  if(knownRemote>0) {
659  Dune::dverb<<rank_<<": publishing "<<knownRemote<<" for index "<<index->global()<< " for processes ";
660 
661  // Update MessageInformation
662  for(auto valid = collIter.begin(); valid != end; ++valid) {
663  ++(infoSend_[valid.process()].publish);
664  (infoSend_[valid.process()].pairs) += knownRemote;
665  Dune::dverb<<valid.process()<<" ";
666  Dune::dverb<<"(publish="<<infoSend_[valid.process()].publish<<", pairs="<<infoSend_[valid.process()].pairs
667  <<") ";
668  }
669  Dune::dverb<<std::endl;
670  }
671  }
672 
673  const auto end = infoSend_.end();
674 
675  // Now determine the buffersizes needed for each neighbour using MPI_Pack_size
676  MessageInformation dummy;
677 
678  auto messageIter= infoSend_.begin();
679  const auto rend = remoteIndices_.end();
680  int neighbour=0;
681 
682  for(auto remote = remoteIndices_.begin(); remote != rend; ++remote, ++neighbour) {
683  MessageInformation* message;
684  MessageInformation recv;
685 
686  if(messageIter != end && messageIter->first==remote->first) {
687  // We want to send message information to that process
688  message = const_cast<MessageInformation*>(&(messageIter->second));
689  ++messageIter;
690  }else
691  // We do not want to send information but the other process might.
692  message = &dummy;
693 
694  sendBufferSizes_[neighbour]=0;
695  int tsize;
696  // The number of indices published
697  MPI_Pack_size(1, MPI_INT,remoteIndices_.communicator(), &tsize);
698  sendBufferSizes_[neighbour] += tsize;
699 
700  for(int i=0; i < message->publish; ++i) {
701  // The global index
702  MPI_Pack_size(1, MPITraits<GlobalIndex>::getType(), remoteIndices_.communicator(), &tsize);
703  sendBufferSizes_[neighbour] += tsize;
704  // The attribute in the local index
705  MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
706  sendBufferSizes_[neighbour] += tsize;
707  // The number of corresponding remote indices
708  MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
709  sendBufferSizes_[neighbour] += tsize;
710  }
711  for(int i=0; i < message->pairs; ++i) {
712  // The process of the remote index
713  MPI_Pack_size(1, MPI_INT, remoteIndices_.communicator(), &tsize);
714  sendBufferSizes_[neighbour] += tsize;
715  // The attribute of the remote index
716  MPI_Pack_size(1, MPI_CHAR, remoteIndices_.communicator(), &tsize);
717  sendBufferSizes_[neighbour] += tsize;
718  }
719 
720  Dune::dverb<<rank_<<": Buffer (neighbour="<<remote->first<<") size is "<< sendBufferSizes_[neighbour]<<" for publish="<<message->publish<<" pairs="<<message->pairs<<std::endl;
721  }
722 
723  }
724 
725  template<typename T>
727  {
728  DefaultNumberer numberer;
729  sync(numberer);
730  }
731 
732  template<typename T>
733  template<typename T1>
734  void IndicesSyncer<T>::sync(T1& numberer)
735  {
736  // The pointers to the local indices in the remote indices
737  // will become invalid due to the resorting of the index set.
738  // Therefore store the corresponding global indices.
739  // Mark all indices as not added
740  const auto end = remoteIndices_.end();
741 
742  // Number of neighbours might change during the syncing.
743  // save the old neighbours
744  std::size_t noOldNeighbours = remoteIndices_.neighbours();
745  int* oldNeighbours = new int[noOldNeighbours];
746  sendBufferSizes_ = new std::size_t[noOldNeighbours];
747  std::size_t neighbourI = 0;
748 
749  for(auto remote = remoteIndices_.begin(); remote != end; ++remote, ++neighbourI) {
750  oldNeighbours[neighbourI] = remote->first;
751 
752  // Make sure we only have one remote index list.
753  assert(remote->second.first==remote->second.second);
754 
755  RemoteIndexList& rList = *(remote->second.first);
756 
757  // Store the corresponding global indices.
758  GlobalIndexList& global = globalMap_[remote->first];
759  BoolList& added = oldMap_[remote->first];
760  auto riEnd = rList.end();
761 
762  for(auto index = rList.begin();
763  index != riEnd; ++index) {
764  global.push_back(std::make_pair(index->localIndexPair().global(),
765  index->localIndexPair().local().attribute()));
766  added.push_back(true);
767  }
768 
769  Iterators iterators(rList, global, added);
770  iteratorsMap_.insert(std::make_pair(remote->first, iterators));
771  assert(checkReset(iteratorsMap_[remote->first], rList,global,added));
772  }
773 
774  // Exchange indices with each neighbour
775  calculateMessageSizes();
776 
777  // Allocate the buffers
778  receiveBufferSize_=1;
779  sendBuffers_ = new char*[noOldNeighbours];
780 
781  for(std::size_t i=0; i<noOldNeighbours; ++i) {
782  sendBuffers_[i] = new char[sendBufferSizes_[i]];
783  receiveBufferSize_ = std::max(receiveBufferSize_, static_cast<int>(sendBufferSizes_[i]));
784  }
785 
786  receiveBuffer_=new char[receiveBufferSize_];
787 
788  indexSet_.beginResize();
789 
790  Dune::dverb<<rank_<<": Neighbours: ";
791 
792  for(std::size_t i = 0; i<noOldNeighbours; ++i)
793  Dune::dverb<<oldNeighbours[i]<<" ";
794 
795  Dune::dverb<<std::endl;
796 
797  MPI_Request* requests = new MPI_Request[noOldNeighbours];
798  MPI_Status* statuses = new MPI_Status[noOldNeighbours];
799 
800  // Pack Message data and start the sends
801  for(std::size_t i = 0; i<noOldNeighbours; ++i)
802  packAndSend(oldNeighbours[i], sendBuffers_[i], sendBufferSizes_[i], requests[i]);
803 
804  // Probe for incoming messages, receive and unpack them
805  for(std::size_t i = 0; i<noOldNeighbours; ++i)
806  recvAndUnpack(numberer);
807  // }else{
808  // recvAndUnpack(oldNeighbours[i], numberer);
809  // packAndSend(oldNeighbours[i]);
810  // }
811  // }
812 
813  delete[] receiveBuffer_;
814 
815  // Wait for the completion of the sends
816  // Wait for completion of sends
817  if(MPI_SUCCESS!=MPI_Waitall(noOldNeighbours, requests, statuses)) {
818  std::cerr<<": MPI_Error occurred while sending message"<<std::endl;
819  for(std::size_t i=0; i< noOldNeighbours; i++)
820  if(MPI_SUCCESS!=statuses[i].MPI_ERROR)
821  std::cerr<<"Destination "<<statuses[i].MPI_SOURCE<<" error code: "<<statuses[i].MPI_ERROR<<std::endl;
822  }
823 
824  delete[] statuses;
825  delete[] requests;
826 
827  for(std::size_t i=0; i<noOldNeighbours; ++i)
828  delete[] sendBuffers_[i];
829 
830  delete[] sendBuffers_;
831  delete[] sendBufferSizes_;
832 
833  // No need for the iterator tuples any more
834  iteratorsMap_.clear();
835 
836  indexSet_.endResize();
837 
838  delete[] oldNeighbours;
839 
840  repairLocalIndexPointers(globalMap_, remoteIndices_, indexSet_);
841 
842  oldMap_.clear();
843  globalMap_.clear();
844 
845  // update the sequence number
846  remoteIndices_.sourceSeqNo_ = remoteIndices_.destSeqNo_ = indexSet_.seqNo();
847  }
848 
849  template<typename T>
850  void IndicesSyncer<T>::packAndSend(int destination, char* buffer, std::size_t bufferSize, MPI_Request& request)
851  {
852  auto iEnd = indexSet_.end();
853  int bpos = 0;
854  int published = 0;
855  int pairs = 0;
856 
857  assert(checkReset());
858 
859  // Pack the number of indices we publish
860  MPI_Pack(&(infoSend_[destination].publish), 1, MPI_INT, buffer, bufferSize, &bpos,
861  remoteIndices_.communicator());
862 
863  for(auto index = indexSet_.begin(); index != iEnd; ++index) {
864  // Search for corresponding remote indices in all iterator tuples
865  auto iteratorsEnd = iteratorsMap_.end();
866 
867  // advance all iterators to a position with global index >= index->global()
868  for(auto iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators) {
869  while(iterators->second.isNotAtEnd() &&
870  iterators->second.globalIndexPair().first < index->global())
871  ++(iterators->second);
872  assert(!iterators->second.isNotAtEnd() || iterators->second.globalIndexPair().first >= index->global());
873  }
874 
875  // Add all remote indices positioned at global which were already present before calling sync
876  // to the message.
877  // Count how many remote indices we will send
878  int indices = 0;
879  bool knownRemote = false; // Is the remote process supposed to know this index?
880 
881  for(auto iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
882  {
883  std::pair<GlobalIndex,Attribute> p;
884  if (iterators->second.isNotAtEnd())
885  {
886  p = iterators->second.globalIndexPair();
887  }
888 
889  if(iterators->second.isNotAtEnd() && iterators->second.isOld()
890  && iterators->second.globalIndexPair().first == index->global()) {
891  indices++;
892  if(destination == iterators->first)
893  knownRemote = true;
894  }
895  }
896 
897  if(!knownRemote)
898  // We do not need to send any indices
899  continue;
900 
901  Dune::dverb<<rank_<<": sending "<<indices<<" for index "<<index->global()<<" to "<<destination<<std::endl;
902 
903 
904  // Pack the global index, the attribute and the number
905  MPI_Pack(const_cast<GlobalIndex*>(&(index->global())), 1, MPITraits<GlobalIndex>::getType(), buffer, bufferSize, &bpos,
906  remoteIndices_.communicator());
907 
908  char attr = index->local().attribute();
909  MPI_Pack(&attr, 1, MPI_CHAR, buffer, bufferSize, &bpos,
910  remoteIndices_.communicator());
911 
912  // Pack the number of remote indices we send.
913  MPI_Pack(&indices, 1, MPI_INT, buffer, bufferSize, &bpos,
914  remoteIndices_.communicator());
915 
916  // Pack the information about the remote indices
917  for(auto iterators = iteratorsMap_.begin(); iteratorsEnd != iterators; ++iterators)
918  if(iterators->second.isNotAtEnd() && iterators->second.isOld()
919  && iterators->second.globalIndexPair().first == index->global()) {
920  int process = iterators->first;
921 
922  ++pairs;
923  assert(pairs <= infoSend_[destination].pairs);
924  MPI_Pack(&process, 1, MPI_INT, buffer, bufferSize, &bpos,
925  remoteIndices_.communicator());
926  char attr2 = iterators->second.remoteIndex().attribute();
927 
928  MPI_Pack(&attr2, 1, MPI_CHAR, buffer, bufferSize, &bpos,
929  remoteIndices_.communicator());
930  --indices;
931  }
932  assert(indices==0);
933  ++published;
934  Dune::dvverb<<" (publish="<<published<<", pairs="<<pairs<<")"<<std::endl;
935  assert(published <= infoSend_[destination].publish);
936  }
937 
938  // Make sure we send all expected entries
939  assert(published == infoSend_[destination].publish);
940  assert(pairs == infoSend_[destination].pairs);
941  resetIteratorsMap();
942 
943  Dune::dverb << rank_<<": Sending message of "<<bpos<<" bytes to "<<destination<<std::endl;
944 
945  MPI_Issend(buffer, bpos, MPI_PACKED, destination, 345, remoteIndices_.communicator(),&request);
946  }
947 
948  template<typename T>
949  inline void IndicesSyncer<T>::insertIntoRemoteIndexList(int process,
950  const std::pair<GlobalIndex,Attribute>& globalPair,
951  char attribute)
952  {
953  Dune::dverb<<"Inserting from "<<process<<" "<<globalPair.first<<", "<<
954  globalPair.second<<" "<<attribute<<std::endl;
955 
956  resetIteratorsMap();
957 
958  // There might be cases where there no remote indices for that process yet
959  typename IteratorsMap::iterator found = iteratorsMap_.find(process);
960 
961  if( found == iteratorsMap_.end() ) {
962  Dune::dverb<<"Discovered new neighbour "<<process<<std::endl;
963  RemoteIndexList* rlist = new RemoteIndexList();
964  remoteIndices_.remoteIndices_.insert(std::make_pair(process,std::make_pair(rlist,rlist)));
965  Iterators iterators = Iterators(*rlist, globalMap_[process], oldMap_[process]);
966  found = iteratorsMap_.insert(std::make_pair(process, iterators)).first;
967  }
968 
969  Iterators& iterators = found->second;
970 
971  // Search for the remote index
972  while(iterators.isNotAtEnd() && iterators.globalIndexPair() < globalPair) {
973  // Increment all iterators
974  ++iterators;
975 
976  }
977 
978  if(iterators.isAtEnd() || iterators.globalIndexPair() != globalPair) {
979  // The entry is not yet known
980  // Insert in the list and do not change the first iterator.
981  iterators.insert(RemoteIndex(Attribute(attribute)),globalPair);
982  return;
983  }
984 
985  // Global indices match
986  bool indexIsThere=false;
987  for(Iterators tmpIterators = iterators;
988  !tmpIterators.isAtEnd() && tmpIterators.globalIndexPair() == globalPair;
989  ++tmpIterators)
990  //entry already exists with the same attribute
991  if(tmpIterators.globalIndexPair().second == attribute) {
992  indexIsThere=true;
993  break;
994  }
995 
996  if(!indexIsThere)
997  // The entry is not yet known
998  // Insert in the list and do not change the first iterator.
999  iterators.insert(RemoteIndex(Attribute(attribute)),globalPair);
1000  }
1001 
1002  template<typename T>
1003  template<typename T1>
1004  void IndicesSyncer<T>::recvAndUnpack(T1& numberer)
1005  {
1006  const ParallelIndexSet& constIndexSet = indexSet_;
1007  auto iEnd = constIndexSet.end();
1008  auto index = constIndexSet.begin();
1009  int bpos = 0;
1010  int publish;
1011 
1012  assert(checkReset());
1013 
1014  MPI_Status status;
1015 
1016  // We have to determine the message size and source before the receive
1017  MPI_Probe(MPI_ANY_SOURCE, 345, remoteIndices_.communicator(), &status);
1018 
1019  int source=status.MPI_SOURCE;
1020  int count;
1021  MPI_Get_count(&status, MPI_PACKED, &count);
1022 
1023  Dune::dvverb<<rank_<<": Receiving message from "<< source<<" with "<<count<<" bytes"<<std::endl;
1024 
1025  if(count>receiveBufferSize_) {
1026  receiveBufferSize_=count;
1027  delete[] receiveBuffer_;
1028  receiveBuffer_ = new char[receiveBufferSize_];
1029  }
1030 
1031  MPI_Recv(receiveBuffer_, count, MPI_PACKED, source, 345, remoteIndices_.communicator(), &status);
1032 
1033  // How many global entries were published?
1034  MPI_Unpack(receiveBuffer_, count, &bpos, &publish, 1, MPI_INT, remoteIndices_.communicator());
1035 
1036  // Now unpack the remote indices and add them.
1037  while(publish>0) {
1038 
1039  // Unpack information about the local index on the source process
1040  GlobalIndex global; // global index of the current entry
1041  char sourceAttribute; // Attribute on the source process
1042  int pairs;
1043 
1044  MPI_Unpack(receiveBuffer_, count, &bpos, &global, 1, MPITraits<GlobalIndex>::getType(),
1045  remoteIndices_.communicator());
1046  MPI_Unpack(receiveBuffer_, count, &bpos, &sourceAttribute, 1, MPI_CHAR,
1047  remoteIndices_.communicator());
1048  MPI_Unpack(receiveBuffer_, count, &bpos, &pairs, 1, MPI_INT,
1049  remoteIndices_.communicator());
1050 
1051  // Insert the entry on the remote process to our
1052  // remote index list
1053  SLList<std::pair<int,Attribute> > sourceAttributeList;
1054  sourceAttributeList.push_back(std::make_pair(source,Attribute(sourceAttribute)));
1055 #ifndef NDEBUG
1056  bool foundSelf = false;
1057 #endif
1058  Attribute myAttribute=Attribute();
1059 
1060  // Unpack the remote indices
1061  for(; pairs>0; --pairs) {
1062  // Unpack the process id that knows the index
1063  int process;
1064  char attribute;
1065  MPI_Unpack(receiveBuffer_, count, &bpos, &process, 1, MPI_INT,
1066  remoteIndices_.communicator());
1067  // Unpack the attribute
1068  MPI_Unpack(receiveBuffer_, count, &bpos, &attribute, 1, MPI_CHAR,
1069  remoteIndices_.communicator());
1070 
1071  if(process==rank_) {
1072 #ifndef NDEBUG
1073  foundSelf=true;
1074 #endif
1075  myAttribute=Attribute(attribute);
1076  // Now we know the local attribute of the global index
1077  //Only add the index if it is unknown.
1078  // Do we know that global index already?
1079  auto pos = std::lower_bound(index, iEnd, IndexPair(global));
1080 
1081  if(pos == iEnd || pos->global() != global) {
1082  // no entry with this global index
1083  indexSet_.add(global,
1084  ParallelLocalIndex<Attribute>(numberer(global),
1085  myAttribute, true));
1086  Dune::dvverb << "Adding "<<global<<" "<<myAttribute<<std::endl;
1087  continue;
1088  }
1089 
1090  // because of above the global indices match. Add only if the attribute is different
1091  bool indexIsThere = false;
1092  index=pos;
1093 
1094  for(; pos->global()==global; ++pos)
1095  if(pos->local().attribute() == myAttribute) {
1096  Dune::dvverb<<"found "<<global<<" "<<myAttribute<<std::endl;
1097  indexIsThere = true;
1098  break;
1099  }
1100 
1101  if(!indexIsThere) {
1102  indexSet_.add(global,
1103  ParallelLocalIndex<Attribute>(numberer(global),
1104  myAttribute, true));
1105  Dune::dvverb << "Adding "<<global<<" "<<myAttribute<<std::endl;
1106  }
1107 
1108  }else{
1109  sourceAttributeList.push_back(std::make_pair(process,Attribute(attribute)));
1110  }
1111  }
1112  assert(foundSelf);
1113  // Insert remote indices
1114  typedef typename SLList<std::pair<int,Attribute> >::const_iterator Iter;
1115  for(Iter i=sourceAttributeList.begin(), end=sourceAttributeList.end();
1116  i!=end; ++i)
1117  insertIntoRemoteIndexList(i->first, std::make_pair(global, myAttribute),
1118  i->second);
1119  --publish;
1120  }
1121 
1122  resetIteratorsMap();
1123  }
1124 
1125  template<typename T>
1126  void IndicesSyncer<T>::resetIteratorsMap(){
1127 
1128  // Reset iterators in all tuples.
1129  const auto remoteEnd = remoteIndices_.remoteIndices_.end();
1130  auto iterators = iteratorsMap_.begin();
1131  auto global = globalMap_.begin();
1132  auto added = oldMap_.begin();
1133 
1134  for(auto remote = remoteIndices_.remoteIndices_.begin();
1135  remote != remoteEnd; ++remote, ++global, ++added, ++iterators) {
1136  iterators->second.reset(*(remote->second.first), global->second, added->second);
1137  }
1138  }
1139 
1140  template<typename T>
1141  bool IndicesSyncer<T>::checkReset(const Iterators& iterators, RemoteIndexList& rList, GlobalIndexList& gList,
1142  BoolList& bList){
1143 
1144  if(std::get<0>(iterators.iterators_) != rList.begin())
1145  return false;
1146  if(std::get<1>(iterators.iterators_) != gList.begin())
1147  return false;
1148  if(std::get<2>(iterators.iterators_) != bList.begin())
1149  return false;
1150  return true;
1151  }
1152 
1153 
1154  template<typename T>
1155  bool IndicesSyncer<T>::checkReset(){
1156 
1157  // Reset iterators in all tuples.
1158  const auto remoteEnd = remoteIndices_.remoteIndices_.end();
1159  auto iterators = iteratorsMap_.begin();
1160  auto global = globalMap_.begin();
1161  auto added = oldMap_.begin();
1162  bool ret = true;
1163 
1164  for(auto remote = remoteIndices_.remoteIndices_.begin();
1165  remote != remoteEnd; ++remote, ++global, ++added, ++iterators) {
1166  if(!checkReset(iterators->second, *(remote->second.first), global->second,
1167  added->second))
1168  ret=false;
1169  }
1170  return ret;
1171  }
1172 }
1173 
1174 #endif
1175 #endif
Provides a map between global and local indices.
Classes describing a distributed indexset.
Implements a singly linked list together with the necessary iterators.
Standard Dune debug streams.
void push_back(const MemberType &item)
Add a new entry to the end of the list.
Definition: sllist.hh:649
SLListIterator< T, A > iterator
The mutable iterator of the list.
Definition: sllist.hh:67
iterator end()
Get an iterator pointing to the end of the list.
Definition: sllist.hh:780
ModifyIterator beginModify()
Get an iterator capable of deleting and inserting elements.
Definition: sllist.hh:793
SLListConstIterator< T, A > const_iterator
The constant iterator of the list.
Definition: sllist.hh:72
SLListModifyIterator< T, A > ModifyIterator
The type of the iterator capable of deletion and insertion.
Definition: sllist.hh:101
iterator begin()
Get an iterator pointing to the first element in the list.
Definition: sllist.hh:768
EnableIfInterOperable< T1, T2, bool >::type operator!=(const ForwardIteratorFacade< T1, V1, R1, D > &lhs, const ForwardIteratorFacade< T2, V2, R2, D > &rhs)
Checks for inequality.
Definition: iteratorfacades.hh:257
EnableIfInterOperable< T1, T2, bool >::type operator<(const RandomAccessIteratorFacade< T1, V1, R1, D > &lhs, const RandomAccessIteratorFacade< T2, V2, R2, D > &rhs)
Comparison operator.
Definition: iteratorfacades.hh:635
EnableIfInterOperable< T1, T2, bool >::type operator==(const ForwardIteratorFacade< T1, V1, R1, D > &lhs, const ForwardIteratorFacade< T2, V2, R2, D > &rhs)
Checks for equality.
Definition: iteratorfacades.hh:235
Dune::RemoteIndices< ParallelIndexSet > RemoteIndices
Type of the remote indices.
Definition: indicessyncer.hh:57
int publish
The number of indices we publish for the other process.
Definition: indicessyncer.hh:124
void repairLocalIndexPointers(std::map< int, SLList< std::pair< typename T::GlobalIndex, typename T::LocalIndex::Attribute >, A > > &globalMap, RemoteIndices< T, A1 > &remoteIndices, const T &indexSet)
Repair the pointers to the local indices in the remote indices.
Definition: indicessyncer.hh:485
ParallelIndexSet::GlobalIndex GlobalIndex
Type of the global index used in the index set.
Definition: indicessyncer.hh:49
int pairs
The number of pairs (attribute and process number) we publish to the neighbour process.
Definition: indicessyncer.hh:129
ParallelIndexSet::LocalIndex::Attribute Attribute
Type of the attribute used in the index set.
Definition: indicessyncer.hh:52
IndicesSyncer(ParallelIndexSet &indexSet, RemoteIndices &remoteIndices)
Constructor.
Definition: indicessyncer.hh:534
T ParallelIndexSet
The type of the index set.
Definition: indicessyncer.hh:43
Attribute & attribute() const
TG GlobalIndex
the type of the global index. This type has to provide at least a operator< for sorting.
Definition: indexset.hh:224
void storeGlobalIndicesOfRemoteIndices(std::map< int, SLList< std::pair< typename T::GlobalIndex, typename T::LocalIndex::Attribute >, A > > &globalMap, const RemoteIndices< T, A1 > &remoteIndices)
Stores the corresponding global indices of the remote index information.
Definition: indicessyncer.hh:459
void sync()
Sync the index set.
Definition: indicessyncer.hh:726
int seqNo() const
Get the internal sequence number.
ParallelIndexSet::IndexPair IndexPair
The type of the index pair.
Definition: indicessyncer.hh:46
MessageInformation()
Definition: indicessyncer.hh:120
DVVerbType dvverb(std::cout)
stream for very verbose output.
Definition: stdstreams.hh:93
DVerbType dverb(std::cout)
Singleton of verbose debug stream.
Definition: stdstreams.hh:114
Dune namespace.
Definition: alignedallocator.hh:11
auto max(const AlignedNumber< T, align > &a, const AlignedNumber< T, align > &b)
Definition: debugalign.hh:412
static MPI_Datatype getType()
Definition: mpitraits.hh:46
A pair consisting of a global and local index.
Definition: indexset.hh:83
Class for recomputing missing indices of a distributed index set.
Definition: indicessyncer.hh:39
An index present on the local process with an additional attribute flag.
Definition: plocalindex.hh:47
The indices present on remote processes.
Definition: remoteindices.hh:187
MPI_Comm communicator() const
Get the mpi communicator used.
Definition: remoteindices.hh:1696
const_iterator end() const
Get an iterator over all remote index lists.
Definition: remoteindices.hh:1529
int neighbours() const
Get the number of processors we share indices with.
Definition: remoteindices.hh:1446
typename std::allocator_traits< A >::template rebind_alloc< RemoteIndex > Allocator
The type of the allocator for the remote index list.
Definition: remoteindices.hh:235
const_iterator begin() const
Get an iterator over all remote index lists.
Definition: remoteindices.hh:1522
Information about an index residing on another processor.
Definition: remoteindices.hh:71
A mutable iterator for the SLList.
Definition: sllist.hh:269
A single linked list.
Definition: sllist.hh:42