Actual source code: ex62.c
1: static char help[] = "Tests `PetscGarbageKeySortedIntersect()`\n\n";
3: #include <petscsys.h>
4: #include <petsc/private/garbagecollector.h>
6: /* This program tests `PetscGarbageKeySortedIntersect(), which is the
7: public (MPI) interface to
8: `PetscErrorCode GarbageKeySortedIntersect_Private()`.
9: Sets are sent packed in arrays, with the first entry as the number of
10: set elements and the sets the remaining elements. This is because the
11: MPI reduction operation must have the call signature:
12: void PetscGarbageKeySortedIntersect(void *inset, void *inoutset, PetscMPIInt *length, MPI_Datatype *dtype)
13: This is a thin wrapper for the private routine:
14: PetscErrorCode GarbageKeySortedIntersect_Private(PetscInt64 seta[], PetscInt *lena, PetscInt64 setb[], PetscInt lenb)
15: Where
16: seta = (PetscInt64 *)inoutset;
17: setb = (PetscInt64 *)inset;
18: And the arguments are passed as:
19: &seta[1], (PetscInt *)&seta[0], &setb[1], (PetscInt)setb[0]
20: */
22: /* Populate a set with upto the first 49 unique Fibonnaci numbers */
23: PetscErrorCode Fibonnaci(PetscInt64 **set, PetscInt n)
24: {
25: PetscInt ii;
26: PetscInt64 fib[] = {1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584,
27: 4181, 6765, 10946, 17711, 28657, 46368, 75025, 121393, 196418, 317811, 514229, 832040, 1346269, 2178309, 3524578, 5702887, 9227465,
28: 14930352, 24157817, 39088169, 63245986, 102334155, 165580141, 267914296, 433494437, 701408733, 1134903170, 1836311903, 2971215073, 4807526976, 7778742049, 12586269025};
30: PetscFunctionBeginUser;
31: PetscAssert(n < 50, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "n must be less than 50");
32: PetscCall(PetscMalloc1(n + 1, set));
33: (*set)[0] = n;
34: for (ii = 0; ii < n; ii++) { (*set)[ii + 1] = fib[ii]; }
35: PetscFunctionReturn(PETSC_SUCCESS);
36: }
38: /* Populate a set with Square numbers */
39: PetscErrorCode Square(PetscInt64 **set, PetscInt n)
40: {
41: PetscInt64 ii;
43: PetscFunctionBeginUser;
44: PetscCall(PetscMalloc1(n + 1, set));
45: (*set)[0] = n;
46: for (ii = 1; ii < n + 1; ii++) { (*set)[ii] = ii * ii; }
47: PetscFunctionReturn(PETSC_SUCCESS);
48: }
50: /* Populate a set with Cube numbers */
51: PetscErrorCode Cube(PetscInt64 **set, PetscInt n)
52: {
53: PetscInt64 ii;
55: PetscFunctionBeginUser;
56: PetscCall(PetscMalloc1(n + 1, set));
57: (*set)[0] = n;
58: for (ii = 1; ii < n + 1; ii++) { (*set)[ii] = ii * ii * ii; }
59: PetscFunctionReturn(PETSC_SUCCESS);
60: }
62: /* Populate a set with numbers to sixth power */
63: PetscErrorCode Sixth(PetscInt64 **set, PetscInt n)
64: {
65: PetscInt64 ii;
67: PetscFunctionBeginUser;
68: PetscCall(PetscMalloc1(n + 1, set));
69: (*set)[0] = n;
70: for (ii = 1; ii < n + 1; ii++) { (*set)[ii] = ii * ii * ii * ii * ii * ii; }
71: PetscFunctionReturn(PETSC_SUCCESS);
72: }
74: /* Print out the contents of a set */
75: PetscErrorCode PrintSet(PetscInt64 *set)
76: {
77: char text[64];
79: PetscFunctionBeginUser;
80: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "["));
81: for (PetscInt64 ii = 1; ii <= set[0]; ii++) {
82: PetscCall(PetscFormatConvert(" %" PetscInt64_FMT ",", text));
83: PetscCall(PetscPrintf(PETSC_COMM_WORLD, text, set[ii]));
84: }
85: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));
86: PetscFunctionReturn(PETSC_SUCCESS);
87: }
89: /* Check set equality */
90: PetscErrorCode AssertSetsEqual(PetscInt64 *set, PetscInt64 *true_set)
91: {
92: PetscInt ii;
94: PetscFunctionBeginUser;
95: PetscAssert(set[0] == true_set[0], PETSC_COMM_WORLD, PETSC_ERR_ARG_INCOMP, "Sets of different sizes");
96: for (ii = 1; ii < set[0] + 1; ii++) PetscAssert(set[ii] == true_set[ii], PETSC_COMM_WORLD, PETSC_ERR_ARG_INCOMP, "Sets are different");
97: PetscFunctionReturn(PETSC_SUCCESS);
98: }
100: /* Tests functionality when two enpty sets are passed */
101: PetscErrorCode test_empty_empty()
102: {
103: PetscInt64 *set_a, *set_b;
104: PetscInt64 truth[] = {0};
105: PetscMPIInt length = 1;
107: PetscFunctionBeginUser;
108: PetscCall(PetscMalloc1(1, &set_a));
109: PetscCall(PetscMalloc1(1, &set_b));
111: set_a[0] = 0;
113: set_b[0] = 0;
115: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
116: PetscCall(PrintSet(set_a));
117: PetscCall(AssertSetsEqual(set_a, truth));
119: PetscCall(PetscFree(set_a));
120: PetscCall(PetscFree(set_b));
121: PetscFunctionReturn(PETSC_SUCCESS);
122: }
124: /* Tests functionality when seta is empty */
125: PetscErrorCode test_a_empty()
126: {
127: PetscInt64 *set_a, *set_b;
128: PetscInt64 truth[] = {0};
129: PetscMPIInt length = 1;
131: PetscFunctionBeginUser;
132: PetscCall(PetscMalloc1(1, &set_a));
133: PetscCall(PetscMalloc1(2, &set_b));
135: set_a[0] = 0;
137: set_b[0] = 1;
138: set_b[1] = 1;
140: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
141: PetscCall(PrintSet(set_a));
142: PetscCall(AssertSetsEqual(set_a, truth));
144: PetscCall(PetscFree(set_a));
145: PetscCall(PetscFree(set_b));
146: PetscFunctionReturn(PETSC_SUCCESS);
147: }
149: /* Tests functionality when setb is empty */
150: PetscErrorCode test_b_empty()
151: {
152: PetscInt64 *set_a, *set_b;
153: PetscInt64 truth[] = {0};
154: PetscMPIInt length = 1;
156: PetscFunctionBeginUser;
157: PetscCall(PetscMalloc1(2, &set_a));
158: PetscCall(PetscMalloc1(1, &set_b));
160: set_a[0] = 1;
161: set_a[1] = 1;
163: set_b[0] = 0;
165: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
166: PetscCall(PrintSet(set_a));
167: PetscCall(AssertSetsEqual(set_a, truth));
169: PetscCall(PetscFree(set_a));
170: PetscCall(PetscFree(set_b));
171: PetscFunctionReturn(PETSC_SUCCESS);
172: }
174: /* Tests functionality when both sets are identical */
175: PetscErrorCode test_identical()
176: {
177: PetscInt64 *set_a, *set_b;
178: PetscInt64 truth[] = {3, 1, 4, 9};
179: PetscMPIInt length = 4;
181: PetscFunctionBeginUser;
182: PetscCall(PetscMalloc1(4, &set_a));
183: PetscCall(PetscMalloc1(4, &set_b));
185: set_a[0] = 3;
186: set_a[1] = 1;
187: set_a[2] = 4;
188: set_a[3] = 9;
190: set_b[0] = 3;
191: set_b[1] = 1;
192: set_b[2] = 4;
193: set_b[3] = 9;
195: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
196: PetscCall(PrintSet(set_a));
197: PetscCall(AssertSetsEqual(set_a, truth));
199: PetscCall(PetscFree(set_a));
200: PetscCall(PetscFree(set_b));
201: PetscFunctionReturn(PETSC_SUCCESS);
202: }
204: /* Tests functionality when sets have no elements in common */
205: PetscErrorCode test_disjoint()
206: {
207: PetscInt64 *set_a, *set_b;
208: PetscInt64 truth[] = {0};
209: PetscMPIInt length = 1;
211: PetscFunctionBeginUser;
212: PetscCall(PetscMalloc1(4, &set_a));
213: PetscCall(PetscMalloc1(4, &set_b));
215: set_a[0] = 3;
216: set_a[1] = 1;
217: set_a[2] = 4;
218: set_a[3] = 9;
220: set_b[0] = 3;
221: set_b[1] = 2;
222: set_b[2] = 6;
223: set_b[3] = 8;
225: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
226: PetscCall(PrintSet(set_a));
227: PetscCall(AssertSetsEqual(set_a, truth));
229: PetscCall(PetscFree(set_a));
230: PetscCall(PetscFree(set_b));
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: /* Tests functionality when sets only have one element in common */
235: PetscErrorCode test_single_common()
236: {
237: PetscInt64 *set_a, *set_b;
238: PetscInt64 truth[] = {1, 4};
239: PetscMPIInt length = 1;
241: PetscFunctionBeginUser;
242: PetscCall(PetscMalloc1(4, &set_a));
243: PetscCall(PetscMalloc1(5, &set_b));
245: set_a[0] = 3;
246: set_a[1] = 1;
247: set_a[2] = 4;
248: set_a[3] = 9;
250: set_b[0] = 3;
251: set_b[1] = 2;
252: set_b[2] = 4;
253: set_b[3] = 6;
254: set_b[4] = 8;
256: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
257: PetscCall(PrintSet(set_a));
258: PetscCall(AssertSetsEqual(set_a, truth));
260: PetscCall(PetscFree(set_a));
261: PetscCall(PetscFree(set_b));
262: PetscFunctionReturn(PETSC_SUCCESS);
263: }
265: /* Specific test case flagged by PETSc issue #1247 */
266: PetscErrorCode test_issue_1247()
267: {
268: PetscInt64 *set_a, *set_b;
269: PetscInt64 truth[] = {0};
270: PetscMPIInt length = 1;
272: PetscFunctionBeginUser;
273: PetscCall(PetscMalloc1(3, &set_a));
274: PetscCall(PetscMalloc1(2, &set_b));
276: set_a[0] = 2;
277: set_a[1] = 2;
278: set_a[2] = 3;
280: set_b[0] = 1;
281: set_b[1] = 1;
283: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
284: PetscCall(PrintSet(set_a));
285: PetscCall(AssertSetsEqual(set_a, truth));
287: PetscCall(PetscFree(set_a));
288: PetscCall(PetscFree(set_b));
289: PetscFunctionReturn(PETSC_SUCCESS);
290: }
292: /* Tests functionality when seta is empty and setb is large */
293: PetscErrorCode test_empty_big()
294: {
295: PetscInt64 *set_a, *set_b;
296: PetscInt64 truth[] = {0};
297: PetscMPIInt length = 1;
299: PetscFunctionBeginUser;
300: PetscCall(PetscMalloc1(1, &set_a));
301: PetscCall(Square(&set_b, 999));
303: set_a[0] = 0;
305: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
306: PetscCall(PrintSet(set_a));
307: PetscCall(AssertSetsEqual(set_a, truth));
309: PetscCall(PetscFree(set_a));
310: PetscCall(PetscFree(set_b));
311: PetscFunctionReturn(PETSC_SUCCESS);
312: }
314: /* Tests functionality when seta is small and setb is large */
315: PetscErrorCode test_small_big()
316: {
317: PetscInt64 *set_a, *set_b;
318: PetscInt64 truth[] = {3, 1, 4, 9};
319: PetscMPIInt length = 1;
321: PetscFunctionBeginUser;
322: PetscCall(PetscMalloc1(5, &set_a));
323: PetscCall(Square(&set_b, 999));
325: set_a[0] = 4;
326: set_a[1] = 1;
327: set_a[2] = 4;
328: set_a[3] = 8;
329: set_a[4] = 9;
331: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
332: PetscCall(PrintSet(set_a));
333: PetscCall(AssertSetsEqual(set_a, truth));
335: PetscCall(PetscFree(set_a));
336: PetscCall(PetscFree(set_b));
337: PetscFunctionReturn(PETSC_SUCCESS);
338: }
340: /* Tests functionality when seta is medium sized and setb is large */
341: PetscErrorCode test_moderate_big()
342: {
343: PetscInt64 *set_a, *set_b;
344: PetscInt64 truth[] = {2, 1, 144};
345: PetscMPIInt length = 1;
347: PetscFunctionBeginUser;
348: PetscCall(Fibonnaci(&set_a, 49));
349: PetscCall(Square(&set_b, 999));
351: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
352: PetscCall(PrintSet(set_a));
353: PetscCall(AssertSetsEqual(set_a, truth));
355: PetscCall(PetscFree(set_a));
356: PetscCall(PetscFree(set_b));
357: PetscFunctionReturn(PETSC_SUCCESS);
358: }
360: /* Tests functionality when seta and setb are large */
361: PetscErrorCode test_big_big()
362: {
363: PetscInt64 *set_a, *set_b;
364: PetscInt64 *truth;
365: PetscMPIInt length = 1;
367: PetscFunctionBeginUser;
368: PetscCall(Cube(&set_a, 999));
369: PetscCall(Square(&set_b, 999));
371: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
372: PetscCall(PrintSet(set_a));
374: PetscCall(Sixth(&truth, 9));
375: PetscCall(AssertSetsEqual(set_a, truth));
377: PetscCall(PetscFree(set_a));
378: PetscCall(PetscFree(set_b));
379: PetscCall(PetscFree(truth));
380: PetscFunctionReturn(PETSC_SUCCESS);
381: }
383: /* Tests functionality when setb is empty and setb is large */
384: PetscErrorCode test_big_empty()
385: {
386: PetscInt64 *set_a, *set_b;
387: PetscInt64 truth[] = {0};
388: PetscMPIInt length = 1;
390: PetscFunctionBeginUser;
391: PetscCall(Cube(&set_a, 999));
392: PetscCall(PetscMalloc1(1, &set_b));
394: set_b[0] = 0;
396: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
397: PetscCall(PrintSet(set_a));
398: PetscCall(AssertSetsEqual(set_a, truth));
400: PetscCall(PetscFree(set_a));
401: PetscCall(PetscFree(set_b));
402: PetscFunctionReturn(PETSC_SUCCESS);
403: }
405: /* Tests functionality when setb is small and setb is large */
406: PetscErrorCode test_big_small()
407: {
408: PetscInt64 *set_a, *set_b;
409: PetscInt64 truth[] = {2, 1, 8};
410: PetscMPIInt length = 1;
412: PetscFunctionBeginUser;
413: PetscCall(Cube(&set_a, 999));
414: PetscCall(PetscMalloc1(5, &set_b));
416: set_b[0] = 4;
417: set_b[1] = 1;
418: set_b[2] = 4;
419: set_b[3] = 8;
420: set_b[4] = 9;
422: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
423: PetscCall(PrintSet(set_a));
424: PetscCall(AssertSetsEqual(set_a, truth));
426: PetscCall(PetscFree(set_a));
427: PetscCall(PetscFree(set_b));
428: PetscFunctionReturn(PETSC_SUCCESS);
429: }
431: /* Tests functionality when setb is medium sized and setb is large */
432: PetscErrorCode test_big_moderate()
433: {
434: PetscInt64 *set_a, *set_b;
435: PetscInt64 truth[] = {2, 1, 8};
436: PetscMPIInt length = 1;
438: PetscFunctionBeginUser;
439: PetscCall(Cube(&set_a, 999));
440: PetscCall(Fibonnaci(&set_b, 49));
442: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
443: PetscCall(PrintSet(set_a));
444: PetscCall(AssertSetsEqual(set_a, truth));
446: PetscCall(PetscFree(set_a));
447: PetscCall(PetscFree(set_b));
448: PetscFunctionReturn(PETSC_SUCCESS);
449: }
451: /* Tests functionality when seta and setb are large, in the opposite
452: order to test_big_big() */
453: PetscErrorCode test_big_big_reversed()
454: {
455: PetscInt64 *set_a, *set_b;
456: PetscInt64 *truth;
457: PetscMPIInt length = 1;
459: PetscFunctionBeginUser;
460: PetscCall(Cube(&set_a, 999));
461: PetscCall(Square(&set_b, 999));
463: PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL);
464: PetscCall(PrintSet(set_a));
466: PetscCall(Sixth(&truth, 9));
467: PetscCall(AssertSetsEqual(set_a, truth));
469: PetscCall(PetscFree(set_a));
470: PetscCall(PetscFree(set_b));
471: PetscCall(PetscFree(truth));
472: PetscFunctionReturn(PETSC_SUCCESS);
473: }
475: /* Main executes the individual tests in a predefined order */
476: int main(int argc, char **argv)
477: {
478: PetscFunctionBeginUser;
479: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
481: /* Small tests */
482: /* Test different edge cases with small sets */
483: PetscCall(test_empty_empty());
484: PetscCall(test_a_empty());
485: PetscCall(test_b_empty());
486: PetscCall(test_identical());
487: PetscCall(test_disjoint());
488: PetscCall(test_single_common());
489: PetscCall(test_issue_1247());
491: /* Big tests */
492: /* Test different edge cases with big sets */
493: PetscCall(test_empty_big());
494: PetscCall(test_small_big());
495: PetscCall(test_moderate_big());
496: PetscCall(test_big_big());
497: PetscCall(test_big_empty());
498: PetscCall(test_big_small());
499: PetscCall(test_big_moderate());
500: PetscCall(test_big_big_reversed());
502: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "ALL PASSED\n"));
503: PetscCall(PetscFinalize());
504: return 0;
505: }
507: /*TEST
509: test:
510: suffix: 0
512: TEST*/