Actual source code: mpiutils.h
1: #ifndef MPIUTILS_H
2: #define MPIUTILS_H
4: #include <petscsys.h>
6: PETSC_EXTERN PetscErrorCode PetscGatherNumberOfMessages_Private(MPI_Comm, const PetscMPIInt[], const PetscInt[], PetscMPIInt *);
7: PETSC_EXTERN PetscErrorCode PetscGatherMessageLengths_Private(MPI_Comm, PetscMPIInt, PetscMPIInt, const PetscInt[], PetscMPIInt **, PetscInt **);
9: #if !defined(PETSC_HAVE_MPI_LARGE_COUNT) /* No matter PetscInt is 32-bit or 64-bit, without MPI large count we always do casting before MPI calls */
10: /* Cast PetscInt <a> to PetscMPIInt <b>, where <a> is likely used for the 'count' argument in MPI routines.
11: It is similar to PetscMPIIntCast() execept that here it returns an MPI error code.
12: */
13: static inline PetscMPIInt PetscMPIIntCast_Internal(PetscInt a, PetscMPIInt *b)
14: {
15: *b = (PetscMPIInt)(a);
16: if (PetscDefined(USE_64BIT_INDICIES) && PetscUnlikely(a > PETSC_MPI_INT_MAX)) return MPI_ERR_COUNT;
17: return MPI_SUCCESS;
18: }
20: static inline PetscMPIInt MPIU_Send(const void *buf, PetscInt count, MPI_Datatype datatype, PetscMPIInt dest, PetscMPIInt tag, MPI_Comm comm)
21: {
22: PetscMPIInt count2;
24: PetscFunctionBegin;
25: PetscCallMPI(PetscMPIIntCast_Internal(count, &count2));
26: PetscCallMPI(MPI_Send(buf, count2, datatype, dest, tag, comm));
27: PetscFunctionReturn(MPI_SUCCESS);
28: }
30: static inline PetscMPIInt MPIU_Send_init(const void *buf, PetscInt count, MPI_Datatype datatype, PetscMPIInt dest, PetscMPIInt tag, MPI_Comm comm, MPI_Request *request)
31: {
32: PetscMPIInt count2;
34: PetscFunctionBegin;
35: PetscCallMPI(PetscMPIIntCast_Internal(count, &count2));
36: PetscCallMPI(MPI_Send_init(buf, count2, datatype, dest, tag, comm, request));
37: PetscFunctionReturn(MPI_SUCCESS);
38: }
40: static inline PetscMPIInt MPIU_Isend(const void *buf, PetscInt count, MPI_Datatype datatype, PetscMPIInt dest, PetscMPIInt tag, MPI_Comm comm, MPI_Request *request)
41: {
42: PetscMPIInt count2;
44: PetscFunctionBegin;
45: PetscCallMPI(PetscMPIIntCast_Internal(count, &count2));
46: PetscCallMPI(MPI_Isend(buf, count2, datatype, dest, tag, comm, request));
47: PetscFunctionReturn(MPI_SUCCESS);
48: }
50: static inline PetscMPIInt MPIU_Recv(void *buf, PetscInt count, MPI_Datatype datatype, PetscMPIInt source, PetscMPIInt tag, MPI_Comm comm, MPI_Status *status)
51: {
52: PetscMPIInt count2;
54: PetscFunctionBegin;
55: PetscCallMPI(PetscMPIIntCast_Internal(count, &count2));
56: PetscCallMPI(MPI_Recv(buf, count2, datatype, source, tag, comm, status));
57: PetscFunctionReturn(MPI_SUCCESS);
58: }
60: static inline PetscMPIInt MPIU_Recv_init(void *buf, PetscInt count, MPI_Datatype datatype, PetscMPIInt source, PetscMPIInt tag, MPI_Comm comm, MPI_Request *request)
61: {
62: PetscMPIInt count2;
64: PetscFunctionBegin;
65: PetscCallMPI(PetscMPIIntCast_Internal(count, &count2));
66: PetscCallMPI(MPI_Recv_init(buf, count2, datatype, source, tag, comm, request));
67: PetscFunctionReturn(MPI_SUCCESS);
68: }
70: static inline PetscMPIInt MPIU_Irecv(void *buf, PetscInt count, MPI_Datatype datatype, PetscMPIInt source, PetscMPIInt tag, MPI_Comm comm, MPI_Request *request)
71: {
72: PetscMPIInt count2;
74: PetscFunctionBegin;
75: PetscCallMPI(PetscMPIIntCast_Internal(count, &count2));
76: PetscCallMPI(MPI_Irecv(buf, count2, datatype, source, tag, comm, request));
77: PetscFunctionReturn(MPI_SUCCESS);
78: }
79: #if defined(PETSC_HAVE_MPI_REDUCE_LOCAL)
80: static inline PetscMPIInt MPIU_Reduce_local(const void *inbuf, void *inoutbuf, PetscInt count, MPI_Datatype datatype, MPI_Op op)
81: {
82: PetscMPIInt count2;
84: PetscFunctionBegin;
85: PetscCallMPI(PetscMPIIntCast_Internal(count, &count2));
86: PetscCallMPI(MPI_Reduce_local(inbuf, inoutbuf, count, datatype, op));
87: PetscFunctionReturn(MPI_SUCCESS);
88: }
89: #endif
91: #elif defined(PETSC_USE_64BIT_INDICES)
92: #define MPIU_Send(buf, count, datatype, dest, tag, comm) MPI_Send_c(buf, count, datatype, dest, tag, comm)
93: #define MPIU_Send_init(buf, count, datatype, dest, tag, comm, request) MPI_Send_init_c(buf, count, datatype, dest, tag, comm, request)
94: #define MPIU_Isend(buf, count, datatype, dest, tag, comm, request) MPI_Isend_c(buf, count, datatype, dest, tag, comm, request)
95: #define MPIU_Recv(buf, count, datatype, source, tag, comm, status) MPI_Recv_c(buf, count, datatype, source, tag, comm, status)
96: #define MPIU_Recv_init(buf, count, datatype, source, tag, comm, request) MPI_Recv_init_c(buf, count, datatype, source, tag, comm, request)
97: #define MPIU_Irecv(buf, count, datatype, source, tag, comm, request) MPI_Irecv_c(buf, count, datatype, source, tag, comm, request)
98: #if defined(PETSC_HAVE_MPI_REDUCE_LOCAL)
99: #define MPIU_Reduce_local(inbuf, inoutbuf, count, datatype, op) MPI_Reduce_local_c(inbuf, inoutbuf, count, datatype, op)
100: #endif
101: #else
102: #define MPIU_Send(buf, count, datatype, dest, tag, comm) MPI_Send(buf, count, datatype, dest, tag, comm)
103: #define MPIU_Send_init(buf, count, datatype, dest, tag, comm, request) MPI_Send_init(buf, count, datatype, dest, tag, comm, request)
104: #define MPIU_Isend(buf, count, datatype, dest, tag, comm, request) MPI_Isend(buf, count, datatype, dest, tag, comm, request)
105: #define MPIU_Recv(buf, count, datatype, source, tag, comm, status) MPI_Recv(buf, count, datatype, source, tag, comm, status)
106: #define MPIU_Recv_init(buf, count, datatype, source, tag, comm, request) MPI_Recv_init(buf, count, datatype, source, tag, comm, request)
107: #define MPIU_Irecv(buf, count, datatype, source, tag, comm, request) MPI_Irecv(buf, count, datatype, source, tag, comm, request)
108: #if defined(PETSC_HAVE_MPI_REDUCE_LOCAL)
109: #define MPIU_Reduce_local(inbuf, inoutbuf, count, datatype, op) MPI_Reduce_local(inbuf, inoutbuf, count, datatype, op)
110: #endif
111: #endif
113: /* These APIs use arrays of MPI_Count/MPI_Aint */
114: #if defined(PETSC_HAVE_MPI_LARGE_COUNT) && defined(PETSC_USE_64BIT_INDICES)
115: #define MPIU_Neighbor_alltoallv(a, b, c, d, e, f, g, h, i) MPI_Neighbor_alltoallv_c(a, b, c, d, e, f, g, h, i)
116: #define MPIU_Ineighbor_alltoallv(a, b, c, d, e, f, g, h, i, j) MPI_Ineighbor_alltoallv_c(a, b, c, d, e, f, g, h, i, j)
117: #else
118: #define MPIU_Neighbor_alltoallv(a, b, c, d, e, f, g, h, i) MPI_Neighbor_alltoallv(a, b, c, d, e, f, g, h, i)
119: #define MPIU_Ineighbor_alltoallv(a, b, c, d, e, f, g, h, i, j) MPI_Ineighbor_alltoallv(a, b, c, d, e, f, g, h, i, j)
120: #endif
122: #endif