dune-common 2.1.1
communicator.hh
Go to the documentation of this file.
00001 // $Id$
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 #include <dune/common/stdstreams.hh>
00010 
00011 #if HAVE_MPI
00012 // MPI header
00013 #include <mpi.h>
00014 
00015 namespace Dune
00016 {
00096   struct SizeOne
00097   {};
00098   
00104   struct VariableSize
00105   {};
00106   
00107 
00113   template<class V>
00114   struct CommPolicy
00115   {
00127     typedef V Type;
00128     
00134     typedef typename V::value_type IndexedType;
00135     
00140     typedef SizeOne IndexedTypeFlag;
00141 
00150     static const void* getAddress(const V& v, int index);
00151     
00157     static int getSize(const V&, int index);
00158   };
00159 
00160   template<class K, int n> class FieldVector;
00161   
00162   template<class B, class A> class VariableBlockVector;
00163   
00164   template<class K, class A, int n>
00165   struct CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >
00166   {
00167     typedef VariableBlockVector<FieldVector<K, n>, A> Type;
00168     
00169     typedef typename Type::B IndexedType;
00170     
00171     typedef VariableSize IndexedTypeFlag;
00172 
00173     static const void* getAddress(const Type& v, int i);
00174     
00175     static int getSize(const Type& v, int i);
00176   };
00177     
00181   class CommunicationError : public IOError
00182   {};
00183 
00187   template<class T>
00188   struct CopyGatherScatter
00189   {
00190     typedef typename CommPolicy<T>::IndexedType IndexedType;
00191     
00192     static const IndexedType& gather(const T& vec, std::size_t i);
00193     
00194     static void scatter(T& vec, const IndexedType& v, std::size_t i);
00195     
00196   };
00197 
00209   template<typename T>
00210   class DatatypeCommunicator : public InterfaceBuilder
00211   {
00212   public:    
00213 
00217     typedef T ParallelIndexSet;
00218 
00222     typedef Dune::RemoteIndices<ParallelIndexSet> RemoteIndices;
00223 
00227     typedef typename RemoteIndices::GlobalIndex GlobalIndex;
00228     
00232     typedef typename RemoteIndices::Attribute Attribute;
00233 
00237     typedef typename RemoteIndices::LocalIndex LocalIndex;
00238     
00242     DatatypeCommunicator();
00243     
00247     ~DatatypeCommunicator();
00248 
00275     template<class T1, class T2, class V>
00276     void build(const RemoteIndices& remoteIndices, const T1& sourceFlags, V& sendData, const T2& destFlags, V& receiveData);
00277     
00281     void forward();
00282     
00286     void backward();
00287    
00291     void free();        
00292   private:
00293     enum { 
00297       commTag_ = 234
00298     };
00299     
00303     const RemoteIndices* remoteIndices_;
00304     
00305     typedef std::map<int,std::pair<MPI_Datatype,MPI_Datatype> > 
00306     MessageTypeMap;
00307     
00311     MessageTypeMap messageTypes;
00312     
00316     void* data_;
00317 
00318     MPI_Request* requests_[2];
00319 
00323     bool created_;
00324     
00328     template<class V, bool FORWARD>
00329     void createRequests(V& sendData, V& receiveData);
00330     
00334     template<class T1, class T2, class V, bool send>
00335     void createDataTypes(const T1& source, const T2& destination, V& data);
00336 
00340     void sendRecv(MPI_Request* req);
00341     
00345     struct IndexedTypeInformation
00346     {
00352       void build(int i)
00353       {
00354         length = new int[i];
00355         displ  = new MPI_Aint[i];
00356         size = i;
00357       }
00358       
00362       void free()
00363       {
00364         delete[] length;
00365         delete[] displ;
00366       }
00368       int* length;
00370       MPI_Aint* displ;
00376       int elements;
00380       int size;
00381     };
00382     
00388     template<class V>
00389     struct MPIDatatypeInformation
00390     {
00395       MPIDatatypeInformation(const V& data): data_(data)
00396       {}
00397       
00403       void reserve(int proc, int size)
00404       {
00405         information_[proc].build(size);
00406       }
00413       void add(int proc, int local)
00414       {
00415         IndexedTypeInformation& info=information_[proc];
00416         assert(info.elements<info.size);
00417         MPI_Address( const_cast<void*>(CommPolicy<V>::getAddress(data_, local)),
00418                     info.displ+info.elements);
00419         info.length[info.elements]=CommPolicy<V>::getSize(data_, local);
00420         info.elements++;
00421       }
00422       
00427       std::map<int,IndexedTypeInformation> information_;
00431       const V& data_;
00432       
00433     };
00434     
00435   };
00436 
00446   class BufferedCommunicator
00447   {
00448     
00449   public:
00453     BufferedCommunicator();
00454     
00461     template<class Data, class Interface>
00462     typename enable_if<is_same<SizeOne,typename CommPolicy<Data>::IndexedTypeFlag>::value, void>::type
00463     build(const Interface& interface);
00464 
00472     template<class Data, class Interface>
00473     void build(const Data& source, const Data& target, const Interface& interface);
00474     
00503     template<class GatherScatter, class Data>
00504     void forward(const Data& source, Data& dest);
00505     
00534     template<class GatherScatter, class Data>
00535     void backward(Data& source, const Data& dest);
00536 
00562     template<class GatherScatter, class Data>
00563     void forward(Data& data);
00564     
00590     template<class GatherScatter, class Data>
00591     void backward(Data& data);
00592     
00596     void free();
00597     
00601     ~BufferedCommunicator();
00602     
00603   private:
00604     
00608     typedef std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >
00609     InterfaceMap;
00610 
00611 
00615     template<class Data, typename IndexedTypeFlag>
00616     struct MessageSizeCalculator
00617     {};
00618     
00623     template<class Data>
00624     struct MessageSizeCalculator<Data,SizeOne>
00625     {
00632       inline int operator()(const InterfaceInformation& info) const;
00641       inline int operator()(const Data& data, const InterfaceInformation& info) const;
00642     };
00643 
00648     template<class Data>
00649     struct MessageSizeCalculator<Data,VariableSize>
00650     {
00659       inline int operator()(const Data& data, const InterfaceInformation& info) const;
00660     };
00661     
00665     template<class Data, class GatherScatter, bool send, typename IndexedTypeFlag>
00666     struct MessageGatherer
00667     {};
00668     
00673     template<class Data, class GatherScatter, bool send>
00674     struct MessageGatherer<Data,GatherScatter,send,SizeOne>
00675     {
00677       typedef typename CommPolicy<Data>::IndexedType Type;
00678       
00683       typedef GatherScatter Gatherer;
00684       
00685       enum{
00691         forward=send
00692       };
00693       
00701       inline void operator()(const InterfaceMap& interface, const Data& data, Type* buffer, size_t bufferSize) const;
00702     };
00703     
00708     template<class Data, class GatherScatter, bool send>
00709     struct MessageGatherer<Data,GatherScatter,send,VariableSize>
00710     {
00712       typedef typename CommPolicy<Data>::IndexedType Type;
00713       
00718       typedef GatherScatter Gatherer;
00719       
00720       enum{
00726         forward=send
00727       };
00728       
00736       inline void operator()(const InterfaceMap& interface, const Data& data, Type* buffer, size_t bufferSize) const;
00737     };
00738 
00742     template<class Data, class GatherScatter, bool send, typename IndexedTypeFlag>
00743     struct MessageScatterer
00744     {};
00745 
00750     template<class Data, class GatherScatter, bool send>
00751     struct MessageScatterer<Data,GatherScatter,send,SizeOne>
00752     {
00754       typedef typename CommPolicy<Data>::IndexedType Type;
00755       
00760       typedef GatherScatter Scatterer;
00761       
00762       enum{
00768         forward=send
00769       };
00770       
00778       inline void operator()(const InterfaceMap& interface, Data& data, Type* buffer, const int& proc) const;
00779     };
00784     template<class Data, class GatherScatter, bool send>
00785     struct MessageScatterer<Data,GatherScatter,send,VariableSize>
00786     {
00788       typedef typename CommPolicy<Data>::IndexedType Type;
00789       
00794       typedef GatherScatter Scatterer;
00795       
00796       enum{
00802         forward=send
00803       };
00804       
00812       inline void operator()(const InterfaceMap& interface, Data& data, Type* buffer, const int& proc) const;
00813     };
00814 
00818     struct MessageInformation
00819     {
00821       MessageInformation()
00822         : start_(0), size_(0)
00823       {}
00824       
00832       MessageInformation(size_t start, size_t size)
00833         :start_(start), size_(size)
00834       {}
00838       size_t start_;
00842       size_t size_;
00843     };
00844 
00851     typedef std::map<int,std::pair<MessageInformation,MessageInformation> >
00852     InformationMap;
00856     InformationMap messageInformation_;
00860     char* buffers_[2];
00864     size_t bufferSize_[2];
00865     
00866     enum{
00870       commTag_
00871     };
00872     
00876     std::map<int,std::pair<InterfaceInformation,InterfaceInformation> > interfaces_;
00877 
00878     MPI_Comm communicator_;
00879 
00883     template<class GatherScatter, bool FORWARD, class Data>
00884     void sendRecv(const Data& source, Data& target);
00885     
00886   };
00887   
00888 #ifndef DOXYGEN
00889 
00890   template<class V>
00891   inline const void* CommPolicy<V>::getAddress(const V& v, int index)
00892   {
00893     return &(v[index]);
00894   }
00895   
00896   template<class V>
00897   inline int CommPolicy<V>::getSize(const V& v, int index)
00898   {
00899     return 1;
00900   }
00901 
00902   template<class K, class A, int n>
00903   inline const void* CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >::getAddress(const Type& v, int index)
00904   {
00905     return &(v[index][0]);
00906   }
00907 
00908   template<class K, class A, int n>
00909   inline int CommPolicy<VariableBlockVector<FieldVector<K, n>, A> >::getSize(const Type& v, int index)
00910   {
00911     return v[index].getsize();
00912   }
00913 
00914   template<class T>
00915   inline const typename CopyGatherScatter<T>::IndexedType& CopyGatherScatter<T>::gather(const T& vec, std::size_t i)
00916   {
00917     return vec[i];
00918   }
00919 
00920   template<class T>
00921   inline void CopyGatherScatter<T>::scatter(T& vec, const IndexedType& v, std::size_t i)
00922   {
00923     vec[i]=v;
00924   }
00925 
00926   template<typename T>
00927   DatatypeCommunicator<T>::DatatypeCommunicator()
00928     : remoteIndices_(0), created_(false)
00929   {
00930     requests_[0]=0;
00931     requests_[1]=0;
00932   }
00933   
00934  
00935 
00936   template<typename T>
00937   DatatypeCommunicator<T>::~DatatypeCommunicator()
00938   {
00939     free();
00940   }
00941   
00942   template<typename T>
00943   template<class T1, class T2, class V>
00944   inline void DatatypeCommunicator<T>::build(const RemoteIndices& remoteIndices,
00945                                                  const T1& source, V& sendData,
00946                                                  const T2& destination, V& receiveData)
00947   {
00948     remoteIndices_ = &remoteIndices;
00949     free();
00950     createDataTypes<T1,T2,V,false>(source,destination, receiveData);
00951     createDataTypes<T1,T2,V,true>(source,destination, sendData);
00952     createRequests<V,true>(sendData, receiveData);
00953     createRequests<V,false>(receiveData, sendData);
00954     created_=true;
00955   }
00956 
00957   template<typename T>
00958   void DatatypeCommunicator<T>::free()
00959   {
00960     if(created_){
00961       delete[] requests_[0];
00962       delete[] requests_[1];
00963       typedef MessageTypeMap::iterator iterator;
00964       typedef MessageTypeMap::const_iterator const_iterator;
00965       
00966       const const_iterator end=messageTypes.end();
00967       
00968       for(iterator process = messageTypes.begin(); process != end; ++process){
00969         MPI_Datatype *type = &(process->second.first);
00970         int finalized=0;
00971 #if MPI_2
00972         MPI_Finalized(&finalized);
00973 #endif
00974         if(*type!=MPI_DATATYPE_NULL && !finalized)
00975           MPI_Type_free(type);
00976         type = &(process->second.second);
00977         if(*type!=MPI_DATATYPE_NULL && !finalized)
00978           MPI_Type_free(type);
00979       }
00980       messageTypes.clear();
00981       created_=false;
00982     }
00983     
00984   }
00985   
00986   template<typename T>
00987   template<class T1, class T2, class V, bool send>
00988   void DatatypeCommunicator<T>::createDataTypes(const T1& sourceFlags, const T2& destFlags, V& data)
00989   {
00990 
00991     MPIDatatypeInformation<V>  dataInfo(data);
00992     this->template buildInterface<RemoteIndices,T1,T2,MPIDatatypeInformation<V>,send>(*remoteIndices_,sourceFlags, destFlags, dataInfo);
00993 
00994     typedef typename RemoteIndices::RemoteIndexMap::const_iterator const_iterator;
00995     const const_iterator end=this->remoteIndices_->end();
00996     
00997     // Allocate MPI_Datatypes and deallocate memory for the type construction.
00998     for(const_iterator process=this->remoteIndices_->begin(); process != end; ++process){
00999       IndexedTypeInformation& info=dataInfo.information_[process->first];
01000       // Shift the displacement
01001       MPI_Aint base;
01002       MPI_Address(const_cast<void *>(CommPolicy<V>::getAddress(data, 0)), &base);
01003       
01004       for(int i=0; i< info.elements; i++){
01005         info.displ[i]-=base;
01006       }
01007       
01008       // Create data type
01009       MPI_Datatype* type = &( send ? messageTypes[process->first].first : messageTypes[process->first].second);
01010       MPI_Type_hindexed(info.elements, info.length, info.displ, 
01011                        MPITraits<typename CommPolicy<V>::IndexedType>::getType(),
01012                        type);
01013       MPI_Type_commit(type);
01014       // Deallocate memory
01015       info.free();
01016     }
01017   }
01018   
01019   template<typename T>
01020   template<class V, bool createForward>
01021   void DatatypeCommunicator<T>::createRequests(V& sendData, V& receiveData)
01022   {
01023     typedef std::map<int,std::pair<MPI_Datatype,MPI_Datatype> >::const_iterator MapIterator;
01024     int rank;
01025     static int index = createForward?1:0;
01026     int noMessages = messageTypes.size();
01027     // allocate request handles
01028     requests_[index] = new MPI_Request[2*noMessages];
01029     const MapIterator end = messageTypes.end();
01030     int request=0;
01031     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01032     
01033     // Set up the requests for receiving first
01034     for(MapIterator process = messageTypes.begin(); process != end; 
01035         ++process, ++request){
01036       MPI_Datatype type = createForward ? process->second.second : process->second.first;
01037       void* address = const_cast<void*>(CommPolicy<V>::getAddress(receiveData,0));
01038       MPI_Recv_init(address, 1, type, process->first, commTag_, this->remoteIndices_->communicator(), requests_[index]+request);
01039     }
01040     
01041     // And now the send requests
01042 
01043     for(MapIterator process = messageTypes.begin(); process != end; 
01044         ++process, ++request){
01045       MPI_Datatype type = createForward ? process->second.first : process->second.second;
01046       void* address =  const_cast<void*>(CommPolicy<V>::getAddress(sendData, 0));
01047       MPI_Ssend_init(address, 1, type, process->first, commTag_, this->remoteIndices_->communicator(), requests_[index]+request);
01048     }
01049   }   
01050         
01051   template<typename T>
01052   void DatatypeCommunicator<T>::forward()
01053   {
01054     sendRecv(requests_[1]);
01055   }
01056   
01057   template<typename T>
01058   void DatatypeCommunicator<T>::backward()
01059   {
01060     sendRecv(requests_[0]);
01061   }
01062 
01063   template<typename T>
01064   void DatatypeCommunicator<T>::sendRecv(MPI_Request* requests)
01065   {
01066     int noMessages = messageTypes.size();
01067     // Start the receive calls first
01068     MPI_Startall(noMessages, requests);
01069     // Now the send calls
01070     MPI_Startall(noMessages, requests+noMessages);
01071     
01072     // Wait for completion of the communication send first then receive
01073     MPI_Status* status=new MPI_Status[2*noMessages];
01074     for(int i=0; i<2*noMessages; i++)
01075       status[i].MPI_ERROR=MPI_SUCCESS;
01076     
01077     int send = MPI_Waitall(noMessages, requests+noMessages, status+noMessages);
01078     int receive = MPI_Waitall(noMessages, requests, status);
01079     
01080     // Error checks
01081     int success=1, globalSuccess=0;
01082     if(send==MPI_ERR_IN_STATUS){
01083       int rank;
01084       MPI_Comm_rank(this->remoteIndices_->communicator(), &rank);
01085       std::cerr<<rank<<": Error in sending :"<<std::endl;
01086       // Search for the error
01087       for(int i=noMessages; i< 2*noMessages; i++)
01088         if(status[i].MPI_ERROR!=MPI_SUCCESS){
01089           char message[300];
01090           int messageLength;
01091           MPI_Error_string(status[i].MPI_ERROR, message, &messageLength);
01092           std::cerr<<" source="<<status[i].MPI_SOURCE<<" message: ";
01093           for(int i=0; i< messageLength; i++)
01094             std::cout<<message[i];
01095         }
01096       std::cerr<<std::endl;
01097       success=0;
01098     }
01099     
01100     if(receive==MPI_ERR_IN_STATUS){
01101       int rank;
01102       MPI_Comm_rank(this->remoteIndices_->communicator(), &rank);
01103       std::cerr<<rank<<": Error in receiving!"<<std::endl;
01104       // Search for the error
01105       for(int i=0; i< noMessages; i++)
01106         if(status[i].MPI_ERROR!=MPI_SUCCESS){
01107           char message[300];
01108           int messageLength;
01109           MPI_Error_string(status[i].MPI_ERROR, message, &messageLength);
01110           std::cerr<<" source="<<status[i].MPI_SOURCE<<" message: ";
01111           for(int i=0; i< messageLength; i++)
01112             std::cerr<<message[i];
01113         }
01114       std::cerr<<std::endl;
01115       success=0;
01116     }
01117 
01118     MPI_Allreduce(&success, &globalSuccess, 1, MPI_INT, MPI_MIN, this->remoteIndices_->communicator());
01119     
01120     delete[] status;
01121 
01122     if(!globalSuccess)
01123       DUNE_THROW(CommunicationError, "A communication error occurred!");
01124       
01125   }
01126   
01127   inline BufferedCommunicator::BufferedCommunicator()
01128   {
01129     buffers_[0]=0;
01130     buffers_[1]=0;
01131     bufferSize_[0]=0;
01132     bufferSize_[1]=0;
01133   }
01134 
01135   template<class Data, class Interface>
01136   typename enable_if<is_same<SizeOne, typename CommPolicy<Data>::IndexedTypeFlag>::value, void>::type
01137   BufferedCommunicator::build(const Interface& interface)
01138   {
01139     interfaces_=interface.interfaces();
01140     communicator_=interface.communicator();
01141     typedef typename std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >
01142       ::const_iterator const_iterator;
01143     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01144     const const_iterator end = interfaces_.end();
01145     int lrank;
01146     MPI_Comm_rank(communicator_, &lrank);
01147     
01148     bufferSize_[0]=0;
01149     bufferSize_[1]=0;
01150 
01151     for(const_iterator interfacePair = interfaces_.begin();
01152         interfacePair != end; ++interfacePair){  
01153       int noSend = MessageSizeCalculator<Data,Flag>()(interfacePair->second.first);
01154       int noRecv = MessageSizeCalculator<Data,Flag>()(interfacePair->second.second);
01155       messageInformation_.insert(std::make_pair(interfacePair->first,
01156                                  std::make_pair(MessageInformation(bufferSize_[0],
01157                                                                    noSend*sizeof(typename CommPolicy<Data>::IndexedType)),
01158                                                 MessageInformation(bufferSize_[1],
01159                                                                    noRecv*sizeof(typename CommPolicy<Data>::IndexedType)))));
01160       bufferSize_[0] += noSend;
01161       bufferSize_[1] += noRecv;
01162     }
01163 
01164     // allocate the buffers
01165     bufferSize_[0] *= sizeof(typename CommPolicy<Data>::IndexedType);
01166     bufferSize_[1] *= sizeof(typename CommPolicy<Data>::IndexedType);
01167     
01168     buffers_[0] = new char[bufferSize_[0]];
01169     buffers_[1] = new char[bufferSize_[1]];    
01170   }  
01171   
01172   template<class Data, class Interface>
01173   void BufferedCommunicator::build(const Data& source, const Data& dest, const Interface& interface)
01174   {
01175     
01176     interfaces_=interface.interfaces();
01177     communicator_=interface.communicator();
01178     typedef typename std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >
01179       ::const_iterator const_iterator;
01180     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01181     const const_iterator end = interfaces_.end();
01182     
01183     bufferSize_[0]=0;
01184     bufferSize_[1]=0;
01185 
01186     for(const_iterator interfacePair = interfaces_.begin();
01187         interfacePair != end; ++interfacePair){      
01188       int noSend = MessageSizeCalculator<Data,Flag>()(source, interfacePair->second.first);
01189       int noRecv = MessageSizeCalculator<Data,Flag>()(dest, interfacePair->second.second);
01190       
01191       messageInformation_.insert(std::make_pair(interfacePair->first,
01192                                                 std::make_pair(MessageInformation(bufferSize_[0],
01193                                                                                   noSend*sizeof(typename CommPolicy<Data>::IndexedType)),
01194                                                                MessageInformation(bufferSize_[1],
01195                                                                                   noRecv*sizeof(typename CommPolicy<Data>::IndexedType)))));
01196       bufferSize_[0] += noSend;
01197       bufferSize_[1] += noRecv;
01198     }
01199 
01200     bufferSize_[0] *= sizeof(typename CommPolicy<Data>::IndexedType);
01201     bufferSize_[1] *= sizeof(typename CommPolicy<Data>::IndexedType);
01202     // allocate the buffers
01203     buffers_[0] = new char[bufferSize_[0]];
01204     buffers_[1] = new char[bufferSize_[1]];
01205   }
01206   
01207   inline void BufferedCommunicator::free()
01208   {
01209       messageInformation_.clear();
01210       if(buffers_[0])
01211         delete[] buffers_[0];
01212       
01213       if(buffers_[1])
01214         delete[] buffers_[1];
01215       buffers_[0]=buffers_[1]=0;
01216   }
01217 
01218   inline BufferedCommunicator::~BufferedCommunicator()
01219   {
01220     free();
01221   }
01222   
01223   template<class Data>
01224   inline int BufferedCommunicator::MessageSizeCalculator<Data,SizeOne>::operator()
01225     (const InterfaceInformation& info) const
01226   {
01227     return info.size();
01228   }
01229 
01230   
01231   template<class Data>
01232   inline int BufferedCommunicator::MessageSizeCalculator<Data,SizeOne>::operator()
01233     (const Data& data, const InterfaceInformation& info) const
01234   {
01235     return operator()(info);
01236   }
01237 
01238   
01239   template<class Data>
01240   inline int BufferedCommunicator::MessageSizeCalculator<Data, VariableSize>::operator()
01241     (const Data& data, const InterfaceInformation& info) const
01242   {
01243     int entries=0;
01244 
01245     for(size_t i=0; i < info.size(); i++)
01246       entries += CommPolicy<Data>::getSize(data,info[i]);
01247 
01248     return entries;
01249   }
01250   
01251   
01252   template<class Data, class GatherScatter, bool FORWARD>
01253   inline void BufferedCommunicator::MessageGatherer<Data,GatherScatter,FORWARD,VariableSize>::operator()(const InterfaceMap& interfaces,const Data& data, Type* buffer, size_t bufferSize) const
01254   {
01255     typedef typename InterfaceMap::const_iterator
01256       const_iterator;
01257 
01258     int rank;
01259     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01260     const const_iterator end = interfaces.end();
01261     size_t index=0;
01262     
01263     for(const_iterator interfacePair = interfaces.begin();
01264         interfacePair != end; ++interfacePair){
01265       int size = forward ? interfacePair->second.first.size() :
01266         interfacePair->second.second.size();
01267       
01268       for(int i=0; i < size; i++){  
01269         int local = forward ? interfacePair->second.first[i] :
01270           interfacePair->second.second[i];
01271         for(std::size_t j=0; j < CommPolicy<Data>::getSize(data, local);j++, index++){
01272 
01273 #ifdef DUNE_ISTL_WITH_CHECKING
01274           assert(bufferSize>=(index+1)*sizeof(typename CommPolicy<Data>::IndexedType));
01275 #endif
01276           buffer[index]=GatherScatter::gather(data, local, j);
01277         }
01278         
01279       }
01280     }
01281         
01282   }
01283 
01284   
01285   template<class Data, class GatherScatter, bool FORWARD>
01286   inline void BufferedCommunicator::MessageGatherer<Data,GatherScatter,FORWARD,SizeOne>::operator()(const InterfaceMap& interfaces, const Data& data, Type* buffer, size_t bufferSize)const
01287   {
01288     typedef typename InterfaceMap::const_iterator
01289       const_iterator;
01290     const const_iterator end = interfaces.end();
01291     size_t index = 0;
01292     
01293     int rank;
01294     MPI_Comm_rank(MPI_COMM_WORLD, &rank);
01295 
01296     for(const_iterator interfacePair = interfaces.begin();
01297         interfacePair != end; ++interfacePair){
01298       size_t size = FORWARD ? interfacePair->second.first.size() :
01299         interfacePair->second.second.size();
01300       
01301       for(size_t i=0; i < size; i++){
01302 
01303 #ifdef DUNE_ISTL_WITH_CHECKING
01304           assert(bufferSize>=(index+1)*sizeof(typename CommPolicy<Data>::IndexedType));
01305 #endif
01306 
01307         buffer[index++] = GatherScatter::gather(data, FORWARD ? interfacePair->second.first[i] :
01308                                                 interfacePair->second.second[i]);
01309       }
01310     }
01311         
01312   }
01313   
01314   
01315   template<class Data, class GatherScatter, bool FORWARD>
01316   inline void BufferedCommunicator::MessageScatterer<Data,GatherScatter,FORWARD,VariableSize>::operator()(const InterfaceMap& interfaces, Data& data, Type* buffer, const int& proc)const
01317   {
01318     typedef typename InterfaceMap::value_type::second_type::first_type Information;
01319     const typename InterfaceMap::const_iterator infoPair = interfaces.find(proc);
01320 
01321     assert(infoPair!=interfaces.end());
01322     
01323     const Information& info = FORWARD ? infoPair->second.second :
01324        infoPair->second.first;
01325 
01326     for(size_t i=0, index=0; i < info.size(); i++){
01327         for(size_t j=0; j < CommPolicy<Data>::getSize(data, info[i]); j++)
01328           GatherScatter::scatter(data, buffer[index++], info[i], j);
01329     }
01330   }
01331 
01332   
01333   template<class Data, class GatherScatter, bool FORWARD>
01334   inline void BufferedCommunicator::MessageScatterer<Data,GatherScatter,FORWARD,SizeOne>::operator()(const InterfaceMap& interfaces, Data& data, Type* buffer, const int& proc)const
01335   {
01336     typedef typename InterfaceMap::value_type::second_type::first_type Information;
01337     const typename InterfaceMap::const_iterator infoPair = interfaces.find(proc);
01338 
01339     assert(infoPair!=interfaces.end());
01340     
01341     const Information& info = FORWARD ? infoPair->second.second :
01342       infoPair->second.first;
01343     
01344     for(size_t i=0; i < info.size(); i++){
01345       GatherScatter::scatter(data, buffer[i], info[i]);
01346     }
01347   }
01348   
01349   
01350   template<class GatherScatter,class Data>
01351   void BufferedCommunicator::forward(Data& data)
01352   {
01353     this->template sendRecv<GatherScatter,true>(data, data);
01354   }
01355   
01356   
01357   template<class GatherScatter, class Data>
01358   void BufferedCommunicator::backward(Data& data)
01359   {
01360     this->template sendRecv<GatherScatter,false>(data, data);
01361   }
01362 
01363   
01364   template<class GatherScatter, class Data>
01365   void BufferedCommunicator::forward(const Data& source, Data& dest)
01366   {
01367     this->template sendRecv<GatherScatter,true>(source, dest);
01368   }
01369   
01370   
01371   template<class GatherScatter, class Data>
01372   void BufferedCommunicator::backward(Data& source, const Data& dest)
01373   {
01374     this->template sendRecv<GatherScatter,false>(dest, source);
01375   }
01376 
01377   
01378   template<class GatherScatter, bool FORWARD, class Data>
01379   void BufferedCommunicator::sendRecv(const Data& source, Data& dest) 
01380   {
01381     int rank, lrank;
01382     
01383     MPI_Comm_rank(MPI_COMM_WORLD,&rank);
01384     MPI_Comm_rank(MPI_COMM_WORLD,&lrank);
01385     
01386     typedef typename CommPolicy<Data>::IndexedType Type;
01387     Type *sendBuffer, *recvBuffer;
01388     size_t sendBufferSize, recvBufferSize;
01389     
01390     if(FORWARD){
01391       sendBuffer = reinterpret_cast<Type*>(buffers_[0]);
01392       sendBufferSize = bufferSize_[0];
01393       recvBuffer = reinterpret_cast<Type*>(buffers_[1]);
01394       recvBufferSize = bufferSize_[1];
01395     }else{
01396       sendBuffer = reinterpret_cast<Type*>(buffers_[1]);
01397       sendBufferSize = bufferSize_[1];
01398       recvBuffer = reinterpret_cast<Type*>(buffers_[0]);
01399       recvBufferSize = bufferSize_[0];
01400     }
01401     typedef typename CommPolicy<Data>::IndexedTypeFlag Flag;
01402 
01403     MessageGatherer<Data,GatherScatter,FORWARD,Flag>()(interfaces_, source, sendBuffer, sendBufferSize);
01404     
01405     MPI_Request* sendRequests = new MPI_Request[messageInformation_.size()];
01406     MPI_Request* recvRequests = new MPI_Request[messageInformation_.size()];
01407     
01408     // Setup receive first
01409     typedef typename InformationMap::const_iterator const_iterator;
01410 
01411     const const_iterator end = messageInformation_.end();
01412     size_t i=0;
01413     int* processMap = new int[messageInformation_.size()];
01414     
01415     for(const_iterator info = messageInformation_.begin(); info != end; ++info, ++i){
01416           processMap[i]=info->first;
01417           if(FORWARD){
01418             assert(info->second.second.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.second.size_ <= recvBufferSize );
01419             Dune::dvverb<<rank<<": receiving "<<info->second.second.size_<<" from "<<info->first<<std::endl;
01420             MPI_Irecv(recvBuffer+info->second.second.start_, info->second.second.size_,
01421                       MPI_BYTE, info->first, commTag_, communicator_,
01422                       recvRequests+i);
01423           }else{
01424             assert(info->second.first.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.first.size_ <= recvBufferSize );
01425             Dune::dvverb<<rank<<": receiving "<<info->second.first.size_<<" to "<<info->first<<std::endl;
01426             MPI_Irecv(recvBuffer+info->second.first.start_, info->second.first.size_,
01427                       MPI_BYTE, info->first, commTag_, communicator_,
01428                       recvRequests+i);
01429           }
01430     }
01431     
01432     // now the send requests
01433     i=0;
01434     for(const_iterator info = messageInformation_.begin(); info != end; ++info, ++i)
01435       if(FORWARD){
01436         assert(info->second.second.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.second.size_ <= recvBufferSize );
01437             Dune::dvverb<<rank<<": sending "<<info->second.first.size_<<" to "<<info->first<<std::endl;
01438         assert(info->second.first.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.first.size_ <= sendBufferSize );
01439         MPI_Issend(sendBuffer+info->second.first.start_, info->second.first.size_,
01440                    MPI_BYTE, info->first, commTag_, communicator_,
01441                    sendRequests+i);
01442       }else{
01443         assert(info->second.second.start_*sizeof(typename CommPolicy<Data>::IndexedType)+info->second.second.size_ <= sendBufferSize );
01444         Dune::dvverb<<rank<<": sending "<<info->second.second.size_<<" to "<<info->first<<std::endl;
01445         MPI_Issend(sendBuffer+info->second.second.start_, info->second.second.size_,
01446                    MPI_BYTE, info->first, commTag_, communicator_,
01447                    sendRequests+i);
01448       }
01449     
01450     // Wait for completion of receive and immediately start scatter
01451     i=0;
01452     int success = 1;
01453     int finished = MPI_UNDEFINED;
01454     MPI_Status status;//[messageInformation_.size()];
01455     //MPI_Waitall(messageInformation_.size(), recvRequests, status);
01456     
01457     for(i=0;i< messageInformation_.size();i++){
01458       status.MPI_ERROR=MPI_SUCCESS;
01459       MPI_Waitany(messageInformation_.size(), recvRequests, &finished, &status);
01460       assert(finished != MPI_UNDEFINED);
01461       
01462       if(status.MPI_ERROR==MPI_SUCCESS){
01463         int& proc = processMap[finished];
01464         typename InformationMap::const_iterator infoIter = messageInformation_.find(proc);
01465         assert(infoIter != messageInformation_.end());
01466 
01467         MessageInformation info = (FORWARD)? infoIter->second.second : infoIter->second.first;
01468         assert(info.start_+info.size_ <= recvBufferSize);
01469 
01470         MessageScatterer<Data,GatherScatter,FORWARD,Flag>()(interfaces_, dest, recvBuffer+info.start_, proc);
01471       }else{
01472         std::cerr<<rank<<": MPI_Error occurred while receiving message from "<<processMap[finished]<<std::endl;
01473         success=0;
01474       }
01475     }
01476 
01477     MPI_Status recvStatus;
01478     
01479     // Wait for completion of sends
01480     for(i=0;i< messageInformation_.size();i++)
01481       if(MPI_SUCCESS!=MPI_Wait(sendRequests+i, &recvStatus)){
01482         std::cerr<<rank<<": MPI_Error occurred while sending message to "<<processMap[finished]<<std::endl;
01483         success=0;
01484       }
01485     /*
01486     int globalSuccess;
01487     MPI_Allreduce(&success, &globalSuccess, 1, MPI_INT, MPI_MIN, interface_->communicator());
01488 
01489     if(!globalSuccess)
01490       DUNE_THROW(CommunicationError, "A communication error occurred!");
01491     */
01492     delete[] processMap;
01493     delete[] sendRequests;
01494     delete[] recvRequests;
01495 
01496   }
01497 
01498 #endif  // DOXYGEN
01499   
01501 }
01502 
01503 #endif
01504 
01505 #endif