mpitraits.hh

Go to the documentation of this file.
00001 // $Id: mpitraits.hh 1108 2009-10-28 10:15:36Z mnolte $
00002 #ifndef DUNE_MPITRAITS_HH
00003 #define DUNE_MPITRAITS_HH
00004 
00005 #if HAVE_MPI
00006 #include"mpi.h"
00007 #endif
00008 
00009 namespace Dune
00010 {
00029   template<typename T>
00030   class MPITraits
00031   {};
00032 
00033 #if HAVE_MPI
00034 
00035   // A Macro for defining traits for the primitive data types
00036 #define ComposeMPITraits(p,m) \
00037   template<> \
00038   struct MPITraits<p>{ \
00039     static inline MPI_Datatype getType(){ \
00040       return m; \
00041     } \
00042   }
00043   
00044   ComposeMPITraits(char, MPI_CHAR);
00045   ComposeMPITraits(unsigned char,MPI_UNSIGNED_CHAR);
00046   ComposeMPITraits(short,MPI_SHORT);
00047   ComposeMPITraits(unsigned short,MPI_UNSIGNED_SHORT);
00048   ComposeMPITraits(int,MPI_INT);
00049   ComposeMPITraits(unsigned int,MPI_UNSIGNED);
00050   ComposeMPITraits(long,MPI_LONG);
00051   ComposeMPITraits(unsigned long,MPI_UNSIGNED_LONG);
00052   ComposeMPITraits(float,MPI_FLOAT);
00053   ComposeMPITraits(double,MPI_DOUBLE);
00054   ComposeMPITraits(long double,MPI_LONG_DOUBLE);
00055 
00056 #undef ComposeMPITraits
00057 
00058   template<class K, int n> class FieldVector;
00059   
00060   template<class K, int n>
00061   struct MPITraits<FieldVector<K,n> >
00062   {
00063     static MPI_Datatype datatype;
00064     static MPI_Datatype vectortype;
00065     
00066     static inline MPI_Datatype getType()
00067     {
00068       if(datatype==MPI_DATATYPE_NULL){
00069         MPI_Type_contiguous(n, MPITraits<K>::getType(), &vectortype);
00070         MPI_Type_commit(&vectortype);
00071         FieldVector<K,n> fvector;
00072         MPI_Aint base;
00073         MPI_Aint displ;
00074         MPI_Address(&fvector, &base);
00075         MPI_Address(&(fvector[0]), &displ);
00076         displ -= base;
00077         int length[1]={1};
00078             
00079         MPI_Type_struct(1, length, &displ, &vectortype, &datatype);
00080         MPI_Type_commit(&datatype);
00081       }
00082       return datatype;
00083     }
00084         
00085   };
00086 
00087   template<class K, int n>
00088   MPI_Datatype MPITraits<FieldVector<K,n> >::datatype = MPI_DATATYPE_NULL;
00089   template<class K, int n>
00090   MPI_Datatype MPITraits<FieldVector<K,n> >::vectortype = {MPI_DATATYPE_NULL};
00091 
00092 
00093   template<int k>
00094   class bigunsignedint;
00095   
00096   template<int k>
00097   struct MPITraits<bigunsignedint<k> >
00098   {
00099     static MPI_Datatype datatype;
00100     static MPI_Datatype vectortype;
00101     
00102     static inline MPI_Datatype getType()
00103     {
00104       if(datatype==MPI_DATATYPE_NULL){
00105         MPI_Type_contiguous(bigunsignedint<k>::n, MPITraits<unsigned short>::getType(),
00106                             &vectortype);
00107         //MPI_Type_commit(&vectortype);
00108         bigunsignedint<k> data;
00109         MPI_Aint base;
00110         MPI_Aint displ;
00111         MPI_Address(&data, &base);
00112         MPI_Address(&(data.digit), &displ);
00113         displ -= base;
00114         int length[1]={1};      
00115         MPI_Type_struct(1, length, &displ, &vectortype, &datatype);
00116         MPI_Type_commit(&datatype);
00117       }
00118       return datatype;
00119     }
00120   };
00121 
00122   
00123   template<int k>
00124   MPI_Datatype MPITraits<bigunsignedint<k> >::datatype = MPI_DATATYPE_NULL;
00125   template<int k>
00126   MPI_Datatype MPITraits<bigunsignedint<k> >::vectortype = MPI_DATATYPE_NULL;
00127 
00128 #endif
00129 
00131 }
00132 
00133 #endif

Generated on Sun Nov 15 22:29:36 2009 for dune-istl by  doxygen 1.5.6