1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
3 * (C) 2001 by Argonne National Laboratory.
4 * See COPYRIGHT in top-level directory.
10 /* Prototypes for picky compilers */
11 void SetData(double *, double *, int, int, int, int, int, int);
12 int CheckData(double *, int, int, int, int, int, int);
14 This is an example of using scatterv to send a matrix from one
15 process to all others, with the matrix stored in Fortran order.
16 Note the use of an explicit UB to enable the sources to overlap.
18 This tests scatterv to make sure that it uses the datatype size
19 and extent correctly. It requires number of processors that
20 can be split with MPI_Dims_create.
24 void SetData(double *sendbuf, double *recvbuf, int nx, int ny,
25 int myrow, int mycol, int nrow, int ncol)
27 int coldim, i, j, m, k;
30 if (myrow == 0 && mycol == 0) {
32 for (j = 0; j < ncol; j++) {
33 for (i = 0; i < nrow; i++) {
34 p = sendbuf + i * nx + j * (ny * coldim);
35 for (m = 0; m < ny; m++) {
36 for (k = 0; k < nx; k++) {
37 p[k] = 1000 * j + 100 * i + m * nx + k;
44 for (i = 0; i < nx * ny; i++)
48 int CheckData(double *recvbuf, int nx, int ny, int myrow, int mycol, int nrow, int expect_no_value)
56 for (m = 0; m < ny; m++) {
57 for (k = 0; k < nx; k++) {
58 /* If expect_no_value is true then we assume that the pre-scatterv
59 * value should remain in the recvbuf for our portion of the array.
60 * This is the case for the root process when using MPI_IN_PLACE. */
64 val = 1000 * mycol + 100 * myrow + m * nx + k;
69 printf("Error in (%d,%d) [%d,%d] location, got %f expected %f\n",
70 m, k, myrow, mycol, p[k], val);
72 else if (errs == 10) {
73 printf("Too many errors; suppressing printing\n");
82 int main(int argc, char **argv)
84 int rank, size, myrow, mycol, nx, ny, stride, cnt, i, j, errs, errs_in_place, tot_errs;
85 double *sendbuf, *recvbuf;
86 MPI_Datatype vec, block, types[2];
91 int dims[2], periods[2], coords[2], lcoords[2];
95 MPI_Init(&argc, &argv);
96 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
97 MPI_Comm_size(MPI_COMM_WORLD, &size);
99 /* Get a 2-d decomposition of the processes */
102 MPI_Dims_create(size, 2, dims);
105 MPI_Cart_create(MPI_COMM_WORLD, 2, dims, periods, 0, &comm2d);
106 MPI_Cart_get(comm2d, 2, dims, periods, coords);
111 printf("Decomposition is [%d x %d]\n", dims[0], dims[1]);
114 /* Get the size of the matrix */
117 stride = nx * dims[0];
119 recvbuf = (double *) malloc(nx * ny * sizeof(double));
121 MPI_Abort(MPI_COMM_WORLD, 1);
124 if (myrow == 0 && mycol == 0) {
125 sendbuf = (double *) malloc(nx * ny * size * sizeof(double));
127 MPI_Abort(MPI_COMM_WORLD, 1);
130 sendcounts = (int *) malloc(size * sizeof(int));
131 scdispls = (int *) malloc(size * sizeof(int));
133 MPI_Type_vector(ny, nx, stride, MPI_DOUBLE, &vec);
139 displs[1] = nx * sizeof(double);
141 MPI_Type_struct(2, blens, displs, types, &block);
143 MPI_Type_commit(&block);
145 /* Set up the transfer */
147 for (i = 0; i < dims[1]; i++) {
148 for (j = 0; j < dims[0]; j++) {
150 /* Using Cart_coords makes sure that ranks (used by
151 * sendrecv) matches the cartesian coordinates (used to
152 * set data in the matrix) */
153 MPI_Cart_coords(comm2d, cnt, 2, lcoords);
154 scdispls[cnt++] = lcoords[0] + lcoords[1] * (dims[0] * ny);
158 SetData(sendbuf, recvbuf, nx, ny, myrow, mycol, dims[0], dims[1]);
159 MPI_Scatterv(sendbuf, sendcounts, scdispls, block, recvbuf, nx * ny, MPI_DOUBLE, 0, comm2d);
160 if ((errs = CheckData(recvbuf, nx, ny, myrow, mycol, dims[0], 0))) {
161 fprintf(stdout, "Failed to transfer data\n");
164 /* once more, but this time passing MPI_IN_PLACE for the root */
165 SetData(sendbuf, recvbuf, nx, ny, myrow, mycol, dims[0], dims[1]);
166 MPI_Scatterv(sendbuf, sendcounts, scdispls, block,
167 (rank == 0 ? MPI_IN_PLACE : recvbuf), nx * ny, MPI_DOUBLE, 0, comm2d);
168 errs_in_place = CheckData(recvbuf, nx, ny, myrow, mycol, dims[0], (rank == 0));
170 fprintf(stdout, "Failed to transfer data (MPI_IN_PLACE)\n");
173 errs += errs_in_place;
174 MPI_Allreduce(&errs, &tot_errs, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD);
177 printf(" No Errors\n");
179 printf("%d errors in use of MPI_SCATTERV\n", tot_errs);
187 MPI_Type_free(&block);
188 MPI_Comm_free(&comm2d);