Actual source code: mpimesg.c
2: #include petsc.h
7: /*@C
8: PetscGatherNumberOfMessages - Computes the number of messages a node expects to receive
10: Collective on MPI_Comm
12: Input Parameters:
13: + comm - Communicator
14: . iflags - an array of integers of length sizeof(comm). A '1' in ilengths[i] represent a
15: message from current node to ith node. Optionally PETSC_NULL
16: - ilengths - Non zero ilengths[i] represent a message to i of length ilengths[i].
17: Optionally PETSC_NULL.
19: Output Parameters:
20: . nrecvs - number of messages received
22: Level: developer
24: Concepts: mpi utility
26: Notes:
27: With this info, the correct message lengths can be determined using
28: PetscGatherMessageLengths()
30: Either iflags or ilengths should be provided. If iflags is not
31: provided (PETSC_NULL) it can be computed from ilengths. If iflags is
32: provided, ilengths is not required.
34: .seealso: PetscGatherMessageLengths()
35: @*/
36: PetscErrorCode PetscGatherNumberOfMessages(MPI_Comm comm,PetscMPIInt *iflags,PetscMPIInt *ilengths,PetscMPIInt *nrecvs)
37: {
38: PetscMPIInt size,rank,*recv_buf,i,*iflags_local,*iflags_localm;
43: MPI_Comm_size(comm,&size);
44: MPI_Comm_rank(comm,&rank);
46: PetscMalloc2(size,PetscMPIInt,&recv_buf,size,PetscMPIInt,&iflags_localm);
48: /* If iflags not provided, compute iflags from ilengths */
49: if (!iflags) {
50: if (!ilengths) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Either iflags or ilengths should be provided");
51: iflags_local = iflags_localm;
52: for (i=0; i<size; i++) {
53: if (ilengths[i]) iflags_local[i] = 1;
54: else iflags_local[i] = 0;
55: }
56: } else {
57: iflags_local = iflags;
58: }
60: /* Post an allreduce to determine the numer of messages the current node will receive */
61: MPI_Allreduce(iflags_local,recv_buf,size,MPI_INT,MPI_SUM,comm);
62: *nrecvs = recv_buf[rank];
64: PetscFree2(recv_buf,iflags_localm);
65: return(0);
66: }
71: /*@C
72: PetscGatherMessageLengths - Computes info about messages that a MPI-node will receive,
73: including (from-id,length) pairs for each message.
75: Collective on MPI_Comm
77: Input Parameters:
78: + comm - Communicator
79: . nsends - number of messages that are to be sent.
80: . nrecvs - number of messages being received
81: - ilengths - an array of integers of length sizeof(comm)
82: a non zero ilengths[i] represent a message to i of length ilengths[i]
85: Output Parameters:
86: + onodes - list of node-ids from which messages are expected
87: - olengths - corresponding message lengths
89: Level: developer
91: Concepts: mpi utility
93: Notes:
94: With this info, the correct MPI_Irecv() can be posted with the correct
95: from-id, with a buffer with the right amount of memory required.
97: The calling function deallocates the memory in onodes and olengths
99: To determine nrecevs, one can use PetscGatherNumberOfMessages()
101: .seealso: PetscGatherNumberOfMessages()
102: @*/
103: PetscErrorCode PetscGatherMessageLengths(MPI_Comm comm,PetscMPIInt nsends,PetscMPIInt nrecvs,PetscMPIInt *ilengths,PetscMPIInt **onodes,PetscMPIInt **olengths)
104: {
106: PetscMPIInt size,tag,i,j;
107: MPI_Request *s_waits,*r_waits;
108: MPI_Status *w_status;
111: MPI_Comm_size(comm,&size);
112: PetscCommGetNewTag(comm,&tag);
114: /* cannot use PetscMalloc3() here because in the call to MPI_Waitall() they MUST be contiguous */
115: PetscMalloc2(nrecvs+nsends,MPI_Request,&r_waits,nrecvs+nsends,MPI_Status,&w_status);
116: s_waits = r_waits+nrecvs;
118: /* Post the Irecv to get the message length-info */
119: PetscMalloc(nrecvs*sizeof(PetscMPIInt),olengths);
120: for (i=0; i<nrecvs; i++) {
121: MPI_Irecv((*olengths)+i,1,MPI_INT,MPI_ANY_SOURCE,tag,comm,r_waits+i);
122: }
124: /* Post the Isends with the message length-info */
125: for (i=0,j=0; i<size; ++i) {
126: if (ilengths[i]) {
127: MPI_Isend(ilengths+i,1,MPI_INT,i,tag,comm,s_waits+j);
128: j++;
129: }
130: }
132: /* Post waits on sends and receivs */
133: MPI_Waitall(nrecvs+nsends,r_waits,w_status);
134:
135: /* Pack up the received data */
136: PetscMalloc(nrecvs*sizeof(PetscMPIInt),onodes);
137: for (i=0; i<nrecvs; ++i) {
138: (*onodes)[i] = w_status[i].MPI_SOURCE;
139: }
140: PetscFree2(r_waits,w_status);
141: return(0);
142: }
145: /*@C
146: PetscGatherMessageLengths2 - Computes info about messages that a MPI-node will receive,
147: including (from-id,length) pairs for each message. Same functionality as PetscGatherMessageLengths()
148: except it takes TWO ilenths and output TWO olengths.
150: Collective on MPI_Comm
152: Input Parameters:
153: + comm - Communicator
154: . nsends - number of messages that are to be sent.
155: . nrecvs - number of messages being received
156: - ilengths1, ilengths2 - array of integers of length sizeof(comm)
157: a non zero ilengths[i] represent a message to i of length ilengths[i]
159: Output Parameters:
160: + onodes - list of node-ids from which messages are expected
161: - olengths1, olengths2 - corresponding message lengths
163: Level: developer
165: Concepts: mpi utility
167: Notes:
168: With this info, the correct MPI_Irecv() can be posted with the correct
169: from-id, with a buffer with the right amount of memory required.
171: The calling function deallocates the memory in onodes and olengths
173: To determine nrecevs, one can use PetscGatherNumberOfMessages()
175: .seealso: PetscGatherMessageLengths() and PetscGatherNumberOfMessages()
176: @*/
177: PetscErrorCode PetscGatherMessageLengths2(MPI_Comm comm,PetscMPIInt nsends,PetscMPIInt nrecvs,PetscMPIInt *ilengths1,PetscMPIInt *ilengths2,PetscMPIInt **onodes,PetscMPIInt **olengths1,PetscMPIInt **olengths2)
178: {
180: PetscMPIInt size,tag,i,j,*buf_s,*buf_r,*buf_j;
181: MPI_Request *s_waits,*r_waits;
182: MPI_Status *w_status;
185: MPI_Comm_size(comm,&size);
186: PetscCommGetNewTag(comm,&tag);
188: /* cannot use PetscMalloc5() because r_waits and s_waits must be contiquous for the call to MPI_Waitall() */
189: PetscMalloc4(nrecvs+nsends,MPI_Request,&r_waits,2*nrecvs,PetscMPIInt,&buf_r,2*nsends,PetscMPIInt,&buf_s,nrecvs+nsends,MPI_Status,&w_status);
190: s_waits = r_waits + nrecvs;
192: /* Post the Irecv to get the message length-info */
193: PetscMalloc((nrecvs+1)*sizeof(PetscMPIInt),olengths1);
194: PetscMalloc((nrecvs+1)*sizeof(PetscMPIInt),olengths2);
195: for (i=0; i<nrecvs; i++) {
196: buf_j = buf_r + (2*i);
197: MPI_Irecv(buf_j,2,MPI_INT,MPI_ANY_SOURCE,tag,comm,r_waits+i);
198: }
200: /* Post the Isends with the message length-info */
201: for (i=0,j=0; i<size; ++i) {
202: if (ilengths1[i]) {
203: buf_j = buf_s + (2*j);
204: buf_j[0] = *(ilengths1+i);
205: buf_j[1] = *(ilengths2+i);
206: MPI_Isend(buf_j,2,MPI_INT,i,tag,comm,s_waits+j);
207: j++;
208: }
209: }
210:
211: /* Post waits on sends and receivs */
212: MPI_Waitall(nrecvs+nsends,r_waits,w_status);
214:
215: /* Pack up the received data */
216: PetscMalloc((nrecvs+1)*sizeof(PetscMPIInt),onodes);
217: for (i=0; i<nrecvs; ++i) {
218: (*onodes)[i] = w_status[i].MPI_SOURCE;
219: buf_j = buf_r + (2*i);
220: (*olengths1)[i] = buf_j[0];
221: (*olengths2)[i] = buf_j[1];
222: }
224: PetscFree4(r_waits,buf_r,buf_s,w_status);
225: return(0);
226: }
228: /*
230: Allocate a bufffer sufficient to hold messages of size specified in olengths.
231: And post Irecvs on these buffers using node info from onodes
232:
233: */
236: PetscErrorCode PetscPostIrecvInt(MPI_Comm comm,PetscMPIInt tag,PetscMPIInt nrecvs,PetscMPIInt *onodes,PetscMPIInt *olengths,PetscInt ***rbuf,MPI_Request **r_waits)
237: {
239: PetscInt len=0,**rbuf_t,i;
240: MPI_Request *r_waits_t;
244: /* compute memory required for recv buffers */
245: for (i=0; i<nrecvs; i++) len += olengths[i]; /* each message length */
246: len *= sizeof(PetscInt);
247: len += (nrecvs+1)*sizeof(PetscInt*); /* Array of pointers for each message */
249: /* allocate memory for recv buffers */
250: PetscMalloc(len,&rbuf_t);
251: rbuf_t[0] = (PetscInt*)(rbuf_t + nrecvs);
252: for (i=1; i<nrecvs; ++i) rbuf_t[i] = rbuf_t[i-1] + olengths[i-1];
254: /* Post the receives */
255: PetscMalloc((nrecvs+1)*sizeof(MPI_Request),&r_waits_t);
256: for (i=0; i<nrecvs; ++i) {
257: MPI_Irecv(rbuf_t[i],olengths[i],MPIU_INT,onodes[i],tag,comm,r_waits_t+i);
258: }
260: *rbuf = rbuf_t;
261: *r_waits = r_waits_t;
262: return(0);
263: }
267: PetscErrorCode PetscPostIrecvScalar(MPI_Comm comm,PetscMPIInt tag,PetscMPIInt nrecvs,PetscMPIInt *onodes,PetscMPIInt *olengths,PetscScalar ***rbuf,MPI_Request **r_waits)
268: {
270: PetscMPIInt len=0,i;
271: PetscScalar **rbuf_t;
272: MPI_Request *r_waits_t;
276: /* compute memory required for recv buffers */
277: for (i=0; i<nrecvs; i++) len += olengths[i]; /* each message length */
278: len *= sizeof(PetscScalar);
279: len += (nrecvs+1)*sizeof(PetscScalar*); /* Array of pointers for each message */
282: /* allocate memory for recv buffers */
283: PetscMalloc(len,&rbuf_t);
284: rbuf_t[0] = (PetscScalar*)(rbuf_t + nrecvs);
285: for (i=1; i<nrecvs; ++i) rbuf_t[i] = rbuf_t[i-1] + olengths[i-1];
287: /* Post the receives */
288: PetscMalloc((nrecvs+1)*sizeof(MPI_Request),&r_waits_t);
289: for (i=0; i<nrecvs; ++i) {
290: MPI_Irecv(rbuf_t[i],olengths[i],MPIU_SCALAR,onodes[i],tag,comm,r_waits_t+i);
291: }
293: *rbuf = rbuf_t;
294: *r_waits = r_waits_t;
295: return(0);
296: }