3 #ifndef DUNE_REMOTEINDICES_HH
4 #define DUNE_REMOTEINDICES_HH
37 template<
typename TG,
typename TA>
39 class MPITraits<IndexPair<TG,ParallelLocalIndex<TA> > >
42 inline static MPI_Datatype getType();
44 static MPI_Datatype type;
48 template<
typename T,
typename A>
51 template<
typename T1,
typename T2>
58 template<
typename T1,
typename T2>
59 std::ostream&
operator<<(std::ostream& os,
const RemoteIndex<T1,T2>& index);
62 template<
typename T,
typename A,
bool mode>
63 class RemoteIndexListModifier;
69 template<
typename T1,
typename T2>
73 friend class IndicesSyncer;
75 template<
typename T,
typename A,
typename A1>
76 friend void repairLocalIndexPointers(std::map<
int,SLList<std::pair<typename T::GlobalIndex, typename T::LocalIndex::Attribute>,A> >&,
80 template<
typename T,
typename A,
bool mode>
81 friend class RemoteIndexListModifier;
88 typedef T1 GlobalIndex;
102 typedef IndexPair<GlobalIndex,ParallelLocalIndex<Attribute> >
109 const Attribute attribute()
const;
116 const PairType& localIndexPair()
const;
129 RemoteIndex(
const T2& attribute,
130 const PairType* local);
138 RemoteIndex(
const T2& attribute);
145 const PairType* localIndex_;
151 template<
class T,
class A>
152 std::ostream&
operator<<(std::ostream& os,
const RemoteIndices<T,A>& indices);
154 class InterfaceBuilder;
156 template<
class T,
class A>
157 class CollectiveIterator;
164 template<
typename T1,
typename T2>
165 class OwnerOverlapCopyCommunication;
184 template<
class T,
class A=std::allocator<RemoteIndex<
typename T::GlobalIndex,
185 typename T::LocalIndex::Attribute> > >
188 friend class InterfaceBuilder;
189 friend class IndicesSyncer<T>;
190 template<
typename T1,
typename A2,
typename A1>
191 friend void repairLocalIndexPointers(std::map<
int,SLList<std::pair<typename T1::GlobalIndex, typename T1::LocalIndex::Attribute>,A2> >&,
192 RemoteIndices<T1,A1>&,
195 template<
class G,
class T1,
class T2>
196 friend void fillIndexSetHoles(
const G& graph, Dune::OwnerOverlapCopyCommunication<T1,T2>& oocomm);
197 friend std::ostream& operator<<<>(std::ostream&,
const RemoteIndices<T>&);
204 typedef T ParallelIndexSet;
208 typedef CollectiveIterator<T,A> CollectiveIteratorT;
224 typedef typename LocalIndex::Attribute Attribute;
229 typedef Dune::RemoteIndex<GlobalIndex,Attribute> RemoteIndex;
235 typedef typename A::template rebind<RemoteIndex>::other Allocator;
242 typedef std::map<int, std::pair<RemoteIndexList*,RemoteIndexList*> >
245 typedef typename RemoteIndexMap::const_iterator const_iterator;
264 inline RemoteIndices(
const ParallelIndexSet& source,
const ParallelIndexSet& destination,
265 const MPI_Comm& comm,
const std::vector<int>& neighbours=std::vector<int>(),
bool includeSelf=
false);
276 void setIncludeSelf(
bool includeSelf);
294 void setIndexSets(
const ParallelIndexSet& source,
const ParallelIndexSet& destination,
295 const MPI_Comm& comm,
const std::vector<int>& neighbours=std::vector<int>());
298 void setNeighbours(
const C& neighbours)
300 neighbourIds.clear();
301 neighbourIds.insert(neighbours.begin(), neighbours.end());
305 const std::set<int>& getNeighbours()
const
324 template<
bool ignorePublic>
336 inline bool isSynced()
const;
341 inline MPI_Comm communicator()
const;
357 template<
bool mode,
bool send>
358 inline RemoteIndexListModifier<T,A,mode> getModifier(
int process);
366 inline const_iterator find(
int proc)
const;
372 inline const_iterator begin()
const;
378 inline const_iterator end()
const;
384 inline CollectiveIteratorT iterator()
const;
395 inline int neighbours()
const;
398 inline const ParallelIndexSet& sourceIndexSet()
const;
401 inline const ParallelIndexSet& destinationIndexSet()
const;
405 RemoteIndices(
const RemoteIndices&) =
delete;
408 const ParallelIndexSet* source_;
411 const ParallelIndexSet* target_;
418 std::set<int> neighbourIds;
421 const static int commTag_=333;
454 typedef IndexPair<GlobalIndex, LocalIndex>
463 RemoteIndexMap remoteIndices_;
475 template<
bool ignorePublic>
476 inline void buildRemote(
bool includeSelf);
483 inline int noPublic(
const ParallelIndexSet& indexSet);
496 template<
bool ignorePublic>
497 inline void packEntries(PairType** myPairs,
const ParallelIndexSet& indexSet,
498 char* p_out, MPI_Datatype type,
int bufferSize,
499 int* position,
int n);
514 inline void unpackIndices(RemoteIndexList& remote,
int remoteEntries,
515 PairType** local,
int localEntries,
char* p_in,
516 MPI_Datatype type,
int* position,
int bufferSize,
519 inline void unpackIndices(RemoteIndexList& send, RemoteIndexList& receive,
520 int remoteEntries, PairType** localSource,
521 int localSourceEntries, PairType** localDest,
522 int localDestEntries,
char* p_in,
523 MPI_Datatype type,
int* position,
int bufferSize);
525 void unpackCreateRemote(
char* p_in, PairType** sourcePairs, PairType** DestPairs,
526 int remoteProc,
int sourcePublish,
int destPublish,
527 int bufferSize,
bool sendTwo,
bool fromOurSelf=
false);
547 template<
class T,
class A,
bool mode>
548 class RemoteIndexListModifier
551 template<
typename T1,
typename A1>
552 friend class RemoteIndices;
555 class InvalidPosition :
public RangeError
573 typedef T ParallelIndexSet;
588 typedef typename LocalIndex::Attribute Attribute;
593 typedef Dune::RemoteIndex<GlobalIndex,Attribute> RemoteIndex;
607 typedef SLListModifyIterator<RemoteIndex,Allocator> ModifyIterator;
612 typedef typename RemoteIndexList::const_iterator ConstIterator;
627 void insert(
const RemoteIndex& index);
644 void insert(
const RemoteIndex& index,
const GlobalIndex& global);
653 bool remove(
const GlobalIndex& global);
667 void repairLocalIndexPointers();
670 RemoteIndexListModifier(
const RemoteIndexListModifier&);
676 RemoteIndexListModifier()
687 RemoteIndexListModifier(
const ParallelIndexSet& indexSet,
688 RemoteIndexList& rList);
690 typedef SLList<GlobalIndex,Allocator> GlobalList;
691 typedef typename GlobalList::ModifyIterator GlobalModifyIterator;
692 RemoteIndexList* rList_;
693 const ParallelIndexSet* indexSet_;
695 ModifyIterator iter_;
696 GlobalModifyIterator giter_;
706 template<
class T,
class A>
707 class CollectiveIterator
713 typedef T ParallelIndexSet;
728 typedef typename LocalIndex::Attribute Attribute;
731 typedef Dune::RemoteIndex<GlobalIndex,Attribute> RemoteIndex;
734 typedef typename A::template rebind<RemoteIndex>::other Allocator;
740 typedef std::map<int,std::pair<
typename RemoteIndexList::const_iterator,
741 const typename RemoteIndexList::const_iterator> >
747 typedef std::map<int, std::pair<RemoteIndexList*,RemoteIndexList*> >
755 inline CollectiveIterator(
const RemoteIndexMap& map_,
bool send);
765 inline void advance(
const GlobalIndex& global);
776 inline void advance(
const GlobalIndex& global,
const Attribute& attribute);
778 CollectiveIterator& operator++();
794 typedef typename Map::iterator RealIterator;
795 typedef typename Map::iterator ConstRealIterator;
799 iterator(
const RealIterator& iter,
const ConstRealIterator& end, GlobalIndex& index)
800 : iter_(iter), end_(end), index_(index), hasAttribute(false)
803 while(iter_!=end_ && iter_->second.first->localIndexPair().global()!=index_)
807 iterator(
const RealIterator& iter,
const ConstRealIterator& end, GlobalIndex index,
809 : iter_(iter), end_(end), index_(index), attribute_(attribute), hasAttribute(true)
812 while(iter_!=end_ && (iter_->second.first->localIndexPair().global()!=index_
813 || iter_->second.first->localIndexPair().local().attribute()!=attribute))
817 iterator(
const iterator& other)
818 : iter_(other.iter_), end_(other.end_), index_(other.index_)
822 iterator& operator++()
826 while(iter_!=end_ && (iter_->second.first->localIndexPair().global()!=index_ ||
828 iter_->second.first->localIndexPair().local().attribute()!=attribute_)))
830 assert(iter_==end_ ||
831 (iter_->second.first->localIndexPair().global()==index_));
832 assert(iter_==end_ || !hasAttribute ||
833 (iter_->second.first->localIndexPair().local().attribute()==attribute_));
840 return *(iter_->second.first);
850 const RemoteIndex* operator->()
const
852 return iter_->second.first.operator->();
858 return other.iter_==iter_;
864 return other.iter_!=iter_;
873 Attribute attribute_;
885 Attribute attribute_;
889 template<
typename TG,
typename TA>
890 MPI_Datatype MPITraits<IndexPair<TG,ParallelLocalIndex<TA> > >::getType()
892 if(type==MPI_DATATYPE_NULL) {
893 int length[2] = {1, 1};
896 MPI_Datatype types[2] = {MPITraits<TG>::getType(),
897 MPITraits<ParallelLocalIndex<TA> >::getType()};
898 IndexPair<TG,ParallelLocalIndex<TA> > rep;
899 MPI_Get_address(&rep, &base);
900 MPI_Get_address(&(rep.global_), &disp[0]);
901 MPI_Get_address(&(rep.local_), &disp[1]);
902 for (MPI_Aint& d : disp)
906 MPI_Type_create_struct(2, length, disp, types, &tmp);
908 MPI_Type_create_resized(tmp, 0,
sizeof(IndexPair<TG,ParallelLocalIndex<TA> >), &type);
909 MPI_Type_commit(&type);
916 template<
typename TG,
typename TA>
917 MPI_Datatype MPITraits<IndexPair<TG,ParallelLocalIndex<TA> > >::type=MPI_DATATYPE_NULL;
919 template<
typename T1,
typename T2>
920 RemoteIndex<T1,T2>::RemoteIndex(
const T2& attribute,
const PairType* local)
921 : localIndex_(local), attribute_(attribute)
924 template<
typename T1,
typename T2>
925 RemoteIndex<T1,T2>::RemoteIndex(
const T2& attribute)
926 : localIndex_(0), attribute_(attribute)
929 template<
typename T1,
typename T2>
930 RemoteIndex<T1,T2>::RemoteIndex()
931 : localIndex_(0), attribute_()
933 template<
typename T1,
typename T2>
936 return localIndex_==ri.localIndex_ && attribute_==ri.attribute;
939 template<
typename T1,
typename T2>
942 return localIndex_!=ri.localIndex_ || attribute_!=ri.attribute_;
945 template<
typename T1,
typename T2>
946 inline const T2 RemoteIndex<T1,T2>::attribute()
const
948 return T2(attribute_);
951 template<
typename T1,
typename T2>
952 inline const IndexPair<T1,ParallelLocalIndex<T2> >& RemoteIndex<T1,T2>::localIndexPair()
const
957 template<
typename T,
typename A>
958 inline RemoteIndices<T,A>::RemoteIndices(
const ParallelIndexSet& source,
959 const ParallelIndexSet& destination,
960 const MPI_Comm& comm,
961 const std::vector<int>& neighbours,
963 : source_(&source), target_(&destination), comm_(comm),
964 sourceSeqNo_(-1), destSeqNo_(-1), publicIgnored(false), firstBuild(true),
965 includeSelf(includeSelf_)
967 setNeighbours(neighbours);
970 template<
typename T,
typename A>
971 void RemoteIndices<T,A>::setIncludeSelf(
bool b)
976 template<
typename T,
typename A>
977 RemoteIndices<T,A>::RemoteIndices()
978 : source_(0), target_(0), sourceSeqNo_(-1),
979 destSeqNo_(-1), publicIgnored(false), firstBuild(true),
983 template<
class T,
typename A>
984 void RemoteIndices<T,A>::setIndexSets(
const ParallelIndexSet& source,
985 const ParallelIndexSet& destination,
986 const MPI_Comm& comm,
987 const std::vector<int>& neighbours)
991 target_ = &destination;
994 setNeighbours(neighbours);
997 template<
typename T,
typename A>
998 const typename RemoteIndices<T,A>::ParallelIndexSet&
999 RemoteIndices<T,A>::sourceIndexSet()
const
1005 template<
typename T,
typename A>
1006 const typename RemoteIndices<T,A>::ParallelIndexSet&
1007 RemoteIndices<T,A>::destinationIndexSet()
const
1013 template<
typename T,
typename A>
1014 RemoteIndices<T,A>::~RemoteIndices()
1019 template<
typename T,
typename A>
1020 template<
bool ignorePublic>
1021 inline void RemoteIndices<T,A>::packEntries(IndexPair<GlobalIndex,LocalIndex>** pairs,
1022 const ParallelIndexSet& indexSet,
1023 char* p_out, MPI_Datatype type,
1025 int *position,
int n)
1029 typedef typename ParallelIndexSet::const_iterator const_iterator;
1030 typedef IndexPair<GlobalIndex,LocalIndex> PairType;
1031 const const_iterator end = indexSet.end();
1035 for(const_iterator index = indexSet.begin(); index != end; ++index)
1036 if(ignorePublic || index->local().isPublic()) {
1038 MPI_Pack(
const_cast<PairType*
>(&(*index)), 1,
1040 p_out, bufferSize, position, comm_);
1041 pairs[i++] =
const_cast<PairType*
>(&(*index));
1047 template<
typename T,
typename A>
1048 inline int RemoteIndices<T,A>::noPublic(
const ParallelIndexSet& indexSet)
1050 typedef typename ParallelIndexSet::const_iterator const_iterator;
1054 const const_iterator end=indexSet.end();
1055 for(const_iterator index=indexSet.begin(); index!=end; ++index)
1056 if(index->local().isPublic())
1064 template<
typename T,
typename A>
1065 inline void RemoteIndices<T,A>::unpackCreateRemote(
char* p_in, PairType** sourcePairs,
1066 PairType** destPairs,
int remoteProc,
1067 int sourcePublish,
int destPublish,
1068 int bufferSize,
bool sendTwo,
1073 int noRemoteSource=-1, noRemoteDest=-1;
1074 char twoIndexSets=0;
1077 MPI_Unpack(p_in, bufferSize, &position, &twoIndexSets, 1, MPI_CHAR, comm_);
1079 MPI_Unpack(p_in, bufferSize, &position, &noRemoteSource, 1, MPI_INT, comm_);
1081 MPI_Unpack(p_in, bufferSize, &position, &noRemoteDest, 1, MPI_INT, comm_);
1085 RemoteIndexList* receive=
new RemoteIndexList();
1087 RemoteIndexList* send=0;
1089 MPI_Datatype type= MPITraits<PairType>::getType();
1093 send =
new RemoteIndexList();
1095 unpackIndices(*send, *receive, noRemoteSource, sourcePairs, sourcePublish,
1096 destPairs, destPublish, p_in, type, &position, bufferSize);
1099 unpackIndices(*receive, noRemoteSource, sourcePairs, sourcePublish,
1100 p_in, type, &position, bufferSize, fromOurSelf);
1105 int oldPos=position;
1107 unpackIndices(*receive, noRemoteSource, destPairs, destPublish,
1108 p_in, type, &position, bufferSize, fromOurSelf);
1113 send =
new RemoteIndexList();
1114 unpackIndices(*send, noRemoteDest, sourcePairs, sourcePublish,
1115 p_in, type, &position, bufferSize, fromOurSelf);
1118 if(receive->empty() && send->empty()) {
1126 remoteIndices_.insert(std::make_pair(remoteProc,
1127 std::make_pair(send,receive)));
1132 template<
typename T,
typename A>
1133 template<
bool ignorePublic>
1134 inline void RemoteIndices<T,A>::buildRemote(
bool includeSelf_)
1138 MPI_Comm_rank(comm_, &rank);
1139 MPI_Comm_size(comm_, &procs);
1143 int sourcePublish, destPublish;
1146 char sendTwo = (source_ != target_);
1148 if(procs==1 && !(sendTwo || includeSelf_))
1152 sourcePublish = (ignorePublic) ? source_->size() : noPublic(*source_);
1155 destPublish = (ignorePublic) ? target_->size() : noPublic(*target_);
1160 int maxPublish, publish=sourcePublish+destPublish;
1163 MPI_Allreduce(&publish, &maxPublish, 1, MPI_INT, MPI_MAX, comm_);
1166 typedef IndexPair<GlobalIndex,LocalIndex> PairType;
1168 PairType** destPairs;
1169 PairType** sourcePairs =
new PairType*[sourcePublish>0 ? sourcePublish : 1];
1172 destPairs =
new PairType*[destPublish>0 ? destPublish : 1];
1174 destPairs=sourcePairs;
1176 char** buffer =
new char*[2];
1183 MPI_Datatype type = MPITraits<PairType>::getType();
1185 MPI_Pack_size(maxPublish, type, comm_,
1187 MPI_Pack_size(1, MPI_INT, comm_,
1189 MPI_Pack_size(1, MPI_CHAR, comm_,
1195 bufferSize += 2 * intSize + charSize;
1197 if(bufferSize<=0) bufferSize=1;
1199 buffer[0] =
new char[bufferSize];
1200 buffer[1] =
new char[bufferSize];
1204 MPI_Pack(&sendTwo, 1, MPI_CHAR, buffer[0], bufferSize, &position,
1208 MPI_Pack(&sourcePublish, 1, MPI_INT, buffer[0], bufferSize, &position,
1210 MPI_Pack(&destPublish, 1, MPI_INT, buffer[0], bufferSize, &position,
1214 packEntries<ignorePublic>(sourcePairs, *source_, buffer[0], type,
1215 bufferSize, &position, sourcePublish);
1218 packEntries<ignorePublic>(destPairs, *target_, buffer[0], type,
1219 bufferSize, &position, destPublish);
1223 if(sendTwo|| includeSelf_)
1224 unpackCreateRemote(buffer[0], sourcePairs, destPairs, rank, sourcePublish,
1225 destPublish, bufferSize, sendTwo, includeSelf_);
1227 neighbourIds.erase(rank);
1229 if(neighbourIds.size()==0)
1231 Dune::dvverb<<rank<<
": Sending messages in a ring"<<std::endl;
1233 for(
int proc=1; proc<procs; proc++) {
1235 char* p_out = buffer[1-(proc%2)];
1236 char* p_in = buffer[proc%2];
1240 MPI_Ssend(p_out, bufferSize, MPI_PACKED, (rank+1)%procs,
1242 MPI_Recv(p_in, bufferSize, MPI_PACKED, (rank+procs-1)%procs,
1243 commTag_, comm_, &status);
1245 MPI_Recv(p_in, bufferSize, MPI_PACKED, (rank+procs-1)%procs,
1246 commTag_, comm_, &status);
1247 MPI_Ssend(p_out, bufferSize, MPI_PACKED, (rank+1)%procs,
1253 int remoteProc = (rank+procs-proc)%procs;
1255 unpackCreateRemote(p_in, sourcePairs, destPairs, remoteProc, sourcePublish,
1256 destPublish, bufferSize, sendTwo);
1263 MPI_Request* requests=
new MPI_Request[neighbourIds.size()];
1264 MPI_Request* req=requests;
1266 typedef typename std::set<int>::size_type size_type;
1267 size_type noNeighbours=neighbourIds.size();
1270 for(std::set<int>::iterator neighbour=neighbourIds.begin();
1271 neighbour!= neighbourIds.end(); ++neighbour) {
1273 MPI_Issend(buffer[0], position , MPI_PACKED, *neighbour, commTag_, comm_, req++);
1278 for(size_type received=0; received <noNeighbours; ++received)
1282 MPI_Probe(MPI_ANY_SOURCE, commTag_, comm_, &status);
1283 int remoteProc=status.MPI_SOURCE;
1285 MPI_Get_count(&status, MPI_PACKED, &size);
1287 MPI_Recv(buffer[1], size, MPI_PACKED, remoteProc,
1288 commTag_, comm_, &status);
1290 unpackCreateRemote(buffer[1], sourcePairs, destPairs, remoteProc, sourcePublish,
1291 destPublish, bufferSize, sendTwo);
1294 MPI_Status* statuses =
new MPI_Status[neighbourIds.size()];
1296 if(MPI_ERR_IN_STATUS==MPI_Waitall(neighbourIds.size(), requests, statuses)) {
1297 for(size_type i=0; i < neighbourIds.size(); ++i)
1298 if(statuses[i].MPI_ERROR!=MPI_SUCCESS) {
1299 std::cerr<<rank<<
": MPI_Error occurred while receiving message."<<std::endl;
1300 MPI_Abort(comm_, 999);
1309 if(destPairs!=sourcePairs)
1312 delete[] sourcePairs;
1318 template<
typename T,
typename A>
1319 inline void RemoteIndices<T,A>::unpackIndices(RemoteIndexList& remote,
1329 if(remoteEntries==0)
1333 MPI_Unpack(p_in, bufferSize, position, &index, 1,
1335 GlobalIndex oldGlobal=index.global();
1336 int n_in=0, localIndex=0;
1339 while(localIndex<localEntries) {
1340 if(local[localIndex]->global()==index.global()) {
1341 int oldLocalIndex=localIndex;
1343 while(localIndex<localEntries &&
1344 local[localIndex]->global()==index.global()) {
1345 if(!fromOurSelf || index.local().attribute() !=
1346 local[localIndex]->local().attribute())
1348 remote.push_back(RemoteIndex(index.local().attribute(),
1349 local[localIndex]));
1354 if((++n_in) < remoteEntries) {
1355 MPI_Unpack(p_in, bufferSize, position, &index, 1,
1357 if(index.global()==oldGlobal)
1359 localIndex=oldLocalIndex;
1361 oldGlobal=index.global();
1369 if (local[localIndex]->global()<index.global()) {
1374 if((++n_in) < remoteEntries) {
1375 MPI_Unpack(p_in, bufferSize, position, &index, 1,
1377 oldGlobal=index.global();
1385 while(++n_in < remoteEntries)
1386 MPI_Unpack(p_in, bufferSize, position, &index, 1,
1391 template<
typename T,
typename A>
1392 inline void RemoteIndices<T,A>::unpackIndices(RemoteIndexList& send,
1393 RemoteIndexList& receive,
1395 PairType** localSource,
1396 int localSourceEntries,
1397 PairType** localDest,
1398 int localDestEntries,
1404 int n_in=0, sourceIndex=0, destIndex=0;
1407 while(n_in<remoteEntries && (sourceIndex<localSourceEntries || destIndex<localDestEntries)) {
1410 MPI_Unpack(p_in, bufferSize, position, &index, 1,
1415 while(sourceIndex<localSourceEntries && localSource[sourceIndex]->global()<index.global())
1418 while(destIndex<localDestEntries && localDest[destIndex]->global()<index.global())
1422 if(sourceIndex<localSourceEntries && localSource[sourceIndex]->global()==index.global())
1423 send.push_back(RemoteIndex(index.local().attribute(),
1424 localSource[sourceIndex]));
1426 if(destIndex < localDestEntries && localDest[destIndex]->global() == index.global())
1427 receive.push_back(RemoteIndex(index.local().attribute(),
1428 localDest[sourceIndex]));
1433 template<
typename T,
typename A>
1434 inline void RemoteIndices<T,A>::free()
1436 typedef typename RemoteIndexMap::iterator Iterator;
1437 Iterator lend = remoteIndices_.end();
1438 for(Iterator lists=remoteIndices_.begin(); lists != lend; ++lists) {
1439 if(lists->second.first==lists->second.second) {
1441 delete lists->second.first;
1443 delete lists->second.first;
1444 delete lists->second.second;
1447 remoteIndices_.clear();
1451 template<
typename T,
typename A>
1452 inline int RemoteIndices<T,A>::neighbours()
const
1454 return remoteIndices_.size();
1457 template<
typename T,
typename A>
1458 template<
bool ignorePublic>
1459 inline void RemoteIndices<T,A>::rebuild()
1463 ignorePublic!=publicIgnored || !
1467 buildRemote<ignorePublic>(includeSelf);
1469 sourceSeqNo_ = source_->seqNo();
1470 destSeqNo_ = target_->seqNo();
1472 publicIgnored=ignorePublic;
1478 template<
typename T,
typename A>
1479 inline bool RemoteIndices<T,A>::isSynced()
const
1481 return sourceSeqNo_==source_->seqNo() && destSeqNo_ ==target_->seqNo();
1484 template<
typename T,
typename A>
1485 template<
bool mode,
bool send>
1486 RemoteIndexListModifier<T,A,mode> RemoteIndices<T,A>::getModifier(
int process)
1492 sourceSeqNo_ = source_->seqNo();
1493 destSeqNo_ = target_->seqNo();
1495 typename RemoteIndexMap::iterator found = remoteIndices_.find(process);
1497 if(found == remoteIndices_.end())
1499 if(source_ != target_)
1500 found = remoteIndices_.insert(found, std::make_pair(process,
1501 std::make_pair(
new RemoteIndexList(),
1502 new RemoteIndexList())));
1504 RemoteIndexList* rlist =
new RemoteIndexList();
1505 found = remoteIndices_.insert(found,
1506 std::make_pair(process,
1507 std::make_pair(rlist, rlist)));
1514 return RemoteIndexListModifier<T,A,mode>(*source_, *(found->second.first));
1516 return RemoteIndexListModifier<T,A,mode>(*target_, *(found->second.second));
1519 template<
typename T,
typename A>
1520 inline typename RemoteIndices<T,A>::const_iterator
1521 RemoteIndices<T,A>::find(
int proc)
const
1523 return remoteIndices_.find(proc);
1526 template<
typename T,
typename A>
1527 inline typename RemoteIndices<T,A>::const_iterator
1528 RemoteIndices<T,A>::begin()
const
1530 return remoteIndices_.begin();
1533 template<
typename T,
typename A>
1534 inline typename RemoteIndices<T,A>::const_iterator
1535 RemoteIndices<T,A>::end()
const
1537 return remoteIndices_.end();
1541 template<
typename T,
typename A>
1544 if(neighbours()!=ri.neighbours())
1547 typedef RemoteIndexList RList;
1548 typedef typename std::map<int,std::pair<RList*,RList*> >::const_iterator const_iterator;
1550 const const_iterator rend = remoteIndices_.end();
1552 for(const_iterator rindex = remoteIndices_.begin(), rindex1=ri.remoteIndices_.begin(); rindex!=rend; ++rindex, ++rindex1) {
1553 if(rindex->first != rindex1->first)
1555 if(*(rindex->second.first) != *(rindex1->second.first))
1557 if(*(rindex->second.second) != *(rindex1->second.second))
1563 template<
class T,
class A,
bool mode>
1564 RemoteIndexListModifier<T,A,mode>::RemoteIndexListModifier(
const ParallelIndexSet& indexSet,
1565 RemoteIndexList& rList)
1566 : rList_(&rList), indexSet_(&indexSet), iter_(rList.beginModify()), end_(rList.end()), first_(true)
1568 if(MODIFYINDEXSET) {
1570 for(ConstIterator iter=iter_; iter != end_; ++iter)
1571 glist_.push_back(iter->localIndexPair().global());
1572 giter_ = glist_.beginModify();
1576 template<
typename T,
typename A,
bool mode>
1577 RemoteIndexListModifier<T,A,mode>::RemoteIndexListModifier(
const RemoteIndexListModifier<T,A,mode>& other)
1578 : rList_(other.rList_), indexSet_(other.indexSet_),
1579 glist_(other.glist_), iter_(other.iter_), giter_(other.giter_), end_(other.end_),
1580 first_(other.first_), last_(other.last_)
1583 template<
typename T,
typename A,
bool mode>
1584 inline void RemoteIndexListModifier<T,A,mode>::repairLocalIndexPointers()
1586 if(MODIFYINDEXSET) {
1588 #ifdef DUNE_ISTL_WITH_CHECKING
1589 if(indexSet_->state()!=
GROUND)
1590 DUNE_THROW(InvalidIndexSetState,
"Index has to be in ground mode for repairing pointers to indices");
1592 typedef typename ParallelIndexSet::const_iterator IndexIterator;
1593 typedef typename GlobalList::const_iterator GlobalIterator;
1594 typedef typename RemoteIndexList::iterator Iterator;
1595 GlobalIterator giter = glist_.begin();
1596 IndexIterator index = indexSet_->begin();
1598 for(Iterator iter=rList_->begin(); iter != end_; ++iter) {
1599 while(index->global()<*giter) {
1601 #ifdef DUNE_ISTL_WITH_CHECKING
1602 if(index == indexSet_->end())
1603 DUNE_THROW(InvalidPosition,
"No such global index in set!");
1607 #ifdef DUNE_ISTL_WITH_CHECKING
1608 if(index->global() != *giter)
1609 DUNE_THROW(InvalidPosition,
"No such global index in set!");
1611 iter->localIndex_ = &(*index);
1616 template<
typename T,
typename A,
bool mode>
1617 inline void RemoteIndexListModifier<T,A,mode>::insert(
const RemoteIndex& index)
1619 static_assert(!mode,
"Not allowed if the mode indicates that new indices"
1620 "might be added to the underlying index set. Use "
1621 "insert(const RemoteIndex&, const GlobalIndex&) instead");
1623 #ifdef DUNE_ISTL_WITH_CHECKING
1624 if(!first_ && index.localIndexPair().global()<last_)
1625 DUNE_THROW(InvalidPosition,
"Modifcation of remote indices have to occur with ascending global index.");
1628 while(iter_ != end_ && iter_->localIndexPair().global() < index.localIndexPair().global()) {
1633 assert(iter_==end_ || iter_->localIndexPair().global() != index.localIndexPair().global());
1634 iter_.insert(index);
1635 last_ = index.localIndexPair().global();
1639 template<
typename T,
typename A,
bool mode>
1640 inline void RemoteIndexListModifier<T,A,mode>::insert(
const RemoteIndex& index,
const GlobalIndex& global)
1642 static_assert(mode,
"Not allowed if the mode indicates that no new indices"
1643 "might be added to the underlying index set. Use "
1644 "insert(const RemoteIndex&) instead");
1645 #ifdef DUNE_ISTL_WITH_CHECKING
1646 if(!first_ && global<last_)
1647 DUNE_THROW(InvalidPosition,
"Modification of remote indices have to occur with ascending global index.");
1650 while(iter_ != end_ && *giter_ < global) {
1656 assert(iter_->localIndexPair().global() != global);
1657 iter_.insert(index);
1658 giter_.insert(global);
1664 template<
typename T,
typename A,
bool mode>
1665 bool RemoteIndexListModifier<T,A,mode>::remove(
const GlobalIndex& global)
1667 #ifdef DUNE_ISTL_WITH_CHECKING
1668 if(!first_ && global<last_)
1669 DUNE_THROW(InvalidPosition,
"Modifcation of remote indices have to occur with ascending global index.");
1674 if(MODIFYINDEXSET) {
1676 while(iter_!=end_ && *giter_< global) {
1680 if(*giter_ == global) {
1686 while(iter_!=end_ && iter_->localIndexPair().global() < global)
1689 if(iter_->localIndexPair().global()==global) {
1700 template<
typename T,
typename A>
1702 inline typename RemoteIndices<T,A>::CollectiveIteratorT RemoteIndices<T,A>::iterator()
const
1704 return CollectiveIterator<T,A>(remoteIndices_, send);
1707 template<
typename T,
typename A>
1708 inline MPI_Comm RemoteIndices<T,A>::communicator()
const
1714 template<
typename T,
typename A>
1715 CollectiveIterator<T,A>::CollectiveIterator(
const RemoteIndexMap& pmap,
bool send)
1717 typedef typename RemoteIndexMap::const_iterator const_iterator;
1719 const const_iterator end=pmap.end();
1720 for(const_iterator process=pmap.begin(); process != end; ++process) {
1721 const RemoteIndexList* list = send ? process->second.first : process->second.second;
1722 typedef typename RemoteIndexList::const_iterator iterator;
1723 map_.insert(std::make_pair(process->first,
1724 std::pair<iterator, const iterator>(list->begin(), list->end())));
1728 template<
typename T,
typename A>
1729 inline void CollectiveIterator<T,A>::advance(
const GlobalIndex& index)
1731 typedef typename Map::iterator iterator;
1732 typedef typename Map::const_iterator const_iterator;
1733 const const_iterator end = map_.end();
1735 for(iterator iter = map_.begin(); iter != end;) {
1737 typename RemoteIndexList::const_iterator current = iter->second.first;
1738 typename RemoteIndexList::const_iterator rend = iter->second.second;
1739 RemoteIndex remoteIndex;
1741 remoteIndex = *current;
1743 while(iter->second.first!=iter->second.second && iter->second.first->localIndexPair().global()<index)
1744 ++(iter->second.first);
1747 if(iter->second.first == iter->second.second)
1757 template<
typename T,
typename A>
1758 inline void CollectiveIterator<T,A>::advance(
const GlobalIndex& index,
1759 const Attribute& attribute)
1761 typedef typename Map::iterator iterator;
1762 typedef typename Map::const_iterator const_iterator;
1763 const const_iterator end = map_.end();
1765 for(iterator iter = map_.begin(); iter != end;) {
1767 typename RemoteIndexList::const_iterator current = iter->second.first;
1768 typename RemoteIndexList::const_iterator rend = iter->second.second;
1769 RemoteIndex remoteIndex;
1771 remoteIndex = *current;
1774 while(iter->second.first!=iter->second.second && iter->second.first->localIndexPair().global()<index)
1775 ++(iter->second.first);
1778 while(iter->second.first!=iter->second.second
1779 && iter->second.first->localIndexPair().global()==index
1780 && iter->second.first->localIndexPair().local().attribute()<attribute)
1781 ++(iter->second.first);
1784 if(iter->second.first == iter->second.second)
1791 attribute_=attribute;
1795 template<
typename T,
typename A>
1796 inline CollectiveIterator<T,A>& CollectiveIterator<T,A>::operator++()
1798 typedef typename Map::iterator iterator;
1799 typedef typename Map::const_iterator const_iterator;
1800 const const_iterator end = map_.end();
1802 for(iterator iter = map_.begin(); iter != end;) {
1804 typename RemoteIndexList::const_iterator current = iter->second.first;
1805 typename RemoteIndexList::const_iterator rend = iter->second.second;
1808 if(iter->second.first->localIndexPair().global()==index_ &&
1809 (noattribute || iter->second.first->localIndexPair().local().attribute() == attribute_))
1810 ++(iter->second.first);
1813 if(iter->second.first == iter->second.second)
1822 template<
typename T,
typename A>
1823 inline bool CollectiveIterator<T,A>::empty()
1825 return map_.empty();
1828 template<
typename T,
typename A>
1829 inline typename CollectiveIterator<T,A>::iterator
1830 CollectiveIterator<T,A>::begin()
1833 return iterator(map_.begin(), map_.end(), index_);
1835 return iterator(map_.begin(), map_.end(), index_,
1839 template<
typename T,
typename A>
1840 inline typename CollectiveIterator<T,A>::iterator
1841 CollectiveIterator<T,A>::end()
1843 return iterator(map_.end(), map_.end(), index_);
1846 template<
typename TG,
typename TA>
1847 inline std::ostream&
operator<<(std::ostream& os,
const RemoteIndex<TG,TA>& index)
1849 os<<
"[global="<<index.localIndexPair().global()<<
", remote attribute="<<index.attribute()<<
" local attribute="<<index.localIndexPair().local().attribute()<<
"]";
1853 template<
typename T,
typename A>
1854 inline std::ostream&
operator<<(std::ostream& os,
const RemoteIndices<T,A>& indices)
1857 MPI_Comm_rank(indices.comm_, &rank);
1859 typedef typename RemoteIndices<T,A>::RemoteIndexList RList;
1860 typedef typename std::map<int,std::pair<RList*,RList*> >::const_iterator const_iterator;
1862 const const_iterator rend = indices.remoteIndices_.end();
1864 for(const_iterator rindex = indices.remoteIndices_.begin(); rindex!=rend; ++rindex) {
1865 os<<rank<<
": Prozess "<<rindex->first<<
":";
1867 if(!rindex->second.first->empty()) {
1870 const typename RList::const_iterator send= rindex->second.first->end();
1872 for(
typename RList::const_iterator index = rindex->second.first->begin();
1873 index != send; ++index)
1877 if(!rindex->second.second->empty()) {
1878 os<<rank<<
": Prozess "<<rindex->first<<
": "<<
"receive: ";
1880 for(
const auto& index : *(rindex->second.second))
1883 os<<std::endl<<std::flush;