communicator.hh

Go to the documentation of this file.
00001 // $Id: communicator.hh 1123 2009-11-13 11:20:55Z 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
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   class BufferedCommunicator
00446   {
00447     
00448   public:
00452     BufferedCommunicator();
00453     
00460     template<class Data, class Interface>
00461     typename enable_if<is_same<SizeOne,typename CommPolicy<Data>::IndexedTypeFlag>::value, void>::type
00462     build(const Interface& interface);
00463 
00471     template<class Data, class Interface>
00472     void build(const Data& source, const Data& target, const Interface& interface);
00473     
00502     template<class GatherScatter, class Data>
00503     void forward(const Data& source, Data& dest);
00504     
00533     template<class GatherScatter, class Data>
00534     void backward(Data& source, const Data& dest);
00535 
00561     template<class GatherScatter, class Data>
00562     void forward(Data& data);
00563     
00589     template<class GatherScatter, class Data>
00590     void backward(Data& data);
00591     
00595     void free();
00596     
00600     ~BufferedCommunicator();
00601     
00602   private:
00603     
00607     typedef std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >
00608     InterfaceMap;
00609 
00610 
00614     template<class Data, typename IndexedTypeFlag>
00615     struct MessageSizeCalculator
00616     {};
00617     
00622     template<class Data>
00623     struct MessageSizeCalculator<Data,SizeOne>
00624     {
00631       inline int operator()(const InterfaceInformation& info) const;
00640       inline int operator()(const Data& data, const InterfaceInformation& info) const;
00641     };
00642 
00647     template<class Data>
00648     struct MessageSizeCalculator<Data,VariableSize>
00649     {
00658       inline int operator()(const Data& data, const InterfaceInformation& info) const;
00659     };
00660     
00664     template<class Data, class GatherScatter, bool send, typename IndexedTypeFlag>
00665     struct MessageGatherer
00666     {};
00667     
00672     template<class Data, class GatherScatter, bool send>
00673     struct MessageGatherer<Data,GatherScatter,send,SizeOne>
00674     {
00676       typedef typename CommPolicy<Data>::IndexedType Type;
00677       
00682       typedef GatherScatter Gatherer;
00683       
00684       enum{
00690         forward=send
00691       };
00692       
00700       inline void operator()(const InterfaceMap& interface, const Data& data, Type* buffer, size_t bufferSize) const;
00701     };
00702     
00707     template<class Data, class GatherScatter, bool send>
00708     struct MessageGatherer<Data,GatherScatter,send,VariableSize>
00709     {
00711       typedef typename CommPolicy<Data>::IndexedType Type;
00712       
00717       typedef GatherScatter Gatherer;
00718       
00719       enum{
00725         forward=send
00726       };
00727       
00735       inline void operator()(const InterfaceMap& interface, const Data& data, Type* buffer, size_t bufferSize) const;
00736     };
00737 
00741     template<class Data, class GatherScatter, bool send, typename IndexedTypeFlag>
00742     struct MessageScatterer
00743     {};
00744 
00749     template<class Data, class GatherScatter, bool send>
00750     struct MessageScatterer<Data,GatherScatter,send,SizeOne>
00751     {
00753       typedef typename CommPolicy<Data>::IndexedType Type;
00754       
00759       typedef GatherScatter Scatterer;
00760       
00761       enum{
00767         forward=send
00768       };
00769       
00777       inline void operator()(const InterfaceMap& interface, Data& data, Type* buffer, const int& proc) const;
00778     };
00783     template<class Data, class GatherScatter, bool send>
00784     struct MessageScatterer<Data,GatherScatter,send,VariableSize>
00785     {
00787       typedef typename CommPolicy<Data>::IndexedType Type;
00788       
00793       typedef GatherScatter Scatterer;
00794       
00795       enum{
00801         forward=send
00802       };
00803       
00811       inline void operator()(const InterfaceMap& interface, Data& data, Type* buffer, const int& proc) const;
00812     };
00813 
00817     struct MessageInformation
00818     {
00820       MessageInformation()
00821         : start_(0), size_(0)
00822       {}
00823       
00831       MessageInformation(size_t start, size_t size)
00832         :start_(start), size_(size)
00833       {}
00837       size_t start_;
00841       size_t size_;
00842     };
00843 
00850     typedef std::map<int,std::pair<MessageInformation,MessageInformation> >
00851     InformationMap;
00855     InformationMap messageInformation_;
00859     char* buffers_[2];
00863     size_t bufferSize_[2];
00864     
00865     enum{
00869       commTag_
00870     };
00871     
00875     std::map<int,std::pair<InterfaceInformation,InterfaceInformation> > interfaces_;
00876 
00877     MPI_Comm communicator_;
00878 
00882     template<class GatherScatter, bool FORWARD, class Data>
00883     void sendRecv(const Data& source, Data& target);
00884     
00885   };
00886   
00887 #ifndef DOXYGEN
00888 
00889   template<class V>
00890   inline const void* CommPolicy<V>::getAddress(const V& v, int index)
00891   {
00892     return &(v[index]);
00893   }
00894   
00895   template<class V>
00896   inline int CommPolicy<V>::getSize(const V& v, int index)
00897   {
00898     return 1;
00899   }
00900 
00901   template<class K, class A, int n>
00902   inline const void* CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >::getAddress(const Type& v, int index)
00903   {
00904     return &(v[index][0]);
00905   }
00906 
00907   template<class K, class A, int n>
00908   inline int CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >::getSize(const Type& v, int index)
00909   {
00910     return v[index].getsize();
00911   }
00912 
00913   template<class T>
00914   inline const typename CopyGatherScatter<T>::IndexedType& CopyGatherScatter<T>::gather(const T& vec, std::size_t i)
00915   {
00916     return vec[i];
00917   }
00918 
00919   template<class T>
00920   inline void CopyGatherScatter<T>::scatter(T& vec, const IndexedType& v, std::size_t i)
00921   {
00922     vec[i]=v;
00923   }
00924 
00925   template<typename T>
00926   DatatypeCommunicator<T>::DatatypeCommunicator()
00927     : remoteIndices_(0), created_(false)
00928   {
00929     requests_[0]=0;
00930     requests_[1]=0;
00931   }
00932   
00933  
00934 
00935   template<typename T>
00936   DatatypeCommunicator<T>::~DatatypeCommunicator()
00937   {
00938     free();
00939   }
00940   
00941   template<typename T>
00942   template<class T1, class T2, class V>
00943   inline void DatatypeCommunicator<T>::build(const RemoteIndices& remoteIndices,
00944                                                  const T1& source, V& sendData,
00945                                                  const T2& destination, V& receiveData)
00946   {
00947     remoteIndices_ = &remoteIndices;
00948     free();
00949     createDataTypes<T1,T2,V,false>(source,destination, receiveData);
00950     createDataTypes<T1,T2,V,true>(source,destination, sendData);
00951     createRequests<V,true>(sendData, receiveData);
00952     createRequests<V,false>(receiveData, sendData);
00953     created_=true;
00954   }
00955 
00956   template<typename T>
00957   void DatatypeCommunicator<T>::free()
00958   {
00959     if(created_){
00960       delete[] requests_[0];
00961       delete[] requests_[1];
00962       typedef MessageTypeMap::iterator iterator;
00963       typedef MessageTypeMap::const_iterator const_iterator;
00964       
00965       const const_iterator end=messageTypes.end();
00966       
00967       for(iterator process = messageTypes.begin(); process != end; ++process){
00968         MPI_Datatype *type = &(process->second.first);
00969         int finalized=0;
00970 #if MPI_2
00971         MPI_Finalized(&finalized);
00972 #endif
00973         if(*type!=MPI_DATATYPE_NULL && !finalized)
00974           MPI_Type_free(type);
00975         type = &(process->second.second);
00976         if(*type!=MPI_DATATYPE_NULL && !finalized)
00977           MPI_Type_free(type);
00978       }
00979       messageTypes.clear();
00980       created_=false;
00981     }
00982     
00983   }
00984   
00985   template<typename T>
00986   template<class T1, class T2, class V, bool send>
00987   void DatatypeCommunicator<T>::createDataTypes(const T1& sourceFlags, const T2& destFlags, V& data)
00988   {
00989 
00990     MPIDatatypeInformation<V>  dataInfo(data);
00991     this->template buildInterface<RemoteIndices,T1,T2,MPIDatatypeInformation<V>,send>(*remoteIndices_,sourceFlags, destFlags, dataInfo);
00992 
00993     typedef typename RemoteIndices::RemoteIndexMap::const_iterator const_iterator;
00994     const const_iterator end=this->remoteIndices_->end();
00995     
00996     // Allocate MPI_Datatypes and deallocate memory for the type construction.
00997     for(const_iterator process=this->remoteIndices_->begin(); process != end; ++process){
00998       IndexedTypeInformation& info=dataInfo.information_[process->first];
00999       // Shift the displacement
01000       MPI_Aint base;
01001       MPI_Address(const_cast<void *>(CommPolicy<V>::getAddress(data, 0)), &base);
01002       
01003       for(int i=0; i< info.elements; i++){
01004         info.displ[i]-=base;
01005       }
01006       
01007       // Create data type
01008       MPI_Datatype* type = &( send ? messageTypes[process->first].first : messageTypes[process->first].second);
01009       MPI_Type_hindexed(info.elements, info.length, info.displ, 
01010                        MPITraits<typename CommPolicy<V>::IndexedType>::getType(),
01011                        type);
01012       MPI_Type_commit(type);
01013       // Deallocate memory
01014       info.free();
01015     }
01016   }
01017   
01018   template<typename T>
01019   template<class V, bool createForward>
01020   void DatatypeCommunicator<T>::createRequests(V& sendData, V& receiveData)
01021   {
01022     typedef std::map<int,std::pair<MPI_Datatype,MPI_Datatype> >::const_iterator MapIterator;
01023     int rank;
01024     static int index = createForward?1:0;
01025     int noMessages = messageTypes.size();
01026     // allocate request handles
01027     requests_[index] = new MPI_Request[2*noMessages];
01028     const MapIterator end = messageTypes.end();
01029     int request=0;
01030     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01031     
01032     // Set up the requests for receiving first
01033     for(MapIterator process = messageTypes.begin(); process != end; 
01034         ++process, ++request){
01035       MPI_Datatype type = createForward ? process->second.second : process->second.first;
01036       void* address = const_cast<void*>(CommPolicy<V>::getAddress(receiveData,0));
01037       MPI_Recv_init(address, 1, type, process->first, commTag_, this->remoteIndices_->communicator(), requests_[index]+request);
01038     }
01039     
01040     // And now the send requests
01041 
01042     for(MapIterator process = messageTypes.begin(); process != end; 
01043         ++process, ++request){
01044       MPI_Datatype type = createForward ? process->second.first : process->second.second;
01045       void* address =  const_cast<void*>(CommPolicy<V>::getAddress(sendData, 0));
01046       MPI_Ssend_init(address, 1, type, process->first, commTag_, this->remoteIndices_->communicator(), requests_[index]+request);
01047     }
01048   }   
01049         
01050   template<typename T>
01051   void DatatypeCommunicator<T>::forward()
01052   {
01053     sendRecv(requests_[1]);
01054   }
01055   
01056   template<typename T>
01057   void DatatypeCommunicator<T>::backward()
01058   {
01059     sendRecv(requests_[0]);
01060   }
01061 
01062   template<typename T>
01063   void DatatypeCommunicator<T>::sendRecv(MPI_Request* requests)
01064   {
01065     int noMessages = messageTypes.size();
01066     // Start the receive calls first
01067     MPI_Startall(noMessages, requests);
01068     // Now the send calls
01069     MPI_Startall(noMessages, requests+noMessages);
01070     
01071     // Wait for completion of the communication send first then receive
01072     MPI_Status* status=new MPI_Status[2*noMessages];
01073     for(int i=0; i<2*noMessages; i++)
01074       status[i].MPI_ERROR=MPI_SUCCESS;
01075     
01076     int send = MPI_Waitall(noMessages, requests+noMessages, status+noMessages);
01077     int receive = MPI_Waitall(noMessages, requests, status);
01078     
01079     // Error checks
01080     int success=1, globalSuccess=0;
01081     if(send==MPI_ERR_IN_STATUS){
01082       int rank;
01083       MPI_Comm_rank(this->remoteIndices_->communicator(), &rank);
01084       std::cerr<<rank<<": Error in sending :"<<std::endl;
01085       // Search for the error
01086       for(int i=noMessages; i< 2*noMessages; i++)
01087         if(status[i].MPI_ERROR!=MPI_SUCCESS){
01088           char message[300];
01089           int messageLength;
01090           MPI_Error_string(status[i].MPI_ERROR, message, &messageLength);
01091           std::cerr<<" source="<<status[i].MPI_SOURCE<<" message: ";
01092           for(int i=0; i< messageLength; i++)
01093             std::cout<<message[i];
01094         }
01095       std::cerr<<std::endl;
01096       success=0;
01097     }
01098     
01099     if(receive==MPI_ERR_IN_STATUS){
01100       int rank;
01101       MPI_Comm_rank(this->remoteIndices_->communicator(), &rank);
01102       std::cerr<<rank<<": Error in receiving!"<<std::endl;
01103       // Search for the error
01104       for(int i=0; i< noMessages; i++)
01105         if(status[i].MPI_ERROR!=MPI_SUCCESS){
01106           char message[300];
01107           int messageLength;
01108           MPI_Error_string(status[i].MPI_ERROR, message, &messageLength);
01109           std::cerr<<" source="<<status[i].MPI_SOURCE<<" message: ";
01110           for(int i=0; i< messageLength; i++)
01111             std::cerr<<message[i];
01112         }
01113       std::cerr<<std::endl;
01114       success=0;
01115     }
01116 
01117     MPI_Allreduce(&success, &globalSuccess, 1, MPI_INT, MPI_MIN, this->remoteIndices_->communicator());
01118     
01119     delete[] status;
01120 
01121     if(!globalSuccess)
01122       DUNE_THROW(CommunicationError, "A communication error occurred!");
01123       
01124   }
01125   
01126   inline BufferedCommunicator::BufferedCommunicator()
01127   {
01128     buffers_[0]=0;
01129     buffers_[1]=0;
01130     bufferSize_[0]=0;
01131     bufferSize_[1]=0;
01132   }
01133 
01134   template<class Data, class Interface>
01135   typename enable_if<is_same<SizeOne, typename CommPolicy<Data>::IndexedTypeFlag>::value, void>::type
01136   BufferedCommunicator::build(const Interface& interface)
01137   {
01138     interfaces_=interface.interfaces();
01139     communicator_=interface.communicator();
01140     typedef typename std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >
01141       ::const_iterator const_iterator;
01142     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01143     const const_iterator end = interfaces_.end();
01144     int lrank;
01145     MPI_Comm_rank(communicator_, &lrank);
01146     
01147     bufferSize_[0]=0;
01148     bufferSize_[1]=0;
01149 
01150     for(const_iterator interfacePair = interfaces_.begin();
01151         interfacePair != end; ++interfacePair){  
01152       int noSend = MessageSizeCalculator<Data,Flag>()(interfacePair->second.first);
01153       int noRecv = MessageSizeCalculator<Data,Flag>()(interfacePair->second.second);
01154       messageInformation_.insert(std::make_pair(interfacePair->first,
01155                                  std::make_pair(MessageInformation(bufferSize_[0],
01156                                                                    noSend*sizeof(typename CommPolicy<Data>::IndexedType)),
01157                                                 MessageInformation(bufferSize_[1],
01158                                                                    noRecv*sizeof(typename CommPolicy<Data>::IndexedType)))));
01159       bufferSize_[0] += noSend;
01160       bufferSize_[1] += noRecv;
01161     }
01162 
01163     // allocate the buffers
01164     bufferSize_[0] *= sizeof(typename CommPolicy<Data>::IndexedType);
01165     bufferSize_[1] *= sizeof(typename CommPolicy<Data>::IndexedType);
01166     
01167     buffers_[0] = new char[bufferSize_[0]];
01168     buffers_[1] = new char[bufferSize_[1]];    
01169   }  
01170   
01171   template<class Data, class Interface>
01172   void BufferedCommunicator::build(const Data& source, const Data& dest, const Interface& interface)
01173   {
01174     
01175     interfaces_=interface.interfaces();
01176     communicator_=interface.communicator();
01177     typedef typename std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >
01178       ::const_iterator const_iterator;
01179     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01180     const const_iterator end = interfaces_.end();
01181     
01182     bufferSize_[0]=0;
01183     bufferSize_[1]=0;
01184 
01185     for(const_iterator interfacePair = interfaces_.begin();
01186         interfacePair != end; ++interfacePair){      
01187       int noSend = MessageSizeCalculator<Data,Flag>()(source, interfacePair->second.first);
01188       int noRecv = MessageSizeCalculator<Data,Flag>()(dest, interfacePair->second.second);
01189       
01190       messageInformation_.insert(std::make_pair(interfacePair->first,
01191                                                 std::make_pair(MessageInformation(bufferSize_[0],
01192                                                                                   noSend*sizeof(typename CommPolicy<Data>::IndexedType)),
01193                                                                MessageInformation(bufferSize_[1],
01194                                                                                   noRecv*sizeof(typename CommPolicy<Data>::IndexedType)))));
01195       bufferSize_[0] += noSend;
01196       bufferSize_[1] += noRecv;
01197     }
01198 
01199     bufferSize_[0] *= sizeof(typename CommPolicy<Data>::IndexedType);
01200     bufferSize_[1] *= sizeof(typename CommPolicy<Data>::IndexedType);
01201     // allocate the buffers
01202     buffers_[0] = new char[bufferSize_[0]];
01203     buffers_[1] = new char[bufferSize_[1]];
01204   }
01205   
01206   inline void BufferedCommunicator::free()
01207   {
01208       messageInformation_.clear();
01209       if(buffers_[0])
01210         delete[] buffers_[0];
01211       
01212       if(buffers_[1])
01213         delete[] buffers_[1];
01214       buffers_[0]=buffers_[1]=0;
01215   }
01216 
01217   inline BufferedCommunicator::~BufferedCommunicator()
01218   {
01219     free();
01220   }
01221   
01222   template<class Data>
01223   inline int BufferedCommunicator::MessageSizeCalculator<Data,SizeOne>::operator()
01224     (const InterfaceInformation& info) const
01225   {
01226     return info.size();
01227   }
01228 
01229   
01230   template<class Data>
01231   inline int BufferedCommunicator::MessageSizeCalculator<Data,SizeOne>::operator()
01232     (const Data& data, const InterfaceInformation& info) const
01233   {
01234     return operator()(info);
01235   }
01236 
01237   
01238   template<class Data>
01239   inline int BufferedCommunicator::MessageSizeCalculator<Data, VariableSize>::operator()
01240     (const Data& data, const InterfaceInformation& info) const
01241   {
01242     int entries=0;
01243 
01244     for(size_t i=0; i < info.size(); i++)
01245       entries += CommPolicy<Data>::getSize(data,info[i]);
01246 
01247     return entries;
01248   }
01249   
01250   
01251   template<class Data, class GatherScatter, bool FORWARD>
01252   inline void BufferedCommunicator::MessageGatherer<Data,GatherScatter,FORWARD,VariableSize>::operator()(const InterfaceMap& interfaces,const Data& data, Type* buffer, size_t bufferSize) const
01253   {
01254     typedef typename InterfaceMap::const_iterator
01255       const_iterator;
01256 
01257     int rank;
01258     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01259     const const_iterator end = interfaces.end();
01260     size_t index=0;
01261     
01262     for(const_iterator interfacePair = interfaces.begin();
01263         interfacePair != end; ++interfacePair){
01264       int size = forward ? interfacePair->second.first.size() :
01265         interfacePair->second.second.size();
01266       
01267       for(int i=0; i < size; i++){  
01268         int local = forward ? interfacePair->second.first[i] :
01269           interfacePair->second.second[i];
01270         for(std::size_t j=0; j < CommPolicy<Data>::getSize(data, local);j++, index++){
01271 
01272 #ifdef DUNE_ISTL_WITH_CHECKING
01273           assert(bufferSize>=(index+1)*sizeof(typename CommPolicy<Data>::IndexedType));
01274 #endif
01275           buffer[index]=GatherScatter::gather(data, local, j);
01276         }
01277         
01278       }
01279     }
01280         
01281   }
01282 
01283   
01284   template<class Data, class GatherScatter, bool FORWARD>
01285   inline void BufferedCommunicator::MessageGatherer<Data,GatherScatter,FORWARD,SizeOne>::operator()(const InterfaceMap& interfaces, const Data& data, Type* buffer, size_t bufferSize)const
01286   {
01287     typedef typename InterfaceMap::const_iterator
01288       const_iterator;
01289     const const_iterator end = interfaces.end();
01290     size_t index = 0;
01291     
01292     int rank;
01293     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01294 
01295     for(const_iterator interfacePair = interfaces.begin();
01296         interfacePair != end; ++interfacePair){
01297       size_t size = FORWARD ? interfacePair->second.first.size() :
01298         interfacePair->second.second.size();
01299       
01300       for(size_t i=0; i < size; i++){
01301 
01302 #ifdef DUNE_ISTL_WITH_CHECKING
01303           assert(bufferSize>=(index+1)*sizeof(typename CommPolicy<Data>::IndexedType));
01304 #endif
01305 
01306         buffer[index++] = GatherScatter::gather(data, FORWARD ? interfacePair->second.first[i] :
01307                                                 interfacePair->second.second[i]);
01308       }
01309     }
01310         
01311   }
01312   
01313   
01314   template<class Data, class GatherScatter, bool FORWARD>
01315   inline void BufferedCommunicator::MessageScatterer<Data,GatherScatter,FORWARD,VariableSize>::operator()(const InterfaceMap& interfaces, Data& data, Type* buffer, const int& proc)const
01316   {
01317     typedef typename InterfaceMap::value_type::second_type::first_type Information;
01318     const typename InterfaceMap::const_iterator infoPair = interfaces.find(proc);
01319 
01320     assert(infoPair!=interfaces.end());
01321     
01322     const Information& info = FORWARD ? infoPair->second.second :
01323        infoPair->second.first;
01324 
01325     for(size_t i=0, index=0; i < info.size(); i++){
01326         for(size_t j=0; j < CommPolicy<Data>::getSize(data, info[i]); j++)
01327           GatherScatter::scatter(data, buffer[index++], info[i], j);
01328     }
01329   }
01330 
01331   
01332   template<class Data, class GatherScatter, bool FORWARD>
01333   inline void BufferedCommunicator::MessageScatterer<Data,GatherScatter,FORWARD,SizeOne>::operator()(const InterfaceMap& interfaces, Data& data, Type* buffer, const int& proc)const
01334   {
01335     typedef typename InterfaceMap::value_type::second_type::first_type Information;
01336     const typename InterfaceMap::const_iterator infoPair = interfaces.find(proc);
01337 
01338     assert(infoPair!=interfaces.end());
01339     
01340     const Information& info = FORWARD ? infoPair->second.second :
01341       infoPair->second.first;
01342     
01343     for(size_t i=0; i < info.size(); i++){
01344       GatherScatter::scatter(data, buffer[i], info[i]);
01345     }
01346   }
01347   
01348   
01349   template<class GatherScatter,class Data>
01350   void BufferedCommunicator::forward(Data& data)
01351   {
01352     this->template sendRecv<GatherScatter,true>(data, data);
01353   }
01354   
01355   
01356   template<class GatherScatter, class Data>
01357   void BufferedCommunicator::backward(Data& data)
01358   {
01359     this->template sendRecv<GatherScatter,false>(data, data);
01360   }
01361 
01362   
01363   template<class GatherScatter, class Data>
01364   void BufferedCommunicator::forward(const Data& source, Data& dest)
01365   {
01366     this->template sendRecv<GatherScatter,true>(source, dest);
01367   }
01368   
01369   
01370   template<class GatherScatter, class Data>
01371   void BufferedCommunicator::backward(Data& source, const Data& dest)
01372   {
01373     this->template sendRecv<GatherScatter,false>(dest, source);
01374   }
01375 
01376   
01377   template<class GatherScatter, bool FORWARD, class Data>
01378   void BufferedCommunicator::sendRecv(const Data& source, Data& dest) 
01379   {
01380     int rank, lrank;
01381     
01382     MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01383     MPI_Comm_rank(MPI_COMM_WORLD,&lrank);
01384     
01385     typedef typename CommPolicy<Data>::IndexedType Type;
01386     Type *sendBuffer, *recvBuffer;
01387     size_t sendBufferSize, recvBufferSize;
01388     
01389     if(FORWARD){
01390       sendBuffer = reinterpret_cast<Type*>(buffers_[0]);
01391       sendBufferSize = bufferSize_[0];
01392       recvBuffer = reinterpret_cast<Type*>(buffers_[1]);
01393       recvBufferSize = bufferSize_[1];
01394     }else{
01395       sendBuffer = reinterpret_cast<Type*>(buffers_[1]);
01396       sendBufferSize = bufferSize_[1];
01397       recvBuffer = reinterpret_cast<Type*>(buffers_[0]);
01398       recvBufferSize = bufferSize_[0];
01399     }
01400     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01401 
01402     MessageGatherer<Data,GatherScatter,FORWARD,Flag>()(interfaces_, source, sendBuffer, sendBufferSize);
01403     
01404     MPI_Request* sendRequests = new MPI_Request[messageInformation_.size()];
01405     MPI_Request* recvRequests = new MPI_Request[messageInformation_.size()];
01406     
01407     // Setup receive first
01408     typedef typename InformationMap::const_iterator const_iterator;
01409 
01410     const const_iterator end = messageInformation_.end();
01411     size_t i=0;
01412     int* processMap = new int[messageInformation_.size()];
01413     
01414     for(const_iterator info = messageInformation_.begin(); info != end; ++info, ++i){
01415           processMap[i]=info->first;
01416           if(FORWARD){
01417             assert(info->second.second.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.second.size_ <= recvBufferSize );
01418             MPI_Irecv(recvBuffer+info->second.second.start_, info->second.second.size_,
01419                       MPI_BYTE, info->first, commTag_, communicator_,
01420                       recvRequests+i);
01421           }else{
01422             assert(info->second.first.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.first.size_ <= recvBufferSize );
01423             MPI_Irecv(recvBuffer+info->second.first.start_, info->second.first.size_,
01424                       MPI_BYTE, info->first, commTag_, communicator_,
01425                       recvRequests+i);
01426           }
01427     }
01428     
01429     // now the send requests
01430     i=0;
01431     for(const_iterator info = messageInformation_.begin(); info != end; ++info, ++i)
01432       if(FORWARD){
01433         assert(info->second.first.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.first.size_ <= sendBufferSize );
01434         MPI_Issend(sendBuffer+info->second.first.start_, info->second.first.size_,
01435                    MPI_BYTE, info->first, commTag_, communicator_,
01436                    sendRequests+i);
01437       }else{
01438         assert(info->second.second.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.second.size_ <= sendBufferSize );
01439         MPI_Issend(sendBuffer+info->second.second.start_, info->second.second.size_,
01440                    MPI_BYTE, info->first, commTag_, communicator_,
01441                    sendRequests+i);
01442       }
01443     
01444     // Wait for completion of receive and immediately start scatter
01445     i=0;
01446     int success = 1;
01447     int finished = MPI_UNDEFINED;
01448     MPI_Status status;//[messageInformation_.size()];
01449     //MPI_Waitall(messageInformation_.size(), recvRequests, status);
01450     
01451     for(i=0;i< messageInformation_.size();i++){
01452       status.MPI_ERROR=MPI_SUCCESS;
01453       MPI_Waitany(messageInformation_.size(), recvRequests, &finished, &status);
01454       assert(finished != MPI_UNDEFINED);
01455       
01456       if(status.MPI_ERROR==MPI_SUCCESS){
01457         int& proc = processMap[finished];
01458         typename InformationMap::const_iterator infoIter = messageInformation_.find(proc);
01459         assert(infoIter != messageInformation_.end());
01460 
01461         MessageInformation info = (FORWARD)? infoIter->second.second : infoIter->second.first;
01462         assert(info.start_+info.size_ <= recvBufferSize);
01463 
01464         MessageScatterer<Data,GatherScatter,FORWARD,Flag>()(interfaces_, dest, recvBuffer+info.start_, proc);
01465       }else{
01466         std::cerr<<rank<<": MPI_Error occurred while receiving message from "<<processMap[finished]<<std::endl;
01467         success=0;
01468       }
01469     }
01470 
01471     MPI_Status recvStatus;
01472     
01473     // Wait for completion of sends
01474     for(i=0;i< messageInformation_.size();i++)
01475       if(MPI_SUCCESS!=MPI_Wait(sendRequests+i, &recvStatus)){
01476         std::cerr<<rank<<": MPI_Error occurred while sending message to "<<processMap[finished]<<std::endl;
01477         success=0;
01478       }
01479     /*
01480     int globalSuccess;
01481     MPI_Allreduce(&success, &globalSuccess, 1, MPI_INT, MPI_MIN, interface_->communicator());
01482 
01483     if(!globalSuccess)
01484       DUNE_THROW(CommunicationError, "A communication error occurred!");
01485     */
01486     delete[] processMap;
01487     delete[] sendRequests;
01488     delete[] recvRequests;
01489 
01490   }
01491 
01492 #endif  // DOXYGEN
01493   
01495 }
01496 
01497 #endif
01498 
01499 #endif
Generated on Sat Apr 24 11:13:45 2010 for dune-istl by  doxygen 1.6.3