interface.hh
Go to the documentation of this file.00001
00002 #ifndef DUNE_INTERFACE_HH
00003 #define DUNE_INTERFACE_HH
00004
00005 #if HAVE_MPI
00006
00007 #include"remoteindices.hh"
00008 #include<dune/common/enumset.hh>
00009
00010 namespace Dune
00011 {
00031 template<typename T>
00032 class InterfaceBuilder
00033 {
00034 public:
00035 class RemotexIndicesStateError : public Exception
00036 {};
00037
00041 typedef T ParallelIndexSet;
00042
00046 typedef Dune::RemoteIndices<ParallelIndexSet> RemoteIndices;
00047
00051 typedef typename RemoteIndices::GlobalIndex GlobalIndex;
00052
00056 typedef typename RemoteIndices::Attribute Attribute;
00057
00058 virtual ~InterfaceBuilder()
00059 {}
00060
00061 protected:
00065 InterfaceBuilder()
00066 {}
00067
00105 template<class T1, class T2, class Op, bool send>
00106 void buildInterface (const RemoteIndices& remoteIndices,
00107 const T1& sourceFlags, const T2& destFlags,
00108 Op& functor) const;
00109 };
00110
00118 class InterfaceInformation
00119 {
00120
00121 public:
00122
00126 size_t size() const
00127 {
00128 return size_;
00129 }
00134 uint32_t& operator[](size_t i)
00135 {
00136 assert(i<size_);
00137 return indices_[i];
00138 }
00143 uint32_t operator[](size_t i) const
00144 {
00145 assert(i<size_);
00146 return indices_[i];
00147 }
00152 void reserve(size_t size)
00153 {
00154 indices_ = new uint32_t[size];
00155 maxSize_ = size;
00156
00157 }
00161 void free()
00162 {
00163 delete[] indices_;
00164 maxSize_ = 0;
00165 size_=0;
00166 indices_=0;
00167 }
00171 void add(uint32_t index)
00172 {
00173 assert(size_<maxSize_);
00174 indices_[size_++]=index;
00175 }
00176
00177 InterfaceInformation()
00178 : size_(0), maxSize_(0), indices_(0)
00179 {}
00180
00181 virtual ~InterfaceInformation()
00182 {
00183 if(indices_!=0){
00184 std::cerr<<"InterfaceInformation::free() should be called before "<<
00185 "destructor!"<<std::endl;
00186 }
00187
00188 }
00189
00190 private:
00194 size_t size_;
00198 size_t maxSize_;
00202 uint32_t* indices_;
00203 };
00204
00216 template<typename T>
00217 class Interface : public InterfaceBuilder<T>
00218 {
00219
00220 public:
00221 typedef InterfaceInformation Information;
00222
00223 typedef std::map<int,std::pair<Information,Information> > InformationMap;
00224
00228 typedef T ParallelIndexSet;
00229
00233 typedef Dune::RemoteIndices<ParallelIndexSet> RemoteIndices;
00234
00235
00239 typedef typename RemoteIndices::GlobalIndex GlobalIndex;
00240
00244 typedef typename RemoteIndices::Attribute Attribute;
00245
00262 template<typename T1, typename T2>
00263 void build(const RemoteIndices& remoteIndices, const T1& sourceFlags,
00264 const T2& destFlags);
00265
00269 void free();
00270
00274 MPI_Comm communicator() const;
00275
00284 const InformationMap& interfaces() const;
00285
00286 Interface()
00287 : interfaces_()
00288 {}
00289
00293 void print() const;
00294
00298 virtual ~Interface();
00299
00300 private:
00304 const RemoteIndices* remoteIndices_;
00305
00313 InformationMap interfaces_;
00314
00315 template<bool send>
00316 class InformationBuilder
00317 {
00318 public:
00319 InformationBuilder(InformationMap& interfaces)
00320 : interfaces_(interfaces)
00321 {}
00322
00323 void reserve(int proc, int size)
00324 {
00325 if(send)
00326 interfaces_[proc].first.reserve(size);
00327 else
00328 interfaces_[proc].second.reserve(size);
00329 }
00330 void add(int proc, int local)
00331 {
00332 if(send){
00333 interfaces_[proc].first.add(local);
00334 }else{
00335 interfaces_[proc].second.add(local);
00336 }
00337 }
00338
00339 private:
00340 InformationMap& interfaces_;
00341 };
00342 };
00343
00344 template<typename T>
00345 template<class T1, class T2, class Op, bool send>
00346 void InterfaceBuilder<T>::buildInterface(const RemoteIndices& remoteIndices, const T1& sourceFlags, const T2& destFlags, Op& interfaceInformation) const
00347 {
00348
00349 if(!remoteIndices.isSynced())
00350 DUNE_THROW(RemotexIndicesStateError,"RemoteIndices is not in sync with the index set. Call RemoteIndices::rebuild first!");
00351
00352 typedef typename RemoteIndices::RemoteIndexMap::const_iterator const_iterator;
00353 typedef typename RemoteIndices::ParallelIndexSet::const_iterator LocalIterator;
00354
00355 const const_iterator end=remoteIndices.end();
00356
00357 int rank;
00358
00359 MPI_Comm_rank(remoteIndices.communicator(), &rank);
00360
00361
00362 for(const_iterator process=remoteIndices.begin(); process != end; ++process){
00363
00364 int size=0;
00365 LocalIterator localIndex = send ? remoteIndices.source_->begin() : remoteIndices.target_->begin();
00366 const LocalIterator localEnd = send ? remoteIndices.source_->end() : remoteIndices.target_->end();
00367 typedef typename RemoteIndices::RemoteIndexList::const_iterator RemoteIterator;
00368 const RemoteIterator remoteEnd = send ? process->second.first->end() :
00369 process->second.second->end();
00370 RemoteIterator remote = send ? process->second.first->begin() : process->second.second->begin();
00371
00372 while(localIndex!=localEnd && remote!=remoteEnd){
00373 if( send ? destFlags.contains(remote->attribute()) :
00374 sourceFlags.contains(remote->attribute())){
00375
00376 while(localIndex->global()<remote->localIndexPair().global()){
00377 localIndex++;
00378 assert(localIndex != localEnd);
00379 }
00380 assert(localIndex->global()==remote->localIndexPair().global());
00381
00382
00383 if( send ? sourceFlags.contains(localIndex->local().attribute()) :
00384 destFlags.contains(localIndex->local().attribute()))
00385 ++size;
00386 }
00387 ++remote;
00388 }
00389 interfaceInformation.reserve(process->first, size);
00390 }
00391
00392
00393
00394 CollectiveIterator<T> remote = remoteIndices.template iterator<send>();
00395 LocalIterator localIndex = send ? remoteIndices.source_->begin() : remoteIndices.target_->begin();
00396 const LocalIterator localEnd = send ? remoteIndices.source_->end() : remoteIndices.target_->end();
00397
00398 while(localIndex!=localEnd && !remote.empty()){
00399 if( send ? sourceFlags.contains(localIndex->local().attribute()) :
00400 destFlags.contains(localIndex->local().attribute()))
00401 {
00402
00403 remote.advance(localIndex->global());
00404
00405 typedef typename CollectiveIterator<T>::iterator ValidIterator;
00406 const ValidIterator end = remote.end();
00407 ValidIterator validEntry = remote.begin();
00408
00409 for(int i=0; validEntry != end; ++i){
00410 if( send ? destFlags.contains(validEntry->attribute()) :
00411 sourceFlags.contains(validEntry->attribute())){
00412
00413 interfaceInformation.add(validEntry.process(),localIndex->local());
00414 }
00415 ++validEntry;
00416 }
00417 }
00418 ++localIndex;
00419 }
00420 }
00421
00422 template<typename T>
00423 inline MPI_Comm Interface<T>::communicator() const
00424 {
00425 return remoteIndices_->communicator();
00426
00427 }
00428
00429 template<typename T>
00430 inline const std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface<T>::interfaces() const
00431 {
00432 return interfaces_;
00433 }
00434
00435 template<typename T>
00436 void Interface<T>::print() const
00437 {
00438 typedef typename InformationMap::const_iterator const_iterator;
00439 const const_iterator end=interfaces_.end();
00440 int rank;
00441 MPI_Comm_rank(communicator(), &rank);
00442
00443 for(const_iterator infoPair=interfaces_.begin(); infoPair!=end; ++infoPair){
00444 {
00445 std::cout<<rank<<": send for process "<<infoPair->first<<": ";
00446 const InterfaceInformation& info(infoPair->second.first);
00447 for(size_t i=0; i < info.size(); i++)
00448 std::cout<<info[i]<<" ";
00449 std::cout<<std::endl;
00450 }{
00451
00452 std::cout<<rank<<": receive for process "<<infoPair->first<<": ";
00453 const InterfaceInformation& info(infoPair->second.second);
00454 for(size_t i=0; i < info.size(); i++)
00455 std::cout<<info[i]<<" ";
00456 std::cout<<std::endl;
00457 }
00458
00459 }
00460 }
00461
00462 template<typename T>
00463 template<typename T1, typename T2>
00464 void Interface<T>::build(const RemoteIndices& remoteIndices, const T1& sourceFlags,
00465 const T2& destFlags)
00466 {
00467 remoteIndices_=&remoteIndices;
00468
00469 assert(interfaces_.empty());
00470
00471
00472 InformationBuilder<true> sendInformation(interfaces_);
00473 this->template buildInterface<T1,T2,InformationBuilder<true>,true>(remoteIndices, sourceFlags,
00474 destFlags, sendInformation);
00475
00476
00477 InformationBuilder<false> recvInformation(interfaces_);
00478 this->template buildInterface<T1,T2,InformationBuilder<false>,false>(remoteIndices,sourceFlags,
00479 destFlags, recvInformation);
00480
00481 }
00482
00483 template<typename T>
00484 void Interface<T>::free()
00485 {
00486 typedef InformationMap::iterator iterator;
00487 typedef InformationMap::const_iterator const_iterator;
00488 const const_iterator end = interfaces_.end();
00489 for(iterator interfacePair = interfaces_.begin(); interfacePair != end; ++interfacePair){
00490 interfacePair->second.first.free();
00491 interfacePair->second.second.free();
00492 }
00493 interfaces_.clear();
00494 }
00495 template<typename T>
00496 Interface<T>::~Interface()
00497 {
00498 free();
00499 }
00501 }
00502
00503 #endif // HAVE_MPI
00504
00505 #endif