1 static char help[] = "Mesh Orientation Tutorial\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petscdmplextransform.h> 5 6 typedef struct { 7 PetscBool genArr; /* Generate all possible cell arrangements */ 8 PetscBool refArr; /* Refine all possible cell arrangements */ 9 PetscBool printTable; /* Print the CAyley table */ 10 PetscInt orntBounds[2]; /* Bounds for the orientation check */ 11 PetscInt numOrnt; /* Number of specific orientations specified, or -1 for all orientations */ 12 PetscInt ornts[48]; /* Specific orientations if specified */ 13 PetscInt initOrnt; /* Initial orientation for starting mesh */ 14 } AppCtx; 15 16 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) { 17 PetscInt n = 2; 18 PetscBool flg; 19 20 PetscFunctionBeginUser; 21 options->genArr = PETSC_FALSE; 22 options->refArr = PETSC_FALSE; 23 options->printTable = PETSC_FALSE; 24 options->orntBounds[0] = PETSC_MIN_INT; 25 options->orntBounds[1] = PETSC_MAX_INT; 26 options->numOrnt = -1; 27 options->initOrnt = 0; 28 29 PetscOptionsBegin(comm, "", "Mesh Orientation Tutorials Options", "DMPLEX"); 30 PetscCall(PetscOptionsBool("-gen_arrangements", "Flag for generating all arrangements of the cell", "ex11.c", options->genArr, &options->genArr, NULL)); 31 PetscCall(PetscOptionsBool("-ref_arrangements", "Flag for refining all arrangements of the cell", "ex11.c", options->refArr, &options->refArr, NULL)); 32 PetscCall(PetscOptionsBool("-print_table", "Print the Cayley table", "ex11.c", options->printTable, &options->printTable, NULL)); 33 PetscCall(PetscOptionsIntArray("-ornt_bounds", "Bounds for orientation checks", "ex11.c", options->orntBounds, &n, NULL)); 34 n = 48; 35 PetscCall(PetscOptionsIntArray("-ornts", "Specific orientations for checks", "ex11.c", options->ornts, &n, &flg)); 36 if (flg) { 37 options->numOrnt = n; 38 PetscCall(PetscSortInt(n, options->ornts)); 39 } 40 PetscCall(PetscOptionsInt("-init_ornt", "Initial orientation for starting mesh", "ex11.c", options->initOrnt, &options->initOrnt, NULL)); 41 PetscOptionsEnd(); 42 PetscFunctionReturn(0); 43 } 44 45 static PetscBool ignoreOrnt(AppCtx *user, PetscInt o) { 46 PetscInt loc; 47 PetscErrorCode ierr; 48 49 if (user->numOrnt < 0) return PETSC_FALSE; 50 ierr = PetscFindInt(o, user->numOrnt, user->ornts, &loc); 51 if (loc < 0 || ierr) return PETSC_TRUE; 52 return PETSC_FALSE; 53 } 54 55 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) { 56 PetscFunctionBeginUser; 57 PetscCall(DMCreate(comm, dm)); 58 PetscCall(DMSetType(*dm, DMPLEX)); 59 PetscCall(DMSetFromOptions(*dm)); 60 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 61 PetscFunctionReturn(0); 62 } 63 64 static PetscErrorCode CheckCellVertices(DM dm, PetscInt cell, PetscInt o) { 65 DMPolytopeType ct; 66 const PetscInt *arrVerts; 67 PetscInt *closure = NULL; 68 PetscInt Ncl, cl, Nv, vStart, vEnd, v; 69 MPI_Comm comm; 70 71 PetscFunctionBeginUser; 72 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm)); 73 PetscCall(DMPlexGetCellType(dm, cell, &ct)); 74 PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd)); 75 PetscCall(DMPlexGetTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure)); 76 for (cl = 0, Nv = 0; cl < Ncl * 2; cl += 2) { 77 const PetscInt vertex = closure[cl]; 78 79 if (vertex < vStart || vertex >= vEnd) continue; 80 closure[Nv++] = vertex; 81 } 82 PetscCheck(Nv == DMPolytopeTypeGetNumVertices(ct), comm, PETSC_ERR_ARG_WRONG, "Cell %" PetscInt_FMT " has %" PetscInt_FMT " vertices != %" PetscInt_FMT " vertices in a %s", cell, Nv, DMPolytopeTypeGetNumVertices(ct), DMPolytopeTypes[ct]); 83 arrVerts = DMPolytopeTypeGetVertexArrangment(ct, o); 84 for (v = 0; v < Nv; ++v) { 85 PetscCheck(closure[v] == arrVerts[v] + vStart, comm, PETSC_ERR_ARG_WRONG, "Cell %" PetscInt_FMT " vertex[%" PetscInt_FMT "]: %" PetscInt_FMT " should be %" PetscInt_FMT " for arrangement %" PetscInt_FMT, cell, v, closure[v], arrVerts[v] + vStart, o); 86 } 87 PetscCall(DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure)); 88 PetscFunctionReturn(0); 89 } 90 91 /* Transform cell with group operation o */ 92 static PetscErrorCode ReorientCell(DM dm, PetscInt cell, PetscInt o, PetscBool swapCoords) { 93 DM cdm; 94 Vec coordinates; 95 PetscScalar *coords, *ccoords = NULL; 96 PetscInt *closure = NULL; 97 PetscInt cdim, d, Nc, Ncl, cl, vStart, vEnd, Nv; 98 99 PetscFunctionBegin; 100 /* Change vertex coordinates so that it plots as we expect */ 101 PetscCall(DMGetCoordinateDM(dm, &cdm)); 102 PetscCall(DMGetCoordinateDim(dm, &cdim)); 103 PetscCall(DMGetCoordinatesLocal(dm, &coordinates)); 104 PetscCall(DMPlexVecGetClosure(cdm, NULL, coordinates, cell, &Nc, &ccoords)); 105 /* Reorient cone */ 106 PetscCall(DMPlexOrientPoint(dm, cell, o)); 107 /* Finish resetting coordinates */ 108 if (swapCoords) { 109 PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd)); 110 PetscCall(VecGetArrayWrite(coordinates, &coords)); 111 PetscCall(DMPlexGetTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure)); 112 for (cl = 0, Nv = 0; cl < Ncl * 2; cl += 2) { 113 const PetscInt vertex = closure[cl]; 114 PetscScalar *vcoords; 115 116 if (vertex < vStart || vertex >= vEnd) continue; 117 PetscCall(DMPlexPointLocalRef(cdm, vertex, coords, &vcoords)); 118 for (d = 0; d < cdim; ++d) vcoords[d] = ccoords[Nv * cdim + d]; 119 ++Nv; 120 } 121 PetscCall(DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure)); 122 PetscCall(VecRestoreArrayWrite(coordinates, &coords)); 123 } 124 PetscCall(DMPlexVecRestoreClosure(cdm, NULL, coordinates, cell, &Nc, &ccoords)); 125 PetscFunctionReturn(0); 126 } 127 128 static PetscErrorCode GenerateArrangments(DM dm, AppCtx *user) { 129 DM odm; 130 DMPolytopeType ct; 131 PetscInt No, o; 132 const char *name; 133 134 PetscFunctionBeginUser; 135 if (!user->genArr) PetscFunctionReturn(0); 136 PetscCall(PetscObjectGetName((PetscObject)dm, &name)); 137 PetscCall(DMPlexGetCellType(dm, 0, &ct)); 138 No = DMPolytopeTypeGetNumArrangments(ct) / 2; 139 for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) { 140 if (ignoreOrnt(user, o)) continue; 141 PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &odm)); 142 PetscCall(ReorientCell(odm, 0, o, PETSC_TRUE)); 143 PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "%s orientation %" PetscInt_FMT "\n", name, o)); 144 PetscCall(DMViewFromOptions(odm, NULL, "-gen_dm_view")); 145 PetscCall(CheckCellVertices(odm, 0, o)); 146 PetscCall(DMDestroy(&odm)); 147 } 148 PetscFunctionReturn(0); 149 } 150 151 static PetscErrorCode VerifyCayleyTable(DM dm, AppCtx *user) { 152 DM dm1, dm2; 153 DMPolytopeType ct; 154 const PetscInt *refcone, *cone; 155 PetscInt No, o1, o2, o3, o4; 156 PetscBool equal; 157 const char *name; 158 159 PetscFunctionBeginUser; 160 if (!user->genArr) PetscFunctionReturn(0); 161 PetscCall(PetscObjectGetName((PetscObject)dm, &name)); 162 PetscCall(DMPlexGetCellType(dm, 0, &ct)); 163 PetscCall(DMPlexGetCone(dm, 0, &refcone)); 164 No = DMPolytopeTypeGetNumArrangments(ct) / 2; 165 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Cayley Table for %s\n", DMPolytopeTypes[ct])); 166 for (o1 = PetscMax(-No, user->orntBounds[0]); o1 < PetscMin(No, user->orntBounds[1]); ++o1) { 167 for (o2 = PetscMax(-No, user->orntBounds[0]); o2 < PetscMin(No, user->orntBounds[1]); ++o2) { 168 PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm1)); 169 PetscCall(DMPlexOrientPoint(dm1, 0, o2)); 170 PetscCall(DMPlexCheckFaces(dm1, 0)); 171 PetscCall(DMPlexOrientPoint(dm1, 0, o1)); 172 PetscCall(DMPlexCheckFaces(dm1, 0)); 173 o3 = DMPolytopeTypeComposeOrientation(ct, o1, o2); 174 /* First verification */ 175 PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm2)); 176 PetscCall(DMPlexOrientPoint(dm2, 0, o3)); 177 PetscCall(DMPlexCheckFaces(dm2, 0)); 178 PetscCall(DMPlexEqual(dm1, dm2, &equal)); 179 if (!equal) { 180 PetscCall(DMViewFromOptions(dm1, NULL, "-error_dm_view")); 181 PetscCall(DMViewFromOptions(dm2, NULL, "-error_dm_view")); 182 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cayley table error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " != %" PetscInt_FMT, DMPolytopeTypes[ct], o1, o2, o3); 183 } 184 /* Second verification */ 185 PetscCall(DMPlexGetCone(dm1, 0, &cone)); 186 PetscCall(DMPolytopeGetOrientation(ct, refcone, cone, &o4)); 187 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", ", o4)); 188 PetscCheck(o3 == o4, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cayley table error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " = %" PetscInt_FMT " != %" PetscInt_FMT, DMPolytopeTypes[ct], o1, o2, o3, o4); 189 PetscCall(DMDestroy(&dm1)); 190 PetscCall(DMDestroy(&dm2)); 191 } 192 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 193 } 194 PetscFunctionReturn(0); 195 } 196 197 static PetscErrorCode VerifyInverse(DM dm, AppCtx *user) { 198 DM dm1, dm2; 199 DMPolytopeType ct; 200 const PetscInt *refcone, *cone; 201 PetscInt No, o, oi, o2; 202 PetscBool equal; 203 const char *name; 204 205 PetscFunctionBeginUser; 206 if (!user->genArr) PetscFunctionReturn(0); 207 PetscCall(PetscObjectGetName((PetscObject)dm, &name)); 208 PetscCall(DMPlexGetCellType(dm, 0, &ct)); 209 PetscCall(DMPlexGetCone(dm, 0, &refcone)); 210 No = DMPolytopeTypeGetNumArrangments(ct) / 2; 211 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Inverse table for %s\n", DMPolytopeTypes[ct])); 212 for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) { 213 if (ignoreOrnt(user, o)) continue; 214 oi = DMPolytopeTypeComposeOrientationInv(ct, 0, o); 215 PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm1)); 216 PetscCall(DMPlexOrientPoint(dm1, 0, o)); 217 PetscCall(DMPlexCheckFaces(dm1, 0)); 218 PetscCall(DMPlexOrientPoint(dm1, 0, oi)); 219 PetscCall(DMPlexCheckFaces(dm1, 0)); 220 /* First verification */ 221 PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm2)); 222 PetscCall(DMPlexEqual(dm1, dm2, &equal)); 223 if (!equal) { 224 PetscCall(DMViewFromOptions(dm1, NULL, "-error_dm_view")); 225 PetscCall(DMViewFromOptions(dm2, NULL, "-error_dm_view")); 226 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inverse error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " != 0", DMPolytopeTypes[ct], o, oi); 227 } 228 /* Second verification */ 229 PetscCall(DMPlexGetCone(dm1, 0, &cone)); 230 PetscCall(DMPolytopeGetOrientation(ct, refcone, cone, &o2)); 231 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", ", oi)); 232 PetscCheck(o2 == 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inverse error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " = %" PetscInt_FMT " != 0", DMPolytopeTypes[ct], o, oi, o2); 233 PetscCall(DMDestroy(&dm1)); 234 PetscCall(DMDestroy(&dm2)); 235 } 236 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 237 PetscFunctionReturn(0); 238 } 239 240 /* Suppose that point p has the same arrangement as o from canonical, compare the subcells to canonical subcells */ 241 static PetscErrorCode CheckSubcells(DM dm, DM odm, PetscInt p, PetscInt o, AppCtx *user) { 242 DMPlexTransform tr, otr; 243 DMPolytopeType ct; 244 DMPolytopeType *rct; 245 const PetscInt *cone, *ornt, *ocone, *oornt; 246 PetscInt *rsize, *rcone, *rornt; 247 PetscInt Nct, n, oi, debug = 0; 248 249 PetscFunctionBeginUser; 250 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr)); 251 PetscCall(DMPlexTransformSetDM(tr, dm)); 252 PetscCall(DMPlexTransformSetFromOptions(tr)); 253 PetscCall(DMPlexTransformSetUp(tr)); 254 255 PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)odm), &otr)); 256 PetscCall(DMPlexTransformSetDM(otr, odm)); 257 PetscCall(DMPlexTransformSetFromOptions(otr)); 258 PetscCall(DMPlexTransformSetUp(otr)); 259 260 PetscCall(DMPlexGetCellType(dm, p, &ct)); 261 PetscCall(DMPlexGetCone(dm, p, &cone)); 262 PetscCall(DMPlexGetConeOrientation(dm, p, &ornt)); 263 PetscCall(DMPlexGetCone(odm, p, &ocone)); 264 PetscCall(DMPlexGetConeOrientation(odm, p, &oornt)); 265 oi = DMPolytopeTypeComposeOrientationInv(ct, 0, o); 266 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Orientation %" PetscInt_FMT "\n", oi)); 267 268 PetscCall(DMPlexTransformCellTransform(tr, ct, p, NULL, &Nct, &rct, &rsize, &rcone, &rornt)); 269 for (n = 0; n < Nct; ++n) { 270 DMPolytopeType ctNew = rct[n]; 271 PetscInt r, ro; 272 273 if (debug) PetscCall(PetscPrintf(PETSC_COMM_SELF, " Checking type %s\n", DMPolytopeTypes[ctNew])); 274 for (r = 0; r < rsize[n]; ++r) { 275 const PetscInt *qcone, *qornt, *oqcone, *oqornt; 276 PetscInt pNew, opNew, oo, pr, fo; 277 PetscBool restore = PETSC_TRUE; 278 279 PetscCall(DMPlexTransformGetTargetPoint(tr, ct, ctNew, p, r, &pNew)); 280 PetscCall(DMPlexTransformGetCone(tr, pNew, &qcone, &qornt)); 281 if (debug) { 282 PetscInt c; 283 284 PetscCall(PetscPrintf(PETSC_COMM_SELF, " Checking replica %" PetscInt_FMT " (%" PetscInt_FMT ")\n Original Cone", r, pNew)); 285 for (c = 0; c < DMPolytopeTypeGetConeSize(ctNew); ++c) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT " (%" PetscInt_FMT ")", qcone[c], qornt[c])); 286 PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 287 } 288 for (ro = 0; ro < rsize[n]; ++ro) { 289 PetscBool found; 290 291 PetscCall(DMPlexTransformGetTargetPoint(otr, ct, ctNew, p, ro, &opNew)); 292 PetscCall(DMPlexTransformGetConeOriented(otr, opNew, o, &oqcone, &oqornt)); 293 PetscCall(DMPolytopeMatchOrientation(ctNew, oqcone, qcone, &oo, &found)); 294 if (found) break; 295 PetscCall(DMPlexTransformRestoreCone(otr, pNew, &oqcone, &oqornt)); 296 } 297 if (debug) { 298 PetscInt c; 299 300 PetscCall(PetscPrintf(PETSC_COMM_SELF, " Checking transform replica %" PetscInt_FMT " (%" PetscInt_FMT ") (%" PetscInt_FMT ")\n Transform Cone", ro, opNew, o)); 301 for (c = 0; c < DMPolytopeTypeGetConeSize(ctNew); ++c) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT " (%" PetscInt_FMT ")", oqcone[c], oqornt[c])); 302 PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 303 PetscCall(PetscPrintf(PETSC_COMM_SELF, " Matched %" PetscInt_FMT "\n", oo)); 304 } 305 if (ro == rsize[n]) { 306 /* The tetrahedron has 3 pairs of opposing edges, and any pair can be connected by the interior segment */ 307 if (ct == DM_POLYTOPE_TETRAHEDRON) { 308 /* The segment in a tetrahedron does not map into itself under the group action */ 309 if (ctNew == DM_POLYTOPE_SEGMENT) { 310 restore = PETSC_FALSE; 311 ro = r; 312 oo = 0; 313 } 314 /* The last four interior faces do not map into themselves under the group action */ 315 if (r > 3 && ctNew == DM_POLYTOPE_TRIANGLE) { 316 restore = PETSC_FALSE; 317 ro = r; 318 oo = 0; 319 } 320 /* The last four interior faces do not map into themselves under the group action */ 321 if (r > 3 && ctNew == DM_POLYTOPE_TETRAHEDRON) { 322 restore = PETSC_FALSE; 323 ro = r; 324 oo = 0; 325 } 326 } 327 PetscCheck(!restore, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Unable to find matching %s %" PetscInt_FMT " orientation for cell orientation %" PetscInt_FMT, DMPolytopeTypes[ctNew], r, o); 328 } 329 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", %" PetscInt_FMT ", ", ro, oo)); 330 PetscCall(DMPlexTransformGetSubcellOrientation(tr, ct, p, oi, ctNew, r, 0, &pr, &fo)); 331 if (!user->printTable) { 332 PetscCheck(pr == ro, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Choose wrong replica %" PetscInt_FMT " != %" PetscInt_FMT, pr, ro); 333 PetscCheck(fo == oo, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Choose wrong orientation %" PetscInt_FMT " != %" PetscInt_FMT, fo, oo); 334 } 335 PetscCall(DMPlexTransformRestoreCone(tr, pNew, &qcone, &qornt)); 336 if (restore) PetscCall(DMPlexTransformRestoreCone(otr, pNew, &oqcone, &oqornt)); 337 } 338 if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n")); 339 } 340 PetscCall(DMPlexTransformDestroy(&tr)); 341 PetscCall(DMPlexTransformDestroy(&otr)); 342 PetscFunctionReturn(0); 343 } 344 345 static PetscErrorCode RefineArrangments(DM dm, AppCtx *user) { 346 DM odm, rdm; 347 DMPolytopeType ct; 348 PetscInt No, o; 349 const char *name; 350 351 PetscFunctionBeginUser; 352 if (!user->refArr) PetscFunctionReturn(0); 353 PetscCall(PetscObjectGetName((PetscObject)dm, &name)); 354 PetscCall(DMPlexGetCellType(dm, 0, &ct)); 355 No = DMPolytopeTypeGetNumArrangments(ct) / 2; 356 for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) { 357 if (ignoreOrnt(user, o)) continue; 358 PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &odm)); 359 if (user->initOrnt) PetscCall(ReorientCell(odm, 0, user->initOrnt, PETSC_FALSE)); 360 PetscCall(ReorientCell(odm, 0, o, PETSC_TRUE)); 361 PetscCall(DMViewFromOptions(odm, NULL, "-orig_dm_view")); 362 PetscCall(DMRefine(odm, MPI_COMM_NULL, &rdm)); 363 PetscCall(DMSetFromOptions(rdm)); 364 PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "%s orientation %" PetscInt_FMT "\n", name, o)); 365 PetscCall(DMViewFromOptions(rdm, NULL, "-ref_dm_view")); 366 PetscCall(CheckSubcells(dm, odm, 0, o, user)); 367 PetscCall(DMDestroy(&odm)); 368 PetscCall(DMDestroy(&rdm)); 369 } 370 PetscFunctionReturn(0); 371 } 372 373 int main(int argc, char **argv) { 374 DM dm; 375 AppCtx user; 376 377 PetscFunctionBeginUser; 378 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 379 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user)); 380 PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm)); 381 if (user.initOrnt) { 382 PetscCall(ReorientCell(dm, 0, user.initOrnt, PETSC_FALSE)); 383 PetscCall(DMViewFromOptions(dm, NULL, "-ornt_dm_view")); 384 } 385 PetscCall(GenerateArrangments(dm, &user)); 386 PetscCall(VerifyCayleyTable(dm, &user)); 387 PetscCall(VerifyInverse(dm, &user)); 388 PetscCall(RefineArrangments(dm, &user)); 389 PetscCall(DMDestroy(&dm)); 390 PetscCall(PetscFinalize()); 391 return 0; 392 } 393 394 /*TEST 395 396 testset: 397 args: -dm_coord_space 0 -dm_plex_reference_cell_domain -gen_arrangements \ 398 -gen_dm_view ::ascii_latex -dm_plex_view_tikzscale 0.5 399 400 test: 401 suffix: segment 402 args: -dm_plex_cell segment \ 403 -dm_plex_view_numbers_depth 1,0 -dm_plex_view_colors_depth 1,0 404 405 test: 406 suffix: triangle 407 args: -dm_plex_cell triangle \ 408 -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 409 410 test: 411 suffix: quadrilateral 412 args: -dm_plex_cell quadrilateral \ 413 -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 414 415 test: 416 suffix: tensor_segment 417 args: -dm_plex_cell tensor_quad \ 418 -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 419 420 test: 421 suffix: tetrahedron 422 args: -dm_plex_cell tetrahedron \ 423 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 424 425 test: 426 suffix: hexahedron 427 args: -dm_plex_cell hexahedron \ 428 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 0.3 429 430 test: 431 suffix: triangular_prism 432 args: -dm_plex_cell triangular_prism \ 433 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 434 435 test: 436 suffix: tensor_triangular_prism 437 args: -dm_plex_cell tensor_triangular_prism \ 438 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 439 440 test: 441 suffix: tensor_quadrilateral_prism 442 args: -dm_plex_cell tensor_quadrilateral_prism \ 443 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 444 445 test: 446 suffix: pyramid 447 args: -dm_plex_cell pyramid \ 448 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 449 450 testset: 451 args: -dm_coord_space 0 -dm_plex_reference_cell_domain -ref_arrangements -dm_plex_check_all \ 452 -ref_dm_view ::ascii_latex -dm_plex_view_tikzscale 1.0 453 454 test: 455 suffix: ref_segment 456 args: -dm_plex_cell segment \ 457 -dm_plex_view_numbers_depth 1,0 -dm_plex_view_colors_depth 1,0 458 459 test: 460 suffix: ref_triangle 461 args: -dm_plex_cell triangle \ 462 -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 463 464 test: 465 suffix: ref_quadrilateral 466 args: -dm_plex_cell quadrilateral \ 467 -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 468 469 test: 470 suffix: ref_tensor_segment 471 args: -dm_plex_cell tensor_quad \ 472 -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 473 474 test: 475 suffix: ref_tetrahedron 476 args: -dm_plex_cell tetrahedron \ 477 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 478 479 test: 480 suffix: ref_hexahedron 481 args: -dm_plex_cell hexahedron \ 482 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 483 484 test: 485 suffix: ref_triangular_prism 486 args: -dm_plex_cell triangular_prism \ 487 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 488 489 test: 490 suffix: ref_tensor_triangular_prism 491 args: -dm_plex_cell tensor_triangular_prism \ 492 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 493 494 test: 495 suffix: ref_tensor_quadrilateral_prism 496 args: -dm_plex_cell tensor_quadrilateral_prism \ 497 -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 498 499 # ToBox should recreate the coordinate space since the cell shape changes 500 testset: 501 args: -dm_coord_space 0 -dm_plex_transform_type refine_tobox -ref_arrangements -dm_plex_check_all 502 503 test: 504 suffix: tobox_triangle 505 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle \ 506 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 507 508 test: 509 suffix: tobox_tensor_segment 510 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad \ 511 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 512 513 test: 514 suffix: tobox_tetrahedron 515 args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron \ 516 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0 517 518 test: 519 suffix: tobox_triangular_prism 520 args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism \ 521 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0 522 523 test: 524 suffix: tobox_tensor_triangular_prism 525 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism \ 526 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0 527 528 test: 529 suffix: tobox_tensor_quadrilateral_prism 530 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism \ 531 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0 532 533 testset: 534 args: -dm_coord_space 0 -dm_plex_transform_type refine_alfeld -ref_arrangements -dm_plex_check_all 535 536 test: 537 suffix: alfeld_triangle 538 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle \ 539 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 540 541 test: 542 suffix: alfeld_tetrahedron 543 args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron \ 544 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0 545 546 testset: 547 args: -dm_plex_transform_type refine_sbr -ref_arrangements -dm_plex_check_all 548 549 # This splits edge 1 of the triangle, and reflects about the added edge 550 test: 551 suffix: sbr_triangle_0 552 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -ornts -2,0 \ 553 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 554 555 # This splits edge 0 of the triangle, and reflects about the added edge 556 test: 557 suffix: sbr_triangle_1 558 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -init_ornt 1 -ornts -3,0 \ 559 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 560 561 # This splits edge 2 of the triangle, and reflects about the added edge 562 test: 563 suffix: sbr_triangle_2 564 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -init_ornt 2 -ornts -1,0 \ 565 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 566 567 # This splits edges 1 and 2 of the triangle 568 test: 569 suffix: sbr_triangle_3 570 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -ornts 0 \ 571 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 572 573 # This splits edges 1 and 0 of the triangle 574 test: 575 suffix: sbr_triangle_4 576 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -ornts 0 \ 577 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 578 579 # This splits edges 0 and 1 of the triangle 580 test: 581 suffix: sbr_triangle_5 582 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -init_ornt 1 -ornts 0 \ 583 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 584 585 # This splits edges 0 and 2 of the triangle 586 test: 587 suffix: sbr_triangle_6 588 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -init_ornt 1 -ornts 0 \ 589 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 590 591 # This splits edges 2 and 0 of the triangle 592 test: 593 suffix: sbr_triangle_7 594 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -init_ornt 2 -ornts 0 \ 595 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 596 597 # This splits edges 2 and 1 of the triangle 598 test: 599 suffix: sbr_triangle_8 600 args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -init_ornt 2 -ornts 0 \ 601 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 602 603 testset: 604 args: -dm_plex_transform_type refine_boundary_layer -dm_plex_transform_bl_splits 2 -ref_arrangements -dm_plex_check_all 605 606 test: 607 suffix: bl_tensor_segment 608 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad \ 609 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0 610 611 # The subcell check is broken because at orientation 3, the internal triangles do not get properly permuted for the check 612 test: 613 suffix: bl_tensor_triangular_prism 614 requires: TODO 615 args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism \ 616 -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0 617 618 TEST*/ 619