Actual source code: ex2.c

  1: static char help[] = "Benchmark DMPlexInterpolate.\n\n";

  3: #include <petscdmplex.h>

  5: int main(int argc, char **argv)
  6: {
  7:   DM            dm, dm2;
  8:   PetscLogStage stage;

 10:   PetscFunctionBeginUser;
 11:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 12:   PetscCall(PetscLogStageRegister("Interpolate", &stage));

 14:   PetscCall(DMCreate(PETSC_COMM_WORLD, &dm));
 15:   PetscCall(DMSetType(dm, DMPLEX));
 16:   PetscCall(DMSetFromOptions(dm));
 17:   PetscCall(DMViewFromOptions(dm, NULL, "-init_dm_view"));

 19:   PetscCall(DMPlexUninterpolate(dm, &dm2));
 20:   PetscCall(DMDestroy(&dm));
 21:   dm = dm2;
 22:   PetscCall(DMViewFromOptions(dm, NULL, "-unint_dm_view"));

 24:   PetscCall(PetscLogStagePush(stage));
 25:   PetscCall(DMPlexInterpolate(dm, &dm2));
 26:   PetscCall(PetscLogStagePop());

 28:   PetscCall(DMViewFromOptions(dm2, NULL, "-interp_dm_view"));

 30:   PetscCall(DMDestroy(&dm2));
 31:   PetscCall(DMDestroy(&dm));
 32:   PetscCall(PetscFinalize());
 33:   return 0;
 34: }

 36: /*TEST

 38:   test:
 39:     suffix: 0
 40:     args: -dm_plex_simplex 0

 42: TEST*/