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 main(int argc, char *argv[]);
22 int parse_args(int argc, char **argv);
23 int struct_negdisp_test(void);
24 int vector_negstride_test(void);
25 int indexed_negdisp_test(void);
26 int struct_struct_test(void);
27 int flatten_test(void);
29 int build_array_section_type(MPI_Aint aext, MPI_Aint astart, MPI_Aint aend, MPI_Datatype *datatype);
31 int main(int argc, char *argv[])
36 MPI_Init(&argc, &argv);
37 parse_args(argc, argv);
39 /* To improve reporting of problems about operations, we
40 change the error handler to errors return */
41 MPI_Comm_set_errhandler( MPI_COMM_WORLD, MPI_ERRORS_RETURN );
43 err = struct_negdisp_test();
44 if (verbose && err) fprintf(stderr, "error in struct_negdisp_test\n");
47 err = vector_negstride_test();
48 if (verbose && err) fprintf(stderr, "error in vector_negstride_test\n");
51 err = indexed_negdisp_test();
52 if (verbose && err) fprintf(stderr, "error in indexed_negdisp_test\n");
55 err = struct_struct_test();
56 if (verbose && err) fprintf(stderr, "error in struct_struct_test\n");
60 if (verbose && err) fprintf(stderr, "error in flatten_test\n");
63 /* print message and exit */
65 fprintf(stderr, "Found %d errors\n", errs);
68 printf(" No Errors\n");
74 /* test uses a struct type that describes data that is contiguous,
75 * but processed in a noncontiguous way.
77 int struct_negdisp_test(void)
80 int sendbuf[6] = { 1, 2, 3, 4, 5, 6 };
81 int recvbuf[6] = { -1, -2, -3, -4, -5, -6 };
82 MPI_Datatype mystruct;
86 MPI_Aint disps[2] = { 0, -1*((int) sizeof(int)) };
87 int blks[2] = { 1, 1, };
88 MPI_Datatype types[2] = { MPI_INT, MPI_INT };
90 err = MPI_Type_struct(2, blks, disps, types, &mystruct);
91 if (err != MPI_SUCCESS) {
94 fprintf(stderr, "MPI_Type_struct returned error\n");
98 MPI_Type_commit(&mystruct);
100 err = MPI_Irecv(recvbuf+1, 4, MPI_INT, 0, 0, MPI_COMM_SELF, &request);
101 if (err != MPI_SUCCESS) {
104 fprintf(stderr, "MPI_Irecv returned error\n");
108 err = MPI_Send(sendbuf+2, 2, mystruct, 0, 0, MPI_COMM_SELF);
109 if (err != MPI_SUCCESS) {
112 fprintf(stderr, "MPI_Send returned error\n");
116 err = MPI_Wait(&request, &status);
117 if (err != MPI_SUCCESS) {
120 fprintf(stderr, "MPI_Wait returned error\n");
125 if (recvbuf[0] != -1) {
128 fprintf(stderr, "recvbuf[0] = %d; should be %d\n", recvbuf[0], -1);
131 if (recvbuf[1] != 3) {
134 fprintf(stderr, "recvbuf[1] = %d; should be %d\n", recvbuf[1], 3);
137 if (recvbuf[2] != 2) {
140 fprintf(stderr, "recvbuf[2] = %d; should be %d\n", recvbuf[2], 2);
143 if (recvbuf[3] != 5) {
146 fprintf(stderr, "recvbuf[3] = %d; should be %d\n", recvbuf[3], 5);
149 if (recvbuf[4] != 4) {
152 fprintf(stderr, "recvbuf[4] = %d; should be %d\n", recvbuf[4], 4);
155 if (recvbuf[5] != -6) {
158 fprintf(stderr, "recvbuf[5] = %d; should be %d\n", recvbuf[5], -6);
162 MPI_Type_free(&mystruct);
167 /* test uses a vector type that describes data that is contiguous,
168 * but processed in a noncontiguous way. this is effectively the
169 * same type as in the struct_negdisp_test above.
171 int vector_negstride_test(void)
174 int sendbuf[6] = { 1, 2, 3, 4, 5, 6 };
175 int recvbuf[6] = { -1, -2, -3, -4, -5, -6 };
176 MPI_Datatype myvector;
180 err = MPI_Type_vector(2, 1, -1, MPI_INT, &myvector);
181 if (err != MPI_SUCCESS) {
184 fprintf(stderr, "MPI_Type_vector returned error\n");
188 MPI_Type_commit(&myvector);
190 err = MPI_Irecv(recvbuf+1, 4, MPI_INT, 0, 0, MPI_COMM_SELF, &request);
191 if (err != MPI_SUCCESS) {
194 fprintf(stderr, "MPI_Irecv returned error\n");
198 err = MPI_Send(sendbuf+2, 2, myvector, 0, 0, MPI_COMM_SELF);
199 if (err != MPI_SUCCESS) {
202 fprintf(stderr, "MPI_Send returned error\n");
206 err = MPI_Wait(&request, &status);
207 if (err != MPI_SUCCESS) {
210 fprintf(stderr, "MPI_Wait returned error\n");
215 if (recvbuf[0] != -1) {
218 fprintf(stderr, "recvbuf[0] = %d; should be %d\n", recvbuf[0], -1);
221 if (recvbuf[1] != 3) {
224 fprintf(stderr, "recvbuf[1] = %d; should be %d\n", recvbuf[1], 3);
227 if (recvbuf[2] != 2) {
230 fprintf(stderr, "recvbuf[2] = %d; should be %d\n", recvbuf[2], 2);
233 if (recvbuf[3] != 5) {
236 fprintf(stderr, "recvbuf[3] = %d; should be %d\n", recvbuf[3], 5);
239 if (recvbuf[4] != 4) {
242 fprintf(stderr, "recvbuf[4] = %d; should be %d\n", recvbuf[4], 4);
245 if (recvbuf[5] != -6) {
248 fprintf(stderr, "recvbuf[5] = %d; should be %d\n", recvbuf[5], -6);
252 MPI_Type_free(&myvector);
257 /* test uses a indexed type that describes data that is contiguous,
258 * but processed in a noncontiguous way. this is effectively the same
259 * type as in the two tests above.
261 int indexed_negdisp_test(void)
264 int sendbuf[6] = { 1, 2, 3, 4, 5, 6 };
265 int recvbuf[6] = { -1, -2, -3, -4, -5, -6 };
266 MPI_Datatype myindexed;
270 int disps[2] = { 0, -1 };
271 int blks[2] = { 1, 1 };
273 err = MPI_Type_indexed(2, blks, disps, MPI_INT, &myindexed);
274 if (err != MPI_SUCCESS) {
277 fprintf(stderr, "MPI_Type_indexed returned error\n");
281 MPI_Type_commit(&myindexed);
283 err = MPI_Irecv(recvbuf+1, 4, MPI_INT, 0, 0, MPI_COMM_SELF, &request);
284 if (err != MPI_SUCCESS) {
287 fprintf(stderr, "MPI_Irecv returned error\n");
291 err = MPI_Send(sendbuf+2, 2, myindexed, 0, 0, MPI_COMM_SELF);
292 if (err != MPI_SUCCESS) {
295 fprintf(stderr, "MPI_Send returned error\n");
299 err = MPI_Wait(&request, &status);
300 if (err != MPI_SUCCESS) {
303 fprintf(stderr, "MPI_Wait returned error\n");
308 if (recvbuf[0] != -1) {
311 fprintf(stderr, "recvbuf[0] = %d; should be %d\n", recvbuf[0], -1);
314 if (recvbuf[1] != 3) {
317 fprintf(stderr, "recvbuf[1] = %d; should be %d\n", recvbuf[1], 3);
320 if (recvbuf[2] != 2) {
323 fprintf(stderr, "recvbuf[2] = %d; should be %d\n", recvbuf[2], 2);
326 if (recvbuf[3] != 5) {
329 fprintf(stderr, "recvbuf[3] = %d; should be %d\n", recvbuf[3], 5);
332 if (recvbuf[4] != 4) {
335 fprintf(stderr, "recvbuf[4] = %d; should be %d\n", recvbuf[4], 4);
338 if (recvbuf[5] != -6) {
341 fprintf(stderr, "recvbuf[5] = %d; should be %d\n", recvbuf[5], -6);
345 MPI_Type_free(&myindexed);
350 #define check_err(fn_name_) \
352 if (err != MPI_SUCCESS) { \
356 char err_str_[MPI_MAX_ERROR_STRING]; \
357 MPI_Error_string(err, err_str_, &len_); \
358 fprintf(stderr, #fn_name_ " failed at line %d, err=%d: %s\n", \
359 __LINE__, err, err_str_); \
363 /* test case from tt#1030 ported to C
365 * Thanks to Matthias Lieber for reporting the bug and providing a good test
367 int struct_struct_test(void)
370 int i, j, dt_size = 0;
375 MPI_Aint displ[COUNT];
377 MPI_Datatype types[COUNT];
378 MPI_Datatype datatype;
380 /* A slight difference from the F90 test: F90 arrays are column-major, C
381 * arrays are row-major. So we invert the order of dimensions. */
384 int array[N][M] = { {-1, -1, -1, -1}, {-1, -1, -1, -1} };
385 int expected[N][M] = { {-1, 1, 2, 5}, {-1, 3, 4, 6} };
387 MPI_Aint astart, aend;
388 MPI_Aint size_exp = 0;
390 /* 1st section selects elements 1 and 2 out of 2nd dimension, complete 1st dim.
391 * should receive the values 1, 2, 3, 4 */
394 err = build_array_section_type(M, astart, aend, &types[0]);
397 if (verbose) fprintf(stderr, "build_array_section_type failed\n");
402 size_exp = size_exp + N * (aend-astart+1) * sizeof(int);
404 /* 2nd section selects last element of 2nd dimension, complete 1st dim.
405 * should receive the values 5, 6 */
408 err = build_array_section_type(M, astart, aend, &types[1]);
411 if (verbose) fprintf(stderr, "build_array_section_type failed\n");
416 size_exp = size_exp + N * (aend-astart+1) * sizeof(int);
419 err = MPI_Type_create_struct(COUNT, blens, displ, types, &datatype);
420 check_err(MPI_Type_create_struct);
421 err = MPI_Type_commit(&datatype);
422 check_err(MPI_Type_commit);
424 err = MPI_Type_size(datatype, &dt_size);
425 check_err(MPI_Type_size);
426 if (dt_size != size_exp) {
428 if (verbose) fprintf(stderr, "unexpected type size\n");
432 /* send the type to ourselves to make sure that the type describes data correctly */
433 for (i = 0; i < (N*M) ; ++i)
434 seq_array[i] = i + 1; /* source values 1..(N*M) */
435 err = MPI_Isend(&seq_array[0], dt_size/sizeof(int), MPI_INT, 0, 42, MPI_COMM_SELF, &req[0]);
436 check_err(MPI_Isend);
437 err = MPI_Irecv(&array[0][0], 1, datatype, 0, 42, MPI_COMM_SELF, &req[1]);
438 check_err(MPI_Irecv);
439 err = MPI_Waitall(2, req, MPI_STATUSES_IGNORE);
440 check_err(MPI_Waitall);
442 /* check against expected */
443 for (i = 0; i < N; ++i) {
444 for (j = 0; j < M; ++j) {
445 if (array[i][j] != expected[i][j]) {
448 fprintf(stderr, "array[%d][%d]=%d, should be %d\n", i, j, array[i][j], expected[i][j]);
453 err = MPI_Type_free(&datatype);
454 check_err(MPI_Type_free);
455 err = MPI_Type_free(&types[0]);
456 check_err(MPI_Type_free);
457 err = MPI_Type_free(&types[1]);
458 check_err(MPI_Type_free);
466 /* create a datatype for a 1D int array subsection
468 - a subsection of the first dimension is defined via astart, aend
469 - indexes are assumed to start with 0, that means:
470 - 0 <= astart <= aend < aext
471 - astart and aend are inclusive
475 aext = 8, astart=2, aend=4 would produce:
477 index | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
478 1D array ###############################
479 datatype LB ########### UB
481 int build_array_section_type(MPI_Aint aext, MPI_Aint astart, MPI_Aint aend, MPI_Datatype *datatype)
485 MPI_Aint displ[COUNT];
487 MPI_Datatype types[COUNT];
489 *datatype = MPI_DATATYPE_NULL;
491 /* lower bound marker */
496 /* subsection starting at astart */
497 displ[1] = astart * sizeof(int);
499 blens[1] = aend - astart + 1;
501 /* upper bound marker */
503 displ[2] = aext * sizeof(int);
506 err = MPI_Type_create_struct(COUNT, blens, displ, types, datatype);
507 if (err != MPI_SUCCESS) {
510 fprintf(stderr, "MPI_Type_create_struct failed, err=%d\n", err);
518 /* start_idx is the "zero" point for the unpack */
519 static int pack_and_check_expected(MPI_Datatype type, const char *name,
520 int start_idx, int size,
521 int *array, int *expected)
526 int *pack_buf = NULL;
529 int sendbuf[8] = {0,1,2,3,4,5,6,7};
531 err = MPI_Type_size(type, &type_size);
532 check_err(MPI_Type_size);
533 assert(sizeof(sendbuf) >= type_size);
535 err = MPI_Pack_size(type_size/sizeof(int), MPI_INT, MPI_COMM_SELF, &pack_size);
536 check_err(MPI_Pack_size);
537 pack_buf = malloc(pack_size);
541 err = MPI_Pack(&sendbuf[0], type_size/sizeof(int), MPI_INT, pack_buf, pack_size, &pos, MPI_COMM_SELF);
544 err = MPI_Unpack(pack_buf, pack_size, &pos, &array[start_idx], 1, type, MPI_COMM_SELF);
545 check_err(MPI_Unpack);
548 /* check against expected */
549 for (i = 0; i < size; ++i) {
550 if (array[i] != expected[i]) {
553 fprintf(stderr, "%s: array[%d]=%d, should be %d\n", name, i, array[i], expected[i]);
560 /* regression for tt#1030, checks for bad offset math in the
561 * blockindexed and indexed dataloop flattening code */
562 int flatten_test(void)
566 /* real indices 0 1 2 3 4 5 6 7 8
567 * indices w/ &array[3] -3 -2 -1 0 1 2 3 4 5 */
568 int array[ARR_SIZE] = {-1,-1,-1,-1,-1,-1,-1,-1,-1};
569 int expected[ARR_SIZE] = {-1, 0, 1,-1, 2,-1, 3,-1, 4};
570 MPI_Datatype idx_type = MPI_DATATYPE_NULL;
571 MPI_Datatype blkidx_type = MPI_DATATYPE_NULL;
572 MPI_Datatype combo = MPI_DATATYPE_NULL;
575 MPI_Aint adispl[COUNT];
577 MPI_Datatype types[COUNT];
579 /* indexed type layout:
581 * 2101 <-- pos (left of 0 is neg)
583 * different blens to prevent optimization into a blockindexed
586 displ[0] = -2; /* elements, puts byte after block end at 0 */
588 displ[1] = 1; /*elements*/
590 err = MPI_Type_indexed(COUNT, blens, displ, MPI_INT, &idx_type);
591 check_err(MPI_Type_indexed);
592 err = MPI_Type_commit(&idx_type);
593 check_err(MPI_Type_commit);
595 /* indexed type layout:
597 * 2101 <-- pos (left of 0 is neg)
601 err = MPI_Type_create_indexed_block(COUNT, 1, displ, MPI_INT, &blkidx_type);
602 check_err(MPI_Type_indexed_block);
603 err = MPI_Type_commit(&blkidx_type);
604 check_err(MPI_Type_commit);
606 /* struct type layout:
607 * II_I_B_B (I=idx_type, B=blkidx_type)
608 * 21012345 <-- pos (left of 0 is neg)
611 adispl[0] = 0; /*bytes*/
615 adispl[1] = 4 * sizeof(int); /* bytes */
616 types[1] = blkidx_type;
618 /* must be a struct in order to trigger flattening code */
619 err = MPI_Type_create_struct(COUNT, blens, adispl, types, &combo);
620 check_err(MPI_Type_indexed);
621 err = MPI_Type_commit(&combo);
622 check_err(MPI_Type_commit);
624 /* pack/unpack with &array[3] */
625 errs += pack_and_check_expected(combo, "combo", 3, ARR_SIZE, array, expected);
627 MPI_Type_free(&combo);
628 MPI_Type_free(&idx_type);
629 MPI_Type_free(&blkidx_type);
636 int parse_args(int argc, char **argv)
641 while ((ret = getopt(argc, argv, "v")) >= 0)
650 if (argc > 1 && strcmp(argv[1], "-v") == 0)