1*c4762a1bSJed Brown static char help[] = "Make a 2D grid of patches and view them\n\n"; 2*c4762a1bSJed Brown 3*c4762a1bSJed Brown /* 4*c4762a1bSJed Brown Serial Test 5*c4762a1bSJed Brown Parallel Test where all zooms are serials 6*c4762a1bSJed Brown Parallel Test where zooms are parallel 7*c4762a1bSJed Brown 8*c4762a1bSJed Brown Return DMPatch from Zoom 9*c4762a1bSJed Brown Override refine from DMPatch to split cells 10*c4762a1bSJed Brown */ 11*c4762a1bSJed Brown #include <petscdmpatch.h> 12*c4762a1bSJed Brown 13*c4762a1bSJed Brown typedef struct { 14*c4762a1bSJed Brown PetscInt debug; /* The debugging level */ 15*c4762a1bSJed Brown PetscInt dim; /* The spatial dimension */ 16*c4762a1bSJed Brown MatStencil patchSize; /* Size of patches */ 17*c4762a1bSJed Brown MatStencil gridSize; /* Size of patch grid */ 18*c4762a1bSJed Brown MatStencil commSize; /* Size of patch comm */ 19*c4762a1bSJed Brown } AppCtx; 20*c4762a1bSJed Brown 21*c4762a1bSJed Brown PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 22*c4762a1bSJed Brown { 23*c4762a1bSJed Brown PetscInt patchSize, commSize, gridSize; 24*c4762a1bSJed Brown PetscErrorCode ierr; 25*c4762a1bSJed Brown 26*c4762a1bSJed Brown PetscFunctionBegin; 27*c4762a1bSJed Brown options->debug = 0; 28*c4762a1bSJed Brown options->dim = 2; 29*c4762a1bSJed Brown patchSize = 0; 30*c4762a1bSJed Brown commSize = 0; 31*c4762a1bSJed Brown gridSize = 0; 32*c4762a1bSJed Brown 33*c4762a1bSJed Brown ierr = PetscOptionsBegin(comm, "", "Patch Test Options", "DMPATCH");CHKERRQ(ierr); 34*c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-debug", "The debugging level", "ex1.c", options->debug, &options->debug, NULL,0);CHKERRQ(ierr); 35*c4762a1bSJed Brown ierr = PetscOptionsRangeInt("-dim", "The spatial dimension", "ex1.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr); 36*c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-patch_size", "The patch size in each dimension", "ex1.c", patchSize, &patchSize, NULL,0);CHKERRQ(ierr); 37*c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-comm_size", "The comm size in each dimension", "ex1.c", commSize, &commSize, NULL,0);CHKERRQ(ierr); 38*c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-grid_size", "The grid size in each dimension", "ex1.c", gridSize, &gridSize, NULL,1);CHKERRQ(ierr); 39*c4762a1bSJed Brown ierr = PetscOptionsEnd(); 40*c4762a1bSJed Brown 41*c4762a1bSJed Brown options->patchSize.i = options->patchSize.j = options->patchSize.k = 1; 42*c4762a1bSJed Brown options->commSize.i = options->commSize.j = options->commSize.k = 1; 43*c4762a1bSJed Brown options->gridSize.i = options->gridSize.j = options->gridSize.k = 1; 44*c4762a1bSJed Brown if (options->dim > 0) {options->patchSize.i = patchSize; options->commSize.i = commSize; options->gridSize.i = gridSize;} 45*c4762a1bSJed Brown if (options->dim > 1) {options->patchSize.j = patchSize; options->commSize.j = commSize; options->gridSize.j = gridSize;} 46*c4762a1bSJed Brown if (options->dim > 2) {options->patchSize.k = patchSize; options->commSize.k = commSize; options->gridSize.k = gridSize;} 47*c4762a1bSJed Brown PetscFunctionReturn(0); 48*c4762a1bSJed Brown }; 49*c4762a1bSJed Brown 50*c4762a1bSJed Brown int main(int argc, char **argv) 51*c4762a1bSJed Brown { 52*c4762a1bSJed Brown DM dm; 53*c4762a1bSJed Brown AppCtx user; /* user-defined work context */ 54*c4762a1bSJed Brown PetscErrorCode ierr; 55*c4762a1bSJed Brown 56*c4762a1bSJed Brown ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 57*c4762a1bSJed Brown ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 58*c4762a1bSJed Brown ierr = DMPatchCreateGrid(PETSC_COMM_WORLD, user.dim, user.patchSize, user.commSize, user.gridSize, &dm);CHKERRQ(ierr); 59*c4762a1bSJed Brown ierr = PetscObjectSetName((PetscObject) dm, "Patch Mesh");CHKERRQ(ierr); 60*c4762a1bSJed Brown ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 61*c4762a1bSJed Brown ierr = DMSetUp(dm);CHKERRQ(ierr); 62*c4762a1bSJed Brown ierr = DMView(dm, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 63*c4762a1bSJed Brown ierr = DMPatchSolve(dm);CHKERRQ(ierr); 64*c4762a1bSJed Brown ierr = DMDestroy(&dm);CHKERRQ(ierr); 65*c4762a1bSJed Brown ierr = PetscFinalize(); 66*c4762a1bSJed Brown return ierr; 67*c4762a1bSJed Brown } 68