Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
d5d9ecf58fde732a0ddda4ad2ebc21a4c5e50698
[simgrid.git] / src / smpi / smpi_pmpi.c
1 /* Copyright (c) 2007-2014. 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   smpi_process_mark_as_initialized();
31 #ifdef HAVE_TRACING
32   int rank = smpi_process_index();
33   TRACE_smpi_init(rank);
34   TRACE_smpi_computing_init(rank);
35   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
36   extra->type = TRACING_INIT;
37   TRACE_smpi_collective_in(rank, -1, __FUNCTION__, extra);
38   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
39 #endif
40   smpi_bench_begin();
41   return MPI_SUCCESS;
42 }
43
44 int PMPI_Finalize(void)
45 {
46   smpi_bench_end();
47 #ifdef HAVE_TRACING
48   int rank = smpi_process_index();
49   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
50   extra->type = TRACING_FINALIZE;
51   TRACE_smpi_collective_in(rank, -1, __FUNCTION__, extra);
52 #endif
53   smpi_process_finalize();
54 #ifdef HAVE_TRACING
55   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
56   TRACE_smpi_finalize(smpi_process_index());
57 #endif
58   smpi_process_destroy();
59   return MPI_SUCCESS;
60 }
61
62 int PMPI_Finalized(int* flag)
63 {
64   *flag=smpi_process_finalized();
65   return MPI_SUCCESS;
66 }
67
68 int PMPI_Get_version (int *version,int *subversion){
69   *version = MPI_VERSION;
70   *subversion= MPI_SUBVERSION;
71   return MPI_SUCCESS;
72 }
73
74 int PMPI_Get_library_version (char *version,int *len){
75   int retval = MPI_SUCCESS;
76   smpi_bench_end();
77   snprintf(version,MPI_MAX_LIBRARY_VERSION_STRING,"SMPI Version %d.%d. Copyright The Simgrid Team 2007-2014",SIMGRID_VERSION_MAJOR,
78           SIMGRID_VERSION_MINOR);
79   *len = strlen(version) > MPI_MAX_LIBRARY_VERSION_STRING ? MPI_MAX_LIBRARY_VERSION_STRING : strlen(version);
80   smpi_bench_begin();
81   return retval;
82 }
83
84 int PMPI_Init_thread(int *argc, char ***argv, int required, int *provided)
85 {
86   if (provided != NULL) {
87     *provided = MPI_THREAD_MULTIPLE;
88   }
89   return MPI_Init(argc, argv);
90 }
91
92 int PMPI_Query_thread(int *provided)
93 {
94   int retval = 0;
95
96   if (provided == NULL) {
97     retval = MPI_ERR_ARG;
98   } else {
99     *provided = MPI_THREAD_MULTIPLE;
100     retval = MPI_SUCCESS;
101   }
102   return retval;
103 }
104
105 int PMPI_Is_thread_main(int *flag)
106 {
107   int retval = 0;
108
109   if (flag == NULL) {
110     retval = MPI_ERR_ARG;
111   } else {
112     *flag = smpi_process_index() == 0;
113     retval = MPI_SUCCESS;
114   }
115   return retval;
116 }
117
118 int PMPI_Abort(MPI_Comm comm, int errorcode)
119 {
120   smpi_bench_end();
121   smpi_process_destroy();
122   // FIXME: should kill all processes in comm instead
123   simcall_process_kill(SIMIX_process_self());
124   return MPI_SUCCESS;
125 }
126
127 double PMPI_Wtime(void)
128 {
129   double time;
130   if (smpi_process_initialized() && !smpi_process_finalized() && !smpi_process_get_sampling()) {
131     smpi_bench_end();
132     time = SIMIX_get_clock();
133     smpi_bench_begin();
134   } else {
135     time = SIMIX_get_clock();
136   }
137   return time;
138 }
139
140 extern double sg_maxmin_precision;
141 double PMPI_Wtick(void)
142 {
143   return sg_maxmin_precision;
144 }
145
146 int PMPI_Address(void *location, MPI_Aint * address)
147 {
148   int retval = 0;
149
150   if (!address) {
151     retval = MPI_ERR_ARG;
152   } else {
153     *address = (MPI_Aint) location;
154     retval = MPI_SUCCESS;
155   }
156   return retval;
157 }
158
159 int PMPI_Get_address(void *location, MPI_Aint * address)
160 {
161   return PMPI_Address(location, address);
162 }
163
164 int PMPI_Type_free(MPI_Datatype * datatype)
165 {
166   int retval = 0;
167   /* Free a predefined datatype is an error according to the standard, and
168      should be checked for */
169   if (*datatype == MPI_DATATYPE_NULL) {
170     retval = MPI_ERR_ARG;
171   } else {
172     smpi_datatype_free(datatype);
173     retval = MPI_SUCCESS;
174   }
175   return retval;
176 }
177
178 int PMPI_Type_size(MPI_Datatype datatype, int *size)
179 {
180   int retval = 0;
181
182   if (datatype == MPI_DATATYPE_NULL) {
183     retval = MPI_ERR_TYPE;
184   } else if (size == NULL) {
185     retval = MPI_ERR_ARG;
186   } else {
187     *size = (int) smpi_datatype_size(datatype);
188     retval = MPI_SUCCESS;
189   }
190   return retval;
191 }
192
193 int PMPI_Type_get_extent(MPI_Datatype datatype, MPI_Aint * lb, MPI_Aint * extent)
194 {
195   int retval = 0;
196
197   if (datatype == MPI_DATATYPE_NULL) {
198     retval = MPI_ERR_TYPE;
199   } else if (lb == NULL || extent == NULL) {
200     retval = MPI_ERR_ARG;
201   } else {
202     retval = smpi_datatype_extent(datatype, lb, extent);
203   }
204   return retval;
205 }
206
207 int PMPI_Type_get_true_extent(MPI_Datatype datatype, MPI_Aint * lb, MPI_Aint * extent)
208 {
209   return PMPI_Type_get_extent(datatype, lb, extent);
210 }
211
212 int PMPI_Type_extent(MPI_Datatype datatype, MPI_Aint * extent)
213 {
214   int retval = 0;
215
216   if (datatype == MPI_DATATYPE_NULL) {
217     retval = MPI_ERR_TYPE;
218   } else if (extent == NULL) {
219     retval = MPI_ERR_ARG;
220   } else {
221     *extent = smpi_datatype_get_extent(datatype);
222     retval = MPI_SUCCESS;
223   }
224   return retval;
225 }
226
227 int PMPI_Type_lb(MPI_Datatype datatype, MPI_Aint * disp)
228 {
229   int retval = 0;
230
231   if (datatype == MPI_DATATYPE_NULL) {
232     retval = MPI_ERR_TYPE;
233   } else if (disp == NULL) {
234     retval = MPI_ERR_ARG;
235   } else {
236     *disp = smpi_datatype_lb(datatype);
237     retval = MPI_SUCCESS;
238   }
239   return retval;
240 }
241
242 int PMPI_Type_ub(MPI_Datatype datatype, MPI_Aint * disp)
243 {
244   int retval = 0;
245
246   if (datatype == MPI_DATATYPE_NULL) {
247     retval = MPI_ERR_TYPE;
248   } else if (disp == NULL) {
249     retval = MPI_ERR_ARG;
250   } else {
251     *disp = smpi_datatype_ub(datatype);
252     retval = MPI_SUCCESS;
253   }
254   return retval;
255 }
256
257 int PMPI_Op_create(MPI_User_function * function, int commute, MPI_Op * op)
258 {
259   int retval = 0;
260
261   if (function == NULL || op == NULL) {
262     retval = MPI_ERR_ARG;
263   } else {
264     *op = smpi_op_new(function, commute);
265     retval = MPI_SUCCESS;
266   }
267   return retval;
268 }
269
270 int PMPI_Op_free(MPI_Op * op)
271 {
272   int retval = 0;
273
274   if (op == NULL) {
275     retval = MPI_ERR_ARG;
276   } else if (*op == MPI_OP_NULL) {
277     retval = MPI_ERR_OP;
278   } else {
279     smpi_op_destroy(*op);
280     *op = MPI_OP_NULL;
281     retval = MPI_SUCCESS;
282   }
283   return retval;
284 }
285
286 int PMPI_Group_free(MPI_Group * group)
287 {
288   int retval = 0;
289
290   if (group == NULL) {
291     retval = MPI_ERR_ARG;
292   } else {
293     smpi_group_destroy(*group);
294     *group = MPI_GROUP_NULL;
295     retval = MPI_SUCCESS;
296   }
297   return retval;
298 }
299
300 int PMPI_Group_size(MPI_Group group, int *size)
301 {
302   int retval = 0;
303
304   if (group == MPI_GROUP_NULL) {
305     retval = MPI_ERR_GROUP;
306   } else if (size == NULL) {
307     retval = MPI_ERR_ARG;
308   } else {
309     *size = smpi_group_size(group);
310     retval = MPI_SUCCESS;
311   }
312   return retval;
313 }
314
315 int PMPI_Group_rank(MPI_Group group, int *rank)
316 {
317   int retval = 0;
318
319   if (group == MPI_GROUP_NULL) {
320     retval = MPI_ERR_GROUP;
321   } else if (rank == NULL) {
322     retval = MPI_ERR_ARG;
323   } else {
324     *rank = smpi_group_rank(group, smpi_process_index());
325     retval = MPI_SUCCESS;
326   }
327   return retval;
328 }
329
330 int PMPI_Group_translate_ranks(MPI_Group group1, int n, int *ranks1,
331                               MPI_Group group2, int *ranks2)
332 {
333   int retval, i, index;
334   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
335     retval = MPI_ERR_GROUP;
336   } else {
337     for (i = 0; i < n; i++) {
338       if(ranks1[i]==MPI_PROC_NULL){
339         ranks2[i]=MPI_PROC_NULL;
340       }else{
341         index = smpi_group_index(group1, ranks1[i]);
342         ranks2[i] = smpi_group_rank(group2, index);
343       }
344     }
345     retval = MPI_SUCCESS;
346   }
347   return retval;
348 }
349
350 int PMPI_Group_compare(MPI_Group group1, MPI_Group group2, int *result)
351 {
352   int retval = 0;
353
354   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
355     retval = MPI_ERR_GROUP;
356   } else if (result == NULL) {
357     retval = MPI_ERR_ARG;
358   } else {
359     *result = smpi_group_compare(group1, group2);
360     retval = MPI_SUCCESS;
361   }
362   return retval;
363 }
364
365 int PMPI_Group_union(MPI_Group group1, MPI_Group group2,
366                     MPI_Group * newgroup)
367 {
368   int retval, i, proc1, proc2, size, size2;
369
370   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
371     retval = MPI_ERR_GROUP;
372   } else if (newgroup == NULL) {
373     retval = MPI_ERR_ARG;
374   } else {
375     size = smpi_group_size(group1);
376     size2 = smpi_group_size(group2);
377     for (i = 0; i < size2; i++) {
378       proc2 = smpi_group_index(group2, i);
379       proc1 = smpi_group_rank(group1, proc2);
380       if (proc1 == MPI_UNDEFINED) {
381         size++;
382       }
383     }
384     if (size == 0) {
385       *newgroup = MPI_GROUP_EMPTY;
386     } else {
387       *newgroup = smpi_group_new(size);
388       size2 = smpi_group_size(group1);
389       for (i = 0; i < size2; i++) {
390         proc1 = smpi_group_index(group1, i);
391         smpi_group_set_mapping(*newgroup, proc1, i);
392       }
393       for (i = size2; i < size; i++) {
394         proc2 = smpi_group_index(group2, i - size2);
395         smpi_group_set_mapping(*newgroup, proc2, i);
396       }
397     }
398     retval = MPI_SUCCESS;
399   }
400   return retval;
401 }
402
403 int PMPI_Group_intersection(MPI_Group group1, MPI_Group group2,
404                            MPI_Group * newgroup)
405 {
406   int retval, i, proc1, proc2, size;
407
408   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
409     retval = MPI_ERR_GROUP;
410   } else if (newgroup == NULL) {
411     retval = MPI_ERR_ARG;
412   } else {
413     size = smpi_group_size(group2);
414     for (i = 0; i < size; i++) {
415       proc2 = smpi_group_index(group2, i);
416       proc1 = smpi_group_rank(group1, proc2);
417       if (proc1 == MPI_UNDEFINED) {
418         size--;
419       }
420     }
421     if (size == 0) {
422       *newgroup = MPI_GROUP_EMPTY;
423     } else {
424       *newgroup = smpi_group_new(size);
425       int j=0;
426       for (i = 0; i < smpi_group_size(group2); i++) {
427         proc2 = smpi_group_index(group2, i);
428         proc1 = smpi_group_rank(group1, proc2);
429         if (proc1 != MPI_UNDEFINED) {
430           smpi_group_set_mapping(*newgroup, proc2, j);
431           j++;
432         }
433       }
434     }
435     retval = MPI_SUCCESS;
436   }
437   return retval;
438 }
439
440 int PMPI_Group_difference(MPI_Group group1, MPI_Group group2, MPI_Group * newgroup)
441 {
442   int retval, i, proc1, proc2, size, size2;
443
444   if (group1 == MPI_GROUP_NULL || group2 == MPI_GROUP_NULL) {
445     retval = MPI_ERR_GROUP;
446   } else if (newgroup == NULL) {
447     retval = MPI_ERR_ARG;
448   } else {
449     size = size2 = smpi_group_size(group1);
450     for (i = 0; i < size2; i++) {
451       proc1 = smpi_group_index(group1, i);
452       proc2 = smpi_group_rank(group2, proc1);
453       if (proc2 != MPI_UNDEFINED) {
454         size--;
455       }
456     }
457     if (size == 0) {
458       *newgroup = MPI_GROUP_EMPTY;
459     } else {
460       *newgroup = smpi_group_new(size);
461       for (i = 0; i < size2; i++) {
462         proc1 = smpi_group_index(group1, i);
463         proc2 = smpi_group_rank(group2, proc1);
464         if (proc2 == MPI_UNDEFINED) {
465           smpi_group_set_mapping(*newgroup, proc1, i);
466         }
467       }
468     }
469     retval = MPI_SUCCESS;
470   }
471   return retval;
472 }
473
474 int PMPI_Group_incl(MPI_Group group, int n, int *ranks, MPI_Group * newgroup)
475 {
476   int retval, i, index;
477
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 if (n == smpi_group_size(group)) {
486       *newgroup = group;
487       if(group!= smpi_comm_group(MPI_COMM_WORLD)
488                 && group != MPI_GROUP_NULL
489                 && group != smpi_comm_group(MPI_COMM_SELF)
490                 && group != MPI_GROUP_EMPTY)
491       smpi_group_use(group);
492     } else {
493       *newgroup = smpi_group_new(n);
494       for (i = 0; i < n; i++) {
495         index = smpi_group_index(group, ranks[i]);
496         smpi_group_set_mapping(*newgroup, index, i);
497       }
498     }
499     retval = MPI_SUCCESS;
500   }
501   return retval;
502 }
503
504 int PMPI_Group_excl(MPI_Group group, int n, int *ranks, MPI_Group * newgroup)
505 {
506   int retval, i, j, newsize, oldsize, index;
507
508   if (group == MPI_GROUP_NULL) {
509     retval = MPI_ERR_GROUP;
510   } else if (newgroup == NULL) {
511     retval = MPI_ERR_ARG;
512   } else {
513     if (n == 0) {
514       *newgroup = group;
515       if(group!= smpi_comm_group(MPI_COMM_WORLD)
516                 && group != MPI_GROUP_NULL
517                 && group != smpi_comm_group(MPI_COMM_SELF)
518                 && group != MPI_GROUP_EMPTY)
519       smpi_group_use(group);
520     } else if (n == smpi_group_size(group)) {
521       *newgroup = MPI_GROUP_EMPTY;
522     } else {
523       oldsize=smpi_group_size(group);
524       newsize = oldsize - n;
525       *newgroup = smpi_group_new(newsize);
526
527       int* to_exclude=xbt_new0(int, smpi_group_size(group));
528       for(i=0; i<oldsize; i++)
529         to_exclude[i]=0;
530       for(i=0; i<n; i++)
531         to_exclude[ranks[i]]=1;
532
533       j=0;
534       for(i=0; i<oldsize; i++){
535         if(to_exclude[i]==0){
536           index = smpi_group_index(group, i);
537           smpi_group_set_mapping(*newgroup, index, j);
538           j++;
539         }
540       }
541
542       xbt_free(to_exclude);
543     }
544     retval = MPI_SUCCESS;
545   }
546   return retval;
547 }
548
549 int PMPI_Group_range_incl(MPI_Group group, int n, int ranges[][3],
550                          MPI_Group * newgroup)
551 {
552   int retval, i, j, rank, size, index;
553
554   if (group == MPI_GROUP_NULL) {
555     retval = MPI_ERR_GROUP;
556   } else if (newgroup == NULL) {
557     retval = MPI_ERR_ARG;
558   } else {
559     if (n == 0) {
560       *newgroup = MPI_GROUP_EMPTY;
561     } else {
562       size = 0;
563       for (i = 0; i < n; i++) {
564         for (rank = ranges[i][0];       /* First */
565              rank >= 0; /* Last */
566               ) {
567           size++;
568
569           rank += ranges[i][2]; /* Stride */
570           if (ranges[i][0]<ranges[i][1]){
571               if(rank > ranges[i][1])
572                 break;
573           }else{
574               if(rank < ranges[i][1])
575                 break;
576           }
577         }
578       }
579
580       *newgroup = smpi_group_new(size);
581       j = 0;
582       for (i = 0; i < n; i++) {
583         for (rank = ranges[i][0];     /* First */
584              rank >= 0; /* Last */
585              ) {
586           index = smpi_group_index(group, rank);
587           smpi_group_set_mapping(*newgroup, index, j);
588           j++;
589           rank += ranges[i][2]; /* Stride */
590           if (ranges[i][0]<ranges[i][1]){
591             if(rank > ranges[i][1])
592               break;
593           }else{
594             if(rank < ranges[i][1])
595               break;
596           }
597         }
598       }
599     }
600     retval = MPI_SUCCESS;
601   }
602   return retval;
603 }
604
605 int PMPI_Group_range_excl(MPI_Group group, int n, int ranges[][3],
606                          MPI_Group * newgroup)
607 {
608   int retval, i, rank, newrank,oldrank, size, index, add;
609
610   if (group == MPI_GROUP_NULL) {
611     retval = MPI_ERR_GROUP;
612   } else if (newgroup == NULL) {
613     retval = MPI_ERR_ARG;
614   } else {
615     if (n == 0) {
616       *newgroup = group;
617       if(group!= smpi_comm_group(MPI_COMM_WORLD)
618                 && group != MPI_GROUP_NULL
619                 && group != smpi_comm_group(MPI_COMM_SELF)
620                 && group != MPI_GROUP_EMPTY)
621       smpi_group_use(group);
622     } else {
623       size = smpi_group_size(group);
624       for (i = 0; i < n; i++) {
625         for (rank = ranges[i][0];       /* First */
626              rank >= 0; /* Last */
627               ) {
628           size--;
629
630           rank += ranges[i][2]; /* Stride */
631           if (ranges[i][0]<ranges[i][1]){
632               if(rank > ranges[i][1])
633                 break;
634           }else{
635               if(rank < ranges[i][1])
636                 break;
637           }
638         }
639       }
640       if (size == 0) {
641         *newgroup = MPI_GROUP_EMPTY;
642       } else {
643         *newgroup = smpi_group_new(size);
644         newrank=0;
645         oldrank=0;
646         while (newrank < size) {
647           add=1;
648           for (i = 0; i < n; i++) {
649             for (rank = ranges[i][0];rank >= 0;){
650               if(rank==oldrank){
651                   add=0;
652                   break;
653               }
654
655               rank += ranges[i][2]; /* Stride */
656
657               if (ranges[i][0]<ranges[i][1]){
658                   if(rank > ranges[i][1])
659                     break;
660               }else{
661                   if(rank < ranges[i][1])
662                     break;
663               }
664             }
665           }
666           if(add==1){
667             index = smpi_group_index(group, oldrank);
668             smpi_group_set_mapping(*newgroup, index, newrank);
669             newrank++;
670           }
671           oldrank++;
672         }
673       }
674     }
675
676     retval = MPI_SUCCESS;
677   }
678   return retval;
679 }
680
681 int PMPI_Comm_rank(MPI_Comm comm, int *rank)
682 {
683   int retval = 0;
684   if (comm == MPI_COMM_NULL) {
685     retval = MPI_ERR_COMM;
686   } else if (rank == NULL) {
687     retval = MPI_ERR_ARG;
688   } else {
689     *rank = smpi_comm_rank(comm);
690     retval = MPI_SUCCESS;
691   }
692   return retval;
693 }
694
695 int PMPI_Comm_size(MPI_Comm comm, int *size)
696 {
697   int retval = 0;
698   if (comm == MPI_COMM_NULL) {
699     retval = MPI_ERR_COMM;
700   } else if (size == NULL) {
701     retval = MPI_ERR_ARG;
702   } else {
703     *size = smpi_comm_size(comm);
704     retval = MPI_SUCCESS;
705   }
706   return retval;
707 }
708
709 int PMPI_Comm_get_name (MPI_Comm comm, char* name, int* len)
710 {
711   int retval = 0;
712
713   if (comm == MPI_COMM_NULL)  {
714     retval = MPI_ERR_COMM;
715   } else if (name == NULL || len == NULL)  {
716     retval = MPI_ERR_ARG;
717   } else {
718     smpi_comm_get_name(comm, name, len);
719     retval = MPI_SUCCESS;
720   }
721   return retval;
722 }
723
724 int PMPI_Comm_group(MPI_Comm comm, MPI_Group * group)
725 {
726   int retval = 0;
727
728   if (comm == MPI_COMM_NULL) {
729     retval = MPI_ERR_COMM;
730   } else if (group == NULL) {
731     retval = MPI_ERR_ARG;
732   } else {
733     *group = smpi_comm_group(comm);
734     if(*group!= smpi_comm_group(MPI_COMM_WORLD)
735               && *group != MPI_GROUP_NULL
736               && *group != smpi_comm_group(MPI_COMM_SELF)
737               && *group != MPI_GROUP_EMPTY)
738     smpi_group_use(*group);
739     retval = MPI_SUCCESS;
740   }
741   return retval;
742 }
743
744 int PMPI_Comm_compare(MPI_Comm comm1, MPI_Comm comm2, int *result)
745 {
746   int retval = 0;
747
748   if (comm1 == MPI_COMM_NULL || comm2 == MPI_COMM_NULL) {
749     retval = MPI_ERR_COMM;
750   } else if (result == NULL) {
751     retval = MPI_ERR_ARG;
752   } else {
753     if (comm1 == comm2) {       /* Same communicators means same groups */
754       *result = MPI_IDENT;
755     } else {
756       *result =
757           smpi_group_compare(smpi_comm_group(comm1),
758                              smpi_comm_group(comm2));
759       if (*result == MPI_IDENT) {
760         *result = MPI_CONGRUENT;
761       }
762     }
763     retval = MPI_SUCCESS;
764   }
765   return retval;
766 }
767
768 int PMPI_Comm_dup(MPI_Comm comm, MPI_Comm * newcomm)
769 {
770   int retval = 0;
771
772   if (comm == MPI_COMM_NULL) {
773     retval = MPI_ERR_COMM;
774   } else if (newcomm == NULL) {
775     retval = MPI_ERR_ARG;
776   } else {
777     *newcomm = smpi_comm_new(smpi_comm_group(comm), smpi_comm_topo(comm));
778     retval = MPI_SUCCESS;
779   }
780   return retval;
781 }
782
783 int PMPI_Comm_create(MPI_Comm comm, MPI_Group group, MPI_Comm * newcomm)
784 {
785   int retval = 0;
786
787   if (comm == MPI_COMM_NULL) {
788     retval = MPI_ERR_COMM;
789   } else if (group == MPI_GROUP_NULL) {
790     retval = MPI_ERR_GROUP;
791   } else if (newcomm == NULL) {
792     retval = MPI_ERR_ARG;
793   } else if(smpi_group_rank(group,smpi_process_index())==MPI_UNDEFINED){
794     *newcomm= MPI_COMM_NULL;
795     retval = MPI_SUCCESS;
796   }else{
797
798     *newcomm = smpi_comm_new(group, NULL);
799     retval = MPI_SUCCESS;
800   }
801   return retval;
802 }
803
804 int PMPI_Comm_free(MPI_Comm * comm)
805 {
806   int retval = 0;
807
808   if (comm == NULL) {
809     retval = MPI_ERR_ARG;
810   } else if (*comm == MPI_COMM_NULL) {
811     retval = MPI_ERR_COMM;
812   } else {
813     smpi_comm_destroy(*comm);
814     *comm = MPI_COMM_NULL;
815     retval = MPI_SUCCESS;
816   }
817   return retval;
818 }
819
820 int PMPI_Comm_disconnect(MPI_Comm * comm)
821 {
822   /* TODO: wait until all communication in comm are done */
823   int retval = 0;
824
825   if (comm == NULL) {
826     retval = MPI_ERR_ARG;
827   } else if (*comm == MPI_COMM_NULL) {
828     retval = MPI_ERR_COMM;
829   } else {
830     smpi_comm_destroy(*comm);
831     *comm = MPI_COMM_NULL;
832     retval = MPI_SUCCESS;
833   }
834   return retval;
835 }
836
837 int PMPI_Comm_split(MPI_Comm comm, int color, int key, MPI_Comm* comm_out)
838 {
839   int retval = 0;
840   smpi_bench_end();
841
842   if (comm_out == NULL) {
843     retval = MPI_ERR_ARG;
844   } else if (comm == MPI_COMM_NULL) {
845     retval = MPI_ERR_COMM;
846   } else {
847     *comm_out = smpi_comm_split(comm, color, key);
848     retval = MPI_SUCCESS;
849   }
850   smpi_bench_begin();
851
852   return retval;
853 }
854
855 int PMPI_Send_init(void *buf, int count, MPI_Datatype datatype, int dst,
856                    int tag, MPI_Comm comm, MPI_Request * request)
857 {
858   int retval = 0;
859
860   smpi_bench_end();
861   if (request == NULL) {
862       retval = MPI_ERR_ARG;
863   } else if (comm == MPI_COMM_NULL) {
864       retval = MPI_ERR_COMM;
865   } else if (!is_datatype_valid(datatype)) {
866       retval = MPI_ERR_TYPE;
867   } else if (dst == MPI_PROC_NULL) {
868       retval = MPI_SUCCESS;
869   } else {
870       *request = smpi_mpi_send_init(buf, count, datatype, dst, tag, comm);
871       retval = MPI_SUCCESS;
872   }
873   smpi_bench_begin();
874   if (retval != MPI_SUCCESS && request)
875     *request = MPI_REQUEST_NULL;
876   return retval;
877 }
878
879 int PMPI_Recv_init(void *buf, int count, MPI_Datatype datatype, int src,
880                    int tag, MPI_Comm comm, MPI_Request * request)
881 {
882   int retval = 0;
883
884   smpi_bench_end();
885   if (request == NULL) {
886     retval = MPI_ERR_ARG;
887   } else if (comm == MPI_COMM_NULL) {
888     retval = MPI_ERR_COMM;
889   } else if (!is_datatype_valid(datatype)) {
890       retval = MPI_ERR_TYPE;
891   } else if (src == MPI_PROC_NULL) {
892     retval = MPI_SUCCESS;
893   } else {
894     *request = smpi_mpi_recv_init(buf, count, datatype, src, tag, comm);
895     retval = MPI_SUCCESS;
896   }
897   smpi_bench_begin();
898   if (retval != MPI_SUCCESS && request)
899     *request = MPI_REQUEST_NULL;
900   return retval;
901 }
902
903 int PMPI_Ssend_init(void* buf, int count, MPI_Datatype datatype,
904                     int dst, int tag, MPI_Comm comm, MPI_Request* request)
905 {
906   int retval = 0;
907
908   smpi_bench_end();
909   if (request == NULL) {
910     retval = MPI_ERR_ARG;
911   } else if (comm == MPI_COMM_NULL) {
912     retval = MPI_ERR_COMM;
913   } else if (!is_datatype_valid(datatype)) {
914       retval = MPI_ERR_TYPE;
915   } else if (dst == MPI_PROC_NULL) {
916     retval = MPI_SUCCESS;
917   } else {
918     *request = smpi_mpi_ssend_init(buf, count, datatype, dst, tag, comm);
919     retval = MPI_SUCCESS;
920   }
921   smpi_bench_begin();
922   if (retval != MPI_SUCCESS && request)
923     *request = MPI_REQUEST_NULL;
924   return retval;
925 }
926
927 int PMPI_Start(MPI_Request * request)
928 {
929   int retval = 0;
930
931   smpi_bench_end();
932   if (request == NULL || *request == MPI_REQUEST_NULL) {
933     retval = MPI_ERR_REQUEST;
934   } else {
935     smpi_mpi_start(*request);
936     retval = MPI_SUCCESS;
937   }
938   smpi_bench_begin();
939   return retval;
940 }
941
942 int PMPI_Startall(int count, MPI_Request * requests)
943 {
944   int retval;
945   int i = 0;
946   smpi_bench_end();
947   if (requests == NULL) {
948     retval = MPI_ERR_ARG;
949   } else {
950     retval = MPI_SUCCESS;
951     for (i = 0 ;  i < count ; i++) {
952       if(requests[i] == MPI_REQUEST_NULL) {
953         retval = MPI_ERR_REQUEST;
954       }
955     }
956     if(retval != MPI_ERR_REQUEST) {
957       smpi_mpi_startall(count, requests);
958     }
959   }
960   smpi_bench_begin();
961   return retval;
962 }
963
964 int PMPI_Request_free(MPI_Request * request)
965 {
966   int retval = 0;
967
968   smpi_bench_end();
969   if (*request == MPI_REQUEST_NULL) {
970     retval = MPI_ERR_ARG;
971   } else {
972     smpi_mpi_request_free(request);
973     retval = MPI_SUCCESS;
974   }
975   smpi_bench_begin();
976   return retval;
977 }
978
979 int PMPI_Irecv(void *buf, int count, MPI_Datatype datatype, int src,
980                int tag, MPI_Comm comm, MPI_Request * request)
981 {
982   int retval = 0;
983
984   smpi_bench_end();
985
986   if (request == NULL) {
987     retval = MPI_ERR_ARG;
988   } else if (comm == MPI_COMM_NULL) {
989     retval = MPI_ERR_COMM;
990   } else if (src == MPI_PROC_NULL) {
991     *request = MPI_REQUEST_NULL;
992     retval = MPI_SUCCESS;
993   } else if (src!=MPI_ANY_SOURCE && (src >= smpi_group_size(smpi_comm_group(comm)) || src <0)){
994     retval = MPI_ERR_COMM;
995   } else if (count < 0) {
996     retval = MPI_ERR_COUNT;
997   } else if (buf==NULL && count > 0) {
998     retval = MPI_ERR_COUNT;
999   } else if (!is_datatype_valid(datatype)) {
1000       retval = MPI_ERR_TYPE;
1001   } else if(tag<0 && tag !=  MPI_ANY_TAG){
1002     retval = MPI_ERR_TAG;
1003   } else {
1004
1005 #ifdef HAVE_TRACING
1006     int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1007     int src_traced = smpi_group_index(smpi_comm_group(comm), src);
1008
1009     instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1010     extra->type = TRACING_IRECV;
1011     extra->send_size = count;
1012     extra->src = src_traced;
1013     extra->dst = rank;
1014     extra->datatype1 = encode_datatype(datatype);
1015     TRACE_smpi_ptp_in(rank, src_traced, rank, __FUNCTION__, extra);
1016 #endif
1017
1018     *request = smpi_mpi_irecv(buf, count, datatype, src, tag, comm);
1019     retval = MPI_SUCCESS;
1020
1021 #ifdef HAVE_TRACING
1022     TRACE_smpi_ptp_out(rank, src_traced, rank, __FUNCTION__);
1023     (*request)->recv = 1;
1024 #endif
1025   }
1026
1027   smpi_bench_begin();
1028   if (retval != MPI_SUCCESS && request)
1029     *request = MPI_REQUEST_NULL;
1030   return retval;
1031 }
1032
1033
1034 int PMPI_Isend(void *buf, int count, MPI_Datatype datatype, int dst,
1035                int tag, MPI_Comm comm, MPI_Request * request)
1036 {
1037   int retval = 0;
1038
1039   smpi_bench_end();
1040   if (request == NULL) {
1041     retval = MPI_ERR_ARG;
1042   } else if (comm == MPI_COMM_NULL) {
1043     retval = MPI_ERR_COMM;
1044   } else if (dst == MPI_PROC_NULL) {
1045     *request = MPI_REQUEST_NULL;
1046     retval = MPI_SUCCESS;
1047   } else if (dst >= smpi_group_size(smpi_comm_group(comm)) || dst <0){
1048     retval = MPI_ERR_RANK;
1049   } else if (count < 0) {
1050     retval = MPI_ERR_COUNT;
1051   } else if (buf==NULL && count > 0) {
1052     retval = MPI_ERR_COUNT;
1053   } else if (!is_datatype_valid(datatype)) {
1054       retval = MPI_ERR_TYPE;
1055   } else if(tag<0 && tag !=  MPI_ANY_TAG){
1056     retval = MPI_ERR_TAG;
1057   } else {
1058
1059 #ifdef HAVE_TRACING
1060     int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1061     int dst_traced = smpi_group_index(smpi_comm_group(comm), dst);
1062
1063     instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1064     extra->type = TRACING_ISEND;
1065     extra->send_size = count;
1066     extra->src = rank;
1067     extra->dst = dst_traced;
1068     extra->datatype1 = encode_datatype(datatype);
1069     TRACE_smpi_ptp_in(rank, rank, dst_traced, __FUNCTION__, extra);
1070     TRACE_smpi_send(rank, rank, dst_traced, count*smpi_datatype_size(datatype));
1071 #endif
1072
1073     *request = smpi_mpi_isend(buf, count, datatype, dst, tag, comm);
1074     retval = MPI_SUCCESS;
1075
1076 #ifdef HAVE_TRACING
1077     TRACE_smpi_ptp_out(rank, rank, dst_traced, __FUNCTION__);
1078     (*request)->send = 1;
1079 #endif
1080   }
1081
1082   smpi_bench_begin();
1083   if (retval != MPI_SUCCESS && request)
1084     *request = MPI_REQUEST_NULL;
1085   return retval;
1086 }
1087
1088 int PMPI_Issend(void* buf, int count, MPI_Datatype datatype,
1089                 int dst, int tag, MPI_Comm comm, MPI_Request* request)
1090 {
1091   int retval = 0;
1092
1093   smpi_bench_end();
1094   if (request == NULL) {
1095     retval = MPI_ERR_ARG;
1096   } else if (comm == MPI_COMM_NULL) {
1097     retval = MPI_ERR_COMM;
1098   } else if (dst == MPI_PROC_NULL) {
1099     *request = MPI_REQUEST_NULL;
1100     retval = MPI_SUCCESS;
1101   } else if (dst >= smpi_group_size(smpi_comm_group(comm)) || dst <0){
1102     retval = MPI_ERR_RANK;
1103   } else if (count < 0) {
1104     retval = MPI_ERR_COUNT;
1105   } else if (buf==NULL && count > 0) {
1106     retval = MPI_ERR_COUNT;
1107   } else if (!is_datatype_valid(datatype)) {
1108       retval = MPI_ERR_TYPE;
1109   } else if(tag<0 && tag !=  MPI_ANY_TAG){
1110     retval = MPI_ERR_TAG;
1111   } else {
1112
1113 #ifdef HAVE_TRACING
1114     int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1115     int dst_traced = smpi_group_index(smpi_comm_group(comm), dst);
1116     instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1117     extra->type = TRACING_ISSEND;
1118     extra->send_size = count;
1119     extra->src = rank;
1120     extra->dst = dst_traced;
1121     extra->datatype1 = encode_datatype(datatype);
1122     TRACE_smpi_ptp_in(rank, rank, dst_traced, __FUNCTION__, extra);
1123     TRACE_smpi_send(rank, rank, dst_traced, count*smpi_datatype_size(datatype));
1124 #endif
1125
1126     *request = smpi_mpi_issend(buf, count, datatype, dst, tag, comm);
1127     retval = MPI_SUCCESS;
1128
1129 #ifdef HAVE_TRACING
1130     TRACE_smpi_ptp_out(rank, rank, dst_traced, __FUNCTION__);
1131     (*request)->send = 1;
1132 #endif
1133   }
1134
1135   smpi_bench_begin();
1136   if (retval != MPI_SUCCESS && request)
1137     *request = MPI_REQUEST_NULL;
1138   return retval;
1139 }
1140
1141 int PMPI_Recv(void *buf, int count, MPI_Datatype datatype, int src, int tag,
1142              MPI_Comm comm, MPI_Status * status)
1143 {
1144   int retval = 0;
1145
1146   smpi_bench_end();
1147   if (comm == MPI_COMM_NULL) {
1148     retval = MPI_ERR_COMM;
1149   } else if (src == MPI_PROC_NULL) {
1150     smpi_empty_status(status);
1151     status->MPI_SOURCE = MPI_PROC_NULL;
1152     retval = MPI_SUCCESS;
1153   } else if (src!=MPI_ANY_SOURCE && (src >= smpi_group_size(smpi_comm_group(comm)) || src <0)){
1154     retval = MPI_ERR_RANK;
1155   } else if (count < 0) {
1156     retval = MPI_ERR_COUNT;
1157   } else if (buf==NULL && count > 0) {
1158     retval = MPI_ERR_COUNT;
1159   } else if (!is_datatype_valid(datatype)) {
1160       retval = MPI_ERR_TYPE;
1161   } else if(tag<0 && tag !=  MPI_ANY_TAG){
1162     retval = MPI_ERR_TAG;
1163   } else {
1164 #ifdef HAVE_TRACING
1165   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1166   int src_traced = smpi_group_index(smpi_comm_group(comm), src);
1167   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1168   extra->type = TRACING_RECV;
1169   extra->send_size = count;
1170   extra->src = src_traced;
1171   extra->dst = rank;
1172   extra->datatype1 = encode_datatype(datatype);
1173   TRACE_smpi_ptp_in(rank, src_traced, rank, __FUNCTION__, extra);
1174 #endif
1175
1176     smpi_mpi_recv(buf, count, datatype, src, tag, comm, status);
1177     retval = MPI_SUCCESS;
1178
1179 #ifdef HAVE_TRACING
1180   //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
1181   if(status!=MPI_STATUS_IGNORE){
1182     src_traced = smpi_group_index(smpi_comm_group(comm), status->MPI_SOURCE);
1183     TRACE_smpi_recv(rank, src_traced, rank);
1184   }
1185   TRACE_smpi_ptp_out(rank, src_traced, rank, __FUNCTION__);
1186 #endif
1187   }
1188
1189   smpi_bench_begin();
1190   return retval;
1191 }
1192
1193 int PMPI_Send(void *buf, int count, MPI_Datatype datatype, int dst, int tag,
1194              MPI_Comm comm)
1195 {
1196   int retval = 0;
1197
1198   smpi_bench_end();
1199
1200   if (comm == MPI_COMM_NULL) {
1201     retval = MPI_ERR_COMM;
1202   } else if (dst == MPI_PROC_NULL) {
1203     retval = MPI_SUCCESS;
1204   } else if (dst >= smpi_group_size(smpi_comm_group(comm)) || dst <0){
1205     retval = MPI_ERR_RANK;
1206   } else if (count < 0) {
1207     retval = MPI_ERR_COUNT;
1208   } else if (buf==NULL && count > 0) {
1209     retval = MPI_ERR_COUNT;
1210   } else if (!is_datatype_valid(datatype)) {
1211       retval = MPI_ERR_TYPE;
1212   } else if(tag<0 && tag !=  MPI_ANY_TAG){
1213     retval = MPI_ERR_TAG;
1214   } else {
1215
1216 #ifdef HAVE_TRACING
1217   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1218   int dst_traced = smpi_group_index(smpi_comm_group(comm), dst);
1219   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1220   extra->type = TRACING_SEND;
1221   extra->send_size = count;
1222   extra->src = rank;
1223   extra->dst = dst_traced;
1224   extra->datatype1 = encode_datatype(datatype);
1225   TRACE_smpi_ptp_in(rank, rank, dst_traced, __FUNCTION__, extra);
1226   TRACE_smpi_send(rank, rank, dst_traced,count*smpi_datatype_size(datatype));
1227 #endif
1228
1229     smpi_mpi_send(buf, count, datatype, dst, tag, comm);
1230     retval = MPI_SUCCESS;
1231
1232 #ifdef HAVE_TRACING
1233   TRACE_smpi_ptp_out(rank, rank, dst_traced, __FUNCTION__);
1234 #endif
1235   }
1236
1237   smpi_bench_begin();
1238   return retval;
1239 }
1240
1241
1242
1243 int PMPI_Ssend(void* buf, int count, MPI_Datatype datatype, int dst, int tag, MPI_Comm comm) {
1244   int retval = 0;
1245
1246    smpi_bench_end();
1247
1248    if (comm == MPI_COMM_NULL) {
1249      retval = MPI_ERR_COMM;
1250    } else if (dst == MPI_PROC_NULL) {
1251      retval = MPI_SUCCESS;
1252    } else if (dst >= smpi_group_size(smpi_comm_group(comm)) || dst <0){
1253      retval = MPI_ERR_RANK;
1254    } else if (count < 0) {
1255      retval = MPI_ERR_COUNT;
1256    } else if (buf==NULL && count > 0) {
1257      retval = MPI_ERR_COUNT;
1258    } else if (!is_datatype_valid(datatype)){
1259      retval = MPI_ERR_TYPE;
1260    } else if(tag<0 && tag !=  MPI_ANY_TAG){
1261      retval = MPI_ERR_TAG;
1262    } else {
1263
1264  #ifdef HAVE_TRACING
1265    int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1266    int dst_traced = smpi_group_index(smpi_comm_group(comm), dst);
1267    instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1268    extra->type = TRACING_SSEND;
1269    extra->send_size = count;
1270    extra->src = rank;
1271    extra->dst = dst_traced;
1272    extra->datatype1 = encode_datatype(datatype);
1273    TRACE_smpi_ptp_in(rank, rank, dst_traced, __FUNCTION__, extra);   TRACE_smpi_send(rank, rank, dst_traced,count*smpi_datatype_size(datatype));
1274  #endif
1275
1276      smpi_mpi_ssend(buf, count, datatype, dst, tag, comm);
1277      retval = MPI_SUCCESS;
1278
1279  #ifdef HAVE_TRACING
1280    TRACE_smpi_ptp_out(rank, rank, dst_traced, __FUNCTION__);
1281  #endif
1282    }
1283
1284    smpi_bench_begin();
1285    return retval;}
1286
1287
1288 int PMPI_Sendrecv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1289                  int dst, int sendtag, void *recvbuf, int recvcount,
1290                  MPI_Datatype recvtype, int src, int recvtag,
1291                  MPI_Comm comm, MPI_Status * status)
1292 {
1293   int retval = 0;
1294
1295   smpi_bench_end();
1296
1297   if (comm == MPI_COMM_NULL) {
1298     retval = MPI_ERR_COMM;
1299   } else if (!is_datatype_valid(sendtype)
1300              || !is_datatype_valid(recvtype)) {
1301     retval = MPI_ERR_TYPE;
1302   } else if (src == MPI_PROC_NULL || dst == MPI_PROC_NULL) {
1303       smpi_empty_status(status);
1304       status->MPI_SOURCE = MPI_PROC_NULL;
1305       retval = MPI_SUCCESS;
1306   }else if (dst >= smpi_group_size(smpi_comm_group(comm)) || dst <0 ||
1307       (src!=MPI_ANY_SOURCE && (src >= smpi_group_size(smpi_comm_group(comm)) || src <0))){
1308     retval = MPI_ERR_RANK;
1309   } else if (sendcount < 0 || recvcount<0) {
1310       retval = MPI_ERR_COUNT;
1311   } else if ((sendbuf==NULL && sendcount > 0)||(recvbuf==NULL && recvcount>0)) {
1312     retval = MPI_ERR_COUNT;
1313   } else if((sendtag<0 && sendtag !=  MPI_ANY_TAG)||(recvtag<0 && recvtag != MPI_ANY_TAG)){
1314     retval = MPI_ERR_TAG;
1315   } else {
1316
1317 #ifdef HAVE_TRACING
1318   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1319   int dst_traced = smpi_group_index(smpi_comm_group(comm), dst);
1320   int src_traced = smpi_group_index(smpi_comm_group(comm), src);
1321   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1322   extra->type = TRACING_SENDRECV;
1323   extra->send_size = sendcount;
1324   extra->recv_size = recvcount;
1325   extra->src = src_traced;
1326   extra->dst = dst_traced;
1327   extra->datatype1 = encode_datatype(sendtype);
1328   extra->datatype2 = encode_datatype(recvtype);
1329
1330   TRACE_smpi_ptp_in(rank, src_traced, dst_traced, __FUNCTION__, extra);
1331   TRACE_smpi_send(rank, rank, dst_traced,sendcount*smpi_datatype_size(sendtype));
1332 #endif
1333
1334
1335     smpi_mpi_sendrecv(sendbuf, sendcount, sendtype, dst, sendtag, recvbuf,
1336                       recvcount, recvtype, src, recvtag, comm, status);
1337     retval = MPI_SUCCESS;
1338
1339 #ifdef HAVE_TRACING
1340   TRACE_smpi_ptp_out(rank, src_traced, dst_traced, __FUNCTION__);
1341   TRACE_smpi_recv(rank, src_traced, rank);
1342 #endif
1343
1344   }
1345
1346   smpi_bench_begin();
1347   return retval;
1348 }
1349
1350 int PMPI_Sendrecv_replace(void *buf, int count, MPI_Datatype datatype,
1351                          int dst, int sendtag, int src, int recvtag,
1352                          MPI_Comm comm, MPI_Status * status)
1353 {
1354   //TODO: suboptimal implementation
1355   void *recvbuf;
1356   int retval = 0;
1357   if (!is_datatype_valid(datatype)) {
1358       retval = MPI_ERR_TYPE;
1359   } else if (count < 0) {
1360       retval = MPI_ERR_COUNT;
1361   } else {
1362     int size = smpi_datatype_get_extent(datatype) * count;
1363     recvbuf = xbt_new0(char, size);
1364     retval =
1365         MPI_Sendrecv(buf, count, datatype, dst, sendtag, recvbuf, count,
1366                      datatype, src, recvtag, comm, status);
1367     if(retval==MPI_SUCCESS){
1368         smpi_datatype_copy(recvbuf, count, datatype, buf, count, datatype);
1369     }
1370     xbt_free(recvbuf);
1371
1372   }
1373   return retval;
1374 }
1375
1376 int PMPI_Test(MPI_Request * request, int *flag, MPI_Status * status)
1377 {
1378   int retval = 0;
1379   smpi_bench_end();
1380   if (request == NULL || flag == NULL) {
1381     retval = MPI_ERR_ARG;
1382   } else if (*request == MPI_REQUEST_NULL) {
1383     *flag= TRUE;
1384     smpi_empty_status(status);
1385     retval = MPI_ERR_REQUEST;
1386   } else {
1387 #ifdef HAVE_TRACING
1388     int rank = request && (*request)->comm != MPI_COMM_NULL
1389       ? smpi_process_index()
1390       : -1;
1391
1392     instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1393     extra->type = TRACING_TEST;
1394     TRACE_smpi_testing_in(rank, extra);
1395 #endif
1396     *flag = smpi_mpi_test(request, status);
1397 #ifdef HAVE_TRACING
1398     TRACE_smpi_testing_out(rank);
1399 #endif
1400     retval = MPI_SUCCESS;
1401   }
1402   smpi_bench_begin();
1403   return retval;
1404 }
1405
1406 int PMPI_Testany(int count, MPI_Request requests[], int *index, int *flag,
1407                 MPI_Status * status)
1408 {
1409   int retval = 0;
1410
1411   smpi_bench_end();
1412   if (index == NULL || flag == NULL) {
1413     retval = MPI_ERR_ARG;
1414   } else {
1415     *flag = smpi_mpi_testany(count, requests, index, status);
1416     retval = MPI_SUCCESS;
1417   }
1418   smpi_bench_begin();
1419   return retval;
1420 }
1421
1422 int PMPI_Testall(int count, MPI_Request* requests, int* flag, MPI_Status* statuses)
1423 {
1424   int retval = 0;
1425
1426   smpi_bench_end();
1427   if (flag == NULL) {
1428     retval = MPI_ERR_ARG;
1429   } else {
1430     *flag = smpi_mpi_testall(count, requests, statuses);
1431     retval = MPI_SUCCESS;
1432   }
1433   smpi_bench_begin();
1434   return retval;
1435 }
1436
1437 int PMPI_Probe(int source, int tag, MPI_Comm comm, MPI_Status* status) {
1438   int retval = 0;
1439   smpi_bench_end();
1440
1441   if (status == NULL) {
1442     retval = MPI_ERR_ARG;
1443   } else if (comm == MPI_COMM_NULL) {
1444     retval = MPI_ERR_COMM;
1445   } else if (source == MPI_PROC_NULL) {
1446     smpi_empty_status(status);
1447     status->MPI_SOURCE = MPI_PROC_NULL;
1448     retval = MPI_SUCCESS;
1449   } else {
1450     smpi_mpi_probe(source, tag, comm, status);
1451     retval = MPI_SUCCESS;
1452   }
1453   smpi_bench_begin();
1454   return retval;
1455 }
1456
1457
1458 int PMPI_Iprobe(int source, int tag, MPI_Comm comm, int* flag, MPI_Status* status) {
1459   int retval = 0;
1460   smpi_bench_end();
1461
1462   if (flag == NULL) {
1463     retval = MPI_ERR_ARG;
1464   } else if (status == NULL) {
1465     retval = MPI_ERR_ARG;
1466   } else if (comm == MPI_COMM_NULL) {
1467     retval = MPI_ERR_COMM;
1468   } else if (source == MPI_PROC_NULL) {
1469     *flag=TRUE;
1470     smpi_empty_status(status);
1471     status->MPI_SOURCE = MPI_PROC_NULL;
1472     retval = MPI_SUCCESS;
1473   } else {
1474     smpi_mpi_iprobe(source, tag, comm, flag, status);
1475     retval = MPI_SUCCESS;
1476   }
1477   smpi_bench_begin();
1478   return retval;
1479 }
1480
1481 int PMPI_Wait(MPI_Request * request, MPI_Status * status)
1482 {
1483   int retval = 0;
1484
1485   smpi_bench_end();
1486
1487   smpi_empty_status(status);
1488
1489   if (request == NULL) {
1490     retval = MPI_ERR_ARG;
1491   } else if (*request == MPI_REQUEST_NULL) {
1492     retval = MPI_ERR_REQUEST;
1493   } else {
1494
1495 #ifdef HAVE_TRACING
1496     int rank = request && (*request)->comm != MPI_COMM_NULL
1497       ? smpi_process_index()
1498       : -1;
1499
1500     int src_traced = (*request)->src;
1501     int dst_traced = (*request)->dst;
1502     MPI_Comm comm = (*request)->comm;
1503     int is_wait_for_receive = (*request)->recv;
1504     instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1505     extra->type = TRACING_WAIT;
1506     TRACE_smpi_ptp_in(rank, src_traced, dst_traced, __FUNCTION__, extra);
1507 #endif
1508
1509     smpi_mpi_wait(request, status);
1510     retval = MPI_SUCCESS;
1511
1512 #ifdef HAVE_TRACING
1513     //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
1514     TRACE_smpi_ptp_out(rank, src_traced, dst_traced, __FUNCTION__);
1515     if (is_wait_for_receive) {
1516       if(src_traced==MPI_ANY_SOURCE)
1517         src_traced = (status!=MPI_STATUS_IGNORE) ?
1518           smpi_group_rank(smpi_comm_group(comm), status->MPI_SOURCE) :
1519           src_traced;
1520       TRACE_smpi_recv(rank, src_traced, dst_traced);
1521     }
1522 #endif
1523
1524   }
1525
1526   smpi_bench_begin();
1527   return retval;
1528 }
1529
1530 int PMPI_Waitany(int count, MPI_Request requests[], int *index, MPI_Status * status)
1531 {
1532   if (index == NULL)
1533     return MPI_ERR_ARG;
1534
1535   smpi_bench_end();
1536 #ifdef HAVE_TRACING
1537   //save requests information for tracing
1538   int i;
1539   int *srcs = xbt_new0(int, count);
1540   int *dsts = xbt_new0(int, count);
1541   int *recvs = xbt_new0(int, count);
1542   MPI_Comm *comms = xbt_new0(MPI_Comm, count);
1543
1544   for (i = 0; i < count; i++) {
1545     MPI_Request req = requests[i];      //already received requests are no longer valid
1546     if (req) {
1547       srcs[i] = req->src;
1548       dsts[i] = req->dst;
1549       recvs[i] = req->recv;
1550       comms[i] = req->comm;
1551     }
1552   }
1553   int rank_traced = smpi_process_index();
1554   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1555   extra->type = TRACING_WAITANY;
1556   extra->send_size=count;
1557   TRACE_smpi_ptp_in(rank_traced, -1, -1, __FUNCTION__,extra);
1558
1559 #endif
1560   *index = smpi_mpi_waitany(count, requests, status);
1561 #ifdef HAVE_TRACING
1562   if(*index!=MPI_UNDEFINED){
1563     int src_traced = srcs[*index];
1564     //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
1565     int dst_traced = dsts[*index];
1566     int is_wait_for_receive = recvs[*index];
1567     if (is_wait_for_receive) {
1568       if(srcs[*index]==MPI_ANY_SOURCE)
1569         src_traced = (status!=MPI_STATUSES_IGNORE) ?
1570                       smpi_group_rank(smpi_comm_group(comms[*index]), status->MPI_SOURCE) :
1571                       srcs[*index];
1572       TRACE_smpi_recv(rank_traced, src_traced, dst_traced);
1573     }
1574     TRACE_smpi_ptp_out(rank_traced, src_traced, dst_traced, __FUNCTION__);
1575     xbt_free(srcs);
1576     xbt_free(dsts);
1577     xbt_free(recvs);
1578     xbt_free(comms);
1579
1580   }
1581 #endif
1582   smpi_bench_begin();
1583   return MPI_SUCCESS;
1584 }
1585
1586 int PMPI_Waitall(int count, MPI_Request requests[], MPI_Status status[])
1587 {
1588
1589   smpi_bench_end();
1590 #ifdef HAVE_TRACING
1591   //save information from requests
1592   int i;
1593   int *srcs = xbt_new0(int, count);
1594   int *dsts = xbt_new0(int, count);
1595   int *recvs = xbt_new0(int, count);
1596   int *valid = xbt_new0(int, count);
1597   MPI_Comm *comms = xbt_new0(MPI_Comm, count);
1598
1599   //int valid_count = 0;
1600   for (i = 0; i < count; i++) {
1601     MPI_Request req = requests[i];
1602     if(req!=MPI_REQUEST_NULL){
1603       srcs[i] = req->src;
1604       dsts[i] = req->dst;
1605       recvs[i] = req->recv;
1606       comms[i] = req->comm;
1607       valid[i]=1;;
1608     }else{
1609       valid[i]=0;
1610     }
1611   }
1612   int rank_traced = smpi_process_index();
1613   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1614   extra->type = TRACING_WAITALL;
1615   extra->send_size=count;
1616   TRACE_smpi_ptp_in(rank_traced, -1, -1, __FUNCTION__,extra);
1617 #endif
1618   int retval = smpi_mpi_waitall(count, requests, status);
1619 #ifdef HAVE_TRACING
1620   for (i = 0; i < count; i++) {
1621     if(valid[i]){
1622     //int src_traced = srcs[*index];
1623     //the src may not have been known at the beginning of the recv (MPI_ANY_SOURCE)
1624       int src_traced = srcs[i];
1625       int dst_traced = dsts[i];
1626       int is_wait_for_receive = recvs[i];
1627       if (is_wait_for_receive) {
1628         if(src_traced==MPI_ANY_SOURCE)
1629         src_traced = (status!=MPI_STATUSES_IGNORE) ?
1630                           smpi_group_rank(smpi_comm_group(comms[i]), status[i].MPI_SOURCE) :
1631                           srcs[i];
1632         TRACE_smpi_recv(rank_traced, src_traced, dst_traced);
1633       }
1634     }
1635   }
1636   TRACE_smpi_ptp_out(rank_traced, -1, -1, __FUNCTION__);
1637   xbt_free(srcs);
1638   xbt_free(dsts);
1639   xbt_free(recvs);
1640   xbt_free(valid);
1641   xbt_free(comms);
1642
1643 #endif
1644   smpi_bench_begin();
1645   return retval;
1646 }
1647
1648 int PMPI_Waitsome(int incount, MPI_Request requests[], int *outcount,
1649                  int *indices, MPI_Status status[])
1650 {
1651   int retval = 0;
1652
1653   smpi_bench_end();
1654   if (outcount == NULL) {
1655     retval = MPI_ERR_ARG;
1656   } else {
1657     *outcount = smpi_mpi_waitsome(incount, requests, indices, status);
1658     retval = MPI_SUCCESS;
1659   }
1660   smpi_bench_begin();
1661   return retval;
1662 }
1663
1664 int PMPI_Testsome(int incount, MPI_Request requests[], int* outcount,
1665                  int* indices, MPI_Status status[])
1666 {
1667   int retval = 0;
1668
1669    smpi_bench_end();
1670    if (outcount == NULL) {
1671      retval = MPI_ERR_ARG;
1672    } else {
1673      *outcount = smpi_mpi_testsome(incount, requests, indices, status);
1674      retval = MPI_SUCCESS;
1675    }
1676    smpi_bench_begin();
1677    return retval;
1678 }
1679
1680
1681 int PMPI_Bcast(void *buf, int count, MPI_Datatype datatype, int root, MPI_Comm comm)
1682 {
1683   int retval = 0;
1684
1685   smpi_bench_end();
1686
1687   if (comm == MPI_COMM_NULL) {
1688     retval = MPI_ERR_COMM;
1689   } else {
1690 #ifdef HAVE_TRACING
1691   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1692   int root_traced = smpi_group_index(smpi_comm_group(comm), root);
1693
1694   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1695   extra->type = TRACING_BCAST;
1696   extra->send_size = count;
1697   extra->root = root_traced;
1698   extra->datatype1 = encode_datatype(datatype);
1699   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__, extra);
1700
1701 #endif
1702     mpi_coll_bcast_fun(buf, count, datatype, root, comm);
1703     retval = MPI_SUCCESS;
1704 #ifdef HAVE_TRACING
1705   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1706 #endif
1707   }
1708
1709   smpi_bench_begin();
1710   return retval;
1711 }
1712
1713 int PMPI_Barrier(MPI_Comm comm)
1714 {
1715   int retval = 0;
1716
1717   smpi_bench_end();
1718
1719   if (comm == MPI_COMM_NULL) {
1720     retval = MPI_ERR_COMM;
1721   } else {
1722 #ifdef HAVE_TRACING
1723   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1724   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1725   extra->type = TRACING_BARRIER;
1726   TRACE_smpi_collective_in(rank, -1, __FUNCTION__, extra);
1727 #endif
1728     mpi_coll_barrier_fun(comm);
1729     retval = MPI_SUCCESS;
1730 #ifdef HAVE_TRACING
1731   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1732 #endif
1733   }
1734
1735   smpi_bench_begin();
1736   return retval;
1737 }
1738
1739 int PMPI_Gather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1740                void *recvbuf, int recvcount, MPI_Datatype recvtype,
1741                int root, MPI_Comm comm)
1742 {
1743   int retval = 0;
1744
1745   smpi_bench_end();
1746
1747   if (comm == MPI_COMM_NULL) {
1748     retval = MPI_ERR_COMM;
1749   } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
1750             ((smpi_comm_rank(comm) == root) && (recvtype == MPI_DATATYPE_NULL))){
1751     retval = MPI_ERR_TYPE;
1752   } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) ||
1753             ((smpi_comm_rank(comm) == root) && (recvcount <0))){
1754     retval = MPI_ERR_COUNT;
1755   } else {
1756
1757     char* sendtmpbuf = (char*) sendbuf;
1758     int sendtmpcount = sendcount;
1759     MPI_Datatype sendtmptype = sendtype;
1760     if( (smpi_comm_rank(comm) == root) && (sendbuf == MPI_IN_PLACE )) {
1761       sendtmpcount=0;
1762       sendtmptype=recvtype;
1763     }
1764 #ifdef HAVE_TRACING
1765   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1766   int root_traced = smpi_group_index(smpi_comm_group(comm), root);
1767   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1768   extra->type = TRACING_GATHER;
1769   extra->send_size = sendtmpcount;
1770   extra->recv_size = recvcount;
1771   extra->root = root_traced;
1772   extra->datatype1 = encode_datatype(sendtmptype);
1773   extra->datatype2 = encode_datatype(recvtype);
1774
1775   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__, extra);
1776 #endif
1777     mpi_coll_gather_fun(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcount,
1778                     recvtype, root, comm);
1779
1780
1781     retval = MPI_SUCCESS;
1782 #ifdef HAVE_TRACING
1783   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1784 #endif
1785   }
1786
1787   smpi_bench_begin();
1788   return retval;
1789 }
1790
1791 int PMPI_Gatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1792                 void *recvbuf, int *recvcounts, int *displs,
1793                 MPI_Datatype recvtype, int root, MPI_Comm comm)
1794 {
1795   int retval = 0;
1796
1797   smpi_bench_end();
1798
1799   if (comm == MPI_COMM_NULL) {
1800     retval = MPI_ERR_COMM;
1801   } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
1802             ((smpi_comm_rank(comm) == root) && (recvtype == MPI_DATATYPE_NULL))){
1803     retval = MPI_ERR_TYPE;
1804   } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
1805     retval = MPI_ERR_COUNT;
1806   } else if (recvcounts == NULL || displs == NULL) {
1807     retval = MPI_ERR_ARG;
1808   } else {
1809     char* sendtmpbuf = (char*) sendbuf;
1810     int sendtmpcount = sendcount;
1811     MPI_Datatype sendtmptype = sendtype;
1812     if( (smpi_comm_rank(comm) == root) && (sendbuf == MPI_IN_PLACE )) {
1813       sendtmpcount=0;
1814       sendtmptype=recvtype;
1815     }
1816
1817 #ifdef HAVE_TRACING
1818   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1819   int root_traced = smpi_group_index(smpi_comm_group(comm), root);
1820   int i=0;
1821   int size = smpi_comm_size(comm);
1822   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1823   extra->type = TRACING_GATHERV;
1824   extra->send_size = sendtmpcount;
1825   extra->recvcounts= xbt_malloc(size*sizeof(int));
1826   for(i=0; i< size; i++)//copy data to avoid bad free
1827     extra->recvcounts[i] = recvcounts[i];
1828   extra->num_processes = size;
1829   extra->root = root_traced;
1830   extra->datatype1 = encode_datatype(sendtmptype);
1831   extra->datatype2 = encode_datatype(recvtype);
1832
1833   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__,extra);
1834 #endif
1835     smpi_mpi_gatherv(sendtmpbuf, sendtmpcount, sendtmptype, recvbuf, recvcounts,
1836                      displs, recvtype, root, comm);
1837     retval = MPI_SUCCESS;
1838 #ifdef HAVE_TRACING
1839   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1840 #endif
1841   }
1842
1843   smpi_bench_begin();
1844   return retval;
1845 }
1846
1847 int PMPI_Allgather(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1848                   void *recvbuf, int recvcount, MPI_Datatype recvtype,
1849                   MPI_Comm comm)
1850 {
1851   int retval = 0;
1852
1853   smpi_bench_end();
1854
1855   if (comm == MPI_COMM_NULL) {
1856     retval = MPI_ERR_COMM;
1857   } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
1858             (recvtype == MPI_DATATYPE_NULL)){
1859     retval = MPI_ERR_TYPE;
1860   } else if ((( sendbuf != MPI_IN_PLACE) && (sendcount <0)) ||
1861             (recvcount <0)){
1862     retval = MPI_ERR_COUNT;
1863   } else {
1864     if(sendbuf == MPI_IN_PLACE) {
1865       sendbuf=((char*)recvbuf)+smpi_datatype_get_extent(recvtype)*recvcount*smpi_comm_rank(comm);
1866       sendcount=recvcount;
1867       sendtype=recvtype;
1868     }
1869 #ifdef HAVE_TRACING
1870   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1871   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1872   extra->type = TRACING_ALLGATHER;
1873   extra->send_size = sendcount;
1874   extra->recv_size = recvcount;
1875   extra->datatype1 = encode_datatype(sendtype);
1876   extra->datatype2 = encode_datatype(recvtype);
1877
1878   TRACE_smpi_collective_in(rank, -1, __FUNCTION__, extra);
1879 #endif
1880     mpi_coll_allgather_fun(sendbuf, sendcount, sendtype, recvbuf, recvcount,
1881                            recvtype, comm);
1882     retval = MPI_SUCCESS;
1883
1884 #ifdef HAVE_TRACING
1885   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1886 #endif
1887   }
1888   smpi_bench_begin();
1889   return retval;
1890 }
1891
1892 int PMPI_Allgatherv(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1893                    void *recvbuf, int *recvcounts, int *displs,
1894                    MPI_Datatype recvtype, MPI_Comm comm)
1895 {
1896   int retval = 0;
1897
1898   smpi_bench_end();
1899
1900   if (comm == MPI_COMM_NULL) {
1901     retval = MPI_ERR_COMM;
1902   } else if ((( sendbuf != MPI_IN_PLACE) && (sendtype == MPI_DATATYPE_NULL)) ||
1903             (recvtype == MPI_DATATYPE_NULL)){
1904     retval = MPI_ERR_TYPE;
1905   } else if (( sendbuf != MPI_IN_PLACE) && (sendcount <0)){
1906     retval = MPI_ERR_COUNT;
1907   } else if (recvcounts == NULL || displs == NULL) {
1908     retval = MPI_ERR_ARG;
1909   } else {
1910
1911     if(sendbuf == MPI_IN_PLACE) {
1912       sendbuf=((char*)recvbuf)+smpi_datatype_get_extent(recvtype)*displs[smpi_comm_rank(comm)];
1913       sendcount=recvcounts[smpi_comm_rank(comm)];
1914       sendtype=recvtype;
1915     }
1916 #ifdef HAVE_TRACING
1917   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1918   int i=0;
1919   int size = smpi_comm_size(comm);
1920   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1921   extra->type = TRACING_ALLGATHERV;
1922   extra->send_size = sendcount;
1923   extra->recvcounts= xbt_malloc(size*sizeof(int));
1924   for(i=0; i< size; i++)//copy data to avoid bad free
1925     extra->recvcounts[i] = recvcounts[i];
1926   extra->num_processes = size;
1927   extra->datatype1 = encode_datatype(sendtype);
1928   extra->datatype2 = encode_datatype(recvtype);
1929
1930   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
1931 #endif
1932     mpi_coll_allgatherv_fun(sendbuf, sendcount, sendtype, recvbuf, recvcounts,
1933                         displs, recvtype, comm);
1934     retval = MPI_SUCCESS;
1935 #ifdef HAVE_TRACING
1936   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
1937 #endif
1938   }
1939
1940   smpi_bench_begin();
1941   return retval;
1942 }
1943
1944 int PMPI_Scatter(void *sendbuf, int sendcount, MPI_Datatype sendtype,
1945                 void *recvbuf, int recvcount, MPI_Datatype recvtype,
1946                 int root, MPI_Comm comm)
1947 {
1948   int retval = 0;
1949
1950   smpi_bench_end();
1951
1952   if (comm == MPI_COMM_NULL) {
1953     retval = MPI_ERR_COMM;
1954   } else if (((smpi_comm_rank(comm)==root) && (sendtype == MPI_DATATYPE_NULL))
1955              || ((recvbuf !=MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
1956     retval = MPI_ERR_TYPE;
1957   } else {
1958
1959     if (recvbuf == MPI_IN_PLACE) {
1960         recvtype=sendtype;
1961         recvcount=sendcount;
1962     }
1963 #ifdef HAVE_TRACING
1964   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
1965   int root_traced = smpi_group_index(smpi_comm_group(comm), root);
1966   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
1967   extra->type = TRACING_SCATTER;
1968   extra->send_size = sendcount;
1969   extra->recv_size= recvcount;
1970   extra->root = root_traced;
1971   extra->datatype1 = encode_datatype(sendtype);
1972   extra->datatype2 = encode_datatype(recvtype);
1973
1974   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__,extra);
1975 #endif
1976     mpi_coll_scatter_fun(sendbuf, sendcount, sendtype, recvbuf, recvcount,
1977                      recvtype, root, comm);
1978     retval = MPI_SUCCESS;
1979 #ifdef HAVE_TRACING
1980   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
1981 #endif
1982   }
1983
1984   smpi_bench_begin();
1985   return retval;
1986 }
1987
1988 int PMPI_Scatterv(void *sendbuf, int *sendcounts, int *displs,
1989                  MPI_Datatype sendtype, void *recvbuf, int recvcount,
1990                  MPI_Datatype recvtype, int root, MPI_Comm comm)
1991 {
1992   int retval = 0;
1993
1994   smpi_bench_end();
1995
1996   if (comm == MPI_COMM_NULL) {
1997     retval = MPI_ERR_COMM;
1998   } else if (sendcounts == NULL || displs == NULL) {
1999     retval = MPI_ERR_ARG;
2000   } else if (((smpi_comm_rank(comm)==root) && (sendtype == MPI_DATATYPE_NULL))
2001              || ((recvbuf !=MPI_IN_PLACE) && (recvtype == MPI_DATATYPE_NULL))) {
2002     retval = MPI_ERR_TYPE;
2003   } else {
2004     if (recvbuf == MPI_IN_PLACE) {
2005         recvtype=sendtype;
2006         recvcount=sendcounts[smpi_comm_rank(comm)];
2007     }
2008 #ifdef HAVE_TRACING
2009   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2010   int root_traced = smpi_group_index(smpi_comm_group(comm), root);
2011   int i=0;
2012   int size = smpi_comm_size(comm);
2013   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2014   extra->type = TRACING_SCATTERV;
2015   extra->recv_size = recvcount;
2016   extra->sendcounts= xbt_malloc(size*sizeof(int));
2017   for(i=0; i< size; i++)//copy data to avoid bad free
2018     extra->sendcounts[i] = sendcounts[i];
2019   extra->num_processes = size;
2020   extra->root = root_traced;
2021   extra->datatype1 = encode_datatype(sendtype);
2022   extra->datatype2 = encode_datatype(recvtype);
2023
2024   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__,extra);
2025
2026 #endif
2027     smpi_mpi_scatterv(sendbuf, sendcounts, displs, sendtype, recvbuf,
2028                       recvcount, recvtype, root, comm);
2029     retval = MPI_SUCCESS;
2030 #ifdef HAVE_TRACING
2031   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
2032 #endif
2033   }
2034
2035   smpi_bench_begin();
2036   return retval;
2037 }
2038
2039 int PMPI_Reduce(void *sendbuf, void *recvbuf, int count,
2040                MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm)
2041 {
2042   int retval = 0;
2043
2044   smpi_bench_end();
2045
2046   if (comm == MPI_COMM_NULL) {
2047     retval = MPI_ERR_COMM;
2048   } else if (!is_datatype_valid(datatype) || op == MPI_OP_NULL) {
2049     retval = MPI_ERR_ARG;
2050   } else {
2051 #ifdef HAVE_TRACING
2052   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2053   int root_traced = smpi_group_index(smpi_comm_group(comm), root);
2054   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2055   extra->type = TRACING_REDUCE;
2056   extra->send_size = count;
2057   extra->datatype1 = encode_datatype(datatype);
2058   extra->root = root_traced;
2059
2060   TRACE_smpi_collective_in(rank, root_traced, __FUNCTION__,extra);
2061 #endif
2062     mpi_coll_reduce_fun(sendbuf, recvbuf, count, datatype, op, root, comm);
2063
2064     retval = MPI_SUCCESS;
2065 #ifdef HAVE_TRACING
2066   TRACE_smpi_collective_out(rank, root_traced, __FUNCTION__);
2067 #endif
2068   }
2069
2070   smpi_bench_begin();
2071   return retval;
2072 }
2073
2074 int PMPI_Reduce_local(void *inbuf, void *inoutbuf, int count,
2075     MPI_Datatype datatype, MPI_Op op){
2076   int retval = 0;
2077
2078     smpi_bench_end();
2079     if (!is_datatype_valid(datatype) || op == MPI_OP_NULL) {
2080       retval = MPI_ERR_ARG;
2081     } else {
2082       smpi_op_apply(op, inbuf, inoutbuf, &count, &datatype);
2083       retval=MPI_SUCCESS;
2084     }
2085     smpi_bench_begin();
2086     return retval;
2087 }
2088
2089 int PMPI_Allreduce(void *sendbuf, void *recvbuf, int count,
2090                   MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
2091 {
2092   int retval = 0;
2093
2094   smpi_bench_end();
2095
2096   if (comm == MPI_COMM_NULL) {
2097     retval = MPI_ERR_COMM;
2098   } else if (!is_datatype_valid(datatype)) {
2099     retval = MPI_ERR_TYPE;
2100   } else if (op == MPI_OP_NULL) {
2101     retval = MPI_ERR_OP;
2102   } else {
2103
2104     char* sendtmpbuf = (char*) sendbuf;
2105     if( sendbuf == MPI_IN_PLACE ) {
2106       sendtmpbuf = (char *)xbt_malloc(count*smpi_datatype_get_extent(datatype));
2107       smpi_datatype_copy(recvbuf, count, datatype,sendtmpbuf, count, datatype);
2108     }
2109 #ifdef HAVE_TRACING
2110   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2111   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2112   extra->type = TRACING_ALLREDUCE;
2113   extra->send_size = count;
2114   extra->datatype1 = encode_datatype(datatype);
2115
2116   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
2117 #endif
2118     mpi_coll_allreduce_fun(sendtmpbuf, recvbuf, count, datatype, op, comm);
2119
2120     if( sendbuf == MPI_IN_PLACE ) {
2121       xbt_free(sendtmpbuf);
2122     }
2123
2124     retval = MPI_SUCCESS;
2125 #ifdef HAVE_TRACING
2126   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
2127 #endif
2128   }
2129
2130   smpi_bench_begin();
2131   return retval;
2132 }
2133
2134 int PMPI_Scan(void *sendbuf, void *recvbuf, int count,
2135              MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
2136 {
2137   int retval = 0;
2138
2139   smpi_bench_end();
2140
2141   if (comm == MPI_COMM_NULL) {
2142     retval = MPI_ERR_COMM;
2143   } else if (!is_datatype_valid(datatype)) {
2144     retval = MPI_ERR_TYPE;
2145   } else if (op == MPI_OP_NULL) {
2146     retval = MPI_ERR_OP;
2147   } else {
2148 #ifdef HAVE_TRACING
2149   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2150   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2151   extra->type = TRACING_SCAN;
2152   extra->send_size = count;
2153   extra->datatype1 = encode_datatype(datatype);
2154
2155   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
2156 #endif
2157     smpi_mpi_scan(sendbuf, recvbuf, count, datatype, op, comm);
2158     retval = MPI_SUCCESS;
2159 #ifdef HAVE_TRACING
2160   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
2161 #endif
2162   }
2163
2164   smpi_bench_begin();
2165   return retval;
2166 }
2167
2168 int PMPI_Exscan(void *sendbuf, void *recvbuf, int count, MPI_Datatype datatype,
2169                 MPI_Op op, MPI_Comm comm){
2170   int retval = 0;
2171
2172   smpi_bench_end();
2173
2174   if (comm == MPI_COMM_NULL) {
2175     retval = MPI_ERR_COMM;
2176   } else if (!is_datatype_valid(datatype)) {
2177     retval = MPI_ERR_TYPE;
2178   } else if (op == MPI_OP_NULL) {
2179     retval = MPI_ERR_OP;
2180   } else {
2181 #ifdef HAVE_TRACING
2182   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2183   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2184   extra->type = TRACING_EXSCAN;
2185   extra->send_size = count;
2186   extra->datatype1 = encode_datatype(datatype);
2187
2188   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
2189 #endif
2190     smpi_mpi_exscan(sendbuf, recvbuf, count, datatype, op, comm);
2191     retval = MPI_SUCCESS;
2192 #ifdef HAVE_TRACING
2193   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
2194 #endif
2195   }
2196
2197   smpi_bench_begin();
2198   return retval;
2199 }
2200
2201 int PMPI_Reduce_scatter(void *sendbuf, void *recvbuf, int *recvcounts,
2202                        MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
2203 {
2204   int retval = 0;
2205   smpi_bench_end();
2206
2207   if (comm == MPI_COMM_NULL) {
2208     retval = MPI_ERR_COMM;
2209   } else if (!is_datatype_valid(datatype)) {
2210     retval = MPI_ERR_TYPE;
2211   } else if (op == MPI_OP_NULL) {
2212     retval = MPI_ERR_OP;
2213   } else if (recvcounts == NULL) {
2214     retval = MPI_ERR_ARG;
2215   } else {
2216 #ifdef HAVE_TRACING
2217   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2218   int i=0;
2219   int size = smpi_comm_size(comm);
2220   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2221   extra->type = TRACING_REDUCE_SCATTER;
2222   extra->send_size = 0;
2223   extra->recvcounts= xbt_malloc(size*sizeof(int));
2224   for(i=0; i< size; i++)//copy data to avoid bad free
2225     extra->recvcounts[i] = recvcounts[i];
2226   extra->num_processes = size;
2227   extra->datatype1 = encode_datatype(datatype);
2228
2229   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
2230 #endif
2231     void* sendtmpbuf=sendbuf;
2232     if(sendbuf==MPI_IN_PLACE){
2233       sendtmpbuf=recvbuf;
2234     }
2235
2236     mpi_coll_reduce_scatter_fun(sendtmpbuf, recvbuf, recvcounts,
2237                        datatype,  op, comm);
2238     retval = MPI_SUCCESS;
2239 #ifdef HAVE_TRACING
2240   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
2241 #endif
2242   }
2243
2244   smpi_bench_begin();
2245   return retval;
2246 }
2247
2248 int PMPI_Reduce_scatter_block(void *sendbuf, void *recvbuf, int recvcount,
2249                        MPI_Datatype datatype, MPI_Op op, MPI_Comm comm)
2250 {
2251   int retval,i;
2252   smpi_bench_end();
2253
2254   if (comm == MPI_COMM_NULL) {
2255     retval = MPI_ERR_COMM;
2256   } else if (!is_datatype_valid(datatype)) {
2257     retval = MPI_ERR_TYPE;
2258   } else if (op == MPI_OP_NULL) {
2259     retval = MPI_ERR_OP;
2260   } else if (recvcount < 0) {
2261     retval = MPI_ERR_ARG;
2262   } else {
2263     int count=smpi_comm_size(comm);
2264
2265 #ifdef HAVE_TRACING
2266   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2267   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2268   extra->type = TRACING_REDUCE_SCATTER;
2269   extra->send_size = 0;
2270   extra->recvcounts= xbt_malloc(count*sizeof(int));
2271   for(i=0; i< count; i++)//copy data to avoid bad free
2272     extra->recvcounts[i] = recvcount;
2273   extra->num_processes = count;
2274   extra->datatype1 = encode_datatype(datatype);
2275
2276   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
2277 #endif
2278     int* recvcounts=(int*)xbt_malloc(count);
2279     for (i=0; i<count;i++)recvcounts[i]=recvcount;
2280     mpi_coll_reduce_scatter_fun(sendbuf, recvbuf, recvcounts,
2281                        datatype,  op, comm);
2282     xbt_free(recvcounts);
2283     retval = MPI_SUCCESS;
2284 #ifdef HAVE_TRACING
2285   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
2286 #endif
2287   }
2288
2289   smpi_bench_begin();
2290   return retval;
2291 }
2292
2293 int PMPI_Alltoall(void *sendbuf, int sendcount, MPI_Datatype sendtype,
2294                  void *recvbuf, int recvcount, MPI_Datatype recvtype,
2295                  MPI_Comm comm)
2296 {
2297   int retval = 0;
2298
2299   smpi_bench_end();
2300
2301   if (comm == MPI_COMM_NULL) {
2302     retval = MPI_ERR_COMM;
2303   } else if (sendtype == MPI_DATATYPE_NULL
2304              || recvtype == MPI_DATATYPE_NULL) {
2305     retval = MPI_ERR_TYPE;
2306   } else {
2307 #ifdef HAVE_TRACING
2308   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2309   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2310   extra->type = TRACING_ALLTOALL;
2311   extra->send_size = sendcount;
2312   extra->recv_size = recvcount;
2313   extra->datatype1 = encode_datatype(sendtype);
2314   extra->datatype2 = encode_datatype(recvtype);
2315
2316   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
2317 #endif
2318     retval = mpi_coll_alltoall_fun(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm);
2319 #ifdef HAVE_TRACING
2320   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
2321 #endif
2322   }
2323
2324   smpi_bench_begin();
2325   return retval;
2326 }
2327
2328 int PMPI_Alltoallv(void *sendbuf, int *sendcounts, int *senddisps,
2329                   MPI_Datatype sendtype, void *recvbuf, int *recvcounts,
2330                   int *recvdisps, MPI_Datatype recvtype, MPI_Comm comm)
2331 {
2332   int retval = 0;
2333
2334   smpi_bench_end();
2335
2336   if (comm == MPI_COMM_NULL) {
2337     retval = MPI_ERR_COMM;
2338   } else if (sendtype == MPI_DATATYPE_NULL
2339              || recvtype == MPI_DATATYPE_NULL) {
2340     retval = MPI_ERR_TYPE;
2341   } else if (sendcounts == NULL || senddisps == NULL || recvcounts == NULL
2342              || recvdisps == NULL) {
2343     retval = MPI_ERR_ARG;
2344   } else {
2345 #ifdef HAVE_TRACING
2346   int rank = comm != MPI_COMM_NULL ? smpi_process_index() : -1;
2347   int i=0;
2348   int size = smpi_comm_size(comm);
2349   instr_extra_data extra = xbt_new0(s_instr_extra_data_t,1);
2350   extra->type = TRACING_ALLTOALLV;
2351   extra->send_size = 0;
2352   extra->recv_size = 0;
2353   extra->recvcounts= xbt_malloc(size*sizeof(int));
2354   extra->sendcounts= xbt_malloc(size*sizeof(int));
2355
2356   for(i=0; i< size; i++){//copy data to avoid bad free
2357     extra->send_size += sendcounts[i];
2358     extra->recv_size += recvcounts[i];
2359
2360     extra->sendcounts[i] = sendcounts[i];
2361     extra->recvcounts[i] = recvcounts[i];
2362   }
2363   extra->num_processes = size;
2364
2365   extra->datatype1 = encode_datatype(sendtype);
2366   extra->datatype2 = encode_datatype(recvtype);
2367
2368   TRACE_smpi_collective_in(rank, -1, __FUNCTION__,extra);
2369 #endif
2370     retval =
2371         mpi_coll_alltoallv_fun(sendbuf, sendcounts, senddisps, sendtype,
2372                                   recvbuf, recvcounts, recvdisps, recvtype,
2373                                   comm);
2374 #ifdef HAVE_TRACING
2375   TRACE_smpi_collective_out(rank, -1, __FUNCTION__);
2376 #endif
2377   }
2378
2379   smpi_bench_begin();
2380   return retval;
2381 }
2382
2383
2384 int PMPI_Get_processor_name(char *name, int *resultlen)
2385 {
2386   int retval = MPI_SUCCESS;
2387
2388   strncpy(name, SIMIX_host_get_name(SIMIX_host_self()),
2389           strlen(SIMIX_host_get_name(SIMIX_host_self())) < MPI_MAX_PROCESSOR_NAME - 1 ?
2390           strlen(SIMIX_host_get_name(SIMIX_host_self())) +1 :
2391           MPI_MAX_PROCESSOR_NAME - 1 );
2392   *resultlen =
2393       strlen(name) >
2394       MPI_MAX_PROCESSOR_NAME ? MPI_MAX_PROCESSOR_NAME : strlen(name);
2395
2396   return retval;
2397 }
2398
2399 int PMPI_Get_count(MPI_Status * status, MPI_Datatype datatype, int *count)
2400 {
2401   int retval = MPI_SUCCESS;
2402   size_t size;
2403
2404   if (status == NULL || count == NULL) {
2405     retval = MPI_ERR_ARG;
2406   } else if (!is_datatype_valid(datatype)) {
2407     retval = MPI_ERR_TYPE;
2408   } else {
2409     size = smpi_datatype_size(datatype);
2410     if (size == 0) {
2411       *count = 0;
2412     } else if (status->count % size != 0) {
2413       retval = MPI_UNDEFINED;
2414     } else {
2415       *count = smpi_mpi_get_count(status, datatype);
2416     }
2417   }
2418   return retval;
2419 }
2420
2421 int PMPI_Type_contiguous(int count, MPI_Datatype old_type, MPI_Datatype* new_type) {
2422   int retval = 0;
2423
2424   if (old_type == MPI_DATATYPE_NULL) {
2425     retval = MPI_ERR_TYPE;
2426   } else if (count<0){
2427     retval = MPI_ERR_COUNT;
2428   } else {
2429     retval = smpi_datatype_contiguous(count, old_type, new_type, 0);
2430   }
2431   return retval;
2432 }
2433
2434 int PMPI_Type_commit(MPI_Datatype* datatype) {
2435   int retval = 0;
2436
2437   if (datatype == NULL || *datatype == MPI_DATATYPE_NULL) {
2438     retval = MPI_ERR_TYPE;
2439   } else {
2440     smpi_datatype_commit(datatype);
2441     retval = MPI_SUCCESS;
2442   }
2443   return retval;
2444 }
2445
2446
2447 int PMPI_Type_vector(int count, int blocklen, int stride, MPI_Datatype old_type, MPI_Datatype* new_type) {
2448   int retval = 0;
2449
2450   if (old_type == MPI_DATATYPE_NULL) {
2451     retval = MPI_ERR_TYPE;
2452   } else if (count<0 || blocklen<0){
2453     retval = MPI_ERR_COUNT;
2454   } else {
2455     retval = smpi_datatype_vector(count, blocklen, stride, old_type, new_type);
2456   }
2457   return retval;
2458 }
2459
2460 int PMPI_Type_hvector(int count, int blocklen, MPI_Aint stride, MPI_Datatype old_type, MPI_Datatype* new_type) {
2461   int retval = 0;
2462
2463   if (old_type == MPI_DATATYPE_NULL) {
2464     retval = MPI_ERR_TYPE;
2465   } else if (count<0 || blocklen<0){
2466     retval = MPI_ERR_COUNT;
2467   } else {
2468     retval = smpi_datatype_hvector(count, blocklen, stride, old_type, new_type);
2469   }
2470   return retval;
2471 }
2472
2473 int PMPI_Type_create_hvector(int count, int blocklen, MPI_Aint stride, MPI_Datatype old_type, MPI_Datatype* new_type) {
2474   return MPI_Type_hvector(count, blocklen, stride, old_type, new_type);
2475 }
2476
2477 int PMPI_Type_indexed(int count, int* blocklens, int* indices, MPI_Datatype old_type, MPI_Datatype* new_type) {
2478   int retval = 0;
2479
2480   if (old_type == MPI_DATATYPE_NULL) {
2481     retval = MPI_ERR_TYPE;
2482   } else if (count<0){
2483     retval = MPI_ERR_COUNT;
2484   } else {
2485     retval = smpi_datatype_indexed(count, blocklens, indices, old_type, new_type);
2486   }
2487   return retval;
2488 }
2489
2490 int PMPI_Type_create_indexed(int count, int* blocklens, int* indices, MPI_Datatype old_type, MPI_Datatype* new_type) {
2491   int retval = 0;
2492
2493   if (old_type == MPI_DATATYPE_NULL) {
2494     retval = MPI_ERR_TYPE;
2495   } else if (count<0){
2496     retval = MPI_ERR_COUNT;
2497   } else {
2498     retval = smpi_datatype_indexed(count, blocklens, indices, old_type, new_type);
2499   }
2500   return retval;
2501 }
2502
2503 int PMPI_Type_create_indexed_block(int count, int blocklength, int* indices, MPI_Datatype old_type, MPI_Datatype* new_type) {
2504   int retval,i;
2505
2506   if (old_type == MPI_DATATYPE_NULL) {
2507     retval = MPI_ERR_TYPE;
2508   } else if (count<0){
2509     retval = MPI_ERR_COUNT;
2510   } else {
2511     int* blocklens=(int*)xbt_malloc(blocklength*count);
2512     for (i=0; i<count;i++)blocklens[i]=blocklength;
2513     retval = smpi_datatype_indexed(count, blocklens, indices, old_type, new_type);
2514     xbt_free(blocklens);
2515   }
2516   return retval;
2517 }
2518
2519
2520 int PMPI_Type_hindexed(int count, int* blocklens, MPI_Aint* indices, MPI_Datatype old_type, MPI_Datatype* new_type) {
2521   int retval = 0;
2522
2523   if (old_type == MPI_DATATYPE_NULL) {
2524     retval = MPI_ERR_TYPE;
2525   } else if (count<0){
2526     retval = MPI_ERR_COUNT;
2527   } else {
2528     retval = smpi_datatype_hindexed(count, blocklens, indices, old_type, new_type);
2529   }
2530   return retval;
2531 }
2532
2533 int PMPI_Type_create_hindexed(int count, int* blocklens, MPI_Aint* indices, MPI_Datatype old_type, MPI_Datatype* new_type) {
2534   return PMPI_Type_hindexed(count, blocklens,indices,old_type,new_type);
2535 }
2536
2537 int PMPI_Type_create_hindexed_block(int count, int blocklength, MPI_Aint* indices, MPI_Datatype old_type, MPI_Datatype* new_type) {
2538   int retval,i;
2539
2540   if (old_type == MPI_DATATYPE_NULL) {
2541     retval = MPI_ERR_TYPE;
2542   } else if (count<0){
2543     retval = MPI_ERR_COUNT;
2544   } else {
2545     int* blocklens=(int*)xbt_malloc(blocklength*count);
2546     for (i=0; i<count;i++)blocklens[i]=blocklength;
2547     retval = smpi_datatype_hindexed(count, blocklens, indices, old_type, new_type);
2548     xbt_free(blocklens);
2549   }
2550   return retval;
2551 }
2552
2553
2554 int PMPI_Type_struct(int count, int* blocklens, MPI_Aint* indices, MPI_Datatype* old_types, MPI_Datatype* new_type) {
2555   int retval = 0;
2556
2557   if (count<0){
2558     retval = MPI_ERR_COUNT;
2559   } else {
2560     retval = smpi_datatype_struct(count, blocklens, indices, old_types, new_type);
2561   }
2562   return retval;
2563 }
2564
2565 int PMPI_Type_create_struct(int count, int* blocklens, MPI_Aint* indices, MPI_Datatype* old_types, MPI_Datatype* new_type) {
2566   return PMPI_Type_struct(count, blocklens, indices, old_types, new_type);
2567 }
2568
2569
2570 int PMPI_Error_class(int errorcode, int* errorclass) {
2571   // assume smpi uses only standard mpi error codes
2572   *errorclass=errorcode;
2573   return MPI_SUCCESS;
2574 }
2575
2576
2577 int PMPI_Initialized(int* flag) {
2578    *flag=smpi_process_initialized();
2579    return MPI_SUCCESS;
2580 }
2581
2582 /* The topo part of MPI_COMM_WORLD should always be NULL. When other topologies
2583  * will be implemented, not only should we check if the topology is NULL, but
2584  * we should check if it is the good topology type (so we have to add a
2585  *  MPIR_Topo_Type field, and replace the MPI_Topology field by an union)*/
2586
2587 int PMPI_Cart_create(MPI_Comm comm_old, int ndims, int* dims, int* periodic, int reorder, MPI_Comm* comm_cart) {
2588   int retval = 0;
2589   smpi_bench_end();
2590   if (comm_old == MPI_COMM_NULL){
2591     return  MPI_ERR_COMM;
2592   }
2593   else if (ndims < 0 ||
2594            (ndims > 0 && (dims == NULL || 
2595                           periodic == NULL)) ||
2596            comm_cart == NULL) {
2597     return MPI_ERR_ARG;
2598   }
2599   retval = smpi_mpi_cart_create(comm_old, ndims, dims, periodic, reorder, comm_cart);
2600
2601   smpi_bench_begin();
2602
2603   return retval;
2604 }
2605
2606 int PMPI_Cart_rank(MPI_Comm comm, int* coords, int* rank) {
2607   if(comm == MPI_COMM_NULL || smpi_comm_topo(comm) == NULL) {
2608     return MPI_ERR_TOPOLOGY;
2609   }
2610   if (coords == NULL) {
2611     return MPI_ERR_ARG;
2612   }
2613   return smpi_mpi_cart_rank(comm, coords, rank);
2614 }
2615
2616 int PMPI_Cart_shift(MPI_Comm comm, int direction, int displ, int* source, int* dest) {
2617   if(comm == MPI_COMM_NULL || smpi_comm_topo(comm) == NULL) {
2618     return MPI_ERR_TOPOLOGY;
2619   }
2620   if (source == NULL || dest == NULL || direction < 0 ) {
2621     return MPI_ERR_ARG;
2622   }
2623   return smpi_mpi_cart_shift(comm, direction, displ, source, dest);
2624 }
2625
2626 int PMPI_Cart_coords(MPI_Comm comm, int rank, int maxdims, int* coords) {
2627   if(comm == MPI_COMM_NULL || smpi_comm_topo(comm) == NULL) {
2628     return MPI_ERR_TOPOLOGY;
2629   }
2630   if (rank < 0 || rank >= smpi_comm_size(comm)) {
2631     return MPI_ERR_RANK;
2632   }
2633   if (maxdims <= 0) {
2634     return MPI_ERR_ARG;
2635   }
2636   if(coords == NULL) {
2637     return MPI_ERR_ARG;
2638   }
2639   return smpi_mpi_cart_coords(comm, rank, maxdims, coords);
2640 }
2641
2642 int PMPI_Cart_get(MPI_Comm comm, int maxdims, int* dims, int* periods, int* coords) {
2643   if(comm == NULL || smpi_comm_topo(comm) == NULL) {
2644     return MPI_ERR_TOPOLOGY;
2645   }
2646   if(maxdims <= 0 || dims == NULL || periods == NULL || coords == NULL) {
2647     return MPI_ERR_ARG;
2648   }
2649   return smpi_mpi_cart_get(comm, maxdims, dims, periods, coords);
2650 }
2651
2652 int PMPI_Cartdim_get(MPI_Comm comm, int* ndims) {
2653   if (comm == MPI_COMM_NULL || smpi_comm_topo(comm) == NULL) {
2654     return MPI_ERR_TOPOLOGY;
2655   }
2656   if (ndims == NULL) {
2657     return MPI_ERR_ARG;
2658   }
2659   return smpi_mpi_cartdim_get(comm, ndims);
2660 }
2661
2662 int PMPI_Dims_create(int nnodes, int ndims, int* dims) {
2663   if(dims == NULL) {
2664     return MPI_ERR_ARG;
2665   }
2666   if (ndims < 1 || nnodes < 1) {
2667     return MPI_ERR_DIMS;
2668   }
2669
2670   return smpi_mpi_dims_create(nnodes, ndims, dims);
2671 }
2672
2673 int PMPI_Cart_sub(MPI_Comm comm, int* remain_dims, MPI_Comm* comm_new) {
2674   if(comm == MPI_COMM_NULL || smpi_comm_topo(comm) == NULL) {
2675     return MPI_ERR_TOPOLOGY;
2676   }
2677   if (comm_new == NULL) {
2678     return MPI_ERR_ARG;
2679   }
2680   return smpi_mpi_cart_sub(comm, remain_dims, comm_new);
2681 }
2682
2683
2684 /* The following calls are not yet implemented and will fail at runtime. */
2685 /* Once implemented, please move them above this notice. */
2686
2687 #define NOT_YET_IMPLEMENTED {                                           \
2688     XBT_WARN("Not yet implemented : %s. Please contact the Simgrid team if support is needed", __FUNCTION__); \
2689     return MPI_SUCCESS;                                                 \
2690   }
2691
2692 int PMPI_Type_create_resized(MPI_Datatype oldtype,MPI_Aint lb, MPI_Aint extent, MPI_Datatype *newtype){
2693   NOT_YET_IMPLEMENTED
2694 }
2695
2696 int PMPI_Type_dup(MPI_Datatype datatype, MPI_Datatype *newtype){
2697   NOT_YET_IMPLEMENTED
2698 }
2699
2700 int PMPI_Type_set_name(MPI_Datatype  datatype, char * name)
2701 {
2702   NOT_YET_IMPLEMENTED
2703 }
2704
2705 int PMPI_Type_get_name(MPI_Datatype  datatype, char * name, int* len)
2706 {
2707   NOT_YET_IMPLEMENTED
2708 }
2709
2710 int PMPI_Pack_size(int incount, MPI_Datatype datatype, MPI_Comm comm, int* size) {
2711   NOT_YET_IMPLEMENTED
2712 }
2713
2714
2715 int PMPI_Cart_map(MPI_Comm comm_old, int ndims, int* dims, int* periods, int* newrank) {
2716   NOT_YET_IMPLEMENTED
2717 }
2718
2719
2720 int PMPI_Graph_create(MPI_Comm comm_old, int nnodes, int* index, int* edges, int reorder, MPI_Comm* comm_graph) {
2721   NOT_YET_IMPLEMENTED
2722 }
2723
2724 int PMPI_Graph_get(MPI_Comm comm, int maxindex, int maxedges, int* index, int* edges) {
2725   NOT_YET_IMPLEMENTED
2726 }
2727
2728 int PMPI_Graph_map(MPI_Comm comm_old, int nnodes, int* index, int* edges, int* newrank) {
2729   NOT_YET_IMPLEMENTED
2730 }
2731
2732 int PMPI_Graph_neighbors(MPI_Comm comm, int rank, int maxneighbors, int* neighbors) {
2733   NOT_YET_IMPLEMENTED
2734 }
2735
2736 int PMPI_Graph_neighbors_count(MPI_Comm comm, int rank, int* nneighbors) {
2737   NOT_YET_IMPLEMENTED
2738 }
2739
2740 int PMPI_Graphdims_get(MPI_Comm comm, int* nnodes, int* nedges) {
2741   NOT_YET_IMPLEMENTED
2742 }
2743
2744 int PMPI_Topo_test(MPI_Comm comm, int* top_type) {
2745   NOT_YET_IMPLEMENTED
2746 }
2747
2748 int PMPI_Errhandler_create(MPI_Handler_function* function, MPI_Errhandler* errhandler) {
2749   NOT_YET_IMPLEMENTED
2750 }
2751
2752 int PMPI_Errhandler_free(MPI_Errhandler* errhandler) {
2753   NOT_YET_IMPLEMENTED
2754 }
2755
2756 int PMPI_Errhandler_get(MPI_Comm comm, MPI_Errhandler* errhandler) {
2757   NOT_YET_IMPLEMENTED
2758 }
2759
2760 int PMPI_Error_string(int errorcode, char* string, int* resultlen) {
2761   NOT_YET_IMPLEMENTED
2762 }
2763
2764 int PMPI_Errhandler_set(MPI_Comm comm, MPI_Errhandler errhandler) {
2765   NOT_YET_IMPLEMENTED
2766 }
2767
2768 int PMPI_Comm_set_errhandler(MPI_Comm comm, MPI_Errhandler errhandler) {
2769   NOT_YET_IMPLEMENTED
2770 }
2771
2772 int PMPI_Comm_get_errhandler(MPI_Comm comm, MPI_Errhandler* errhandler) {
2773   NOT_YET_IMPLEMENTED
2774 }
2775
2776 int PMPI_Cancel(MPI_Request* request) {
2777   NOT_YET_IMPLEMENTED
2778 }
2779
2780 int PMPI_Buffer_attach(void* buffer, int size) {
2781   NOT_YET_IMPLEMENTED
2782 }
2783
2784 int PMPI_Buffer_detach(void* buffer, int* size) {
2785   NOT_YET_IMPLEMENTED
2786 }
2787
2788 int PMPI_Comm_test_inter(MPI_Comm comm, int* flag) {
2789   NOT_YET_IMPLEMENTED
2790 }
2791
2792 int PMPI_Comm_get_attr (MPI_Comm comm, int comm_keyval, void *attribute_val, int *flag)
2793 {
2794   NOT_YET_IMPLEMENTED
2795 }
2796
2797 int PMPI_Comm_set_attr (MPI_Comm comm, int comm_keyval, void *attribute_val)
2798 {
2799   NOT_YET_IMPLEMENTED
2800 }
2801
2802 int PMPI_Comm_delete_attr (MPI_Comm comm, int comm_keyval)
2803 {
2804   NOT_YET_IMPLEMENTED
2805 }
2806
2807 int PMPI_Comm_create_keyval(MPI_Comm_copy_attr_function* copy_fn, MPI_Comm_delete_attr_function* delete_fn, int* keyval, void* extra_state)
2808 {
2809   NOT_YET_IMPLEMENTED
2810 }
2811
2812 int PMPI_Comm_free_keyval(int* keyval) {
2813   NOT_YET_IMPLEMENTED
2814 }
2815
2816 int PMPI_Pcontrol(const int level )
2817 {
2818   NOT_YET_IMPLEMENTED
2819 }
2820
2821 int PMPI_Unpack(void* inbuf, int insize, int* position, void* outbuf, int outcount, MPI_Datatype type, MPI_Comm comm) {
2822   NOT_YET_IMPLEMENTED
2823 }
2824
2825 int PMPI_Type_get_attr (MPI_Datatype type, int type_keyval, void *attribute_val, int* flag)
2826 {
2827   NOT_YET_IMPLEMENTED
2828 }
2829
2830 int PMPI_Type_set_attr (MPI_Datatype type, int type_keyval, void *attribute_val)
2831 {
2832   NOT_YET_IMPLEMENTED
2833 }
2834
2835 int PMPI_Type_delete_attr (MPI_Datatype type, int comm_keyval)
2836 {
2837   NOT_YET_IMPLEMENTED
2838 }
2839
2840 int PMPI_Type_create_keyval(MPI_Type_copy_attr_function* copy_fn, MPI_Type_delete_attr_function* delete_fn, int* keyval, void* extra_state)
2841 {
2842   NOT_YET_IMPLEMENTED
2843 }
2844
2845 int PMPI_Type_free_keyval(int* keyval) {
2846   NOT_YET_IMPLEMENTED
2847 }
2848
2849 int PMPI_Intercomm_create(MPI_Comm local_comm, int local_leader, MPI_Comm peer_comm, int remote_leader, int tag, MPI_Comm* comm_out) {
2850   NOT_YET_IMPLEMENTED
2851 }
2852
2853 int PMPI_Intercomm_merge(MPI_Comm comm, int high, MPI_Comm* comm_out) {
2854   NOT_YET_IMPLEMENTED
2855 }
2856
2857 int PMPI_Bsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm) {
2858   NOT_YET_IMPLEMENTED
2859 }
2860
2861 int PMPI_Bsend_init(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
2862   NOT_YET_IMPLEMENTED
2863 }
2864
2865 int PMPI_Ibsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
2866   NOT_YET_IMPLEMENTED
2867 }
2868
2869 int PMPI_Comm_remote_group(MPI_Comm comm, MPI_Group* group) {
2870   NOT_YET_IMPLEMENTED
2871 }
2872
2873 int PMPI_Comm_remote_size(MPI_Comm comm, int* size) {
2874   NOT_YET_IMPLEMENTED
2875 }
2876
2877 int PMPI_Attr_delete(MPI_Comm comm, int keyval) {
2878   NOT_YET_IMPLEMENTED
2879 }
2880
2881 int PMPI_Attr_get(MPI_Comm comm, int keyval, void* attr_value, int* flag) {
2882   NOT_YET_IMPLEMENTED
2883 }
2884
2885 int PMPI_Attr_put(MPI_Comm comm, int keyval, void* attr_value) {
2886   NOT_YET_IMPLEMENTED
2887 }
2888
2889 int PMPI_Rsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm) {
2890   NOT_YET_IMPLEMENTED
2891 }
2892
2893 int PMPI_Rsend_init(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
2894   NOT_YET_IMPLEMENTED
2895 }
2896
2897 int PMPI_Irsend(void* buf, int count, MPI_Datatype datatype, int dest, int tag, MPI_Comm comm, MPI_Request* request) {
2898   NOT_YET_IMPLEMENTED
2899 }
2900
2901 int PMPI_Keyval_create(MPI_Copy_function* copy_fn, MPI_Delete_function* delete_fn, int* keyval, void* extra_state) {
2902   NOT_YET_IMPLEMENTED
2903 }
2904
2905 int PMPI_Keyval_free(int* keyval) {
2906   NOT_YET_IMPLEMENTED
2907 }
2908
2909 int PMPI_Test_cancelled(MPI_Status* status, int* flag) {
2910   NOT_YET_IMPLEMENTED
2911 }
2912
2913 int PMPI_Pack(void* inbuf, int incount, MPI_Datatype type, void* outbuf, int outcount, int* position, MPI_Comm comm) {
2914   NOT_YET_IMPLEMENTED
2915 }
2916
2917 int PMPI_Pack_external_size(char *datarep, int incount, MPI_Datatype datatype, MPI_Aint *size){
2918   NOT_YET_IMPLEMENTED
2919 }
2920
2921 int PMPI_Pack_external(char *datarep, void *inbuf, int incount, MPI_Datatype datatype, void *outbuf, MPI_Aint outcount, MPI_Aint *position){
2922   NOT_YET_IMPLEMENTED
2923 }
2924
2925 int PMPI_Unpack_external( char *datarep, void *inbuf, MPI_Aint insize, MPI_Aint *position, void *outbuf, int outcount, MPI_Datatype datatype){
2926   NOT_YET_IMPLEMENTED
2927 }
2928
2929 int PMPI_Get_elements(MPI_Status* status, MPI_Datatype datatype, int* elements) {
2930   NOT_YET_IMPLEMENTED
2931 }
2932
2933 int PMPI_Win_fence( int assert,  MPI_Win win){
2934   NOT_YET_IMPLEMENTED
2935 }
2936
2937 int PMPI_Win_free( MPI_Win* win){
2938   NOT_YET_IMPLEMENTED
2939 }
2940
2941 int PMPI_Win_create( void *base, MPI_Aint size, int disp_unit, MPI_Info info, MPI_Comm comm, MPI_Win *win){
2942   NOT_YET_IMPLEMENTED
2943 }
2944
2945 int PMPI_Info_create( MPI_Info *info){
2946   NOT_YET_IMPLEMENTED
2947 }
2948
2949 int PMPI_Info_set( MPI_Info info, char *key, char *value){
2950   NOT_YET_IMPLEMENTED
2951 }
2952
2953 int PMPI_Info_free( MPI_Info *info){
2954   NOT_YET_IMPLEMENTED
2955 }
2956
2957 int PMPI_Get( void *origin_addr, int origin_count, MPI_Datatype origin_datatype, int target_rank,
2958               MPI_Aint target_disp, int target_count, MPI_Datatype target_datatype, MPI_Win win){
2959   NOT_YET_IMPLEMENTED
2960 }
2961
2962 int PMPI_Type_get_envelope( MPI_Datatype datatype, int *num_integers,
2963                             int *num_addresses, int *num_datatypes, int *combiner){
2964   NOT_YET_IMPLEMENTED
2965 }
2966
2967 int PMPI_Type_get_contents(MPI_Datatype datatype, int max_integers, int max_addresses,
2968                            int max_datatypes, int* array_of_integers, MPI_Aint* array_of_addresses,
2969                            MPI_Datatype* array_of_datatypes){
2970   NOT_YET_IMPLEMENTED
2971 }
2972
2973 int PMPI_Type_create_darray(int size, int rank, int ndims, int* array_of_gsizes,
2974                             int* array_of_distribs, int* array_of_dargs, int* array_of_psizes,
2975                             int order, MPI_Datatype oldtype, MPI_Datatype *newtype) {
2976   NOT_YET_IMPLEMENTED
2977 }
2978
2979 int PMPI_Type_create_subarray(int ndims,int *array_of_sizes, int *array_of_subsizes, int *array_of_starts, int order, MPI_Datatype oldtype, MPI_Datatype *newtype){
2980   NOT_YET_IMPLEMENTED
2981 }
2982
2983 int PMPI_Type_match_size(int typeclass,int size,MPI_Datatype *datatype){
2984   NOT_YET_IMPLEMENTED
2985 }
2986
2987 int PMPI_Alltoallw( void *sendbuf, int *sendcnts, int *sdispls, MPI_Datatype *sendtypes,
2988                     void *recvbuf, int *recvcnts, int *rdispls, MPI_Datatype *recvtypes,
2989                     MPI_Comm comm){
2990   NOT_YET_IMPLEMENTED
2991 }
2992
2993 int PMPI_Comm_set_name(MPI_Comm comm, char* name){
2994   NOT_YET_IMPLEMENTED
2995 }
2996
2997 int PMPI_Comm_dup_with_info(MPI_Comm comm, MPI_Info info, MPI_Comm * newcomm){
2998   NOT_YET_IMPLEMENTED
2999 }
3000
3001 int PMPI_Comm_split_type(MPI_Comm comm, int split_type, int key, MPI_Info info, MPI_Comm *newcomm){
3002   NOT_YET_IMPLEMENTED
3003 }
3004
3005 int PMPI_Comm_set_info (MPI_Comm comm, MPI_Info info){
3006   NOT_YET_IMPLEMENTED
3007 }
3008
3009 int PMPI_Comm_get_info (MPI_Comm comm, MPI_Info* info){
3010   NOT_YET_IMPLEMENTED
3011 }
3012
3013 int PMPI_Info_get(MPI_Info info,char *key,int valuelen, char *value, int *flag){
3014   NOT_YET_IMPLEMENTED
3015 }
3016
3017 int PMPI_Comm_create_errhandler( MPI_Comm_errhandler_fn *function, MPI_Errhandler *errhandler){
3018   NOT_YET_IMPLEMENTED
3019 }
3020
3021 int PMPI_Add_error_class( int *errorclass){
3022   NOT_YET_IMPLEMENTED
3023 }
3024
3025 int PMPI_Add_error_code(  int errorclass, int *errorcode){
3026   NOT_YET_IMPLEMENTED
3027 }
3028
3029 int PMPI_Add_error_string( int errorcode, char *string){
3030   NOT_YET_IMPLEMENTED
3031 }
3032
3033 int PMPI_Comm_call_errhandler(MPI_Comm comm,int errorcode){
3034   NOT_YET_IMPLEMENTED
3035 }
3036
3037 int PMPI_Info_dup(MPI_Info info, MPI_Info *newinfo){
3038   NOT_YET_IMPLEMENTED
3039 }
3040
3041 int PMPI_Info_delete(MPI_Info info, char *key){
3042   NOT_YET_IMPLEMENTED
3043 }
3044
3045 int PMPI_Info_get_nkeys( MPI_Info info, int *nkeys){
3046   NOT_YET_IMPLEMENTED
3047 }
3048
3049 int PMPI_Info_get_nthkey( MPI_Info info, int n, char *key){
3050   NOT_YET_IMPLEMENTED
3051 }
3052
3053 int PMPI_Info_get_valuelen( MPI_Info info, char *key, int *valuelen, int *flag){
3054   NOT_YET_IMPLEMENTED
3055 }
3056
3057 int PMPI_Request_get_status( MPI_Request request, int *flag, MPI_Status *status){
3058   NOT_YET_IMPLEMENTED
3059 }
3060
3061 int PMPI_Grequest_start( MPI_Grequest_query_function *query_fn, MPI_Grequest_free_function *free_fn, MPI_Grequest_cancel_function *cancel_fn, void *extra_state, MPI_Request *request){
3062   NOT_YET_IMPLEMENTED
3063 }
3064
3065 int PMPI_Grequest_complete( MPI_Request request){
3066   NOT_YET_IMPLEMENTED
3067 }
3068
3069 int PMPI_Status_set_cancelled(MPI_Status *status,int flag){
3070   NOT_YET_IMPLEMENTED
3071 }
3072
3073 int PMPI_Status_set_elements( MPI_Status *status, MPI_Datatype datatype, int count){
3074   NOT_YET_IMPLEMENTED
3075 }
3076
3077 int PMPI_Comm_connect( char *port_name, MPI_Info info, int root, MPI_Comm comm, MPI_Comm *newcomm){
3078   NOT_YET_IMPLEMENTED
3079 }
3080
3081 int PMPI_Publish_name( char *service_name, MPI_Info info, char *port_name){
3082   NOT_YET_IMPLEMENTED
3083 }
3084
3085 int PMPI_Unpublish_name( char *service_name, MPI_Info info, char *port_name){
3086   NOT_YET_IMPLEMENTED
3087 }
3088
3089 int PMPI_Lookup_name( char *service_name, MPI_Info info, char *port_name){
3090   NOT_YET_IMPLEMENTED
3091 }
3092
3093 int PMPI_Comm_join( int fd, MPI_Comm *intercomm){
3094   NOT_YET_IMPLEMENTED
3095 }
3096
3097 int PMPI_Open_port( MPI_Info info, char *port_name){
3098   NOT_YET_IMPLEMENTED
3099 }
3100
3101 int PMPI_Close_port(char *port_name){
3102   NOT_YET_IMPLEMENTED
3103 }
3104
3105 int PMPI_Comm_accept( char *port_name, MPI_Info info, int root, MPI_Comm comm, MPI_Comm *newcomm){
3106   NOT_YET_IMPLEMENTED
3107 }
3108
3109 int PMPI_Comm_spawn( char *command, char **argv, int maxprocs, MPI_Info info, int root, MPI_Comm comm, MPI_Comm *intercomm, int* array_of_errcodes){
3110   NOT_YET_IMPLEMENTED
3111 }
3112
3113 int PMPI_Comm_spawn_multiple( int count, char **array_of_commands, char*** array_of_argv,
3114                               int* array_of_maxprocs, MPI_Info* array_of_info, int root,
3115                               MPI_Comm comm, MPI_Comm *intercomm, int* array_of_errcodes){
3116   NOT_YET_IMPLEMENTED
3117 }
3118
3119 int PMPI_Comm_get_parent( MPI_Comm *parent){
3120   NOT_YET_IMPLEMENTED
3121 }