Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
Plug some easy memleaks
[simgrid.git] / src / smpi / smpi_base.c
1 #include "private.h"
2 #include "xbt/time.h"
3
4 XBT_LOG_NEW_DEFAULT_SUBCATEGORY(smpi_base, smpi,
5                                 "Logging specific to SMPI (base)");
6 XBT_LOG_EXTERNAL_CATEGORY(smpi_base);
7 XBT_LOG_EXTERNAL_CATEGORY(smpi_bench);
8 XBT_LOG_EXTERNAL_CATEGORY(smpi_kernel);
9 XBT_LOG_EXTERNAL_CATEGORY(smpi_mpi);
10 XBT_LOG_EXTERNAL_CATEGORY(smpi_receiver);
11 XBT_LOG_EXTERNAL_CATEGORY(smpi_sender);
12 XBT_LOG_EXTERNAL_CATEGORY(smpi_util);
13
14 smpi_mpi_global_t smpi_mpi_global = NULL;
15
16 /**
17  * Operations of MPI_OP : implemented=land,sum,min,max
18  **/
19 void smpi_mpi_land_func(void *a, void *b, int *length,
20                         MPI_Datatype * datatype);
21
22 void smpi_mpi_land_func(void *a, void *b, int *length,
23                         MPI_Datatype * datatype)
24 {
25   int i;
26   if (*datatype == smpi_mpi_global->mpi_int) {
27     int *x = a, *y = b;
28     for (i = 0; i < *length; i++) {
29       y[i] = x[i] && y[i];
30     }
31   }
32 }
33
34 /**
35  * sum two vectors element-wise
36  *
37  * @param a the first vectors
38  * @param b the second vectors
39  * @return the second vector is modified and contains the element-wise sums
40  **/
41 void smpi_mpi_sum_func(void *a, void *b, int *length,
42                        MPI_Datatype * datatype);
43
44 void smpi_mpi_sum_func(void *a, void *b, int *length, MPI_Datatype * datatype)
45 {
46           int i;
47           if (*datatype == smpi_mpi_global->mpi_byte) {
48                                 char *x = a, *y = b;
49                                 for (i = 0; i < *length; i++) {
50                                           y[i] = x[i] + y[i];
51                                 }
52           } else if (*datatype == smpi_mpi_global->mpi_int) {
53                                 int *x = a, *y = b;
54                                 for (i = 0; i < *length; i++) {
55                                           y[i] = x[i] + y[i];
56                                 }
57           } else if (*datatype == smpi_mpi_global->mpi_float) {
58                                 float *x = a, *y = b;
59                                 for (i = 0; i < *length; i++) {
60                                           y[i] = x[i] + y[i];
61                                 }
62           } else if (*datatype == smpi_mpi_global->mpi_double) {
63                                 double *x = a, *y = b;
64                                 for (i = 0; i < *length; i++) {
65                                           y[i] = x[i] + y[i];
66                                 }
67           }
68 }
69 /**
70  * compute the min of two vectors element-wise
71  **/
72 void smpi_mpi_min_func(void *a, void *b, int *length, MPI_Datatype * datatype);
73
74 void smpi_mpi_min_func(void *a, void *b, int *length, MPI_Datatype * datatype)
75 {
76           int i;
77           if (*datatype == smpi_mpi_global->mpi_byte) {
78                                 char *x = a, *y = b;
79                                 for (i = 0; i < *length; i++) {
80                                           y[i] = x[i] < y[i] ? x[i] : y[i];
81                                 }
82           } else {
83           if (*datatype == smpi_mpi_global->mpi_int) {
84                                 int *x = a, *y = b;
85                                 for (i = 0; i < *length; i++) {
86                                           y[i] = x[i] < y[i] ? x[i] : y[i];
87                                 }
88           } else {
89           if (*datatype == smpi_mpi_global->mpi_float) {
90                                 float *x = a, *y = b;
91                                 for (i = 0; i < *length; i++) {
92                                           y[i] = x[i] < y[i] ? x[i] : y[i];
93                                 }
94           } else {
95           if (*datatype == smpi_mpi_global->mpi_double) {
96                                 double *x = a, *y = b;
97                                 for (i = 0; i < *length; i++) {
98                                           y[i] = x[i] < y[i] ? x[i] : y[i];
99                                 }
100
101           }}}}
102 }
103 /**
104  * compute the max of two vectors element-wise
105  **/
106 void smpi_mpi_max_func(void *a, void *b, int *length, MPI_Datatype * datatype);
107
108 void smpi_mpi_max_func(void *a, void *b, int *length, MPI_Datatype * datatype)
109 {
110           int i;
111           if (*datatype == smpi_mpi_global->mpi_byte) {
112                                 char *x = a, *y = b;
113                                 for (i = 0; i < *length; i++) {
114                                           y[i] = x[i] > y[i] ? x[i] : y[i];
115                                 }
116           } else if (*datatype == smpi_mpi_global->mpi_int) {
117                                 int *x = a, *y = b;
118                                 for (i = 0; i < *length; i++) {
119                                           y[i] = x[i] > y[i] ? x[i] : y[i];
120                                 }
121           } else if (*datatype == smpi_mpi_global->mpi_float) {
122                                 float *x = a, *y = b;
123                                 for (i = 0; i < *length; i++) {
124                                           y[i] = x[i] > y[i] ? x[i] : y[i];
125                                 }
126           } else if (*datatype == smpi_mpi_global->mpi_double) {
127                                 double *x = a, *y = b;
128                                 for (i = 0; i < *length; i++) {
129                                           y[i] = x[i] > y[i] ? x[i] : y[i];
130                                 }
131
132           }
133 }
134
135
136
137
138 /**
139  * tell the MPI rank of the calling process (from its SIMIX process id)
140  **/
141 int smpi_mpi_comm_rank(smpi_mpi_communicator_t comm)
142 {
143   return comm->index_to_rank_map[smpi_process_index()];
144 }
145
146 void smpi_process_init(int *argc, char***argv)
147 {
148   smpi_process_data_t pdata;
149
150   // initialize some local variables
151
152   pdata = xbt_new(s_smpi_process_data_t, 1);
153   SIMIX_process_set_data(SIMIX_process_self(),pdata);
154
155   /* get rank from command line, and remove it from argv */
156   pdata->index = atoi( (*argv)[1] );
157   DEBUG1("I'm rank %d",pdata->index);
158   if (*argc>2) {
159           memmove((*argv)[1],(*argv)[2], sizeof(char*)* (*argc-2));
160           (*argv)[ (*argc)-1] = NULL;
161   }
162   (*argc)--;
163
164   pdata->mutex = SIMIX_mutex_init();
165   pdata->cond = SIMIX_cond_init();
166   pdata->finalize = 0;
167
168   pdata->pending_recv_request_queue = xbt_fifo_new();
169   pdata->pending_send_request_queue = xbt_fifo_new();
170   pdata->received_message_queue = xbt_fifo_new();
171
172   pdata->main = SIMIX_process_self();
173   pdata->sender = SIMIX_process_create("smpi_sender",
174           smpi_sender, pdata,
175           SIMIX_host_get_name(SIMIX_host_self()), 0, NULL,
176           /*props */ NULL);
177   pdata->receiver = SIMIX_process_create("smpi_receiver",
178           smpi_receiver, pdata,
179           SIMIX_host_get_name(SIMIX_host_self()), 0, NULL,
180           /*props */ NULL);
181
182   smpi_global->main_processes[pdata->index] = SIMIX_process_self();
183   return;
184 }
185
186 void smpi_process_finalize()
187 {
188   smpi_process_data_t pdata =  SIMIX_process_get_data(SIMIX_process_self());
189
190   pdata->finalize = 2; /* Tell sender and receiver to quit */
191   SIMIX_process_resume(pdata->sender);
192   SIMIX_process_resume(pdata->receiver);
193   while (pdata->finalize>0) { /* wait until it's done */
194           SIMIX_cond_wait(pdata->cond,pdata->mutex);
195   }
196
197   SIMIX_mutex_destroy(pdata->mutex);
198   SIMIX_cond_destroy(pdata->cond);
199   xbt_fifo_free(pdata->pending_recv_request_queue);
200   xbt_fifo_free(pdata->pending_send_request_queue);
201   xbt_fifo_free(pdata->received_message_queue);
202   xbt_free(pdata);
203 }
204
205 int smpi_mpi_barrier(smpi_mpi_communicator_t comm)
206 {
207
208   SIMIX_mutex_lock(comm->barrier_mutex);
209   ++comm->barrier_count;
210   if (comm->barrier_count > comm->size) {       // only happens on second barrier...
211     comm->barrier_count = 0;
212   } else if (comm->barrier_count == comm->size) {
213     SIMIX_cond_broadcast(comm->barrier_cond);
214   }
215   while (comm->barrier_count < comm->size) {
216     SIMIX_cond_wait(comm->barrier_cond, comm->barrier_mutex);
217   }
218   SIMIX_mutex_unlock(comm->barrier_mutex);
219
220   return MPI_SUCCESS;
221 }
222
223 int smpi_mpi_isend(smpi_mpi_request_t request)
224 {
225         smpi_process_data_t pdata =  SIMIX_process_get_data(SIMIX_process_self());
226   int retval = MPI_SUCCESS;
227
228   if (NULL == request) {
229     retval = MPI_ERR_INTERN;
230   } else {
231     xbt_fifo_push(pdata->pending_send_request_queue, request);
232     SIMIX_process_resume(pdata->sender);
233   }
234
235   return retval;
236 }
237
238 int smpi_mpi_irecv(smpi_mpi_request_t request)
239 {
240   int retval = MPI_SUCCESS;
241   smpi_process_data_t pdata =  SIMIX_process_get_data(SIMIX_process_self());
242
243   if (NULL == request) {
244     retval = MPI_ERR_INTERN;
245   } else {
246     xbt_fifo_push(pdata->pending_recv_request_queue, request);
247
248     if (SIMIX_process_is_suspended(pdata->receiver)) {
249       SIMIX_process_resume(pdata->receiver);
250     }
251   }
252
253   return retval;
254 }
255
256 int smpi_mpi_wait(smpi_mpi_request_t request, smpi_mpi_status_t * status)
257 {
258   int retval = MPI_SUCCESS;
259
260   if (NULL == request) {
261     retval = MPI_ERR_INTERN;
262   } else {
263     SIMIX_mutex_lock(request->mutex);
264     while (!request->completed) {
265       SIMIX_cond_wait(request->cond, request->mutex);
266     }
267     if (NULL != status) {
268       status->MPI_SOURCE = request->src;
269       status->MPI_TAG = request->tag;
270       status->MPI_ERROR = MPI_SUCCESS;
271     }
272     SIMIX_mutex_unlock(request->mutex);
273   }
274
275   return retval;
276 }
277
278 int smpi_mpi_waitall(int count, smpi_mpi_request_t requests[], smpi_mpi_status_t status[]) {
279         int cpt;
280         int index;
281         int retval;
282         smpi_mpi_status_t stat;
283
284         for (cpt=0; cpt<count;cpt++) {
285                 retval = smpi_mpi_waitany(count,requests, &index,&stat);
286                 if (retval != MPI_SUCCESS)
287                         return retval;
288                 memcpy(&(status[index]),&stat,sizeof(stat));
289         }
290         return MPI_SUCCESS;
291 }
292
293 int smpi_mpi_waitany(int count, smpi_mpi_request_t *requests, int *index, smpi_mpi_status_t *status) {
294           int cpt;
295
296           *index = MPI_UNDEFINED;
297           if (NULL == requests) {
298             return MPI_ERR_INTERN;
299           }
300           /* First check if one of them is already done */
301           for (cpt=0;cpt<count;cpt++) {
302                   if (requests[cpt]->completed && !requests[cpt]->consumed) { /* got ya */
303                           *index=cpt;
304                           goto found_request;
305                   }
306           }
307           /* If none found, block */
308           /* FIXME: should use a SIMIX_cond_waitany, when implemented. For now, block on the first one */
309           while (1) {
310                   for (cpt=0;cpt<count;cpt++) {
311                           if (!requests[cpt]->completed) { /* this one is not done, wait on it */
312                                   while (!requests[cpt]->completed)
313                                       SIMIX_cond_wait(requests[cpt]->cond, requests[cpt]->mutex);
314
315                                   *index=cpt;
316                                   goto found_request;
317                           }
318                   }
319                   if (cpt == count) /* they are all done. Damn user */
320                           return MPI_ERR_REQUEST;
321           }
322
323           found_request:
324           requests[*index]->consumed = 1;
325
326            if (NULL != status) {
327               status->MPI_SOURCE = requests[*index]->src;
328               status->MPI_TAG = requests[*index]->tag;
329               status->MPI_ERROR = MPI_SUCCESS;
330             }
331           return MPI_SUCCESS;
332
333 }