1c4762a1bSJed Brown static char help[] = "Define a simple field over the mesh\n\n";
2c4762a1bSJed Brown
3c4762a1bSJed Brown #include <petscdmplex.h>
4c4762a1bSJed Brown
main(int argc,char ** argv)5d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv)
6d71ae5a4SJacob Faibussowitsch {
730602db0SMatthew G. Knepley DM dm;
8c4762a1bSJed Brown Vec u;
9c4762a1bSJed Brown PetscSection section;
10*3124aba7SMatthew G. Knepley PetscInt dim, numFields, numBC;
11c4762a1bSJed Brown PetscInt numComp[3];
12c4762a1bSJed Brown PetscInt numDof[12];
13c4762a1bSJed Brown PetscInt bcField[1];
14c4762a1bSJed Brown IS bcPointIS[1];
15c4762a1bSJed Brown
16327415f7SBarry Smith PetscFunctionBeginUser;
179566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, NULL, help));
18c4762a1bSJed Brown /* Create a mesh */
199566063dSJacob Faibussowitsch PetscCall(DMCreate(PETSC_COMM_WORLD, &dm));
209566063dSJacob Faibussowitsch PetscCall(DMSetType(dm, DMPLEX));
219566063dSJacob Faibussowitsch PetscCall(DMSetFromOptions(dm));
229566063dSJacob Faibussowitsch PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
239566063dSJacob Faibussowitsch PetscCall(DMGetDimension(dm, &dim));
24c4762a1bSJed Brown /* Create a scalar field u, a vector field v, and a surface vector field w */
25c4762a1bSJed Brown numFields = 3;
26c4762a1bSJed Brown numComp[0] = 1;
27c4762a1bSJed Brown numComp[1] = dim;
28c4762a1bSJed Brown numComp[2] = dim - 1;
29*3124aba7SMatthew G. Knepley for (PetscInt i = 0; i < numFields * (dim + 1); ++i) numDof[i] = 0;
30c4762a1bSJed Brown /* Let u be defined on vertices */
31c4762a1bSJed Brown numDof[0 * (dim + 1) + 0] = 1;
32c4762a1bSJed Brown /* Let v be defined on cells */
33c4762a1bSJed Brown numDof[1 * (dim + 1) + dim] = dim;
34c4762a1bSJed Brown /* Let w be defined on faces */
35c4762a1bSJed Brown numDof[2 * (dim + 1) + dim - 1] = dim - 1;
36c4762a1bSJed Brown /* Setup boundary conditions */
37c4762a1bSJed Brown numBC = 1;
38c4762a1bSJed Brown /* Prescribe a Dirichlet condition on u on the boundary
39c4762a1bSJed Brown Label "marker" is made by the mesh creation routine */
40c4762a1bSJed Brown bcField[0] = 0;
419566063dSJacob Faibussowitsch PetscCall(DMGetStratumIS(dm, "marker", 1, &bcPointIS[0]));
42c4762a1bSJed Brown /* Create a PetscSection with this data layout */
439566063dSJacob Faibussowitsch PetscCall(DMSetNumFields(dm, numFields));
449566063dSJacob Faibussowitsch PetscCall(DMPlexCreateSection(dm, NULL, numComp, numDof, numBC, bcField, NULL, bcPointIS, NULL, §ion));
459566063dSJacob Faibussowitsch PetscCall(ISDestroy(&bcPointIS[0]));
46c4762a1bSJed Brown /* Name the Field variables */
479566063dSJacob Faibussowitsch PetscCall(PetscSectionSetFieldName(section, 0, "u"));
489566063dSJacob Faibussowitsch PetscCall(PetscSectionSetFieldName(section, 1, "v"));
499566063dSJacob Faibussowitsch PetscCall(PetscSectionSetFieldName(section, 2, "w"));
500cae8135SMatthew G. Knepley PetscCall(PetscSectionViewFromOptions(section, NULL, "-mysection_view"));
51c4762a1bSJed Brown /* Tell the DM to use this data layout */
529566063dSJacob Faibussowitsch PetscCall(DMSetLocalSection(dm, section));
53c4762a1bSJed Brown /* Create a Vec with this layout and view it */
549566063dSJacob Faibussowitsch PetscCall(DMGetGlobalVector(dm, &u));
55*3124aba7SMatthew G. Knepley PetscCall(PetscObjectSetName((PetscObject)u, "Solution"));
56*3124aba7SMatthew G. Knepley PetscCall(VecSet(u, 1.));
57*3124aba7SMatthew G. Knepley PetscCall(VecViewFromOptions(u, NULL, "-vec_view"));
589566063dSJacob Faibussowitsch PetscCall(DMRestoreGlobalVector(dm, &u));
59c4762a1bSJed Brown /* Cleanup */
609566063dSJacob Faibussowitsch PetscCall(PetscSectionDestroy(§ion));
619566063dSJacob Faibussowitsch PetscCall(DMDestroy(&dm));
629566063dSJacob Faibussowitsch PetscCall(PetscFinalize());
63b122ec5aSJacob Faibussowitsch return 0;
64c4762a1bSJed Brown }
65c4762a1bSJed Brown
66c4762a1bSJed Brown /*TEST
67c4762a1bSJed Brown
68c4762a1bSJed Brown test:
69c4762a1bSJed Brown suffix: 0
70c4762a1bSJed Brown requires: triangle
71*3124aba7SMatthew G. Knepley args: -mysection_view -vec_view vtk:sol.vtu -info :~sys,mat
720cae8135SMatthew G. Knepley
73c4762a1bSJed Brown test:
74c4762a1bSJed Brown suffix: 1
75c4762a1bSJed Brown requires: ctetgen
760cae8135SMatthew G. Knepley args: -dm_plex_dim 3 -mysection_view -info :~sys,mat
770cae8135SMatthew G. Knepley
780cae8135SMatthew G. Knepley test:
790cae8135SMatthew G. Knepley suffix: filter_0
800cae8135SMatthew G. Knepley args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3 -dm_plex_label_lower 0,1,2,3,4,6 \
810cae8135SMatthew G. Knepley -dm_refine 1 -dm_plex_transform_type transform_filter -dm_plex_transform_active lower \
820cae8135SMatthew G. Knepley -dm_view
830cae8135SMatthew G. Knepley
840cae8135SMatthew G. Knepley test:
850cae8135SMatthew G. Knepley suffix: submesh_0
860cae8135SMatthew G. Knepley requires: triangle
870cae8135SMatthew G. Knepley args: -dm_plex_label_middle 9,12,15,23,31 -dm_plex_submesh middle -dm_view
88c4762a1bSJed Brown
89c4762a1bSJed Brown TEST*/
90