Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
Give LD_LIBRARY_PATH to smpirun.
[simgrid.git] / src / smpi / smpi_mpi.c
1 /* Copyright (c) 2007, 2008, 2009, 2010. The SimGrid Team.
2  * All rights reserved.                                                     */
3
4 /* This program is free software; you can redistribute it and/or modify it
5   * under the terms of the license (GNU LGPL) which comes with this package. */
6
7 #include "private.h"
8 #include "smpi_coll_private.h"
9 #include "smpi_mpi_dt_private.h"
10
11 XBT_LOG_NEW_DEFAULT_SUBCATEGORY(smpi_mpi, smpi,
12                                 "Logging specific to SMPI (mpi)");
13
14 /* MPI User level calls */
15
16 int MPI_Init(int* argc, char*** argv) {
17   smpi_process_init(argc, argv);
18   smpi_bench_begin(-1, NULL);
19   return MPI_SUCCESS;
20 }
21
22 int MPI_Finalize(void) {
23   smpi_bench_end(-1, NULL);
24   smpi_process_destroy();
25   return MPI_SUCCESS;
26 }
27
28 int MPI_Init_thread(int* argc, char*** argv, int required, int* provided) {
29   if(provided != NULL) {
30     *provided = MPI_THREAD_MULTIPLE;
31   }
32   return MPI_Init(argc, argv);
33 }
34
35 int MPI_Query_thread(int* provided) {
36   int retval;
37
38   smpi_bench_end(-1, NULL);
39   if(provided == NULL) {
40     retval = MPI_ERR_ARG;
41   } else {
42     *provided = MPI_THREAD_MULTIPLE;
43     retval = MPI_SUCCESS;
44   }
45   smpi_bench_begin(-1, NULL);
46   return retval;
47 }
48
49 int MPI_Is_thread_main(int* flag) {
50   int retval;
51
52   smpi_bench_end(-1, NULL);
53   if(flag == NULL) {
54     retval = MPI_ERR_ARG;
55   } else {
56     *flag = smpi_process_index() == 0;
57     retval = MPI_SUCCESS;
58   }
59   smpi_bench_begin(-1, NULL);
60   return retval;
61 }
62
63 int MPI_Abort(MPI_Comm comm, int errorcode) {
64   smpi_bench_end(-1, NULL);
65   smpi_process_destroy();
66   // FIXME: should kill all processes in comm instead
67   SIMIX_process_kill(SIMIX_process_self());
68   return MPI_SUCCESS;
69 }
70
71 double MPI_Wtime(void) {
72   double time;
73
74   smpi_bench_end(-1, NULL);
75   time = SIMIX_get_clock();
76   smpi_bench_begin(-1, NULL);
77   return time;
78 }
79
80 int MPI_Address(void *location, MPI_Aint *address) {
81   int retval;
82
83   smpi_bench_end(-1, NULL);
84   if(!address) {
85     retval = MPI_ERR_ARG;
86   } else {
87     *address = (MPI_Aint)location;
88   }
89   smpi_bench_begin(-1, NULL);
90   return retval;
91 }
92
93 int MPI_Type_free(MPI_Datatype* datatype) {
94   int retval;
95
96   smpi_bench_end(-1, NULL);
97   if(!datatype) {
98     retval = MPI_ERR_ARG;
99   } else {
100     // FIXME: always fail for now
101     retval = MPI_ERR_TYPE;
102   }
103   smpi_bench_begin(-1, NULL);
104   return retval;
105 }
106
107 int MPI_Type_size(MPI_Datatype datatype, size_t* size) {
108   int retval;
109
110   smpi_bench_end(-1, NULL);
111   if(datatype == MPI_DATATYPE_NULL) {
112     retval = MPI_ERR_TYPE;
113   } else if(size == NULL) {
114     retval = MPI_ERR_ARG;
115   } else {
116     *size = smpi_datatype_size(datatype);
117     retval = MPI_SUCCESS;
118   }
119   smpi_bench_begin(-1, NULL);
120   return retval;
121 }
122
123 int MPI_Type_get_extent(MPI_Datatype datatype, MPI_Aint* lb, MPI_Aint* extent) {
124   int retval;
125
126   smpi_bench_end(-1, NULL);
127   if(datatype == MPI_DATATYPE_NULL) {
128     retval = MPI_ERR_TYPE;
129   } else if(lb == NULL || extent == NULL) {
130     retval = MPI_ERR_ARG;
131   } else {
132     retval = smpi_datatype_extent(datatype, lb, extent);
133   }
134   smpi_bench_begin(-1, NULL);
135   return retval;
136 }
137
138 int MPI_Type_extent(MPI_Datatype datatype, MPI_Aint* extent) {
139   int retval;
140   MPI_Aint dummy;
141
142   smpi_bench_end(-1, NULL);
143   if(datatype == MPI_DATATYPE_NULL) {
144     retval = MPI_ERR_TYPE;
145   } else if(extent == NULL) {
146     retval = MPI_ERR_ARG;
147   } else {
148     retval = smpi_datatype_extent(datatype, &dummy, extent);
149   }
150   smpi_bench_begin(-1, NULL);
151   return retval;
152 }
153
154 int MPI_Type_lb(MPI_Datatype datatype, MPI_Aint* disp) {
155   int retval;
156
157   smpi_bench_end(-1, NULL);
158   if(datatype == MPI_DATATYPE_NULL) {
159     retval = MPI_ERR_TYPE;
160   } else if(disp == NULL) {
161     retval = MPI_ERR_ARG;
162   } else {
163     *disp = smpi_datatype_lb(datatype);
164     retval = MPI_SUCCESS;
165   }
166   smpi_bench_begin(-1, NULL);
167   return retval;
168 }
169
170 int MPI_Type_ub(MPI_Datatype datatype, MPI_Aint* disp) {
171   int retval;
172
173   smpi_bench_end(-1, NULL);
174   if(datatype == MPI_DATATYPE_NULL) {
175     retval = MPI_ERR_TYPE;
176   } else if(disp == NULL) {
177     retval = MPI_ERR_ARG;
178   } else {
179     *disp = smpi_datatype_ub(datatype);
180     retval = MPI_SUCCESS;
181   }
182   smpi_bench_begin(-1, NULL);
183   return retval;
184 }
185
186 int MPI_Op_create(MPI_User_function* function, int commute, MPI_Op* op) {
187   int retval;
188
189   smpi_bench_end(-1, NULL);
190   if(function == NULL || op == NULL) {
191     retval = MPI_ERR_ARG;
192   } else {
193     *op = smpi_op_new(function, commute);
194     retval = MPI_SUCCESS;
195   }
196   smpi_bench_begin(-1, NULL);
197   return retval;
198 }
199
200 int MPI_Op_free(MPI_Op* op) {
201   int retval;
202
203   smpi_bench_end(-1, NULL);
204   if(op == NULL) {
205     retval = MPI_ERR_ARG;
206   } else if(*op == MPI_OP_NULL) {
207     retval = MPI_ERR_OP;
208   } else {
209     smpi_op_destroy(*op);
210     *op = MPI_OP_NULL;
211     retval = MPI_SUCCESS;
212   }
213   smpi_bench_begin(-1, NULL);
214   return retval;
215 }
216
217 int MPI_Group_free(MPI_Group *group) {
218   int retval;
219
220   smpi_bench_end(-1, NULL);
221   if(group == NULL) {
222     retval = MPI_ERR_ARG;
223   } else {
224     smpi_group_destroy(*group);
225     *group = MPI_GROUP_NULL;
226     retval = MPI_SUCCESS;
227   }
228   smpi_bench_begin(-1, NULL);
229   return retval;
230 }
231
232 int MPI_Group_size(MPI_Group group, int* size) {
233   int retval;
234
235   smpi_bench_end(-1, NULL);
236   if(group == MPI_GROUP_NULL) {
237     retval = MPI_ERR_GROUP;
238   } else if(size == NULL) {
239     retval = MPI_ERR_ARG;
240   } else {
241     *size = smpi_group_size(group);
242     retval = MPI_SUCCESS;
243   }
244   smpi_bench_begin(-1, NULL);
245   return retval;
246 }
247
248 int MPI_Group_rank(MPI_Group group, int* rank) {
249   int retval;
250
251   smpi_bench_end(-1, NULL);
252   if(group == MPI_GROUP_NULL) {
253     retval = MPI_ERR_GROUP;
254   } else if(rank == NULL) {
255     retval = MPI_ERR_ARG;
256   } else {
257     *rank = smpi_group_rank(group, smpi_process_index());
258     retval = MPI_SUCCESS;
259   }
260   smpi_bench_begin(-1, NULL);
261   return retval;
262 }
263
264 int MPI_Group_translate_ranks (MPI_Group group1, int n, int* ranks1, MPI_Group group2, int* ranks2) {
265   int retval, i, index;
266
267   smpi_bench_end(-1, NULL);
268   if(group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
269     retval = MPI_ERR_GROUP;
270   } else {
271     for(i = 0; i < n; i++) {
272       index = smpi_group_index(group1, ranks1[i]);
273       ranks2[i] = smpi_group_rank(group2, index);
274     }
275     retval = MPI_SUCCESS;
276   }
277   smpi_bench_begin(-1, NULL);
278   return retval;
279 }
280
281 int MPI_Group_compare(MPI_Group group1, MPI_Group group2, int* result) {
282   int retval;
283
284   smpi_bench_end(-1, NULL);
285   if(group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
286     retval = MPI_ERR_GROUP;
287   } else if(result == NULL) {
288     retval = MPI_ERR_ARG;
289   } else {
290     *result = smpi_group_compare(group1, group2);
291     retval = MPI_SUCCESS;
292   }
293   smpi_bench_begin(-1, NULL);
294   return retval;
295 }
296
297 int MPI_Group_union(MPI_Group group1, MPI_Group group2, MPI_Group* newgroup) {
298   int retval, i, proc1, proc2, size, size2;
299
300   smpi_bench_end(-1, NULL);
301   if(group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
302     retval = MPI_ERR_GROUP;
303   } else if(newgroup == NULL) {
304     retval = MPI_ERR_ARG;
305   } else {
306     size = smpi_group_size(group1);
307     size2 = smpi_group_size(group2);
308     for(i = 0; i < size2; i++) {
309       proc2 = smpi_group_index(group2, i);
310       proc1 = smpi_group_rank(group1, proc2);
311       if(proc1 == MPI_UNDEFINED) {
312         size++;
313       }
314     }
315     if(size == 0) {
316       *newgroup = MPI_GROUP_EMPTY;
317     } else {
318       *newgroup = smpi_group_new(size);
319       size2 = smpi_group_size(group1);
320       for(i = 0; i < size2; i++) {
321         proc1 = smpi_group_index(group1, i);
322         smpi_group_set_mapping(*newgroup, proc1, i);
323       }
324       for(i = size2; i < size; i++) {
325         proc2 = smpi_group_index(group2, i - size2);
326         smpi_group_set_mapping(*newgroup, proc2, i);
327       }
328     }
329     smpi_group_use(*newgroup);
330     retval = MPI_SUCCESS;
331   }
332   smpi_bench_begin(-1, NULL);
333   return retval;
334 }
335
336 int MPI_Group_intersection(MPI_Group group1, MPI_Group group2, MPI_Group* newgroup) {
337    int retval, i, proc1, proc2, size, size2;
338
339   smpi_bench_end(-1, NULL);
340   if(group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
341     retval = MPI_ERR_GROUP;
342   } else if(newgroup == NULL) {
343     retval = MPI_ERR_ARG;
344   } else {
345     size = smpi_group_size(group1);
346     size2 = smpi_group_size(group2);
347     for(i = 0; i < size2; i++) {
348       proc2 = smpi_group_index(group2, i);
349       proc1 = smpi_group_rank(group1, proc2);
350       if(proc1 == MPI_UNDEFINED) {
351         size--;
352       }
353     }
354     if(size == 0) {
355       *newgroup = MPI_GROUP_EMPTY;
356     } else {
357       *newgroup = smpi_group_new(size);
358       size2 = smpi_group_size(group1);
359       for(i = 0; i < size2; i++) {
360         proc1 = smpi_group_index(group1, i);
361         proc2 = smpi_group_rank(group2, proc1);
362         if(proc2 != MPI_UNDEFINED) {
363           smpi_group_set_mapping(*newgroup, proc1, i);
364         }
365       }
366     }
367     smpi_group_use(*newgroup);
368     retval = MPI_SUCCESS;
369   }
370   smpi_bench_begin(-1, NULL);
371   return retval;
372 }
373
374 int MPI_Group_difference(MPI_Group group1, MPI_Group group2, MPI_Group* newgroup) {
375   int retval, i, proc1, proc2, size, size2;
376
377   smpi_bench_end(-1, NULL);
378   if(group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
379     retval = MPI_ERR_GROUP;
380   } else if(newgroup == NULL) {
381     retval = MPI_ERR_ARG;
382   } else {
383     size = size2 = smpi_group_size(group1);
384     for(i = 0; i < size2; i++) {
385       proc1 = smpi_group_index(group1, i);
386       proc2 = smpi_group_rank(group2, proc1);
387       if(proc2 != MPI_UNDEFINED) {
388         size--;
389       }
390     }
391     if(size == 0) {
392       *newgroup = MPI_GROUP_EMPTY;
393     } else {
394       *newgroup = smpi_group_new(size);
395       for(i = 0; i < size2; i++) {
396         proc1 = smpi_group_index(group1, i);
397         proc2 = smpi_group_rank(group2, proc1);
398         if(proc2 == MPI_UNDEFINED) {
399           smpi_group_set_mapping(*newgroup, proc1, i);
400         }
401       }
402     }
403     smpi_group_use(*newgroup);
404     retval = MPI_SUCCESS;
405   }
406   smpi_bench_begin(-1, NULL);
407   return retval;
408 }
409
410 int MPI_Group_incl(MPI_Group group, int n, int* ranks, MPI_Group* newgroup) {
411   int retval, i, index;
412
413   smpi_bench_end(-1, NULL);
414   if(group == MPI_GROUP_NULL) {
415     retval = MPI_ERR_GROUP;
416   } else if(newgroup == NULL) {
417     retval = MPI_ERR_ARG;
418   } else {
419     if(n == 0) {
420       *newgroup = MPI_GROUP_EMPTY;
421     } else if(n == smpi_group_size(group)) {
422       *newgroup = group;
423     } else {
424       *newgroup = smpi_group_new(n);
425       for(i = 0; i < n; i++) {
426         index = smpi_group_index(group, ranks[i]);
427         smpi_group_set_mapping(*newgroup, index, i);
428       }
429     }
430     smpi_group_use(*newgroup);
431     retval = MPI_SUCCESS;
432   }
433   smpi_bench_begin(-1, NULL);
434   return retval;
435 }
436
437 int MPI_Group_excl(MPI_Group group, int n, int* ranks, MPI_Group* newgroup) {
438   int retval, i, size, rank, index;
439
440   smpi_bench_end(-1, NULL);
441   if(group == MPI_GROUP_NULL) {
442     retval = MPI_ERR_GROUP;
443   } else if(newgroup == NULL) {
444     retval = MPI_ERR_ARG;
445   } else {
446     if(n == 0) {
447       *newgroup = group;
448     } else if(n == smpi_group_size(group)) {
449       *newgroup = MPI_GROUP_EMPTY;
450     } else {
451       size = smpi_group_size(group) - n;
452       *newgroup = smpi_group_new(size);
453       rank = 0;
454       while(rank < size) {
455         for(i = 0; i < n; i++) {
456           if(ranks[i] == rank) {
457             break;
458           }
459         }
460         if(i >= n) {
461           index = smpi_group_index(group, rank);
462           smpi_group_set_mapping(*newgroup, index, rank);
463           rank++;
464         }
465       }
466     }
467     smpi_group_use(*newgroup);
468     retval = MPI_SUCCESS;
469   }
470   smpi_bench_begin(-1, NULL);
471   return retval;
472 }
473
474 int MPI_Group_range_incl(MPI_Group group, int n, int ranges[][3], MPI_Group* newgroup) {
475   int retval, i, j, rank, size, index;
476
477   smpi_bench_end(-1, NULL);
478   if(group == MPI_GROUP_NULL) {
479     retval = MPI_ERR_GROUP;
480   } else if(newgroup == NULL) {
481     retval = MPI_ERR_ARG;
482   } else {
483     if(n == 0) {
484       *newgroup = MPI_GROUP_EMPTY;
485     } else {
486       size = 0;
487       for(i = 0; i < n; i++) {
488         for(rank = ranges[i][0]; /* First */
489             rank >= 0 && rank <= ranges[i][1]; /* Last */
490             rank += ranges[i][2] /* Stride */) {
491           size++;
492         }
493       }
494       if(size == smpi_group_size(group)) {
495         *newgroup = group;
496       } else {
497         *newgroup = smpi_group_new(size);
498         j = 0;
499         for(i = 0; i < n; i++) {
500           for(rank = ranges[i][0]; /* First */
501               rank >= 0 && rank <= ranges[i][1]; /* Last */
502               rank += ranges[i][2] /* Stride */) {
503             index = smpi_group_index(group, rank);
504             smpi_group_set_mapping(*newgroup, index, j);
505             j++;
506           }
507         }
508       }
509     }
510     smpi_group_use(*newgroup);
511     retval = MPI_SUCCESS;
512   }
513   smpi_bench_begin(-1, NULL);
514   return retval;
515 }
516
517 int MPI_Group_range_excl(MPI_Group group, int n, int ranges[][3], MPI_Group* newgroup) {
518   int retval, i, newrank, rank, size, index, add;
519
520   smpi_bench_end(-1, NULL);
521   if(group == MPI_GROUP_NULL) {
522     retval = MPI_ERR_GROUP;
523   } else if(newgroup == NULL) {
524     retval = MPI_ERR_ARG;
525   } else {
526     if(n == 0) {
527       *newgroup = group;
528     } else {
529       size = smpi_group_size(group);
530       for(i = 0; i < n; i++) {
531         for(rank = ranges[i][0]; /* First */
532             rank >= 0 && rank <= ranges[i][1]; /* Last */
533             rank += ranges[i][2] /* Stride */) {
534           size--;
535         }
536       }
537       if(size == 0) {
538         *newgroup = MPI_GROUP_EMPTY;
539       } else {
540         *newgroup = smpi_group_new(size);
541         newrank = 0;
542         while(newrank < size) {
543           for(i = 0; i < n; i++) {
544             add = 1;
545             for(rank = ranges[i][0]; /* First */
546                 rank >= 0 && rank <= ranges[i][1]; /* Last */
547                 rank += ranges[i][2] /* Stride */) {
548               if(rank == newrank) {
549                 add = 0;
550                 break;
551               }
552             }
553             if(add == 1) {
554               index = smpi_group_index(group, newrank);
555               smpi_group_set_mapping(*newgroup, index, newrank);
556             }
557           }
558         }
559       }
560     }
561     smpi_group_use(*newgroup);
562     retval = MPI_SUCCESS;
563   }
564   smpi_bench_begin(-1, NULL);
565   return retval;
566 }
567
568 int MPI_Comm_rank(MPI_Comm comm, int* rank) {
569   int retval;
570
571   smpi_bench_end(-1, NULL);
572   if(comm == MPI_COMM_NULL) {
573     retval = MPI_ERR_COMM;
574   } else {
575     *rank = smpi_comm_rank(comm);
576     retval = MPI_SUCCESS;
577   }
578   smpi_bench_begin(-1, NULL);
579   return retval;
580 }
581
582 int MPI_Comm_size(MPI_Comm comm, int* size) {
583   int retval;
584
585   smpi_bench_end(-1, NULL);
586   if(comm == MPI_COMM_NULL) {
587     retval = MPI_ERR_COMM;
588   } else if(size == NULL) {
589     retval = MPI_ERR_ARG;
590   } else {
591     *size = smpi_comm_size(comm);
592     retval = MPI_SUCCESS;
593   }
594   smpi_bench_begin(-1, NULL);
595   return retval;
596 }
597
598 int MPI_Comm_group(MPI_Comm comm, MPI_Group* group) {
599   int retval;
600
601   smpi_bench_end(-1, NULL);
602   if(comm == MPI_COMM_NULL) {
603     retval = MPI_ERR_COMM;
604   } else if(group == NULL) {
605     retval = MPI_ERR_ARG;
606   } else {
607     *group = smpi_comm_group(comm);
608     retval = MPI_SUCCESS;
609   }
610   smpi_bench_begin(-1, NULL);
611   return retval;
612 }
613
614 int MPI_Comm_compare(MPI_Comm comm1, MPI_Comm comm2, int* result) {
615   int retval;
616
617   smpi_bench_end(-1, NULL);
618   if(comm1 == MPI_COMM_NULL || comm2 == MPI_COMM_NULL) {
619     retval = MPI_ERR_COMM;
620   } else if(result == NULL) {
621     retval = MPI_ERR_ARG;
622   } else {
623     if(comm1 == comm2) { /* Same communicators means same groups */
624       *result = MPI_IDENT;
625     } else {
626       *result = smpi_group_compare(smpi_comm_group(comm1), smpi_comm_group(comm2));
627       if(*result == MPI_IDENT) {
628         *result = MPI_CONGRUENT;
629       }
630     }
631     retval = MPI_SUCCESS;
632   }
633   smpi_bench_begin(-1, NULL);
634   return retval;
635 }
636
637 int MPI_Comm_dup(MPI_Comm comm, MPI_Comm* newcomm) {
638   int retval;
639
640   smpi_bench_end(-1, NULL);
641   if(comm == MPI_COMM_NULL) {
642     retval = MPI_ERR_COMM;
643   } else if(newcomm == NULL) {
644     retval = MPI_ERR_ARG;
645   } else {
646     *newcomm = smpi_comm_new(smpi_comm_group(comm));
647     retval = MPI_SUCCESS;
648   }
649   smpi_bench_begin(-1, NULL);
650   return retval;
651 }
652
653 int MPI_Comm_create(MPI_Comm comm, MPI_Group group, MPI_Comm* newcomm) {
654   int retval;
655
656   smpi_bench_end(-1, NULL);
657   if(comm == MPI_COMM_NULL) {
658     retval = MPI_ERR_COMM;
659   } else if(group == MPI_GROUP_NULL) {
660     retval = MPI_ERR_GROUP;
661   } else if(newcomm == NULL) {
662     retval = MPI_ERR_ARG;
663   } else {
664     *newcomm = smpi_comm_new(group);
665     retval = MPI_SUCCESS;
666   }
667   smpi_bench_begin(-1, NULL);
668   return retval;
669 }
670
671 int MPI_Comm_free(MPI_Comm* comm) {
672   int retval;
673
674   smpi_bench_end(-1, NULL);
675   if(comm == NULL) {
676     retval = MPI_ERR_ARG;
677   } else if(*comm == MPI_COMM_NULL) {
678     retval = MPI_ERR_COMM;
679   } else {
680     smpi_comm_destroy(*comm);
681     *comm = MPI_COMM_NULL;
682     retval = MPI_SUCCESS;
683   }
684   smpi_bench_begin(-1, NULL);
685   return retval;
686 }
687
688 int MPI_Send_init(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request* request) {
689   int retval;
690   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
691
692   smpi_bench_end(rank, "Send_init");
693   if(request == NULL) {
694     retval = MPI_ERR_ARG;
695   } else if (comm == MPI_COMM_NULL) {
696     retval = MPI_ERR_COMM;
697   } else {
698     *request = smpi_mpi_send_init(buf, count, datatype, dst, tag, comm);
699     retval = MPI_SUCCESS;
700   }
701   smpi_bench_begin(rank, "Send_init");
702   return retval;
703 }
704
705 int MPI_Recv_init(void* buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Request* request) {
706   int retval;
707   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
708
709   smpi_bench_end(rank, "Recv_init");
710   if(request == NULL) {
711     retval = MPI_ERR_ARG;
712   } else if (comm == MPI_COMM_NULL) {
713     retval = MPI_ERR_COMM;
714   } else {
715     *request = smpi_mpi_recv_init(buf, count, datatype, src, tag, comm);
716     retval = MPI_SUCCESS;
717   }
718   smpi_bench_begin(rank, "Recv_init");
719   return retval;
720 }
721
722 int MPI_Start(MPI_Request* request) {
723   int retval;
724   MPI_Comm comm = (*request)->comm;
725   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
726
727   smpi_bench_end(rank, "Start");
728   if(request == NULL) {
729     retval = MPI_ERR_ARG;
730   } else {
731     smpi_mpi_start(*request);
732     retval = MPI_SUCCESS;
733   }
734   smpi_bench_begin(rank, "Start");
735   return retval;
736 }
737
738 int MPI_Startall(int count, MPI_Request* requests) {
739   int retval;
740   MPI_Comm comm = count > 0 && requests ? requests[0]->comm : MPI_COMM_NULL;
741   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
742
743   smpi_bench_end(rank, "Startall");
744   if(requests == NULL) {
745     retval = MPI_ERR_ARG;
746   } else {
747     smpi_mpi_startall(count, requests);
748     retval = MPI_SUCCESS;
749   }
750   smpi_bench_begin(rank, "Startall");
751   return retval;
752 }
753
754 int MPI_Request_free(MPI_Request* request) {
755   int retval;
756   MPI_Comm comm = (*request)->comm;
757   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
758
759   smpi_bench_end(rank, "Request_free");
760   if(request == NULL) {
761     retval = MPI_ERR_ARG;
762   } else {
763     smpi_mpi_request_free(request);
764     retval = MPI_SUCCESS;
765   }
766   smpi_bench_begin(rank, "Request_free");
767   return retval;
768 }
769
770 int MPI_Irecv(void* buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Request* request) {
771   int retval;
772   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
773
774   smpi_bench_end(rank, "Irecv");
775   if(request == NULL) {
776     retval = MPI_ERR_ARG;
777   } else if (comm == MPI_COMM_NULL) {
778     retval = MPI_ERR_COMM;
779   } else {
780     *request = smpi_mpi_irecv(buf, count, datatype, src, tag, comm);
781     retval = MPI_SUCCESS;
782   }
783   smpi_bench_begin(rank, "Irecv");
784   return retval;
785 }
786
787 int MPI_Isend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm, MPI_Request* request) {
788   int retval;
789   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
790
791   smpi_bench_end(rank, "Isend");
792   if(request == NULL) {
793     retval = MPI_ERR_ARG;
794   } else if (comm == MPI_COMM_NULL) {
795     retval = MPI_ERR_COMM;
796   } else {
797     *request = smpi_mpi_isend(buf, count, datatype, dst, tag, comm);
798     retval = MPI_SUCCESS;
799   }
800   smpi_bench_begin(rank, "Isend");
801   return retval;
802 }
803
804 int MPI_Recv(void* buf, int count, MPI_Datatype datatype, int src, int tag, MPI_Comm comm, MPI_Status* status) {
805   int retval;
806   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
807
808   smpi_bench_end(rank, "Recv");
809   if (comm == MPI_COMM_NULL) {
810     retval = MPI_ERR_COMM;
811   } else {
812     smpi_mpi_recv(buf, count, datatype, src, tag, comm, status);
813     retval = MPI_SUCCESS;
814   }
815   smpi_bench_begin(rank, "Recv");
816   return retval;
817 }
818
819 int MPI_Send(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm) {
820   int retval;
821   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
822
823   smpi_bench_end(rank, "Send");
824   if (comm == MPI_COMM_NULL) {
825     retval = MPI_ERR_COMM;
826   } else {
827     smpi_mpi_send(buf, count, datatype, dst, tag, comm);
828     retval = MPI_SUCCESS;
829   }
830   smpi_bench_begin(rank, "Send");
831   return retval;
832 }
833
834 int MPI_Sendrecv(void* sendbuf, int sendcount, MPI_Datatype sendtype, int dst, int sendtag, void* recvbuf, int recvcount, MPI_Datatype recvtype, int src, int recvtag, MPI_Comm comm, MPI_Status* status) {
835   int retval;
836   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
837
838   smpi_bench_end(rank, "Sendrecv");
839   if (comm == MPI_COMM_NULL) {
840     retval = MPI_ERR_COMM;
841   } else if (sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
842     retval = MPI_ERR_TYPE;
843   } else {
844     smpi_mpi_sendrecv(sendbuf, sendcount, sendtype, dst, sendtag, recvbuf, recvcount, recvtype, src, recvtag, comm, status);
845     retval = MPI_SUCCESS;
846   }
847   smpi_bench_begin(rank, "Sendrecv");
848   return retval;
849 }
850
851 int MPI_Sendrecv_replace(void* buf, int count, MPI_Datatype datatype, int dst, int sendtag, int src, int recvtag, MPI_Comm comm, MPI_Status* status) {
852   //TODO: suboptimal implementation
853   void* recvbuf;
854   int retval, size;
855
856   size = smpi_datatype_size(datatype) * count;
857   recvbuf = xbt_new(char, size);
858   retval = MPI_Sendrecv(buf, count, datatype, dst, sendtag, recvbuf, count, datatype, src, recvtag, comm, status);
859   memcpy(buf, recvbuf, size * sizeof(char));
860   xbt_free(recvbuf);
861   return retval;
862 }
863
864 int MPI_Test(MPI_Request* request, int* flag, MPI_Status* status) {
865   int retval;
866   int rank = request && (*request)->comm != MPI_COMM_NULL
867              ? smpi_comm_rank((*request)->comm)
868              : -1;
869
870   smpi_bench_end(rank, "Test");
871   if(request == NULL || flag == NULL) {
872     retval = MPI_ERR_ARG;
873   } else if(*request == MPI_REQUEST_NULL) {
874     retval = MPI_ERR_REQUEST;
875   } else {
876     *flag = smpi_mpi_test(request, status);
877     retval = MPI_SUCCESS;
878   }
879   smpi_bench_begin(rank, "Test");
880   return retval;
881 }
882
883 int MPI_Testany(int count, MPI_Request requests[], int* index, int* flag, MPI_Status* status) {
884   int retval;
885
886   smpi_bench_end(-1, NULL); //FIXME
887   if(index == NULL || flag == NULL) {
888     retval = MPI_ERR_ARG;
889   } else {
890     *flag = smpi_mpi_testany(count, requests, index, status);
891     retval = MPI_SUCCESS;
892   }
893   smpi_bench_begin(-1, NULL);
894   return retval;
895 }
896
897 int MPI_Wait(MPI_Request* request, MPI_Status* status) {
898   int retval;
899   int rank = request && (*request)->comm != MPI_COMM_NULL
900              ? smpi_comm_rank((*request)->comm)
901              : -1;
902
903   smpi_bench_end(rank, "Wait");
904   if(request == NULL) {
905     retval = MPI_ERR_ARG;
906   } else if(*request == MPI_REQUEST_NULL) {
907     retval = MPI_ERR_REQUEST;
908   } else {
909     smpi_mpi_wait(request, status);
910     retval = MPI_SUCCESS;
911   }
912   smpi_bench_begin(rank, "Wait");
913   return retval;
914 }
915
916 int MPI_Waitany(int count, MPI_Request requests[], int* index, MPI_Status* status) {
917   int retval;
918
919   smpi_bench_end(-1, NULL); //FIXME
920   if(index == NULL) {
921     retval = MPI_ERR_ARG;
922   } else {
923     *index = smpi_mpi_waitany(count, requests, status);
924     retval = MPI_SUCCESS;
925   }
926   smpi_bench_begin(-1, NULL);
927   return retval;
928 }
929
930 int MPI_Waitall(int count, MPI_Request requests[],  MPI_Status status[]) {
931   smpi_bench_end(-1, NULL); //FIXME
932   smpi_mpi_waitall(count, requests, status);
933   smpi_bench_begin(-1, NULL);
934   return MPI_SUCCESS;
935 }
936
937 int MPI_Waitsome(int incount, MPI_Request requests[], int* outcount, int* indices, MPI_Status status[]) {
938   int retval;
939
940   smpi_bench_end(-1, NULL); //FIXME
941   if(outcount == NULL || indices == NULL) {
942     retval = MPI_ERR_ARG;
943   } else {
944     *outcount = smpi_mpi_waitsome(incount, requests, indices, status);
945     retval = MPI_SUCCESS;
946   }
947   smpi_bench_begin(-1, NULL);
948   return retval;
949 }
950
951 int MPI_Bcast(void* buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm) {
952   int retval;
953   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
954
955   smpi_bench_end(rank, "Bcast");
956   if(comm == MPI_COMM_NULL) {
957     retval = MPI_ERR_COMM;
958   } else {
959     smpi_mpi_bcast(buf, count, datatype, root, comm);
960     retval = MPI_SUCCESS;
961   }
962   smpi_bench_begin(rank, "Bcast");
963   return retval;
964 }
965
966 int MPI_Barrier(MPI_Comm comm) {
967   int retval;
968   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
969
970   smpi_bench_end(rank, "Barrier");
971   if(comm == MPI_COMM_NULL) {
972     retval = MPI_ERR_COMM;
973   } else {
974     smpi_mpi_barrier(comm);
975     retval = MPI_SUCCESS;
976   }
977   smpi_bench_begin(rank, "Barrier");
978   return retval;
979 }
980
981 int MPI_Gather(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm) {
982   int retval;
983   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
984
985   smpi_bench_end(rank, "Gather");
986   if(comm == MPI_COMM_NULL) {
987     retval = MPI_ERR_COMM;
988   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
989     retval = MPI_ERR_TYPE;
990   } else {
991     smpi_mpi_gather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
992     retval = MPI_SUCCESS;
993   }
994   smpi_bench_begin(rank, "Gather");
995   return retval;
996 }
997
998 int MPI_Gatherv(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int* recvcounts, int* displs, MPI_Datatype recvtype, int root, MPI_Comm comm) {
999   int retval;
1000   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1001
1002   smpi_bench_end(rank, "Gatherv");
1003   if(comm == MPI_COMM_NULL) {
1004     retval = MPI_ERR_COMM;
1005   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
1006     retval = MPI_ERR_TYPE;
1007   } else if(recvcounts == NULL || displs == NULL) {
1008     retval = MPI_ERR_ARG;
1009   } else {
1010     smpi_mpi_gatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, root, comm);
1011     retval = MPI_SUCCESS;
1012   }
1013   smpi_bench_begin(rank, "Gatherv");
1014   return retval;
1015 }
1016
1017 int MPI_Allgather(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm) {
1018   int retval;
1019   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1020
1021   smpi_bench_end(rank, "Allgather");
1022   if(comm == MPI_COMM_NULL) {
1023     retval = MPI_ERR_COMM;
1024   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
1025     retval = MPI_ERR_TYPE;
1026   } else {
1027     smpi_mpi_allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
1028     retval = MPI_SUCCESS;
1029   }
1030   smpi_bench_begin(rank, "Allgather");
1031   return retval;
1032 }
1033
1034 int MPI_Allgatherv(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int* recvcounts, int* displs, MPI_Datatype recvtype, MPI_Comm comm) {
1035   int retval;
1036   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1037
1038   smpi_bench_end(rank, "Allgatherv");
1039   if(comm == MPI_COMM_NULL) {
1040     retval = MPI_ERR_COMM;
1041   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
1042     retval = MPI_ERR_TYPE;
1043   } else if(recvcounts == NULL || displs == NULL) {
1044     retval = MPI_ERR_ARG;
1045   } else {
1046     smpi_mpi_allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts, displs, recvtype, comm);
1047     retval = MPI_SUCCESS;
1048   }
1049   smpi_bench_begin(rank, "Allgatherv");
1050   return retval;
1051 }
1052
1053 int MPI_Scatter(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm) {
1054   int retval;
1055   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1056
1057   smpi_bench_end(rank, "Scatter");
1058   if(comm == MPI_COMM_NULL) {
1059     retval = MPI_ERR_COMM;
1060   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
1061     retval = MPI_ERR_TYPE;
1062   } else {
1063     smpi_mpi_scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, root, comm);
1064     retval = MPI_SUCCESS;
1065   }
1066   smpi_bench_begin(rank, "Scatter");
1067   return retval;
1068 }
1069
1070 int MPI_Scatterv(void* sendbuf, int* sendcounts, int* displs, MPI_Datatype sendtype, void* recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm) {
1071   int retval;
1072   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1073
1074   smpi_bench_end(rank, "Scatterv");
1075   if(comm == MPI_COMM_NULL) {
1076     retval = MPI_ERR_COMM;
1077   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
1078     retval = MPI_ERR_TYPE;
1079   } else if(sendcounts == NULL || displs == NULL) {
1080     retval = MPI_ERR_ARG;
1081   } else {
1082     smpi_mpi_scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf, recvcount, recvtype, root, comm);
1083     retval = MPI_SUCCESS;
1084   }
1085   smpi_bench_begin(rank, "Scatterv");
1086   return retval;
1087 }
1088
1089 int MPI_Reduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm) {
1090   int retval;
1091   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1092
1093   smpi_bench_end(rank, "Reduce");
1094   if(comm == MPI_COMM_NULL) {
1095     retval = MPI_ERR_COMM;
1096   } else if(datatype == MPI_DATATYPE_NULL || op == MPI_OP_NULL) {
1097     retval = MPI_ERR_ARG;
1098   } else {
1099     smpi_mpi_reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
1100     retval = MPI_SUCCESS;
1101   }
1102   smpi_bench_begin(rank, "Reduce");
1103   return retval;
1104 }
1105
1106 int MPI_Allreduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) {
1107   int retval;
1108   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1109
1110   smpi_bench_end(rank, "Allreduce");
1111   if(comm == MPI_COMM_NULL) {
1112     retval = MPI_ERR_COMM;
1113   } else if(datatype == MPI_DATATYPE_NULL) {
1114     retval = MPI_ERR_TYPE;
1115   } else if(op == MPI_OP_NULL) {
1116     retval = MPI_ERR_OP;
1117   } else {
1118     smpi_mpi_allreduce(sendbuf, recvbuf, count, datatype, op, comm);
1119     retval = MPI_SUCCESS;
1120   }
1121   smpi_bench_begin(rank, "Allreduce");
1122   return retval;
1123 }
1124
1125 int MPI_Scan(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) {
1126   int retval;
1127   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1128
1129   smpi_bench_end(rank, "Scan");
1130   if(comm == MPI_COMM_NULL) {
1131     retval = MPI_ERR_COMM;
1132   } else if(datatype == MPI_DATATYPE_NULL) {
1133     retval = MPI_ERR_TYPE;
1134   } else if(op == MPI_OP_NULL) {
1135     retval = MPI_ERR_OP;
1136   } else {
1137     smpi_mpi_scan(sendbuf, recvbuf, count, datatype, op, comm);
1138     retval = MPI_SUCCESS;
1139   }
1140   smpi_bench_begin(rank, "Scan");
1141   return retval;
1142 }
1143
1144 int MPI_Reduce_scatter(void* sendbuf, void* recvbuf, int* recvcounts, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) {
1145   int retval, i, size, count;
1146   int* displs;
1147   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1148
1149   smpi_bench_end(rank, "Reduce_scatter");
1150   if(comm == MPI_COMM_NULL) {
1151     retval = MPI_ERR_COMM;
1152   } else if(datatype == MPI_DATATYPE_NULL) {
1153     retval = MPI_ERR_TYPE;
1154   } else if(op == MPI_OP_NULL) {
1155     retval = MPI_ERR_OP;
1156   } else if(recvcounts == NULL) {
1157     retval = MPI_ERR_ARG;
1158   } else {
1159     /* arbitrarily choose root as rank 0 */
1160     /* TODO: faster direct implementation ? */
1161     size = smpi_comm_size(comm);
1162     count = 0;
1163     displs = xbt_new(int, size);
1164     for(i = 0; i < size; i++) {
1165       count += recvcounts[i];
1166       displs[i] = 0;
1167     }
1168     smpi_mpi_reduce(sendbuf, recvbuf, count, datatype, op, 0, comm);
1169     smpi_mpi_scatterv(recvbuf, recvcounts, displs, datatype, recvbuf, recvcounts[rank], datatype, 0, comm);
1170     xbt_free(displs);
1171     retval = MPI_SUCCESS;
1172   }
1173   smpi_bench_begin(rank, "Reduce_scatter");
1174   return retval;
1175 }
1176
1177 int MPI_Alltoall(void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount, MPI_Datatype recvtype, MPI_Comm comm) {
1178   int retval, size, sendsize;
1179   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1180
1181   smpi_bench_end(rank, "Alltoall");
1182   if(comm == MPI_COMM_NULL) {
1183     retval = MPI_ERR_COMM;
1184   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
1185     retval = MPI_ERR_TYPE;
1186   } else {
1187     size = smpi_comm_size(comm);
1188     sendsize = smpi_datatype_size(sendtype) * sendcount;
1189     if(sendsize < 200 && size > 12) {
1190       retval = smpi_coll_tuned_alltoall_bruck(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
1191     } else if(sendsize < 3000) {
1192       retval = smpi_coll_tuned_alltoall_basic_linear(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
1193     } else {
1194       retval = smpi_coll_tuned_alltoall_pairwise(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
1195     }
1196   }
1197   smpi_bench_begin(rank, "Alltoall");
1198   return retval;
1199 }
1200
1201 int MPI_Alltoallv(void* sendbuf, int* sendcounts, int* senddisps, MPI_Datatype sendtype, void* recvbuf, int *recvcounts, int* recvdisps, MPI_Datatype recvtype, MPI_Comm comm) {
1202   int retval;
1203   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1204
1205   smpi_bench_end(rank, "Alltoallv");
1206   if(comm == MPI_COMM_NULL) {
1207     retval = MPI_ERR_COMM;
1208   } else if(sendtype == MPI_DATATYPE_NULL || recvtype == MPI_DATATYPE_NULL) {
1209     retval = MPI_ERR_TYPE;
1210   } else if(sendcounts == NULL || senddisps == NULL || recvcounts == NULL || recvdisps == NULL) {
1211     retval = MPI_ERR_ARG;
1212   } else {
1213     retval = smpi_coll_basic_alltoallv(sendbuf, sendcounts, senddisps, sendtype, recvbuf, recvcounts, recvdisps, recvtype, comm);
1214   }
1215   smpi_bench_begin(rank, "Alltoallv");
1216   return retval;
1217 }
1218
1219
1220 int MPI_Get_processor_name( char *name, int *resultlen ) {
1221   int retval = MPI_SUCCESS;
1222   smpi_bench_end(-1, NULL);
1223   strncpy( name , SIMIX_host_get_name(SIMIX_host_self()), MPI_MAX_PROCESSOR_NAME-1);
1224   *resultlen= strlen(name) > MPI_MAX_PROCESSOR_NAME ? MPI_MAX_PROCESSOR_NAME : strlen(name);
1225
1226   smpi_bench_begin(-1, NULL);
1227   return retval;
1228 }
1229