Logo AND Algorithmique Numérique Distribuée

Public GIT Repository
Fix dead assignments.
[simgrid.git] / src / smpi / colls / alltoall-3dmesh.c
index cae32a2..6f43a42 100644 (file)
@@ -1,4 +1,10 @@
-#include "colls.h"
+/* Copyright (c) 2013-2014. The SimGrid Team.
+ * All rights reserved.                                                     */
+
+/* This program is free software; you can redistribute it and/or modify it
+ * under the terms of the license (GNU LGPL) which comes with this package. */
+
+#include "colls_private.h"
 #include <math.h>
 
 /*****************************************************************************
@@ -50,16 +56,16 @@ int smpi_coll_tuned_alltoall_3dmesh(void *send_buff, int send_count,
   MPI_Status status, *statuses;
   int i, j, src, dst, rank, num_procs, num_reqs, X, Y, Z, block_size, count;
   int my_z, two_dsize, my_row_base, my_col_base, my_z_base, src_row_base;
-  int src_z_base, send_offset, recv_offset, tag = 1, failure = 0, success = 1;
+  int src_z_base, send_offset, recv_offset, tag = COLL_TAG_ALLTOALL;
 
   char *tmp_buff1, *tmp_buff2;
 
-  MPI_Comm_rank(comm, &rank);
-  MPI_Comm_size(comm, &num_procs);
-  MPI_Type_extent(send_type, &extent);
+  rank = smpi_comm_rank(comm);
+  num_procs = smpi_comm_size(comm);
+  extent = smpi_datatype_get_extent(send_type);
 
   if (!alltoall_check_is_3dmesh(num_procs, &X, &Y, &Z))
-    return failure;
+    return MPI_ERR_OTHER;
 
   num_reqs = X;
   if (Y > X)
@@ -76,33 +82,17 @@ int smpi_coll_tuned_alltoall_3dmesh(void *send_buff, int send_count,
 
   block_size = extent * send_count;
 
-  tmp_buff1 = (char *) malloc(block_size * num_procs * two_dsize);
-  if (!tmp_buff1) {
-    printf("alltoall-3Dmesh:97: cannot allocate memory\n");
-    MPI_Finalize();
-    exit(failure);
-  }
-
-  tmp_buff2 = (char *) malloc(block_size * two_dsize);
-  if (!tmp_buff2) {
-    printf("alltoall-3Dmesh:105: cannot allocate memory\n");
-    MPI_Finalize();
-    exit(failure);
-  }
+  tmp_buff1 = (char *) xbt_malloc(block_size * num_procs * two_dsize);
+  tmp_buff2 = (char *) xbt_malloc(block_size * two_dsize);
 
-  statuses = (MPI_Status *) malloc(num_reqs * sizeof(MPI_Status));
-  reqs = (MPI_Request *) malloc(num_reqs * sizeof(MPI_Request));
-  if (!reqs) {
-    printf("alltoall-3Dmesh:113: cannot allocate memory\n");
-    MPI_Finalize();
-    exit(failure);
-  }
+  statuses = (MPI_Status *) xbt_malloc(num_reqs * sizeof(MPI_Status));
+  reqs = (MPI_Request *) xbt_malloc(num_reqs * sizeof(MPI_Request));
 
   req_ptr = reqs;
 
-  send_offset = recv_offset = (rank % two_dsize) * block_size * num_procs;
+  recv_offset = (rank % two_dsize) * block_size * num_procs;
 
-  MPI_Sendrecv(send_buff, send_count * num_procs, send_type, rank, tag,
+  smpi_mpi_sendrecv(send_buff, send_count * num_procs, send_type, rank, tag,
                tmp_buff1 + recv_offset, num_procs * recv_count,
                recv_type, rank, tag, comm, &status);
 
@@ -113,18 +103,17 @@ int smpi_coll_tuned_alltoall_3dmesh(void *send_buff, int send_count,
     if (src == rank)
       continue;
     recv_offset = (src % two_dsize) * block_size * num_procs;
-    MPI_Irecv(tmp_buff1 + recv_offset, count, recv_type, src, tag, comm,
-              req_ptr++);
+    *(req_ptr++) = smpi_mpi_irecv(tmp_buff1 + recv_offset, count, recv_type, src, tag, comm);
   }
 
   for (i = 0; i < Y; i++) {
     dst = i + my_row_base;
     if (dst == rank)
       continue;
-    MPI_Send(send_buff, count, send_type, dst, tag, comm);
+    smpi_mpi_send(send_buff, count, send_type, dst, tag, comm);
   }
 
-  MPI_Waitall(Y - 1, reqs, statuses);
+  smpi_mpi_waitall(Y - 1, reqs, statuses);
   req_ptr = reqs;
 
 
@@ -136,8 +125,8 @@ int smpi_coll_tuned_alltoall_3dmesh(void *send_buff, int send_count,
     src_row_base = (src / X) * X;
 
     recv_offset = (src_row_base % two_dsize) * block_size * num_procs;
-    MPI_Irecv(tmp_buff1 + recv_offset, recv_count * num_procs * Y,
-              recv_type, src, tag, comm, req_ptr++);
+    *(req_ptr++) = smpi_mpi_irecv(tmp_buff1 + recv_offset, recv_count * num_procs * Y,
+              recv_type, src, tag, comm);
   }
 
   send_offset = (my_row_base % two_dsize) * block_size * num_procs;
@@ -145,17 +134,17 @@ int smpi_coll_tuned_alltoall_3dmesh(void *send_buff, int send_count,
     dst = (i * Y + my_col_base);
     if (dst == rank)
       continue;
-    MPI_Send(tmp_buff1 + send_offset, send_count * num_procs * Y, send_type,
+    smpi_mpi_send(tmp_buff1 + send_offset, send_count * num_procs * Y, send_type,
              dst, tag, comm);
   }
 
-  MPI_Waitall(X - 1, reqs, statuses);
+  smpi_mpi_waitall(X - 1, reqs, statuses);
   req_ptr = reqs;
 
   for (i = 0; i < two_dsize; i++) {
     send_offset = (rank * block_size) + (i * block_size * num_procs);
     recv_offset = (my_z_base * block_size) + (i * block_size);
-    MPI_Sendrecv(tmp_buff1 + send_offset, send_count, send_type, rank, tag,
+    smpi_mpi_sendrecv(tmp_buff1 + send_offset, send_count, send_type, rank, tag,
                  (char *) recv_buff + recv_offset, recv_count, recv_type,
                  rank, tag, comm, &status);
   }
@@ -166,8 +155,8 @@ int smpi_coll_tuned_alltoall_3dmesh(void *send_buff, int send_count,
 
     recv_offset = (src_z_base * block_size);
 
-    MPI_Irecv((char *) recv_buff + recv_offset, recv_count * two_dsize,
-              recv_type, src, tag, comm, req_ptr++);
+    *(req_ptr++) = smpi_mpi_irecv((char *) recv_buff + recv_offset, recv_count * two_dsize,
+              recv_type, src, tag, comm);
   }
 
   for (i = 1; i < Z; i++) {
@@ -176,22 +165,22 @@ int smpi_coll_tuned_alltoall_3dmesh(void *send_buff, int send_count,
     recv_offset = 0;
     for (j = 0; j < two_dsize; j++) {
       send_offset = (dst + j * num_procs) * block_size;
-      MPI_Sendrecv(tmp_buff1 + send_offset, send_count, send_type,
+      smpi_mpi_sendrecv(tmp_buff1 + send_offset, send_count, send_type,
                    rank, tag, tmp_buff2 + recv_offset, recv_count,
                    recv_type, rank, tag, comm, &status);
 
       recv_offset += block_size;
     }
 
-    MPI_Send(tmp_buff2, send_count * two_dsize, send_type, dst, tag, comm);
+    smpi_mpi_send(tmp_buff2, send_count * two_dsize, send_type, dst, tag, comm);
 
   }
 
-  MPI_Waitall(Z - 1, reqs, statuses);
+  smpi_mpi_waitall(Z - 1, reqs, statuses);
 
   free(reqs);
   free(statuses);
   free(tmp_buff1);
   free(tmp_buff2);
-  return success;
+  return MPI_SUCCESS;
 }