1c4762a1bSJed Brown static char help[] = "Orient a mesh in parallel\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscdmplex.h> 4c4762a1bSJed Brown 5c4762a1bSJed Brown typedef struct { 6c4762a1bSJed Brown /* Domain and mesh definition */ 7c4762a1bSJed Brown PetscInt dim; /* The topological mesh dimension */ 8c4762a1bSJed Brown PetscBool cellSimplex; /* Use simplices or hexes */ 9c4762a1bSJed Brown char filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */ 10c4762a1bSJed Brown PetscBool testPartition; /* Use a fixed partitioning for testing */ 11c4762a1bSJed Brown PetscInt testNum; /* Labels the different test partitions */ 12c4762a1bSJed Brown } AppCtx; 13c4762a1bSJed Brown 14c4762a1bSJed Brown static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 15c4762a1bSJed Brown { 16c4762a1bSJed Brown PetscErrorCode ierr; 17c4762a1bSJed Brown 18c4762a1bSJed Brown PetscFunctionBeginUser; 19c4762a1bSJed Brown options->dim = 2; 20c4762a1bSJed Brown options->cellSimplex = PETSC_TRUE; 21c4762a1bSJed Brown options->filename[0] = '\0'; 22c4762a1bSJed Brown options->testPartition = PETSC_TRUE; 23c4762a1bSJed Brown options->testNum = 0; 24c4762a1bSJed Brown 25c4762a1bSJed Brown ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr); 26c4762a1bSJed Brown ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex13.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr); 27c4762a1bSJed Brown ierr = PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex13.c", options->cellSimplex, &options->cellSimplex, NULL);CHKERRQ(ierr); 28*589a23caSBarry Smith ierr = PetscOptionsString("-filename", "The mesh file", "ex13.c", options->filename, options->filename, sizeof(options->filename), NULL);CHKERRQ(ierr); 29c4762a1bSJed Brown ierr = PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex13.c", options->testPartition, &options->testPartition, NULL);CHKERRQ(ierr); 30c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-test_num", "The test partition number", "ex13.c", options->testNum, &options->testNum, NULL,0);CHKERRQ(ierr); 31c4762a1bSJed Brown ierr = PetscOptionsEnd(); 32c4762a1bSJed Brown PetscFunctionReturn(0); 33c4762a1bSJed Brown } 34c4762a1bSJed Brown 35c4762a1bSJed Brown static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 36c4762a1bSJed Brown { 37c4762a1bSJed Brown DM dmDist = NULL; 38c4762a1bSJed Brown PetscInt dim = user->dim; 39c4762a1bSJed Brown PetscBool cellSimplex = user->cellSimplex; 40c4762a1bSJed Brown const char *filename = user->filename; 41c4762a1bSJed Brown size_t len; 42c4762a1bSJed Brown PetscErrorCode ierr; 43c4762a1bSJed Brown 44c4762a1bSJed Brown PetscFunctionBeginUser; 45c4762a1bSJed Brown ierr = PetscStrlen(filename, &len);CHKERRQ(ierr); 46c4762a1bSJed Brown if (len) {ierr = DMPlexCreateFromFile(comm, filename, PETSC_TRUE, dm);CHKERRQ(ierr);} 47c4762a1bSJed Brown else {ierr = DMPlexCreateBoxMesh(comm, dim, cellSimplex, NULL, NULL, NULL, NULL, PETSC_TRUE, dm);CHKERRQ(ierr);} 48c4762a1bSJed Brown if (user->testPartition) { 49c4762a1bSJed Brown PetscPartitioner part; 50c4762a1bSJed Brown PetscInt *sizes = NULL; 51c4762a1bSJed Brown PetscInt *points = NULL; 52c4762a1bSJed Brown PetscMPIInt rank, size; 53c4762a1bSJed Brown 54c4762a1bSJed Brown ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 55c4762a1bSJed Brown ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr); 56c4762a1bSJed Brown if (!rank) { 57c4762a1bSJed Brown if (dim == 2 && cellSimplex && size == 2) { 58c4762a1bSJed Brown switch (user->testNum) { 59c4762a1bSJed Brown case 0: { 60c4762a1bSJed Brown PetscInt triSizes_p2[2] = {4, 4}; 61c4762a1bSJed Brown PetscInt triPoints_p2[8] = {3, 5, 6, 7, 0, 1, 2, 4}; 62c4762a1bSJed Brown 63c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 8, &points);CHKERRQ(ierr); 64c4762a1bSJed Brown ierr = PetscArraycpy(sizes, triSizes_p2, 2);CHKERRQ(ierr); 65c4762a1bSJed Brown ierr = PetscArraycpy(points, triPoints_p2, 8);CHKERRQ(ierr); 66c4762a1bSJed Brown break;} 67c4762a1bSJed Brown case 1: { 68c4762a1bSJed Brown PetscInt triSizes_p2[2] = {6, 2}; 69c4762a1bSJed Brown PetscInt triPoints_p2[8] = {1, 2, 3, 4, 6, 7, 0, 5}; 70c4762a1bSJed Brown 71c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 8, &points);CHKERRQ(ierr); 72c4762a1bSJed Brown ierr = PetscArraycpy(sizes, triSizes_p2, 2);CHKERRQ(ierr); 73c4762a1bSJed Brown ierr = PetscArraycpy(points, triPoints_p2, 8);CHKERRQ(ierr); 74c4762a1bSJed Brown break;} 75c4762a1bSJed Brown default: 76c4762a1bSJed Brown SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %d for triangular mesh on 2 procs", user->testNum); 77c4762a1bSJed Brown } 78c4762a1bSJed Brown } else if (dim == 2 && cellSimplex && size == 3) { 79c4762a1bSJed Brown PetscInt triSizes_p3[3] = {3, 3, 2}; 80c4762a1bSJed Brown PetscInt triPoints_p3[8] = {1, 2, 4, 3, 6, 7, 0, 5}; 81c4762a1bSJed Brown 82c4762a1bSJed Brown ierr = PetscMalloc2(3, &sizes, 8, &points);CHKERRQ(ierr); 83c4762a1bSJed Brown ierr = PetscArraycpy(sizes, triSizes_p3, 3);CHKERRQ(ierr); 84c4762a1bSJed Brown ierr = PetscArraycpy(points, triPoints_p3, 8);CHKERRQ(ierr); 85c4762a1bSJed Brown } else if (dim == 2 && !cellSimplex && size == 2) { 86c4762a1bSJed Brown PetscInt quadSizes_p2[2] = {2, 2}; 87c4762a1bSJed Brown PetscInt quadPoints_p2[4] = {2, 3, 0, 1}; 88c4762a1bSJed Brown 89c4762a1bSJed Brown ierr = PetscMalloc2(2, &sizes, 4, &points);CHKERRQ(ierr); 90c4762a1bSJed Brown ierr = PetscArraycpy(sizes, quadSizes_p2, 2);CHKERRQ(ierr); 91c4762a1bSJed Brown ierr = PetscArraycpy(points, quadPoints_p2, 4);CHKERRQ(ierr); 92c4762a1bSJed Brown } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition"); 93c4762a1bSJed Brown } 94c4762a1bSJed Brown ierr = DMPlexGetPartitioner(*dm, &part);CHKERRQ(ierr); 95c4762a1bSJed Brown ierr = PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);CHKERRQ(ierr); 96c4762a1bSJed Brown ierr = PetscPartitionerShellSetPartition(part, size, sizes, points);CHKERRQ(ierr); 97c4762a1bSJed Brown ierr = PetscFree2(sizes, points);CHKERRQ(ierr); 98c4762a1bSJed Brown } 99c4762a1bSJed Brown ierr = DMPlexDistribute(*dm, 0, NULL, &dmDist);CHKERRQ(ierr); 100c4762a1bSJed Brown if (dmDist) { 101c4762a1bSJed Brown ierr = DMDestroy(dm);CHKERRQ(ierr); 102c4762a1bSJed Brown *dm = dmDist; 103c4762a1bSJed Brown } 104c4762a1bSJed Brown ierr = PetscObjectSetName((PetscObject) *dm, cellSimplex ? "Simplicial Mesh" : "Tensor Product Mesh");CHKERRQ(ierr); 105c4762a1bSJed Brown ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 106c4762a1bSJed Brown PetscFunctionReturn(0); 107c4762a1bSJed Brown } 108c4762a1bSJed Brown 109c4762a1bSJed Brown static PetscErrorCode ScrambleOrientation(DM dm, AppCtx *user) 110c4762a1bSJed Brown { 111c4762a1bSJed Brown PetscInt h, cStart, cEnd, c; 112c4762a1bSJed Brown PetscErrorCode ierr; 113c4762a1bSJed Brown 114c4762a1bSJed Brown PetscFunctionBeginUser; 115c4762a1bSJed Brown ierr = DMPlexGetVTKCellHeight(dm, &h);CHKERRQ(ierr); 116c4762a1bSJed Brown ierr = DMPlexGetHeightStratum(dm, h, &cStart, &cEnd);CHKERRQ(ierr); 117c4762a1bSJed Brown for (c = cStart; c < cEnd; ++c) { 118c4762a1bSJed Brown /* Could use PetscRand instead */ 119c4762a1bSJed Brown if (c%2) {ierr = DMPlexReverseCell(dm, c);CHKERRQ(ierr);} 120c4762a1bSJed Brown } 121c4762a1bSJed Brown PetscFunctionReturn(0); 122c4762a1bSJed Brown } 123c4762a1bSJed Brown 124c4762a1bSJed Brown static PetscErrorCode TestOrientation(DM dm, AppCtx *user) 125c4762a1bSJed Brown { 126c4762a1bSJed Brown PetscErrorCode ierr; 127c4762a1bSJed Brown 128c4762a1bSJed Brown PetscFunctionBeginUser; 129c4762a1bSJed Brown ierr = ScrambleOrientation(dm, user);CHKERRQ(ierr); 130c4762a1bSJed Brown ierr = DMPlexOrient(dm);CHKERRQ(ierr); 131c4762a1bSJed Brown ierr = DMViewFromOptions(dm, NULL, "-oriented_dm_view");CHKERRQ(ierr); 132c4762a1bSJed Brown PetscFunctionReturn(0); 133c4762a1bSJed Brown } 134c4762a1bSJed Brown 135c4762a1bSJed Brown int main(int argc, char **argv) 136c4762a1bSJed Brown { 137c4762a1bSJed Brown DM dm; 138c4762a1bSJed Brown AppCtx user; /* user-defined work context */ 139c4762a1bSJed Brown PetscErrorCode ierr; 140c4762a1bSJed Brown 141c4762a1bSJed Brown ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 142c4762a1bSJed Brown ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 143c4762a1bSJed Brown ierr = CreateMesh(PETSC_COMM_WORLD, &user, &dm);CHKERRQ(ierr); 144c4762a1bSJed Brown ierr = TestOrientation(dm, &user);CHKERRQ(ierr); 145c4762a1bSJed Brown ierr = DMDestroy(&dm);CHKERRQ(ierr); 146c4762a1bSJed Brown ierr = PetscFinalize(); 147c4762a1bSJed Brown return ierr; 148c4762a1bSJed Brown } 149c4762a1bSJed Brown 150c4762a1bSJed Brown /*TEST 151c4762a1bSJed Brown test: 152c4762a1bSJed Brown suffix: 0 153c4762a1bSJed Brown requires: triangle 154c4762a1bSJed Brown args: -test_partition 0 -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view 155c4762a1bSJed Brown test: 156c4762a1bSJed Brown suffix: 1 157c4762a1bSJed Brown requires: triangle 158c4762a1bSJed Brown nsize: 2 159c4762a1bSJed Brown args: -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view 160c4762a1bSJed Brown test: 161c4762a1bSJed Brown suffix: 2 162c4762a1bSJed Brown requires: triangle 163c4762a1bSJed Brown nsize: 2 164c4762a1bSJed Brown args: -test_num 1 -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view 165c4762a1bSJed Brown test: 166c4762a1bSJed Brown suffix: 3 167c4762a1bSJed Brown requires: triangle 168c4762a1bSJed Brown nsize: 3 169c4762a1bSJed Brown args: -dm_view ascii::ascii_info_detail -oriented_dm_view ascii::ascii_info_detail -orientation_view -orientation_view_synchronized 170c4762a1bSJed Brown 171c4762a1bSJed Brown TEST*/ 172