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