DUNE-FEM (unstable)

mpitraits.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3// SPDX-FileCopyrightInfo: Copyright © DUNE Project contributors, see file LICENSE.md in module root
4// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
5#ifndef DUNE_COMMON_PARALLEL_MPITRAITS_HH
6#define DUNE_COMMON_PARALLEL_MPITRAITS_HH
7
18#if HAVE_MPI
19
20#include <complex>
21#include <cstddef>
22#include <cstdint>
23#include <type_traits>
24#include <utility>
25
26#include <mpi.h>
27
28namespace Dune
29{
39 template<typename T>
40 struct MPITraits
41 {
42 private:
43 MPITraits(){}
44 MPITraits(const MPITraits&){}
45 static MPI_Datatype datatype;
46 static MPI_Datatype vectortype;
47 public:
48 static inline MPI_Datatype getType()
49 {
50 if(datatype==MPI_DATATYPE_NULL) {
51 MPI_Type_contiguous(sizeof(T),MPI_BYTE,&datatype);
52 MPI_Type_commit(&datatype);
53 }
54 return datatype;
55 }
56 static constexpr bool is_intrinsic = false;
57 };
58 template<class T>
59 MPI_Datatype MPITraits<T>::datatype = MPI_DATATYPE_NULL;
60
61#ifndef DOXYGEN
62
63 // A Macro for defining traits for the primitive data types
64#define ComposeMPITraits(p,m) \
65 template<> \
66 struct MPITraits<p>{ \
67 static inline MPI_Datatype getType(){ \
68 return m; \
69 } \
70 static constexpr bool is_intrinsic = true; \
71 }
72
73 ComposeMPITraits(char, MPI_CHAR);
74 ComposeMPITraits(unsigned char,MPI_UNSIGNED_CHAR);
75 ComposeMPITraits(short,MPI_SHORT);
76 ComposeMPITraits(unsigned short,MPI_UNSIGNED_SHORT);
77 ComposeMPITraits(int,MPI_INT);
78 ComposeMPITraits(unsigned int,MPI_UNSIGNED);
79 ComposeMPITraits(long,MPI_LONG);
80 ComposeMPITraits(unsigned long,MPI_UNSIGNED_LONG);
81 ComposeMPITraits(float,MPI_FLOAT);
82 ComposeMPITraits(double,MPI_DOUBLE);
83 ComposeMPITraits(long double,MPI_LONG_DOUBLE);
84 ComposeMPITraits(std::complex<double>, MPI_CXX_DOUBLE_COMPLEX);
85 ComposeMPITraits(std::complex<long double>, MPI_CXX_LONG_DOUBLE_COMPLEX);
86 ComposeMPITraits(std::complex<float>, MPI_CXX_FLOAT_COMPLEX);
87
88
89#undef ComposeMPITraits
90
91 template<class K, int n> class FieldVector;
92
93 template<class K, int n>
94 struct MPITraits<FieldVector<K,n> >
95 {
96 static MPI_Datatype datatype;
97 static MPI_Datatype vectortype;
98
99 static inline MPI_Datatype getType()
100 {
101 if(datatype==MPI_DATATYPE_NULL) {
102 MPI_Type_contiguous(n, MPITraits<K>::getType(), &vectortype);
103 MPI_Type_commit(&vectortype);
104 FieldVector<K,n> fvector;
105 MPI_Aint base;
106 MPI_Aint displ;
107 MPI_Get_address(&fvector, &base);
108 MPI_Get_address(&(fvector[0]), &displ);
109 displ -= base;
110 int length[1]={1};
111
112 MPI_Type_create_struct(1, length, &displ, &vectortype, &datatype);
113 MPI_Type_commit(&datatype);
114 }
115 return datatype;
116 }
117
118 };
119
120 template<class K, int n>
121 MPI_Datatype MPITraits<FieldVector<K,n> >::datatype = MPI_DATATYPE_NULL;
122 template<class K, int n>
123 MPI_Datatype MPITraits<FieldVector<K,n> >::vectortype = {MPI_DATATYPE_NULL};
124
125
126 template<int k>
127 class bigunsignedint;
128
129 template<int k>
130 struct MPITraits<bigunsignedint<k> >
131 {
132 static MPI_Datatype datatype;
133 static MPI_Datatype vectortype;
134
135 static inline MPI_Datatype getType()
136 {
137 if(datatype==MPI_DATATYPE_NULL) {
138 MPI_Type_contiguous(bigunsignedint<k>::n, MPITraits<std::uint16_t>::getType(),
139 &vectortype);
140 //MPI_Type_commit(&vectortype);
141 bigunsignedint<k> data;
142 MPI_Aint base;
143 MPI_Aint displ;
144 MPI_Get_address(&data, &base);
145 MPI_Get_address(&(data.digit), &displ);
146 displ -= base;
147 int length[1]={1};
148 MPI_Type_create_struct(1, length, &displ, &vectortype, &datatype);
149 MPI_Type_commit(&datatype);
150 }
151 return datatype;
152 }
153 };
154}
155
156namespace Dune
157{
158 template<int k>
159 MPI_Datatype MPITraits<bigunsignedint<k> >::datatype = MPI_DATATYPE_NULL;
160 template<int k>
161 MPI_Datatype MPITraits<bigunsignedint<k> >::vectortype = MPI_DATATYPE_NULL;
162
163 template<typename T1, typename T2>
164 struct MPITraits<std::pair<T1,T2 > >
165 {
166 public:
167 inline static MPI_Datatype getType();
168 private:
169 static MPI_Datatype type;
170 };
171 template<typename T1, typename T2>
172 MPI_Datatype MPITraits<std::pair<T1,T2> >::getType()
173 {
174 if(type==MPI_DATATYPE_NULL) {
175 int length[2] = {1, 1};
176 MPI_Aint disp[2];
177 MPI_Datatype types[2] = {MPITraits<T1>::getType(),
178 MPITraits<T2>::getType()};
179
180 using Pair = std::pair<T1, T2>;
181 static_assert(std::is_standard_layout<Pair>::value, "offsetof() is only defined for standard layout types");
182 disp[0] = offsetof(Pair, first);
183 disp[1] = offsetof(Pair, second);
184
185 MPI_Datatype tmp;
186 MPI_Type_create_struct(2, length, disp, types, &tmp);
187
188 MPI_Type_create_resized(tmp, 0, sizeof(Pair), &type);
189 MPI_Type_commit(&type);
190
191 MPI_Type_free(&tmp);
192 }
193 return type;
194 }
195
196 template<typename T1, typename T2>
197 MPI_Datatype MPITraits<std::pair<T1,T2> >::type=MPI_DATATYPE_NULL;
198
199#endif // !DOXYGEN
200
201} // namespace Dune
202
203#endif // HAVE_MPI
204
207#endif // DUNE_COMMON_PARALLEL_MPITRAITS_HH
vector space out of a tensor product of fields.
Definition: fvector.hh:91
Dune namespace.
Definition: alignedallocator.hh:13
STL namespace.
A traits class describing the mapping of types onto MPI_Datatypes.
Definition: mpitraits.hh:41
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)