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