DOLFIN
DOLFIN C++ interface
Loading...
Searching...
No Matches
MPI.h
1// Copyright (C) 2007-2014 Magnus Vikstrøm and Garth N. Wells
2//
3// This file is part of DOLFIN.
4//
5// DOLFIN is free software: you can redistribute it and/or modify
6// it under the terms of the GNU Lesser General Public License as published by
7// the Free Software Foundation, either version 3 of the License, or
8// (at your option) any later version.
9//
10// DOLFIN is distributed in the hope that it will be useful,
11// but WITHOUT ANY WARRANTY; without even the implied warranty of
12// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13// GNU Lesser General Public License for more details.
14//
15// You should have received a copy of the GNU Lesser General Public License
16// along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
17//
18// Modified by Ola Skavhaug 2008-2009
19// Modified by Anders Logg 2008-2011
20// Modified by Niclas Jansson 2009
21// Modified by Joachim B Haga 2012
22// Modified by Martin Sandve Alnes 2014
23
24#ifndef __MPI_DOLFIN_WRAPPER_H
25#define __MPI_DOLFIN_WRAPPER_H
26
27#include <cstdint>
28#include <iostream>
29#include <numeric>
30#include <type_traits>
31#include <utility>
32#include <vector>
33
34#ifdef HAS_MPI
35#define MPICH_IGNORE_CXX_SEEK 1
36#include <mpi.h>
37#endif
38
39#include <dolfin/log/log.h>
40#include <dolfin/log/Table.h>
41
42#ifndef HAS_MPI
43typedef int MPI_Comm;
44#define MPI_COMM_WORLD 2
45#define MPI_COMM_SELF 1
46#define MPI_COMM_NULL 0
47#endif
48
49namespace dolfin
50{
51
52 #ifdef HAS_MPI
53 class MPIInfo
54 {
55 public:
56 MPIInfo();
57 ~MPIInfo();
58 MPI_Info& operator*();
59 private:
60
61 MPI_Info info;
62
63 };
64 #endif
65
69
70 class MPI
71 {
72 public:
73
74 // Create a duplicate MPI communicator and manage lifetime of the
75 // communicator
76 class Comm
77 {
78 public:
79
81 Comm(MPI_Comm comm);
82
83 // The disabling of the below is turned off because the SWIG
84 // docstring generator fails on it.
85
86 // Disable default constructor
87 //Comm() = default;
88
89 // Disable copy constructor
90 //Comm(const Comm& comm) = delete;
91
92 // Disable assignment operator
93 //void operator=(Comm const &comm) = delete;
94
96 ~Comm();
97
99 void free();
100
103 void reset(MPI_Comm comm);
104
106 unsigned int rank() const;
107
111 unsigned int size() const;
112
114 void barrier() const;
115
117 MPI_Comm comm() const;
118
119 private:
120
121 // MPI communicator
122 MPI_Comm _comm;
123 };
124
126 static unsigned int rank(MPI_Comm comm);
127
130 static unsigned int size(MPI_Comm comm);
131
134 static bool is_broadcaster(MPI_Comm comm);
135
138 static bool is_receiver(MPI_Comm comm);
139
141 static void barrier(MPI_Comm comm);
142
145 template<typename T>
146 static void all_to_all(MPI_Comm comm,
147 std::vector<std::vector<T>>& in_values,
148 std::vector<std::vector<T>>& out_values);
149
150
153 template<typename T>
154 static void all_to_all(MPI_Comm comm,
155 std::vector<std::vector<T>>& in_values,
156 std::vector<T>& out_values);
157
159 template<typename T>
160 static void broadcast(MPI_Comm comm, std::vector<T>& value,
161 unsigned int broadcaster=0);
162
164 template<typename T>
165 static void broadcast(MPI_Comm comm, T& value,
166 unsigned int broadcaster=0);
167
169 template<typename T>
170 static void scatter(MPI_Comm comm,
171 const std::vector<std::vector<T>>& in_values,
172 std::vector<T>& out_value,
173 unsigned int sending_process=0);
174
176 template<typename T>
177 static void scatter(MPI_Comm comm,
178 const std::vector<T>& in_values,
179 T& out_value, unsigned int sending_process=0);
180
182 template<typename T>
183 static void gather(MPI_Comm comm, const std::vector<T>& in_values,
184 std::vector<T>& out_values,
185 unsigned int receiving_process=0);
186
188 static void gather(MPI_Comm comm, const std::string& in_values,
189 std::vector<std::string>& out_values,
190 unsigned int receiving_process=0);
191
194 template<typename T>
195 static void all_gather(MPI_Comm comm,
196 const std::vector<T>& in_values,
197 std::vector<T>& out_values);
198
200 template<typename T>
201 static void all_gather(MPI_Comm comm,
202 const std::vector<T>& in_values,
203 std::vector<std::vector<T>>& out_values);
204
206 template<typename T>
207 static void all_gather(MPI_Comm comm, const T in_value,
208 std::vector<T>& out_values);
209
212 static void all_gather(MPI_Comm comm, const std::string& in_values,
213 std::vector<std::string>& out_values);
214
216 template<typename T> static T max(MPI_Comm comm, const T& value);
217
219 template<typename T> static T min(MPI_Comm comm, const T& value);
220
222 template<typename T> static T sum(MPI_Comm comm, const T& value);
223
225 template<typename T> static T avg(MPI_Comm comm, const T& value);
226
228 template<typename T, typename X> static
229 T all_reduce(MPI_Comm comm, const T& value, X op);
230
233 static std::size_t global_offset(MPI_Comm comm,
234 std::size_t range, bool exclusive);
235
237 template<typename T>
238 static void send_recv(MPI_Comm comm,
239 const std::vector<T>& send_value,
240 unsigned int dest, int send_tag,
241 std::vector<T>& recv_value,
242 unsigned int source, int recv_tag);
243
245 template<typename T>
246 static void send_recv(MPI_Comm comm,
247 const std::vector<T>& send_value, unsigned int dest,
248 std::vector<T>& recv_value, unsigned int source);
249
252 static std::pair<std::int64_t, std::int64_t>
253 local_range(MPI_Comm comm, std::int64_t N);
254
257 static std::pair<std::int64_t, std::int64_t>
258 local_range(MPI_Comm comm, int process, std::int64_t N);
259
262 static std::pair<std::int64_t, std::int64_t>
263 compute_local_range(int process, std::int64_t N, int size);
264
266 static unsigned int index_owner(MPI_Comm comm,
267 std::size_t index, std::size_t N);
268
269 #ifdef HAS_MPI
272 static MPI_Op MPI_AVG();
273 #endif
274
275 private:
276
277 #ifndef HAS_MPI
278 static void error_no_mpi(const char *where)
279 {
280 dolfin_error("MPI.h",
281 where,
282 "DOLFIN has been configured without MPI support");
283 }
284 #endif
285
286 #ifdef HAS_MPI
287 // Return MPI data type
288 template<typename T>
289 struct dependent_false : std::false_type {};
290 template<typename T> static MPI_Datatype mpi_type()
291 {
292 static_assert(dependent_false<T>::value, "Unknown MPI type");
293 dolfin_error("MPI.h",
294 "perform MPI operation",
295 "MPI data type unknown");
296 return MPI_CHAR;
297 }
298 #endif
299
300 #ifdef HAS_MPI
301 // Maps some MPI_Op values to string
302 static std::map<MPI_Op, std::string> operation_map;
303 #endif
304
305 };
306
307 #ifdef HAS_MPI
308 // Specialisations for MPI_Datatypes
309 template<> inline MPI_Datatype MPI::mpi_type<float>() { return MPI_FLOAT; }
310 template<> inline MPI_Datatype MPI::mpi_type<double>() { return MPI_DOUBLE; }
311 template<> inline MPI_Datatype MPI::mpi_type<short int>() { return MPI_SHORT; }
312 template<> inline MPI_Datatype MPI::mpi_type<int>() { return MPI_INT; }
313 template<> inline MPI_Datatype MPI::mpi_type<long int>() { return MPI_LONG; }
314 template<> inline MPI_Datatype MPI::mpi_type<unsigned int>() { return MPI_UNSIGNED; }
315 template<> inline MPI_Datatype MPI::mpi_type<unsigned long int>() { return MPI_UNSIGNED_LONG; }
316 template<> inline MPI_Datatype MPI::mpi_type<long long>() { return MPI_LONG_LONG; }
317 #endif
318 //---------------------------------------------------------------------------
319 template<typename T>
320 void dolfin::MPI::broadcast(MPI_Comm comm, std::vector<T>& value,
321 unsigned int broadcaster)
322 {
323 #ifdef HAS_MPI
324 // Broadcast cast size
325 std::size_t bsize = value.size();
326 MPI_Bcast(&bsize, 1, mpi_type<std::size_t>(), broadcaster, comm);
327
328 // Broadcast
329 value.resize(bsize);
330 MPI_Bcast(const_cast<T*>(value.data()), bsize, mpi_type<T>(),
331 broadcaster, comm);
332 #endif
333 }
334 //---------------------------------------------------------------------------
335 template<typename T>
336 void dolfin::MPI::broadcast(MPI_Comm comm, T& value,
337 unsigned int broadcaster)
338 {
339 #ifdef HAS_MPI
340 MPI_Bcast(&value, 1, mpi_type<T>(), broadcaster, comm);
341 #endif
342 }
343 //---------------------------------------------------------------------------
344 template<typename T>
345 void dolfin::MPI::all_to_all(MPI_Comm comm,
346 std::vector<std::vector<T>>& in_values,
347 std::vector<std::vector<T>>& out_values)
348 {
349 #ifdef HAS_MPI
350 const std::size_t comm_size = MPI::size(comm);
351
352 // Data size per destination
353 dolfin_assert(in_values.size() == comm_size);
354 std::vector<int> data_size_send(comm_size);
355 std::vector<int> data_offset_send(comm_size + 1, 0);
356 for (std::size_t p = 0; p < comm_size; ++p)
357 {
358 data_size_send[p] = in_values[p].size();
359 data_offset_send[p + 1] = data_offset_send[p] + data_size_send[p];
360 }
361
362 // Get received data sizes
363 std::vector<int> data_size_recv(comm_size);
364 MPI_Alltoall(data_size_send.data(), 1, mpi_type<int>(),
365 data_size_recv.data(), 1, mpi_type<int>(),
366 comm);
367
368 // Pack data and build receive offset
369 std::vector<int> data_offset_recv(comm_size + 1 ,0);
370 std::vector<T> data_send(data_offset_send[comm_size]);
371 for (std::size_t p = 0; p < comm_size; ++p)
372 {
373 data_offset_recv[p + 1] = data_offset_recv[p] + data_size_recv[p];
374 std::copy(in_values[p].begin(), in_values[p].end(),
375 data_send.begin() + data_offset_send[p]);
376 }
377
378 // Send/receive data
379 std::vector<T> data_recv(data_offset_recv[comm_size]);
380 MPI_Alltoallv(data_send.data(), data_size_send.data(),
381 data_offset_send.data(), mpi_type<T>(),
382 data_recv.data(), data_size_recv.data(),
383 data_offset_recv.data(), mpi_type<T>(), comm);
384
385 // Repack data
386 out_values.resize(comm_size);
387 for (std::size_t p = 0; p < comm_size; ++p)
388 {
389 out_values[p].resize(data_size_recv[p]);
390 std::copy(data_recv.begin() + data_offset_recv[p],
391 data_recv.begin() + data_offset_recv[p + 1],
392 out_values[p].begin());
393 }
394 #else
395 dolfin_assert(in_values.size() == 1);
396 out_values = in_values;
397 #endif
398 }
399 //---------------------------------------------------------------------------
400 template<typename T>
401 void dolfin::MPI::all_to_all(MPI_Comm comm,
402 std::vector<std::vector<T>>& in_values,
403 std::vector<T>& out_values)
404 {
405 #ifdef HAS_MPI
406 const std::size_t comm_size = MPI::size(comm);
407
408 // Data size per destination
409 dolfin_assert(in_values.size() == comm_size);
410 std::vector<int> data_size_send(comm_size);
411 std::vector<int> data_offset_send(comm_size + 1, 0);
412 for (std::size_t p = 0; p < comm_size; ++p)
413 {
414 data_size_send[p] = in_values[p].size();
415 data_offset_send[p + 1] = data_offset_send[p] + data_size_send[p];
416 }
417
418 // Get received data sizes
419 std::vector<int> data_size_recv(comm_size);
420 MPI_Alltoall(data_size_send.data(), 1, mpi_type<int>(),
421 data_size_recv.data(), 1, mpi_type<int>(),
422 comm);
423
424 // Pack data and build receive offset
425 std::vector<int> data_offset_recv(comm_size + 1 ,0);
426 std::vector<T> data_send(data_offset_send[comm_size]);
427 for (std::size_t p = 0; p < comm_size; ++p)
428 {
429 data_offset_recv[p + 1] = data_offset_recv[p] + data_size_recv[p];
430 std::copy(in_values[p].begin(), in_values[p].end(),
431 data_send.begin() + data_offset_send[p]);
432 }
433
434 // Send/receive data
435 out_values.resize(data_offset_recv[comm_size]);
436 MPI_Alltoallv(data_send.data(), data_size_send.data(),
437 data_offset_send.data(), mpi_type<T>(),
438 out_values.data(), data_size_recv.data(),
439 data_offset_recv.data(), mpi_type<T>(), comm);
440
441 #else
442 dolfin_assert(in_values.size() == 1);
443 out_values = in_values[0];
444 #endif
445 }
446 //---------------------------------------------------------------------------
447#ifndef DOXYGEN_IGNORE
448 template<> inline
449 void dolfin::MPI::all_to_all(MPI_Comm comm,
450 std::vector<std::vector<bool>>& in_values,
451 std::vector<std::vector<bool>>& out_values)
452 {
453 #ifdef HAS_MPI
454 // Copy to short int
455 std::vector<std::vector<short int>> send(in_values.size());
456 for (std::size_t i = 0; i < in_values.size(); ++i)
457 send[i].assign(in_values[i].begin(), in_values[i].end());
458
459 // Communicate data
460 std::vector<std::vector<short int>> recv;
461 all_to_all(comm, send, recv);
462
463 // Copy back to bool
464 out_values.resize(recv.size());
465 for (std::size_t i = 0; i < recv.size(); ++i)
466 out_values[i].assign(recv[i].begin(), recv[i].end());
467 #else
468 dolfin_assert(in_values.size() == 1);
469 out_values = in_values;
470 #endif
471 }
472
473 template<> inline
474 void dolfin::MPI::all_to_all(MPI_Comm comm,
475 std::vector<std::vector<bool>>& in_values,
476 std::vector<bool>& out_values)
477 {
478 #ifdef HAS_MPI
479 // Copy to short int
480 std::vector<std::vector<short int>> send(in_values.size());
481 for (std::size_t i = 0; i < in_values.size(); ++i)
482 send[i].assign(in_values[i].begin(), in_values[i].end());
483
484 // Communicate data
485 std::vector<short int> recv;
486 all_to_all(comm, send, recv);
487
488 // Copy back to bool
489 out_values.assign(recv.begin(), recv.end());
490 #else
491 dolfin_assert(in_values.size() == 1);
492 out_values = in_values[0];
493 #endif
494 }
495
496#endif
497 //---------------------------------------------------------------------------
498 template<typename T>
499 void dolfin::MPI::scatter(MPI_Comm comm,
500 const std::vector<std::vector<T>>& in_values,
501 std::vector<T>& out_value,
502 unsigned int sending_process)
503 {
504 #ifdef HAS_MPI
505
506 // Scatter number of values to each process
507 const std::size_t comm_size = MPI::size(comm);
508 std::vector<int> all_num_values;
509 if (MPI::rank(comm) == sending_process)
510 {
511 dolfin_assert(in_values.size() == comm_size);
512 all_num_values.resize(comm_size);
513 for (std::size_t i = 0; i < comm_size; ++i)
514 all_num_values[i] = in_values[i].size();
515 }
516 int my_num_values = 0;
517 scatter(comm, all_num_values, my_num_values, sending_process);
518
519 // Prepare send buffer and offsets
520 std::vector<T> sendbuf;
521 std::vector<int> offsets;
522 if (MPI::rank(comm) == sending_process)
523 {
524 // Build offsets
525 offsets.resize(comm_size + 1, 0);
526 for (std::size_t i = 1; i <= comm_size; ++i)
527 offsets[i] = offsets[i - 1] + all_num_values[i - 1];
528
529 // Allocate send buffer and fill
530 const std::size_t n = std::accumulate(all_num_values.begin(),
531 all_num_values.end(), 0);
532 sendbuf.resize(n);
533 for (std::size_t p = 0; p < in_values.size(); ++p)
534 {
535 std::copy(in_values[p].begin(), in_values[p].end(),
536 sendbuf.begin() + offsets[p]);
537 }
538 }
539
540 // Scatter
541 out_value.resize(my_num_values);
542 MPI_Scatterv(const_cast<T*>(sendbuf.data()), all_num_values.data(),
543 offsets.data(), mpi_type<T>(),
544 out_value.data(), my_num_values,
545 mpi_type<T>(), sending_process, comm);
546 #else
547 dolfin_assert(sending_process == 0);
548 dolfin_assert(in_values.size() == 1);
549 out_value = in_values[0];
550 #endif
551 }
552 //---------------------------------------------------------------------------
553#ifndef DOXYGEN_IGNORE
554 template<> inline
555 void dolfin::MPI::scatter(MPI_Comm comm,
556 const std::vector<std::vector<bool>>& in_values,
557 std::vector<bool>& out_value,
558 unsigned int sending_process)
559 {
560 #ifdef HAS_MPI
561 // Copy data
562 std::vector<std::vector<short int>> in(in_values.size());
563 for (std::size_t i = 0; i < in_values.size(); ++i)
564 in[i] = std::vector<short int>(in_values[i].begin(), in_values[i].end());
565
566 std::vector<short int> out;
567 scatter(comm, in, out, sending_process);
568
569 out_value.resize(out.size());
570 std::copy(out.begin(), out.end(), out_value.begin());
571 #else
572 dolfin_assert(sending_process == 0);
573 dolfin_assert(in_values.size() == 1);
574 out_value = in_values[0];
575 #endif
576 }
577#endif
578 //---------------------------------------------------------------------------
579 template<typename T>
580 void dolfin::MPI::scatter(MPI_Comm comm,
581 const std::vector<T>& in_values,
582 T& out_value, unsigned int sending_process)
583 {
584 #ifdef HAS_MPI
585 if (MPI::rank(comm) == sending_process)
586 dolfin_assert(in_values.size() == MPI::size(comm));
587
588 MPI_Scatter(const_cast<T*>(in_values.data()), 1, mpi_type<T>(),
589 &out_value, 1, mpi_type<T>(), sending_process, comm);
590 #else
591 dolfin_assert(sending_process == 0);
592 dolfin_assert(in_values.size() == 1);
593 out_value = in_values[0];
594 #endif
595 }
596 //---------------------------------------------------------------------------
597 template<typename T>
598 void dolfin::MPI::gather(MPI_Comm comm,
599 const std::vector<T>& in_values,
600 std::vector<T>& out_values,
601 unsigned int receiving_process)
602 {
603 #ifdef HAS_MPI
604 const std::size_t comm_size = MPI::size(comm);
605
606 // Get data size on each process
607 std::vector<int> pcounts(comm_size);
608 const int local_size = in_values.size();
609 MPI_Gather(const_cast<int*>(&local_size), 1, mpi_type<int>(),
610 pcounts.data(), 1, mpi_type<int>(),
611 receiving_process, comm);
612
613 // Build offsets
614 std::vector<int> offsets(comm_size + 1, 0);
615 for (std::size_t i = 1; i <= comm_size; ++i)
616 offsets[i] = offsets[i - 1] + pcounts[i - 1];
617
618 const std::size_t n = std::accumulate(pcounts.begin(), pcounts.end(), 0);
619 out_values.resize(n);
620 MPI_Gatherv(const_cast<T*>(in_values.data()), in_values.size(),
621 mpi_type<T>(),
622 out_values.data(), pcounts.data(), offsets.data(),
623 mpi_type<T>(), receiving_process, comm);
624 #else
625 out_values = in_values;
626 #endif
627 }
628 //---------------------------------------------------------------------------
629 inline void dolfin::MPI::gather(MPI_Comm comm,
630 const std::string& in_values,
631 std::vector<std::string>& out_values,
632 unsigned int receiving_process)
633 {
634 #ifdef HAS_MPI
635 const std::size_t comm_size = MPI::size(comm);
636
637 // Get data size on each process
638 std::vector<int> pcounts(comm_size);
639 int local_size = in_values.size();
640 MPI_Gather(&local_size, 1, MPI_INT,
641 pcounts.data(), 1,MPI_INT,
642 receiving_process, comm);
643
644 // Build offsets
645 std::vector<int> offsets(comm_size + 1, 0);
646 for (std::size_t i = 1; i <= comm_size; ++i)
647 offsets[i] = offsets[i - 1] + pcounts[i - 1];
648
649 // Gather
650 const std::size_t n = std::accumulate(pcounts.begin(), pcounts.end(), 0);
651 std::vector<char> _out(n);
652 MPI_Gatherv(const_cast<char*>(in_values.data()), in_values.size(),
653 MPI_CHAR,
654 _out.data(), pcounts.data(), offsets.data(),
655 MPI_CHAR, receiving_process, comm);
656
657 // Rebuild
658 out_values.resize(comm_size);
659 for (std::size_t p = 0; p < comm_size; ++p)
660 {
661 out_values[p] = std::string(_out.begin() + offsets[p],
662 _out.begin() + offsets[p + 1]);
663 }
664 #else
665 out_values.clear();
666 out_values.push_back(in_values);
667 #endif
668 }
669 //---------------------------------------------------------------------------
670 inline void dolfin::MPI::all_gather(MPI_Comm comm,
671 const std::string& in_values,
672 std::vector<std::string>& out_values)
673 {
674 #ifdef HAS_MPI
675 const std::size_t comm_size = MPI::size(comm);
676
677 // Get data size on each process
678 std::vector<int> pcounts(comm_size);
679 int local_size = in_values.size();
680 MPI_Allgather(&local_size, 1, MPI_INT, pcounts.data(), 1, MPI_INT, comm);
681
682 // Build offsets
683 std::vector<int> offsets(comm_size + 1, 0);
684 for (std::size_t i = 1; i <= comm_size; ++i)
685 offsets[i] = offsets[i - 1] + pcounts[i - 1];
686
687 // Gather
688 const std::size_t n = std::accumulate(pcounts.begin(), pcounts.end(), 0);
689 std::vector<char> _out(n);
690 MPI_Allgatherv(const_cast<char*>(in_values.data()), in_values.size(),
691 MPI_CHAR, _out.data(), pcounts.data(), offsets.data(),
692 MPI_CHAR, comm);
693
694 // Rebuild
695 out_values.resize(comm_size);
696 for (std::size_t p = 0; p < comm_size; ++p)
697 {
698 out_values[p] = std::string(_out.begin() + offsets[p],
699 _out.begin() + offsets[p + 1]);
700 }
701 #else
702 out_values.clear();
703 out_values.push_back(in_values);
704 #endif
705 }
706 //-------------------------------------------------------------------------
707 template<typename T>
708 void dolfin::MPI::all_gather(MPI_Comm comm,
709 const std::vector<T>& in_values,
710 std::vector<T>& out_values)
711 {
712 #ifdef HAS_MPI
713 out_values.resize(in_values.size()*MPI::size(comm));
714 MPI_Allgather(const_cast<T*>(in_values.data()), in_values.size(),
715 mpi_type<T>(),
716 out_values.data(), in_values.size(), mpi_type<T>(),
717 comm);
718#else
719 out_values = in_values;
720 #endif
721 }
722 //---------------------------------------------------------------------------
723 template<typename T>
724 void dolfin::MPI::all_gather(MPI_Comm comm,
725 const std::vector<T>& in_values,
726 std::vector<std::vector<T>>& out_values)
727 {
728 #ifdef HAS_MPI
729 const std::size_t comm_size = MPI::size(comm);
730
731 // Get data size on each process
732 std::vector<int> pcounts;
733 const int local_size = in_values.size();
734 MPI::all_gather(comm, local_size, pcounts);
735 dolfin_assert(pcounts.size() == comm_size);
736
737 // Build offsets
738 std::vector<int> offsets(comm_size + 1, 0);
739 for (std::size_t i = 1; i <= comm_size; ++i)
740 offsets[i] = offsets[i - 1] + pcounts[i - 1];
741
742 // Gather data
743 const std::size_t n = std::accumulate(pcounts.begin(), pcounts.end(), 0);
744 std::vector<T> recvbuf(n);
745 MPI_Allgatherv(const_cast<T*>(in_values.data()), in_values.size(),
746 mpi_type<T>(),
747 recvbuf.data(), pcounts.data(), offsets.data(),
748 mpi_type<T>(), comm);
749
750 // Repack data
751 out_values.resize(comm_size);
752 for (std::size_t p = 0; p < comm_size; ++p)
753 {
754 out_values[p].resize(pcounts[p]);
755 for (int i = 0; i < pcounts[p]; ++i)
756 out_values[p][i] = recvbuf[offsets[p] + i];
757 }
758 #else
759 out_values.clear();
760 out_values.push_back(in_values);
761 #endif
762 }
763 //---------------------------------------------------------------------------
764 template<typename T>
765 void dolfin::MPI::all_gather(MPI_Comm comm, const T in_value,
766 std::vector<T>& out_values)
767 {
768 #ifdef HAS_MPI
769 out_values.resize(MPI::size(comm));
770 MPI_Allgather(const_cast<T*>(&in_value), 1, mpi_type<T>(),
771 out_values.data(), 1, mpi_type<T>(), comm);
772 #else
773 out_values.clear();
774 out_values.push_back(in_value);
775 #endif
776 }
777 //---------------------------------------------------------------------------
778 template<typename T, typename X>
779 T dolfin::MPI::all_reduce(MPI_Comm comm, const T& value, X op)
780 {
781 #ifdef HAS_MPI
782 T out;
783 MPI_Allreduce(const_cast<T*>(&value), &out, 1, mpi_type<T>(), op, comm);
784 return out;
785 #else
786 return value;
787 #endif
788 }
789 //---------------------------------------------------------------------------
790 template<typename T> T dolfin::MPI::max(MPI_Comm comm, const T& value)
791 {
792 #ifdef HAS_MPI
793 // Enforce cast to MPI_Op; this is needed because template dispatch may
794 // not recognize this is possible, e.g. C-enum to unsigned int in SGI MPT
795 MPI_Op op = static_cast<MPI_Op>(MPI_MAX);
796 return all_reduce(comm, value, op);
797 #else
798 return value;
799 #endif
800 }
801 //---------------------------------------------------------------------------
802 template<typename T> T dolfin::MPI::min(MPI_Comm comm, const T& value)
803 {
804 #ifdef HAS_MPI
805 // Enforce cast to MPI_Op; this is needed because template dispatch may
806 // not recognize this is possible, e.g. C-enum to unsigned int in SGI MPT
807 MPI_Op op = static_cast<MPI_Op>(MPI_MIN);
808 return all_reduce(comm, value, op);
809 #else
810 return value;
811 #endif
812 }
813 //---------------------------------------------------------------------------
814 template<typename T> T dolfin::MPI::sum(MPI_Comm comm, const T& value)
815 {
816 #ifdef HAS_MPI
817 // Enforce cast to MPI_Op; this is needed because template dispatch may
818 // not recognize this is possible, e.g. C-enum to unsigned int in SGI MPT
819 MPI_Op op = static_cast<MPI_Op>(MPI_SUM);
820 return all_reduce(comm, value, op);
821 #else
822 return value;
823 #endif
824 }
825 //---------------------------------------------------------------------------
826 template<typename T> T dolfin::MPI::avg(MPI_Comm comm, const T& value)
827 {
828 #ifdef HAS_MPI
829 dolfin_error("MPI.h",
830 "perform average reduction",
831 "Not implemented for this type");
832 #else
833 return value;
834 #endif
835 }
836 //---------------------------------------------------------------------------
837 template<typename T>
838 void dolfin::MPI::send_recv(MPI_Comm comm,
839 const std::vector<T>& send_value,
840 unsigned int dest, int send_tag,
841 std::vector<T>& recv_value,
842 unsigned int source, int recv_tag)
843 {
844 #ifdef HAS_MPI
845 std::size_t send_size = send_value.size();
846 std::size_t recv_size = 0;
847 MPI_Status mpi_status;
848 MPI_Sendrecv(&send_size, 1, mpi_type<std::size_t>(), dest, send_tag,
849 &recv_size, 1, mpi_type<std::size_t>(), source, recv_tag,
850 comm, &mpi_status);
851
852 recv_value.resize(recv_size);
853 MPI_Sendrecv(const_cast<T*>(send_value.data()), send_value.size(),
854 mpi_type<T>(), dest, send_tag,
855 recv_value.data(), recv_size, mpi_type<T>(),
856 source, recv_tag, comm, &mpi_status);
857 #else
858 dolfin_error("MPI.h",
859 "call MPI::send_recv",
860 "DOLFIN has been configured without MPI support");
861 #endif
862 }
863 //---------------------------------------------------------------------------
864 template<typename T>
865 void dolfin::MPI::send_recv(MPI_Comm comm,
866 const std::vector<T>& send_value,
867 unsigned int dest,
868 std::vector<T>& recv_value,
869 unsigned int source)
870 {
871 MPI::send_recv(comm, send_value, dest, 0, recv_value, source, 0);
872 }
873 //---------------------------------------------------------------------------
874 // Specialization for dolfin::log::Table class
875 // NOTE: This function is not trully "all_reduce", it reduces to rank 0
876 // and returns zero Table on other ranks.
877 #ifdef HAS_MPI
878 template<>
879 Table dolfin::MPI::all_reduce(MPI_Comm, const Table&, MPI_Op);
880 #endif
881 //---------------------------------------------------------------------------
882 // Specialization for dolfin::log::Table class
883 #ifdef HAS_MPI
884 template<>
885 Table dolfin::MPI::avg(MPI_Comm, const Table&);
886 #endif
887 //---------------------------------------------------------------------------
888}
889
890#endif
Definition MPI.h:54
Definition MPI.h:77
Comm(MPI_Comm comm)
Duplicate communicator and wrap duplicate.
Definition MPI.cpp:31
MPI_Comm comm() const
Return the underlying MPI_Comm object.
Definition MPI.cpp:117
unsigned int size() const
Definition MPI.cpp:73
unsigned int rank() const
Return process rank for the communicator.
Definition MPI.cpp:68
~Comm()
Destructor (frees wrapped communicator)
Definition MPI.cpp:51
void free()
Free (destroy) communicator. Calls function 'MPI_Comm_free'.
Definition MPI.cpp:56
void barrier() const
Set a barrier (synchronization point)
Definition MPI.cpp:84
void reset(MPI_Comm comm)
Definition MPI.cpp:91
Definition MPI.h:71
static MPI_Op MPI_AVG()
Definition MPI.cpp:359
static T all_reduce(MPI_Comm comm, const T &value, X op)
All reduce.
Definition MPI.h:779
static std::size_t global_offset(MPI_Comm comm, std::size_t range, bool exclusive)
Definition MPI.cpp:185
static void barrier(MPI_Comm comm)
Set a barrier (synchronization point)
Definition MPI.cpp:178
static void all_to_all(MPI_Comm comm, std::vector< std::vector< T > > &in_values, std::vector< std::vector< T > > &out_values)
Definition MPI.h:345
static T avg(MPI_Comm comm, const T &value)
Return average across comm; implemented only for T == Table.
Definition MPI.h:826
static void broadcast(MPI_Comm comm, std::vector< T > &value, unsigned int broadcaster=0)
Broadcast vector of value from broadcaster to all processes.
Definition MPI.h:320
static T min(MPI_Comm comm, const T &value)
Return global min value.
Definition MPI.h:802
static T max(MPI_Comm comm, const T &value)
Return global max value.
Definition MPI.h:790
static unsigned int size(MPI_Comm comm)
Definition MPI.cpp:154
static bool is_receiver(MPI_Comm comm)
Definition MPI.cpp:172
static T sum(MPI_Comm comm, const T &value)
Sum values and return sum.
Definition MPI.h:814
static unsigned int rank(MPI_Comm comm)
Return process rank for the communicator.
Definition MPI.cpp:143
static std::pair< std::int64_t, std::int64_t > local_range(MPI_Comm comm, std::int64_t N)
Definition MPI.cpp:201
static unsigned int index_owner(MPI_Comm comm, std::size_t index, std::size_t N)
Return which process owns index (inverse of local_range)
Definition MPI.cpp:230
static void gather(MPI_Comm comm, const std::vector< T > &in_values, std::vector< T > &out_values, unsigned int receiving_process=0)
Gather values on one process.
Definition MPI.h:598
static std::pair< std::int64_t, std::int64_t > compute_local_range(int process, std::int64_t N, int size)
Definition MPI.cpp:213
static void send_recv(MPI_Comm comm, const std::vector< T > &send_value, unsigned int dest, int send_tag, std::vector< T > &recv_value, unsigned int source, int recv_tag)
Send-receive data between processes (blocking)
Definition MPI.h:838
static bool is_broadcaster(MPI_Comm comm)
Definition MPI.cpp:166
static void scatter(MPI_Comm comm, const std::vector< std::vector< T > > &in_values, std::vector< T > &out_value, unsigned int sending_process=0)
Scatter vector in_values[i] to process i.
Definition MPI.h:499
static void all_gather(MPI_Comm comm, const std::vector< T > &in_values, std::vector< T > &out_values)
Definition MPI.h:708
Definition Table.h:50
Definition adapt.h:30
void begin(std::string msg,...)
Begin task (increase indentation level)
Definition log.cpp:153
void dolfin_error(std::string location, std::string task, std::string reason,...)
Definition log.cpp:129
void end()
End task (decrease indentation level)
Definition log.cpp:168
void assign(std::shared_ptr< Function > receiving_func, std::shared_ptr< const Function > assigning_func)
Definition assign.cpp:27