Communicator.inc.hpp 4.03 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
#pragma once

#if HAVE_MPI

namespace AMDiS { namespace Mpi {

// send mpi datatype
template <class Data, REQUIRES_(not Concepts::RawPointer<Data>)>
void Communicator::send(Data const& data, int to, Tag tag) const
{
  MPI_Send(to_void_ptr(&data), 1, type_to_mpi<Data>(), to, tag.value, comm_);
}


// send array of mpi datatypes
template <class Data>
void Communicator::send(Data const* data, std::size_t size, int to, Tag tag) const
{
  MPI_Send(to_void_ptr(data), int(size), type_to_mpi<Data>(), to, tag.value, comm_);
}


template <class T>
void Communicator::send(std::vector<T> const& vec, int to, Tag tag) const
{
  MPI_Send(to_void_ptr(vec.data()), int(vec.size()), type_to_mpi<T>(), to, tag.value, comm_);
}


// -------------------------------------------------------------------------------------


// send mpi datatype (non-blocking)
template <class Data, REQUIRES_(not Concepts::RawPointer<Data>)>
Request Communicator::isend(Data const& data, int to, Tag tag) const
{
  MPI_Request request;
  MPI_Isend(to_void_ptr(&data), 1, type_to_mpi<Data>(), to, tag.value, comm_, &request);
  return {request};
}


// send array of mpi datatypes (non-blocking)
template <class Data>
Request Communicator::isend(Data const* data, std::size_t size, int to, Tag tag) const
{
  MPI_Request request;
  MPI_Isend(to_void_ptr(data), size, type_to_mpi<Data>(), to, tag.value, comm_, &request);
  return {request};
}


template <class T>
Request Communicator::isend(std::vector<T> const& vec, int to, Tag tag) const
{
  MPI_Request request;
  MPI_Isend(to_void_ptr(vec.data()), int(vec.size()), type_to_mpi<T>(), to, tag.value, comm_, &request);
  return {request};
}

// -------------------------------------------------------------------------------------

// receive mpi datatype
template <class Data, REQUIRES_(not Concepts::RawPointer<Data>)>
MPI_Status Communicator::recv(Data& data, int from, Tag tag) const
{
  MPI_Status status;
  MPI_Recv(&data, 1, type_to_mpi<Data>(), from, tag.value, comm_, &status);
  return status;
}


// receive array of mpi datatypes
template <class Data>
MPI_Status Communicator::recv(Data* data, std::size_t size, int from, Tag tag) const
{
  MPI_Status status;
  MPI_Recv(data, size, type_to_mpi<Data>(), from, tag.value, comm_, &status);
  return status;
}


// receive array of mpi datatypes
template <class T>
MPI_Status Communicator::recv(std::vector<T>& vec, int from, Tag tag) const
{
  MPI_Status status;
  MPI_Probe(from, tag.value, comm_, &status);

  int size = 0;
  MPI_Get_count(&status, type_to_mpi<T>(), &size);
  int min_size = std::max(size,1);

  vec.resize(min_size);
  MPI_Recv(vec.data(), min_size, type_to_mpi<T>(), from, tag.value, comm_, MPI_STATUS_IGNORE);
  if (size != min_size)
    vec.resize(size);
  return status;
}


// -------------------------------------------------------------------------------------

// receive mpi datatype
template <class Data, REQUIRES_(not Concepts::RawPointer<Data>)>
Request Communicator::irecv(Data& data, int from, Tag tag) const
{
  MPI_Request request;
  MPI_Irecv(&data, 1, type_to_mpi<Data>(), from, tag.value, comm_, &request);
  return {request};
}


// receive array of mpi datatypes
template <class Data>
Request Communicator::irecv(Data* data, std::size_t size, int from, Tag tag) const
{
  MPI_Request request;
  MPI_Irecv(data, size, type_to_mpi<Data>(), from, tag.value, comm_, &request);
  return {request};
}


template <class T>
Request Communicator::irecv(std::vector<T>& vec, int from, Tag tag) const
{
  return Request{ RecvDynamicSize(from,tag.value,comm_,
    [comm=comm_,&vec](MPI_Status status) -> MPI_Request
    {
      int size = 0;
      MPI_Get_count(&status, type_to_mpi<T>(), &size);
      int min_size = std::max(size,1);

      vec.resize(min_size);
      MPI_Request req;
      MPI_Irecv(vec.data(), min_size, type_to_mpi<T>(), status.MPI_SOURCE, status.MPI_TAG, comm, &req);
      return req;
    },
    [&vec](MPI_Status status)
    {
      int size = 0;
      MPI_Get_count(&status, type_to_mpi<T>(), &size);
      vec.resize(size);
    }) };
}

}} // end namespace AMDiS::Mpi

#endif // HAVE_MPI