Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
add option to generate states for code outside smpi to allow computation timing ...
[simgrid.git] / src / smpi / smpi_pmpi.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_mpi_dt_private.h"
9
10 XBT_LOG_NEW_DEFAULT_SUBCATEGORY(smpi_pmpi, smpi,
11                                 "Logging specific to SMPI (pmpi)");
12
13 #ifdef HAVE_TRACING
14 //this function need to be here because of the calls to smpi_bench
15 void TRACE_smpi_set_category(const char *category)
16 {
17   //need to end bench otherwise categories for execution tasks are wrong
18   smpi_bench_end();
19   TRACE_internal_smpi_set_category (category);
20   //begin bench after changing process's category
21   smpi_bench_begin();
22 }
23 #endif
24
25 /* PMPI User level calls */
26
27 int PMPI_Init(int *argc, char ***argv)
28 {
29   smpi_process_init(argc, argv);
30 #ifdef HAVE_TRACING
31   int rank = smpi_process_index();
32   TRACE_smpi_init(rank);
33
34   TRACE_smpi_computing_init(rank);
35 #endif
36   smpi_bench_begin();
37   return MPI_SUCCESS;
38 }
39
40 int PMPI_Finalize(void)
41 {
42   smpi_process_finalize();
43   smpi_bench_end();
44 #ifdef HAVE_TRACING
45   int rank = smpi_process_index();
46   TRACE_smpi_computing_out(rank);
47   TRACE_smpi_finalize(smpi_process_index());
48 #endif
49   smpi_process_destroy();
50   return MPI_SUCCESS;
51 }
52
53 int PMPI_Init_thread(int *argc, char ***argv, int required, int *provided)
54 {
55   if (provided != NULL) {
56     *provided = MPI_THREAD_MULTIPLE;
57   }
58   return MPI_Init(argc, argv);
59 }
60
61 int PMPI_Query_thread(int *provided)
62 {
63   int retval;
64
65   smpi_bench_end();
66   if (provided == NULL) {
67     retval = MPI_ERR_ARG;
68   } else {
69     *provided = MPI_THREAD_MULTIPLE;
70     retval = MPI_SUCCESS;
71   }
72   smpi_bench_begin();
73   return retval;
74 }
75
76 int PMPI_Is_thread_main(int *flag)
77 {
78   int retval;
79
80   smpi_bench_end();
81   if (flag == NULL) {
82     retval = MPI_ERR_ARG;
83   } else {
84     *flag = smpi_process_index() == 0;
85     retval = MPI_SUCCESS;
86   }
87   smpi_bench_begin();
88   return retval;
89 }
90
91 int PMPI_Abort(MPI_Comm comm, int errorcode)
92 {
93   smpi_bench_end();
94   smpi_process_destroy();
95 #ifdef HAVE_TRACING
96   int rank = smpi_process_index();
97   TRACE_smpi_computing_out(rank);
98 #endif
99   // FIXME: should kill all processes in comm instead
100   simcall_process_kill(SIMIX_process_self());
101   return MPI_SUCCESS;
102 }
103
104 double PMPI_Wtime(void)
105 {
106   double time;
107
108   smpi_bench_end();
109   time = SIMIX_get_clock();
110   smpi_bench_begin();
111   return time;
112 }
113
114 int PMPI_Address(void *location, MPI_Aint * address)
115 {
116   int retval;
117
118   smpi_bench_end();
119   if (!address) {
120     retval = MPI_ERR_ARG;
121   } else {
122     *address = (MPI_Aint) location;
123   }
124   smpi_bench_begin();
125   return retval;
126 }
127
128 int PMPI_Type_free(MPI_Datatype * datatype)
129 {
130   int retval;
131
132   smpi_bench_end();
133   if (!datatype) {
134     retval = MPI_ERR_ARG;
135   } else {
136     // FIXME: always fail for now
137     retval = MPI_ERR_TYPE;
138   }
139   smpi_bench_begin();
140   return retval;
141 }
142
143 int PMPI_Type_size(MPI_Datatype datatype, int *size)
144 {
145   int retval;
146
147   smpi_bench_end();
148   if (datatype == MPI_DATATYPE_NULL) {
149     retval = MPI_ERR_TYPE;
150   } else if (size == NULL) {
151     retval = MPI_ERR_ARG;
152   } else {
153     *size = (int) smpi_datatype_size(datatype);
154     retval = MPI_SUCCESS;
155   }
156   smpi_bench_begin();
157   return retval;
158 }
159
160 int PMPI_Type_get_extent(MPI_Datatype datatype, MPI_Aint * lb, MPI_Aint * extent)
161 {
162   int retval;
163
164   smpi_bench_end();
165   if (datatype == MPI_DATATYPE_NULL) {
166     retval = MPI_ERR_TYPE;
167   } else if (lb == NULL || extent == NULL) {
168     retval = MPI_ERR_ARG;
169   } else {
170     retval = smpi_datatype_extent(datatype, lb, extent);
171   }
172   smpi_bench_begin();
173   return retval;
174 }
175
176 int PMPI_Type_extent(MPI_Datatype datatype, MPI_Aint * extent)
177 {
178   int retval;
179   MPI_Aint dummy;
180
181   smpi_bench_end();
182   if (datatype == MPI_DATATYPE_NULL) {
183     retval = MPI_ERR_TYPE;
184   } else if (extent == NULL) {
185     retval = MPI_ERR_ARG;
186   } else {
187     retval = smpi_datatype_extent(datatype, &dummy, extent);
188   }
189   smpi_bench_begin();
190   return retval;
191 }
192
193 int PMPI_Type_lb(MPI_Datatype datatype, MPI_Aint * disp)
194 {
195   int retval;
196
197   smpi_bench_end();
198   if (datatype == MPI_DATATYPE_NULL) {
199     retval = MPI_ERR_TYPE;
200   } else if (disp == NULL) {
201     retval = MPI_ERR_ARG;
202   } else {
203     *disp = smpi_datatype_lb(datatype);
204     retval = MPI_SUCCESS;
205   }
206   smpi_bench_begin();
207   return retval;
208 }
209
210 int PMPI_Type_ub(MPI_Datatype datatype, MPI_Aint * disp)
211 {
212   int retval;
213
214   smpi_bench_end();
215   if (datatype == MPI_DATATYPE_NULL) {
216     retval = MPI_ERR_TYPE;
217   } else if (disp == NULL) {
218     retval = MPI_ERR_ARG;
219   } else {
220     *disp = smpi_datatype_ub(datatype);
221     retval = MPI_SUCCESS;
222   }
223   smpi_bench_begin();
224   return retval;
225 }
226
227 int PMPI_Op_create(MPI_User_function * function, int commute, MPI_Op * op)
228 {
229   int retval;
230
231   smpi_bench_end();
232   if (function == NULL || op == NULL) {
233     retval = MPI_ERR_ARG;
234   } else {
235     *op = smpi_op_new(function, commute);
236     retval = MPI_SUCCESS;
237   }
238   smpi_bench_begin();
239   return retval;
240 }
241
242 int PMPI_Op_free(MPI_Op * op)
243 {
244   int retval;
245
246   smpi_bench_end();
247   if (op == NULL) {
248     retval = MPI_ERR_ARG;
249   } else if (*op == MPI_OP_NULL) {
250     retval = MPI_ERR_OP;
251   } else {
252     smpi_op_destroy(*op);
253     *op = MPI_OP_NULL;
254     retval = MPI_SUCCESS;
255   }
256   smpi_bench_begin();
257   return retval;
258 }
259
260 int PMPI_Group_free(MPI_Group * group)
261 {
262   int retval;
263
264   smpi_bench_end();
265   if (group == NULL) {
266     retval = MPI_ERR_ARG;
267   } else {
268     smpi_group_destroy(*group);
269     *group = MPI_GROUP_NULL;
270     retval = MPI_SUCCESS;
271   }
272   smpi_bench_begin();
273   return retval;
274 }
275
276 int PMPI_Group_size(MPI_Group group, int *size)
277 {
278   int retval;
279
280   smpi_bench_end();
281   if (group == MPI_GROUP_NULL) {
282     retval = MPI_ERR_GROUP;
283   } else if (size == NULL) {
284     retval = MPI_ERR_ARG;
285   } else {
286     *size = smpi_group_size(group);
287     retval = MPI_SUCCESS;
288   }
289   smpi_bench_begin();
290   return retval;
291 }
292
293 int PMPI_Group_rank(MPI_Group group, int *rank)
294 {
295   int retval;
296
297   smpi_bench_end();
298   if (group == MPI_GROUP_NULL) {
299     retval = MPI_ERR_GROUP;
300   } else if (rank == NULL) {
301     retval = MPI_ERR_ARG;
302   } else {
303     *rank = smpi_group_rank(group, smpi_process_index());
304     retval = MPI_SUCCESS;
305   }
306   smpi_bench_begin();
307   return retval;
308 }
309
310 int PMPI_Group_translate_ranks(MPI_Group group1, int n, int *ranks1,
311                               MPI_Group group2, int *ranks2)
312 {
313   int retval, i, index;
314
315   smpi_bench_end();
316   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
317     retval = MPI_ERR_GROUP;
318   } else {
319     for (i = 0; i < n; i++) {
320       index = smpi_group_index(group1, ranks1[i]);
321       ranks2[i] = smpi_group_rank(group2, index);
322     }
323     retval = MPI_SUCCESS;
324   }
325   smpi_bench_begin();
326   return retval;
327 }
328
329 int PMPI_Group_compare(MPI_Group group1, MPI_Group group2, int *result)
330 {
331   int retval;
332
333   smpi_bench_end();
334   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
335     retval = MPI_ERR_GROUP;
336   } else if (result == NULL) {
337     retval = MPI_ERR_ARG;
338   } else {
339     *result = smpi_group_compare(group1, group2);
340     retval = MPI_SUCCESS;
341   }
342   smpi_bench_begin();
343   return retval;
344 }
345
346 int PMPI_Group_union(MPI_Group group1, MPI_Group group2,
347                     MPI_Group * newgroup)
348 {
349   int retval, i, proc1, proc2, size, size2;
350
351   smpi_bench_end();
352   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
353     retval = MPI_ERR_GROUP;
354   } else if (newgroup == NULL) {
355     retval = MPI_ERR_ARG;
356   } else {
357     size = smpi_group_size(group1);
358     size2 = smpi_group_size(group2);
359     for (i = 0; i < size2; i++) {
360       proc2 = smpi_group_index(group2, i);
361       proc1 = smpi_group_rank(group1, proc2);
362       if (proc1 == MPI_UNDEFINED) {
363         size++;
364       }
365     }
366     if (size == 0) {
367       *newgroup = MPI_GROUP_EMPTY;
368     } else {
369       *newgroup = smpi_group_new(size);
370       size2 = smpi_group_size(group1);
371       for (i = 0; i < size2; i++) {
372         proc1 = smpi_group_index(group1, i);
373         smpi_group_set_mapping(*newgroup, proc1, i);
374       }
375       for (i = size2; i < size; i++) {
376         proc2 = smpi_group_index(group2, i - size2);
377         smpi_group_set_mapping(*newgroup, proc2, i);
378       }
379     }
380     smpi_group_use(*newgroup);
381     retval = MPI_SUCCESS;
382   }
383   smpi_bench_begin();
384   return retval;
385 }
386
387 int PMPI_Group_intersection(MPI_Group group1, MPI_Group group2,
388                            MPI_Group * newgroup)
389 {
390   int retval, i, proc1, proc2, size, size2;
391
392   smpi_bench_end();
393   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
394     retval = MPI_ERR_GROUP;
395   } else if (newgroup == NULL) {
396     retval = MPI_ERR_ARG;
397   } else {
398     size = smpi_group_size(group1);
399     size2 = smpi_group_size(group2);
400     for (i = 0; i < size2; i++) {
401       proc2 = smpi_group_index(group2, i);
402       proc1 = smpi_group_rank(group1, proc2);
403       if (proc1 == MPI_UNDEFINED) {
404         size--;
405       }
406     }
407     if (size == 0) {
408       *newgroup = MPI_GROUP_EMPTY;
409     } else {
410       *newgroup = smpi_group_new(size);
411       size2 = smpi_group_size(group1);
412       for (i = 0; i < size2; i++) {
413         proc1 = smpi_group_index(group1, i);
414         proc2 = smpi_group_rank(group2, proc1);
415         if (proc2 != MPI_UNDEFINED) {
416           smpi_group_set_mapping(*newgroup, proc1, i);
417         }
418       }
419     }
420     smpi_group_use(*newgroup);
421     retval = MPI_SUCCESS;
422   }
423   smpi_bench_begin();
424   return retval;
425 }
426
427 int PMPI_Group_difference(MPI_Group group1, MPI_Group group2, MPI_Group * newgroup)
428 {
429   int retval, i, proc1, proc2, size, size2;
430
431   smpi_bench_end();
432   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
433     retval = MPI_ERR_GROUP;
434   } else if (newgroup == NULL) {
435     retval = MPI_ERR_ARG;
436   } else {
437     size = size2 = smpi_group_size(group1);
438     for (i = 0; i < size2; i++) {
439       proc1 = smpi_group_index(group1, i);
440       proc2 = smpi_group_rank(group2, proc1);
441       if (proc2 != MPI_UNDEFINED) {
442         size--;
443       }
444     }
445     if (size == 0) {
446       *newgroup = MPI_GROUP_EMPTY;
447     } else {
448       *newgroup = smpi_group_new(size);
449       for (i = 0; i < size2; i++) {
450         proc1 = smpi_group_index(group1, i);
451         proc2 = smpi_group_rank(group2, proc1);
452         if (proc2 == MPI_UNDEFINED) {
453           smpi_group_set_mapping(*newgroup, proc1, i);
454         }
455       }
456     }
457     smpi_group_use(*newgroup);
458     retval = MPI_SUCCESS;
459   }
460   smpi_bench_begin();
461   return retval;
462 }
463
464 int PMPI_Group_incl(MPI_Group group, int n, int *ranks, MPI_Group * newgroup)
465 {
466   int retval, i, index;
467
468   smpi_bench_end();
469   if (group == MPI_GROUP_NULL) {
470     retval = MPI_ERR_GROUP;
471   } else if (newgroup == NULL) {
472     retval = MPI_ERR_ARG;
473   } else {
474     if (n == 0) {
475       *newgroup = MPI_GROUP_EMPTY;
476     } else if (n == smpi_group_size(group)) {
477       *newgroup = group;
478     } else {
479       *newgroup = smpi_group_new(n);
480       for (i = 0; i < n; i++) {
481         index = smpi_group_index(group, ranks[i]);
482         smpi_group_set_mapping(*newgroup, index, i);
483       }
484     }
485     smpi_group_use(*newgroup);
486     retval = MPI_SUCCESS;
487   }
488   smpi_bench_begin();
489   return retval;
490 }
491
492 int PMPI_Group_excl(MPI_Group group, int n, int *ranks, MPI_Group * newgroup)
493 {
494   int retval, i, size, rank, index;
495
496   smpi_bench_end();
497   if (group == MPI_GROUP_NULL) {
498     retval = MPI_ERR_GROUP;
499   } else if (newgroup == NULL) {
500     retval = MPI_ERR_ARG;
501   } else {
502     if (n == 0) {
503       *newgroup = group;
504     } else if (n == smpi_group_size(group)) {
505       *newgroup = MPI_GROUP_EMPTY;
506     } else {
507       size = smpi_group_size(group) - n;
508       *newgroup = smpi_group_new(size);
509       rank = 0;
510       while (rank < size) {
511         for (i = 0; i < n; i++) {
512           if (ranks[i] == rank) {
513             break;
514           }
515         }
516         if (i >= n) {
517           index = smpi_group_index(group, rank);
518           smpi_group_set_mapping(*newgroup, index, rank);
519           rank++;
520         }
521       }
522     }
523     smpi_group_use(*newgroup);
524     retval = MPI_SUCCESS;
525   }
526   smpi_bench_begin();
527   return retval;
528 }
529
530 int PMPI_Group_range_incl(MPI_Group group, int n, int ranges[][3],
531                          MPI_Group * newgroup)
532 {
533   int retval, i, j, rank, size, index;
534
535   smpi_bench_end();
536   if (group == MPI_GROUP_NULL) {
537     retval = MPI_ERR_GROUP;
538   } else if (newgroup == NULL) {
539     retval = MPI_ERR_ARG;
540   } else {
541     if (n == 0) {
542       *newgroup = MPI_GROUP_EMPTY;
543     } else {
544       size = 0;
545       for (i = 0; i < n; i++) {
546         for (rank = ranges[i][0];       /* First */
547              rank >= 0 && rank <= ranges[i][1]; /* Last */
548              rank += ranges[i][2] /* Stride */ ) {
549           size++;
550         }
551       }
552       if (size == smpi_group_size(group)) {
553         *newgroup = group;
554       } else {
555         *newgroup = smpi_group_new(size);
556         j = 0;
557         for (i = 0; i < n; i++) {
558           for (rank = ranges[i][0];     /* First */
559                rank >= 0 && rank <= ranges[i][1];       /* Last */
560                rank += ranges[i][2] /* Stride */ ) {
561             index = smpi_group_index(group, rank);
562             smpi_group_set_mapping(*newgroup, index, j);
563             j++;
564           }
565         }
566       }
567     }
568     smpi_group_use(*newgroup);
569     retval = MPI_SUCCESS;
570   }
571   smpi_bench_begin();
572   return retval;
573 }
574
575 int PMPI_Group_range_excl(MPI_Group group, int n, int ranges[][3],
576                          MPI_Group * newgroup)
577 {
578   int retval, i, newrank, rank, size, index, add;
579
580   smpi_bench_end();
581   if (group == MPI_GROUP_NULL) {
582     retval = MPI_ERR_GROUP;
583   } else if (newgroup == NULL) {
584     retval = MPI_ERR_ARG;
585   } else {
586     if (n == 0) {
587       *newgroup = group;
588     } else {
589       size = smpi_group_size(group);
590       for (i = 0; i < n; i++) {
591         for (rank = ranges[i][0];       /* First */
592              rank >= 0 && rank <= ranges[i][1]; /* Last */
593              rank += ranges[i][2] /* Stride */ ) {
594           size--;
595         }
596       }
597       if (size == 0) {
598         *newgroup = MPI_GROUP_EMPTY;
599       } else {
600         *newgroup = smpi_group_new(size);
601         newrank = 0;
602         while (newrank < size) {
603           for (i = 0; i < n; i++) {
604             add = 1;
605             for (rank = ranges[i][0];   /* First */
606                  rank >= 0 && rank <= ranges[i][1];     /* Last */
607                  rank += ranges[i][2] /* Stride */ ) {
608               if (rank == newrank) {
609                 add = 0;
610                 break;
611               }
612             }
613             if (add == 1) {
614               index = smpi_group_index(group, newrank);
615               smpi_group_set_mapping(*newgroup, index, newrank);
616             }
617           }
618         }
619       }
620     }
621     smpi_group_use(*newgroup);
622     retval = MPI_SUCCESS;
623   }
624   smpi_bench_begin();
625   return retval;
626 }
627
628 int PMPI_Comm_rank(MPI_Comm comm, int *rank)
629 {
630   int retval;
631
632   smpi_bench_end();
633   if (comm == MPI_COMM_NULL) {
634     retval = MPI_ERR_COMM;
635   } else if (rank == NULL) {
636     retval = MPI_ERR_ARG;
637   } else {
638     *rank = smpi_comm_rank(comm);
639     retval = MPI_SUCCESS;
640   }
641   smpi_bench_begin();
642   return retval;
643 }
644
645 int PMPI_Comm_size(MPI_Comm comm, int *size)
646 {
647   int retval;
648
649   smpi_bench_end();
650   if (comm == MPI_COMM_NULL) {
651     retval = MPI_ERR_COMM;
652   } else if (size == NULL) {
653     retval = MPI_ERR_ARG;
654   } else {
655     *size = smpi_comm_size(comm);
656     retval = MPI_SUCCESS;
657   }
658   smpi_bench_begin();
659   return retval;
660 }
661
662 int PMPI_Comm_get_name (MPI_Comm comm, char* name, int* len)
663 {
664   int retval;
665
666   smpi_bench_end();
667   if (comm == MPI_COMM_NULL)  {
668     retval = MPI_ERR_COMM;
669   } else if (name == NULL || len == NULL)  {
670     retval = MPI_ERR_ARG;
671   } else {
672     smpi_comm_get_name(comm, name, len);
673     retval = MPI_SUCCESS;
674   }
675   smpi_bench_begin();
676   return retval;
677 }
678
679 int PMPI_Comm_group(MPI_Comm comm, MPI_Group * group)
680 {
681   int retval;
682
683   smpi_bench_end();
684   if (comm == MPI_COMM_NULL) {
685     retval = MPI_ERR_COMM;
686   } else if (group == NULL) {
687     retval = MPI_ERR_ARG;
688   } else {
689     *group = smpi_comm_group(comm);
690     retval = MPI_SUCCESS;
691   }
692   smpi_bench_begin();
693   return retval;
694 }
695
696 int PMPI_Comm_compare(MPI_Comm comm1, MPI_Comm comm2, int *result)
697 {
698   int retval;
699
700   smpi_bench_end();
701   if (comm1 == MPI_COMM_NULL || comm2 == MPI_COMM_NULL) {
702     retval = MPI_ERR_COMM;
703   } else if (result == NULL) {
704     retval = MPI_ERR_ARG;
705   } else {
706     if (comm1 == comm2) {       /* Same communicators means same groups */
707       *result = MPI_IDENT;
708     } else {
709       *result =
710           smpi_group_compare(smpi_comm_group(comm1),
711                              smpi_comm_group(comm2));
712       if (*result == MPI_IDENT) {
713         *result = MPI_CONGRUENT;
714       }
715     }
716     retval = MPI_SUCCESS;
717   }
718   smpi_bench_begin();
719   return retval;
720 }
721
722 int PMPI_Comm_dup(MPI_Comm comm, MPI_Comm * newcomm)
723 {
724   int retval;
725
726   smpi_bench_end();
727   if (comm == MPI_COMM_NULL) {
728     retval = MPI_ERR_COMM;
729   } else if (newcomm == NULL) {
730     retval = MPI_ERR_ARG;
731   } else {
732     *newcomm = smpi_comm_new(smpi_comm_group(comm));
733     retval = MPI_SUCCESS;
734   }
735   smpi_bench_begin();
736   return retval;
737 }
738
739 int PMPI_Comm_create(MPI_Comm comm, MPI_Group group, MPI_Comm * newcomm)
740 {
741   int retval;
742
743   smpi_bench_end();
744   if (comm == MPI_COMM_NULL) {
745     retval = MPI_ERR_COMM;
746   } else if (group == MPI_GROUP_NULL) {
747     retval = MPI_ERR_GROUP;
748   } else if (newcomm == NULL) {
749     retval = MPI_ERR_ARG;
750   } else {
751     *newcomm = smpi_comm_new(group);
752     retval = MPI_SUCCESS;
753   }
754   smpi_bench_begin();
755   return retval;
756 }
757
758 int PMPI_Comm_free(MPI_Comm * comm)
759 {
760   int retval;
761
762   smpi_bench_end();
763   if (comm == NULL) {
764     retval = MPI_ERR_ARG;
765   } else if (*comm == MPI_COMM_NULL) {
766     retval = MPI_ERR_COMM;
767   } else {
768     smpi_comm_destroy(*comm);
769     *comm = MPI_COMM_NULL;
770     retval = MPI_SUCCESS;
771   }
772   smpi_bench_begin();
773   return retval;
774 }
775
776 int PMPI_Comm_disconnect(MPI_Comm * comm)
777 {
778   /* TODO: wait until all communication in comm are done */
779   int retval;
780
781   smpi_bench_end();
782   if (comm == NULL) {
783     retval = MPI_ERR_ARG;
784   } else if (*comm == MPI_COMM_NULL) {
785     retval = MPI_ERR_COMM;
786   } else {
787     smpi_comm_destroy(*comm);
788     *comm = MPI_COMM_NULL;
789     retval = MPI_SUCCESS;
790   }
791   smpi_bench_begin();
792   return retval;
793 }
794
795 int PMPI_Comm_split(MPI_Comm comm, int color, int key, MPI_Comm* comm_out)
796 {
797   int retval;
798
799   smpi_bench_end();
800   if (comm_out == NULL) {
801     retval = MPI_ERR_ARG;
802   } else if (comm == MPI_COMM_NULL) {
803     retval = MPI_ERR_COMM;
804   } else {
805     *comm_out = smpi_comm_split(comm, color, key);
806     retval = MPI_SUCCESS;
807   }
808   smpi_bench_begin();
809   return retval;
810 }
811
812 int PMPI_Send_init(void *buf, int count, MPI_Datatype datatype, int dst,
813                   int tag, MPI_Comm comm, MPI_Request * request)
814 {
815   int retval;
816
817   smpi_bench_end();
818   if (request == NULL) {
819     retval = MPI_ERR_ARG;
820   } else if (comm == MPI_COMM_NULL) {
821     retval = MPI_ERR_COMM;
822   } else {
823     *request = smpi_mpi_send_init(buf, count, datatype, dst, tag, comm);
824     retval = MPI_SUCCESS;
825   }
826   smpi_bench_begin();
827   return retval;
828 }
829
830 int PMPI_Recv_init(void *buf, int count, MPI_Datatype datatype, int src,
831                   int tag, MPI_Comm comm, MPI_Request * request)
832 {
833   int retval;
834
835   smpi_bench_end();
836   if (request == NULL) {
837     retval = MPI_ERR_ARG;
838   } else if (comm == MPI_COMM_NULL) {
839     retval = MPI_ERR_COMM;
840   } else {
841     *request = smpi_mpi_recv_init(buf, count, datatype, src, tag, comm);
842     retval = MPI_SUCCESS;
843   }
844   smpi_bench_begin();
845   return retval;
846 }
847
848 int PMPI_Start(MPI_Request * request)
849 {
850   int retval;
851
852   smpi_bench_end();
853   if (request == NULL) {
854     retval = MPI_ERR_ARG;
855   } else {
856     smpi_mpi_start(*request);
857     retval = MPI_SUCCESS;
858   }
859   smpi_bench_begin();
860   return retval;
861 }
862
863 int PMPI_Startall(int count, MPI_Request * requests)
864 {
865   int retval;
866
867   smpi_bench_end();
868   if (requests == NULL) {
869     retval = MPI_ERR_ARG;
870   } else {
871     smpi_mpi_startall(count, requests);
872     retval = MPI_SUCCESS;
873   }
874   smpi_bench_begin();
875   return retval;
876 }
877
878 int PMPI_Request_free(MPI_Request * request)
879 {
880   int retval;
881
882   smpi_bench_end();
883   if (request == NULL) {
884     retval = MPI_ERR_ARG;
885   } else {
886     smpi_mpi_request_free(request);
887     retval = MPI_SUCCESS;
888   }
889   smpi_bench_begin();
890   return retval;
891 }
892
893 int PMPI_Irecv(void *buf, int count, MPI_Datatype datatype, int src,
894               int tag, MPI_Comm comm, MPI_Request * request)
895 {
896   int retval;
897
898   smpi_bench_end();
899 #ifdef HAVE_TRACING
900   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
901   int src_traced = smpi_group_rank(smpi_comm_group(comm), src);
902   TRACE_smpi_ptp_in(rank, src_traced, rank, __FUNCTION__);
903 #endif
904   if (request == NULL) {
905     retval = MPI_ERR_ARG;
906   } else if (comm == MPI_COMM_NULL) {
907     retval = MPI_ERR_COMM;
908   } else {
909     *request = smpi_mpi_irecv(buf, count, datatype, src, tag, comm);
910     retval = MPI_SUCCESS;
911   }
912 #ifdef HAVE_TRACING
913   TRACE_smpi_ptp_out(rank, src_traced, rank, __FUNCTION__);
914   (*request)->recv = 1;
915 #endif
916   smpi_bench_begin();
917   return retval;
918 }
919
920 int PMPI_Isend(void *buf, int count, MPI_Datatype datatype, int dst,
921               int tag, MPI_Comm comm, MPI_Request * request)
922 {
923   int retval;
924
925   smpi_bench_end();
926 #ifdef HAVE_TRACING
927   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
928   TRACE_smpi_computing_out(rank);
929   int dst_traced = smpi_group_rank(smpi_comm_group(comm), dst);
930   TRACE_smpi_ptp_in(rank, rank, dst_traced, __FUNCTION__);
931   TRACE_smpi_send(rank, rank, dst_traced);
932 #endif
933   if (request == NULL) {
934     retval = MPI_ERR_ARG;
935   } else if (comm == MPI_COMM_NULL) {
936     retval = MPI_ERR_COMM;
937   } else {
938     *request = smpi_mpi_isend(buf, count, datatype, dst, tag, comm);
939     retval = MPI_SUCCESS;
940   }
941 #ifdef HAVE_TRACING
942   TRACE_smpi_ptp_out(rank, rank, dst_traced, __FUNCTION__);
943   (*request)->send = 1;
944   TRACE_smpi_computing_in(rank);
945 #endif
946   smpi_bench_begin();
947   return retval;
948 }
949
950 int PMPI_Recv(void *buf, int count, MPI_Datatype datatype, int src, int tag,
951              MPI_Comm comm, MPI_Status * status)
952 {
953   int retval;
954
955   smpi_bench_end();
956 #ifdef HAVE_TRACING
957   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
958   int src_traced = smpi_group_rank(smpi_comm_group(comm), src);
959   TRACE_smpi_computing_out(rank);
960
961   TRACE_smpi_ptp_in(rank, src_traced, rank, __FUNCTION__);
962 #endif
963   if (comm == MPI_COMM_NULL) {
964     retval = MPI_ERR_COMM;
965   } else {
966     smpi_mpi_recv(buf, count, datatype, src, tag, comm, status);
967     retval = MPI_SUCCESS;
968   }
969 #ifdef HAVE_TRACING
970   TRACE_smpi_ptp_out(rank, src_traced, rank, __FUNCTION__);
971   TRACE_smpi_recv(rank, src_traced, rank);
972   TRACE_smpi_computing_in(rank);
973 #endif
974   smpi_bench_begin();
975   return retval;
976 }
977
978 int PMPI_Send(void *buf, int count, MPI_Datatype datatype, int dst, int tag,
979              MPI_Comm comm)
980 {
981   int retval;
982
983   smpi_bench_end();
984 #ifdef HAVE_TRACING
985   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
986   TRACE_smpi_computing_out(rank);
987   int dst_traced = smpi_group_rank(smpi_comm_group(comm), dst);
988   TRACE_smpi_ptp_in(rank, rank, dst_traced, __FUNCTION__);
989   TRACE_smpi_send(rank, rank, dst_traced);
990 #endif
991   if (comm == MPI_COMM_NULL) {
992     retval = MPI_ERR_COMM;
993   } else {
994     smpi_mpi_send(buf, count, datatype, dst, tag, comm);
995     retval = MPI_SUCCESS;
996   }
997 #ifdef HAVE_TRACING
998   TRACE_smpi_ptp_out(rank, rank, dst_traced, __FUNCTION__);
999   TRACE_smpi_computing_in(rank);
1000 #endif
1001   smpi_bench_begin();
1002   return retval;
1003 }
1004
1005 int PMPI_Sendrecv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1006                  int dst, int sendtag, void *recvbuf, int recvcount,
1007                  MPI_Datatype recvtype, int src, int recvtag,
1008                  MPI_Comm comm, MPI_Status * status)
1009 {
1010   int retval;
1011
1012   smpi_bench_end();
1013 #ifdef HAVE_TRACING
1014   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1015   TRACE_smpi_computing_out(rank);
1016   int dst_traced = smpi_group_rank(smpi_comm_group(comm), dst);
1017   int src_traced = smpi_group_rank(smpi_comm_group(comm), src);
1018   TRACE_smpi_ptp_in(rank, src_traced, dst_traced, __FUNCTION__);
1019   TRACE_smpi_send(rank, rank, dst_traced);
1020   TRACE_smpi_send(rank, src_traced, rank);
1021 #endif
1022   if (comm == MPI_COMM_NULL) {
1023     retval = MPI_ERR_COMM;
1024   } else if (sendtype == MPI_DATATYPE_NULL
1025              || recvtype == MPI_DATATYPE_NULL) {
1026     retval = MPI_ERR_TYPE;
1027   } else {
1028     smpi_mpi_sendrecv(sendbuf, sendcount, sendtype, dst, sendtag, recvbuf,
1029                       recvcount, recvtype, src, recvtag, comm, status);
1030     retval = MPI_SUCCESS;
1031   }
1032 #ifdef HAVE_TRACING
1033   TRACE_smpi_ptp_out(rank, src_traced, dst_traced, __FUNCTION__);
1034   TRACE_smpi_recv(rank, rank, dst_traced);
1035   TRACE_smpi_recv(rank, src_traced, rank);
1036   TRACE_smpi_computing_in(rank);
1037
1038 #endif
1039   smpi_bench_begin();
1040   return retval;
1041 }
1042
1043 int PMPI_Sendrecv_replace(void *buf, int count, MPI_Datatype datatype,
1044                          int dst, int sendtag, int src, int recvtag,
1045                          MPI_Comm comm, MPI_Status * status)
1046 {
1047   //TODO: suboptimal implementation
1048   void *recvbuf;
1049   int retval, size;
1050
1051   size = smpi_datatype_size(datatype) * count;
1052   recvbuf = xbt_new(char, size);
1053   retval =
1054       MPI_Sendrecv(buf, count, datatype, dst, sendtag, recvbuf, count,
1055                    datatype, src, recvtag, comm, status);
1056   memcpy(buf, recvbuf, size * sizeof(char));
1057   xbt_free(recvbuf);
1058   return retval;
1059 }
1060
1061 int PMPI_Test(MPI_Request * request, int *flag, MPI_Status * status)
1062 {
1063   int retval;
1064
1065   smpi_bench_end();
1066   if (request == NULL || flag == NULL) {
1067     retval = MPI_ERR_ARG;
1068   } else if (*request == MPI_REQUEST_NULL) {
1069     retval = MPI_ERR_REQUEST;
1070   } else {
1071     *flag = smpi_mpi_test(request, status);
1072     retval = MPI_SUCCESS;
1073   }
1074   smpi_bench_begin();
1075   return retval;
1076 }
1077
1078 int PMPI_Testany(int count, MPI_Request requests[], int *index, int *flag,
1079                 MPI_Status * status)
1080 {
1081   int retval;
1082
1083   smpi_bench_end();
1084   if (index == NULL || flag == NULL) {
1085     retval = MPI_ERR_ARG;
1086   } else {
1087     *flag = smpi_mpi_testany(count, requests, index, status);
1088     retval = MPI_SUCCESS;
1089   }
1090   smpi_bench_begin();
1091   return retval;
1092 }
1093
1094 int PMPI_Wait(MPI_Request * request, MPI_Status * status)
1095 {
1096   int retval;
1097
1098   smpi_bench_end();
1099 #ifdef HAVE_TRACING
1100   int rank = request && (*request)->comm != MPI_COMM_NULL
1101       ? smpi_comm_rank((*request)->comm)
1102       : -1;
1103   TRACE_smpi_computing_out(rank);
1104
1105   MPI_Group group = smpi_comm_group((*request)->comm);
1106   int src_traced = smpi_group_rank(group, (*request)->src);
1107   int dst_traced = smpi_group_rank(group, (*request)->dst);
1108   int is_wait_for_receive = (*request)->recv;
1109   TRACE_smpi_ptp_in(rank, src_traced, dst_traced, __FUNCTION__);
1110 #endif
1111   if (request == NULL) {
1112     retval = MPI_ERR_ARG;
1113   } else if (*request == MPI_REQUEST_NULL) {
1114     retval = MPI_ERR_REQUEST;
1115   } else {
1116     smpi_mpi_wait(request, status);
1117     retval = MPI_SUCCESS;
1118   }
1119 #ifdef HAVE_TRACING
1120   TRACE_smpi_ptp_out(rank, src_traced, dst_traced, __FUNCTION__);
1121   if (is_wait_for_receive) {
1122     TRACE_smpi_recv(rank, src_traced, dst_traced);
1123   }
1124   TRACE_smpi_computing_in(rank);
1125
1126 #endif
1127   smpi_bench_begin();
1128   return retval;
1129 }
1130
1131 int PMPI_Waitany(int count, MPI_Request requests[], int *index, MPI_Status * status)
1132 {
1133   int retval;
1134
1135   smpi_bench_end();
1136 #ifdef HAVE_TRACING
1137   //save requests information for tracing
1138   int i;
1139   xbt_dynar_t srcs = xbt_dynar_new(sizeof(int), NULL);
1140   xbt_dynar_t dsts = xbt_dynar_new(sizeof(int), NULL);
1141   xbt_dynar_t recvs = xbt_dynar_new(sizeof(int), NULL);
1142   for (i = 0; i < count; i++) {
1143     MPI_Request req = requests[i];      //already received requests are no longer valid
1144     if (req) {
1145       int *asrc = xbt_new(int, 1);
1146       int *adst = xbt_new(int, 1);
1147       int *arecv = xbt_new(int, 1);
1148       *asrc = req->src;
1149       *adst = req->dst;
1150       *arecv = req->recv;
1151       xbt_dynar_insert_at(srcs, i, asrc);
1152       xbt_dynar_insert_at(dsts, i, adst);
1153       xbt_dynar_insert_at(recvs, i, arecv);
1154       xbt_free(asrc);
1155       xbt_free(adst);
1156       xbt_free(arecv);
1157     } else {
1158       int *t = xbt_new(int, 1);
1159       xbt_dynar_insert_at(srcs, i, t);
1160       xbt_dynar_insert_at(dsts, i, t);
1161       xbt_dynar_insert_at(recvs, i, t);
1162       xbt_free(t);
1163     }
1164   }
1165   int rank_traced = smpi_comm_rank(MPI_COMM_WORLD);
1166   TRACE_smpi_computing_out(rank_traced);
1167
1168   TRACE_smpi_ptp_in(rank_traced, -1, -1, __FUNCTION__);
1169
1170 #endif
1171   if (index == NULL) {
1172     retval = MPI_ERR_ARG;
1173   } else {
1174     *index = smpi_mpi_waitany(count, requests, status);
1175     retval = MPI_SUCCESS;
1176   }
1177 #ifdef HAVE_TRACING
1178   int src_traced, dst_traced, is_wait_for_receive;
1179   xbt_dynar_get_cpy(srcs, *index, &src_traced);
1180   xbt_dynar_get_cpy(dsts, *index, &dst_traced);
1181   xbt_dynar_get_cpy(recvs, *index, &is_wait_for_receive);
1182   if (is_wait_for_receive) {
1183     TRACE_smpi_recv(rank_traced, src_traced, dst_traced);
1184   }
1185   TRACE_smpi_ptp_out(rank_traced, src_traced, dst_traced, __FUNCTION__);
1186   //clean-up of dynars
1187   xbt_dynar_free(&srcs);
1188   xbt_dynar_free(&dsts);
1189   xbt_dynar_free(&recvs);
1190   TRACE_smpi_computing_in(rank_traced);
1191
1192 #endif
1193   smpi_bench_begin();
1194   return retval;
1195 }
1196
1197 int PMPI_Waitall(int count, MPI_Request requests[], MPI_Status status[])
1198 {
1199
1200   smpi_bench_end();
1201 #ifdef HAVE_TRACING
1202   //save information from requests
1203   int i;
1204   xbt_dynar_t srcs = xbt_dynar_new(sizeof(int), NULL);
1205   xbt_dynar_t dsts = xbt_dynar_new(sizeof(int), NULL);
1206   xbt_dynar_t recvs = xbt_dynar_new(sizeof(int), NULL);
1207   for (i = 0; i < count; i++) {
1208     MPI_Request req = requests[i];      //all req should be valid in Waitall
1209     int *asrc = xbt_new(int, 1);
1210     int *adst = xbt_new(int, 1);
1211     int *arecv = xbt_new(int, 1);
1212     *asrc = req->src;
1213     *adst = req->dst;
1214     *arecv = req->recv;
1215     xbt_dynar_insert_at(srcs, i, asrc);
1216     xbt_dynar_insert_at(dsts, i, adst);
1217     xbt_dynar_insert_at(recvs, i, arecv);
1218     xbt_free(asrc);
1219     xbt_free(adst);
1220     xbt_free(arecv);
1221   }
1222   int rank_traced = smpi_comm_rank (MPI_COMM_WORLD);
1223   TRACE_smpi_computing_out(rank_traced);
1224
1225   TRACE_smpi_ptp_in(rank_traced, -1, -1, __FUNCTION__);
1226 #endif
1227   smpi_mpi_waitall(count, requests, status);
1228 #ifdef HAVE_TRACING
1229   for (i = 0; i < count; i++) {
1230     int src_traced, dst_traced, is_wait_for_receive;
1231     xbt_dynar_get_cpy(srcs, i, &src_traced);
1232     xbt_dynar_get_cpy(dsts, i, &dst_traced);
1233     xbt_dynar_get_cpy(recvs, i, &is_wait_for_receive);
1234     if (is_wait_for_receive) {
1235       TRACE_smpi_recv(rank_traced, src_traced, dst_traced);
1236     }
1237   }
1238   TRACE_smpi_ptp_out(rank_traced, -1, -1, __FUNCTION__);
1239   //clean-up of dynars
1240   xbt_dynar_free(&srcs);
1241   xbt_dynar_free(&dsts);
1242   xbt_dynar_free(&recvs);
1243   TRACE_smpi_computing_in(rank_traced);
1244 #endif
1245   smpi_bench_begin();
1246   return MPI_SUCCESS;
1247 }
1248
1249 int PMPI_Waitsome(int incount, MPI_Request requests[], int *outcount,
1250                  int *indices, MPI_Status status[])
1251 {
1252   int retval;
1253
1254   smpi_bench_end();
1255   if (outcount == NULL || indices == NULL) {
1256     retval = MPI_ERR_ARG;
1257   } else {
1258     *outcount = smpi_mpi_waitsome(incount, requests, indices, status);
1259     retval = MPI_SUCCESS;
1260   }
1261   smpi_bench_begin();
1262   return retval;
1263 }
1264
1265 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
1266 {
1267   int retval;
1268
1269   smpi_bench_end();
1270 #ifdef HAVE_TRACING
1271   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1272   TRACE_smpi_computing_out(rank);
1273   int root_traced = smpi_group_rank(smpi_comm_group(comm), root);
1274   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__);
1275 #endif
1276   if (comm == MPI_COMM_NULL) {
1277     retval = MPI_ERR_COMM;
1278   } else {
1279     smpi_mpi_bcast(buf, count, datatype, root, comm);
1280     retval = MPI_SUCCESS;
1281   }
1282 #ifdef HAVE_TRACING
1283   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1284   TRACE_smpi_computing_in(rank);
1285 #endif
1286   smpi_bench_begin();
1287   return retval;
1288 }
1289
1290 int PMPI_Barrier(MPI_Comm comm)
1291 {
1292   int retval;
1293
1294   smpi_bench_end();
1295 #ifdef HAVE_TRACING
1296   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1297   TRACE_smpi_computing_out(rank);
1298   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1299 #endif
1300   if (comm == MPI_COMM_NULL) {
1301     retval = MPI_ERR_COMM;
1302   } else {
1303     smpi_mpi_barrier(comm);
1304     retval = MPI_SUCCESS;
1305   }
1306 #ifdef HAVE_TRACING
1307   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1308   TRACE_smpi_computing_in(rank);
1309 #endif
1310   smpi_bench_begin();
1311   return retval;
1312 }
1313
1314 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1315                void *recvbuf, int recvcount, MPI_Datatype recvtype,
1316                int root, MPI_Comm comm)
1317 {
1318   int retval;
1319
1320   smpi_bench_end();
1321 #ifdef HAVE_TRACING
1322   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1323   TRACE_smpi_computing_out(rank);
1324   int root_traced = smpi_group_rank(smpi_comm_group(comm), root);
1325   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__);
1326 #endif
1327   if (comm == MPI_COMM_NULL) {
1328     retval = MPI_ERR_COMM;
1329   } else if (sendtype == MPI_DATATYPE_NULL
1330              || recvtype == MPI_DATATYPE_NULL) {
1331     retval = MPI_ERR_TYPE;
1332   } else {
1333     smpi_mpi_gather(sendbuf, sendcount, sendtype, recvbuf, recvcount,
1334                     recvtype, root, comm);
1335     retval = MPI_SUCCESS;
1336   }
1337 #ifdef HAVE_TRACING
1338   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1339   TRACE_smpi_computing_in(rank);
1340 #endif
1341   smpi_bench_begin();
1342   return retval;
1343 }
1344
1345 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1346                 void *recvbuf, int *recvcounts, int *displs,
1347                 MPI_Datatype recvtype, int root, MPI_Comm comm)
1348 {
1349   int retval;
1350
1351   smpi_bench_end();
1352 #ifdef HAVE_TRACING
1353   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1354   TRACE_smpi_computing_out(rank);
1355   int root_traced = smpi_group_rank(smpi_comm_group(comm), root);
1356   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__);
1357 #endif
1358   if (comm == MPI_COMM_NULL) {
1359     retval = MPI_ERR_COMM;
1360   } else if (sendtype == MPI_DATATYPE_NULL
1361              || recvtype == MPI_DATATYPE_NULL) {
1362     retval = MPI_ERR_TYPE;
1363   } else if (recvcounts == NULL || displs == NULL) {
1364     retval = MPI_ERR_ARG;
1365   } else {
1366     smpi_mpi_gatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts,
1367                      displs, recvtype, root, comm);
1368     retval = MPI_SUCCESS;
1369   }
1370 #ifdef HAVE_TRACING
1371   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1372   TRACE_smpi_computing_in(rank);
1373 #endif
1374   smpi_bench_begin();
1375   return retval;
1376 }
1377
1378 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1379                   void *recvbuf, int recvcount, MPI_Datatype recvtype,
1380                   MPI_Comm comm)
1381 {
1382   int retval;
1383
1384   smpi_bench_end();
1385 #ifdef HAVE_TRACING
1386   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1387   TRACE_smpi_computing_out(rank);
1388   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1389 #endif
1390   if (comm == MPI_COMM_NULL) {
1391     retval = MPI_ERR_COMM;
1392   } else if (sendtype == MPI_DATATYPE_NULL
1393              || recvtype == MPI_DATATYPE_NULL) {
1394     retval = MPI_ERR_TYPE;
1395   } else {
1396     smpi_mpi_allgather(sendbuf, sendcount, sendtype, recvbuf, recvcount,
1397                        recvtype, comm);
1398     retval = MPI_SUCCESS;
1399   }
1400 #ifdef HAVE_TRACING
1401   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1402 #endif
1403   smpi_bench_begin();
1404   return retval;
1405 }
1406
1407 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1408                    void *recvbuf, int *recvcounts, int *displs,
1409                    MPI_Datatype recvtype, MPI_Comm comm)
1410 {
1411   int retval;
1412
1413   smpi_bench_end();
1414 #ifdef HAVE_TRACING
1415   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1416   TRACE_smpi_computing_out(rank);
1417   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1418 #endif
1419   if (comm == MPI_COMM_NULL) {
1420     retval = MPI_ERR_COMM;
1421   } else if (sendtype == MPI_DATATYPE_NULL
1422              || recvtype == MPI_DATATYPE_NULL) {
1423     retval = MPI_ERR_TYPE;
1424   } else if (recvcounts == NULL || displs == NULL) {
1425     retval = MPI_ERR_ARG;
1426   } else {
1427     smpi_mpi_allgatherv(sendbuf, sendcount, sendtype, recvbuf, recvcounts,
1428                         displs, recvtype, comm);
1429     retval = MPI_SUCCESS;
1430   }
1431 #ifdef HAVE_TRACING
1432   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1433   TRACE_smpi_computing_in(rank);
1434 #endif
1435   smpi_bench_begin();
1436   return retval;
1437 }
1438
1439 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1440                 void *recvbuf, int recvcount, MPI_Datatype recvtype,
1441                 int root, MPI_Comm comm)
1442 {
1443   int retval;
1444
1445   smpi_bench_end();
1446 #ifdef HAVE_TRACING
1447   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1448   TRACE_smpi_computing_out(rank);
1449   int root_traced = smpi_group_rank(smpi_comm_group(comm), root);
1450   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__);
1451 #endif
1452   if (comm == MPI_COMM_NULL) {
1453     retval = MPI_ERR_COMM;
1454   } else if (sendtype == MPI_DATATYPE_NULL
1455              || recvtype == MPI_DATATYPE_NULL) {
1456     retval = MPI_ERR_TYPE;
1457   } else {
1458     smpi_mpi_scatter(sendbuf, sendcount, sendtype, recvbuf, recvcount,
1459                      recvtype, root, comm);
1460     retval = MPI_SUCCESS;
1461   }
1462 #ifdef HAVE_TRACING
1463   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1464   TRACE_smpi_computing_in(rank);
1465 #endif
1466   smpi_bench_begin();
1467   return retval;
1468 }
1469
1470 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
1471                  MPI_Datatype sendtype, void *recvbuf, int recvcount,
1472                  MPI_Datatype recvtype, int root, MPI_Comm comm)
1473 {
1474   int retval;
1475
1476   smpi_bench_end();
1477 #ifdef HAVE_TRACING
1478   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1479   TRACE_smpi_computing_out(rank);
1480   int root_traced = smpi_group_rank(smpi_comm_group(comm), root);
1481   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__);
1482 #endif
1483   if (comm == MPI_COMM_NULL) {
1484     retval = MPI_ERR_COMM;
1485   } else if (sendtype == MPI_DATATYPE_NULL
1486              || recvtype == MPI_DATATYPE_NULL) {
1487     retval = MPI_ERR_TYPE;
1488   } else if (sendcounts == NULL || displs == NULL) {
1489     retval = MPI_ERR_ARG;
1490   } else {
1491     smpi_mpi_scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf,
1492                       recvcount, recvtype, root, comm);
1493     retval = MPI_SUCCESS;
1494   }
1495 #ifdef HAVE_TRACING
1496   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1497   TRACE_smpi_computing_in(rank);
1498 #endif
1499   smpi_bench_begin();
1500   return retval;
1501 }
1502
1503 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count,
1504                MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
1505 {
1506   int retval;
1507
1508   smpi_bench_end();
1509 #ifdef HAVE_TRACING
1510   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1511   TRACE_smpi_computing_out(rank);
1512   int root_traced = smpi_group_rank(smpi_comm_group(comm), root);
1513   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__);
1514 #endif
1515   if (comm == MPI_COMM_NULL) {
1516     retval = MPI_ERR_COMM;
1517   } else if (datatype == MPI_DATATYPE_NULL || op == MPI_OP_NULL) {
1518     retval = MPI_ERR_ARG;
1519   } else {
1520     smpi_mpi_reduce(sendbuf, recvbuf, count, datatype, op, root, comm);
1521     retval = MPI_SUCCESS;
1522   }
1523 #ifdef HAVE_TRACING
1524   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1525   TRACE_smpi_computing_in(rank);
1526 #endif
1527   smpi_bench_begin();
1528   return retval;
1529 }
1530
1531 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count,
1532                   MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
1533 {
1534   int retval;
1535
1536   smpi_bench_end();
1537 #ifdef HAVE_TRACING
1538   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1539   TRACE_smpi_computing_out(rank);
1540   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1541 #endif
1542   if (comm == MPI_COMM_NULL) {
1543     retval = MPI_ERR_COMM;
1544   } else if (datatype == MPI_DATATYPE_NULL) {
1545     retval = MPI_ERR_TYPE;
1546   } else if (op == MPI_OP_NULL) {
1547     retval = MPI_ERR_OP;
1548   } else {
1549     smpi_mpi_allreduce(sendbuf, recvbuf, count, datatype, op, comm);
1550     retval = MPI_SUCCESS;
1551   }
1552 #ifdef HAVE_TRACING
1553   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1554   TRACE_smpi_computing_in(rank);
1555 #endif
1556   smpi_bench_begin();
1557   return retval;
1558 }
1559
1560 int PMPI_Scan(void *sendbuf, void *recvbuf, int count,
1561              MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
1562 {
1563   int retval;
1564
1565   smpi_bench_end();
1566 #ifdef HAVE_TRACING
1567   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1568   TRACE_smpi_computing_out(rank);
1569   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1570 #endif
1571   if (comm == MPI_COMM_NULL) {
1572     retval = MPI_ERR_COMM;
1573   } else if (datatype == MPI_DATATYPE_NULL) {
1574     retval = MPI_ERR_TYPE;
1575   } else if (op == MPI_OP_NULL) {
1576     retval = MPI_ERR_OP;
1577   } else {
1578     smpi_mpi_scan(sendbuf, recvbuf, count, datatype, op, comm);
1579     retval = MPI_SUCCESS;
1580   }
1581 #ifdef HAVE_TRACING
1582   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1583   TRACE_smpi_computing_in(rank);
1584 #endif
1585   smpi_bench_begin();
1586   return retval;
1587 }
1588
1589 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts,
1590                        MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
1591 {
1592   int retval, i, size, count;
1593   int *displs;
1594   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1595
1596   smpi_bench_end();
1597 #ifdef HAVE_TRACING
1598   TRACE_smpi_computing_out(rank);
1599   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1600 #endif
1601   if (comm == MPI_COMM_NULL) {
1602     retval = MPI_ERR_COMM;
1603   } else if (datatype == MPI_DATATYPE_NULL) {
1604     retval = MPI_ERR_TYPE;
1605   } else if (op == MPI_OP_NULL) {
1606     retval = MPI_ERR_OP;
1607   } else if (recvcounts == NULL) {
1608     retval = MPI_ERR_ARG;
1609   } else {
1610     /* arbitrarily choose root as rank 0 */
1611     /* TODO: faster direct implementation ? */
1612     size = smpi_comm_size(comm);
1613     count = 0;
1614     displs = xbt_new(int, size);
1615     for (i = 0; i < size; i++) {
1616       count += recvcounts[i];
1617       displs[i] = 0;
1618     }
1619     smpi_mpi_reduce(sendbuf, recvbuf, count, datatype, op, 0, comm);
1620     smpi_mpi_scatterv(recvbuf, recvcounts, displs, datatype, recvbuf,
1621                       recvcounts[rank], datatype, 0, comm);
1622     xbt_free(displs);
1623     retval = MPI_SUCCESS;
1624   }
1625 #ifdef HAVE_TRACING
1626   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1627   TRACE_smpi_computing_in(rank);
1628 #endif
1629   smpi_bench_begin();
1630   return retval;
1631 }
1632
1633 int PMPI_Alltoall(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1634                  void *recvbuf, int recvcount, MPI_Datatype recvtype,
1635                  MPI_Comm comm)
1636 {
1637   int retval, size, sendsize;
1638
1639   smpi_bench_end();
1640 #ifdef HAVE_TRACING
1641   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1642   TRACE_smpi_computing_out(rank);
1643   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1644 #endif
1645   if (comm == MPI_COMM_NULL) {
1646     retval = MPI_ERR_COMM;
1647   } else if (sendtype == MPI_DATATYPE_NULL
1648              || recvtype == MPI_DATATYPE_NULL) {
1649     retval = MPI_ERR_TYPE;
1650   } else {
1651     size = smpi_comm_size(comm);
1652     sendsize = smpi_datatype_size(sendtype) * sendcount;
1653     if (sendsize < 200 && size > 12) {
1654       retval =
1655           smpi_coll_tuned_alltoall_bruck(sendbuf, sendcount, sendtype,
1656                                          recvbuf, recvcount, recvtype,
1657                                          comm);
1658     } else if (sendsize < 3000) {
1659       retval =
1660           smpi_coll_tuned_alltoall_basic_linear(sendbuf, sendcount,
1661                                                 sendtype, recvbuf,
1662                                                 recvcount, recvtype, comm);
1663     } else {
1664       retval =
1665           smpi_coll_tuned_alltoall_pairwise(sendbuf, sendcount, sendtype,
1666                                             recvbuf, recvcount, recvtype,
1667                                             comm);
1668     }
1669   }
1670 #ifdef HAVE_TRACING
1671   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1672   TRACE_smpi_computing_in(rank);
1673 #endif
1674   smpi_bench_begin();
1675   return retval;
1676 }
1677
1678 int PMPI_Alltoallv(void *sendbuf, int *sendcounts, int *senddisps,
1679                   MPI_Datatype sendtype, void *recvbuf, int *recvcounts,
1680                   int *recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
1681 {
1682   int retval;
1683
1684   smpi_bench_end();
1685 #ifdef HAVE_TRACING
1686   int rank = comm != MPI_COMM_NULL ? smpi_comm_rank(comm) : -1;
1687   TRACE_smpi_computing_out(rank);
1688   TRACE_smpi_collective_in(rank, -1, __FUNCTION__);
1689 #endif
1690   if (comm == MPI_COMM_NULL) {
1691     retval = MPI_ERR_COMM;
1692   } else if (sendtype == MPI_DATATYPE_NULL
1693              || recvtype == MPI_DATATYPE_NULL) {
1694     retval = MPI_ERR_TYPE;
1695   } else if (sendcounts == NULL || senddisps == NULL || recvcounts == NULL
1696              || recvdisps == NULL) {
1697     retval = MPI_ERR_ARG;
1698   } else {
1699     retval =
1700         smpi_coll_basic_alltoallv(sendbuf, sendcounts, senddisps, sendtype,
1701                                   recvbuf, recvcounts, recvdisps, recvtype,
1702                                   comm);
1703   }
1704 #ifdef HAVE_TRACING
1705   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1706   TRACE_smpi_computing_in(rank);
1707 #endif
1708   smpi_bench_begin();
1709   return retval;
1710 }
1711
1712
1713 int PMPI_Get_processor_name(char *name, int *resultlen)
1714 {
1715   int retval = MPI_SUCCESS;
1716
1717   smpi_bench_end();
1718   strncpy(name, SIMIX_host_get_name(SIMIX_host_self()),
1719           MPI_MAX_PROCESSOR_NAME - 1);
1720   *resultlen =
1721       strlen(name) >
1722       MPI_MAX_PROCESSOR_NAME ? MPI_MAX_PROCESSOR_NAME : strlen(name);
1723
1724   smpi_bench_begin();
1725   return retval;
1726 }
1727
1728 int PMPI_Get_count(MPI_Status * status, MPI_Datatype datatype, int *count)
1729 {
1730   int retval = MPI_SUCCESS;
1731   size_t size;
1732
1733   smpi_bench_end();
1734   if (status == NULL || count == NULL) {
1735     retval = MPI_ERR_ARG;
1736   } else if (datatype == MPI_DATATYPE_NULL) {
1737     retval = MPI_ERR_TYPE;
1738   } else {
1739     size = smpi_datatype_size(datatype);
1740     if (size == 0) {
1741       *count = 0;
1742     } else if (status->count % size != 0) {
1743       retval = MPI_UNDEFINED;
1744     } else {
1745       *count = smpi_mpi_get_count(status, datatype);
1746     }
1747   }
1748   smpi_bench_begin();
1749   return retval;
1750 }
1751
1752 /* The following calls are not yet implemented and will fail at runtime. */
1753 /* Once implemented, please move them above this notice. */
1754
1755 static int not_yet_implemented(void) {
1756    xbt_die("Not yet implemented");
1757    return MPI_ERR_OTHER;
1758 }
1759
1760 int PMPI_Pack_size(int incount, MPI_Datatype datatype, MPI_Comm comm, int* size) {
1761    return not_yet_implemented();
1762 }
1763
1764 int PMPI_Cart_coords(MPI_Comm comm, int rank, int maxdims, int* coords) {
1765    return not_yet_implemented();
1766 }
1767
1768 int PMPI_Cart_create(MPI_Comm comm_old, int ndims, int* dims, int* periods, int reorder, MPI_Comm* comm_cart) {
1769    return not_yet_implemented();
1770 }
1771
1772 int PMPI_Cart_get(MPI_Comm comm, int maxdims, int* dims, int* periods, int* coords) {
1773    return not_yet_implemented();
1774 }
1775
1776 int PMPI_Cart_map(MPI_Comm comm_old, int ndims, int* dims, int* periods, int* newrank) {
1777    return not_yet_implemented();
1778 }
1779
1780 int PMPI_Cart_rank(MPI_Comm comm, int* coords, int* rank) {
1781    return not_yet_implemented();
1782 }
1783
1784 int PMPI_Cart_shift(MPI_Comm comm, int direction, int displ, int* source, int* dest) {
1785    return not_yet_implemented();
1786 }
1787
1788 int PMPI_Cart_sub(MPI_Comm comm, int* remain_dims, MPI_Comm* comm_new) {
1789    return not_yet_implemented();
1790 }
1791
1792 int PMPI_Cartdim_get(MPI_Comm comm, int* ndims) {
1793    return not_yet_implemented();
1794 }
1795
1796 int PMPI_Graph_create(MPI_Comm comm_old, int nnodes, int* index, int* edges, int reorder, MPI_Comm* comm_graph) {
1797    return not_yet_implemented();
1798 }
1799
1800 int PMPI_Graph_get(MPI_Comm comm, int maxindex, int maxedges, int* index, int* edges) {
1801    return not_yet_implemented();
1802 }
1803
1804 int PMPI_Graph_map(MPI_Comm comm_old, int nnodes, int* index, int* edges, int* newrank) {
1805    return not_yet_implemented();
1806 }
1807
1808 int PMPI_Graph_neighbors(MPI_Comm comm, int rank, int maxneighbors, int* neighbors) {
1809    return not_yet_implemented();
1810 }
1811
1812 int PMPI_Graph_neighbors_count(MPI_Comm comm, int rank, int* nneighbors) {
1813    return not_yet_implemented();
1814 }
1815
1816 int PMPI_Graphdims_get(MPI_Comm comm, int* nnodes, int* nedges) {
1817    return not_yet_implemented();
1818 }
1819
1820 int PMPI_Topo_test(MPI_Comm comm, int* top_type) {
1821    return not_yet_implemented();
1822 }
1823
1824 int PMPI_Error_class(int errorcode, int* errorclass) {
1825    return not_yet_implemented();
1826 }
1827
1828 int PMPI_Errhandler_create(MPI_Handler_function* function, MPI_Errhandler* errhandler) {
1829    return not_yet_implemented();
1830 }
1831
1832 int PMPI_Errhandler_free(MPI_Errhandler* errhandler) {
1833    return not_yet_implemented();
1834 }
1835
1836 int PMPI_Errhandler_get(MPI_Comm comm, MPI_Errhandler* errhandler) {
1837    return not_yet_implemented();
1838 }
1839
1840 int PMPI_Error_string(int errorcode, char* string, int* resultlen) {
1841    return not_yet_implemented();
1842 }
1843
1844 int PMPI_Errhandler_set(MPI_Comm comm, MPI_Errhandler errhandler) {
1845    return not_yet_implemented();
1846 }
1847
1848 int PMPI_Type_contiguous(int count, MPI_Datatype old_type, MPI_Datatype* newtype) {
1849    return not_yet_implemented();
1850 }
1851
1852 int PMPI_Cancel(MPI_Request* request) {
1853    return not_yet_implemented();
1854 }
1855
1856 int PMPI_Buffer_attach(void* buffer, int size) {
1857    return not_yet_implemented();
1858 }
1859
1860 int PMPI_Buffer_detach(void* buffer, int* size) {
1861    return not_yet_implemented();
1862 }
1863
1864 int PMPI_Testsome(int incount, MPI_Request* requests, int* outcount, int* indices, MPI_Status* statuses) {
1865    return not_yet_implemented();
1866 }
1867
1868 int PMPI_Comm_test_inter(MPI_Comm comm, int* flag) {
1869    return not_yet_implemented();
1870 }
1871
1872 int PMPI_Unpack(void* inbuf, int insize, int* position, void* outbuf, int outcount, MPI_Datatype type, MPI_Comm comm) {
1873    return not_yet_implemented();
1874 }
1875
1876 int PMPI_Type_commit(MPI_Datatype* datatype) {
1877    return not_yet_implemented();
1878 }
1879
1880 int PMPI_Type_hindexed(int count, int* blocklens, MPI_Aint* indices, MPI_Datatype old_type, MPI_Datatype* newtype) {
1881    return not_yet_implemented();
1882 }
1883
1884 int PMPI_Type_hvector(int count, int blocklen, MPI_Aint stride, MPI_Datatype old_type, MPI_Datatype* newtype) {
1885    return not_yet_implemented();
1886 }
1887
1888 int PMPI_Type_indexed(int count, int* blocklens, int* indices, MPI_Datatype old_type, MPI_Datatype* newtype) {
1889    return not_yet_implemented();
1890 }
1891
1892 int PMPI_Type_struct(int count, int* blocklens, MPI_Aint* indices, MPI_Datatype* old_types, MPI_Datatype* newtype) {
1893    return not_yet_implemented();
1894 }
1895
1896 int PMPI_Type_vector(int count, int blocklen, int stride, MPI_Datatype old_type, MPI_Datatype* newtype) {
1897    return not_yet_implemented();
1898 }
1899
1900 int PMPI_Ssend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm) {
1901    return not_yet_implemented();
1902 }
1903
1904 int PMPI_Ssend_init(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
1905    return not_yet_implemented();
1906 }
1907
1908 int PMPI_Intercomm_create(MPI_Comm local_comm, int local_leader, MPI_Comm peer_comm, int remote_leader, int tag, MPI_Comm* comm_out) {
1909    return not_yet_implemented();
1910 }
1911
1912 int PMPI_Intercomm_merge(MPI_Comm comm, int high, MPI_Comm* comm_out) {
1913    return not_yet_implemented();
1914 }
1915
1916 int PMPI_Bsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm) {
1917    return not_yet_implemented();
1918 }
1919
1920 int PMPI_Bsend_init(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
1921    return not_yet_implemented();
1922 }
1923
1924 int PMPI_Ibsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
1925    return not_yet_implemented();
1926 }
1927
1928 int PMPI_Comm_remote_group(MPI_Comm comm, MPI_Group* group) {
1929    return not_yet_implemented();
1930 }
1931
1932 int PMPI_Comm_remote_size(MPI_Comm comm, int* size) {
1933    return not_yet_implemented();
1934 }
1935
1936 int PMPI_Issend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
1937    return not_yet_implemented();
1938 }
1939
1940 int PMPI_Probe(int source, int tag, MPI_Comm comm, MPI_Status* status) {
1941    return not_yet_implemented();
1942 }
1943
1944 int PMPI_Attr_delete(MPI_Comm comm, int keyval) {
1945    return not_yet_implemented();
1946 }
1947
1948 int PMPI_Attr_get(MPI_Comm comm, int keyval, void* attr_value, int* flag) {
1949    return not_yet_implemented();
1950 }
1951
1952 int PMPI_Attr_put(MPI_Comm comm, int keyval, void* attr_value) {
1953    return not_yet_implemented();
1954 }
1955
1956 int PMPI_Rsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm) {
1957    return not_yet_implemented();
1958 }
1959
1960 int PMPI_Rsend_init(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
1961    return not_yet_implemented();
1962 }
1963
1964 int PMPI_Irsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
1965    return not_yet_implemented();
1966 }
1967
1968 int PMPI_Keyval_create(MPI_Copy_function* copy_fn, MPI_Delete_function* delete_fn, int* keyval, void* extra_state) {
1969    return not_yet_implemented();
1970 }
1971
1972 int PMPI_Keyval_free(int* keyval) {
1973    return not_yet_implemented();
1974 }
1975
1976 int PMPI_Test_cancelled(MPI_Status* status, int* flag) {
1977    return not_yet_implemented();
1978 }
1979
1980 int PMPI_Pack(void* inbuf, int incount, MPI_Datatype type, void* outbuf, int outcount, int* position, MPI_Comm comm) {
1981    return not_yet_implemented();
1982 }
1983
1984 int PMPI_Testall(int count, MPI_Request* requests, int* flag, MPI_Status* statuses) {
1985    return not_yet_implemented();
1986 }
1987
1988 int PMPI_Get_elements(MPI_Status* status, MPI_Datatype datatype, int* elements) {
1989    return not_yet_implemented();
1990 }
1991
1992 int PMPI_Dims_create(int nnodes, int ndims, int* dims) {
1993    return not_yet_implemented();
1994 }
1995
1996 int PMPI_Iprobe(int source, int tag, MPI_Comm comm, int* flag, MPI_Status* status) {
1997    return not_yet_implemented();
1998 }
1999
2000 int PMPI_Initialized(int* flag) {
2001    return not_yet_implemented();
2002 }