1 static char help[] = "Create a Plex sphere from quads and create a P1 section\n\n"; 2 3 #include <petscdmplex.h> 4 5 typedef struct { 6 PetscInt dim; /* Topological problem dimension */ 7 PetscBool simplex; /* Mesh with simplices */ 8 } AppCtx; 9 10 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 11 { 12 PetscErrorCode ierr; 13 14 PetscFunctionBeginUser; 15 options->dim = 2; 16 options->simplex = PETSC_FALSE; 17 18 ierr = PetscOptionsBegin(comm, "", "Sphere Mesh Options", "DMPLEX");CHKERRQ(ierr); 19 ierr = PetscOptionsRangeInt("-dim", "Problem dimension", "ex7.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr); 20 ierr = PetscOptionsBool("-simplex", "Use simplices, or tensor product cells", "ex7.c", options->simplex, &options->simplex, NULL);CHKERRQ(ierr); 21 ierr = PetscOptionsEnd(); 22 PetscFunctionReturn(0); 23 } 24 25 static PetscErrorCode SetupSection(DM dm) 26 { 27 PetscSection s; 28 PetscInt vStart, vEnd, v; 29 PetscErrorCode ierr; 30 31 PetscFunctionBeginUser; 32 ierr = DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);CHKERRQ(ierr); 33 ierr = PetscSectionCreate(PetscObjectComm((PetscObject) dm), &s);CHKERRQ(ierr); 34 ierr = PetscSectionSetNumFields(s, 1);CHKERRQ(ierr); 35 ierr = PetscSectionSetFieldComponents(s, 0, 1);CHKERRQ(ierr); 36 ierr = PetscSectionSetChart(s, vStart, vEnd);CHKERRQ(ierr); 37 for (v = vStart; v < vEnd; ++v) { 38 ierr = PetscSectionSetDof(s, v, 1);CHKERRQ(ierr); 39 ierr = PetscSectionSetFieldDof(s, v, 0, 1);CHKERRQ(ierr); 40 } 41 ierr = PetscSectionSetUp(s);CHKERRQ(ierr); 42 ierr = DMSetLocalSection(dm, s);CHKERRQ(ierr); 43 ierr = PetscSectionDestroy(&s);CHKERRQ(ierr); 44 PetscFunctionReturn(0); 45 } 46 47 int main(int argc, char **argv) 48 { 49 DM dm; 50 Vec u; 51 AppCtx ctx; 52 PetscErrorCode ierr; 53 54 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 55 ierr = ProcessOptions(PETSC_COMM_WORLD, &ctx);CHKERRQ(ierr); 56 ierr = DMPlexCreateSphereMesh(PETSC_COMM_WORLD, ctx.dim, ctx.simplex, 1.0, &dm);CHKERRQ(ierr); 57 ierr = PetscObjectSetName((PetscObject) dm, "Sphere");CHKERRQ(ierr); 58 /* Distribute mesh over processes */ 59 { 60 DM dmDist = NULL; 61 PetscPartitioner part; 62 ierr = DMPlexGetPartitioner(dm, &part);CHKERRQ(ierr); 63 ierr = PetscPartitionerSetFromOptions(part);CHKERRQ(ierr); 64 ierr = DMPlexDistribute(dm, 0, NULL, &dmDist);CHKERRQ(ierr); 65 if (dmDist) { 66 ierr = DMDestroy(&dm);CHKERRQ(ierr); 67 dm = dmDist; 68 } 69 } 70 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 71 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 72 ierr = SetupSection(dm);CHKERRQ(ierr); 73 ierr = DMGetGlobalVector(dm, &u);CHKERRQ(ierr); 74 ierr = VecSet(u, 2);CHKERRQ(ierr); 75 ierr = VecViewFromOptions(u, NULL, "-vec_view");CHKERRQ(ierr); 76 ierr = DMRestoreGlobalVector(dm, &u);CHKERRQ(ierr); 77 ierr = DMDestroy(&dm);CHKERRQ(ierr); 78 ierr = PetscFinalize(); 79 return ierr; 80 } 81 82 /*TEST 83 84 test: 85 suffix: 2d_quad 86 requires: !__float128 87 args: -dm_view 88 89 test: 90 suffix: 2d_quad_parallel 91 requires: !__float128 92 args: -dm_view -petscpartitioner_type simple 93 nsize: 2 94 95 test: 96 suffix: 2d_tri 97 requires: !__float128 98 args: -simplex -dm_view 99 100 test: 101 suffix: 2d_tri_parallel 102 requires: !__float128 103 args: -simplex -dm_view -petscpartitioner_type simple 104 nsize: 2 105 106 test: 107 suffix: 3d_tri 108 requires: !__float128 109 args: -dim 3 -simplex -dm_view 110 111 test: 112 suffix: 3d_tri_parallel 113 requires: !__float128 114 args: -dim 3 -simplex -dm_view -petscpartitioner_type simple 115 nsize: 2 116 117 TEST*/ 118