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