Actual source code: ex40.c
1: static const char help[] = "Tests for Plex transforms, including regular refinement";
3: #include <petscdmplex.h>
4: #include <petscsf.h>
6: #include <petsc/private/dmpleximpl.h>
8: static PetscErrorCode LabelPoints(DM dm)
9: {
10: DMLabel label;
11: PetscInt pStart, pEnd, p;
12: PetscBool flg = PETSC_FALSE;
14: PetscFunctionBegin;
15: PetscCall(PetscOptionsGetBool(NULL, NULL, "-label_mesh", &flg, NULL));
16: if (!flg) PetscFunctionReturn(PETSC_SUCCESS);
17: PetscCall(DMCreateLabel(dm, "test"));
18: PetscCall(DMGetLabel(dm, "test", &label));
19: PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
20: for (p = pStart; p < pEnd; ++p) PetscCall(DMLabelSetValue(label, p, p));
21: PetscFunctionReturn(PETSC_SUCCESS);
22: }
24: static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm)
25: {
26: PetscFunctionBegin;
27: PetscCall(DMCreate(comm, dm));
28: PetscCall(DMSetType(*dm, DMPLEX));
29: PetscCall(DMSetFromOptions(*dm));
30: PetscCall(LabelPoints(*dm));
31: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, "post_label_"));
32: PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
33: PetscCall(DMSetFromOptions(*dm));
34: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)*dm, NULL));
35: PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
36: PetscFunctionReturn(PETSC_SUCCESS);
37: }
39: int main(int argc, char **argv)
40: {
41: DM dm;
43: PetscFunctionBeginUser;
44: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
45: PetscCall(CreateMesh(PETSC_COMM_WORLD, &dm));
46: PetscCall(DMDestroy(&dm));
47: PetscCall(PetscFinalize());
48: return 0;
49: }
51: /*TEST
52: test:
53: suffix: ref_seg
54: args: -dm_plex_reference_cell_domain -dm_plex_cell segment -dm_refine 1 -dm_plex_check_all
56: test:
57: suffix: ref_tri
58: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2 -dm_plex_check_all
60: test:
61: suffix: box_tri
62: requires: triangle parmetis
63: nsize: {{1 3 5}}
64: args: -dm_plex_box_faces 3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
66: test:
67: suffix: ref_quad
68: args: -dm_plex_reference_cell_domain -dm_plex_cell quadrilateral -dm_refine 2 -dm_plex_check_all
70: test:
71: suffix: box_quad
72: nsize: {{1 3 5}}
73: requires: parmetis
74: args: -dm_plex_box_faces 3,3 -dm_plex_simplex 0 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
76: test:
77: suffix: ref_tet
78: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2 -dm_plex_check_all
80: test:
81: suffix: box_tet
82: requires: ctetgen parmetis
83: nsize: {{1 3 5}}
84: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
86: test:
87: suffix: ref_hex
88: args: -dm_plex_reference_cell_domain -dm_plex_cell hexahedron -dm_refine 2 -dm_plex_check_all
90: test:
91: suffix: box_hex
92: requires: parmetis
93: nsize: {{1 3 5}}
94: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_plex_simplex 0 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
96: test:
97: suffix: ref_trip
98: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2 -dm_plex_check_all
100: test:
101: suffix: ref_tquad
102: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad -dm_refine 2 -dm_plex_check_all
104: test:
105: suffix: ref_ttrip
106: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2 -dm_plex_check_all
108: test:
109: suffix: ref_tquadp
110: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2 -dm_plex_check_all
112: test:
113: suffix: ref_pyramid
114: args: -dm_plex_reference_cell_domain -dm_plex_cell pyramid -dm_refine 2 -dm_plex_check_all
116: testset:
117: args: -dm_coord_space 0 -dm_plex_transform_type refine_tobox -dm_plex_check_all
119: test:
120: suffix: ref_tri_tobox
121: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2
123: test:
124: suffix: box_tri_tobox
125: requires: triangle parmetis
126: nsize: {{1 3 5}}
127: args: -dm_plex_box_faces 3,3 -dm_refine 2 -petscpartitioner_type parmetis
129: test:
130: suffix: ref_tet_tobox
131: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2
133: test:
134: suffix: box_tet_tobox
135: requires: ctetgen parmetis
136: nsize: {{1 3 5}}
137: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -petscpartitioner_type parmetis
139: test:
140: suffix: ref_trip_tobox
141: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2
143: test:
144: suffix: ref_ttrip_tobox
145: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2
147: test:
148: suffix: ref_tquadp_tobox
149: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2
151: testset:
152: args: -dm_coord_space 0 -label_mesh -post_label_dm_extrude 2 -post_label_dm_plex_check_all -dm_view ::ascii_info_detail
154: test:
155: suffix: extrude_quad
156: args: -dm_plex_simplex 0
158: TEST*/