interface.hh

Go to the documentation of this file.
00001 // $Id: interface.hh 1184 2010-03-12 13:10:59Z mdroh_01 $
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   class InterfaceBuilder
00032   {
00033   public: 
00034     class RemotexIndicesStateError : public Exception
00035     {};
00036         
00037     virtual ~InterfaceBuilder()
00038     {}
00039 
00040   protected:
00044     InterfaceBuilder()
00045     {}
00046     
00084     template<class R, class T1, class T2, class Op, bool send>
00085     void buildInterface (const R& remoteIndices, 
00086                          const T1& sourceFlags, const T2& destFlags,
00087                          Op& functor) const;
00088   };
00089   
00097   class InterfaceInformation
00098   {
00099     
00100   public:
00101     
00105     size_t size() const
00106     {
00107       return size_;
00108     }
00113     std::size_t& operator[](size_t i)
00114     {
00115       assert(i<size_);
00116       return indices_[i];
00117     }
00122     std::size_t operator[](size_t i) const
00123     {
00124       assert(i<size_);
00125       return indices_[i];
00126     }
00131     void reserve(size_t size)
00132     {
00133       indices_ = new std::size_t[size];
00134       maxSize_ = size;
00135       
00136     }
00140     void free()
00141     {
00142       if(indices_)
00143         delete[] indices_;
00144       maxSize_ = 0;
00145       size_=0;
00146       indices_=0;
00147     }
00151     void add(std::size_t index)
00152     {
00153       assert(size_<maxSize_);
00154       indices_[size_++]=index;
00155     }
00156     
00157     InterfaceInformation() 
00158       : size_(0), maxSize_(0), indices_(0)
00159     {}
00160     
00161     virtual ~InterfaceInformation()
00162     {
00163     }
00164 
00165     bool operator!=(const InterfaceInformation& o) const
00166     {
00167       return !operator==(o);
00168     }
00169     
00170     bool operator==(const InterfaceInformation& o) const
00171     {
00172       if(size_!=o.size_)
00173         return false;
00174       for(std::size_t i=0; i< size_; ++i)
00175         if(indices_[i]!=o.indices_[i])
00176           return false;
00177       return true;
00178     }
00179     
00180   private:
00184     size_t size_;
00188     size_t maxSize_;
00192     std::size_t* indices_;
00193   };
00194 
00206   class Interface : public InterfaceBuilder
00207   {
00208     
00209   public:
00210     typedef InterfaceInformation Information;
00211 
00212     typedef std::map<int,std::pair<Information,Information> > InformationMap;
00213 
00230     template<typename R, typename T1, typename T2>
00231     void build(const R& remoteIndices, const T1& sourceFlags, 
00232                const T2& destFlags);
00233 
00237     void free();
00238 
00242     MPI_Comm communicator() const;
00243 
00252     const InformationMap& interfaces() const;
00253     
00254     Interface(MPI_Comm comm)
00255       : communicator_(comm), interfaces_()
00256     {}
00257   
00258     Interface()
00259       : communicator_(MPI_COMM_NULL), interfaces_()
00260     {}
00261     
00265     void print() const;
00266 
00267     bool operator!=(const Interface& o) const
00268     {
00269       return ! operator==(o);
00270     }
00271     
00272     bool operator==(const Interface& o) const
00273     {
00274       if(communicator_!=o.communicator_)
00275         return false;
00276       if(interfaces_.size()!=o.interfaces_.size())
00277         return false;
00278       typedef InformationMap::const_iterator MIter;
00279       
00280       for(MIter m=interfaces_.begin(), om=o.interfaces_.begin();
00281           m!=interfaces_.end(); ++m, ++om)
00282         {
00283           if(om->first!=m->first)
00284             return false;
00285           if(om->second.first!=om->second.first)
00286             return false;
00287           if(om->second.second!=om->second.second)
00288             return false;
00289         }
00290       return true;
00291     }
00292       
00296     virtual ~Interface();
00297 
00298     void strip();
00299   protected:
00300     
00309     InformationMap& interfaces();
00310 
00312     MPI_Comm communicator_;
00313         
00314   private:
00322     InformationMap interfaces_;
00323     
00324     template<bool send>
00325     class InformationBuilder
00326     {
00327     public:
00328       InformationBuilder(InformationMap& interfaces)
00329         : interfaces_(interfaces)
00330       {}
00331       
00332       void reserve(int proc, int size)
00333       {
00334         if(send)
00335           interfaces_[proc].first.reserve(size);
00336         else
00337           interfaces_[proc].second.reserve(size);
00338       }
00339       void add(int proc, std::size_t local)
00340       {
00341         if(send){
00342           interfaces_[proc].first.add(local);
00343         }else{
00344           interfaces_[proc].second.add(local);
00345         }
00346       }
00347       
00348     private:
00349       InformationMap& interfaces_;
00350     };
00351   };
00352   
00353   template<class R, class T1, class T2, class Op, bool send>
00354   void InterfaceBuilder::buildInterface(const R& remoteIndices, const T1& sourceFlags, const T2& destFlags, Op& interfaceInformation) const
00355   {
00356 
00357     if(!remoteIndices.isSynced())
00358       DUNE_THROW(RemotexIndicesStateError,"RemoteIndices is not in sync with the index set. Call RemoteIndices::rebuild first!");
00359     // Allocate the memory for the data type construction.
00360     typedef R RemoteIndices;
00361     typedef typename RemoteIndices::RemoteIndexMap::const_iterator const_iterator;
00362     typedef typename RemoteIndices::ParallelIndexSet::const_iterator LocalIterator;
00363 
00364     const const_iterator end=remoteIndices.end();
00365     
00366     int rank;
00367     
00368     MPI_Comm_rank(remoteIndices.communicator(), &rank);
00369     
00370     // Allocate memory for the type construction.
00371     for(const_iterator process=remoteIndices.begin(); process != end; ++process){
00372       // Messure the number of indices send to the remote process first
00373       int size=0;
00374       LocalIterator localIndex = send ? remoteIndices.source_->begin() : remoteIndices.target_->begin();
00375       const LocalIterator localEnd = send ?  remoteIndices.source_->end() : remoteIndices.target_->end();
00376       typedef typename RemoteIndices::RemoteIndexList::const_iterator RemoteIterator;
00377       const RemoteIterator remoteEnd = send ? process->second.first->end() : 
00378         process->second.second->end();
00379       RemoteIterator remote = send ? process->second.first->begin() : process->second.second->begin();
00380         
00381       while(localIndex!=localEnd && remote!=remoteEnd){
00382           if( send ?  destFlags.contains(remote->attribute()) :
00383               sourceFlags.contains(remote->attribute())){
00384             // search for the matching local index
00385             while(localIndex->global()<remote->localIndexPair().global()){
00386               localIndex++;
00387               assert(localIndex != localEnd); // Should never happen
00388             }
00389             assert(localIndex->global()==remote->localIndexPair().global());
00390             
00391             // do we send the index?
00392             if( send ? sourceFlags.contains(localIndex->local().attribute()) :
00393                 destFlags.contains(localIndex->local().attribute()))
00394               ++size;
00395           }
00396           ++remote;
00397       }
00398       interfaceInformation.reserve(process->first, size);
00399     }
00400 
00401     // compare the local and remote indices and set up the types
00402     
00403     typedef typename RemoteIndices::CollectiveIteratorT  CIter;
00404     CIter remote = remoteIndices.template iterator<send>();
00405     LocalIterator localIndex = send ? remoteIndices.source_->begin() : remoteIndices.target_->begin();
00406     const LocalIterator localEnd = send ?  remoteIndices.source_->end() : remoteIndices.target_->end();
00407     
00408     while(localIndex!=localEnd && !remote.empty()){
00409       if( send ? sourceFlags.contains(localIndex->local().attribute()) :
00410           destFlags.contains(localIndex->local().attribute()))
00411       {
00412         // search for matching remote indices
00413         remote.advance(localIndex->global());
00414         // Iterate over the list that are positioned at global
00415         typedef typename CIter::iterator ValidIterator;
00416         const ValidIterator end = remote.end();
00417         ValidIterator validEntry = remote.begin();
00418 
00419         for(int i=0; validEntry != end; ++i){
00420           if( send ?  destFlags.contains(validEntry->attribute()) :
00421               sourceFlags.contains(validEntry->attribute())){
00422             // We will receive data for this index
00423             interfaceInformation.add(validEntry.process(),localIndex->local());
00424           }
00425           ++validEntry;
00426         }
00427       }
00428       ++localIndex;
00429     }
00430   }
00431   
00432   inline MPI_Comm Interface::communicator() const
00433   {
00434     return communicator_;
00435     
00436   }
00437   
00438   
00439   inline const std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface::interfaces() const
00440   {
00441     return interfaces_;
00442   }
00443 
00444   inline std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface::interfaces()
00445   {
00446     return interfaces_;
00447   }
00448 
00449   inline void Interface::print() const
00450   {
00451     typedef InformationMap::const_iterator const_iterator;
00452     const const_iterator end=interfaces_.end();
00453     int rank;
00454     MPI_Comm_rank(communicator(), &rank);
00455     
00456     for(const_iterator infoPair=interfaces_.begin(); infoPair!=end; ++infoPair){
00457       {
00458         std::cout<<rank<<": send for process "<<infoPair->first<<": ";
00459         const InterfaceInformation& info(infoPair->second.first);
00460         for(size_t i=0; i < info.size(); i++)
00461           std::cout<<info[i]<<" ";
00462         std::cout<<std::endl;
00463       }{
00464         
00465         std::cout<<rank<<": receive for process "<<infoPair->first<<": ";
00466         const InterfaceInformation& info(infoPair->second.second);
00467         for(size_t i=0; i < info.size(); i++)
00468           std::cout<<info[i]<<" ";
00469         std::cout<<std::endl;
00470       }
00471       
00472     }
00473   }
00474   
00475   template<typename R, typename T1, typename T2>
00476   inline void Interface::build(const R& remoteIndices, const T1& sourceFlags, 
00477                           const T2& destFlags)
00478   {
00479     communicator_=remoteIndices.communicator();
00480     
00481     assert(interfaces_.empty());
00482 
00483     // Build the send interface
00484     InformationBuilder<true> sendInformation(interfaces_);
00485     this->template buildInterface<R,T1,T2,InformationBuilder<true>,true>(remoteIndices, sourceFlags, 
00486                                                               destFlags, sendInformation);
00487     
00488     // Build the receive interface
00489     InformationBuilder<false> recvInformation(interfaces_);
00490     this->template buildInterface<R,T1,T2,InformationBuilder<false>,false>(remoteIndices,sourceFlags, 
00491                                                                 destFlags, recvInformation);
00492     strip();
00493   }
00494   inline void Interface::strip()
00495   {
00496     typedef InformationMap::iterator const_iterator;
00497     for(const_iterator interfacePair = interfaces_.begin(); interfacePair != interfaces_.end();)
00498       if(interfacePair->second.first.size()==0 && interfacePair->second.second.size()==0){
00499         interfacePair->second.first.free();
00500         interfacePair->second.second.free();
00501         const_iterator toerase=interfacePair++;
00502         interfaces_.erase(toerase);
00503       }else
00504         ++interfacePair;
00505   }
00506   
00507   inline void Interface::free()
00508   {
00509     typedef InformationMap::iterator iterator;
00510     typedef InformationMap::const_iterator const_iterator;
00511     const const_iterator end = interfaces_.end();
00512     for(iterator interfacePair = interfaces_.begin(); interfacePair != end; ++interfacePair){
00513       interfacePair->second.first.free();
00514       interfacePair->second.second.free();
00515     }
00516     interfaces_.clear();
00517   }
00518 
00519   inline Interface::~Interface()
00520   {
00521     free();
00522   }
00524 }
00525 
00526 namespace std
00527 {
00528   inline ostream& operator<<(ostream& os, const Dune::Interface& interface)
00529   {
00530     typedef Dune::Interface::InformationMap InfoMap;
00531     typedef InfoMap::const_iterator Iter;
00532     for(Iter i=interface.interfaces().begin(), end = interface.interfaces().end();
00533         i!=end; ++i)
00534       {
00535         os<<i->first<<": [ source=[";
00536         for(std::size_t j=0; j < i->second.first.size(); ++j)
00537           os<<i->second.first[j]<<" ";
00538         os<<"] size="<<i->second.first.size()<<", target=[";
00539         for(std::size_t j=0; j < i->second.second.size(); ++j)
00540           os<<i->second.second[j]<<" ";
00541         os<<"] size="<<i->second.second.size()<<"\n";
00542       }
00543     return os;
00544   }
00545 }// end namespace std
00546 #endif // HAVE_MPI
00547 
00548 #endif
Generated on Sat Apr 24 11:13:46 2010 for dune-istl by  doxygen 1.6.3