1 static char help[] = "Create a box mesh with DMMoab and test defining a tag on the mesh\n\n"; 2 3 #include <petscdmmoab.h> 4 5 typedef struct { 6 DM dm; /* DM implementation using the MOAB interface */ 7 PetscBool debug; /* The debugging level */ 8 PetscLogEvent createMeshEvent; 9 /* Domain and mesh definition */ 10 PetscInt dim; /* The topological mesh dimension */ 11 PetscInt nele; /* Elements in each dimension */ 12 PetscInt degree; /* Degree of refinement */ 13 PetscBool simplex; /* Use simplex elements */ 14 PetscInt nlevels; /* Number of levels in mesh hierarchy */ 15 PetscInt nghost; /* Number of ghost layers in the mesh */ 16 char input_file[PETSC_MAX_PATH_LEN]; /* Import mesh from file */ 17 char output_file[PETSC_MAX_PATH_LEN]; /* Output mesh file name */ 18 PetscBool write_output; /* Write output mesh and data to file */ 19 } AppCtx; 20 21 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 22 { 23 PetscErrorCode ierr; 24 25 PetscFunctionBegin; 26 options->debug = PETSC_FALSE; 27 options->nlevels = 1; 28 options->nghost = 1; 29 options->dim = 2; 30 options->nele = 5; 31 options->degree = 2; 32 options->simplex = PETSC_FALSE; 33 options->write_output = PETSC_FALSE; 34 options->input_file[0] = '\0'; 35 ierr = PetscStrcpy(options->output_file,"ex3.h5m");CHKERRQ(ierr); 36 37 ierr = PetscOptionsBegin(comm, "", "Uniform Mesh Refinement Options", "DMMOAB");CHKERRQ(ierr); 38 ierr = PetscOptionsBool("-debug", "Enable debug messages", "ex2.cxx", options->debug, &options->debug, NULL);CHKERRQ(ierr); 39 ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex3.cxx", options->dim, &options->dim, NULL,0,3);CHKERRQ(ierr); 40 ierr = PetscOptionsBoundedInt("-n", "The number of elements in each dimension", "ex3.cxx", options->nele, &options->nele, NULL,1);CHKERRQ(ierr); 41 ierr = PetscOptionsBoundedInt("-levels", "Number of levels in the hierarchy", "ex3.cxx", options->nlevels, &options->nlevels, NULL,0);CHKERRQ(ierr); 42 ierr = PetscOptionsBoundedInt("-degree", "Number of degrees at each level of refinement", "ex3.cxx", options->degree, &options->degree, NULL,0);CHKERRQ(ierr); 43 ierr = PetscOptionsBoundedInt("-ghost", "Number of ghost layers in the mesh", "ex3.cxx", options->nghost, &options->nghost, NULL,0);CHKERRQ(ierr); 44 ierr = PetscOptionsBool("-simplex", "Create simplices instead of tensor product elements", "ex3.cxx", options->simplex, &options->simplex, NULL);CHKERRQ(ierr); 45 ierr = PetscOptionsString("-input", "The input mesh file", "ex3.cxx", options->input_file, options->input_file, PETSC_MAX_PATH_LEN, NULL);CHKERRQ(ierr); 46 ierr = PetscOptionsString("-io", "Write out the mesh and solution that is defined on it (Default H5M format)", "ex3.cxx", options->output_file, options->output_file, PETSC_MAX_PATH_LEN, &options->write_output);CHKERRQ(ierr); 47 ierr = PetscOptionsEnd(); 48 49 ierr = PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);CHKERRQ(ierr); 50 PetscFunctionReturn(0); 51 }; 52 53 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user) 54 { 55 size_t len; 56 PetscMPIInt rank; 57 PetscErrorCode ierr; 58 59 PetscFunctionBegin; 60 ierr = PetscLogEventBegin(user->createMeshEvent,0,0,0,0);CHKERRQ(ierr); 61 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 62 ierr = PetscStrlen(user->input_file, &len);CHKERRQ(ierr); 63 if (len) { 64 if (user->debug) PetscPrintf(comm, "Loading mesh from file: %s and creating the coarse level DM object.\n",user->input_file); 65 ierr = DMMoabLoadFromFile(comm, user->dim, user->nghost, user->input_file, "", &user->dm);CHKERRQ(ierr); 66 } 67 else { 68 if (user->debug) { 69 PetscPrintf(comm, "Creating a %D-dimensional structured %s mesh of %Dx%Dx%D in memory and creating a DM object.\n",user->dim,(user->simplex?"simplex":"regular"),user->nele,user->nele,user->nele); 70 } 71 ierr = DMMoabCreateBoxMesh(comm, user->dim, user->simplex, NULL, user->nele, user->nghost, &user->dm);CHKERRQ(ierr); 72 } 73 74 ierr = PetscObjectSetName((PetscObject)user->dm, "Coarse Mesh");CHKERRQ(ierr); 75 ierr = PetscLogEventEnd(user->createMeshEvent,0,0,0,0);CHKERRQ(ierr); 76 PetscFunctionReturn(0); 77 } 78 79 int main(int argc, char **argv) 80 { 81 AppCtx user; /* user-defined work context */ 82 MPI_Comm comm; 83 PetscInt i; 84 Mat R; 85 DM *dmhierarchy; 86 PetscInt *degrees; 87 PetscBool createR; 88 PetscErrorCode ierr; 89 90 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 91 comm = PETSC_COMM_WORLD; 92 createR = PETSC_FALSE; 93 94 ierr = ProcessOptions(comm, &user);CHKERRQ(ierr); 95 ierr = CreateMesh(comm, &user);CHKERRQ(ierr); 96 ierr = DMSetFromOptions(user.dm);CHKERRQ(ierr); 97 98 /* SetUp the data structures for DMMOAB */ 99 ierr = DMSetUp(user.dm);CHKERRQ(ierr); 100 101 ierr = PetscMalloc(sizeof(DM)*(user.nlevels+1),&dmhierarchy); 102 for (i=0; i<=user.nlevels; i++) dmhierarchy[i] = NULL; 103 104 // coarsest grid = 0 105 // finest grid = nlevels 106 dmhierarchy[0] = user.dm; 107 PetscObjectReference((PetscObject)user.dm); 108 109 if (user.nlevels) { 110 ierr = PetscMalloc1(user.nlevels, °rees);CHKERRQ(ierr); 111 for (i=0; i < user.nlevels; i++) degrees[i] = user.degree; 112 if (user.debug) PetscPrintf(comm, "Generate the MOAB mesh hierarchy with %D levels.\n", user.nlevels); 113 ierr = DMMoabGenerateHierarchy(user.dm,user.nlevels,degrees);CHKERRQ(ierr); 114 115 PetscBool usehierarchy=PETSC_FALSE; 116 if (usehierarchy) { 117 ierr = DMRefineHierarchy(user.dm,user.nlevels,&dmhierarchy[1]);CHKERRQ(ierr); 118 } 119 else { 120 if (user.debug) { 121 PetscPrintf(PETSC_COMM_WORLD, "Level %D\n", 0); 122 ierr = DMView(user.dm, 0);CHKERRQ(ierr); 123 } 124 for (i=1; i<=user.nlevels; i++) { 125 if (user.debug) PetscPrintf(PETSC_COMM_WORLD, "Level %D\n", i); 126 ierr = DMRefine(dmhierarchy[i-1],MPI_COMM_NULL,&dmhierarchy[i]);CHKERRQ(ierr); 127 if (createR) { 128 ierr = DMCreateInterpolation(dmhierarchy[i-1],dmhierarchy[i],&R,NULL);CHKERRQ(ierr); 129 } 130 if (user.debug) { 131 ierr = DMView(dmhierarchy[i], 0);CHKERRQ(ierr); 132 if (createR) { 133 ierr = MatView(R,0);CHKERRQ(ierr); 134 } 135 } 136 /* Solvers could now set operator "R" to the multigrid PC object for level i 137 PCMGSetInterpolation(pc,i,R) 138 */ 139 if (createR) { 140 ierr = MatDestroy(&R);CHKERRQ(ierr); 141 } 142 } 143 } 144 } 145 146 if (user.write_output) { 147 if (user.debug) PetscPrintf(comm, "Output mesh hierarchy to file: %s.\n",user.output_file); 148 ierr = DMMoabOutput(dmhierarchy[user.nlevels],(const char*)user.output_file,"");CHKERRQ(ierr); 149 } 150 151 for (i=0; i<=user.nlevels; i++) { 152 ierr = DMDestroy(&dmhierarchy[i]);CHKERRQ(ierr); 153 } 154 ierr = PetscFree(degrees);CHKERRQ(ierr); 155 ierr = PetscFree(dmhierarchy);CHKERRQ(ierr); 156 ierr = DMDestroy(&user.dm);CHKERRQ(ierr); 157 ierr = PetscFinalize(); 158 return 0; 159 } 160 161 /*TEST 162 163 build: 164 requires: moab 165 166 test: 167 args: -debug -n 2 -dim 2 -levels 2 -simplex 168 filter: grep -v "DM_0x" 169 170 test: 171 args: -debug -n 2 -dim 3 -levels 2 172 filter: grep -v "DM_0x" 173 suffix: 1_2 174 175 test: 176 args: -debug -n 2 -dim 3 -ghost 1 -levels 2 177 filter: grep -v "DM_0x" 178 nsize: 2 179 suffix: 2_1 180 181 TEST*/ 182