Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
Add missing dependency for smpi_traced_simple.
[simgrid.git] / examples / smpi / scatter.c
index 9647aee..e4a6303 100644 (file)
 #include <stdio.h>
 #include <mpi.h>
 
-int ibm_test(int rank, int size)  {
+static int ibm_test(int rank, int size)
+{
 
 #define MAXLEN  10000
 
-         int success = 1;
-         int root,i,j,k;
-         int out[MAXLEN*64];
-         int in[MAXLEN];
+  int success = 1;
+  int root = 0;
+  int i, j, k;
+  int *out;
+  int *in;
+
 
+  out = malloc(MAXLEN * 64 * sizeof(int));
+  in = malloc(MAXLEN * sizeof(int));
 
-         for(j=1,root=0;j<=MAXLEN;j*=10,root=(root+1)%size)  {
-                   if(rank == root)
-                               for(i=0;i<j*size;i++)  out[i] = i;
+  for (j = 1; j <= MAXLEN; j *= 10) {
+    root = (root + 1) % size;
+    if (rank == root)
+      for (i = 0; i < j * size; i++)
+        out[i] = i;
 
-                   MPI_Scatter(out,j,MPI_INT,in,j,MPI_INT,root,MPI_COMM_WORLD);
+    MPI_Scatter(out, j, MPI_INT, in, j, MPI_INT, root, MPI_COMM_WORLD);
 
-                   for(k=0;k<j;k++) {
-                               if(in[k] != k+rank*j) {
-                                         fprintf(stderr,"task %d bad answer (%d) at index %d k of %d (should be %d)", rank,in[k],k,j,(k+rank*j));
-                                         return( 0 ); 
-                               }
-                   }
-         }
-        return( success );
-         MPI_Barrier( MPI_COMM_WORLD );
+    for (k = 0; k < j; k++) {
+      if (in[k] != k + rank * j) {
+        fprintf(stderr,
+                "task %d bad answer (%d) at index %d k of %d (should be %d)",
+                rank, in[k], k, j, (k + rank * j));
+        return (0);
+      }
+    }
+  }
+  free(out);
+  free(in);
+  MPI_Barrier(MPI_COMM_WORLD);
+  return (success);
 }
 
 /**
  * small test: the root sends a single distinct double to other processes
  **/
-int small_test( int rank, int size ) {
-  int success=1;
-  int sendcount=size;
-  int recvcount=1;
+static int small_test(int rank, int size)
+{
+  int success = 1;
+  int retval;
+  int sendcount = 1;            // one double to each process
+  int recvcount = 1;
   int i;
-  double *sndbuf =NULL;
+  double *sndbuf = NULL;
   double rcvd;
-  int root=0; // arbitrary choice 
+  int root = 0;                 // arbitrary choice 
 
   // on root, initialize sendbuf
-  if (root == rank )  {
-           sndbuf = malloc( sendcount*sizeof(double));
-           for (i=0; i<sendcount;i++) {
-                       sndbuf[i] = (double) rank;
-           }
+  if (root == rank) {
+    sndbuf = malloc(size * sizeof(double));
+    for (i = 0; i < size; i++) {
+      sndbuf[i] = (double) i;
+    }
   }
 
-   MPI_Scatter(sndbuf, sendcount, MPI_DOUBLE, &rcvd,recvcount,MPI_DOUBLE,root,MPI_COMM_WORLD);
-
-  return(success);
+  retval =
+      MPI_Scatter(sndbuf, sendcount, MPI_DOUBLE, &rcvd, recvcount,
+                  MPI_DOUBLE, root, MPI_COMM_WORLD);
+  if (retval != MPI_SUCCESS) {
+    fprintf(stderr, "(%s:%d) MPI_Scatter() returned retval=%d\n", __FILE__,
+            __LINE__, retval);
+    return 0;
+  }
+  // verification
+  if ((double) rank != rcvd) {
+    fprintf(stderr, "[%d] has %lf instead of %d\n", rank, rcvd, rank);
+    success = 0;
+  }
+  return (success);
 }
 
 
@@ -87,20 +111,25 @@ int main(int argc, char **argv)
   MPI_Comm_size(MPI_COMM_WORLD, &size);
   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
 
-
-
+  /* test 1 */
   if (0 == rank)
     printf("** Small Test Result: ... \n");
-  small_test( rank, size );
+  if (!small_test(rank, size))
+    printf("\t[%d] failed.\n", rank);
+  else
+    printf("\t[%d] ok.\n", rank);
+
 
+  MPI_Barrier(MPI_COMM_WORLD);
 
+  /* test 2 */
   if (0 == rank)
     printf("** IBM Test Result: ... \n");
-/*  if (!ibm_test(rank, size))
+  if (!ibm_test(rank, size))
     printf("\t[%d] failed.\n", rank);
   else
     printf("\t[%d] ok.\n", rank);
-*/
+
   MPI_Finalize();
   return 0;
 }