1 /* -*- Mode: C; c-basic-offset:4 ; indent-tabs-mode:nil ; -*- */
3 * (C) 2010 by Argonne National Laboratory.
4 * See COPYRIGHT in top-level directory.
12 * Tests for lock contention, including special cases within the MPICH code
13 * (any MPI implementation should pass these tests; in the MPICH case, our
14 * coverage analysis showed that the lockcontention.c test was not covering
15 * all cases, and in fact, this test revealed a bug in the code).
17 * In all of these tests, each process writes (or accesses) the values
18 * rank + i*size_of_world for NELM times.
20 * This test strives to avoid operations not strictly permitted by MPI RMA,
21 * for example, it doesn't target the same locations with multiple put/get
22 * calls in the same access epoch.
27 #define MAX_ERRS_REPORT 10
30 * Each process writes data into the rmabuf on the process with target rank
31 * trank. The final result in rmabuf are the consecutive integers starting
32 * from 0. Each process, however, does not write a consecutive block.
33 * Instead, they write these locations:
36 * for j=0,...,NBLOCK-1
37 * j + NBLOCK * (rank + i * wsize)
39 * The value written is the location.
41 * In many cases, multiple RMA operations are needed. Where these must not
42 * overlap, the above pattern is replicated at NBLOCK*NELM*wsize.
43 * (NBLOCK is either 1 or NBLOCK in the code below, depending on use)
46 static int toterrs = 0;
48 int testValues(int, int, int, int *, const char *);
50 int main(int argc, char *argv[])
52 int rank, wsize, i, j, cnt;
53 int *rmabuf, *localbuf, *localbuf2, *vals;
58 MTest_Init(&argc, &argv);
59 MPI_Comm_size(MPI_COMM_WORLD, &wsize);
60 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
63 fprintf(stderr, "Run this program with at least 2 processes\n");
64 MPI_Abort(MPI_COMM_WORLD, 1);
67 windowsize = (2 * NBLOCK + 2) * NELM * wsize;
68 rmabuf = (int *) malloc(windowsize * sizeof(int));
69 localbuf = (int *) malloc(NELM * sizeof(int));
70 localbuf2 = (int *) malloc(NELM * NBLOCK * sizeof(int));
71 vals = (int *) malloc(NELM * sizeof(int));
74 * Initialize the buffers
76 for (i = 0; i < NELM; i++) {
77 localbuf[i] = rank + i * wsize;
80 for (i = 0; i < NELM; i++) {
81 for (j = 0; j < NBLOCK; j++) {
82 localbuf2[cnt++] = j + NBLOCK * (rank + i * wsize);
85 for (i = 0; i < windowsize; i++) {
89 /* Create the window */
90 MPI_Win_create(rmabuf, windowsize * sizeof(int), sizeof(int), MPI_INFO_NULL,
91 MPI_COMM_WORLD, &win);
93 /* Multiple puts, with contention at trank */
94 MPI_Barrier(MPI_COMM_WORLD);
95 for (i = 0; i < NELM; i++) {
96 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
97 MPI_Put(&localbuf[i], 1, MPI_INT, trank, rank + i * wsize, 1, MPI_INT, win);
98 MPI_Put(&localbuf[i], 1, MPI_INT, trank, rank + (i + NELM) * wsize, 1, MPI_INT, win);
99 MPI_Win_unlock(trank, win);
101 MPI_Barrier(MPI_COMM_WORLD);
103 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
104 toterrs += testValues(1, NELM, wsize, rmabuf, "Multiple puts (1)");
105 toterrs += testValues(1, NELM, wsize, rmabuf + wsize * NELM, "Multiple puts (2)");
106 MPI_Win_unlock(trank, win);
108 MPI_Barrier(MPI_COMM_WORLD);
109 /* Reinit the rmabuf */
110 for (i = 0; i < windowsize; i++) {
113 MPI_Barrier(MPI_COMM_WORLD);
115 /* Single put with contention */
117 for (i = 0; i < NELM; i++) {
118 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
119 MPI_Put(&localbuf[i], 1, MPI_INT, trank, rank + i * wsize, 1, MPI_INT, win);
120 MPI_Win_unlock(trank, win);
122 MPI_Barrier(MPI_COMM_WORLD);
124 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
125 toterrs += testValues(1, NELM, wsize, rmabuf, "Single put");
126 MPI_Win_unlock(trank, win);
129 /* Reinit the rmabuf */
130 for (i = 0; i < windowsize; i++) {
133 /* Longer puts with contention at trank */
134 MPI_Barrier(MPI_COMM_WORLD);
135 for (i = 0; i < NELM; i++) {
136 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
138 MPI_Put(&localbuf2[i * NBLOCK], NBLOCK, MPI_INT, trank,
139 NELM * wsize + NBLOCK * (rank + i * wsize), NBLOCK, MPI_INT, win);
140 MPI_Put(&localbuf2[i * NBLOCK], NBLOCK, MPI_INT, trank,
141 NELM * wsize + NBLOCK * (rank + (i + NELM) * wsize), NBLOCK, MPI_INT, win);
143 MPI_Put(&localbuf[i], 1, MPI_INT, trank, rank + i * wsize, 1, MPI_INT, win);
144 MPI_Win_unlock(trank, win);
146 MPI_Barrier(MPI_COMM_WORLD);
148 /* For simplicity in testing, set the values that rank==trank
150 for (i = 0; i < NELM; i++) {
151 for (j = 0; j < NBLOCK; j++) {
152 rmabuf[NELM * wsize + NBLOCK * (trank + i * wsize) + j] =
153 j + NBLOCK * (trank + i * wsize);
154 rmabuf[NELM * wsize + NBLOCK * (trank + (i + NELM) * wsize) + j] =
155 j + NBLOCK * (trank + i * wsize);
158 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
159 toterrs += testValues(1, NELM, wsize, rmabuf, "Long puts (1)");
160 toterrs += testValues(NBLOCK, NELM, wsize, rmabuf + NELM * wsize, "Long puts(2)");
161 toterrs += testValues(NBLOCK, NELM, wsize,
162 rmabuf + NELM * wsize * (1 + NBLOCK), "Long puts(3)");
163 MPI_Win_unlock(trank, win);
166 /* Reinit the rmabuf */
167 for (i = 0; i < windowsize; i++) {
170 for (i = 0; i < NELM; i++)
173 /* Put mixed with Get */
174 MPI_Barrier(MPI_COMM_WORLD);
175 for (i = 0; i < NELM; i++) {
176 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
178 MPI_Put(&localbuf2[i], NBLOCK, MPI_INT, trank,
179 NELM * wsize + NBLOCK * (rank + i * wsize), NBLOCK, MPI_INT, win);
180 MPI_Put(&localbuf[i], 1, MPI_INT, trank, rank + i * wsize, 1, MPI_INT, win);
183 MPI_Get(&vals[i], 1, MPI_INT, trank, i, 1, MPI_INT, win);
185 MPI_Win_unlock(trank, win);
187 MPI_Barrier(MPI_COMM_WORLD);
189 /* Just test the Get */
190 for (i = 0; i < wsize; i++) {
194 if (toterrs < MAX_ERRS_REPORT) {
195 printf("put/get: vals[%d] = %d, expected -1\n", i, vals[i]);
199 else if (vals[i] != i && vals[i] != -1) {
201 if (toterrs < MAX_ERRS_REPORT) {
202 printf("put/get: vals[%d] = %d, expected -1 or %d\n", i, vals[i], i);
208 /* Contention only with get */
209 for (i = 0; i < windowsize; i++) {
212 for (i = 0; i < NELM; i++)
215 MPI_Barrier(MPI_COMM_WORLD);
216 for (i = 0; i < NELM; i++) {
217 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
218 MPI_Get(&vals[i], 1, MPI_INT, trank, i, 1, MPI_INT, win);
219 MPI_Win_unlock(trank, win);
221 MPI_Barrier(MPI_COMM_WORLD);
223 for (i = 0; i < NELM; i++) {
226 if (toterrs < MAX_ERRS_REPORT) {
227 printf("single get: vals[%d] = %d, expected %d\n", i, vals[i], -i);
233 /* Contention with accumulate */
234 MPI_Barrier(MPI_COMM_WORLD);
235 for (i = 0; i < NELM * wsize; i++) {
238 MPI_Barrier(MPI_COMM_WORLD);
239 for (i = 0; i < NELM; i++) {
240 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
241 MPI_Accumulate(&localbuf[i], 1, MPI_INT, trank, rank + i * wsize, 1, MPI_INT, MPI_SUM, win);
242 MPI_Accumulate(&localbuf[i], 1, MPI_INT, trank, rank + i * wsize, 1, MPI_INT, MPI_SUM, win);
243 MPI_Win_unlock(trank, win);
245 MPI_Barrier(MPI_COMM_WORLD);
247 MPI_Win_lock(MPI_LOCK_EXCLUSIVE, trank, 0, win);
248 for (i = 0; i < NELM * wsize; i++) {
249 if (rmabuf[i] != 2 * i) {
251 if (toterrs < MAX_ERRS_REPORT) {
252 printf("2 accumulate: rmabuf[%d] = %d, expected %d\n", i, rmabuf[i], 2 * i);
256 MPI_Win_unlock(trank, win);
266 MTest_Finalize(toterrs);
271 /* Test the values in the rmabuf against the expected values. Return the
273 int testValues(int nb, int nelm, int wsize, int *rmabuf, const char *msg)
277 for (i = 0; i < nb * nelm * wsize; i++) {
278 if (rmabuf[i] != i) {
279 if (toterrs + errs < MAX_ERRS_REPORT) {
280 printf("%s:rmabuf[%d] = %d expected %d\n", msg, i, rmabuf[i], i);