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