1 static char help[] = "Test PetscSFFCompose against some corner cases \n\n"; 2 3 #include <petscsf.h> 4 5 int main(int argc, char **argv) { 6 PetscMPIInt size; 7 PetscSF sfA0, sfA1, sfA2, sfB; 8 PetscInt nroots, nleaves; 9 PetscInt *ilocalA0, *ilocalA1, *ilocalA2, *ilocalB; 10 PetscSFNode *iremoteA0, *iremoteA1, *iremoteA2, *iremoteB; 11 12 PetscFunctionBeginUser; 13 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 14 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 15 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for one MPI process"); 16 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA0)); 17 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA1)); 18 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA2)); 19 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfB)); 20 /* sfA0 */ 21 nroots = 1; 22 nleaves = 0; 23 PetscCall(PetscMalloc1(nleaves, &ilocalA0)); 24 PetscCall(PetscMalloc1(nleaves, &iremoteA0)); 25 PetscCall(PetscSFSetGraph(sfA0, nroots, nleaves, ilocalA0, PETSC_OWN_POINTER, iremoteA0, PETSC_OWN_POINTER)); 26 PetscCall(PetscSFSetUp(sfA0)); 27 PetscCall(PetscObjectSetName((PetscObject)sfA0, "sfA0")); 28 PetscCall(PetscSFView(sfA0, NULL)); 29 /* sfA1 */ 30 nroots = 1; 31 nleaves = 1; 32 PetscCall(PetscMalloc1(nleaves, &ilocalA1)); 33 PetscCall(PetscMalloc1(nleaves, &iremoteA1)); 34 ilocalA1[0] = 1; 35 iremoteA1[0].rank = 0; 36 iremoteA1[0].index = 0; 37 PetscCall(PetscSFSetGraph(sfA1, nroots, nleaves, ilocalA1, PETSC_OWN_POINTER, iremoteA1, PETSC_OWN_POINTER)); 38 PetscCall(PetscSFSetUp(sfA1)); 39 PetscCall(PetscObjectSetName((PetscObject)sfA1, "sfA1")); 40 PetscCall(PetscSFView(sfA1, NULL)); 41 /* sfA2 */ 42 nroots = 1; 43 nleaves = 1; 44 PetscCall(PetscMalloc1(nleaves, &ilocalA2)); 45 PetscCall(PetscMalloc1(nleaves, &iremoteA2)); 46 ilocalA2[0] = 0; 47 iremoteA2[0].rank = 0; 48 iremoteA2[0].index = 0; 49 PetscCall(PetscSFSetGraph(sfA2, nroots, nleaves, ilocalA2, PETSC_OWN_POINTER, iremoteA2, PETSC_OWN_POINTER)); 50 PetscCall(PetscSFSetUp(sfA2)); 51 PetscCall(PetscObjectSetName((PetscObject)sfA2, "sfA2")); 52 PetscCall(PetscSFView(sfA2, NULL)); 53 /* sfB */ 54 nroots = 2; 55 nleaves = 2; 56 PetscCall(PetscMalloc1(nleaves, &ilocalB)); 57 PetscCall(PetscMalloc1(nleaves, &iremoteB)); 58 ilocalB[0] = 100; 59 iremoteB[0].rank = 0; 60 iremoteB[0].index = 0; 61 ilocalB[1] = 101; 62 iremoteB[1].rank = 0; 63 iremoteB[1].index = 1; 64 PetscCall(PetscSFSetGraph(sfB, nroots, nleaves, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER)); 65 PetscCall(PetscSFSetUp(sfB)); 66 PetscCall(PetscObjectSetName((PetscObject)sfB, "sfB")); 67 PetscCall(PetscSFView(sfB, NULL)); 68 /* Test 0 */ 69 { 70 PetscSF sfC; 71 72 PetscCall(PetscSFCompose(sfA0, sfB, &sfC)); 73 PetscCall(PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA0, sfB)")); 74 PetscCall(PetscSFView(sfC, NULL)); 75 PetscCall(PetscSFDestroy(&sfC)); 76 } 77 /* Test 1 */ 78 { 79 PetscSF sfC; 80 81 PetscCall(PetscSFCompose(sfA1, sfB, &sfC)); 82 PetscCall(PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA1, sfB)")); 83 PetscCall(PetscSFView(sfC, NULL)); 84 PetscCall(PetscSFDestroy(&sfC)); 85 } 86 /* Test 2 */ 87 { 88 PetscSF sfC; 89 90 PetscCall(PetscSFCompose(sfA2, sfB, &sfC)); 91 PetscCall(PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA2, sfB)")); 92 PetscCall(PetscSFView(sfC, NULL)); 93 PetscCall(PetscSFDestroy(&sfC)); 94 } 95 PetscCall(PetscSFDestroy(&sfA0)); 96 PetscCall(PetscSFDestroy(&sfA1)); 97 PetscCall(PetscSFDestroy(&sfA2)); 98 PetscCall(PetscSFDestroy(&sfB)); 99 PetscCall(PetscFinalize()); 100 return 0; 101 } 102 103 /*TEST 104 105 test: 106 suffix: 0 107 108 TEST*/ 109