1 static char help[] = "Performance tests for DMPlex query operations\n\n"; 2 3 #include <petscdmplex.h> 4 5 typedef struct { 6 PetscInt dim; /* The topological mesh dimension */ 7 PetscBool cellSimplex; /* Flag for simplices */ 8 PetscInt n; /* The number of faces per dimension for mesh */ 9 } AppCtx; 10 11 static PetscErrorCode ProcessOptions(AppCtx *options) 12 { 13 PetscErrorCode ierr; 14 15 PetscFunctionBegin; 16 options->dim = 2; 17 options->cellSimplex = PETSC_TRUE; 18 options->n = 2; 19 20 ierr = PetscOptionsBegin(PETSC_COMM_SELF, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr); 21 ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex21.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr); 22 ierr = PetscOptionsBool("-cellSimplex", "Flag for simplices", "ex21.c", options->cellSimplex, &options->cellSimplex, NULL);CHKERRQ(ierr); 23 ierr = PetscOptionsBoundedInt("-n", "The number of faces per dimension", "ex21.c", options->n, &options->n, NULL,0);CHKERRQ(ierr); 24 ierr = PetscOptionsEnd();CHKERRQ(ierr); 25 PetscFunctionReturn(0); 26 } 27 28 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 29 { 30 PetscInt dim = user->dim; 31 PetscInt n = user->n; 32 PetscBool cellSimplex = user->cellSimplex; 33 PetscInt cells[3]; 34 PetscErrorCode ierr; 35 36 PetscFunctionBegin; 37 cells[0] = cells[1] = cells[2] = n; 38 ierr = DMPlexCreateBoxMesh(comm, dim, cellSimplex, cells, NULL, NULL, NULL, PETSC_FALSE, dm);CHKERRQ(ierr); 39 { 40 DM ddm = NULL; 41 42 ierr = DMPlexDistribute(*dm, 0, NULL, &ddm);CHKERRQ(ierr); 43 if (ddm) { 44 ierr = DMDestroy(dm);CHKERRQ(ierr); 45 *dm = ddm; 46 } 47 } 48 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 49 ierr = DMViewFromOptions(*dm, NULL, "-orig_dm_view");CHKERRQ(ierr); 50 PetscFunctionReturn(0); 51 } 52 53 static PetscErrorCode TestInterpolate(DM dm, AppCtx *user) 54 { 55 DM idm; 56 PetscErrorCode ierr; 57 58 PetscFunctionBegin; 59 ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr); 60 ierr = DMViewFromOptions(idm, NULL, "-interp_dm_view");CHKERRQ(ierr); 61 ierr = DMDestroy(&idm);CHKERRQ(ierr); 62 PetscFunctionReturn(0); 63 } 64 65 int main(int argc, char **argv) 66 { 67 DM dm; 68 AppCtx user; 69 PetscErrorCode ierr; 70 71 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 72 ierr = ProcessOptions(&user);CHKERRQ(ierr); 73 ierr = PetscLogDefaultBegin();CHKERRQ(ierr); 74 ierr = CreateMesh(PETSC_COMM_WORLD, &user, &dm);CHKERRQ(ierr); 75 ierr = TestInterpolate(dm, &user);CHKERRQ(ierr); 76 ierr = DMDestroy(&dm);CHKERRQ(ierr); 77 ierr = PetscFinalize(); 78 return ierr; 79 } 80 81 /*TEST 82 83 test: 84 suffix: 0 85 requires: triangle 86 args: -dim 2 -n 100 87 test: 88 suffix: 1 89 requires: ctetgen 90 args: -dim 3 -n 20 91 92 TEST*/ 93