communicator.hh

Go to the documentation of this file.
00001 // $Id: communicator.hh 982 2009-01-27 17:00:32Z mblatt $
00002 #ifndef DUNE_COMMUNICATOR
00003 #define DUNE_COMMUNICATOR
00004 
00005 #include "remoteindices.hh"
00006 #include "interface.hh"
00007 #include <dune/common/exceptions.hh>
00008 #include <dune/common/typetraits.hh>
00009 
00010 #if HAVE_MPI
00011 // MPI header
00012 #include <mpi.h>
00013 
00014 namespace Dune
00015 {
00095   struct SizeOne
00096   {};
00097   
00103   struct VariableSize
00104   {};
00105   
00106 
00112   template<class V>
00113   struct CommPolicy
00114   {
00126     typedef V Type;
00127     
00133     typedef typename V::value_type IndexedType;
00134     
00139     typedef SizeOne IndexedTypeFlag;
00140 
00149     static const void* getAddress(const V& v, int index);
00150     
00156     static int getSize(const V&, int index);
00157   };
00158 
00159   template<class K, int n> class FieldVector;
00160   
00161   template<class B, class A> class VariableBlockVector;
00162   
00163   template<class K, class A, int n>
00164   struct CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >
00165   {
00166     typedef VariableBlockVector<FieldVector<K, n>, A> Type;
00167     
00168     typedef typename Type::B IndexedType;
00169     
00170     typedef VariableSize IndexedTypeFlag;
00171 
00172     static const void* getAddress(const Type& v, int i);
00173     
00174     static int getSize(const Type& v, int i);
00175   };
00176     
00180   class CommunicationError : public IOError
00181   {};
00182 
00186   template<class T>
00187   struct CopyGatherScatter
00188   {
00189     typedef typename CommPolicy<T>::IndexedType IndexedType;
00190     
00191     static const IndexedType& gather(const T& vec, std::size_t i);
00192     
00193     static void scatter(T& vec, const IndexedType& v, std::size_t i);
00194     
00195   };
00196 
00208   template<typename T>
00209   class DatatypeCommunicator : public InterfaceBuilder<T>
00210   {
00211   public:    
00212 
00216     typedef T ParallelIndexSet;
00217 
00221     typedef Dune::RemoteIndices<ParallelIndexSet> RemoteIndices;
00222 
00226     typedef typename RemoteIndices::GlobalIndex GlobalIndex;
00227     
00231     typedef typename RemoteIndices::Attribute Attribute;
00232 
00236     typedef typename RemoteIndices::LocalIndex LocalIndex;
00237     
00241     DatatypeCommunicator();
00242     
00246     ~DatatypeCommunicator();
00247 
00274     template<class T1, class T2, class V>
00275     void build(const RemoteIndices& remoteIndices, const T1& sourceFlags, V& sendData, const T2& destFlags, V& receiveData);
00276     
00280     void forward();
00281     
00285     void backward();
00286    
00290     void free();        
00291   private:
00292     enum { 
00296       commTag_ = 234
00297     };
00298     
00302     const RemoteIndices* remoteIndices_;
00303     
00304     typedef std::map<int,std::pair<MPI_Datatype,MPI_Datatype> > 
00305     MessageTypeMap;
00306     
00310     MessageTypeMap messageTypes;
00311     
00315     void* data_;
00316 
00317     MPI_Request* requests_[2];
00318 
00322     bool created_;
00323     
00327     template<class V, bool FORWARD>
00328     void createRequests(V& sendData, V& receiveData);
00329     
00333     template<class T1, class T2, class V, bool send>
00334     void createDataTypes(const T1& source, const T2& destination, V& data);
00335 
00339     void sendRecv(MPI_Request* req);
00340     
00344     struct IndexedTypeInformation
00345     {
00351       void build(int i)
00352       {
00353         length = new int[i];
00354         displ  = new MPI_Aint[i];
00355         size = i;
00356       }
00357       
00361       void free()
00362       {
00363         delete[] length;
00364         delete[] displ;
00365       }
00367       int* length;
00369       MPI_Aint* displ;
00375       int elements;
00379       int size;
00380     };
00381     
00387     template<class V>
00388     struct MPIDatatypeInformation
00389     {
00394       MPIDatatypeInformation(const V& data): data_(data)
00395       {}
00396       
00402       void reserve(int proc, int size)
00403       {
00404         information_[proc].build(size);
00405       }
00412       void add(int proc, int local)
00413       {
00414         IndexedTypeInformation& info=information_[proc];
00415         assert(info.elements<info.size);
00416         MPI_Address( const_cast<void*>(CommPolicy<V>::getAddress(data_, local)),
00417                     info.displ+info.elements);
00418         info.length[info.elements]=CommPolicy<V>::getSize(data_, local);
00419         info.elements++;
00420       }
00421       
00426       std::map<int,IndexedTypeInformation> information_;
00430       const V& data_;
00431       
00432     };
00433     
00434   };
00435 
00445   template<typename T>
00446   class BufferedCommunicator
00447   {
00448   public:
00453     typedef T PIndexSet;
00454     
00458     typedef Dune::Interface<PIndexSet> Interface;
00459     
00463     typedef typename Interface::GlobalIndex GlobalIndex;
00464     
00468     typedef typename Interface::Attribute Attribute;
00469     
00473     BufferedCommunicator();
00474     
00481     template<class Data>
00482     typename EnableIf<is_same<SizeOne,typename CommPolicy<Data>::IndexedTypeFlag>::value, void>::Type
00483     build(const Interface& interface);
00484 
00492     template<class Data>
00493     void build(const Data& source, const Data& target, const Interface& interface);
00494     
00523     template<class GatherScatter, class Data>
00524     void forward(const Data& source, Data& dest);
00525     
00554     template<class GatherScatter, class Data>
00555     void backward(Data& source, const Data& dest);
00556 
00582     template<class GatherScatter, class Data>
00583     void forward(Data& data);
00584     
00610     template<class GatherScatter, class Data>
00611     void backward(Data& data);
00612     
00616     void free();
00617     
00621     ~BufferedCommunicator();
00622     
00623   private:
00624 
00628     template<class Data, typename IndexedTypeFlag>
00629     struct MessageSizeCalculator
00630     {};
00631     
00636     template<class Data>
00637     struct MessageSizeCalculator<Data,SizeOne>
00638     {
00645       inline int operator()(const InterfaceInformation& info) const;
00654       inline int operator()(const Data& data, const InterfaceInformation& info) const;
00655     };
00656 
00661     template<class Data>
00662     struct MessageSizeCalculator<Data,VariableSize>
00663     {
00672       inline int operator()(const Data& data, const InterfaceInformation& info) const;
00673     };
00674     
00678     template<class Data, class GatherScatter, bool send, typename IndexedTypeFlag>
00679     struct MessageGatherer
00680     {};
00681     
00686     template<class Data, class GatherScatter, bool send>
00687     struct MessageGatherer<Data,GatherScatter,send,SizeOne>
00688     {
00690       typedef typename CommPolicy<Data>::IndexedType Type;
00691       
00696       typedef GatherScatter Gatherer;
00697       
00698       enum{
00704         forward=send
00705       };
00706       
00714       inline void operator()(const Interface& interface, const Data& data, Type* buffer, size_t bufferSize) const;
00715     };
00716 
00721     template<class Data, class GatherScatter, bool send>
00722     struct MessageGatherer<Data,GatherScatter,send,VariableSize>
00723     {
00725       typedef typename CommPolicy<Data>::IndexedType Type;
00726       
00731       typedef GatherScatter Gatherer;
00732       
00733       enum{
00739         forward=send
00740       };
00741       
00749       inline void operator()(const Interface& interface, const Data& data, Type* buffer, size_t bufferSize) const;
00750     };
00751 
00755     template<class Data, class GatherScatter, bool send, typename IndexedTypeFlag>
00756     struct MessageScatterer
00757     {};
00758 
00763     template<class Data, class GatherScatter, bool send>
00764     struct MessageScatterer<Data,GatherScatter,send,SizeOne>
00765     {
00767       typedef typename CommPolicy<Data>::IndexedType Type;
00768       
00773       typedef GatherScatter Scatterer;
00774       
00775       enum{
00781         forward=send
00782       };
00783       
00791       inline void operator()(const Interface& interface, Data& data, Type* buffer, const int& proc) const;
00792     };
00797     template<class Data, class GatherScatter, bool send>
00798     struct MessageScatterer<Data,GatherScatter,send,VariableSize>
00799     {
00801       typedef typename CommPolicy<Data>::IndexedType Type;
00802       
00807       typedef GatherScatter Scatterer;
00808       
00809       enum{
00815         forward=send
00816       };
00817       
00825       inline void operator()(const Interface& interface, Data& data, Type* buffer, const int& proc) const;
00826     };
00827 
00831     struct MessageInformation
00832     {
00834       MessageInformation()
00835         : start_(0), size_(0)
00836       {}
00837       
00845       MessageInformation(size_t start, size_t size)
00846         :start_(start), size_(size)
00847       {}
00851       size_t start_;
00855       size_t size_;
00856     };
00857 
00864     typedef std::map<int,std::pair<MessageInformation,MessageInformation> >
00865     InformationMap;
00869     InformationMap messageInformation_;
00873     char* buffers_[2];
00877     size_t bufferSize_[2];
00878     
00879     enum{
00883       commTag_
00884     };
00885     
00889     const Interface* interface_;
00890 
00894     template<class GatherScatter, bool FORWARD, class Data>
00895     void sendRecv(const Data& source, Data& target);
00896     
00897   };
00898   
00899   template<class V>
00900   inline const void* CommPolicy<V>::getAddress(const V& v, int index)
00901   {
00902     return &(v[index]);
00903   }
00904   
00905   template<class V>
00906   inline int CommPolicy<V>::getSize(const V& v, int index)
00907   {
00908     return 1;
00909   }
00910 
00911   template<class K, class A, int n>
00912   inline const void* CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >::getAddress(const Type& v, int index)
00913   {
00914     return &(v[index][0]);
00915   }
00916 
00917   template<class K, class A, int n>
00918   inline int CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >::getSize(const Type& v, int index)
00919   {
00920     return v[index].getsize();
00921   }
00922 
00923   
00924   template<class T>
00925   inline const typename CopyGatherScatter<T>::IndexedType& CopyGatherScatter<T>::gather(const T& vec, std::size_t i)
00926   {
00927     return vec[i];
00928     
00929   }
00930   
00931   template<class T>
00932   inline void CopyGatherScatter<T>::scatter(T& vec, const IndexedType& v, std::size_t i)
00933   {
00934     vec[i]=v;
00935     
00936   }
00937   template<typename T>
00938   DatatypeCommunicator<T>::DatatypeCommunicator()
00939     : remoteIndices_(0), created_(false)
00940   {
00941     requests_[0]=0;
00942     requests_[1]=0;
00943   }
00944   
00945  
00946 
00947   template<typename T>
00948   DatatypeCommunicator<T>::~DatatypeCommunicator()
00949   {
00950     free();
00951   }
00952   
00953   template<typename T>
00954   template<class T1, class T2, class V>
00955   inline void DatatypeCommunicator<T>::build(const RemoteIndices& remoteIndices,
00956                                                  const T1& source, V& sendData,
00957                                                  const T2& destination, V& receiveData)
00958   {
00959     remoteIndices_ = &remoteIndices;
00960     free();
00961     createDataTypes<T1,T2,V,false>(source,destination, receiveData);
00962     createDataTypes<T1,T2,V,true>(source,destination, sendData);
00963     createRequests<V,true>(sendData, receiveData);
00964     createRequests<V,false>(receiveData, sendData);
00965     created_=true;
00966   }
00967 
00968   template<typename T>
00969   void DatatypeCommunicator<T>::free()
00970   {
00971     if(created_){
00972       delete[] requests_[0];
00973       delete[] requests_[1];
00974       typedef MessageTypeMap::iterator iterator;
00975       typedef MessageTypeMap::const_iterator const_iterator;
00976       
00977       const const_iterator end=messageTypes.end();
00978       
00979       for(iterator process = messageTypes.begin(); process != end; ++process){
00980         MPI_Datatype *type = &(process->second.first);
00981         int finalized=0;
00982 #if MPI_2
00983         MPI_Finalized(&finalized);
00984 #endif
00985         if(*type!=MPI_DATATYPE_NULL && !finalized)
00986           MPI_Type_free(type);
00987         type = &(process->second.second);
00988         if(*type!=MPI_DATATYPE_NULL && !finalized)
00989           MPI_Type_free(type);
00990       }
00991       messageTypes.clear();
00992       created_=false;
00993     }
00994     
00995   }
00996   
00997   template<typename T>
00998   template<class T1, class T2, class V, bool send>
00999   void DatatypeCommunicator<T>::createDataTypes(const T1& sourceFlags, const T2& destFlags, V& data)
01000   {
01001 
01002     MPIDatatypeInformation<V>  dataInfo(data);
01003     this->template buildInterface<T1,T2,MPIDatatypeInformation<V>,send>(*remoteIndices_,sourceFlags, destFlags, dataInfo);
01004 
01005     typedef typename RemoteIndices::RemoteIndexMap::const_iterator const_iterator;
01006     const const_iterator end=this->remoteIndices_->end();
01007     
01008     // Allocate MPI_Datatypes and deallocate memory for the type construction.
01009     for(const_iterator process=this->remoteIndices_->begin(); process != end; ++process){
01010       IndexedTypeInformation& info=dataInfo.information_[process->first];
01011       // Shift the displacement
01012       MPI_Aint base;
01013       MPI_Address(const_cast<void *>(CommPolicy<V>::getAddress(data, 0)), &base);
01014       
01015       for(int i=0; i< info.elements; i++){
01016         info.displ[i]-=base;
01017       }
01018       
01019       // Create data type
01020       MPI_Datatype* type = &( send ? messageTypes[process->first].first : messageTypes[process->first].second);
01021       MPI_Type_hindexed(info.elements, info.length, info.displ, 
01022                        MPITraits<typename CommPolicy<V>::IndexedType>::getType(),
01023                        type);
01024       MPI_Type_commit(type);
01025       // Deallocate memory
01026       info.free();
01027     }
01028   }
01029   
01030   template<typename T>
01031   template<class V, bool createForward>
01032   void DatatypeCommunicator<T>::createRequests(V& sendData, V& receiveData)
01033   {
01034     typedef std::map<int,std::pair<MPI_Datatype,MPI_Datatype> >::const_iterator MapIterator;
01035     int rank;
01036     static int index = createForward?1:0;
01037     int noMessages = messageTypes.size();
01038     // allocate request handles
01039     requests_[index] = new MPI_Request[2*noMessages];
01040     const MapIterator end = messageTypes.end();
01041     int request=0;
01042     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01043     
01044     // Set up the requests for receiving first
01045     for(MapIterator process = messageTypes.begin(); process != end; 
01046         ++process, ++request){
01047       MPI_Datatype type = createForward ? process->second.second : process->second.first;
01048       void* address = const_cast<void*>(CommPolicy<V>::getAddress(receiveData,0));
01049       MPI_Recv_init(address, 1, type, process->first, commTag_, this->remoteIndices_->communicator(), requests_[index]+request);
01050     }
01051     
01052     // And now the send requests
01053 
01054     for(MapIterator process = messageTypes.begin(); process != end; 
01055         ++process, ++request){
01056       MPI_Datatype type = createForward ? process->second.first : process->second.second;
01057       void* address =  const_cast<void*>(CommPolicy<V>::getAddress(sendData, 0));
01058       MPI_Ssend_init(address, 1, type, process->first, commTag_, this->remoteIndices_->communicator(), requests_[index]+request);
01059     }
01060   }   
01061         
01062   template<typename T>
01063   void DatatypeCommunicator<T>::forward()
01064   {
01065     sendRecv(requests_[1]);
01066   }
01067   
01068   template<typename T>
01069   void DatatypeCommunicator<T>::backward()
01070   {
01071     sendRecv(requests_[0]);
01072   }
01073 
01074   template<typename T>
01075   void DatatypeCommunicator<T>::sendRecv(MPI_Request* requests)
01076   {
01077     int noMessages = messageTypes.size();
01078     // Start the receive calls first
01079     MPI_Startall(noMessages, requests);
01080     // Now the send calls
01081     MPI_Startall(noMessages, requests+noMessages);
01082     
01083     // Wait for completion of the communication send first then receive
01084     MPI_Status* status=new MPI_Status[2*noMessages];
01085     for(int i=0; i<2*noMessages; i++)
01086       status[i].MPI_ERROR=MPI_SUCCESS;
01087     
01088     int send = MPI_Waitall(noMessages, requests+noMessages, status+noMessages);
01089     int receive = MPI_Waitall(noMessages, requests, status);
01090     
01091     // Error checks
01092     int success=1, globalSuccess=0;
01093     if(send==MPI_ERR_IN_STATUS){
01094       int rank;
01095       MPI_Comm_rank(this->remoteIndices_->communicator(), &rank);
01096       std::cerr<<rank<<": Error in sending :"<<std::endl;
01097       // Search for the error
01098       for(int i=noMessages; i< 2*noMessages; i++)
01099         if(status[i].MPI_ERROR!=MPI_SUCCESS){
01100           char message[300];
01101           int messageLength;
01102           MPI_Error_string(status[i].MPI_ERROR, message, &messageLength);
01103           std::cerr<<" source="<<status[i].MPI_SOURCE<<" message: ";
01104           for(int i=0; i< messageLength; i++)
01105             std::cout<<message[i];
01106         }
01107       std::cerr<<std::endl;
01108       success=0;
01109     }
01110     
01111     if(receive==MPI_ERR_IN_STATUS){
01112       int rank;
01113       MPI_Comm_rank(this->remoteIndices_->communicator(), &rank);
01114       std::cerr<<rank<<": Error in receiving!"<<std::endl;
01115       // Search for the error
01116       for(int i=0; i< noMessages; i++)
01117         if(status[i].MPI_ERROR!=MPI_SUCCESS){
01118           char message[300];
01119           int messageLength;
01120           MPI_Error_string(status[i].MPI_ERROR, message, &messageLength);
01121           std::cerr<<" source="<<status[i].MPI_SOURCE<<" message: ";
01122           for(int i=0; i< messageLength; i++)
01123             std::cerr<<message[i];
01124         }
01125       std::cerr<<std::endl;
01126       success=0;
01127     }
01128 
01129     MPI_Allreduce(&success, &globalSuccess, 1, MPI_INT, MPI_MIN, this->remoteIndices_->communicator());
01130     
01131     delete[] status;
01132 
01133     if(!globalSuccess)
01134       DUNE_THROW(CommunicationError, "A communication error occurred!");
01135       
01136   }
01137   
01138 
01139 
01140   template<typename T>
01141   BufferedCommunicator<T>::BufferedCommunicator()
01142     : interface_(0)
01143   {
01144     buffers_[0]=0;
01145     buffers_[1]=0;
01146     bufferSize_[0]=0;
01147     bufferSize_[1]=0;
01148   }
01149   
01150   
01151   template<typename T>
01152   template<class Data>
01153   typename EnableIf<is_same<SizeOne, typename CommPolicy<Data>::IndexedTypeFlag>::value, void>::Type
01154   BufferedCommunicator<T>::build(const Interface& interface)
01155   {
01156     typedef typename Interface::InformationMap::const_iterator const_iterator;
01157     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01158     const const_iterator end = interface.interfaces().end();
01159     int lrank;
01160     MPI_Comm_rank(interface.communicator(), &lrank);
01161     
01162     bufferSize_[0]=0;
01163     bufferSize_[1]=0;
01164 
01165     for(const_iterator interfacePair = interface.interfaces().begin();
01166         interfacePair != end; ++interfacePair){  
01167       int noSend = MessageSizeCalculator<Data,Flag>()(interfacePair->second.first);
01168       int noRecv = MessageSizeCalculator<Data,Flag>()(interfacePair->second.second);
01169       messageInformation_.insert(std::make_pair(interfacePair->first,
01170                                  std::make_pair(MessageInformation(bufferSize_[0],
01171                                                                    noSend*sizeof(typename CommPolicy<Data>::IndexedType)),
01172                                                 MessageInformation(bufferSize_[1],
01173                                                                    noRecv*sizeof(typename CommPolicy<Data>::IndexedType)))));
01174       bufferSize_[0] += noSend;
01175       bufferSize_[1] += noRecv;
01176     }
01177 
01178     // allocate the buffers
01179     bufferSize_[0] *= sizeof(typename CommPolicy<Data>::IndexedType);
01180     bufferSize_[1] *= sizeof(typename CommPolicy<Data>::IndexedType);
01181     
01182     buffers_[0] = new char[bufferSize_[0]];
01183     buffers_[1] = new char[bufferSize_[1]];
01184     interface_ = &interface;
01185     
01186   }  
01187   
01188   template<typename T>
01189   template<class Data>
01190   void BufferedCommunicator<T>::build(const Data& source, const Data& dest, const Interface& interface)
01191   {
01192     typedef typename Interface::InformationMap::const_iterator const_iterator;
01193     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01194     const const_iterator end = interface.interfaces().end();
01195     
01196     bufferSize_[0]=0;
01197     bufferSize_[1]=0;
01198 
01199     for(const_iterator interfacePair = interface.interfaces().begin();
01200         interfacePair != end; ++interfacePair){      
01201       int noSend = MessageSizeCalculator<Data,Flag>()(source, interfacePair->second.first);
01202       int noRecv = MessageSizeCalculator<Data,Flag>()(dest, interfacePair->second.second);
01203       
01204       messageInformation_.insert(std::make_pair(interfacePair->first,
01205                                                 std::make_pair(MessageInformation(bufferSize_[0],
01206                                                                                   noSend*sizeof(typename CommPolicy<Data>::IndexedType)),
01207                                                                MessageInformation(bufferSize_[1],
01208                                                                                   noRecv*sizeof(typename CommPolicy<Data>::IndexedType)))));
01209       bufferSize_[0] += noSend;
01210       bufferSize_[1] += noRecv;
01211     }
01212 
01213     bufferSize_[0] *= sizeof(typename CommPolicy<Data>::IndexedType);
01214     bufferSize_[1] *= sizeof(typename CommPolicy<Data>::IndexedType);
01215     // allocate the buffers
01216     buffers_[0] = new char[bufferSize_[0]];
01217     buffers_[1] = new char[bufferSize_[1]];
01218     interface_ = &interface;
01219   }
01220   
01221   template<typename T>
01222   void BufferedCommunicator<T>::free()
01223   {
01224     if(interface_!=0){
01225       messageInformation_.clear();
01226       delete[] buffers_[0];
01227       delete[] buffers_[1];
01228       interface_=0;
01229     }
01230   }
01231 
01232   template<typename T>
01233   BufferedCommunicator<T>::~BufferedCommunicator()
01234   {
01235     free();
01236   }
01237   
01238   template<typename T>
01239   template<class Data>
01240   inline int BufferedCommunicator<T>::MessageSizeCalculator<Data,SizeOne>::operator()
01241     (const InterfaceInformation& info) const
01242   {
01243     return info.size();
01244   }
01245 
01246   template<typename T>
01247   template<class Data>
01248   inline int BufferedCommunicator<T>::MessageSizeCalculator<Data,SizeOne>::operator()
01249     (const Data& data, const InterfaceInformation& info) const
01250   {
01251     return operator()(info);
01252   }
01253 
01254   template<typename T>
01255   template<class Data>
01256   inline int BufferedCommunicator<T>::MessageSizeCalculator<Data, VariableSize>::operator()
01257     (const Data& data, const InterfaceInformation& info) const
01258   {
01259     int entries=0;
01260 
01261     for(int i=0; i <  info.size(); i++)
01262       entries += CommPolicy<Data>::getSize(data,info[i]);
01263 
01264     return entries;
01265   }
01266   
01267   template<typename T>
01268   template<class Data, class GatherScatter, bool FORWARD>
01269   inline void BufferedCommunicator<T>::MessageGatherer<Data,GatherScatter,FORWARD,VariableSize>::operator()(const Interface& interface,const Data& data, Type* buffer, size_t bufferSize) const
01270   {
01271     typedef typename Interface::InformationMap::const_iterator
01272       const_iterator;
01273 
01274     int rank;
01275     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01276     const const_iterator end = interface.interfaces().end();
01277     size_t index=0;
01278     
01279     for(const_iterator interfacePair = interface.interfaces().begin();
01280         interfacePair != end; ++interfacePair){
01281       int size = forward ? interfacePair->second.first.size() :
01282         interfacePair->second.second.size();
01283       int proc = interfacePair->first();
01284       
01285       for(int i=0; i < size; i++){  
01286         int local = forward ? interfacePair->second->first[i] :
01287           interfacePair->second->second[i];
01288         for(int j=0; j < CommPolicy<Data>::getSize(data, local);j++, index++){
01289 
01290 #ifdef DUNE_ISTL_WITH_CHECKING
01291           assert(bufferSize>=(index+1)*sizeof(typename CommPolicy<Data>::IndexedType));
01292 #endif
01293           buffer[index]=GatherScatter::gather(data, local, j);
01294         }
01295         
01296       }
01297     }
01298         
01299   }
01300 
01301   template<typename T>
01302   template<class Data, class GatherScatter, bool FORWARD>
01303   inline void BufferedCommunicator<T>::MessageGatherer<Data,GatherScatter,FORWARD,SizeOne>::operator()(const Interface& interface, const Data& data, Type* buffer, size_t bufferSize)const
01304   {
01305     typedef typename Interface::InformationMap::const_iterator
01306       const_iterator;
01307     const const_iterator end = interface.interfaces().end();
01308     size_t index = 0;
01309     
01310     int rank;
01311     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01312 
01313     for(const_iterator interfacePair = interface.interfaces().begin();
01314         interfacePair != end; ++interfacePair){
01315       size_t size = FORWARD ? interfacePair->second.first.size() :
01316         interfacePair->second.second.size();
01317       
01318       for(size_t i=0; i < size; i++){
01319 
01320 #ifdef DUNE_ISTL_WITH_CHECKING
01321           assert(bufferSize>=(index+1)*sizeof(typename CommPolicy<Data>::IndexedType));
01322 #endif
01323 
01324         buffer[index++] = GatherScatter::gather(data, FORWARD ? interfacePair->second.first[i] :
01325                                                 interfacePair->second.second[i]);
01326       }
01327     }
01328         
01329   }
01330   
01331   template<typename T>
01332   template<class Data, class GatherScatter, bool FORWARD>
01333   inline void BufferedCommunicator<T>::MessageScatterer<Data,GatherScatter,FORWARD,VariableSize>::operator()(const Interface& interface, Data& data, Type* buffer, const int& proc)const
01334   {
01335     typedef typename Interface::Information Information;
01336     const typename Interface::InformationMap::const_iterator infoPair = interface.interfaces().find(proc);
01337 
01338     assert(infoPair!=interface.interfaces().end());
01339     
01340     const Information& info = FORWARD ? infoPair->second.second :
01341        infoPair->second.first;
01342 
01343     for(int i=0, index=0; i < info.size(); i++){
01344         for(int j=0; j < CommPolicy<Data>::getSize(data, info[i]); j++)
01345           GatherScatter::scatter(data, buffer[index++], info[i], j);
01346     }
01347   }
01348 
01349   template<typename T>
01350   template<class Data, class GatherScatter, bool FORWARD>
01351   inline void BufferedCommunicator<T>::MessageScatterer<Data,GatherScatter,FORWARD,SizeOne>::operator()(const Interface& interface, Data& data, Type* buffer, const int& proc)const
01352   {
01353     typedef typename Interface::Information Information;
01354     const typename Interface::InformationMap::const_iterator infoPair = interface.interfaces().find(proc);
01355 
01356     assert(infoPair!=interface.interfaces().end());
01357     
01358     const Information& info = FORWARD ? infoPair->second.second :
01359       infoPair->second.first;
01360     
01361     for(size_t i=0; i < info.size(); i++){
01362       GatherScatter::scatter(data, buffer[i], info[i]);
01363     }
01364   }
01365   
01366   template<typename T>
01367   template<class GatherScatter,class Data>
01368   void BufferedCommunicator<T>::forward(Data& data)
01369   {
01370     this->template sendRecv<GatherScatter,true>(data, data);
01371   }
01372   
01373   template<typename T>
01374   template<class GatherScatter, class Data>
01375   void BufferedCommunicator<T>::backward(Data& data)
01376   {
01377     this->template sendRecv<GatherScatter,false>(data, data);
01378   }
01379 
01380   template<typename T>
01381   template<class GatherScatter, class Data>
01382   void BufferedCommunicator<T>::forward(const Data& source, Data& dest)
01383   {
01384     this->template sendRecv<GatherScatter,true>(source, dest);
01385   }
01386   
01387   template<typename T>
01388   template<class GatherScatter, class Data>
01389   void BufferedCommunicator<T>::backward(Data& source, const Data& dest)
01390   {
01391     this->template sendRecv<GatherScatter,false>(dest, source);
01392   }
01393 
01394   template<typename T>
01395   template<class GatherScatter, bool FORWARD, class Data>
01396   void BufferedCommunicator<T>::sendRecv(const Data& source, Data& dest) 
01397   {
01398     int rank, lrank;
01399     
01400     MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01401     MPI_Comm_rank(MPI_COMM_WORLD,&lrank);
01402     
01403     typedef typename CommPolicy<Data>::IndexedType Type;
01404     Type *sendBuffer, *recvBuffer;
01405     size_t sendBufferSize, recvBufferSize;
01406     
01407     if(FORWARD){
01408       sendBuffer = reinterpret_cast<Type*>(buffers_[0]);
01409       sendBufferSize = bufferSize_[0];
01410       recvBuffer = reinterpret_cast<Type*>(buffers_[1]);
01411       recvBufferSize = bufferSize_[1];
01412     }else{
01413       sendBuffer = reinterpret_cast<Type*>(buffers_[1]);
01414       sendBufferSize = bufferSize_[1];
01415       recvBuffer = reinterpret_cast<Type*>(buffers_[0]);
01416       recvBufferSize = bufferSize_[0];
01417     }
01418     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01419 
01420     MessageGatherer<Data,GatherScatter,FORWARD,Flag>()(*interface_, source, sendBuffer, sendBufferSize);
01421     
01422     MPI_Request* sendRequests = new MPI_Request[messageInformation_.size()];
01423     MPI_Request* recvRequests = new MPI_Request[messageInformation_.size()];
01424     
01425     // Setup receive first
01426     typedef typename InformationMap::const_iterator const_iterator;
01427 
01428     const const_iterator end = messageInformation_.end();
01429     size_t i=0;
01430     int* processMap = new int[messageInformation_.size()];
01431     
01432     for(const_iterator info = messageInformation_.begin(); info != end; ++info, ++i){
01433           processMap[i]=info->first;
01434           if(FORWARD){
01435             assert(info->second.second.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.second.size_ <= recvBufferSize );
01436             MPI_Irecv(recvBuffer+info->second.second.start_, info->second.second.size_,
01437                       MPI_BYTE, info->first, commTag_, interface_->communicator(),
01438                       recvRequests+i);
01439           }else{
01440             assert(info->second.first.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.first.size_ <= recvBufferSize );
01441             MPI_Irecv(recvBuffer+info->second.first.start_, info->second.first.size_,
01442                       MPI_BYTE, info->first, commTag_, interface_->communicator(),
01443                       recvRequests+i);
01444           }
01445     }
01446     
01447     // now the send requests
01448     i=0;
01449     for(const_iterator info = messageInformation_.begin(); info != end; ++info, ++i)
01450       if(FORWARD){
01451         assert(info->second.first.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.first.size_ <= sendBufferSize );
01452         MPI_Issend(sendBuffer+info->second.first.start_, info->second.first.size_,
01453                    MPI_BYTE, info->first, commTag_, interface_->communicator(),
01454                    sendRequests+i);
01455       }else{
01456         assert(info->second.second.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.second.size_ <= sendBufferSize );
01457         MPI_Issend(sendBuffer+info->second.second.start_, info->second.second.size_,
01458                    MPI_BYTE, info->first, commTag_, interface_->communicator(),
01459                    sendRequests+i);
01460       }
01461     
01462     // Wait for completion of receive and immediately start scatter
01463     i=0;
01464     int success = 1;
01465     int finished = MPI_UNDEFINED;
01466     MPI_Status status;//[messageInformation_.size()];
01467     //MPI_Waitall(messageInformation_.size(), recvRequests, status);
01468     
01469     for(i=0;i< messageInformation_.size();i++){
01470       status.MPI_ERROR=MPI_SUCCESS;
01471       MPI_Waitany(messageInformation_.size(), recvRequests, &finished, &status);
01472       assert(finished != MPI_UNDEFINED);
01473       
01474       if(status.MPI_ERROR==MPI_SUCCESS){
01475         int& proc = processMap[finished];
01476         typename InformationMap::const_iterator infoIter = messageInformation_.find(proc);
01477         assert(infoIter != messageInformation_.end());
01478 
01479         MessageInformation info = (FORWARD)? infoIter->second.second : infoIter->second.first;
01480         assert(info.start_+info.size_ <= recvBufferSize);
01481 
01482         MessageScatterer<Data,GatherScatter,FORWARD,Flag>()(*interface_, dest, recvBuffer+info.start_, proc);
01483       }else{
01484         std::cerr<<rank<<": MPI_Error occurred while receiving message from "<<processMap[finished]<<std::endl;
01485         success=0;
01486       }
01487     }
01488 
01489     MPI_Status recvStatus;
01490     
01491     // Wait for completion of sends
01492     for(i=0;i< messageInformation_.size();i++)
01493       if(MPI_SUCCESS!=MPI_Wait(sendRequests+i, &recvStatus)){
01494         std::cerr<<rank<<": MPI_Error occurred while sending message to "<<processMap[finished]<<std::endl;
01495         success=0;
01496       }
01497     /*
01498     int globalSuccess;
01499     MPI_Allreduce(&success, &globalSuccess, 1, MPI_INT, MPI_MIN, interface_->communicator());
01500 
01501     if(!globalSuccess)
01502       DUNE_THROW(CommunicationError, "A communication error occurred!");
01503     */
01504     delete[] processMap;
01505     delete[] sendRequests;
01506     delete[] recvRequests;
01507 
01508   }
01509 
01511 }
01512 
01513 #endif
01514 
01515 #endif

Generated on Thu Apr 2 10:40:13 2009 for dune-istl by  doxygen 1.5.6