1 static char help[] = "Tests adaptive refinement using DMForest, and uses HDF5.\n\n"; 2 3 #include <petscdmforest.h> 4 #include <petscdmplex.h> 5 #include <petscviewerhdf5.h> 6 7 int main (int argc, char **argv) 8 { 9 DM base, forest, plex; 10 PetscSection s; 11 PetscViewer viewer; 12 Vec g = NULL, g2 = NULL; 13 PetscReal nrm; 14 PetscBool adapt = PETSC_FALSE, userSection = PETSC_FALSE; 15 PetscInt vStart, vEnd, v, i; 16 PetscErrorCode ierr; 17 18 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 19 ierr = PetscOptionsGetBool(NULL, NULL, "-adapt", &adapt, NULL);CHKERRQ(ierr); 20 ierr = PetscOptionsGetBool(NULL, NULL, "-user_section", &userSection, NULL);CHKERRQ(ierr); 21 22 /* Create a base DMPlex mesh */ 23 ierr = DMCreate(PETSC_COMM_WORLD, &base);CHKERRQ(ierr); 24 ierr = DMSetType(base, DMPLEX);CHKERRQ(ierr); 25 ierr = DMSetFromOptions(base);CHKERRQ(ierr); 26 ierr = DMViewFromOptions(base, NULL, "-dm_view");CHKERRQ(ierr); 27 28 /* Covert Plex mesh to Forest and destroy base */ 29 ierr = DMCreate(PETSC_COMM_WORLD, &forest);CHKERRQ(ierr); 30 ierr = DMSetType(forest, DMP4EST);CHKERRQ(ierr); 31 ierr = DMForestSetBaseDM(forest, base);CHKERRQ(ierr); 32 ierr = DMSetUp(forest);CHKERRQ(ierr); 33 ierr = DMDestroy(&base);CHKERRQ(ierr); 34 ierr = DMViewFromOptions(forest, NULL, "-dm_view");CHKERRQ(ierr); 35 36 if (adapt) { 37 /* Adaptively refine the cell 0 of the mesh */ 38 for (i = 0; i < 3; ++i) { 39 DM postforest; 40 DMLabel adaptLabel = NULL; 41 42 ierr = DMLabelCreate(PETSC_COMM_SELF, "adapt", &adaptLabel);CHKERRQ(ierr); 43 ierr = DMLabelSetValue(adaptLabel, 0, DM_ADAPT_REFINE);CHKERRQ(ierr); 44 ierr = DMForestTemplate(forest, PETSC_COMM_WORLD, &postforest);CHKERRQ(ierr); 45 ierr = DMForestSetAdaptivityLabel(postforest, adaptLabel);CHKERRQ(ierr); 46 ierr = DMLabelDestroy(&adaptLabel);CHKERRQ(ierr); 47 ierr = DMSetUp(postforest);CHKERRQ(ierr); 48 ierr = DMDestroy(&forest);CHKERRQ(ierr); 49 forest = postforest; 50 } 51 } else { 52 /* Adaptively refine all cells of the mesh */ 53 PetscInt cStart, cEnd, c; 54 55 for (i = 0; i < 3; ++i) { 56 DM postforest; 57 DMLabel adaptLabel = NULL; 58 59 ierr = DMLabelCreate(PETSC_COMM_SELF, "adapt", &adaptLabel);CHKERRQ(ierr); 60 61 ierr = DMForestGetCellChart(forest, &cStart, &cEnd);CHKERRQ(ierr); 62 for (c = cStart; c < cEnd; ++c) { 63 ierr = DMLabelSetValue(adaptLabel, c, DM_ADAPT_REFINE);CHKERRQ(ierr); 64 } 65 66 ierr = DMForestTemplate(forest, PETSC_COMM_WORLD, &postforest);CHKERRQ(ierr); 67 ierr = DMForestSetAdaptivityLabel(postforest, adaptLabel);CHKERRQ(ierr); 68 ierr = DMLabelDestroy(&adaptLabel);CHKERRQ(ierr); 69 ierr = DMSetUp(postforest);CHKERRQ(ierr); 70 ierr = DMDestroy(&forest);CHKERRQ(ierr); 71 forest = postforest; 72 } 73 } 74 ierr = DMViewFromOptions(forest, NULL, "-dm_view");CHKERRQ(ierr); 75 76 /* Setup the section*/ 77 if (userSection) { 78 ierr = DMConvert(forest, DMPLEX, &plex);CHKERRQ(ierr); 79 ierr = DMViewFromOptions(plex, NULL, "-dm_view");CHKERRQ(ierr); 80 ierr = DMPlexGetDepthStratum(plex, 0, &vStart, &vEnd);CHKERRQ(ierr); 81 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD, "Vertices [%D, %D)\n", vStart, vEnd);CHKERRQ(ierr); 82 ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD, NULL);CHKERRQ(ierr); 83 ierr = PetscSectionCreate(PetscObjectComm((PetscObject) forest), &s);CHKERRQ(ierr); 84 ierr = PetscSectionSetNumFields(s, 1);CHKERRQ(ierr); 85 ierr = PetscSectionSetChart(s, vStart, vEnd);CHKERRQ(ierr); 86 for (v = vStart; v < vEnd; ++v) { 87 ierr = PetscSectionSetDof(s, v, 1);CHKERRQ(ierr); 88 ierr = PetscSectionSetFieldDof(s, v, 0, 1);CHKERRQ(ierr); 89 } 90 ierr = PetscSectionSetUp(s);CHKERRQ(ierr); 91 ierr = DMSetLocalSection(forest, s);CHKERRQ(ierr); 92 ierr = PetscObjectViewFromOptions((PetscObject) s, NULL, "-my_section_view");CHKERRQ(ierr); 93 ierr = PetscSectionDestroy(&s);CHKERRQ(ierr); 94 ierr = DMDestroy(&plex);CHKERRQ(ierr); 95 } else { 96 PetscFE fe; 97 PetscInt dim; 98 99 ierr = DMGetDimension(forest, &dim);CHKERRQ(ierr); 100 ierr = PetscFECreateLagrange(PETSC_COMM_SELF, dim, 1, PETSC_FALSE, 1, PETSC_DETERMINE, &fe);CHKERRQ(ierr); 101 ierr = DMAddField(forest, NULL, (PetscObject) fe);CHKERRQ(ierr); 102 ierr = PetscFEDestroy(&fe);CHKERRQ(ierr); 103 ierr = DMCreateDS(forest);CHKERRQ(ierr); 104 } 105 106 /* Create the global vector*/ 107 ierr = DMCreateGlobalVector(forest, &g);CHKERRQ(ierr); 108 ierr = PetscObjectSetName((PetscObject) g, "g");CHKERRQ(ierr); 109 ierr = VecSet(g, 1.0);CHKERRQ(ierr); 110 111 /* Test global to local*/ 112 Vec l; 113 ierr = DMCreateLocalVector(forest, &l);CHKERRQ(ierr); 114 ierr = VecZeroEntries(l);CHKERRQ(ierr); 115 ierr = DMGlobalToLocal(forest, g, INSERT_VALUES, l);CHKERRQ(ierr); 116 ierr = VecZeroEntries(g);CHKERRQ(ierr); 117 ierr = DMLocalToGlobal(forest, l, INSERT_VALUES, g);CHKERRQ(ierr); 118 ierr = VecDestroy(&l);CHKERRQ(ierr); 119 120 /* Save a vector*/ 121 ierr = PetscViewerHDF5Open(PETSC_COMM_WORLD, "forestHDF.h5", FILE_MODE_WRITE, &viewer);CHKERRQ(ierr); 122 ierr = VecView(g, viewer);CHKERRQ(ierr); 123 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 124 125 /* Load another vector to load into*/ 126 ierr = DMCreateGlobalVector(forest, &g2);CHKERRQ(ierr); 127 ierr = PetscObjectSetName((PetscObject) g2, "g");CHKERRQ(ierr); 128 ierr = VecZeroEntries(g2);CHKERRQ(ierr); 129 130 /* Load a vector*/ 131 ierr = PetscViewerHDF5Open(PETSC_COMM_WORLD, "forestHDF.h5", FILE_MODE_READ, &viewer);CHKERRQ(ierr); 132 ierr = VecLoad(g2, viewer);CHKERRQ(ierr); 133 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 134 135 /* Check if the data is the same*/ 136 ierr = VecAXPY(g2, -1.0, g);CHKERRQ(ierr); 137 ierr = VecNorm(g2, NORM_INFINITY, &nrm);CHKERRQ(ierr); 138 if (PetscAbsReal(nrm) > PETSC_SMALL) SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Invalid difference norm %g", (double) nrm); 139 140 ierr = VecDestroy(&g);CHKERRQ(ierr); 141 ierr = VecDestroy(&g2);CHKERRQ(ierr); 142 ierr = DMDestroy(&forest);CHKERRQ(ierr); 143 ierr = PetscFinalize(); 144 return ierr; 145 } 146 147 /*TEST 148 149 build: 150 requires: hdf5 p4est 151 152 test: 153 suffix: 0 154 nsize: {{1 2 5}} 155 args: -adapt -dm_plex_simplex 0 -dm_plex_box_faces 2,2 156 157 TEST*/ 158