1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
4 * (C) 2003 by Argonne National Laboratory.
5 * See COPYRIGHT in top-level directory.
13 /* This example sends contiguous data and receives a vector on some nodes
14 and contiguous data on others. There is some evidence that some
15 MPI implementations do not check recvcount on the root process; this
16 test checks for that case
19 int main( int argc, char **argv )
22 double *vecin, *vecout, ivalue;
23 int root, i, n, stride, errs = 0;
27 MTest_Init( &argc, &argv );
29 MPI_Comm_size( MPI_COMM_WORLD, &size );
30 MPI_Comm_rank( MPI_COMM_WORLD, &rank );
34 /* Note that vecout really needs to be only (n-1)*stride+1 doubles, but
35 this is easier and allows a little extra room if there is a bug */
36 vecout = (double *)malloc( n * stride * sizeof(double) );
37 vecin = (double *)malloc( n * size * sizeof(double) );
39 MPI_Type_vector( n, 1, stride, MPI_DOUBLE, &vec );
40 MPI_Type_commit( &vec );
41 MPI_Type_extent( vec, &vextent );
42 if (vextent != ((n-1)*(MPI_Aint)stride + 1) * sizeof(double) ) {
44 printf( "Vector extent is %ld, should be %ld\n",
45 (long) vextent, (long)(((n-1)*stride+1)*sizeof(double)) );
47 /* Note that the exted of type vector is from the first to the
48 last element, not n*stride.
49 E.g., with n=1, the extent is a single double */
51 for (i=0; i<n*size; i++) vecin[i] = (double)i;
52 for (root=0; root<size; root++) {
53 for (i=0; i<n*stride; i++) vecout[i] = -1.0;
55 /* Receive into a vector */
56 MPI_Scatter( vecin, n, MPI_DOUBLE, vecout, 1, vec,
57 root, MPI_COMM_WORLD );
60 if (vecout[i*stride] != ivalue) {
62 printf( "[%d] Expected %f but found %f for vecout[%d] on root\n",
63 rank, ivalue, vecout[i*stride], i *stride );
68 /* Receive into contiguous data */
69 MPI_Scatter( NULL, -1, MPI_DATATYPE_NULL, vecout, n, MPI_DOUBLE,
70 root, MPI_COMM_WORLD );
72 ivalue = rank * n + i;
73 if (vecout[i] != ivalue) {
74 printf( "[%d] Expected %f but found %f for vecout[%d]\n",
75 rank, ivalue, vecout[i], i );
82 MTest_Finalize( errs );
83 MPI_Type_free( &vec );