OpenLB 1.7
Loading...
Searching...
No Matches
mpiManager.cpp
Go to the documentation of this file.
1/* This file is part of the OpenLB library
2 *
3 * Copyright (C) 2007 The OpenLB project
4 *
5 * This program is free software; you can redistribute it and/or
6 * modify it under the terms of the GNU General Public License
7 * as published by the Free Software Foundation; either version 2
8 * of the License, or (at your option) any later version.
9 *
10 * This program 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 General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public
16 * License along with this program; if not, write to the Free
17 * Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
18 * Boston, MA 02110-1301, USA.
19*/
20
21#include "mpiManager.h"
22
23#include <unistd.h>
24
25namespace olb {
26
27namespace singleton {
28
30{
31 static MpiManager instance;
32 return instance;
33}
34
35#ifdef PARALLEL_MODE_MPI
36
37MpiManager::MpiManager() : ok(false), clout(std::cout,"MpiManager")
38{ }
39
40MpiManager::~MpiManager()
41{
42 if (ok) {
43 MPI_Finalize();
44 ok = false;
45 }
46}
47
48void MpiManager::init(int *argc, char ***argv, bool verbose)
49{
50 int ok0{};
51 MPI_Initialized(&ok0);
52 if (ok0) {
53 return;
54 }
55 int ok1 = MPI_Init(argc, argv);
56 int ok2 = MPI_Comm_rank(MPI_COMM_WORLD, &taskId);
57 int ok3 = MPI_Comm_size(MPI_COMM_WORLD, &numTasks);
58 int ok4 = MPI_Comm_set_errhandler(MPI_COMM_WORLD, MPI_ERRORS_ARE_FATAL);
59 ok = (ok1 == MPI_SUCCESS && ok2 == MPI_SUCCESS && ok3 == MPI_SUCCESS && ok4 == MPI_SUCCESS);
60 if (verbose) {
61 clout << "Sucessfully initialized, numThreads=" << getSize() << std::endl;
62 }
63}
64
66{
67 return numTasks;
68}
69
71{
72 return taskId;
73}
74
76{
77 return 0;
78}
79
81{
82 return bossId() == getRank();
83}
84
85double MpiManager::getTime() const
86{
87 if (!ok) {
88 return 0.;
89 }
90 return MPI_Wtime();
91}
92
93void MpiManager::barrier(MPI_Comm comm)
94{
95 if (!ok) {
96 return;
97 }
98 MPI_Barrier(comm);
99}
100
101void MpiManager::synchronizeIO(unsigned tDelay, MPI_Comm comm)
102{
103 usleep(tDelay);
104 barrier(comm);
105}
106
107template <>
108void MpiManager::send<bool>(bool *buf, int count, int dest, int tag, MPI_Comm comm)
109{
110 if (!ok) {
111 return;
112 }
113 MPI_Send(static_cast<void*>(buf), count, MPI_BYTE, dest, tag, comm);
114}
115
116template <>
117void MpiManager::send<char>(char *buf, int count, int dest, int tag, MPI_Comm comm)
118{
119 if (!ok) {
120 return;
121 }
122 MPI_Send(static_cast<void*>(buf), count, MPI_CHAR, dest, tag, comm);
123}
124
125template <>
126void MpiManager::send<std::uint8_t>(std::uint8_t *buf, int count, int dest, int tag, MPI_Comm comm)
127{
128 if (!ok) {
129 return;
130 }
131 MPI_Send(static_cast<void*>(buf), count, MPI_BYTE, dest, tag, comm);
132}
133
134template <>
135void MpiManager::send<int>(int *buf, int count, int dest, int tag, MPI_Comm comm)
136{
137 if (!ok) {
138 return;
139 }
140 MPI_Send(static_cast<void*>(buf), count, MPI_INT, dest, tag, comm);
141}
142
143template <>
144void MpiManager::send<float>(float *buf, int count, int dest, int tag, MPI_Comm comm)
145{
146 if (!ok) {
147 return;
148 }
149 MPI_Send(static_cast<void*>(buf), count, MPI_FLOAT, dest, tag, comm);
150}
151
152template <>
153void MpiManager::send<double>(double *buf, int count, int dest, int tag, MPI_Comm comm)
154{
155 if (!ok) {
156 return;
157 }
158 MPI_Send(static_cast<void*>(buf), count, MPI_DOUBLE, dest, tag, comm);
159}
160
161template <>
162void MpiManager::sendInit<double>(double *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
163{
164 if (ok) {
165 MPI_Send_init(buf, count, MPI_DOUBLE, dest, tag, comm, request);
166 }
167}
168
169template <>
170void MpiManager::sendInit<std::size_t>(std::size_t *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
171{
172 if (ok) {
173 MPI_Send_init(buf, count, MPI_UNSIGNED_LONG, dest, tag, comm, request);
174 }
175}
176
177template <>
178void MpiManager::sendInit<std::uint32_t>(std::uint32_t *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
179{
180 if (ok) {
181 MPI_Send_init(buf, count, MPI_UNSIGNED, dest, tag, comm, request);
182 }
183}
184
185template <>
186void MpiManager::sendInit<std::uint8_t>(std::uint8_t *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
187{
188 if (ok) {
189 MPI_Send_init(buf, count, MPI_BYTE, dest, tag, comm, request);
190 }
191}
192
193template <>
194void MpiManager::sendInit<int>(int *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
195{
196 if (ok) {
197 MPI_Send_init(buf, count, MPI_INT, dest, tag, comm, request);
198 }
199}
200
201template <>
202void MpiManager::sendInit<bool>(bool *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
203{
204 if (ok) {
205 MPI_Send_init(static_cast<void*>(buf), count, MPI_BYTE, dest, tag, comm, request);
206 }
207}
208
209template <>
210void MpiManager::iSend<bool>
211(bool *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
212{
213 if (ok) {
214 MPI_Isend(static_cast<void*>(buf), count, MPI_BYTE, dest, tag, comm, request);
215 }
216}
217
218template <>
219void MpiManager::iSend<char>
220(char *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
221{
222 if (ok) {
223 MPI_Isend(static_cast<void*>(buf), count, MPI_CHAR, dest, tag, comm, request);
224 }
225}
226
227template <>
228void MpiManager::iSend<std::uint8_t>
229(std::uint8_t *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
230{
231 if (ok) {
232 MPI_Isend(static_cast<void*>(buf), count, MPI_BYTE, dest, tag, comm, request);
233 }
234}
235
236template <>
237void MpiManager::iSend<int>
238(int *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
239{
240 if (ok) {
241 MPI_Isend(static_cast<void*>(buf), count, MPI_INT, dest, tag, comm, request);
242 }
243}
244
245template <>
246void MpiManager::iSend<std::size_t>
247(std::size_t *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
248{
249 if (ok) {
250 MPI_Isend(static_cast<void*>(buf), count, MPI_UNSIGNED_LONG, dest, tag, comm, request);
251 }
252}
253
254template <>
255void MpiManager::iSend<std::uint32_t>
256(std::uint32_t *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
257{
258 if (ok) {
259 MPI_Isend(static_cast<void*>(buf), count, MPI_UNSIGNED, dest, tag, comm, request);
260 }
261}
262
263template <>
264void MpiManager::iSend<float>
265(float *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
266{
267 if (ok) {
268 MPI_Isend(static_cast<void*>(buf), count, MPI_FLOAT, dest, tag, comm, request);
269 }
270}
271
272template <>
273void MpiManager::iSend<double>
274(double *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
275{
276 if (ok) {
277 MPI_Isend(static_cast<void*>(buf), count, MPI_DOUBLE, dest, tag, comm, request);
278 }
279}
280
281template <>
282void MpiManager::iSend<long double>
283(long double *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
284{
285 if (ok) {
286 MPI_Isend(static_cast<void*>(buf), count, MPI_LONG_DOUBLE, dest, tag, comm, request);
287 }
288}
289
290
291template <>
292void MpiManager::ibSend<bool>
293(bool *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
294{
295 if (ok) {
296 MPI_Ibsend(static_cast<void*>(buf), count, MPI_BYTE, dest, tag, comm, request);
297 }
298}
299
300template <>
301void MpiManager::ibSend<char>
302(char *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
303{
304 if (ok) {
305 MPI_Ibsend(static_cast<void*>(buf), count, MPI_CHAR, dest, tag, comm, request);
306 }
307}
308
309template <>
310void MpiManager::ibSend<int>
311(int *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
312{
313 if (ok) {
314 MPI_Ibsend(static_cast<void*>(buf), count, MPI_INT, dest, tag, comm, request);
315 }
316}
317
318template <>
319void MpiManager::ibSend<float>
320(float *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
321{
322 if (ok) {
323 MPI_Ibsend(static_cast<void*>(buf), count, MPI_FLOAT, dest, tag, comm, request);
324 }
325}
326
327template <>
328void MpiManager::ibSend<double>
329(double *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
330{
331 if (ok) {
332 MPI_Ibsend(static_cast<void*>(buf), count, MPI_DOUBLE, dest, tag, comm, request);
333 }
334}
335
336std::size_t MpiManager::probeReceiveSize(int source, MPI_Datatype type, int tag, MPI_Comm comm)
337{
338 MPI_Status status;
339 if (MPI_Probe(source, tag, comm, &status) == MPI_SUCCESS) {
340 int requestSize;
341 MPI_Get_count(&status, type, &requestSize);
342 if (requestSize == MPI_UNDEFINED) {
343 throw std::runtime_error("MPI_UNDEFINED in probeReceiveSize(" + std::to_string(source) + "," + std::to_string(tag) + ")" + " ranks " + std::to_string(source) + " -> " + std::to_string(singleton::mpi().getRank()));
344 }
345 return requestSize;
346 } else {
347 throw std::runtime_error("MPI_Probe failed in probeReceiveSize");
348 }
349}
350
351template <>
352std::size_t MpiManager::probeReceiveSize<std::uint32_t>(int source, int tag, MPI_Comm comm)
353{
354 return probeReceiveSize(source, MPI_UNSIGNED, tag, comm);
355}
356
357template <>
358std::size_t MpiManager::probeReceiveSize<std::uint64_t>(int source, int tag, MPI_Comm comm)
359{
360 return probeReceiveSize(source, MPI_UNSIGNED_LONG, tag, comm);
361}
362
363template <>
364void MpiManager::receive<bool>(bool *buf, int count, int source, int tag, MPI_Comm comm)
365{
366 if (!ok) {
367 return;
368 }
369 MPI_Status status;
370 MPI_Recv(static_cast<void*>(buf), count, MPI_BYTE, source, tag, comm, &status);
371}
372
373
374template <>
375void MpiManager::receive<char>(char *buf, int count, int source, int tag, MPI_Comm comm)
376{
377 if (!ok) {
378 return;
379 }
380 MPI_Status status;
381 MPI_Recv(static_cast<void*>(buf), count, MPI_CHAR, source, tag, comm, &status);
382}
383
384template <>
385void MpiManager::receive<std::uint8_t>(std::uint8_t *buf, int count, int source, int tag, MPI_Comm comm)
386{
387 if (!ok) {
388 return;
389 }
390 MPI_Status status;
391 MPI_Recv(static_cast<std::uint8_t*>(buf), count, MPI_BYTE, source, tag, comm, &status);
392}
393
394template <>
395void MpiManager::receive<int>(int *buf, int count, int source, int tag, MPI_Comm comm)
396{
397 if (!ok) {
398 return;
399 }
400 MPI_Status status;
401 MPI_Recv(static_cast<void*>(buf), count, MPI_INT, source, tag, comm, &status);
402}
403
404template <>
405void MpiManager::receive<std::size_t>(std::size_t *buf, int count, int source, int tag, MPI_Comm comm)
406{
407 if (!ok) {
408 return;
409 }
410 MPI_Status status;
411 MPI_Recv(static_cast<void*>(buf), count, MPI_UNSIGNED_LONG, source, tag, comm, &status);
412}
413
414template <>
415void MpiManager::receive<std::uint32_t>(std::uint32_t *buf, int count, int source, int tag, MPI_Comm comm)
416{
417 if (!ok) {
418 return;
419 }
420 MPI_Status status;
421 MPI_Recv(static_cast<void*>(buf), count, MPI_UNSIGNED, source, tag, comm, &status);
422}
423
424template <>
425void MpiManager::receive<float>(float *buf, int count, int source, int tag, MPI_Comm comm)
426{
427 if (!ok) {
428 return;
429 }
430 MPI_Status status;
431 MPI_Recv(static_cast<void*>(buf), count, MPI_FLOAT, source, tag, comm, &status);
432}
433
434template <>
435void MpiManager::receive<double>(double *buf, int count, int source, int tag, MPI_Comm comm)
436{
437 if (!ok) {
438 return;
439 }
440 MPI_Status status;
441 MPI_Recv(static_cast<void*>(buf), count, MPI_DOUBLE, source, tag, comm, &status);
442}
443
444template <>
445void MpiManager::receive<long double>(long double *buf, int count, int source, int tag, MPI_Comm comm)
446{
447 if (!ok) {
448 return;
449 }
450 MPI_Status status;
451 MPI_Recv(static_cast<void*>(buf), count, MPI_LONG_DOUBLE, source, tag, comm, &status);
452}
453
454template <>
455void MpiManager::sendToMaster<bool>(bool* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
456{
457 if (!ok) {
458 return;
459 }
460 if (iAmRoot && !isMainProcessor()) {
461 send(sendBuf, sendCount, 0);
462 }
463 if (isMainProcessor() && !iAmRoot) {
464 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
465 }
466}
467
468template <>
469void MpiManager::sendToMaster<char>(char* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
470{
471 if (!ok) {
472 return;
473 }
474 if (iAmRoot && !isMainProcessor()) {
475 send(sendBuf, sendCount, 0);
476 }
477 if (isMainProcessor() && !iAmRoot) {
478 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
479 }
480}
481
482template <>
483void MpiManager::sendToMaster<int>(int* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
484{
485 if (!ok) {
486 return;
487 }
488 if (iAmRoot && !isMainProcessor()) {
489 send(sendBuf, sendCount, 0);
490 }
491 if (isMainProcessor() && !iAmRoot) {
492 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
493 }
494}
495
496template <>
497void MpiManager::sendToMaster<float>(float* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
498{
499 if (!ok) {
500 return;
501 }
502 if (iAmRoot && !isMainProcessor()) {
503 send(sendBuf, sendCount, 0);
504 }
505 if (isMainProcessor() && !iAmRoot) {
506 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
507 }
508}
509
510template <>
511void MpiManager::sendToMaster<double>(double* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
512{
513 if (!ok) {
514 return;
515 }
516 if (iAmRoot && !isMainProcessor()) {
517 send(sendBuf, sendCount, 0);
518 }
519 if (isMainProcessor() && !iAmRoot) {
520 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
521 }
522}
523
524template <>
525void MpiManager::recvInit<double>(double *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
526{
527 if (ok) {
528 MPI_Recv_init(buf, count, MPI_DOUBLE, dest, tag, comm, request);
529 }
530}
531
532template <>
533void MpiManager::recvInit<int>(int *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
534{
535 if (ok) {
536 MPI_Recv_init(buf, count, MPI_INT, dest, tag, comm, request);
537 }
538}
539
540template <>
541void MpiManager::recvInit<std::uint8_t>(std::uint8_t *buf, int count, int dest, MPI_Request* request, int tag, MPI_Comm comm)
542{
543 if (ok) {
544 MPI_Recv_init(buf, count, MPI_BYTE, dest, tag, comm, request);
545 }
546}
547
548template <>
549void MpiManager::iRecv<bool>(bool *buf, int count, int source, MPI_Request* request, int tag, MPI_Comm comm)
550{
551 if (ok) {
552 MPI_Irecv(static_cast<void*>(buf), count, MPI_BYTE, source, tag, comm, request);
553 }
554}
555
556template <>
557void MpiManager::iRecv<char>(char *buf, int count, int source, MPI_Request* request, int tag, MPI_Comm comm)
558{
559 if (ok) {
560 MPI_Irecv(static_cast<void*>(buf), count, MPI_CHAR, source, tag, comm, request);
561 }
562}
563
564template <>
565void MpiManager::iRecv<int>(int *buf, int count, int source, MPI_Request* request, int tag, MPI_Comm comm)
566{
567 if (ok) {
568 MPI_Irecv(static_cast<void*>(buf), count, MPI_INT, source, tag, comm, request);
569 }
570}
571
572template <>
573void MpiManager::iRecv<float>(float *buf, int count, int source, MPI_Request* request, int tag, MPI_Comm comm)
574{
575 if (ok) {
576 MPI_Irecv(static_cast<void*>(buf), count, MPI_FLOAT, source, tag, comm, request);
577 }
578}
579
580template <>
581void MpiManager::iRecv<double>(double *buf, int count, int source, MPI_Request* request, int tag, MPI_Comm comm)
582{
583 if (ok) {
584 MPI_Irecv(static_cast<void*>(buf), count, MPI_DOUBLE, source, tag, comm, request);
585 }
586}
587
588template <>
589void MpiManager::sendRecv<bool>
590(bool *sendBuf, bool *recvBuf, int count, int dest, int source, int tag, MPI_Comm comm)
591{
592 if (!ok) {
593 return;
594 }
595 MPI_Status status;
596 MPI_Sendrecv(static_cast<void*>(sendBuf),
597 count,
598 MPI_BYTE, dest, tag,
599 static_cast<void*>(recvBuf),
600 count,
601 MPI_BYTE, source, tag, comm, &status);
602}
603
604template <>
605void MpiManager::sendRecv<char>
606(char *sendBuf, char *recvBuf, int count, int dest, int source, int tag, MPI_Comm comm)
607{
608 if (!ok) {
609 return;
610 }
611 MPI_Status status;
612 MPI_Sendrecv(static_cast<void*>(sendBuf),
613 count,
614 MPI_CHAR, dest, tag,
615 static_cast<void*>(recvBuf),
616 count,
617 MPI_CHAR, source, tag, comm, &status);
618}
619
620template <>
621void MpiManager::sendRecv<int>
622(int *sendBuf, int *recvBuf, int count, int dest, int source, int tag, MPI_Comm comm)
623{
624 if (!ok) {
625 return;
626 }
627 MPI_Status status;
628 MPI_Sendrecv(static_cast<void*>(sendBuf),
629 count,
630 MPI_INT, dest, tag,
631 static_cast<void*>(recvBuf),
632 count,
633 MPI_INT, source, tag, comm, &status);
634}
635
636template <>
637void MpiManager::sendRecv<float>
638(float *sendBuf, float *recvBuf, int count, int dest, int source, int tag, MPI_Comm comm)
639{
640 if (!ok) {
641 return;
642 }
643 MPI_Status status;
644 MPI_Sendrecv(static_cast<void*>(sendBuf),
645 count,
646 MPI_FLOAT, dest, tag,
647 static_cast<void*>(recvBuf),
648 count,
649 MPI_FLOAT, source, tag, comm, &status);
650}
651
652template <>
653void MpiManager::sendRecv<long>
654(long *sendBuf, long *recvBuf, int count, int dest, int source, int tag, MPI_Comm comm)
655{
656 if (!ok) {
657 return;
658 }
659 MPI_Status status;
660 MPI_Sendrecv(static_cast<void*>(sendBuf),
661 count,
662 MPI_LONG, dest, tag,
663 static_cast<void*>(recvBuf),
664 count,
665 MPI_LONG, source, tag, comm, &status);
666}
667
668template <>
669void MpiManager::sendRecv<double>
670(double *sendBuf, double *recvBuf, int count, int dest, int source, int tag, MPI_Comm comm)
671{
672 if (!ok) {
673 return;
674 }
675 MPI_Status status;
676 MPI_Sendrecv(static_cast<void*>(sendBuf),
677 count,
678 MPI_DOUBLE, dest, tag,
679 static_cast<void*>(recvBuf),
680 count,
681 MPI_DOUBLE, source, tag, comm, &status);
682}
683
684template <>
685void MpiManager::sendRecv<long double>
686(long double *sendBuf, long double *recvBuf, int count, int dest, int source, int tag, MPI_Comm comm)
687{
688 if (!ok) {
689 return;
690 }
691 MPI_Status status;
692 MPI_Sendrecv(static_cast<void*>(sendBuf),
693 count,
694 MPI_LONG_DOUBLE, dest, tag,
695 static_cast<void*>(recvBuf),
696 count,
697 MPI_LONG_DOUBLE, source, tag, comm, &status);
698}
699
700template <>
701void MpiManager::scatterv<bool>(bool* sendBuf, int* sendCounts, int* displs,
702 bool* recvBuf, int recvCount, int root, MPI_Comm comm)
703{
704 if (!ok) {
705 return;
706 }
707 MPI_Scatterv(static_cast<void*>(sendBuf),
708 sendCounts, displs, MPI_BYTE,
709 static_cast<void*>(recvBuf),
710 recvCount, MPI_BYTE, root, comm);
711}
712
713template <>
714void MpiManager::scatterv<char>(char* sendBuf, int* sendCounts, int* displs,
715 char* recvBuf, int recvCount, int root, MPI_Comm comm)
716{
717 if (!ok) {
718 return;
719 }
720 MPI_Scatterv(static_cast<void*>(sendBuf),
721 sendCounts, displs, MPI_CHAR,
722 static_cast<void*>(recvBuf),
723 recvCount, MPI_CHAR, root, comm);
724}
725
726template <>
727void MpiManager::scatterv<int>(int *sendBuf, int* sendCounts, int* displs,
728 int* recvBuf, int recvCount, int root, MPI_Comm comm)
729{
730 if (!ok) {
731 return;
732 }
733 MPI_Scatterv(static_cast<void*>(sendBuf),
734 sendCounts, displs, MPI_INT,
735 static_cast<void*>(recvBuf),
736 recvCount, MPI_INT, root, comm);
737}
738
739template <>
740void MpiManager::scatterv<float>(float *sendBuf, int* sendCounts, int* displs,
741 float* recvBuf, int recvCount, int root, MPI_Comm comm)
742{
743 if (!ok) {
744 return;
745 }
746 MPI_Scatterv(static_cast<void*>(sendBuf),
747 sendCounts, displs, MPI_FLOAT,
748 static_cast<void*>(recvBuf),
749 recvCount, MPI_FLOAT, root, comm);
750}
751
752template <>
753void MpiManager::scatterv<double>(double *sendBuf, int* sendCounts, int* displs,
754 double* recvBuf, int recvCount, int root, MPI_Comm comm)
755{
756 if (!ok) {
757 return;
758 }
759 MPI_Scatterv(static_cast<void*>(sendBuf),
760 sendCounts, displs, MPI_DOUBLE,
761 static_cast<void*>(recvBuf),
762 recvCount, MPI_DOUBLE, root, comm);
763}
764
765template <>
766void MpiManager::gather<int>(int* sendBuf, int sendCount,
767 int* recvBuf, int recvCount,
768 int root, MPI_Comm comm)
769{
770 if (!ok) {
771 return;
772 }
773 MPI_Gather(static_cast<void*>(sendBuf), sendCount, MPI_INT,
774 static_cast<void*>(recvBuf), recvCount, MPI_INT,
775 root, comm);
776}
777
778template <>
779void MpiManager::gatherv<bool>(bool* sendBuf, int sendCount,
780 bool* recvBuf, int* recvCounts, int* displs,
781 int root, MPI_Comm comm)
782{
783 if (!ok) {
784 return;
785 }
786 MPI_Gatherv(static_cast<void*>(sendBuf), sendCount, MPI_BYTE,
787 static_cast<void*>(recvBuf), recvCounts, displs, MPI_BYTE,
788 root, comm);
789}
790
791template <>
792void MpiManager::gatherv<char>(char* sendBuf, int sendCount,
793 char* recvBuf, int* recvCounts, int* displs,
794 int root, MPI_Comm comm)
795{
796 if (!ok) {
797 return;
798 }
799 MPI_Gatherv(static_cast<void*>(sendBuf), sendCount, MPI_CHAR,
800 static_cast<void*>(recvBuf), recvCounts, displs, MPI_CHAR,
801 root, comm);
802}
803
804template <>
805void MpiManager::gatherv<int>(int* sendBuf, int sendCount,
806 int* recvBuf, int* recvCounts, int* displs,
807 int root, MPI_Comm comm)
808{
809 if (!ok) {
810 return;
811 }
812 MPI_Gatherv(static_cast<void*>(sendBuf), sendCount, MPI_INT,
813 static_cast<void*>(recvBuf), recvCounts, displs, MPI_INT,
814 root, comm);
815}
816
817template <>
818void MpiManager::gatherv<float>(float* sendBuf, int sendCount,
819 float* recvBuf, int* recvCounts, int* displs,
820 int root, MPI_Comm comm)
821{
822 if (!ok) {
823 return;
824 }
825 MPI_Gatherv(static_cast<void*>(sendBuf), sendCount, MPI_FLOAT,
826 static_cast<void*>(recvBuf), recvCounts, displs, MPI_FLOAT,
827 root, comm);
828}
829
830template <>
831void MpiManager::gatherv<double>(double* sendBuf, int sendCount,
832 double* recvBuf, int* recvCounts, int* displs,
833 int root, MPI_Comm comm)
834{
835 if (!ok) {
836 return;
837 }
838 MPI_Gatherv(static_cast<void*>(sendBuf), sendCount, MPI_DOUBLE,
839 static_cast<void*>(recvBuf), recvCounts, displs, MPI_DOUBLE,
840 root, comm);
841}
842
843template <>
844void MpiManager::gatherv<std::size_t>(std::size_t* sendBuf, int sendCount,
845 std::size_t* recvBuf, int* recvCounts, int* displs,
846 int root, MPI_Comm comm)
847{
848 if (!ok) {
849 return;
850 }
851 MPI_Gatherv(static_cast<void*>(sendBuf), sendCount, MPI_UNSIGNED_LONG,
852 static_cast<void*>(recvBuf), recvCounts, displs, MPI_UNSIGNED_LONG,
853 root, comm);
854}
855
856template <>
857void MpiManager::bCast<bool>(bool* sendBuf, int sendCount, int root, MPI_Comm comm)
858{
859 if (!ok) {
860 return;
861 }
862 MPI_Bcast(static_cast<void*>(sendBuf),
863 sendCount, MPI_BYTE, root, comm);
864}
865
866template <>
867void MpiManager::bCast<char>(char* sendBuf, int sendCount, int root, MPI_Comm comm)
868{
869 if (!ok) {
870 return;
871 }
872 MPI_Bcast(static_cast<void*>(sendBuf),
873 sendCount, MPI_CHAR, root, comm);
874}
875
876template <>
877void MpiManager::bCast<unsigned char>(unsigned char* sendBuf, int sendCount, int root, MPI_Comm comm)
878{
879 if (!ok) {
880 return;
881 }
882 MPI_Bcast(static_cast<void*>(sendBuf),
883 sendCount, MPI_UNSIGNED_CHAR, root, comm);
884}
885
886template <>
887void MpiManager::bCast<int>(int* sendBuf, int sendCount, int root, MPI_Comm comm)
888{
889 if (!ok) {
890 return;
891 }
892 MPI_Bcast(static_cast<void*>(sendBuf),
893 sendCount, MPI_INT, root, comm);
894}
895
896template <>
897void MpiManager::bCast<unsigned long>(unsigned long* sendBuf, int sendCount, int root, MPI_Comm comm)
898{
899 if (!ok) {
900 return;
901 }
902 MPI_Bcast(static_cast<void*>(sendBuf),
903 sendCount, MPI_UNSIGNED_LONG, root, comm);
904}
905
906template <>
907void MpiManager::bCast<float>(float* sendBuf, int sendCount, int root, MPI_Comm comm)
908{
909 if (!ok) {
910 return;
911 }
912 MPI_Bcast(static_cast<void*>(sendBuf),
913 sendCount, MPI_FLOAT, root, comm);
914}
915
916template <>
917void MpiManager::bCast<double>(double* sendBuf, int sendCount, int root, MPI_Comm comm)
918{
919 if (!ok) {
920 return;
921 }
922 MPI_Bcast(static_cast<void*>(sendBuf),
923 sendCount, MPI_DOUBLE, root, comm);
924}
925
926
927template <>
928void MpiManager::bCast<std::string>(std::string* sendBuf, int sendCount, int root, MPI_Comm comm)
929{
930 if (!ok) {
931 return;
932 }
933 int length = (int) sendBuf->size();
934 MPI_Bcast(static_cast<void*>(&length), 1, MPI_INT, root, comm);
935 char* buffer = new char[length+1];
936 if (getRank()==root) {
937 std::copy(sendBuf->c_str(), sendBuf->c_str()+length+1, buffer);
938 }
939 MPI_Bcast(static_cast<void*>(buffer), length+1, MPI_CHAR, root, comm);
940 if (getRank()!=root) {
941 *sendBuf = buffer;
942 }
943 delete [] buffer;
944}
945
946template <>
947void MpiManager::bCast<bool>(bool& sendVal, int root, MPI_Comm comm)
948{
949 if (!ok) {
950 return;
951 }
952 MPI_Bcast(&sendVal, 1, MPI_BYTE, root, comm);
953}
954template <>
955void MpiManager::bCast<char>(char& sendVal, int root, MPI_Comm comm)
956{
957 if (!ok) {
958 return;
959 }
960 MPI_Bcast(&sendVal, 1, MPI_CHAR, root, comm);
961}
962
963template <>
964void MpiManager::bCast<unsigned char>(unsigned char& sendVal, int root, MPI_Comm comm)
965{
966 if (!ok) {
967 return;
968 }
969 MPI_Bcast(&sendVal, 1, MPI_UNSIGNED_CHAR, root, comm);
970}
971
972template <>
973void MpiManager::bCast<int>(int& sendVal, int root, MPI_Comm comm)
974{
975 if (!ok) {
976 return;
977 }
978 MPI_Bcast(&sendVal, 1, MPI_INT, root, comm);
979}
980
981template <>
982void MpiManager::bCast<unsigned long>(unsigned long& sendVal, int root, MPI_Comm comm)
983{
984 if (!ok) {
985 return;
986 }
987 MPI_Bcast(&sendVal, 1, MPI_UNSIGNED_LONG, root, comm);
988}
989
990template <>
991void MpiManager::bCast<float>(float& sendVal, int root, MPI_Comm comm)
992{
993 if (!ok) {
994 return;
995 }
996 MPI_Bcast(&sendVal, 1, MPI_FLOAT, root, comm);
997}
998
999template <>
1000void MpiManager::bCast<double>(double& sendVal, int root, MPI_Comm comm)
1001{
1002 if (!ok) {
1003 return;
1004 }
1005 MPI_Bcast(&sendVal, 1, MPI_DOUBLE, root, comm);
1006}
1007
1008template <>
1009void MpiManager::bCastThroughMaster<bool>(bool* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
1010{
1011 if (!ok) {
1012 return;
1013 }
1014 if (iAmRoot && !isMainProcessor()) {
1015 send(sendBuf, sendCount, 0);
1016 }
1017 if (isMainProcessor() && !iAmRoot) {
1018 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
1019 }
1020 bCast(sendBuf, sendCount, 0);
1021}
1022
1023template <>
1024void MpiManager::bCastThroughMaster<char>(char* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
1025{
1026 if (!ok) {
1027 return;
1028 }
1029 if (iAmRoot && !isMainProcessor()) {
1030 send(sendBuf, sendCount, 0);
1031 }
1032 if (isMainProcessor() && !iAmRoot) {
1033 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
1034 }
1035 bCast(sendBuf, sendCount, 0);
1036}
1037
1038template <>
1039void MpiManager::bCastThroughMaster<int>(int* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
1040{
1041 if (!ok) {
1042 return;
1043 }
1044 if (iAmRoot && !isMainProcessor()) {
1045 send(sendBuf, sendCount, 0);
1046 }
1047 if (isMainProcessor() && !iAmRoot) {
1048 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
1049 }
1050 bCast(sendBuf, sendCount, 0);
1051}
1052
1053template <>
1054void MpiManager::bCastThroughMaster<float>(float* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
1055{
1056 if (!ok) {
1057 return;
1058 }
1059 if (iAmRoot && !isMainProcessor()) {
1060 send(sendBuf, sendCount, 0);
1061 }
1062 if (isMainProcessor() && !iAmRoot) {
1063 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
1064 }
1065 bCast(sendBuf, sendCount, 0);
1066}
1067
1068template <>
1069void MpiManager::bCastThroughMaster<double>(double* sendBuf, int sendCount, bool iAmRoot, MPI_Comm comm)
1070{
1071 if (!ok) {
1072 return;
1073 }
1074 if (iAmRoot && !isMainProcessor()) {
1075 send(sendBuf, sendCount, 0);
1076 }
1077 if (isMainProcessor() && !iAmRoot) {
1078 receive(sendBuf, sendCount, MPI_ANY_SOURCE);
1079 }
1080 bCast(sendBuf, sendCount, 0);
1081}
1082
1083template <>
1084void MpiManager::reduce<bool>(bool& sendVal, bool& recvVal, MPI_Op op, int root, MPI_Comm comm)
1085{
1086 if (!ok) {
1087 return;
1088 }
1089 MPI_Reduce(static_cast<void*>(&sendVal),
1090 static_cast<void*>(&recvVal), 1, MPI_BYTE, op, root, comm);
1091}
1092
1093template <>
1094void MpiManager::reduce<char>(char& sendVal, char& recvVal, MPI_Op op, int root, MPI_Comm comm)
1095{
1096 if (!ok) {
1097 return;
1098 }
1099 MPI_Reduce(static_cast<void*>(&sendVal),
1100 static_cast<void*>(&recvVal), 1, MPI_CHAR, op, root, comm);
1101}
1102
1103template <>
1104void MpiManager::reduce<int>(int& sendVal, int& recvVal, MPI_Op op, int root, MPI_Comm comm)
1105{
1106 if (!ok) {
1107 return;
1108 }
1109 MPI_Reduce(static_cast<void*>(&sendVal),
1110 static_cast<void*>(&recvVal), 1, MPI_INT, op, root, comm);
1111}
1112
1113template <>
1114void MpiManager::reduce<float>(float& sendVal, float& recvVal, MPI_Op op, int root, MPI_Comm comm)
1115{
1116 if (!ok) {
1117 return;
1118 }
1119 MPI_Reduce(static_cast<void*>(&sendVal),
1120 static_cast<void*>(&recvVal), 1, MPI_FLOAT, op, root, comm);
1121}
1122
1123template <>
1124void MpiManager::reduce<double>(double& sendVal, double& recvVal, MPI_Op op, int root, MPI_Comm comm)
1125{
1126 if (!ok) {
1127 return;
1128 }
1129 MPI_Reduce(static_cast<void*>(&sendVal),
1130 static_cast<void*>(&recvVal), 1, MPI_DOUBLE, op, root, comm);
1131}
1132
1133template <>
1134void MpiManager::reduce<std::size_t>(std::size_t& sendVal, std::size_t& recvVal, MPI_Op op, int root, MPI_Comm comm)
1135{
1136 if (!ok) {
1137 return;
1138 }
1139 MPI_Reduce(static_cast<void*>(&sendVal),
1140 static_cast<void*>(&recvVal), 1, MPI_UNSIGNED_LONG, op, root, comm);
1141}
1142
1143template <>
1144void MpiManager::reduceVect<char>(std::vector<char>& sendVal, std::vector<char>& recvVal,
1145 MPI_Op op, int root, MPI_Comm comm)
1146{
1147 if (!ok) {
1148 return;
1149 }
1150 MPI_Reduce(static_cast<void*>(&(sendVal[0])),
1151 static_cast<void*>(&(recvVal[0])),
1152 sendVal.size(), MPI_CHAR, op, root, comm);
1153}
1154
1155template <>
1156void MpiManager::reduceVect<int>(std::vector<int>& sendVal, std::vector<int>& recvVal,
1157 MPI_Op op, int root, MPI_Comm comm)
1158{
1159 if (!ok) {
1160 return;
1161 }
1162 MPI_Reduce(static_cast<void*>(&(sendVal[0])),
1163 static_cast<void*>(&(recvVal[0])),
1164 sendVal.size(), MPI_INT, op, root, comm);
1165}
1166
1167template <>
1168void MpiManager::reduceVect<float>(std::vector<float>& sendVal, std::vector<float>& recvVal,
1169 MPI_Op op, int root, MPI_Comm comm)
1170{
1171 if (!ok) {
1172 return;
1173 }
1174 MPI_Reduce(static_cast<void*>(&(sendVal[0])),
1175 static_cast<void*>(&(recvVal[0])),
1176 sendVal.size(), MPI_FLOAT, op, root, comm);
1177}
1178
1179template <>
1180void MpiManager::reduceVect<double>(std::vector<double>& sendVal, std::vector<double>& recvVal,
1181 MPI_Op op, int root, MPI_Comm comm)
1182{
1183 if (!ok) {
1184 return;
1185 }
1186 MPI_Reduce(static_cast<void*>(&(sendVal[0])),
1187 static_cast<void*>(&(recvVal[0])),
1188 sendVal.size(), MPI_DOUBLE, op, root, comm);
1189}
1190
1191template <>
1192void MpiManager::reduceAndBcast<bool>(bool& reductVal, MPI_Op op, int root, MPI_Comm comm)
1193{
1194 if (!ok) {
1195 return;
1196 }
1197 char recvVal;
1198 MPI_Reduce(static_cast<void*>(&reductVal), static_cast<void*>(&recvVal), 1, MPI_BYTE, op, root, comm);
1199 reductVal = recvVal;
1200 MPI_Bcast(&reductVal, 1, MPI_BYTE, root, comm);
1201
1202}
1203
1204template <>
1205void MpiManager::reduceAndBcast<char>(char& reductVal, MPI_Op op, int root, MPI_Comm comm)
1206{
1207 if (!ok) {
1208 return;
1209 }
1210 char recvVal;
1211 MPI_Reduce(&reductVal, &recvVal, 1, MPI_CHAR, op, root, comm);
1212 reductVal = recvVal;
1213 MPI_Bcast(&reductVal, 1, MPI_CHAR, root, comm);
1214
1215}
1216
1217template <>
1218void MpiManager::reduceAndBcast<int>(int& reductVal, MPI_Op op, int root, MPI_Comm comm)
1219{
1220 if (!ok) {
1221 return;
1222 }
1223 int recvVal;
1224 MPI_Reduce(&reductVal, &recvVal, 1, MPI_INT, op, root, comm);
1225 reductVal = recvVal;
1226 MPI_Bcast(&reductVal, 1, MPI_INT, root, comm);
1227
1228}
1229
1230template <>
1231void MpiManager::reduceAndBcast<float>(float& reductVal, MPI_Op op, int root, MPI_Comm comm)
1232{
1233 if (!ok) {
1234 return;
1235 }
1236 float recvVal;
1237 MPI_Reduce(&reductVal, &recvVal, 1, MPI_FLOAT, op, root, comm);
1238 reductVal = recvVal;
1239 MPI_Bcast(&reductVal, 1, MPI_FLOAT, root, comm);
1240
1241}
1242
1243template <>
1244void MpiManager::reduceAndBcast<double>(double& reductVal, MPI_Op op, int root, MPI_Comm comm)
1245{
1246 if (!ok) {
1247 return;
1248 }
1249 double recvVal;
1250 MPI_Reduce(&reductVal, &recvVal, 1, MPI_DOUBLE, op, root, comm);
1251 reductVal = recvVal;
1252 MPI_Bcast(&reductVal, 1, MPI_DOUBLE, root, comm);
1253
1254}
1255
1256template <>
1257void MpiManager::reduceAndBcast<long double>(long double& reductVal, MPI_Op op, int root, MPI_Comm comm)
1258{
1259 if (!ok) {
1260 return;
1261 }
1262 long double recvVal;
1263 MPI_Reduce(&reductVal, &recvVal, 1, MPI_LONG_DOUBLE, op, root, comm);
1264 reductVal = recvVal;
1265 MPI_Bcast(&reductVal, 1, MPI_LONG_DOUBLE, root, comm);
1266
1267}
1268
1269template <>
1270void MpiManager::reduceAndBcast<long>(long& reductVal, MPI_Op op, int root, MPI_Comm comm)
1271{
1272 if (!ok) {
1273 return;
1274 }
1275 long recvVal;
1276 MPI_Reduce(&reductVal, &recvVal, 1, MPI_LONG, op, root, comm);
1277 reductVal = recvVal;
1278 MPI_Bcast(&reductVal, 1, MPI_LONG, root, comm);
1279
1280}
1281
1282template <>
1283void MpiManager::reduceAndBcast<unsigned long>(unsigned long& reductVal, MPI_Op op, int root, MPI_Comm comm)
1284{
1285 if (!ok) {
1286 return;
1287 }
1288 unsigned long recvVal;
1289 MPI_Reduce(&reductVal, &recvVal, 1, MPI_UNSIGNED_LONG, op, root, comm);
1290 reductVal = recvVal;
1291 MPI_Bcast(&reductVal, 1, MPI_UNSIGNED_LONG, root, comm);
1292
1293}
1294
1295void MpiManager::wait(MPI_Request* request, MPI_Status* status)
1296{
1297 if (!ok) {
1298 return;
1299 }
1300 MPI_Wait(request, status);
1301}
1302
1304{
1305 if (!ok || mpiNbHelper.get_size() == 0) {
1306 return;
1307 }
1308 MPI_Waitall(mpiNbHelper.get_size(), mpiNbHelper.get_mpiRequest(), mpiNbHelper.get_mpiStatus());
1309}
1310
1311
1315
1317{
1318 std::swap(_size, rhs._size);
1319 std::swap(_mpiRequest, rhs._mpiRequest);
1320 std::swap(_mpiStatus, rhs._mpiStatus);
1321}
1322
1324{
1325 free();
1326 _mpiRequest.reset(new MPI_Request[n] { });
1327 _mpiStatus.reset(new MPI_Status[n] { });
1328 _size = n;
1329}
1330
1332{
1333 _size = 0;
1334}
1335
1337{
1338 return _size;
1339}
1340
1342{
1343 OLB_PRECONDITION(size_t(i) < _size);
1344 return &_mpiRequest[i];
1345}
1346
1348{
1349 OLB_PRECONDITION(size_t(i) < _size);
1350 return &_mpiStatus[i];
1351}
1352
1354{
1355 MPI_Start(get_mpiRequest(i));
1356}
1357
1359{
1360 MPI_Wait(get_mpiRequest(i), get_mpiStatus(i));
1361}
1362
1364{
1365 int done;
1366 MPI_Test(get_mpiRequest(i), &done, MPI_STATUS_IGNORE);
1367 return done;
1368}
1369
1370#endif // PARALLEL_MODE_MPI
1371
1372} // namespace singleton
1373
1374} // namespace olb
Wrapper functions that simplify the use of MPI.
Definition mpiManager.h:90
void wait(MPI_Request *request, MPI_Status *status)
Complete a non-blocking MPI operation.
void send(T *buf, int count, int dest, int tag=0, MPI_Comm comm=MPI_COMM_WORLD)
Sends data at *buf, blocking.
void bCast(T *sendBuf, int sendCount, int root=0, MPI_Comm comm=MPI_COMM_WORLD)
Broadcast data from one processor to multiple processors.
int getSize() const
Returns the number of processes.
double getTime() const
Returns universal MPI-time in seconds.
bool isMainProcessor() const
Tells whether current processor is main processor.
void synchronizeIO(unsigned tDelay=100, MPI_Comm comm=MPI_COMM_WORLD)
Synchronizes the processes and wait to ensure correct cout order.
int getRank() const
Returns the process ID.
std::size_t probeReceiveSize(int source, MPI_Datatype type, int tag=0, MPI_Comm comm=MPI_COMM_WORLD)
Probe size of incoming message.
void waitAll(MpiNonBlockingHelper &mpiNbHelper)
Complete a series of non-blocking MPI operations.
void barrier(MPI_Comm comm=MPI_COMM_WORLD)
Synchronizes the processes.
void init(int *argc, char ***argv, bool verbose=true)
Initializes the mpi manager.
int bossId() const
Returns process ID of main processor.
void receive(T *buf, int count, int source, int tag=0, MPI_Comm comm=MPI_COMM_WORLD)
Receives data at *buf, blocking.
Helper class for non blocking MPI communication.
Definition mpiManager.h:51
void allocate(unsigned i)
Allocates memory.
MPI_Status * get_mpiStatus(int i=0) const
Get the specified status object.
MPI_Request * get_mpiRequest(int i=0) const
Get the specified request object.
void swap(MpiNonBlockingHelper &rhs)
Swap method.
unsigned get_size() const
Returns the size of the vector _mpiRequest/_mpiStatus.
Wrapper functions that simplify the use of MPI.
MpiManager & mpi()
Top level namespace for all of OpenLB.
#define OLB_PRECONDITION(COND)
Definition olbDebug.h:46