1 static char help[] = "Tests `PetscGarbageKeySortedIntersect()`\n\n"; 2 3 #include <petscsys.h> 4 #include <petsc/private/garbagecollector.h> 5 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 */ 21 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}; 29 30 PetscFunctionBeginUser; 31 PetscAssert((n < 50), PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "n must be less than 50\n"); 32 PetscCall(PetscMalloc1(n + 1, set)); 33 (*set)[0] = (PetscInt64)n; 34 for (ii = 0; ii < n; ii++) { (*set)[ii + 1] = fib[ii]; } 35 PetscFunctionReturn(PETSC_SUCCESS); 36 } 37 38 /* Populate a set with Square numbers */ 39 PetscErrorCode Square(PetscInt64 **set, PetscInt n) 40 { 41 PetscInt64 ii; 42 43 PetscFunctionBeginUser; 44 PetscCall(PetscMalloc1(n + 1, set)); 45 (*set)[0] = (PetscInt64)n; 46 for (ii = 1; ii < n + 1; ii++) { (*set)[ii] = ii * ii; } 47 PetscFunctionReturn(PETSC_SUCCESS); 48 } 49 50 /* Populate a set with Cube numbers */ 51 PetscErrorCode Cube(PetscInt64 **set, PetscInt n) 52 { 53 PetscInt64 ii; 54 55 PetscFunctionBeginUser; 56 PetscCall(PetscMalloc1(n + 1, set)); 57 (*set)[0] = (PetscInt64)n; 58 for (ii = 1; ii < n + 1; ii++) { (*set)[ii] = ii * ii * ii; } 59 PetscFunctionReturn(PETSC_SUCCESS); 60 } 61 62 /* Populate a set with numbers to sixth power */ 63 PetscErrorCode Sixth(PetscInt64 **set, PetscInt n) 64 { 65 PetscInt64 ii; 66 67 PetscFunctionBeginUser; 68 PetscCall(PetscMalloc1(n + 1, set)); 69 (*set)[0] = (PetscInt64)n; 70 for (ii = 1; ii < n + 1; ii++) { (*set)[ii] = ii * ii * ii * ii * ii * ii; } 71 PetscFunctionReturn(PETSC_SUCCESS); 72 } 73 74 /* Print out the contents of a set */ 75 PetscErrorCode PrintSet(PetscInt64 *set) 76 { 77 char text[64]; 78 PetscInt ii; 79 80 PetscFunctionBeginUser; 81 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "[")); 82 for (ii = 1; ii <= (PetscInt)set[0]; ii++) { 83 PetscCall(PetscFormatConvert(" %" PetscInt64_FMT ",", text)); 84 PetscCall(PetscPrintf(PETSC_COMM_WORLD, text, set[ii])); 85 } 86 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n")); 87 PetscFunctionReturn(PETSC_SUCCESS); 88 } 89 90 /* Check set equality */ 91 PetscErrorCode AssertSetsEqual(PetscInt64 *set, PetscInt64 *true_set) 92 { 93 PetscInt ii; 94 95 PetscFunctionBeginUser; 96 PetscAssert((set[0] == true_set[0]), PETSC_COMM_WORLD, PETSC_ERR_ARG_INCOMP, "Sets of different sizes"); 97 for (ii = 1; ii < set[0] + 1; ii++) PetscAssert((set[ii] == true_set[ii]), PETSC_COMM_WORLD, PETSC_ERR_ARG_INCOMP, "Sets are different"); 98 PetscFunctionReturn(PETSC_SUCCESS); 99 } 100 101 /* Tests functionality when two enpty sets are passed */ 102 PetscErrorCode test_empty_empty() 103 { 104 PetscInt64 *set_a, *set_b; 105 PetscInt64 truth[] = {0}; 106 PetscMPIInt length = 1; 107 108 PetscFunctionBeginUser; 109 110 PetscCall(PetscMalloc1(1, &set_a)); 111 PetscCall(PetscMalloc1(1, &set_b)); 112 113 set_a[0] = 0; 114 115 set_b[0] = 0; 116 117 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 118 PetscCall(PrintSet(set_a)); 119 PetscCall(AssertSetsEqual(set_a, truth)); 120 121 PetscCall(PetscFree(set_a)); 122 PetscCall(PetscFree(set_b)); 123 124 PetscFunctionReturn(PETSC_SUCCESS); 125 } 126 127 /* Tests functionality when seta is empty */ 128 PetscErrorCode test_a_empty() 129 { 130 PetscInt64 *set_a, *set_b; 131 PetscInt64 truth[] = {0}; 132 PetscMPIInt length = 1; 133 134 PetscFunctionBeginUser; 135 136 PetscCall(PetscMalloc1(1, &set_a)); 137 PetscCall(PetscMalloc1(2, &set_b)); 138 139 set_a[0] = 0; 140 141 set_b[0] = 1; 142 set_b[1] = 1; 143 144 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 145 PetscCall(PrintSet(set_a)); 146 PetscCall(AssertSetsEqual(set_a, truth)); 147 148 PetscCall(PetscFree(set_a)); 149 PetscCall(PetscFree(set_b)); 150 151 PetscFunctionReturn(PETSC_SUCCESS); 152 } 153 154 /* Tests functionality when setb is empty */ 155 PetscErrorCode test_b_empty() 156 { 157 PetscInt64 *set_a, *set_b; 158 PetscInt64 truth[] = {0}; 159 PetscMPIInt length = 1; 160 161 PetscFunctionBeginUser; 162 163 PetscCall(PetscMalloc1(2, &set_a)); 164 PetscCall(PetscMalloc1(1, &set_b)); 165 166 set_a[0] = 1; 167 set_a[1] = 1; 168 169 set_b[0] = 0; 170 171 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 172 PetscCall(PrintSet(set_a)); 173 PetscCall(AssertSetsEqual(set_a, truth)); 174 175 PetscCall(PetscFree(set_a)); 176 PetscCall(PetscFree(set_b)); 177 178 PetscFunctionReturn(PETSC_SUCCESS); 179 } 180 181 /* Tests functionality when both sets are identical */ 182 PetscErrorCode test_identical() 183 { 184 PetscInt64 *set_a, *set_b; 185 PetscInt64 truth[] = {3, 1, 4, 9}; 186 PetscMPIInt length = 4; 187 188 PetscFunctionBeginUser; 189 190 PetscCall(PetscMalloc1(4, &set_a)); 191 PetscCall(PetscMalloc1(4, &set_b)); 192 193 set_a[0] = 3; 194 set_a[1] = 1; 195 set_a[2] = 4; 196 set_a[3] = 9; 197 198 set_b[0] = 3; 199 set_b[1] = 1; 200 set_b[2] = 4; 201 set_b[3] = 9; 202 203 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 204 PetscCall(PrintSet(set_a)); 205 PetscCall(AssertSetsEqual(set_a, truth)); 206 207 PetscCall(PetscFree(set_a)); 208 PetscCall(PetscFree(set_b)); 209 210 PetscFunctionReturn(PETSC_SUCCESS); 211 } 212 213 /* Tests functionality when sets have no elements in common */ 214 PetscErrorCode test_disjoint() 215 { 216 PetscInt64 *set_a, *set_b; 217 PetscInt64 truth[] = {0}; 218 PetscMPIInt length = 1; 219 220 PetscFunctionBeginUser; 221 222 PetscCall(PetscMalloc1(4, &set_a)); 223 PetscCall(PetscMalloc1(4, &set_b)); 224 225 set_a[0] = 3; 226 set_a[1] = 1; 227 set_a[2] = 4; 228 set_a[3] = 9; 229 230 set_b[0] = 3; 231 set_b[1] = 2; 232 set_b[2] = 6; 233 set_b[3] = 8; 234 235 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 236 PetscCall(PrintSet(set_a)); 237 PetscCall(AssertSetsEqual(set_a, truth)); 238 239 PetscCall(PetscFree(set_a)); 240 PetscCall(PetscFree(set_b)); 241 242 PetscFunctionReturn(PETSC_SUCCESS); 243 } 244 245 /* Tests functionality when sets only have one element in common */ 246 PetscErrorCode test_single_common() 247 { 248 PetscInt64 *set_a, *set_b; 249 PetscInt64 truth[] = {1, 4}; 250 PetscMPIInt length = 1; 251 252 PetscFunctionBeginUser; 253 254 PetscCall(PetscMalloc1(4, &set_a)); 255 PetscCall(PetscMalloc1(5, &set_b)); 256 257 set_a[0] = 3; 258 set_a[1] = 1; 259 set_a[2] = 4; 260 set_a[3] = 9; 261 262 set_b[0] = 3; 263 set_b[1] = 2; 264 set_b[2] = 4; 265 set_b[3] = 6; 266 set_b[4] = 8; 267 268 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 269 PetscCall(PrintSet(set_a)); 270 PetscCall(AssertSetsEqual(set_a, truth)); 271 272 PetscCall(PetscFree(set_a)); 273 PetscCall(PetscFree(set_b)); 274 275 PetscFunctionReturn(PETSC_SUCCESS); 276 } 277 278 /* Specific test case flagged by PETSc issue #1247 */ 279 PetscErrorCode test_issue_1247() 280 { 281 PetscInt64 *set_a, *set_b; 282 PetscInt64 truth[] = {0}; 283 PetscMPIInt length = 1; 284 285 PetscFunctionBeginUser; 286 287 PetscCall(PetscMalloc1(3, &set_a)); 288 PetscCall(PetscMalloc1(2, &set_b)); 289 290 set_a[0] = 2; 291 set_a[1] = 2; 292 set_a[2] = 3; 293 294 set_b[0] = 1; 295 set_b[1] = 1; 296 297 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 298 PetscCall(PrintSet(set_a)); 299 PetscCall(AssertSetsEqual(set_a, truth)); 300 301 PetscCall(PetscFree(set_a)); 302 PetscCall(PetscFree(set_b)); 303 304 PetscFunctionReturn(PETSC_SUCCESS); 305 } 306 307 /* Tests functionality when seta is empty and setb is large */ 308 PetscErrorCode test_empty_big() 309 { 310 PetscInt64 *set_a, *set_b; 311 PetscInt64 truth[] = {0}; 312 PetscMPIInt length = 1; 313 314 PetscFunctionBeginUser; 315 316 PetscCall(PetscMalloc1(1, &set_a)); 317 PetscCall(Square(&set_b, 999)); 318 319 set_a[0] = 0; 320 321 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 322 PetscCall(PrintSet(set_a)); 323 PetscCall(AssertSetsEqual(set_a, truth)); 324 325 PetscCall(PetscFree(set_a)); 326 PetscCall(PetscFree(set_b)); 327 328 PetscFunctionReturn(PETSC_SUCCESS); 329 } 330 331 /* Tests functionality when seta is small and setb is large */ 332 PetscErrorCode test_small_big() 333 { 334 PetscInt64 *set_a, *set_b; 335 PetscInt64 truth[] = {3, 1, 4, 9}; 336 PetscMPIInt length = 1; 337 338 PetscFunctionBeginUser; 339 340 PetscCall(PetscMalloc1(5, &set_a)); 341 PetscCall(Square(&set_b, 999)); 342 343 set_a[0] = 4; 344 set_a[1] = 1; 345 set_a[2] = 4; 346 set_a[3] = 8; 347 set_a[4] = 9; 348 349 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 350 PetscCall(PrintSet(set_a)); 351 PetscCall(AssertSetsEqual(set_a, truth)); 352 353 PetscCall(PetscFree(set_a)); 354 PetscCall(PetscFree(set_b)); 355 356 PetscFunctionReturn(PETSC_SUCCESS); 357 } 358 359 /* Tests functionality when seta is medium sized and setb is large */ 360 PetscErrorCode test_moderate_big() 361 { 362 PetscInt64 *set_a, *set_b; 363 PetscInt64 truth[] = {2, 1, 144}; 364 PetscMPIInt length = 1; 365 366 PetscFunctionBeginUser; 367 368 PetscCall(Fibonnaci(&set_a, 49)); 369 PetscCall(Square(&set_b, 999)); 370 371 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 372 PetscCall(PrintSet(set_a)); 373 PetscCall(AssertSetsEqual(set_a, truth)); 374 375 PetscCall(PetscFree(set_a)); 376 PetscCall(PetscFree(set_b)); 377 378 PetscFunctionReturn(PETSC_SUCCESS); 379 } 380 381 /* Tests functionality when seta and setb are large */ 382 PetscErrorCode test_big_big() 383 { 384 PetscInt64 *set_a, *set_b; 385 PetscInt64 *truth; 386 PetscMPIInt length = 1; 387 388 PetscFunctionBeginUser; 389 390 PetscCall(Cube(&set_a, 999)); 391 PetscCall(Square(&set_b, 999)); 392 393 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 394 PetscCall(PrintSet(set_a)); 395 396 PetscCall(Sixth(&truth, 9)); 397 PetscCall(AssertSetsEqual(set_a, truth)); 398 399 PetscCall(PetscFree(set_a)); 400 PetscCall(PetscFree(set_b)); 401 PetscCall(PetscFree(truth)); 402 403 PetscFunctionReturn(PETSC_SUCCESS); 404 } 405 406 /* Tests functionality when setb is empty and setb is large */ 407 PetscErrorCode test_big_empty() 408 { 409 PetscInt64 *set_a, *set_b; 410 PetscInt64 truth[] = {0}; 411 PetscMPIInt length = 1; 412 413 PetscFunctionBeginUser; 414 415 PetscCall(Cube(&set_a, 999)); 416 PetscCall(PetscMalloc1(1, &set_b)); 417 418 set_b[0] = 0; 419 420 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 421 PetscCall(PrintSet(set_a)); 422 PetscCall(AssertSetsEqual(set_a, truth)); 423 424 PetscCall(PetscFree(set_a)); 425 PetscCall(PetscFree(set_b)); 426 427 PetscFunctionReturn(PETSC_SUCCESS); 428 } 429 430 /* Tests functionality when setb is small and setb is large */ 431 PetscErrorCode test_big_small() 432 { 433 PetscInt64 *set_a, *set_b; 434 PetscInt64 truth[] = {2, 1, 8}; 435 PetscMPIInt length = 1; 436 437 PetscFunctionBeginUser; 438 439 PetscCall(Cube(&set_a, 999)); 440 PetscCall(PetscMalloc1(5, &set_b)); 441 442 set_b[0] = 4; 443 set_b[1] = 1; 444 set_b[2] = 4; 445 set_b[3] = 8; 446 set_b[4] = 9; 447 448 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 449 PetscCall(PrintSet(set_a)); 450 PetscCall(AssertSetsEqual(set_a, truth)); 451 452 PetscCall(PetscFree(set_a)); 453 PetscCall(PetscFree(set_b)); 454 455 PetscFunctionReturn(PETSC_SUCCESS); 456 } 457 458 /* Tests functionality when setb is medium sized and setb is large */ 459 PetscErrorCode test_big_moderate() 460 { 461 PetscInt64 *set_a, *set_b; 462 PetscInt64 truth[] = {2, 1, 8}; 463 PetscMPIInt length = 1; 464 465 PetscFunctionBeginUser; 466 467 PetscCall(Cube(&set_a, 999)); 468 PetscCall(Fibonnaci(&set_b, 49)); 469 470 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 471 PetscCall(PrintSet(set_a)); 472 PetscCall(AssertSetsEqual(set_a, truth)); 473 474 PetscCall(PetscFree(set_a)); 475 PetscCall(PetscFree(set_b)); 476 477 PetscFunctionReturn(PETSC_SUCCESS); 478 } 479 480 /* Tests functionality when seta and setb are large, in the opposite 481 order to test_big_big() */ 482 PetscErrorCode test_big_big_reversed() 483 { 484 PetscInt64 *set_a, *set_b; 485 PetscInt64 *truth; 486 PetscMPIInt length = 1; 487 488 PetscFunctionBeginUser; 489 490 PetscCall(Cube(&set_a, 999)); 491 PetscCall(Square(&set_b, 999)); 492 493 PetscGarbageKeySortedIntersect((void *)set_b, (void *)set_a, &length, NULL); 494 PetscCall(PrintSet(set_a)); 495 496 PetscCall(Sixth(&truth, 9)); 497 PetscCall(AssertSetsEqual(set_a, truth)); 498 499 PetscCall(PetscFree(set_a)); 500 PetscCall(PetscFree(set_b)); 501 PetscCall(PetscFree(truth)); 502 503 PetscFunctionReturn(PETSC_SUCCESS); 504 } 505 506 /* Main executes the individual tests in a predefined order */ 507 int main(int argc, char **argv) 508 { 509 PetscFunctionBeginUser; 510 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 511 512 /* Small tests */ 513 /* Test different edge cases with small sets */ 514 PetscCall(test_empty_empty()); 515 PetscCall(test_a_empty()); 516 PetscCall(test_b_empty()); 517 PetscCall(test_identical()); 518 PetscCall(test_disjoint()); 519 PetscCall(test_single_common()); 520 PetscCall(test_issue_1247()); 521 522 /* Big tests */ 523 /* Test different edge cases with big sets */ 524 PetscCall(test_empty_big()); 525 PetscCall(test_small_big()); 526 PetscCall(test_moderate_big()); 527 PetscCall(test_big_big()); 528 PetscCall(test_big_empty()); 529 PetscCall(test_big_small()); 530 PetscCall(test_big_moderate()); 531 PetscCall(test_big_big_reversed()); 532 533 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "ALL PASSED\n")); 534 PetscCall(PetscFinalize()); 535 return 0; 536 } 537 538 /*TEST 539 540 test: 541 suffix: 0 542 543 TEST*/ 544