1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
3 * (C) 2012 by Argonne National Laboratory.
4 * See COPYRIGHT in top-level directory.
19 #define PUT_VAL 0xdcba97
20 #define ACC_VAL 10771134
23 * Additional tests for lock contention. These are designed to exercise
24 * some of the optimizations within MPICH, but all are valid MPI programs.
25 * Tests structure includes
26 * lock local (must happen at this time since application can use load
27 * store after the lock)
28 * send message to partner
32 * Provide a delay so that
33 * the partner will see the
36 * lock // Note: this may block
37 * rma operations (see below)
40 * unlock send back to partner
41 * receive from partner
42 * check for correct data
44 * The delay may be implemented as a ring of message communication; this
45 * is likely to automatically scale the time to what is needed
48 /* Define a datatype to be used with */
52 /* Define long RMA ops size */
57 void RMATest( int i, MPI_Win win, int master, int *srcbuf, int srcbufsize, int *getbuf, int getbufsize );
58 int RMACheck( int i, int *buf, MPI_Aint bufsize );
59 int RMACheckGet( int i, MPI_Win win, int *getbuf, MPI_Aint getsize);
60 void RMATestInit( int i, int *buf, MPI_Aint bufsize );
62 int main( int argc, char *argv[] )
66 int *rmabuffer=0, *getbuf=0;
67 MPI_Aint bufsize=0, getbufsize=0;
68 int master, partner, next, wrank, wsize, i;
69 int ntest = LAST_TEST;
72 MTest_Init( &argc, &argv );
74 /* Determine who is responsible for each part of the test */
75 MPI_Comm_rank( MPI_COMM_WORLD, &wrank );
76 MPI_Comm_size( MPI_COMM_WORLD, &wsize );
78 fprintf( stderr, "This test requires at least 3 processes\n" );
79 MPI_Abort( MPI_COMM_WORLD, 1 );
85 if (next == partner) next++;
88 if (next == partner) next++;
91 /* Determine the last test to run (by default, run them all) */
92 for (i=1; i<argc; i++) {
93 if (strcmp( "-ntest", argv[i] ) == 0) {
96 ntest = atoi( argv[i] );
99 fprintf( stderr, "Missing value for -ntest\n" );
100 MPI_Abort( MPI_COMM_WORLD, 1 );
105 MPI_Type_vector( veccount, 1, stride, MPI_INT, &vectype );
106 MPI_Type_commit( &vectype );
108 /* Create the RMA window */
110 if (wrank == master) {
112 MPI_Alloc_mem( bufsize*sizeof(int), MPI_INFO_NULL, &rmabuffer );
114 else if (wrank == partner) {
115 getbufsize = RMA_SIZE;
116 getbuf = (int *)malloc( getbufsize*sizeof(int) );
118 fprintf( stderr, "Unable to allocated %d bytes for getbuf\n",
120 MPI_Abort( MPI_COMM_WORLD, 1 );
123 srcbuf = malloc(RMA_SIZE*sizeof(*srcbuf));
126 MPI_Win_create( rmabuffer, bufsize, sizeof(int), MPI_INFO_NULL,
127 MPI_COMM_WORLD, &win );
129 /* Run a sequence of tests */
130 for (i=0; i<=ntest; i++) {
131 if (wrank == master) {
132 MTestPrintfMsg( 0, "Test %d\n", i );
133 /* Because this lock is local, it must return only when the
135 MPI_Win_lock( MPI_LOCK_EXCLUSIVE, 0, master, win );
136 RMATestInit( i, rmabuffer, bufsize );
137 MPI_Send( MPI_BOTTOM, 0, MPI_INT, partner, i, MPI_COMM_WORLD );
138 MPI_Send( MPI_BOTTOM, 0, MPI_INT, next, i, MPI_COMM_WORLD );
139 MPI_Recv( MPI_BOTTOM, 0, MPI_INT, MPI_ANY_SOURCE, i,
140 MPI_COMM_WORLD, MPI_STATUS_IGNORE );
141 MPI_Win_unlock( master, win );
142 MPI_Recv( MPI_BOTTOM, 0, MPI_INT, partner, i, MPI_COMM_WORLD,
144 errs += RMACheck( i, rmabuffer, bufsize );
146 else if (wrank == partner) {
147 MPI_Recv( MPI_BOTTOM, 0, MPI_INT, master, i, MPI_COMM_WORLD,
149 MPI_Win_lock( MPI_LOCK_EXCLUSIVE, 0, master, win );
150 RMATest( i, win, master, srcbuf, RMA_SIZE, getbuf, getbufsize );
151 MPI_Win_unlock( master, win );
152 errs += RMACheckGet( i, win, getbuf, getbufsize );
153 MPI_Send( MPI_BOTTOM, 0, MPI_INT, master, i, MPI_COMM_WORLD );
156 MPI_Recv( MPI_BOTTOM, 0, MPI_INT, MPI_ANY_SOURCE, i,
157 MPI_COMM_WORLD, MPI_STATUS_IGNORE );
158 MPI_Send( MPI_BOTTOM, 0, MPI_INT, next, i, MPI_COMM_WORLD );
163 MPI_Free_mem( rmabuffer );
168 MPI_Win_free( &win );
169 MPI_Type_free( &vectype );
171 MTest_Finalize( errs );
173 return MTestReturnValue( errs );
176 /* Perform the tests.
178 * The srcbuf must be passed in because the buffer must remain valid
179 * until the subsequent unlock call. */
180 void RMATest( int i, MPI_Win win, int master, int *srcbuf, int srcbufsize, int *getbuf, int getbufsize )
183 int *source = srcbuf;
184 assert(srcbufsize == RMA_SIZE);
186 for (j=0; j<srcbufsize; j++) source[j] = -j;
189 case 0: /* Single short put (1 word at OFFSET_1) */
191 MPI_Put( source, 1, MPI_INT, master, OFFSET_1, 1, MPI_INT, win );
193 case 1: /* Single short accumulate (1 word of value 17 at OFFSET_2) */
195 MPI_Accumulate( source, 1, MPI_INT, master,
196 OFFSET_2, 1, MPI_INT, MPI_SUM, win );
198 case 2: /* Single short get (1 word at OFFSET_3) */
200 MPI_Get( getbuf, 1, MPI_INT, master, OFFSET_3, 1, MPI_INT, win );
202 case 3: /* Datatype single put (strided put) */
203 for (j=0; j<veccount; j++) {
204 source[j*stride] = PUT_VAL + j;
206 MPI_Put( source, 1, vectype, master, OFFSET_1, 1, vectype, win );
208 case 4: /* Datatype single accumulate (strided acc) */
209 for (j=0; j<veccount; j++) {
210 source[j*stride] = ACC_VAL + j;
212 MPI_Accumulate( source, 1, vectype, master,
213 OFFSET_2, 1, vectype, MPI_SUM, win );
215 case 5: /* Datatype single get (strided get) */
216 for (j=0; j<veccount; j++) {
219 MPI_Get( getbuf, 1, vectype, master,
220 OFFSET_3, 1, vectype, win );
222 case 6: /* a few small puts (like strided put, but 1 word at a time) */
223 for (j=0; j<veccount; j++) {
224 source[j*stride] = PUT_VAL + j;
226 for (j=0; j<veccount; j++) {
227 MPI_Put( source + j*stride, 1, MPI_INT, master,
228 OFFSET_1+j*stride, 1, MPI_INT, win );
231 case 7: /* a few small accumulates (like strided acc, but 1 word at a time )*/
232 for (j=0; j<veccount; j++) {
233 source[j*stride] = ACC_VAL + j;
235 for (j=0; j<veccount; j++) {
236 MPI_Accumulate( source + j*stride, 1, MPI_INT, master,
237 OFFSET_2+j*stride, 1, MPI_INT, MPI_SUM, win );
240 case 8: /* a few small gets (like strided get, but 1 word at a time) */
241 for (j=0; j<veccount; j++) {
242 getbuf[j*stride] = -j;
244 for (j=0; j<veccount; j++) {
245 MPI_Get( getbuf + j*stride, 1, MPI_INT, master,
246 OFFSET_3+j*stride, 1, MPI_INT, win );
249 case 9: /* Single long put (OFFSET_1) */
250 for (j=0; j<longcount; j++) source[j] = j;
251 MPI_Put( source, longcount, MPI_INT, master,
252 OFFSET_1, longcount, MPI_INT, win );
254 case 10: /* Single long accumulate (OFFSET_2) */
255 for (j=0; j<longcount; j++) source[j] = j;
256 MPI_Accumulate( source, longcount, MPI_INT, master,
257 OFFSET_2, longcount, MPI_INT, MPI_SUM, win );
259 case 11: /* Single long get (OFFSET_3) */
260 for (j=0; j<longcount; j++) getbuf[j] = -j;
261 MPI_Get( getbuf, longcount, MPI_INT, master,
262 OFFSET_3, longcount, MPI_INT, win );
264 case 12: /* a few long puts (start at OFFSET_1, medcount ) */
265 for (j=0; j<mednum; j++) {
266 for (k=0; k<medcount; k++) {
267 source[j*medcount+k] = j*2*medcount+k;
269 MPI_Put( source + j*medcount, medcount, MPI_INT, master,
270 OFFSET_1 + j*2*medcount, medcount, MPI_INT, win );
273 case 13: /* a few long accumulates (start at OFFSET_2, medcount) */
274 for (j=0; j<mednum; j++) {
275 for (k=0; k<medcount; k++) {
276 source[j*medcount+k] = ACC_VAL + j*2*medcount+k;
278 MPI_Accumulate( source + j*medcount, medcount, MPI_INT, master,
279 OFFSET_2 + j*2*medcount, medcount, MPI_INT,
283 case 14: /* a few long gets (start at OFFSET_3, medcount) */
284 for (j=0; j<mednum; j++) {
285 for (k=0; k<medcount; k++) {
286 getbuf[j*medcount+k] = -(j*medcount+k);
288 MPI_Get( getbuf + j*medcount, medcount, MPI_INT, master,
289 OFFSET_3 + j*2*medcount, medcount, MPI_INT, win );
295 int RMACheck( int i, int *buf, MPI_Aint bufsize )
301 case 0: /* Single short put (1 word at OFFSET_1) */
302 if (buf[OFFSET_1] != PUT_VAL) {
304 printf( "case 0: value is %d should be %d\n",
305 buf[OFFSET_1], PUT_VAL );
308 case 1: /* Single short accumulate (1 word of value 17 at OFFSET_2) */
309 if (buf[OFFSET_2] != ACC_VAL + OFFSET_2) {
311 printf( "case 1: value is %d should be %d\n",
312 buf[OFFSET_2], ACC_VAL + OFFSET_2 );
315 case 2: /* Single short get (1 word at OFFSET_3) */
316 /* See RMACheckGet */
318 case 3: /* Datatype single put (strided put) */
319 case 6: /* a few small puts (like strided put, but 1 word at a time) */
320 /* FIXME: The conditional and increment are reversed below. This looks
321 * like a bug, and currently prevents the following test from running. */
322 for (j=0; j++; j<veccount) {
323 if (buf[j*stride] != PUT_VAL + j) {
325 printf( "case %d: value is %d should be %d\n", i,
326 buf[j*stride], PUT_VAL+j );
330 case 4: /* Datatype single accumulate (strided acc) */
331 case 7: /* a few small accumulates (like strided acc, but 1 word at a time )*/
332 /* FIXME: The conditional and increment are reversed below. This looks
333 * like a bug, and currently prevents the following test from running. */
334 for (j=0; j++; j<veccount) {
335 if (buf[j*stride] != ACC_VAL + j + OFFSET_2 + j*stride) {
337 printf( "case %d: value is %d should be %d\n", i,
338 buf[j*stride], ACC_VAL+j+OFFSET_2+j*stride );
342 case 5: /* Datatype single get (strided get) */
343 case 8: /* a few small gets (like strided get, but 1 word at a time) */
344 /* See RMACheckGet */
346 case 9: /* Single long put (OFFSET_1) */
347 for (j=0; j<longcount; j++) {
348 if (buf[OFFSET_1+j] != j) {
350 printf( "case 9: value is %d should be %d\n",
351 buf[OFFSET_1+j], OFFSET_1 + j );
355 case 10: /* Single long accumulate (OFFSET_2) */
356 for (j=0; j<longcount; j++) {
357 if (buf[OFFSET_2+j] != OFFSET_2 + j + j) {
359 printf( "case 10: value is %d should be %d\n",
360 buf[OFFSET_2+j], OFFSET_2 + j + j );
364 case 11: /* Single long get (OFFSET_3) */
365 /* See RMACheckGet */
367 case 12: /* a few long puts (start at OFFSET_1, medcount ) */
368 for (j=0; j<mednum; j++) {
369 for (k=0; k<medcount; k++) {
370 if (buf[OFFSET_1 + j*2*medcount + k] !=
373 printf( "case 12: value is %d should be %d\n",
374 buf[OFFSET_1+j*2*medcount + k], j*2*medcount + k );
379 case 13: /* a few long accumulates (start at OFFSET_2, medcount) */
380 for (j=0; j<mednum; j++) {
381 for (k=0; k<medcount; k++) {
382 if (buf[OFFSET_2 + j*2*medcount + k] !=
383 OFFSET_2 + 2*j*2*medcount+2*k + ACC_VAL ) {
385 printf( "case 13: value is %d should be %d\n",
386 buf[OFFSET_2+j*2*medcount + k],
387 OFFSET_2 + 2*j*2*medcount + k +ACC_VAL);
392 case 14: /* a few long gets (start at OFFSET_3, medcount) */
393 /* See RMACheckGet */
396 fprintf( stderr, "Unrecognized case %d\n", i );
403 int RMACheckGet( int i, MPI_Win win, int *getbuf, MPI_Aint getsize)
410 case 0: /* Single short put (1 word at OFFSET_1) */
412 case 1: /* Single short accumulate (1 word of value 17 at OFFSET_2) */
414 case 2: /* Single short get (1 word at OFFSET_3) */
415 if (getbuf[0] != OFFSET_3) {
417 printf( "case 2: value is %d should be %d\n",
418 getbuf[0], OFFSET_3 );
421 case 3: /* Datatype single put (strided put) */
423 case 4: /* Datatype single accumulate (strided acc) */
425 case 5: /* Datatype single get (strided get) */
426 case 8: /* a few small gets (like strided get, but 1 word at a time) */
427 for (j=0; j<veccount; j++) {
428 if (getbuf[j*stride] != OFFSET_3 + j*stride) {
430 printf( "case %d: value is %d should be %d\n", i,
431 getbuf[j*stride], OFFSET_3 + j*stride );
436 case 6: /* a few small puts (like strided put, but 1 word at a time) */
438 case 7: /* a few small accumulates (like strided acc, but 1 word at a time )*/
440 case 9: /* Single long put (OFFSET_1) */
442 case 10: /* Single long accumulate (OFFSET_2) */
444 case 11: /* Single long get (OFFSET_3) */
445 for (j=0; j<longcount; j++) {
446 if (getbuf[j] != OFFSET_3 + j) {
448 printf( "case 11: value is %d should be %d\n",
449 getbuf[j], OFFSET_3 + j );
453 case 12: /* a few long puts (start at OFFSET_1, medcount ) */
455 case 13: /* a few long accumulates (start at OFFSET_2, medcount) */
457 case 14: /* a few long gets (start at OFFSET_3, medcount) */
458 for (j=0; j<mednum; j++) {
459 for (k=0; k<medcount; k++) {
460 if (getbuf[j*medcount + k] !=
461 OFFSET_3 + j*2*medcount+k ) {
463 printf( "case 14: buf[%d] value is %d should be %d\n",
465 getbuf[j*medcount + k],
466 OFFSET_3 + j*2*medcount + k );
472 fprintf( stderr, "Unrecognized case %d\n", i );
480 void RMATestInit( int i, int *buf, MPI_Aint bufsize )
483 for (j=0; j<bufsize; j++) {