27 #include <queso/MpiComm.h>
28 #include <queso/Environment.h>
30 #ifdef QUESO_HAS_TRILINOS
31 #include <Epetra_MpiComm.h>
32 #include <Epetra_SerialComm.h>
41 #ifdef QUESO_HAS_TRILINOS
42 m_epetraComm( new Epetra_MpiComm(inputRawComm) ),
44 m_rawComm (inputRawComm),
50 int mpiRC = MPI_Comm_rank(inputRawComm,&
m_worldRank);
53 mpiRC = MPI_Comm_rank(inputRawComm,&
m_myPid);
56 mpiRC = MPI_Comm_size(inputRawComm,&
m_numProc);
68 #ifdef QUESO_HAS_TRILINOS
69 m_epetraComm( new Epetra_SerialComm() ),
71 m_rawComm (RawValue_MPI_COMM_SELF),
82 #ifdef QUESO_HAS_TRILINOS
93 #ifdef QUESO_HAS_TRILINOS
113 #ifdef QUESO_HAS_TRILINOS
120 #endif // QUESO_HAS_MPI
126 #ifdef QUESO_HAS_TRILINOS
135 #ifdef QUESO_HAS_TRILINOS
148 int mpiRC = MPI_Allreduce(sendbuf, recvbuf, count, datatype, op,
m_rawComm);
154 template <
typename T>
157 const char* whereMsg,
const char* whatMsg)
const
161 T * sendbuf_noconst =
const_cast<T *
>(sendbuf);
167 size_t dataTypeSize =
sizeof(T);
168 size_t dataTotal = dataTypeSize*count;
169 std::memcpy(recvbuf, sendbuf, dataTotal);
176 #ifdef QUESO_HAS_TRILINOS
195 int mpiRC = MPI_Bcast(buffer, count, datatype, root,
m_rawComm);
206 const char* whereMsg,
const char* whatMsg)
const
214 int mpiRC = MPI_Gather(sendbuf, sendcnt, sendtype,
215 recvbuf, recvcount, recvtype,
222 template <
typename T>
225 int root,
const char* whereMsg,
const char* whatMsg)
const
232 T * sendbuf_noconst =
const_cast<T *
>(sendbuf);
233 int mpiRC = MPI_Gather(sendbuf_noconst, sendcnt,
StandardType<T>(sendbuf),
240 size_t dataTypeSize =
sizeof(T);
241 size_t sendTotal = dataTypeSize*sendcnt;
242 size_t recvTotal = dataTypeSize*recvcount;
243 if (sendTotal != recvTotal) {
244 std::cerr <<
"MpiCommClass::Gather()"
245 <<
": sendTotal != recvTotal"
249 std::memcpy(recvbuf, sendbuf, sendTotal);
258 const char* whereMsg,
const char* whatMsg)
const
266 int mpiRC = MPI_Gatherv(sendbuf, sendcnt, sendtype,
267 recvbuf, recvcnts, displs, recvtype,
274 template <
typename T>
277 int * displs,
int root,
const char * whereMsg,
278 const char * whatMsg)
const
285 T * sendbuf_noconst =
const_cast<T *
>(sendbuf);
286 int mpiRC = MPI_Gatherv(sendbuf_noconst, sendcnt,
StandardType<T>(sendbuf),
287 recvbuf, recvcnts, displs,
293 size_t dataTypeSize =
sizeof(T);
294 size_t sendTotal = dataTypeSize*sendcnt;
295 size_t recvTotal = dataTypeSize*recvcnts[0];
296 if (sendTotal != recvTotal) {
297 std::cerr <<
"MpiCommClass::Gatherv()"
298 <<
": sendTotal != recvTotal"
302 std::memcpy(recvbuf, sendbuf, sendTotal);
309 const char* whereMsg,
const char* whatMsg)
const
313 int mpiRC = MPI_Recv(buf, count, datatype, source, tag,
m_rawComm, status);
322 const char* whereMsg,
const char* whatMsg)
const
326 int mpiRC = MPI_Send(buf, count, datatype, dest, tag,
m_rawComm);
337 for (
int i = 0; i < this->
NumProc(); ++i) {
338 if (i == this->
MyPID()) {
358 #ifdef QUESO_HAS_TRILINOS
370 #ifdef QUESO_HAS_TRILINOS
384 template void MpiComm::Allreduce<int>(
const int *,
390 template void MpiComm::Allreduce<char>(
const char *,
396 template void MpiComm::Allreduce<unsigned int>(
const unsigned int *,
402 template void MpiComm::Allreduce<double>(
const double *,
409 template void MpiComm::Gather<int>(
const int * sendbuf,
414 const char * whereMsg,
415 const char * whatMsg)
const;
417 template void MpiComm::Gather<char>(
const char * sendbuf,
422 const char * whereMsg,
423 const char * whatMsg)
const;
425 template void MpiComm::Gather<unsigned int>(
const unsigned int * sendbuf,
427 unsigned int * recvbuf,
430 const char * whereMsg,
431 const char * whatMsg)
const;
433 template void MpiComm::Gather<double>(
const double * sendbuf,
438 const char * whereMsg,
439 const char * whatMsg)
const;
441 template void MpiComm::Gatherv<int>(
const int * sendbuf,
447 const char * whereMsg,
448 const char * whatMsg)
const;
450 template void MpiComm::Gatherv<char>(
const char * sendbuf,
456 const char * whereMsg,
457 const char * whatMsg)
const;
459 template void MpiComm::Gatherv<unsigned int>(
const unsigned int * sendbuf,
461 unsigned int * recvbuf,
465 const char * whereMsg,
466 const char * whatMsg)
const;
468 template void MpiComm::Gatherv<double>(
const double * sendbuf,
474 const char * whereMsg,
475 const char * whatMsg)
const;
MPI_Comm RawType_MPI_Comm
const Epetra_Comm & epetraMpiComm() const
Extract MPI Communicator from a Epetra_MpiComm object.
MpiComm()
Default Constructor.
The QUESO MPI Communicator Class.
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.
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 rank of the MPI process in QUESO's full communicator.
int inter0Rank() const
Returns the process inter0 rank.
int m_myPid
Process ID of this process.
RawType_MPI_Comm m_rawComm
Embedded wrapped opaque MPI_Comm object.
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.
void Barrier() const
Pause every process in *this communicator until all the processes reach this point.
MPI_Datatype RawType_MPI_Datatype
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.
int subRank() const
Returns the rank of the MPI process in the sub-communicator subComm()
MPI_Status RawType_MPI_Status
int MyPID() const
Return my process ID.
Epetra_Comm * m_epetraComm
int NumProc() const
Returns total number of processes.
MonteCarloSGOptions::MonteCarloSGOptions(const BaseEnvironment &env, const char *prefix, const McOptionsValues &alternativeOptionsValues queso_require_equal_to_msg)(m_env.optionsInputFileName(), std::string(""), std::string("this constructor is incompatible with the existence of an options input file"))
int m_worldRank
World rank.
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.
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.
MpiComm & operator=(const MpiComm &rhs)
Assignment operator.
unsigned int subId() const
Access function to the number of each sub-environment Id: m_subId.
This (virtual) class sets up the environment underlying the use of the QUESO library by an executable...
unsigned int syncVerbosity() const
Access function to private attribute m_syncVerbosity.