1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 3 /*@C 4 DMPlexInvertCell - Flips cell orientations since Plex stores some of them internally with outward normals. 5 6 Input Parameters: 7 + cellType - The cell type 8 - cone - The incoming cone 9 10 Output Parameter: 11 . cone - The inverted cone (in-place) 12 13 Level: developer 14 15 .seealso: `DMPlexGenerate()` 16 @*/ 17 PetscErrorCode DMPlexInvertCell(DMPolytopeType cellType, PetscInt cone[]) { 18 #define SWAPCONE(cone, i, j) \ 19 do { \ 20 PetscInt _cone_tmp; \ 21 _cone_tmp = (cone)[i]; \ 22 (cone)[i] = (cone)[j]; \ 23 (cone)[j] = _cone_tmp; \ 24 } while (0) 25 26 PetscFunctionBegin; 27 switch (cellType) { 28 case DM_POLYTOPE_POINT: break; 29 case DM_POLYTOPE_SEGMENT: break; 30 case DM_POLYTOPE_POINT_PRISM_TENSOR: break; 31 case DM_POLYTOPE_TRIANGLE: break; 32 case DM_POLYTOPE_QUADRILATERAL: break; 33 case DM_POLYTOPE_SEG_PRISM_TENSOR: SWAPCONE(cone, 2, 3); break; 34 case DM_POLYTOPE_TETRAHEDRON: SWAPCONE(cone, 0, 1); break; 35 case DM_POLYTOPE_HEXAHEDRON: SWAPCONE(cone, 1, 3); break; 36 case DM_POLYTOPE_TRI_PRISM: SWAPCONE(cone, 1, 2); break; 37 case DM_POLYTOPE_TRI_PRISM_TENSOR: break; 38 case DM_POLYTOPE_QUAD_PRISM_TENSOR: break; 39 case DM_POLYTOPE_PYRAMID: SWAPCONE(cone, 1, 3); break; 40 default: break; 41 } 42 PetscFunctionReturn(0); 43 #undef SWAPCONE 44 } 45 46 /*@C 47 DMPlexReorderCell - Flips cell orientations since Plex stores some of them internally with outward normals. 48 49 Input Parameters: 50 + dm - The DMPlex object 51 . cell - The cell 52 - cone - The incoming cone 53 54 Output Parameter: 55 . cone - The reordered cone (in-place) 56 57 Level: developer 58 59 .seealso: `DMPlexGenerate()` 60 @*/ 61 PetscErrorCode DMPlexReorderCell(DM dm, PetscInt cell, PetscInt cone[]) { 62 DMPolytopeType cellType; 63 64 PetscFunctionBegin; 65 PetscCall(DMPlexGetCellType(dm, cell, &cellType)); 66 PetscCall(DMPlexInvertCell(cellType, cone)); 67 PetscFunctionReturn(0); 68 } 69 70 /*@C 71 DMPlexTriangleSetOptions - Set the options used for the Triangle mesh generator 72 73 Not Collective 74 75 Inputs Parameters: 76 + dm - The DMPlex object 77 - opts - The command line options 78 79 Level: developer 80 81 .seealso: `DMPlexTetgenSetOptions()`, `DMPlexGenerate()` 82 @*/ 83 PetscErrorCode DMPlexTriangleSetOptions(DM dm, const char *opts) { 84 DM_Plex *mesh = (DM_Plex *)dm->data; 85 86 PetscFunctionBegin; 87 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 88 PetscValidCharPointer(opts, 2); 89 PetscCall(PetscFree(mesh->triangleOpts)); 90 PetscCall(PetscStrallocpy(opts, &mesh->triangleOpts)); 91 PetscFunctionReturn(0); 92 } 93 94 /*@C 95 DMPlexTetgenSetOptions - Set the options used for the Tetgen mesh generator 96 97 Not Collective 98 99 Inputs Parameters: 100 + dm - The DMPlex object 101 - opts - The command line options 102 103 Level: developer 104 105 .seealso: `DMPlexTriangleSetOptions()`, `DMPlexGenerate()` 106 @*/ 107 PetscErrorCode DMPlexTetgenSetOptions(DM dm, const char *opts) { 108 DM_Plex *mesh = (DM_Plex *)dm->data; 109 110 PetscFunctionBegin; 111 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 112 PetscValidCharPointer(opts, 2); 113 PetscCall(PetscFree(mesh->tetgenOpts)); 114 PetscCall(PetscStrallocpy(opts, &mesh->tetgenOpts)); 115 PetscFunctionReturn(0); 116 } 117 118 /*@C 119 DMPlexGenerate - Generates a mesh. 120 121 Not Collective 122 123 Input Parameters: 124 + boundary - The DMPlex boundary object 125 . name - The mesh generation package name 126 - interpolate - Flag to create intermediate mesh elements 127 128 Output Parameter: 129 . mesh - The DMPlex object 130 131 Options Database: 132 + -dm_plex_generate <name> - package to generate mesh, for example, triangle, ctetgen or tetgen 133 - -dm_generator <name> - package to generate mesh, for example, triangle, ctetgen or tetgen 134 135 Level: intermediate 136 137 .seealso: `DMPlexCreate()`, `DMRefine()` 138 @*/ 139 PetscErrorCode DMPlexGenerate(DM boundary, const char name[], PetscBool interpolate, DM *mesh) { 140 DMGeneratorFunctionList fl; 141 char genname[PETSC_MAX_PATH_LEN]; 142 const char *suggestions; 143 PetscInt dim; 144 PetscBool flg; 145 146 PetscFunctionBegin; 147 PetscValidHeaderSpecific(boundary, DM_CLASSID, 1); 148 PetscValidLogicalCollectiveBool(boundary, interpolate, 3); 149 PetscCall(DMGetDimension(boundary, &dim)); 150 PetscCall(PetscOptionsGetString(((PetscObject)boundary)->options, ((PetscObject)boundary)->prefix, "-dm_generator", genname, sizeof(genname), &flg)); 151 if (flg) name = genname; 152 else { 153 PetscCall(PetscOptionsGetString(((PetscObject)boundary)->options, ((PetscObject)boundary)->prefix, "-dm_plex_generate", genname, sizeof(genname), &flg)); 154 if (flg) name = genname; 155 } 156 157 fl = DMGenerateList; 158 if (name) { 159 while (fl) { 160 PetscCall(PetscStrcmp(fl->name, name, &flg)); 161 if (flg) { 162 PetscCall((*fl->generate)(boundary, interpolate, mesh)); 163 PetscFunctionReturn(0); 164 } 165 fl = fl->next; 166 } 167 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Grid generator %s not registered; you may need to add --download-%s to your ./configure options", name, name); 168 } else { 169 while (fl) { 170 if (boundary->dim == fl->dim) { 171 PetscCall((*fl->generate)(boundary, interpolate, mesh)); 172 PetscFunctionReturn(0); 173 } 174 fl = fl->next; 175 } 176 suggestions = ""; 177 if (boundary->dim + 1 == 2) suggestions = " You may need to add --download-triangle to your ./configure options"; 178 else if (boundary->dim + 1 == 3) suggestions = " You may need to add --download-ctetgen or --download-tetgen in your ./configure options"; 179 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No grid generator of dimension %" PetscInt_FMT " registered%s", boundary->dim + 1, suggestions); 180 } 181 } 182