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
55: output_file: output/empty.out
57: test:
58: suffix: ref_tri
59: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2 -dm_plex_check_all
60: output_file: output/empty.out
62: test:
63: suffix: box_tri
64: requires: triangle parmetis
65: nsize: {{1 3 5}}
66: args: -dm_plex_box_faces 3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
67: output_file: output/empty.out
69: test:
70: suffix: ref_quad
71: args: -dm_plex_reference_cell_domain -dm_plex_cell quadrilateral -dm_refine 2 -dm_plex_check_all
72: output_file: output/empty.out
74: test:
75: suffix: box_quad
76: nsize: {{1 3 5}}
77: requires: parmetis
78: args: -dm_plex_box_faces 3,3 -dm_plex_simplex 0 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
79: output_file: output/empty.out
81: test:
82: suffix: box_quad_label
83: args: -dm_plex_box_faces 3,3 -dm_plex_simplex 0 -dm_refine 2 -dm_plex_check_all -dm_plex_transform_label_match_strata {{0 1}separate output} -dm_view
85: test:
86: suffix: ref_tet
87: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2 -dm_plex_check_all
88: output_file: output/empty.out
90: test:
91: suffix: box_tet
92: requires: ctetgen parmetis
93: nsize: {{1 3 5}}
94: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
95: output_file: output/empty.out
97: test:
98: suffix: ref_hex
99: args: -dm_plex_reference_cell_domain -dm_plex_cell hexahedron -dm_refine 2 -dm_plex_check_all
100: output_file: output/empty.out
102: test:
103: suffix: box_hex
104: requires: parmetis
105: nsize: {{1 3 5}}
106: 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
107: output_file: output/empty.out
109: test:
110: suffix: ref_trip
111: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2 -dm_plex_check_all
112: output_file: output/empty.out
114: test:
115: suffix: ref_tquad
116: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad -dm_refine 2 -dm_plex_check_all
117: output_file: output/empty.out
119: test:
120: suffix: ref_ttrip
121: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2 -dm_plex_check_all
122: output_file: output/empty.out
124: test:
125: suffix: ref_tquadp
126: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2 -dm_plex_check_all
127: output_file: output/empty.out
129: test:
130: suffix: ref_pyramid
131: args: -dm_plex_reference_cell_domain -dm_plex_cell pyramid -dm_refine 2 -dm_plex_check_all
132: output_file: output/empty.out
134: testset:
135: args: -dm_coord_space 0 -dm_plex_transform_type refine_tobox -dm_plex_check_all
137: test:
138: suffix: ref_tri_tobox
139: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2
140: output_file: output/empty.out
142: test:
143: suffix: box_tri_tobox
144: requires: triangle parmetis
145: nsize: {{1 3 5}}
146: args: -dm_plex_box_faces 3,3 -dm_refine 2 -petscpartitioner_type parmetis
147: output_file: output/empty.out
149: test:
150: suffix: ref_tet_tobox
151: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2
152: output_file: output/empty.out
154: test:
155: suffix: box_tet_tobox
156: requires: ctetgen parmetis
157: nsize: {{1 3 5}}
158: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -petscpartitioner_type parmetis
159: output_file: output/empty.out
161: test:
162: suffix: ref_trip_tobox
163: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2
164: output_file: output/empty.out
166: test:
167: suffix: ref_ttrip_tobox
168: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2
169: output_file: output/empty.out
171: test:
172: suffix: ref_tquadp_tobox
173: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2
174: output_file: output/empty.out
176: testset:
177: args: -dm_coord_space 0 -label_mesh -post_label_dm_extrude 2 -post_label_dm_plex_check_all -dm_view ::ascii_info_detail
179: test:
180: suffix: extrude_quad
181: args: -dm_plex_simplex 0
183: TEST*/