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.
14 The default behavior of the test routines should be to briefly indicate
15 the cause of any errors - in this test, that means that verbose needs
16 to be set. Verbose should turn on output that is independent of error
19 static int verbose = 1;
21 int parse_args(int argc, char **argv);
22 int struct_negdisp_test(void);
23 int vector_negstride_test(void);
24 int indexed_negdisp_test(void);
25 int struct_struct_test(void);
26 int flatten_test(void);
28 int build_array_section_type(MPI_Aint aext, MPI_Aint astart, MPI_Aint aend, MPI_Datatype *datatype);
30 int main(int argc, char *argv[])
35 MPI_Init(&argc, &argv);
36 parse_args(argc, argv);
38 /* To improve reporting of problems about operations, we
39 change the error handler to errors return */
40 MPI_Comm_set_errhandler( MPI_COMM_WORLD, MPI_ERRORS_RETURN );
42 err = struct_negdisp_test();
43 if (verbose && err) fprintf(stderr, "error in struct_negdisp_test\n");
46 err = vector_negstride_test();
47 if (verbose && err) fprintf(stderr, "error in vector_negstride_test\n");
50 err = indexed_negdisp_test();
51 if (verbose && err) fprintf(stderr, "error in indexed_negdisp_test\n");
54 err = struct_struct_test();
55 if (verbose && err) fprintf(stderr, "error in struct_struct_test\n");
59 if (verbose && err) fprintf(stderr, "error in flatten_test\n");
62 /* print message and exit */
64 fprintf(stderr, "Found %d errors\n", errs);
67 printf(" No Errors\n");
73 /* test uses a struct type that describes data that is contiguous,
74 * but processed in a noncontiguous way.
76 int struct_negdisp_test(void)
79 int sendbuf[6] = { 1, 2, 3, 4, 5, 6 };
80 int recvbuf[6] = { -1, -2, -3, -4, -5, -6 };
81 MPI_Datatype mystruct;
85 MPI_Aint disps[2] = { 0, -1*((int) sizeof(int)) };
86 int blks[2] = { 1, 1, };
87 MPI_Datatype types[2] = { MPI_INT, MPI_INT };
89 err = MPI_Type_struct(2, blks, disps, types, &mystruct);
90 if (err != MPI_SUCCESS) {
93 fprintf(stderr, "MPI_Type_struct returned error\n");
97 MPI_Type_commit(&mystruct);
99 err = MPI_Irecv(recvbuf+1, 4, MPI_INT, 0, 0, MPI_COMM_SELF, &request);
100 if (err != MPI_SUCCESS) {
103 fprintf(stderr, "MPI_Irecv returned error\n");
107 err = MPI_Send(sendbuf+2, 2, mystruct, 0, 0, MPI_COMM_SELF);
108 if (err != MPI_SUCCESS) {
111 fprintf(stderr, "MPI_Send returned error\n");
115 err = MPI_Wait(&request, &status);
116 if (err != MPI_SUCCESS) {
119 fprintf(stderr, "MPI_Wait returned error\n");
124 if (recvbuf[0] != -1) {
127 fprintf(stderr, "recvbuf[0] = %d; should be %d\n", recvbuf[0], -1);
130 if (recvbuf[1] != 3) {
133 fprintf(stderr, "recvbuf[1] = %d; should be %d\n", recvbuf[1], 3);
136 if (recvbuf[2] != 2) {
139 fprintf(stderr, "recvbuf[2] = %d; should be %d\n", recvbuf[2], 2);
142 if (recvbuf[3] != 5) {
145 fprintf(stderr, "recvbuf[3] = %d; should be %d\n", recvbuf[3], 5);
148 if (recvbuf[4] != 4) {
151 fprintf(stderr, "recvbuf[4] = %d; should be %d\n", recvbuf[4], 4);
154 if (recvbuf[5] != -6) {
157 fprintf(stderr, "recvbuf[5] = %d; should be %d\n", recvbuf[5], -6);
161 MPI_Type_free(&mystruct);
166 /* test uses a vector type that describes data that is contiguous,
167 * but processed in a noncontiguous way. this is effectively the
168 * same type as in the struct_negdisp_test above.
170 int vector_negstride_test(void)
173 int sendbuf[6] = { 1, 2, 3, 4, 5, 6 };
174 int recvbuf[6] = { -1, -2, -3, -4, -5, -6 };
175 MPI_Datatype myvector;
179 err = MPI_Type_vector(2, 1, -1, MPI_INT, &myvector);
180 if (err != MPI_SUCCESS) {
183 fprintf(stderr, "MPI_Type_vector returned error\n");
187 MPI_Type_commit(&myvector);
189 err = MPI_Irecv(recvbuf+1, 4, MPI_INT, 0, 0, MPI_COMM_SELF, &request);
190 if (err != MPI_SUCCESS) {
193 fprintf(stderr, "MPI_Irecv returned error\n");
197 err = MPI_Send(sendbuf+2, 2, myvector, 0, 0, MPI_COMM_SELF);
198 if (err != MPI_SUCCESS) {
201 fprintf(stderr, "MPI_Send returned error\n");
205 err = MPI_Wait(&request, &status);
206 if (err != MPI_SUCCESS) {
209 fprintf(stderr, "MPI_Wait returned error\n");
214 if (recvbuf[0] != -1) {
217 fprintf(stderr, "recvbuf[0] = %d; should be %d\n", recvbuf[0], -1);
220 if (recvbuf[1] != 3) {
223 fprintf(stderr, "recvbuf[1] = %d; should be %d\n", recvbuf[1], 3);
226 if (recvbuf[2] != 2) {
229 fprintf(stderr, "recvbuf[2] = %d; should be %d\n", recvbuf[2], 2);
232 if (recvbuf[3] != 5) {
235 fprintf(stderr, "recvbuf[3] = %d; should be %d\n", recvbuf[3], 5);
238 if (recvbuf[4] != 4) {
241 fprintf(stderr, "recvbuf[4] = %d; should be %d\n", recvbuf[4], 4);
244 if (recvbuf[5] != -6) {
247 fprintf(stderr, "recvbuf[5] = %d; should be %d\n", recvbuf[5], -6);
251 MPI_Type_free(&myvector);
256 /* test uses a indexed type that describes data that is contiguous,
257 * but processed in a noncontiguous way. this is effectively the same
258 * type as in the two tests above.
260 int indexed_negdisp_test(void)
263 int sendbuf[6] = { 1, 2, 3, 4, 5, 6 };
264 int recvbuf[6] = { -1, -2, -3, -4, -5, -6 };
265 MPI_Datatype myindexed;
269 int disps[2] = { 0, -1 };
270 int blks[2] = { 1, 1 };
272 err = MPI_Type_indexed(2, blks, disps, MPI_INT, &myindexed);
273 if (err != MPI_SUCCESS) {
276 fprintf(stderr, "MPI_Type_indexed returned error\n");
280 MPI_Type_commit(&myindexed);
282 err = MPI_Irecv(recvbuf+1, 4, MPI_INT, 0, 0, MPI_COMM_SELF, &request);
283 if (err != MPI_SUCCESS) {
286 fprintf(stderr, "MPI_Irecv returned error\n");
290 err = MPI_Send(sendbuf+2, 2, myindexed, 0, 0, MPI_COMM_SELF);
291 if (err != MPI_SUCCESS) {
294 fprintf(stderr, "MPI_Send returned error\n");
298 err = MPI_Wait(&request, &status);
299 if (err != MPI_SUCCESS) {
302 fprintf(stderr, "MPI_Wait returned error\n");
307 if (recvbuf[0] != -1) {
310 fprintf(stderr, "recvbuf[0] = %d; should be %d\n", recvbuf[0], -1);
313 if (recvbuf[1] != 3) {
316 fprintf(stderr, "recvbuf[1] = %d; should be %d\n", recvbuf[1], 3);
319 if (recvbuf[2] != 2) {
322 fprintf(stderr, "recvbuf[2] = %d; should be %d\n", recvbuf[2], 2);
325 if (recvbuf[3] != 5) {
328 fprintf(stderr, "recvbuf[3] = %d; should be %d\n", recvbuf[3], 5);
331 if (recvbuf[4] != 4) {
334 fprintf(stderr, "recvbuf[4] = %d; should be %d\n", recvbuf[4], 4);
337 if (recvbuf[5] != -6) {
340 fprintf(stderr, "recvbuf[5] = %d; should be %d\n", recvbuf[5], -6);
344 MPI_Type_free(&myindexed);
349 #define check_err(fn_name_) \
351 if (err != MPI_SUCCESS) { \
355 char err_str_[MPI_MAX_ERROR_STRING]; \
356 MPI_Error_string(err, err_str_, &len_); \
357 fprintf(stderr, #fn_name_ " failed at line %d, err=%d: %s\n", \
358 __LINE__, err, err_str_); \
362 /* test case from tt#1030 ported to C
364 * Thanks to Matthias Lieber for reporting the bug and providing a good test
366 int struct_struct_test(void)
369 int i, j, dt_size = 0;
374 MPI_Aint displ[COUNT];
376 MPI_Datatype types[COUNT];
377 MPI_Datatype datatype;
379 /* A slight difference from the F90 test: F90 arrays are column-major, C
380 * arrays are row-major. So we invert the order of dimensions. */
383 int array[N][M] = { {-1, -1, -1, -1}, {-1, -1, -1, -1} };
384 int expected[N][M] = { {-1, 1, 2, 5}, {-1, 3, 4, 6} };
386 MPI_Aint astart, aend;
387 MPI_Aint size_exp = 0;
389 /* 1st section selects elements 1 and 2 out of 2nd dimension, complete 1st dim.
390 * should receive the values 1, 2, 3, 4 */
393 err = build_array_section_type(M, astart, aend, &types[0]);
396 if (verbose) fprintf(stderr, "build_array_section_type failed\n");
401 size_exp = size_exp + N * (aend-astart+1) * sizeof(int);
403 /* 2nd section selects last element of 2nd dimension, complete 1st dim.
404 * should receive the values 5, 6 */
407 err = build_array_section_type(M, astart, aend, &types[1]);
410 if (verbose) fprintf(stderr, "build_array_section_type failed\n");
415 size_exp = size_exp + N * (aend-astart+1) * sizeof(int);
418 err = MPI_Type_create_struct(COUNT, blens, displ, types, &datatype);
419 check_err(MPI_Type_create_struct);
420 err = MPI_Type_commit(&datatype);
421 check_err(MPI_Type_commit);
423 err = MPI_Type_size(datatype, &dt_size);
424 check_err(MPI_Type_size);
425 if (dt_size != size_exp) {
427 if (verbose) fprintf(stderr, "unexpected type size\n");
431 /* send the type to ourselves to make sure that the type describes data correctly */
432 for (i = 0; i < (N*M) ; ++i)
433 seq_array[i] = i + 1; /* source values 1..(N*M) */
434 err = MPI_Isend(&seq_array[0], dt_size/sizeof(int), MPI_INT, 0, 42, MPI_COMM_SELF, &req[0]);
435 check_err(MPI_Isend);
436 err = MPI_Irecv(&array[0][0], 1, datatype, 0, 42, MPI_COMM_SELF, &req[1]);
437 check_err(MPI_Irecv);
438 err = MPI_Waitall(2, req, MPI_STATUSES_IGNORE);
439 check_err(MPI_Waitall);
441 /* check against expected */
442 for (i = 0; i < N; ++i) {
443 for (j = 0; j < M; ++j) {
444 if (array[i][j] != expected[i][j]) {
447 fprintf(stderr, "array[%d][%d]=%d, should be %d\n", i, j, array[i][j], expected[i][j]);
452 err = MPI_Type_free(&datatype);
453 check_err(MPI_Type_free);
454 err = MPI_Type_free(&types[0]);
455 check_err(MPI_Type_free);
456 err = MPI_Type_free(&types[1]);
457 check_err(MPI_Type_free);
465 /* create a datatype for a 1D int array subsection
467 - a subsection of the first dimension is defined via astart, aend
468 - indexes are assumed to start with 0, that means:
469 - 0 <= astart <= aend < aext
470 - astart and aend are inclusive
474 aext = 8, astart=2, aend=4 would produce:
476 index | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
477 1D array ###############################
478 datatype LB ########### UB
480 int build_array_section_type(MPI_Aint aext, MPI_Aint astart, MPI_Aint aend, MPI_Datatype *datatype)
484 MPI_Aint displ[COUNT];
486 MPI_Datatype types[COUNT];
488 *datatype = MPI_DATATYPE_NULL;
490 /* lower bound marker */
495 /* subsection starting at astart */
496 displ[1] = astart * sizeof(int);
498 blens[1] = aend - astart + 1;
500 /* upper bound marker */
502 displ[2] = aext * sizeof(int);
505 err = MPI_Type_create_struct(COUNT, blens, displ, types, datatype);
506 if (err != MPI_SUCCESS) {
509 fprintf(stderr, "MPI_Type_create_struct failed, err=%d\n", err);
517 /* start_idx is the "zero" point for the unpack */
518 static int pack_and_check_expected(MPI_Datatype type, const char *name,
519 int start_idx, int size,
520 int *array, int *expected)
525 int *pack_buf = NULL;
528 int sendbuf[8] = {0,1,2,3,4,5,6,7};
530 err = MPI_Type_size(type, &type_size);
531 check_err(MPI_Type_size);
532 assert(sizeof(sendbuf) >= type_size);
534 err = MPI_Pack_size(type_size/sizeof(int), MPI_INT, MPI_COMM_SELF, &pack_size);
535 check_err(MPI_Pack_size);
536 pack_buf = malloc(pack_size);
540 err = MPI_Pack(&sendbuf[0], type_size/sizeof(int), MPI_INT, pack_buf, pack_size, &pos, MPI_COMM_SELF);
543 err = MPI_Unpack(pack_buf, pack_size, &pos, &array[start_idx], 1, type, MPI_COMM_SELF);
544 check_err(MPI_Unpack);
547 /* check against expected */
548 for (i = 0; i < size; ++i) {
549 if (array[i] != expected[i]) {
552 fprintf(stderr, "%s: array[%d]=%d, should be %d\n", name, i, array[i], expected[i]);
559 /* regression for tt#1030, checks for bad offset math in the
560 * blockindexed and indexed dataloop flattening code */
561 int flatten_test(void)
565 /* real indices 0 1 2 3 4 5 6 7 8
566 * indices w/ &array[3] -3 -2 -1 0 1 2 3 4 5 */
567 int array[ARR_SIZE] = {-1,-1,-1,-1,-1,-1,-1,-1,-1};
568 int expected[ARR_SIZE] = {-1, 0, 1,-1, 2,-1, 3,-1, 4};
569 MPI_Datatype idx_type = MPI_DATATYPE_NULL;
570 MPI_Datatype blkidx_type = MPI_DATATYPE_NULL;
571 MPI_Datatype combo = MPI_DATATYPE_NULL;
574 MPI_Aint adispl[COUNT];
576 MPI_Datatype types[COUNT];
578 /* indexed type layout:
580 * 2101 <-- pos (left of 0 is neg)
582 * different blens to prevent optimization into a blockindexed
585 displ[0] = -2; /* elements, puts byte after block end at 0 */
587 displ[1] = 1; /*elements*/
589 err = MPI_Type_indexed(COUNT, blens, displ, MPI_INT, &idx_type);
590 check_err(MPI_Type_indexed);
591 err = MPI_Type_commit(&idx_type);
592 check_err(MPI_Type_commit);
594 /* indexed type layout:
596 * 2101 <-- pos (left of 0 is neg)
600 err = MPI_Type_create_indexed_block(COUNT, 1, displ, MPI_INT, &blkidx_type);
601 check_err(MPI_Type_indexed_block);
602 err = MPI_Type_commit(&blkidx_type);
603 check_err(MPI_Type_commit);
605 /* struct type layout:
606 * II_I_B_B (I=idx_type, B=blkidx_type)
607 * 21012345 <-- pos (left of 0 is neg)
610 adispl[0] = 0; /*bytes*/
614 adispl[1] = 4 * sizeof(int); /* bytes */
615 types[1] = blkidx_type;
617 /* must be a struct in order to trigger flattening code */
618 err = MPI_Type_create_struct(COUNT, blens, adispl, types, &combo);
619 check_err(MPI_Type_indexed);
620 err = MPI_Type_commit(&combo);
621 check_err(MPI_Type_commit);
623 /* pack/unpack with &array[3] */
624 errs += pack_and_check_expected(combo, "combo", 3, ARR_SIZE, array, expected);
626 MPI_Type_free(&combo);
627 MPI_Type_free(&idx_type);
628 MPI_Type_free(&blkidx_type);
635 int parse_args(int argc, char **argv)
640 while ((ret = getopt(argc, argv, "v")) >= 0)
649 if (argc > 1 && strcmp(argv[1], "-v") == 0)