27 #include <queso/MpiComm.h>
28 #include <queso/Environment.h>
36 #ifdef QUESO_HAS_TRILINOS
37 m_epetraMpiComm( new Epetra_MpiComm(inputRawComm) ),
39 m_rawComm (inputRawComm),
45 int mpiRC = MPI_Comm_rank(inputRawComm,&
m_worldRank);
48 mpiRC = MPI_Comm_rank(inputRawComm,&
m_myPid);
51 mpiRC = MPI_Comm_size(inputRawComm,&
m_numProc);
63 #ifdef QUESO_HAS_TRILINOS
64 m_epetraMpiComm( new Epetra_MpiComm(MPI_COMM_SELF) ),
77 #ifdef QUESO_HAS_TRILINOS
88 #ifdef QUESO_HAS_TRILINOS
89 delete m_epetraMpiComm;
90 m_epetraMpiComm = NULL;
107 #ifdef QUESO_HAS_TRILINOS
108 return m_epetraMpiComm->Comm();
116 #ifdef QUESO_HAS_TRILINOS
117 return m_epetraMpiComm->MyPID();
125 #ifdef QUESO_HAS_TRILINOS
126 return m_epetraMpiComm->NumProc();
138 int mpiRC = MPI_Allreduce(sendbuf, recvbuf, count, datatype, op,
m_rawComm);
144 template <
typename T>
147 const char* whereMsg,
const char* whatMsg)
const
151 T * sendbuf_noconst =
const_cast<T *
>(sendbuf);
157 size_t dataTypeSize =
sizeof(T);
158 size_t dataTotal = dataTypeSize*count;
159 std::memcpy(recvbuf, sendbuf, dataTotal);
166 #ifdef QUESO_HAS_TRILINOS
167 return m_epetraMpiComm->Barrier();
185 int mpiRC = MPI_Bcast(buffer, count, datatype, root,
m_rawComm);
196 const char* whereMsg,
const char* whatMsg)
const
204 int mpiRC = MPI_Gather(sendbuf, sendcnt, sendtype,
205 recvbuf, recvcount, recvtype,
212 template <
typename T>
215 int root,
const char* whereMsg,
const char* whatMsg)
const
222 T * sendbuf_noconst =
const_cast<T *
>(sendbuf);
223 int mpiRC = MPI_Gather(sendbuf_noconst, sendcnt,
StandardType<T>(sendbuf),
230 size_t dataTypeSize =
sizeof(T);
231 size_t sendTotal = dataTypeSize*sendcnt;
232 size_t recvTotal = dataTypeSize*recvcount;
233 if (sendTotal != recvTotal) {
234 std::cerr <<
"MpiCommClass::Gather()"
235 <<
": sendTotal != recvTotal"
239 std::memcpy(recvbuf, sendbuf, sendTotal);
248 const char* whereMsg,
const char* whatMsg)
const
256 int mpiRC = MPI_Gatherv(sendbuf, sendcnt, sendtype,
257 recvbuf, recvcnts, displs, recvtype,
264 template <
typename T>
267 int * displs,
int root,
const char * whereMsg,
268 const char * whatMsg)
const
275 T * sendbuf_noconst =
const_cast<T *
>(sendbuf);
276 int mpiRC = MPI_Gatherv(sendbuf_noconst, sendcnt,
StandardType<T>(sendbuf),
277 recvbuf, recvcnts, displs,
283 size_t dataTypeSize =
sizeof(T);
284 size_t sendTotal = dataTypeSize*sendcnt;
285 size_t recvTotal = dataTypeSize*recvcnts[0];
286 if (sendTotal != recvTotal) {
287 std::cerr <<
"MpiCommClass::Gatherv()"
288 <<
": sendTotal != recvTotal"
292 std::memcpy(recvbuf, sendbuf, sendTotal);
299 const char* whereMsg,
const char* whatMsg)
const
303 int mpiRC = MPI_Recv(buf, count, datatype, source, tag,
m_rawComm, status);
312 const char* whereMsg,
const char* whatMsg)
const
316 int mpiRC = MPI_Send(buf, count, datatype, dest, tag,
m_rawComm);
327 for (
int i = 0; i < this->
NumProc(); ++i) {
328 if (i == this->
MyPID()) {
348 #ifdef QUESO_HAS_TRILINOS
349 const Epetra_MpiComm&
350 MpiComm::epetraMpiComm()
const
352 return *m_epetraMpiComm;
360 #ifdef QUESO_HAS_TRILINOS
361 delete m_epetraMpiComm;
362 m_epetraMpiComm =
new Epetra_MpiComm(*src.m_epetraMpiComm);
374 template void MpiComm::Allreduce<int>(
const int *,
380 template void MpiComm::Allreduce<char>(
const char *,
386 template void MpiComm::Allreduce<unsigned int>(
const unsigned int *,
392 template void MpiComm::Allreduce<double>(
const double *,
399 template void MpiComm::Gather<int>(
const int * sendbuf,
404 const char * whereMsg,
405 const char * whatMsg)
const;
407 template void MpiComm::Gather<char>(
const char * sendbuf,
412 const char * whereMsg,
413 const char * whatMsg)
const;
415 template void MpiComm::Gather<unsigned int>(
const unsigned int * sendbuf,
417 unsigned int * recvbuf,
420 const char * whereMsg,
421 const char * whatMsg)
const;
423 template void MpiComm::Gather<double>(
const double * sendbuf,
428 const char * whereMsg,
429 const char * whatMsg)
const;
431 template void MpiComm::Gatherv<int>(
const int * sendbuf,
437 const char * whereMsg,
438 const char * whatMsg)
const;
440 template void MpiComm::Gatherv<char>(
const char * sendbuf,
446 const char * whereMsg,
447 const char * whatMsg)
const;
449 template void MpiComm::Gatherv<unsigned int>(
const unsigned int * sendbuf,
451 unsigned int * recvbuf,
455 const char * whereMsg,
456 const char * whatMsg)
const;
458 template void MpiComm::Gatherv<double>(
const double * sendbuf,
464 const char * whereMsg,
465 const char * whatMsg)
const;
MpiComm()
Default Constructor.
int m_worldRank
World rank.
void Gather(void *sendbuf, int sendcnt, RawType_MPI_Datatype sendtype, void *recvbuf, int recvcount, RawType_MPI_Datatype recvtype, int root, const char *whereMsg, const char *whatMsg) const
Gather values from each process to collect on all processes.
int subRank() const
Access function for sub-rank.
RawType_MPI_Comm Comm() const
Extract MPI Communicator from a MpiComm object.
void Bcast(void *buffer, int count, RawType_MPI_Datatype datatype, int root, const char *whereMsg, const char *whatMsg) const
Broadcast values from the root process to the slave processes.
void syncPrintDebugMsg(const char *msg, unsigned int msgVerbosity, unsigned int numUSecs) const
Synchronizes all the processes and print debug message.
void copy(const MpiComm &src)
Copies from an existing MpiComm instance.
void Recv(void *buf, int count, RawType_MPI_Datatype datatype, int source, int tag, RawType_MPI_Status *status, const char *whereMsg, const char *whatMsg) const
Blocking receive of data from this process to another process.
int fullRank() const
Returns the process full rank.
#define queso_deprecated()
unsigned int syncVerbosity() const
Access function to private attribute m_syncVerbosity.
int m_myPid
Process ID of this process.
int inter0Rank() const
Returns the process inter0 rank.
#define queso_require_equal_to_msg(expr1, expr2, msg)
void Barrier() const
Pause every process in *this communicator until all the processes reach this point.
#define RawValue_MPI_COMM_SELF
This (virtual) class sets up the environment underlying the use of the QUESO library by an executable...
unsigned int subId() const
Access function to the number of each sub-environment Id: m_subId.
int MyPID() const
Return my process ID.
int NumProc() const
Returns total number of processes.
RawType_MPI_Comm m_rawComm
Embedded wrapped opaque MPI_Comm object.
void Gatherv(void *sendbuf, int sendcnt, RawType_MPI_Datatype sendtype, void *recvbuf, int *recvcnts, int *displs, RawType_MPI_Datatype recvtype, int root, const char *whereMsg, const char *whatMsg) const
Gathers into specified locations from all processes in a group.
void Allreduce(void *sendbuf, void *recvbuf, int count, RawType_MPI_Datatype datatype, RawType_MPI_Op op, const char *whereMsg, const char *whatMsg) const
Combines values from all processes and distributes the result back to all processes.
MpiComm & operator=(const MpiComm &rhs)
Assignment operator.
The QUESO MPI Communicator Class.
const BaseEnvironment & m_env
void Send(void *buf, int count, RawType_MPI_Datatype datatype, int dest, int tag, const char *whereMsg, const char *whatMsg) const
Possibly blocking send of data from this process to another process.