1 static char help[] = "Create a Plex sphere from quads and create a P1 section\n\n";
2
3 #include <petscdmplex.h>
4
SetupSection(DM dm)5 static PetscErrorCode SetupSection(DM dm)
6 {
7 PetscSection s;
8 PetscInt vStart, vEnd, v;
9
10 PetscFunctionBeginUser;
11 PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd));
12 PetscCall(PetscSectionCreate(PetscObjectComm((PetscObject)dm), &s));
13 PetscCall(PetscSectionSetNumFields(s, 1));
14 PetscCall(PetscSectionSetFieldComponents(s, 0, 1));
15 PetscCall(PetscSectionSetChart(s, vStart, vEnd));
16 for (v = vStart; v < vEnd; ++v) {
17 PetscCall(PetscSectionSetDof(s, v, 1));
18 PetscCall(PetscSectionSetFieldDof(s, v, 0, 1));
19 }
20 PetscCall(PetscSectionSetUp(s));
21 PetscCall(DMSetLocalSection(dm, s));
22 PetscCall(PetscSectionDestroy(&s));
23 PetscFunctionReturn(PETSC_SUCCESS);
24 }
25
main(int argc,char ** argv)26 int main(int argc, char **argv)
27 {
28 DM dm;
29 Vec u;
30
31 PetscFunctionBeginUser;
32 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
33 PetscCall(DMCreate(PETSC_COMM_WORLD, &dm));
34 PetscCall(DMSetType(dm, DMPLEX));
35 PetscCall(DMSetFromOptions(dm));
36 PetscCall(PetscObjectSetName((PetscObject)dm, "Sphere"));
37 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
38
39 PetscCall(SetupSection(dm));
40 PetscCall(DMGetGlobalVector(dm, &u));
41 PetscCall(VecSet(u, 2));
42 PetscCall(VecViewFromOptions(u, NULL, "-vec_view"));
43 PetscCall(DMRestoreGlobalVector(dm, &u));
44 PetscCall(DMDestroy(&dm));
45 PetscCall(PetscFinalize());
46 return 0;
47 }
48
49 /*TEST
50
51 testset:
52 requires: !__float128
53 args: -dm_plex_shape sphere -dm_view
54
55 test:
56 suffix: 2d_quad
57 args: -dm_plex_simplex 0
58
59 test:
60 suffix: 2d_tri
61 args:
62
63 test:
64 suffix: 3d_tri
65 args: -dm_plex_dim 3
66
67 testset:
68 requires: !__float128
69 args: -dm_plex_shape sphere -petscpartitioner_type simple -dm_view
70
71 test:
72 suffix: 2d_quad_parallel
73 nsize: 2
74 args: -dm_plex_simplex 0
75
76 test:
77 suffix: 2d_tri_parallel
78 nsize: 2
79
80 test:
81 suffix: 3d_tri_parallel
82 nsize: 2
83 args: -dm_plex_dim 3
84
85 TEST*/
86