1c4762a1bSJed Brown static char help[] = "Tests for creation of hybrid meshes\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown /* TODO 4c4762a1bSJed Brown - Propagate hybridSize with distribution 5c4762a1bSJed Brown - Test with multiple fault segments 6c4762a1bSJed Brown - Test with embedded fault 7c4762a1bSJed Brown - Test with multiple faults 8c4762a1bSJed Brown - Move over all PyLith tests? 9c4762a1bSJed Brown */ 10c4762a1bSJed Brown 11c4762a1bSJed Brown #include <petscdmplex.h> 12c4762a1bSJed Brown /* List of test meshes 13c4762a1bSJed Brown 14c4762a1bSJed Brown Triangle 15c4762a1bSJed Brown -------- 16c4762a1bSJed Brown Test 0: 17c4762a1bSJed Brown Two triangles sharing a face 18c4762a1bSJed Brown 19c4762a1bSJed Brown 4 20c4762a1bSJed Brown / | \ 21c4762a1bSJed Brown 8 | 9 22c4762a1bSJed Brown / | \ 23c4762a1bSJed Brown 2 0 7 1 5 24c4762a1bSJed Brown \ | / 25c4762a1bSJed Brown 6 | 10 26c4762a1bSJed Brown \ | / 27c4762a1bSJed Brown 3 28c4762a1bSJed Brown 29c4762a1bSJed Brown should become two triangles separated by a zero-volume cell with 4 vertices 30c4762a1bSJed Brown 31c4762a1bSJed Brown 5--16--8 4--12--6 3 32c4762a1bSJed Brown / | | \ / | | | \ 33c4762a1bSJed Brown 11 | | 12 9 | | | 4 34c4762a1bSJed Brown / | | \ / | | | \ 35c4762a1bSJed Brown 3 0 10 2 14 1 6 2 0 8 1 10 6 0 1 36c4762a1bSJed Brown \ | | / \ | | | / 37c4762a1bSJed Brown 9 | | 13 7 | | | 5 38c4762a1bSJed Brown \ | | / \ | | | / 39c4762a1bSJed Brown 4--15--7 3--11--5 2 40c4762a1bSJed Brown 41c4762a1bSJed Brown Test 1: 42c4762a1bSJed Brown Four triangles sharing two faces which are oriented against each other 43c4762a1bSJed Brown 44c4762a1bSJed Brown 9 45c4762a1bSJed Brown / \ 46c4762a1bSJed Brown / \ 47c4762a1bSJed Brown 17 2 16 48c4762a1bSJed Brown / \ 49c4762a1bSJed Brown / \ 50c4762a1bSJed Brown 8-----15----5 51c4762a1bSJed Brown \ /|\ 52c4762a1bSJed Brown \ / | \ 53c4762a1bSJed Brown 18 3 12 | 14 54c4762a1bSJed Brown \ / | \ 55c4762a1bSJed Brown \ / | \ 56c4762a1bSJed Brown 4 0 11 1 7 57c4762a1bSJed Brown \ | / 58c4762a1bSJed Brown \ | / 59c4762a1bSJed Brown 10 | 13 60c4762a1bSJed Brown \ | / 61c4762a1bSJed Brown \|/ 62c4762a1bSJed Brown 6 63c4762a1bSJed Brown 64c4762a1bSJed Brown Fault mesh 65c4762a1bSJed Brown 66c4762a1bSJed Brown 0 --> 0 67c4762a1bSJed Brown 1 --> 1 68c4762a1bSJed Brown 2 --> 2 69c4762a1bSJed Brown 3 --> 3 70c4762a1bSJed Brown 4 --> 5 71c4762a1bSJed Brown 5 --> 6 72c4762a1bSJed Brown 6 --> 8 73c4762a1bSJed Brown 7 --> 11 74c4762a1bSJed Brown 8 --> 15 75c4762a1bSJed Brown 76c4762a1bSJed Brown 2 77c4762a1bSJed Brown | 78c4762a1bSJed Brown 6----8----4 79c4762a1bSJed Brown | | 80c4762a1bSJed Brown 3 | 81c4762a1bSJed Brown 0-7-1 82c4762a1bSJed Brown | 83c4762a1bSJed Brown | 84c4762a1bSJed Brown 5 85c4762a1bSJed Brown 86c4762a1bSJed Brown should become four triangles separated by two zero-volume cells with 4 vertices 87c4762a1bSJed Brown 88c4762a1bSJed Brown 11 89c4762a1bSJed Brown / \ 90c4762a1bSJed Brown / \ 91c4762a1bSJed Brown / \ 92c4762a1bSJed Brown 22 2 21 93c4762a1bSJed Brown / \ 94c4762a1bSJed Brown / \ 95c4762a1bSJed Brown 10-----20------7 96c4762a1bSJed Brown 28 | 5 26/ \ 97c4762a1bSJed Brown 14----25----12 \ 98c4762a1bSJed Brown \ /| |\ 99c4762a1bSJed Brown \ / | | \ 100c4762a1bSJed Brown 23 3 17 | | 19 101c4762a1bSJed Brown \ / | | \ 102c4762a1bSJed Brown \ / | | \ 103c4762a1bSJed Brown 6 0 24 4 16 1 9 104c4762a1bSJed Brown \ | | / 105c4762a1bSJed Brown \ | | / 106c4762a1bSJed Brown 15 | | 18 107c4762a1bSJed Brown \ | | / 108c4762a1bSJed Brown \| |/ 109c4762a1bSJed Brown 13---8 110c4762a1bSJed Brown 27 111c4762a1bSJed Brown 112c4762a1bSJed Brown Tetrahedron 113c4762a1bSJed Brown ----------- 114c4762a1bSJed Brown Test 0: 115c4762a1bSJed Brown Two tets sharing a face 116c4762a1bSJed Brown 117c4762a1bSJed Brown cell 5 _______ cell 118c4762a1bSJed Brown 0 / | \ \ 1 119c4762a1bSJed Brown 16 | 18 22 120c4762a1bSJed Brown /8 19 10\ \ 121c4762a1bSJed Brown 2-15-|----4--21--6 122c4762a1bSJed Brown \ 9| 7 / / 123c4762a1bSJed Brown 14 | 17 20 124c4762a1bSJed Brown \ | / / 125c4762a1bSJed Brown 3------- 126c4762a1bSJed Brown 127c4762a1bSJed Brown should become two tetrahedrons separated by a zero-volume cell with 3 faces/3 edges/6 vertices 128c4762a1bSJed Brown 129c4762a1bSJed Brown cell 6 ___36___10______ cell 130c4762a1bSJed Brown 0 / | \ |\ \ 1 131c4762a1bSJed Brown 24 | 26 | 32 30 132c4762a1bSJed Brown /12 27 14\ 33 \ \ 133c4762a1bSJed Brown 3-23-|----5--35-|---9--29--7 134c4762a1bSJed Brown \ 13| 11/ |18 / / 135c4762a1bSJed Brown 22 | 25 | 31 28 136c4762a1bSJed Brown \ | / |/ / 137c4762a1bSJed Brown 4----34----8------ 138c4762a1bSJed Brown cell 2 139c4762a1bSJed Brown 140c4762a1bSJed Brown In parallel, 141c4762a1bSJed Brown 142c4762a1bSJed Brown cell 5 ___28____8 4______ cell 143c4762a1bSJed Brown 0 / | \ |\ |\ \ 0 144c4762a1bSJed Brown 19 | 21 | 24 | 13 6 11 145c4762a1bSJed Brown /10 22 12\ 25 \ |8 \ \ 146c4762a1bSJed Brown 2-18-|----4--27-|---7 14 3--10--1 147c4762a1bSJed Brown \ 11| 9 / |13 / | / / 148c4762a1bSJed Brown 17 | 20 | 23 | 12 5 9 149c4762a1bSJed Brown \ | / |/ |/ / 150c4762a1bSJed Brown 3----26----6 2------ 151c4762a1bSJed Brown cell 1 152c4762a1bSJed Brown 153c4762a1bSJed Brown Test 1: 154c4762a1bSJed Brown Four tets sharing two faces 155c4762a1bSJed Brown 156c4762a1bSJed Brown Cells: 0-3,4-5 157c4762a1bSJed Brown Vertices: 6-15 158c4762a1bSJed Brown Faces: 16-29,30-34 159c4762a1bSJed Brown Edges: 35-52,53-56 160c4762a1bSJed Brown 161c4762a1bSJed Brown Quadrilateral 162c4762a1bSJed Brown ------------- 163c4762a1bSJed Brown Test 0: 164c4762a1bSJed Brown Two quads sharing a face 165c4762a1bSJed Brown 166c4762a1bSJed Brown 5--10---4--14---7 167c4762a1bSJed Brown | | | 168c4762a1bSJed Brown 11 0 9 1 13 169c4762a1bSJed Brown | | | 170c4762a1bSJed Brown 2---8---3--12---6 171c4762a1bSJed Brown 172c4762a1bSJed Brown should become two quads separated by a zero-volume cell with 4 vertices 173c4762a1bSJed Brown 174c4762a1bSJed Brown 6--13---5-20-10--17---8 5--10---4-14--7 4---7---2 175c4762a1bSJed Brown | | | | | | | | | 176c4762a1bSJed Brown 14 0 12 2 18 1 16 11 0 9 1 12 8 0 6 177c4762a1bSJed Brown | | | | | | | | | 178c4762a1bSJed Brown 3--11---4-19--9--15---7 2---8---3-13--6 3---5---1 179c4762a1bSJed Brown 180c4762a1bSJed Brown Test 1: 181c4762a1bSJed Brown 182c4762a1bSJed Brown Original mesh with 9 cells, 183c4762a1bSJed Brown 184c4762a1bSJed Brown 9 ----10 ----11 ----12 185c4762a1bSJed Brown | | | | 186c4762a1bSJed Brown | | | | 187c4762a1bSJed Brown | | | | 188c4762a1bSJed Brown | | | | 189c4762a1bSJed Brown 13 ----14 ----15 ----16 190c4762a1bSJed Brown | | | | 191c4762a1bSJed Brown | | | | 192c4762a1bSJed Brown | | | | 193c4762a1bSJed Brown | | | | 194c4762a1bSJed Brown 17 ----18 ----19 ----20 195c4762a1bSJed Brown | | | | 196c4762a1bSJed Brown | | | | 197c4762a1bSJed Brown | | | | 198c4762a1bSJed Brown | | | | 199c4762a1bSJed Brown 21 ----22 ----23 ----24 200c4762a1bSJed Brown 201c4762a1bSJed Brown After first fault, 202c4762a1bSJed Brown 203c4762a1bSJed Brown 12 ----13 ----14-28 ----15 204c4762a1bSJed Brown | | | | | 205c4762a1bSJed Brown | 0 | 1 | 9| 2 | 206c4762a1bSJed Brown | | | | | 207c4762a1bSJed Brown | | | | | 208c4762a1bSJed Brown 16 ----17 ----18-29 ----19 209c4762a1bSJed Brown | | | | | 210c4762a1bSJed Brown | 3 | 4 |10| 5 | 211c4762a1bSJed Brown | | | | | 212c4762a1bSJed Brown | | | | | 213c4762a1bSJed Brown 20 ----21-----22-30 ----23 214c4762a1bSJed Brown | | | \--11- | 215c4762a1bSJed Brown | 6 | 7 | 8 | 216c4762a1bSJed Brown | | | | 217c4762a1bSJed Brown | | | | 218c4762a1bSJed Brown 24 ----25 ----26--------27 219c4762a1bSJed Brown 220c4762a1bSJed Brown After second fault, 221c4762a1bSJed Brown 222c4762a1bSJed Brown 14 ----15 ----16-30 ----17 223c4762a1bSJed Brown | | | | | 224c4762a1bSJed Brown | 0 | 1 | 9| 2 | 225c4762a1bSJed Brown | | | | | 226c4762a1bSJed Brown | | | | | 227c4762a1bSJed Brown 18 ----19 ----20-31 ----21 228c4762a1bSJed Brown | | | | | 229c4762a1bSJed Brown | 3 | 4 |10| 5 | 230c4762a1bSJed Brown | | | | | 231c4762a1bSJed Brown | | | | | 232c4762a1bSJed Brown 33 ----34-----24-32 ----25 233c4762a1bSJed Brown | 12 | 13 / | \-11-- | 234c4762a1bSJed Brown 22 ----23---/ | | 235c4762a1bSJed Brown | | 7 | 8 | 236c4762a1bSJed Brown | 6 | | | 237c4762a1bSJed Brown | | | | 238c4762a1bSJed Brown | | | | 239c4762a1bSJed Brown 26 ----27 ----28--------29 240c4762a1bSJed Brown 241c4762a1bSJed Brown Hexahedron 242c4762a1bSJed Brown ---------- 243c4762a1bSJed Brown Test 0: 244c4762a1bSJed Brown Two hexes sharing a face 245c4762a1bSJed Brown 246c4762a1bSJed Brown cell 9-----31------8-----42------13 cell 247c4762a1bSJed Brown 0 /| /| /| 1 248c4762a1bSJed Brown 32 | 15 30| 21 41| 249c4762a1bSJed Brown / | / | / | 250c4762a1bSJed Brown 6-----29------7-----40------12 | 251c4762a1bSJed Brown | | 18 | | 24 | | 252c4762a1bSJed Brown | 36 | 35 | 44 253c4762a1bSJed Brown |19 | |17 | |23 | 254c4762a1bSJed Brown 33 | 16 34 | 22 43 | 255c4762a1bSJed Brown | 5-----27--|---4-----39--|---11 256c4762a1bSJed Brown | / | / | / 257c4762a1bSJed Brown | 28 14 | 26 20 | 38 258c4762a1bSJed Brown |/ |/ |/ 259c4762a1bSJed Brown 2-----25------3-----37------10 260c4762a1bSJed Brown 261c4762a1bSJed Brown should become two hexes separated by a zero-volume cell with 8 vertices 262c4762a1bSJed Brown 263c4762a1bSJed Brown cell 2 264c4762a1bSJed Brown cell 10-----41------9-----62------18----52------14 cell 265c4762a1bSJed Brown 0 /| /| /| /| 1 266c4762a1bSJed Brown 42 | 20 40| 32 56| 26 51| 267c4762a1bSJed Brown / | / | / | / | 268c4762a1bSJed Brown 7-----39------8-----61------17--|-50------13 | 269c4762a1bSJed Brown | | 23 | | | | 29 | | 270c4762a1bSJed Brown | 46 | 45 | 58 | 54 271c4762a1bSJed Brown |24 | |22 | |30 | |28 | 272c4762a1bSJed Brown 43 | 21 44 | 57 | 27 53 | 273c4762a1bSJed Brown | 6-----37--|---5-----60--|---16----49--|---12 274c4762a1bSJed Brown | / | / | / | / 275c4762a1bSJed Brown | 38 19 | 36 31 | 55 25 | 48 276c4762a1bSJed Brown |/ |/ |/ |/ 277c4762a1bSJed Brown 3-----35------4-----59------15----47------11 278c4762a1bSJed Brown 279c4762a1bSJed Brown In parallel, 280c4762a1bSJed Brown 281c4762a1bSJed Brown cell 2 282c4762a1bSJed Brown cell 9-----31------8-----44------13 8----20------4 cell 283c4762a1bSJed Brown 0 /| /| /| /| /| 1 284c4762a1bSJed Brown 32 | 15 30| 22 38| 24 | 10 19| 285c4762a1bSJed Brown / | / | / | / | / | 286c4762a1bSJed Brown 6-----29------7-----43------12 | 7----18------3 | 287c4762a1bSJed Brown | | 18 | | | | | | 13 | | 288c4762a1bSJed Brown | 36 | 35 | 40 | 26 | 22 289c4762a1bSJed Brown |19 | |17 | |20 | |14 | |12 | 290c4762a1bSJed Brown 33 | 16 34 | 39 | 25 | 11 21 | 291c4762a1bSJed Brown | 5-----27--|---4-----42--|---11 | 6----17--|---2 292c4762a1bSJed Brown | / | / | / | / | / 293c4762a1bSJed Brown | 28 14 | 26 21 | 37 |23 9 | 16 294c4762a1bSJed Brown |/ |/ |/ |/ |/ 295c4762a1bSJed Brown 2-----25------3-----41------10 5----15------1 296c4762a1bSJed Brown 297c4762a1bSJed Brown Test 1: 298c4762a1bSJed Brown 299c4762a1bSJed Brown */ 300c4762a1bSJed Brown 301c4762a1bSJed Brown typedef struct { 302c4762a1bSJed Brown PetscInt debug; /* The debugging level */ 303c4762a1bSJed Brown PetscInt dim; /* The topological mesh dimension */ 304c4762a1bSJed Brown PetscBool cellSimplex; /* Use simplices or hexes */ 305c4762a1bSJed Brown PetscBool testPartition; /* Use a fixed partitioning for testing */ 306c4762a1bSJed Brown PetscInt testNum; /* The particular mesh to test */ 307c4762a1bSJed Brown } AppCtx; 308c4762a1bSJed Brown 309c4762a1bSJed Brown PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 310c4762a1bSJed Brown { 311c4762a1bSJed Brown PetscErrorCode ierr; 312c4762a1bSJed Brown 313c4762a1bSJed Brown PetscFunctionBegin; 314c4762a1bSJed Brown options->debug = 0; 315c4762a1bSJed Brown options->dim = 2; 316c4762a1bSJed Brown options->cellSimplex = PETSC_TRUE; 317c4762a1bSJed Brown options->testPartition = PETSC_TRUE; 318c4762a1bSJed Brown options->testNum = 0; 319c4762a1bSJed Brown 320c4762a1bSJed Brown ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr); 321c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-debug", "The debugging level", "ex5.c", options->debug, &options->debug, NULL,0);CHKERRQ(ierr); 322c4762a1bSJed Brown ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex5.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr); 323c4762a1bSJed Brown ierr = PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex5.c", options->cellSimplex, &options->cellSimplex, NULL);CHKERRQ(ierr); 324c4762a1bSJed Brown ierr = PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex5.c", options->testPartition, &options->testPartition, NULL);CHKERRQ(ierr); 325c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-test_num", "The particular mesh to test", "ex5.c", options->testNum, &options->testNum, NULL,0);CHKERRQ(ierr); 326c4762a1bSJed Brown ierr = PetscOptionsEnd(); 327c4762a1bSJed Brown PetscFunctionReturn(0); 328c4762a1bSJed Brown } 329c4762a1bSJed Brown 330c4762a1bSJed Brown PetscErrorCode CreateSimplex_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 331c4762a1bSJed Brown { 332c4762a1bSJed Brown DM idm; 333c4762a1bSJed Brown PetscInt p; 334c4762a1bSJed Brown PetscMPIInt rank; 335c4762a1bSJed Brown PetscErrorCode ierr; 336c4762a1bSJed Brown 337c4762a1bSJed Brown PetscFunctionBegin; 338c4762a1bSJed Brown ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 339c4762a1bSJed Brown if (!rank) { 340c4762a1bSJed Brown switch (testNum) { 341c4762a1bSJed Brown case 0: 342c4762a1bSJed Brown { 343c4762a1bSJed Brown PetscInt numPoints[2] = {4, 2}; 344c4762a1bSJed Brown PetscInt coneSize[6] = {3, 3, 0, 0, 0, 0}; 345c4762a1bSJed Brown PetscInt cones[6] = {2, 3, 4, 5, 4, 3}; 346c4762a1bSJed Brown PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0}; 347c4762a1bSJed Brown PetscScalar vertexCoords[8] = {-0.5, 0.5, 0.0, 0.0, 0.0, 1.0, 0.5, 0.5}; 348c4762a1bSJed Brown PetscInt markerPoints[8] = {2, 1, 3, 1, 4, 1, 5, 1}; 349c4762a1bSJed Brown PetscInt faultPoints[2] = {3, 4}; 350c4762a1bSJed Brown 351c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 352c4762a1bSJed Brown for (p = 0; p < 4; ++p) {ierr = DMSetLabelValue(*dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);} 353c4762a1bSJed Brown for (p = 0; p < 2; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 354c4762a1bSJed Brown } 355c4762a1bSJed Brown break; 356c4762a1bSJed Brown case 1: 357c4762a1bSJed Brown { 358c4762a1bSJed Brown PetscInt numPoints[2] = {6, 4}; 359c4762a1bSJed Brown PetscInt coneSize[10] = {3, 3, 3, 3, 0, 0, 0, 0, 0, 0}; 360c4762a1bSJed Brown PetscInt cones[12] = {4, 6, 5, 5, 6, 7, 8, 5, 9, 8, 4, 5}; 361c4762a1bSJed Brown PetscInt coneOrientations[12] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 362c4762a1bSJed Brown PetscScalar vertexCoords[12] = {-1.0, 0.0, 0.0, 1.0, 0.0, -1.0, 1.0, 0.0, -2.0, 1.0, -1.0, 2.0}; 363c4762a1bSJed Brown PetscInt markerPoints[6] = {4, 1, 6, 1, 8, 1}; 364c4762a1bSJed Brown PetscInt faultPoints[3] = {5, 6, 8}; 365c4762a1bSJed Brown 366c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 367c4762a1bSJed Brown for (p = 0; p < 3; ++p) {ierr = DMSetLabelValue(*dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);} 368c4762a1bSJed Brown for (p = 0; p < 3; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 369c4762a1bSJed Brown } 370c4762a1bSJed Brown break; 371c4762a1bSJed Brown default: 372c4762a1bSJed Brown SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 373c4762a1bSJed Brown } 374c4762a1bSJed Brown } else { 375c4762a1bSJed Brown PetscInt numPoints[3] = {0, 0, 0}; 376c4762a1bSJed Brown 377c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 378c4762a1bSJed Brown ierr = DMCreateLabel(*dm, "fault");CHKERRQ(ierr); 379c4762a1bSJed Brown } 380c4762a1bSJed Brown ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 381c4762a1bSJed Brown ierr = DMViewFromOptions(idm, NULL, "-in_dm_view");CHKERRQ(ierr); 382c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 383c4762a1bSJed Brown *dm = idm; 384c4762a1bSJed Brown PetscFunctionReturn(0); 385c4762a1bSJed Brown } 386c4762a1bSJed Brown 387c4762a1bSJed Brown PetscErrorCode CreateSimplex_3D(MPI_Comm comm, AppCtx *user, DM dm) 388c4762a1bSJed Brown { 389c4762a1bSJed Brown PetscInt depth = 3, testNum = user->testNum, p; 390c4762a1bSJed Brown PetscMPIInt rank; 391c4762a1bSJed Brown PetscErrorCode ierr; 392c4762a1bSJed Brown 393c4762a1bSJed Brown PetscFunctionBegin; 394c4762a1bSJed Brown ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 395c4762a1bSJed Brown if (!rank) { 396c4762a1bSJed Brown switch (testNum) { 397c4762a1bSJed Brown case 0: 398c4762a1bSJed Brown { 399c4762a1bSJed Brown PetscInt numPoints[4] = {5, 7, 9, 2}; 400c4762a1bSJed Brown PetscInt coneSize[23] = {4, 4, 0, 0, 0, 0, 0, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2}; 401c4762a1bSJed Brown PetscInt cones[47] = {7, 8, 9, 10, 11, 10, 13, 12, 15, 17, 14, 16, 18, 15, 14, 19, 16, 17, 18, 19, 17, 21, 20, 18, 22, 21, 22, 19, 20, 2, 3, 2, 4, 2, 5, 3, 4, 4, 5, 5, 3, 3, 6, 4, 6, 5, 6}; 402c4762a1bSJed Brown PetscInt coneOrientations[47] = {0, 0, 0, 0, 0, -3, 2, 2, 0, -2, -2, 0, -2, -2, 0, -2, -2, 0, 0, 0, 0, 0, -2, 0, 0, -2, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 403c4762a1bSJed Brown PetscScalar vertexCoords[15] = {0.0, 0.0, -0.5, 0.0, -0.5, 0.0, 1.0, 0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.5}; 404c4762a1bSJed Brown PetscInt markerPoints[20] = {2, 1, 3, 1, 4, 1, 5, 1, 14, 1, 15, 1, 16, 1, 17, 1, 18, 1, 19, 1}; 405c4762a1bSJed Brown PetscInt faultPoints[3] = {3, 4, 5}; 406c4762a1bSJed Brown 407c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 408c4762a1bSJed Brown for (p = 0; p < 10; ++p) { 409c4762a1bSJed Brown ierr = DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr); 410c4762a1bSJed Brown } 411c4762a1bSJed Brown for(p = 0; p < 3; ++p) { 412c4762a1bSJed Brown ierr = DMSetLabelValue(dm, "fault", faultPoints[p], 1);CHKERRQ(ierr); 413c4762a1bSJed Brown } 414c4762a1bSJed Brown } 415c4762a1bSJed Brown break; 416c4762a1bSJed Brown case 1: 417c4762a1bSJed Brown { 418c4762a1bSJed Brown PetscInt numPoints[4] = {6, 13, 12, 4}; 419c4762a1bSJed Brown PetscInt coneSize[35] = {4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2}; 420c4762a1bSJed Brown PetscInt cones[78] = {10, 11, 12, 13, 10, 15, 16, 14, 17, 18, 14, 19, 20, 13, 19, 21, 22, 23, 24, 25, 26, 22, 24, 27, 25, 23, 26, 27, 28, 29, 23, 24, 30, 28, 22, 29, 30, 31, 32, 28, 29, 33, 31, 32, 33, 23, 26, 34, 33, 34, 27, 32, 6, 5, 5, 7, 7, 6, 6, 4, 4, 5, 7, 4, 7, 9, 9, 5, 6, 9, 9, 8, 8, 7, 5, 8, 4, 8}; 421c4762a1bSJed Brown PetscInt coneOrientations[78] = { 0, 0, 0, 0, -3, 1, 0, 2, 0, 0, -1, 0, 0, -1, -2, 0, 0, 0, 0, 0, 0, -2, -2, 0, -2, -2, -2, -2, 0, 0, 0, 0, 0, -2, 0, -2, -2, 0, 0, 0, 0, 0, -2, -2, -2, 0, -2, 0, -2, -2, -2, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 422c4762a1bSJed Brown PetscScalar vertexCoords[18] = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0, 0.0, -1.0, 1.0, 0.0, 0.0}; 423c4762a1bSJed Brown PetscInt markerPoints[14] = {5, 1, 6, 1, 7, 1, 10, 1, 22, 1, 23, 1, 24, 1}; 424c4762a1bSJed Brown PetscInt faultPoints[4] = {5, 6, 7, 8}; 425c4762a1bSJed Brown 426c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 427c4762a1bSJed Brown for (p = 0; p < 7; ++p) { 428c4762a1bSJed Brown ierr = DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr); 429c4762a1bSJed Brown } 430c4762a1bSJed Brown for(p = 0; p < 4; ++p) { 431c4762a1bSJed Brown ierr = DMSetLabelValue(dm, "fault", faultPoints[p], 1);CHKERRQ(ierr); 432c4762a1bSJed Brown } 433c4762a1bSJed Brown } 434c4762a1bSJed Brown break; 435c4762a1bSJed Brown default: 436c4762a1bSJed Brown SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 437c4762a1bSJed Brown } 438c4762a1bSJed Brown } else { 439c4762a1bSJed Brown PetscInt numPoints[4] = {0, 0, 0, 0}; 440c4762a1bSJed Brown 441c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 442c4762a1bSJed Brown ierr = DMCreateLabel(dm, "fault");CHKERRQ(ierr); 443c4762a1bSJed Brown } 444c4762a1bSJed Brown PetscFunctionReturn(0); 445c4762a1bSJed Brown } 446c4762a1bSJed Brown 447c4762a1bSJed Brown PetscErrorCode CreateQuad_2D(MPI_Comm comm, PetscInt testNum, DM *dm) 448c4762a1bSJed Brown { 449c4762a1bSJed Brown DM idm; 450c4762a1bSJed Brown PetscInt p; 451c4762a1bSJed Brown PetscMPIInt rank; 452c4762a1bSJed Brown PetscErrorCode ierr; 453c4762a1bSJed Brown 454c4762a1bSJed Brown PetscFunctionBegin; 455c4762a1bSJed Brown ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 456c4762a1bSJed Brown if (!rank) { 457c4762a1bSJed Brown switch (testNum) { 458c4762a1bSJed Brown case 0: 459c4762a1bSJed Brown case 2: 460c4762a1bSJed Brown { 461c4762a1bSJed Brown PetscInt numPoints[2] = {6, 2}; 462c4762a1bSJed Brown PetscInt coneSize[8] = {4, 4, 0, 0, 0, 0, 0, 0}; 463c4762a1bSJed Brown PetscInt cones[8] = {2, 3, 4, 5, 3, 6, 7, 4}; 464c4762a1bSJed Brown PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0}; 465c4762a1bSJed Brown PetscScalar vertexCoords[12] = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0, 0.5, 1.0}; 466c4762a1bSJed Brown PetscInt markerPoints[12] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1}; 467c4762a1bSJed Brown PetscInt faultPoints[2] = {3, 4}; 468c4762a1bSJed Brown 469c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 470c4762a1bSJed Brown for (p = 0; p < 6; ++p) {ierr = DMSetLabelValue(*dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);} 471c4762a1bSJed Brown if (testNum == 0) for (p = 0; p < 2; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 472c4762a1bSJed Brown if (testNum == 2) for (p = 0; p < 2; ++p) {ierr = DMSetLabelValue(*dm, "pfault", faultPoints[p], 1);CHKERRQ(ierr);} 473c4762a1bSJed Brown } 474c4762a1bSJed Brown break; 475c4762a1bSJed Brown case 1: 476c4762a1bSJed Brown { 477c4762a1bSJed Brown PetscInt numPoints[2] = {16, 9}; 478c4762a1bSJed Brown PetscInt coneSize[25] = {4, 4, 4, 4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}; 479c4762a1bSJed Brown PetscInt cones[36] = {9, 13, 14, 10, 480c4762a1bSJed Brown 10, 14, 15, 11, 481c4762a1bSJed Brown 11, 15, 16, 12, 482c4762a1bSJed Brown 13, 17, 18, 14, 483c4762a1bSJed Brown 14, 18, 19, 15, 484c4762a1bSJed Brown 15, 19, 20, 16, 485c4762a1bSJed Brown 17, 21, 22, 18, 486c4762a1bSJed Brown 18, 22, 23, 19, 487c4762a1bSJed Brown 19, 23, 24, 20}; 488c4762a1bSJed Brown PetscInt coneOrientations[36] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 489c4762a1bSJed Brown PetscScalar vertexCoords[32] = {-3.0, 3.0, -1.0, 3.0, 1.0, 3.0, 3.0, 3.0, -3.0, 1.0, -1.0, 1.0, 1.0, 1.0, 3.0, 1.0, 490c4762a1bSJed Brown -3.0, -1.0, -1.0, -1.0, 1.0, -1.0, 3.0, -1.0, -3.0, -3.0, -1.0, -3.0, 1.0, -3.0, 3.0, -3.0}; 491c4762a1bSJed Brown PetscInt faultPoints[3] = {11, 15, 19}; 492c4762a1bSJed Brown PetscInt fault2Points[2] = {17, 18}; 493c4762a1bSJed Brown 494c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 495c4762a1bSJed Brown for (p = 0; p < 3; ++p) {ierr = DMSetLabelValue(*dm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 496c4762a1bSJed Brown for (p = 0; p < 2; ++p) {ierr = DMSetLabelValue(*dm, "fault2", fault2Points[p], 1);CHKERRQ(ierr);} 497c4762a1bSJed Brown } 498c4762a1bSJed Brown break; 499c4762a1bSJed Brown default: 500c4762a1bSJed Brown SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 501c4762a1bSJed Brown } 502c4762a1bSJed Brown } else { 503c4762a1bSJed Brown PetscInt numPoints[3] = {0, 0, 0}; 504c4762a1bSJed Brown 505c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 506c4762a1bSJed Brown if (testNum == 2) {ierr = DMCreateLabel(*dm, "pfault");CHKERRQ(ierr);} 507c4762a1bSJed Brown else {ierr = DMCreateLabel(*dm, "fault");CHKERRQ(ierr);} 508c4762a1bSJed Brown } 509c4762a1bSJed Brown ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 510c4762a1bSJed Brown ierr = DMViewFromOptions(idm, NULL, "-in_dm_view");CHKERRQ(ierr); 511c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 512c4762a1bSJed Brown *dm = idm; 513c4762a1bSJed Brown PetscFunctionReturn(0); 514c4762a1bSJed Brown } 515c4762a1bSJed Brown 516c4762a1bSJed Brown PetscErrorCode CreateHex_3D(MPI_Comm comm, PetscInt testNum, DM *dm) 517c4762a1bSJed Brown { 518c4762a1bSJed Brown DM idm; 519c4762a1bSJed Brown PetscInt depth = 3, p; 520c4762a1bSJed Brown PetscMPIInt rank; 521c4762a1bSJed Brown PetscErrorCode ierr; 522c4762a1bSJed Brown 523c4762a1bSJed Brown PetscFunctionBegin; 524c4762a1bSJed Brown ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 525c4762a1bSJed Brown if (!rank) { 526c4762a1bSJed Brown switch (testNum) { 527c4762a1bSJed Brown case 0: 528c4762a1bSJed Brown { 529c4762a1bSJed Brown PetscInt numPoints[2] = {12, 2}; 530c4762a1bSJed Brown PetscInt coneSize[14] = {8,8, 0,0,0,0,0,0,0,0,0,0,0,0}; 531c4762a1bSJed Brown PetscInt cones[16] = {2,5,4,3,6,7,8,9, 3,4,11,10,7,12,13,8}; 532c4762a1bSJed Brown PetscInt coneOrientations[16] = {0,0,0,0,0,0,0,0, 0,0, 0,0 ,0, 0, 0,0}; 533c4762a1bSJed Brown PetscScalar vertexCoords[36] = {-0.5,0.0,0.0, 0.0,0.0,0.0, 0.0,1.0,0.0, -0.5,1.0,0.0, 534c4762a1bSJed Brown -0.5,0.0,1.0, 0.0,0.0,1.0, 0.0,1.0,1.0, -0.5,1.0,1.0, 535c4762a1bSJed Brown 0.5,0.0,0.0, 0.5,1.0,0.0, 0.5,0.0,1.0, 0.5,1.0,1.0}; 536c4762a1bSJed Brown PetscInt markerPoints[52] = {2,1,3,1,4,1,5,1,6,1,7,1,8,1,9,1}; 537c4762a1bSJed Brown PetscInt faultPoints[4] = {3, 4, 7, 8}; 538c4762a1bSJed Brown 539c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 540c4762a1bSJed Brown ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 541c4762a1bSJed Brown for(p = 0; p < 8; ++p) {ierr = DMSetLabelValue(idm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);} 542c4762a1bSJed Brown for(p = 0; p < 4; ++p) {ierr = DMSetLabelValue(idm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 543c4762a1bSJed Brown } 544c4762a1bSJed Brown break; 545c4762a1bSJed Brown case 1: 546c4762a1bSJed Brown { 547c4762a1bSJed Brown /* Cell Adjacency Graph: 548c4762a1bSJed Brown 0 -- { 8, 13, 21, 24} --> 1 549c4762a1bSJed Brown 0 -- {20, 21, 23, 24} --> 5 F 550c4762a1bSJed Brown 1 -- {10, 15, 21, 24} --> 2 551c4762a1bSJed Brown 1 -- {13, 14, 15, 24} --> 6 552c4762a1bSJed Brown 2 -- {21, 22, 24, 25} --> 4 F 553c4762a1bSJed Brown 3 -- {21, 24, 30, 35} --> 4 554c4762a1bSJed Brown 3 -- {21, 24, 28, 33} --> 5 555c4762a1bSJed Brown */ 556c4762a1bSJed Brown PetscInt numPoints[2] = {30, 7}; 557c4762a1bSJed Brown PetscInt coneSize[37] = {8,8,8,8,8,8,8, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 558c4762a1bSJed Brown PetscInt cones[56] = { 8, 21, 20, 7, 13, 12, 23, 24, 559c4762a1bSJed Brown 14, 15, 10, 9, 13, 8, 21, 24, 560c4762a1bSJed Brown 15, 16, 11, 10, 24, 21, 22, 25, 561c4762a1bSJed Brown 30, 29, 28, 21, 35, 24, 33, 34, 562c4762a1bSJed Brown 24, 21, 30, 35, 25, 36, 31, 22, 563c4762a1bSJed Brown 27, 20, 21, 28, 32, 33, 24, 23, 564c4762a1bSJed Brown 15, 24, 13, 14, 19, 18, 17, 26}; 565c4762a1bSJed Brown PetscInt coneOrientations[56] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 566c4762a1bSJed Brown PetscScalar vertexCoords[90] = {-2.0, -2.0, -2.0, -2.0, -1.0, -2.0, -3.0, 0.0, -2.0, -2.0, 1.0, -2.0, -2.0, 2.0, -2.0, -2.0, -2.0, 0.0, 567c4762a1bSJed Brown -2.0, -1.0, 0.0, -3.0, 0.0, 0.0, -2.0, 1.0, 0.0, -2.0, 2.0, 0.0, -2.0, -1.0, 2.0, -3.0, 0.0, 2.0, 568c4762a1bSJed Brown -2.0, 1.0, 2.0, 0.0, -2.0, -2.0, 0.0, 0.0, -2.0, 0.0, 2.0, -2.0, 0.0, -2.0, 0.0, 0.0, 0.0, 0.0, 569c4762a1bSJed Brown 0.0, 2.0, 0.0, 0.0, 0.0, 2.0, 2.0, -2.0, -2.0, 2.0, -1.0, -2.0, 3.0, 0.0, -2.0, 2.0, 1.0, -2.0, 570c4762a1bSJed Brown 2.0, 2.0, -2.0, 2.0, -2.0, 0.0, 2.0, -1.0, 0.0, 3.0, 0.0, 0.0, 2.0, 1.0, 0.0, 2.0, 2.0, 0.0}; 571c4762a1bSJed Brown PetscInt faultPoints[6] = {20, 21, 22, 23, 24, 25}; 572c4762a1bSJed Brown 573c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 574c4762a1bSJed Brown ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 575c4762a1bSJed Brown for(p = 0; p < 6; ++p) {ierr = DMSetLabelValue(idm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 576c4762a1bSJed Brown } 577c4762a1bSJed Brown break; 578c4762a1bSJed Brown case 2: 579c4762a1bSJed Brown { 580c4762a1bSJed Brown /* Buried fault edge */ 581c4762a1bSJed Brown PetscInt numPoints[2] = {18, 4}; 582c4762a1bSJed Brown PetscInt coneSize[22] = {8,8,8,8, 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 583c4762a1bSJed Brown PetscInt cones[32] = { 4, 5, 8, 7, 13, 16, 17, 14, 584c4762a1bSJed Brown 5, 6, 9, 8, 14, 17, 18, 15, 585c4762a1bSJed Brown 7, 8, 11, 10, 16, 19, 20, 17, 586c4762a1bSJed Brown 8, 9, 12, 11, 17, 20, 21, 18}; 587c4762a1bSJed Brown PetscInt coneOrientations[32] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; 588c4762a1bSJed Brown PetscScalar vertexCoords[54] = {-2.0, -2.0, 0.0, -2.0, 0.0, 0.0, -2.0, 2.0, 0.0, 0.0, -2.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 589c4762a1bSJed Brown 2.0, -2.0, 0.0, 2.0, 0.0, 0.0, 2.0, 2.0, 0.0, -2.0, -2.0, 2.0, -2.0, 0.0, 2.0, -2.0, 2.0, 2.0, 590c4762a1bSJed Brown 0.0, -2.0, 2.0, 0.0, 0.0, 2.0, 0.0, 2.0, 2.0, 2.0, -2.0, 2.0, 2.0, 0.0, 2.0, 2.0, 2.0, 2.0}; 591c4762a1bSJed Brown PetscInt faultPoints[4] = {7, 8, 16, 17}; 592c4762a1bSJed Brown 593c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr); 594c4762a1bSJed Brown ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 595c4762a1bSJed Brown for(p = 0; p < 4; ++p) {ierr = DMSetLabelValue(idm, "fault", faultPoints[p], 1);CHKERRQ(ierr);} 596c4762a1bSJed Brown } 597c4762a1bSJed Brown break; 598c4762a1bSJed Brown default: SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum); 599c4762a1bSJed Brown } 600c4762a1bSJed Brown } else { 601c4762a1bSJed Brown PetscInt numPoints[4] = {0, 0, 0, 0}; 602c4762a1bSJed Brown 603c4762a1bSJed Brown ierr = DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 604c4762a1bSJed Brown ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr); 605c4762a1bSJed Brown ierr = DMCreateLabel(idm, "fault");CHKERRQ(ierr); 606c4762a1bSJed Brown } 607c4762a1bSJed Brown ierr = DMViewFromOptions(idm, NULL, "-in_dm_view");CHKERRQ(ierr); 608c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 609c4762a1bSJed Brown *dm = idm; 610c4762a1bSJed Brown PetscFunctionReturn(0); 611c4762a1bSJed Brown } 612c4762a1bSJed Brown 613c4762a1bSJed Brown PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 614c4762a1bSJed Brown { 615c4762a1bSJed Brown PetscInt dim = user->dim; 616c4762a1bSJed Brown PetscBool cellSimplex = user->cellSimplex, hasFault, hasFault2, hasParallelFault; 617c4762a1bSJed Brown PetscMPIInt rank, size; 618c4762a1bSJed Brown PetscErrorCode ierr; 619c4762a1bSJed Brown 620c4762a1bSJed Brown PetscFunctionBegin; 621c4762a1bSJed Brown ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 622c4762a1bSJed Brown ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr); 623c4762a1bSJed Brown ierr = DMCreate(comm, dm);CHKERRQ(ierr); 624c4762a1bSJed Brown ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 625c4762a1bSJed Brown ierr = DMSetDimension(*dm, dim);CHKERRQ(ierr); 626c4762a1bSJed Brown switch (dim) { 627c4762a1bSJed Brown case 2: 628c4762a1bSJed Brown if (cellSimplex) { 629c4762a1bSJed Brown ierr = CreateSimplex_2D(comm, user->testNum, dm);CHKERRQ(ierr); 630c4762a1bSJed Brown } else { 631c4762a1bSJed Brown ierr = CreateQuad_2D(comm, user->testNum, dm);CHKERRQ(ierr); 632c4762a1bSJed Brown } 633c4762a1bSJed Brown break; 634c4762a1bSJed Brown case 3: 635c4762a1bSJed Brown if (cellSimplex) { 636c4762a1bSJed Brown ierr = CreateSimplex_3D(comm, user, *dm);CHKERRQ(ierr); 637c4762a1bSJed Brown } else { 638c4762a1bSJed Brown ierr = CreateHex_3D(comm, user->testNum, dm);CHKERRQ(ierr); 639c4762a1bSJed Brown } 640c4762a1bSJed Brown break; 641c4762a1bSJed Brown default: 642c4762a1bSJed Brown SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make hybrid meshes for dimension %d", dim); 643c4762a1bSJed Brown } 644c4762a1bSJed Brown ierr = PetscObjectSetOptionsPrefix((PetscObject) *dm, "orig_");CHKERRQ(ierr); 645c4762a1bSJed Brown ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 646c4762a1bSJed Brown ierr = DMHasLabel(*dm, "fault", &hasFault);CHKERRQ(ierr); 647c4762a1bSJed Brown if (hasFault) { 648c4762a1bSJed Brown DM dmHybrid = NULL; 649c4762a1bSJed Brown DMLabel faultLabel, faultBdLabel, hybridLabel; 650c4762a1bSJed Brown 651c4762a1bSJed Brown ierr = DMGetLabel(*dm, "fault", &faultLabel);CHKERRQ(ierr); 652c4762a1bSJed Brown ierr = DMGetLabel(*dm, "faultBd", &faultBdLabel);CHKERRQ(ierr); 653c4762a1bSJed Brown ierr = DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, &hybridLabel, NULL, NULL, &dmHybrid);CHKERRQ(ierr); 654c4762a1bSJed Brown ierr = DMLabelView(hybridLabel, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 655c4762a1bSJed Brown ierr = DMLabelDestroy(&hybridLabel);CHKERRQ(ierr); 656c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 657c4762a1bSJed Brown *dm = dmHybrid; 658c4762a1bSJed Brown } 659c4762a1bSJed Brown ierr = DMHasLabel(*dm, "fault2", &hasFault2);CHKERRQ(ierr); 660c4762a1bSJed Brown if (hasFault2) { 661c4762a1bSJed Brown DM dmHybrid = NULL; 662c4762a1bSJed Brown DMLabel faultLabel, faultBdLabel, hybridLabel; 663c4762a1bSJed Brown 664c4762a1bSJed Brown ierr = PetscObjectSetOptionsPrefix((PetscObject) *dm, "faulted_");CHKERRQ(ierr); 665c4762a1bSJed Brown ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 666c4762a1bSJed Brown ierr = DMViewFromOptions(*dm, NULL, "-faulted_dm_view");CHKERRQ(ierr); 667c4762a1bSJed Brown ierr = DMGetLabel(*dm, "fault2", &faultLabel);CHKERRQ(ierr); 668c4762a1bSJed Brown ierr = DMGetLabel(*dm, "fault2Bd", &faultBdLabel);CHKERRQ(ierr); 669c4762a1bSJed Brown ierr = DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, &hybridLabel, NULL, NULL, &dmHybrid);CHKERRQ(ierr); 670c4762a1bSJed Brown ierr = DMLabelView(hybridLabel, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 671c4762a1bSJed Brown ierr = DMLabelDestroy(&hybridLabel);CHKERRQ(ierr); 672c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 673c4762a1bSJed Brown *dm = dmHybrid; 674c4762a1bSJed Brown } 675c4762a1bSJed Brown if (user->testPartition && size > 1) { 676c4762a1bSJed Brown PetscPartitioner part; 677c4762a1bSJed Brown PetscInt *sizes = NULL; 678c4762a1bSJed Brown PetscInt *points = NULL; 679c4762a1bSJed Brown 680c4762a1bSJed Brown if (!rank) { 681c4762a1bSJed Brown if (dim == 2 && cellSimplex && size == 2) { 682c4762a1bSJed Brown switch (user->testNum) { 683c4762a1bSJed Brown case 0: { 684c4762a1bSJed Brown PetscInt triSizes_p2[2] = {1, 2}; 685c4762a1bSJed Brown PetscInt triPoints_p2[3] = {0, 1, 2}; 686c4762a1bSJed Brown 687c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 3, &points);CHKERRQ(ierr); 688c4762a1bSJed Brown ierr = PetscArraycpy(sizes, triSizes_p2, 2);CHKERRQ(ierr); 689c4762a1bSJed Brown ierr = PetscArraycpy(points, triPoints_p2, 3);CHKERRQ(ierr);break;} 690c4762a1bSJed Brown default: 691c4762a1bSJed Brown SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular mesh on 2 procs", user->testNum); 692c4762a1bSJed Brown } 693c4762a1bSJed Brown } else if (dim == 2 && !cellSimplex && size == 2) { 694c4762a1bSJed Brown switch (user->testNum) { 695c4762a1bSJed Brown case 0: { 696c4762a1bSJed Brown PetscInt quadSizes_p2[2] = {1, 2}; 697c4762a1bSJed Brown PetscInt quadPoints_p2[3] = {0, 1, 2}; 698c4762a1bSJed Brown 699c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 3, &points);CHKERRQ(ierr); 700c4762a1bSJed Brown ierr = PetscArraycpy(sizes, quadSizes_p2, 2);CHKERRQ(ierr); 701c4762a1bSJed Brown ierr = PetscArraycpy(points, quadPoints_p2, 3);CHKERRQ(ierr);break;} 702c4762a1bSJed Brown case 2: { 703c4762a1bSJed Brown PetscInt quadSizes_p2[2] = {1, 1}; 704c4762a1bSJed Brown PetscInt quadPoints_p2[2] = {0, 1}; 705c4762a1bSJed Brown 706c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 2, &points);CHKERRQ(ierr); 707c4762a1bSJed Brown ierr = PetscArraycpy(sizes, quadSizes_p2, 2);CHKERRQ(ierr); 708c4762a1bSJed Brown ierr = PetscArraycpy(points, quadPoints_p2, 2);CHKERRQ(ierr);break;} 709c4762a1bSJed Brown default: 710c4762a1bSJed Brown SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for quadrilateral mesh on 2 procs", user->testNum); 711c4762a1bSJed Brown } 712c4762a1bSJed Brown } else if (dim == 3 && cellSimplex && size == 2) { 713c4762a1bSJed Brown switch (user->testNum) { 714c4762a1bSJed Brown case 0: { 715c4762a1bSJed Brown PetscInt tetSizes_p2[2] = {1, 2}; 716c4762a1bSJed Brown PetscInt tetPoints_p2[3] = {0, 1, 2}; 717c4762a1bSJed Brown 718c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 3, &points);CHKERRQ(ierr); 719c4762a1bSJed Brown ierr = PetscArraycpy(sizes, tetSizes_p2, 2);CHKERRQ(ierr); 720c4762a1bSJed Brown ierr = PetscArraycpy(points, tetPoints_p2, 3);CHKERRQ(ierr);break;} 721c4762a1bSJed Brown default: 722c4762a1bSJed Brown SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular mesh on 2 procs", user->testNum); 723c4762a1bSJed Brown } 724c4762a1bSJed Brown } else if (dim == 3 && !cellSimplex && size == 2) { 725c4762a1bSJed Brown switch (user->testNum) { 726c4762a1bSJed Brown case 0: { 727c4762a1bSJed Brown PetscInt hexSizes_p2[2] = {1, 2}; 728c4762a1bSJed Brown PetscInt hexPoints_p2[3] = {0, 1, 2}; 729c4762a1bSJed Brown 730c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 3, &points);CHKERRQ(ierr); 731c4762a1bSJed Brown ierr = PetscArraycpy(sizes, hexSizes_p2, 2);CHKERRQ(ierr); 732c4762a1bSJed Brown ierr = PetscArraycpy(points, hexPoints_p2, 3);CHKERRQ(ierr);break;} 733c4762a1bSJed Brown default: 734c4762a1bSJed Brown SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular mesh on 2 procs", user->testNum); 735c4762a1bSJed Brown } 736c4762a1bSJed Brown } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition"); 737c4762a1bSJed Brown } 738c4762a1bSJed Brown ierr = DMPlexGetPartitioner(*dm, &part);CHKERRQ(ierr); 739c4762a1bSJed Brown ierr = PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);CHKERRQ(ierr); 740c4762a1bSJed Brown ierr = PetscPartitionerShellSetPartition(part, size, sizes, points);CHKERRQ(ierr); 741c4762a1bSJed Brown ierr = PetscFree2(sizes, points);CHKERRQ(ierr); 742c4762a1bSJed Brown } 743c4762a1bSJed Brown { 744c4762a1bSJed Brown DM pdm = NULL; 745c4762a1bSJed Brown 746c4762a1bSJed Brown /* Distribute mesh over processes */ 747c4762a1bSJed Brown ierr = DMPlexDistribute(*dm, 0, NULL, &pdm);CHKERRQ(ierr); 748c4762a1bSJed Brown if (pdm) { 749c4762a1bSJed Brown ierr = DMViewFromOptions(pdm, NULL, "-dm_view");CHKERRQ(ierr); 750c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 751c4762a1bSJed Brown *dm = pdm; 752c4762a1bSJed Brown } 753c4762a1bSJed Brown } 754c4762a1bSJed Brown ierr = DMHasLabel(*dm, "pfault", &hasParallelFault);CHKERRQ(ierr); 755c4762a1bSJed Brown if (hasParallelFault) { 756c4762a1bSJed Brown DM dmHybrid = NULL; 757c4762a1bSJed Brown DMLabel faultLabel, faultBdLabel, hybridLabel; 758c4762a1bSJed Brown 759c4762a1bSJed Brown ierr = DMGetLabel(*dm, "pfault", &faultLabel);CHKERRQ(ierr); 760c4762a1bSJed Brown ierr = DMGetLabel(*dm, "pfaultBd", &faultBdLabel);CHKERRQ(ierr); 761c4762a1bSJed Brown ierr = DMPlexCreateHybridMesh(*dm, faultLabel, faultBdLabel, &hybridLabel, NULL, NULL, &dmHybrid);CHKERRQ(ierr); 762c4762a1bSJed Brown ierr = DMLabelView(hybridLabel, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 763c4762a1bSJed Brown ierr = DMLabelDestroy(&hybridLabel);CHKERRQ(ierr); 764c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 765c4762a1bSJed Brown *dm = dmHybrid; 766c4762a1bSJed Brown } 767c4762a1bSJed Brown ierr = PetscObjectSetName((PetscObject) *dm, "Hybrid Mesh");CHKERRQ(ierr); 768c4762a1bSJed Brown ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 769c4762a1bSJed Brown ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 770c4762a1bSJed Brown PetscFunctionReturn(0); 771c4762a1bSJed Brown } 772c4762a1bSJed Brown 773c4762a1bSJed Brown int main(int argc, char **argv) 774c4762a1bSJed Brown { 775c4762a1bSJed Brown DM dm; 776c4762a1bSJed Brown AppCtx user; /* user-defined work context */ 777c4762a1bSJed Brown PetscErrorCode ierr; 778c4762a1bSJed Brown 779c4762a1bSJed Brown ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 780c4762a1bSJed Brown ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 781c4762a1bSJed Brown ierr = CreateMesh(PETSC_COMM_WORLD, &user, &dm);CHKERRQ(ierr); 782c4762a1bSJed Brown ierr = DMDestroy(&dm);CHKERRQ(ierr); 783c4762a1bSJed Brown ierr = PetscFinalize(); 784c4762a1bSJed Brown return ierr; 785c4762a1bSJed Brown } 786c4762a1bSJed Brown 787c4762a1bSJed Brown /*TEST 788c4762a1bSJed Brown testset: 789*54fcfd0cSMatthew G. Knepley args: -orig_dm_plex_check_all -dm_plex_check_all 790c4762a1bSJed Brown # 2D Simplex 791c4762a1bSJed Brown test: 792c4762a1bSJed Brown suffix: tri_0 793c4762a1bSJed Brown args: -dim 2 794c4762a1bSJed Brown test: 795c4762a1bSJed Brown suffix: tri_1 796c4762a1bSJed Brown nsize: 2 797c4762a1bSJed Brown args: -dim 2 798c4762a1bSJed Brown test: 799c4762a1bSJed Brown suffix: tri_t1_0 800c4762a1bSJed Brown args: -dim 2 -test_num 1 801c4762a1bSJed Brown # 3D Simplex 802c4762a1bSJed Brown test: 803c4762a1bSJed Brown suffix: tet_0 804c4762a1bSJed Brown args: -dim 3 805c4762a1bSJed Brown test: 806c4762a1bSJed Brown suffix: tet_1 807c4762a1bSJed Brown nsize: 2 808c4762a1bSJed Brown args: -dim 3 809c4762a1bSJed Brown test: 810c4762a1bSJed Brown suffix: tet_t1_0 811c4762a1bSJed Brown args: -dim 3 -test_num 1 812c4762a1bSJed Brown 813c4762a1bSJed Brown testset: 814*54fcfd0cSMatthew G. Knepley args: -orig_dm_plex_check_all -dm_plex_check_all 815c4762a1bSJed Brown # 2D Quads 816c4762a1bSJed Brown test: 817c4762a1bSJed Brown suffix: quad_0 818c4762a1bSJed Brown args: -dim 2 -cell_simplex 0 819c4762a1bSJed Brown test: 820c4762a1bSJed Brown suffix: quad_1 821c4762a1bSJed Brown nsize: 2 822c4762a1bSJed Brown args: -dim 2 -cell_simplex 0 823c4762a1bSJed Brown test: 824c4762a1bSJed Brown suffix: quad_t1_0 825*54fcfd0cSMatthew G. Knepley args: -dim 2 -cell_simplex 0 -test_num 1 -faulted_dm_plex_check_all 826c4762a1bSJed Brown # 3D Hex 827c4762a1bSJed Brown test: 828c4762a1bSJed Brown suffix: hex_0 829c4762a1bSJed Brown args: -dim 3 -cell_simplex 0 830c4762a1bSJed Brown test: 831c4762a1bSJed Brown suffix: hex_1 832c4762a1bSJed Brown nsize: 2 833c4762a1bSJed Brown args: -dim 3 -cell_simplex 0 834c4762a1bSJed Brown test: 835c4762a1bSJed Brown suffix: hex_t1_0 836c4762a1bSJed Brown args: -dim 3 -cell_simplex 0 -test_num 1 837c4762a1bSJed Brown test: 838c4762a1bSJed Brown suffix: hex_t2_0 839c4762a1bSJed Brown args: -dim 3 -cell_simplex 0 -test_num 2 840c4762a1bSJed Brown 841c4762a1bSJed Brown TEST*/ 842