Actual source code: ex10.c

  1: static char help[] = "Test PetscSFFCompose against some corner cases \n\n";

  3: #include <petscsf.h>

  5: int main(int argc, char **argv)
  6: {
  7:   PetscMPIInt  size;
  8:   PetscSF      sfA0, sfA1, sfA2, sfB;
  9:   PetscInt     nroots, nleaves;
 10:   PetscInt    *ilocalA0, *ilocalA1, *ilocalA2, *ilocalB;
 11:   PetscSFNode *iremoteA0, *iremoteA1, *iremoteA2, *iremoteB;

 14:   PetscInitialize(&argc, &argv, NULL, help);
 15:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 17:   PetscSFCreate(PETSC_COMM_WORLD, &sfA0);
 18:   PetscSFCreate(PETSC_COMM_WORLD, &sfA1);
 19:   PetscSFCreate(PETSC_COMM_WORLD, &sfA2);
 20:   PetscSFCreate(PETSC_COMM_WORLD, &sfB);
 21:   /* sfA0 */
 22:   nroots  = 1;
 23:   nleaves = 0;
 24:   PetscMalloc1(nleaves, &ilocalA0);
 25:   PetscMalloc1(nleaves, &iremoteA0);
 26:   PetscSFSetGraph(sfA0, nroots, nleaves, ilocalA0, PETSC_OWN_POINTER, iremoteA0, PETSC_OWN_POINTER);
 27:   PetscSFSetUp(sfA0);
 28:   PetscObjectSetName((PetscObject)sfA0, "sfA0");
 29:   PetscSFView(sfA0, NULL);
 30:   /* sfA1 */
 31:   nroots  = 1;
 32:   nleaves = 1;
 33:   PetscMalloc1(nleaves, &ilocalA1);
 34:   PetscMalloc1(nleaves, &iremoteA1);
 35:   ilocalA1[0]        = 1;
 36:   iremoteA1[0].rank  = 0;
 37:   iremoteA1[0].index = 0;
 38:   PetscSFSetGraph(sfA1, nroots, nleaves, ilocalA1, PETSC_OWN_POINTER, iremoteA1, PETSC_OWN_POINTER);
 39:   PetscSFSetUp(sfA1);
 40:   PetscObjectSetName((PetscObject)sfA1, "sfA1");
 41:   PetscSFView(sfA1, NULL);
 42:   /* sfA2 */
 43:   nroots  = 1;
 44:   nleaves = 1;
 45:   PetscMalloc1(nleaves, &ilocalA2);
 46:   PetscMalloc1(nleaves, &iremoteA2);
 47:   ilocalA2[0]        = 0;
 48:   iremoteA2[0].rank  = 0;
 49:   iremoteA2[0].index = 0;
 50:   PetscSFSetGraph(sfA2, nroots, nleaves, ilocalA2, PETSC_OWN_POINTER, iremoteA2, PETSC_OWN_POINTER);
 51:   PetscSFSetUp(sfA2);
 52:   PetscObjectSetName((PetscObject)sfA2, "sfA2");
 53:   PetscSFView(sfA2, NULL);
 54:   /* sfB */
 55:   nroots  = 2;
 56:   nleaves = 2;
 57:   PetscMalloc1(nleaves, &ilocalB);
 58:   PetscMalloc1(nleaves, &iremoteB);
 59:   ilocalB[0]        = 100;
 60:   iremoteB[0].rank  = 0;
 61:   iremoteB[0].index = 0;
 62:   ilocalB[1]        = 101;
 63:   iremoteB[1].rank  = 0;
 64:   iremoteB[1].index = 1;
 65:   PetscSFSetGraph(sfB, nroots, nleaves, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER);
 66:   PetscSFSetUp(sfB);
 67:   PetscObjectSetName((PetscObject)sfB, "sfB");
 68:   PetscSFView(sfB, NULL);
 69:   /* Test 0 */
 70:   {
 71:     PetscSF sfC;

 73:     PetscSFCompose(sfA0, sfB, &sfC);
 74:     PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA0, sfB)");
 75:     PetscSFView(sfC, NULL);
 76:     PetscSFDestroy(&sfC);
 77:   }
 78:   /* Test 1 */
 79:   {
 80:     PetscSF sfC;

 82:     PetscSFCompose(sfA1, sfB, &sfC);
 83:     PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA1, sfB)");
 84:     PetscSFView(sfC, NULL);
 85:     PetscSFDestroy(&sfC);
 86:   }
 87:   /* Test 2 */
 88:   {
 89:     PetscSF sfC;

 91:     PetscSFCompose(sfA2, sfB, &sfC);
 92:     PetscObjectSetName((PetscObject)sfC, "PetscSFCompose(sfA2, sfB)");
 93:     PetscSFView(sfC, NULL);
 94:     PetscSFDestroy(&sfC);
 95:   }
 96:   PetscSFDestroy(&sfA0);
 97:   PetscSFDestroy(&sfA1);
 98:   PetscSFDestroy(&sfA2);
 99:   PetscSFDestroy(&sfB);
100:   PetscFinalize();
101:   return 0;
102: }

104: /*TEST

106:    test:
107:      suffix: 0

109: TEST*/