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