interface.hh
Go to the documentation of this file.00001
00002 #ifndef DUNE_INTERFACE_HH
00003 #define DUNE_INTERFACE_HH
00004
00005 #include"remoteindices.hh"
00006 #include<dune/common/enumset.hh>
00007 namespace Dune
00008 {
00028 template<typename T>
00029 class InterfaceBuilder
00030 {
00031 public:
00032
00036 typedef T ParallelIndexSet;
00037
00041 typedef Dune :: RemoteIndices<ParallelIndexSet> RemoteIndices;
00042
00046 typedef typename RemoteIndices::GlobalIndex GlobalIndex;
00047
00051 typedef typename RemoteIndices::Attribute Attribute;
00052
00053 virtual ~InterfaceBuilder()
00054 {}
00055
00056 protected:
00060 InterfaceBuilder()
00061 {}
00062
00100 template<class T1, class T2, class Op, bool send>
00101 void buildInterface (const RemoteIndices& remoteIndices,
00102 const T1& sourceFlags, const T2& destFlags,
00103 Op& functor) const;
00104 };
00105
00113 class InterfaceInformation
00114 {
00115
00116 public:
00117
00121 size_t size() const
00122 {
00123 return size_;
00124 }
00129 uint32_t& operator[](size_t i)
00130 {
00131 assert(i<size_);
00132 return indices_[i];
00133 }
00138 uint32_t operator[](size_t i) const
00139 {
00140 assert(i<size_);
00141 return indices_[i];
00142 }
00147 void reserve(size_t size)
00148 {
00149 indices_ = new uint32_t[size];
00150 maxSize_ = size;
00151
00152 }
00156 void free()
00157 {
00158 delete[] indices_;
00159 maxSize_ = 0;
00160 size_=0;
00161 indices_=0;
00162 }
00166 void add(uint32_t index)
00167 {
00168 assert(size_<maxSize_);
00169 indices_[size_++]=index;
00170 }
00171
00172 InterfaceInformation()
00173 : size_(0), maxSize_(0), indices_(0)
00174 {}
00175
00176 virtual ~InterfaceInformation()
00177 {
00178 if(indices_!=0){
00179 std::cerr<<"InterfaceInformation::free() should be called before "<<
00180 "destructor!"<<std::endl;
00181 }
00182
00183 }
00184
00185 private:
00189 size_t size_;
00193 size_t maxSize_;
00197 uint32_t* indices_;
00198 };
00199
00211 template<typename T>
00212 class Interface : public InterfaceBuilder<T>
00213 {
00214
00215 public:
00216 typedef InterfaceInformation Information;
00217
00218 typedef std::map<int,std::pair<Information,Information> > InformationMap;
00219
00223 typedef T ParallelIndexSet;
00224
00228 typedef Dune :: RemoteIndices<ParallelIndexSet> RemoteIndices;
00229
00230
00234 typedef typename RemoteIndices::GlobalIndex GlobalIndex;
00235
00239 typedef typename RemoteIndices::Attribute Attribute;
00240
00257 template<typename T1, typename T2>
00258 void build(const RemoteIndices& remoteIndices, const T1& sourceFlags,
00259 const T2& destFlags);
00260
00264 void free();
00265
00269 MPI_Comm communicator() const;
00270
00279 const InformationMap& interfaces() const;
00280
00281 Interface()
00282 : interfaces_()
00283 {}
00284
00288 void print() const;
00289
00293 virtual ~Interface();
00294
00295 private:
00299 const RemoteIndices* remoteIndices_;
00300
00308 InformationMap interfaces_;
00309
00310 template<bool send>
00311 class InformationBuilder
00312 {
00313 public:
00314 InformationBuilder(InformationMap& interfaces)
00315 : interfaces_(interfaces)
00316 {}
00317
00318 void reserve(int proc, int size)
00319 {
00320 if(send)
00321 interfaces_[proc].first.reserve(size);
00322 else
00323 interfaces_[proc].second.reserve(size);
00324 }
00325 void add(int proc, int local)
00326 {
00327 if(send){
00328 interfaces_[proc].first.add(local);
00329 }else{
00330 interfaces_[proc].second.add(local);
00331 }
00332 }
00333
00334 private:
00335 InformationMap& interfaces_;
00336 };
00337 };
00338
00339 template<typename T>
00340 template<class T1, class T2, class Op, bool send>
00341 void InterfaceBuilder<T>::buildInterface(const RemoteIndices& remoteIndices, const T1& sourceFlags, const T2& destFlags, Op& interfaceInformation) const
00342 {
00343
00344 typedef typename RemoteIndices::RemoteIndexMap::const_iterator const_iterator;
00345 typedef typename RemoteIndices::ParallelIndexSet::const_iterator LocalIterator;
00346
00347 const const_iterator end=remoteIndices.end();
00348
00349 int rank;
00350
00351 MPI_Comm_rank(remoteIndices.communicator(), &rank);
00352
00353
00354 for(const_iterator process=remoteIndices.begin(); process != end; ++process){
00355
00356 int size=0;
00357 LocalIterator localIndex = send ? remoteIndices.source_->begin() : remoteIndices.target_->begin();
00358 const LocalIterator localEnd = send ? remoteIndices.source_->end() : remoteIndices.target_->end();
00359 typedef typename RemoteIndices::RemoteIndexList::const_iterator RemoteIterator;
00360 const RemoteIterator remoteEnd = send ? process->second.first->end() :
00361 process->second.second->end();
00362 RemoteIterator remote = send ? process->second.first->begin() : process->second.second->begin();
00363
00364 while(localIndex!=localEnd && remote!=remoteEnd){
00365 if( send ? destFlags.contains(remote->attribute()) :
00366 sourceFlags.contains(remote->attribute())){
00367
00368 while(localIndex->global()<remote->localIndexPair().global()){
00369 localIndex++;
00370 assert(localIndex != localEnd);
00371 }
00372 assert(localIndex->global()==remote->localIndexPair().global());
00373
00374
00375 if( send ? sourceFlags.contains(localIndex->local().attribute()) :
00376 destFlags.contains(localIndex->local().attribute()))
00377 ++size;
00378 }
00379 ++remote;
00380 }
00381 interfaceInformation.reserve(process->first, size);
00382 }
00383
00384
00385
00386 CollectiveIterator<T> remote = remoteIndices.template iterator<send>();
00387 LocalIterator localIndex = send ? remoteIndices.source_->begin() : remoteIndices.target_->begin();
00388 const LocalIterator localEnd = send ? remoteIndices.source_->end() : remoteIndices.target_->end();
00389
00390 while(localIndex!=localEnd && !remote.empty()){
00391 if( send ? sourceFlags.contains(localIndex->local().attribute()) :
00392 destFlags.contains(localIndex->local().attribute()))
00393 {
00394
00395 remote.advance(localIndex->global());
00396
00397 typedef typename CollectiveIterator<T>::iterator ValidIterator;
00398 const ValidIterator end = remote.end();
00399 ValidIterator validEntry = remote.begin();
00400
00401 for(int i=0; validEntry != end; ++i){
00402 if( send ? destFlags.contains(validEntry->attribute()) :
00403 sourceFlags.contains(validEntry->attribute())){
00404
00405 interfaceInformation.add(validEntry.process(),localIndex->local());
00406 }
00407 ++validEntry;
00408 }
00409 }
00410 ++localIndex;
00411 }
00412 }
00413
00414 template<typename T>
00415 inline MPI_Comm Interface<T>::communicator() const
00416 {
00417 return remoteIndices_->communicator();
00418
00419 }
00420
00421 template<typename T>
00422 inline const std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface<T>::interfaces() const
00423 {
00424 return interfaces_;
00425 }
00426
00427 template<typename T>
00428 void Interface<T>::print() const
00429 {
00430 typedef typename InformationMap::const_iterator const_iterator;
00431 const const_iterator end=interfaces_.end();
00432 int rank;
00433 MPI_Comm_rank(communicator(), &rank);
00434
00435 for(const_iterator infoPair=interfaces_.begin(); infoPair!=end; ++infoPair){
00436 {
00437 std::cout<<rank<<": send for process "<<infoPair->first<<": ";
00438 const InterfaceInformation& info(infoPair->second.first);
00439 for(size_t i=0; i < info.size(); i++)
00440 std::cout<<info[i]<<" ";
00441 std::cout<<std::endl;
00442 }{
00443
00444 std::cout<<rank<<": receive for process "<<infoPair->first<<": ";
00445 const InterfaceInformation& info(infoPair->second.second);
00446 for(size_t i=0; i < info.size(); i++)
00447 std::cout<<info[i]<<" ";
00448 std::cout<<std::endl;
00449 }
00450
00451 }
00452 }
00453
00454 template<typename T>
00455 template<typename T1, typename T2>
00456 void Interface<T>::build(const RemoteIndices& remoteIndices, const T1& sourceFlags,
00457 const T2& destFlags)
00458 {
00459 remoteIndices_=&remoteIndices;
00460
00461 assert(interfaces_.empty());
00462
00463
00464 InformationBuilder<true> sendInformation(interfaces_);
00465 this->template buildInterface<T1,T2,InformationBuilder<true>,true>(remoteIndices, sourceFlags,
00466 destFlags, sendInformation);
00467
00468
00469 InformationBuilder<false> recvInformation(interfaces_);
00470 this->template buildInterface<T1,T2,InformationBuilder<false>,false>(remoteIndices,sourceFlags,
00471 destFlags, recvInformation);
00472
00473 }
00474
00475 template<typename T>
00476 void Interface<T>::free()
00477 {
00478 typedef InformationMap::iterator iterator;
00479 typedef InformationMap::const_iterator const_iterator;
00480 const const_iterator end = interfaces_.end();
00481 for(iterator interfacePair = interfaces_.begin(); interfacePair != end; ++interfacePair){
00482 interfacePair->second.first.free();
00483 interfacePair->second.second.free();
00484 }
00485 interfaces_.clear();
00486 }
00487 template<typename T>
00488 Interface<T>::~Interface()
00489 {
00490 free();
00491 }
00493 }
00494
00495 #endif