source: CONFIG_DEVT/IPSLCM6.5_work_ENSEMBLES/oasis3-mct/lib/mct/mpi-serial/send.c

Last change on this file was 5725, checked in by aclsce, 3 years ago

Added new oasis3-MCT version to be used to handle ensembles simulations with XIOS.

File size: 4.9 KB
Line 
1
2#include "mpiP.h"
3
4
5/*
6 * SENDING
7 *
8 */
9
10
11
12static int mpi_match_recv(void *r, void *tag)
13{
14  return( ((Req *)r)->tag == MPI_ANY_TAG ||
15          ((Req *)r)->tag == *((int *)tag) );
16}
17
18
19/*
20 *
21 */
22
23
24
25FC_FUNC( mpi_isend , MPI_ISEND )(void *buf, int *count, int *datatype,
26   int *dest, int *tag, int *comm, int *req, int *ierror)
27{
28
29  *ierror=MPI_Isend(buf,*count,*datatype,*dest,*tag,
30                    *comm, (void *)req);
31
32}
33
34
35
36int MPI_Isend(void *buf, int count, MPI_Datatype datatype,
37              int dest, int tag, MPI_Comm comm, MPI_Request *request)
38{
39  pListitem match;
40  Comm *mycomm;
41  Req *rreq, *sreq;
42
43  mycomm=mpi_handle_to_ptr(comm);   /* (Comm *)comm; */
44
45#ifdef INFO
46  fflush(stdout);
47  fprintf(stderr,"MPI_Isend: Comm=%d  tag=%d  count=%d type=%d\n",
48         mycomm->num,tag,count,datatype);
49#endif
50
51  if (dest!=0 && dest!=MPI_PROC_NULL)
52    {
53      fprintf(stderr,"MPI_Isend: send to %d\n",dest);
54      abort();
55    }
56
57  mpi_alloc_handle(request,(void **) &sreq);
58
59
60  if (dest==MPI_PROC_NULL)
61    {
62      sreq->complete=1;
63      return(MPI_SUCCESS);
64    }
65
66  if ( match=AP_list_search_func(mycomm->recvlist,mpi_match_recv,&tag) )
67    {
68      rreq=(Req *)AP_listitem_data(match);
69      AP_list_delete_item(mycomm->recvlist,match);
70
71//      memcpy(rreq->buf,buf,count * datatype);
72      copy_data2(buf, count, datatype, rreq->buf, count, datatype);
73      rreq->complete=1;
74      rreq->source=0;
75      rreq->tag=tag;                    /* in case rreq->tag was MPI_ANY_TAG */
76
77      sreq->complete=1;
78
79#ifdef DEBUG
80      printf("Completion(send) value=%d tag=%d\n",
81             *((int *)buf),rreq->tag);
82#endif
83
84      return(MPI_SUCCESS);
85    }
86
87  sreq->buf=buf;
88  sreq->tag=tag;
89  sreq->complete=0;
90  sreq->listitem=AP_list_append(mycomm->sendlist,sreq);
91
92#ifdef INFO
93  print_list(mycomm->sendlist,"sendlist for comm ",mycomm->num);
94#endif
95
96  return(MPI_SUCCESS);
97}
98
99
100/*********/
101
102
103FC_FUNC(mpi_send, MPI_SEND) ( void *buf, int *count, int *datatype,
104                                int *dest, int *tag, int *comm, int *ierror)
105{
106  *ierror=MPI_Send(buf, *count, *datatype, *dest, *tag, *comm);
107}
108
109
110
111int MPI_Send(void* buf, int count, MPI_Datatype datatype,
112             int dest, int tag, MPI_Comm comm)
113{
114  MPI_Request request;
115  MPI_Status status;
116
117#ifdef INFO
118  fflush(stdout);
119  fprintf(stderr,"MPI_Send: ");
120#endif
121
122  MPI_Isend(buf,count,datatype,dest,tag,comm,&request);
123  MPI_Wait(&request,&status);
124
125
126  return(MPI_SUCCESS);
127}
128
129
130
131
132/*********/
133
134
135FC_FUNC(mpi_ssend, MPI_SSEND) ( void *buf, int *count, int *datatype,
136                                  int *dest, int *tag, int *comm, int *ierror)
137{
138  *ierror=MPI_Send(buf, *count, *datatype, *dest, *tag, *comm);
139}
140
141
142
143int MPI_Ssend(void* buf, int count, MPI_Datatype datatype,
144              int dest, int tag, MPI_Comm comm)
145{
146  return(MPI_Send(buf,count,datatype,dest,tag,comm));
147}
148
149
150
151/*********/
152
153
154FC_FUNC(mpi_rsend, MPI_RSEND) ( void *buf, int *count, int *datatype,
155                                  int *dest, int *tag, int *comm, int *ierror)
156{
157  *ierror=MPI_Send(buf, *count, *datatype, *dest, *tag, *comm);
158}
159
160
161
162int MPI_Rsend(void* buf, int count, MPI_Datatype datatype,
163              int dest, int tag, MPI_Comm comm)
164{
165  return(MPI_Send(buf,count,datatype,dest,tag,comm));
166}
167
168
169
170
171/*********/
172
173
174
175FC_FUNC( mpi_irsend , MPI_IRSEND )(void *buf, int *count, int *datatype,
176   int *dest, int *tag, int *comm, int *req, int *ierror)
177{
178
179  *ierror=MPI_Irsend(buf,*count,*datatype,*dest,*tag,
180                    *comm, (void *)req);
181
182}
183
184
185
186int MPI_Irsend(void *buf, int count, MPI_Datatype datatype,
187               int dest, int tag, MPI_Comm comm, MPI_Request *request)
188{
189  MPI_Status status;
190  Req *req;
191
192
193  MPI_Isend(buf,count,datatype,dest,tag,comm,request);
194
195  /* Ready mode implied a receive must already be posted,
196   * so the Isend should have completed already.
197   * Can't use MPI_Test here for the error check because
198   * it would clear the request prematurely.
199   */
200
201  req=mpi_handle_to_ptr(*request);
202  if ( !req->complete )
203    {
204      fprintf(stderr,"MPI_Irsend: no matching receive found\n");
205      abort();
206    }
207
208
209  return(MPI_SUCCESS);
210}
211
212
213
214
215/*********/
216
217
218FC_FUNC(mpi_sendrecv, MPI_SENDRECV) (
219     void *sendbuf, int *sendcount, int *sendtype, int *dest, int *sendtag,
220     void *recvbuf, int *recvcount, int *recvtype, int *source, int *recvtag,
221     int *comm, int *status,
222     int *ierror)
223{
224  *ierror=MPI_Sendrecv(sendbuf, *sendcount, *sendtype, *dest, *sendtag,
225                       recvbuf, *recvcount, *recvtype, *source, *recvtag,
226                       *comm, mpi_c_status(status));
227}
228
229
230
231int MPI_Sendrecv(void* sendbuf, int sendcount, MPI_Datatype sendtype,
232                 int dest, int sendtag,
233                 void *recvbuf, int recvcount, MPI_Datatype recvtype,
234                 int source, int recvtag,
235                 MPI_Comm comm, MPI_Status *status)
236{
237  MPI_Request request;
238
239
240  MPI_Irecv(recvbuf, recvcount, recvtype, source, recvtag, comm, &request);
241
242  MPI_Send(sendbuf, sendcount, sendtype, dest, sendtag, comm);
243
244  MPI_Wait(&request,status);
245
246
247  return(MPI_SUCCESS);
248}
249
250
251
Note: See TracBrowser for help on using the repository browser.