Actual source code: mpimesg.c
1: /*$Id: mpimesg.c,v 1.14 2001/08/07 03:02:06 balay Exp $*/
3: #include petsc.h
6: /*@C
7: PetscGatherNumberOfMessages - Computes the number of messages a node expects to receive
9: Collective on MPI_Comm
11: Input Parameters:
12: + comm - Communicator
13: . iflags - an array of integers of length sizeof(comm). A '1' in ilengths[i] represent a
14: message from current node to ith node. Optionally PETSC_NULL
15: - ilengths - Non zero ilengths[i] represent a message to i of length ilengths[i].
16: Optionally PETSC_NULL.
18: Output Parameters:
19: . nrecvs - number of messages received
21: Level: developer
23: Concepts: mpi utility
25: Notes:
26: With this info, the correct message lengths can be determined using
27: PetscGatherMessageLengths()
29: Either iflags or ilengths should be provided. If iflags is not
30: provided (PETSC_NULL) it can be computed from ilengths. If iflags is
31: provided, ilengths is not required.
33: .seealso: PetscGatherMessageLengths()
34: @*/
35: int PetscGatherNumberOfMessages(MPI_Comm comm,int *iflags,int *ilengths,int *nrecvs)
36: {
37: int *recv_buf,size,rank,i,ierr,*iflags_local;
41: MPI_Comm_size(comm,&size);
42: MPI_Comm_rank(comm,&rank);
45: /* If iflags not provided, compute iflags from ilengths */
46: if (!iflags) {
47: if (!ilengths) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Either iflags or ilengths should be provided");
48: PetscMalloc(size*sizeof(int),&iflags_local);
49: for (i=0; i<size; i++) {
50: if (ilengths[i]) iflags_local[i] = 1;
51: else iflags_local[i] = 0;
52: }
53: } else {
54: iflags_local = iflags;
55: }
57: PetscMalloc(size*sizeof(int),&recv_buf);
59: /* Now post an allreduce to determine the numer of messages the current node will receive */
60: ierr = MPI_Allreduce(iflags_local,recv_buf,size,MPI_INT,MPI_SUM,comm);
61: *nrecvs = recv_buf[rank];
63: if (!iflags) {
64: PetscFree(iflags_local);
65: }
66: PetscFree(recv_buf);
67:
68: return(0);
69: }
72: /*@C
73: PetscGatherMessageLengths - Computes info about messages that a MPI-node will receive,
74: including (from-id,length) pairs for each message.
76: Collective on MPI_Comm
78: Input Parameters:
79: + comm - Communicator
80: . nsends - number of messages that are to be sent.
81: . nrecvs - number of messages being received
82: - ilengths - an array of integers of length sizeof(comm)
83: a non zero ilengths[i] represent a message to i of length ilengths[i]
86: Output Parameters:
87: + onodes - list of node-ids from which messages are expected
88: - olengths - corresponding message lengths
90: Level: developer
92: Concepts: mpi utility
94: Notes:
95: With this info, the correct MPI_Irecv() can be posted with the correct
96: from-id, with a buffer with the right amount of memory required.
98: The calling function deallocates the memory in onodes and olengths
100: To determine nrecevs, one can use PetscGatherNumberOfMessages()
102: .seealso: PetscGatherNumberOfMessages()
103: @*/
104: int PetscGatherMessageLengths(MPI_Comm comm,int nsends,int nrecvs,int *ilengths,int **onodes,int **olengths)
105: {
106: int size,i,j,tag,ierr;
107: MPI_Request *s_waits,*r_waits;
108: MPI_Status *w_status;
112: MPI_Comm_size(comm,&size);
113: PetscCommGetNewTag(comm,&tag);
115: PetscMalloc((nrecvs+nsends+1)*sizeof(MPI_Request),&r_waits);
116: s_waits = r_waits + nrecvs;
118: /* Now post the Irecv to get the message length-info */
119: PetscMalloc((nrecvs+1)*sizeof(int),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: /* Now post the Isends with the message lenght-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: }
131:
132: /* Now post waits on sends and receivs */
133: PetscMalloc((nrecvs+nsends+1)*sizeof(MPI_Status),&w_status);
134: MPI_Waitall(nrecvs+nsends,r_waits,w_status);
136:
137: /* Now pack up the received data */
138: PetscMalloc((nrecvs+1)*sizeof(int),onodes);
139: for (i=0; i<nrecvs; ++i) {
140: (*onodes)[i] = w_status[i].MPI_SOURCE;
141: }
143: PetscFree(r_waits);
144: PetscFree(w_status);
145:
146: return(0);
147: }
149: /*
151: Allocate a bufffer sufficient to hold messages of size specified in olengths.
152: And post Irecvs on these buffers using node info from onodes
153:
154: */
155: int PetscPostIrecvInt(MPI_Comm comm,int tag,int nrecvs,int *onodes,int *olengths,int ***rbuf,MPI_Request **r_waits)
156: {
157: int len=0,**rbuf_t,i,ierr;
158: MPI_Request *r_waits_t;
162: /* compute memory required for recv buffers */
163: for (i=0; i<nrecvs; i++) len += olengths[i]; /* each message length */
164: len *= sizeof(int);
165: len += (nrecvs+1)*sizeof(int*); /* Array of pointers for each message */
167: /* allocate memory for recv buffers */
168: ierr = PetscMalloc(len,&rbuf_t);
169: rbuf_t[0] = (int*)(rbuf_t + nrecvs);
170: for (i=1; i<nrecvs; ++i) rbuf_t[i] = rbuf_t[i-1] + olengths[i-1];
172: /* Post the receives */
173: PetscMalloc((nrecvs+1)*sizeof(MPI_Request),&r_waits_t);
174: for (i=0; i<nrecvs; ++i) {
175: MPI_Irecv(rbuf_t[i],olengths[i],MPI_INT,onodes[i],tag,comm,r_waits_t+i);
176: }
178: *rbuf = rbuf_t;
179: *r_waits = r_waits_t;
180: return(0);
181: }
183: int PetscPostIrecvScalar(MPI_Comm comm,int tag,int nrecvs,int *onodes,int *olengths,PetscScalar ***rbuf,MPI_Request **r_waits)
184: {
185: int len=0,i,ierr;
186: PetscScalar **rbuf_t;
187: MPI_Request *r_waits_t;
191: /* compute memory required for recv buffers */
192: for (i=0; i<nrecvs; i++) len += olengths[i]; /* each message length */
193: len *= sizeof(PetscScalar);
194: len += (nrecvs+1)*sizeof(PetscScalar*); /* Array of pointers for each message */
197: /* allocate memory for recv buffers */
198: ierr = PetscMalloc(len,&rbuf_t);
199: rbuf_t[0] = (PetscScalar*)(rbuf_t + nrecvs);
200: for (i=1; i<nrecvs; ++i) rbuf_t[i] = rbuf_t[i-1] + olengths[i-1];
202: /* Post the receives */
203: PetscMalloc((nrecvs+1)*sizeof(MPI_Request),&r_waits_t);
204: for (i=0; i<nrecvs; ++i) {
205: MPI_Irecv(rbuf_t[i],olengths[i],MPIU_SCALAR,onodes[i],tag,comm,r_waits_t+i);
206: }
208: *rbuf = rbuf_t;
209: *r_waits = r_waits_t;
210: return(0);
211: }