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: PetscOptionsGetBool(NULL, NULL, "-label_mesh", &flg, NULL);
15: if (!flg) return 0;
16: DMCreateLabel(dm, "test");
17: DMGetLabel(dm, "test", &label);
18: DMPlexGetChart(dm, &pStart, &pEnd);
19: for (p = pStart; p < pEnd; ++p) DMLabelSetValue(label, p, p);
20: return 0;
21: }
23: static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm)
24: {
25: DMCreate(comm, dm);
26: DMSetType(*dm, DMPLEX);
27: DMSetFromOptions(*dm);
28: LabelPoints(*dm);
29: PetscObjectSetOptionsPrefix((PetscObject)*dm, "post_label_");
30: DMSetFromOptions(*dm);
31: PetscObjectSetOptionsPrefix((PetscObject)*dm, NULL);
32: DMViewFromOptions(*dm, NULL, "-dm_view");
33: return 0;
34: }
36: int main(int argc, char **argv)
37: {
38: DM dm;
41: PetscInitialize(&argc, &argv, NULL, help);
42: CreateMesh(PETSC_COMM_WORLD, &dm);
43: DMDestroy(&dm);
44: PetscFinalize();
45: return 0;
46: }
48: /*TEST
49: test:
50: suffix: ref_seg
51: args: -dm_plex_reference_cell_domain -dm_plex_cell segment -dm_refine 1 -dm_plex_check_all
53: test:
54: suffix: ref_tri
55: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2 -dm_plex_check_all
57: test:
58: suffix: box_tri
59: requires: triangle parmetis
60: nsize: {{1 3 5}}
61: args: -dm_plex_box_faces 3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
63: test:
64: suffix: ref_quad
65: args: -dm_plex_reference_cell_domain -dm_plex_cell quadrilateral -dm_refine 2 -dm_plex_check_all
67: test:
68: suffix: box_quad
69: nsize: {{1 3 5}}
70: requires: parmetis
71: args: -dm_plex_box_faces 3,3 -dm_plex_simplex 0 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
73: test:
74: suffix: ref_tet
75: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2 -dm_plex_check_all
77: test:
78: suffix: box_tet
79: requires: ctetgen parmetis
80: nsize: {{1 3 5}}
81: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -dm_plex_check_all -petscpartitioner_type parmetis
83: test:
84: suffix: ref_hex
85: args: -dm_plex_reference_cell_domain -dm_plex_cell hexahedron -dm_refine 2 -dm_plex_check_all
87: test:
88: suffix: box_hex
89: requires: parmetis
90: nsize: {{1 3 5}}
91: 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
93: test:
94: suffix: ref_trip
95: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2 -dm_plex_check_all
97: test:
98: suffix: ref_tquad
99: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad -dm_refine 2 -dm_plex_check_all
101: test:
102: suffix: ref_ttrip
103: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2 -dm_plex_check_all
105: test:
106: suffix: ref_tquadp
107: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2 -dm_plex_check_all
109: test:
110: suffix: ref_pyramid
111: args: -dm_plex_reference_cell_domain -dm_plex_cell pyramid -dm_refine 2 -dm_plex_check_all
113: testset:
114: args: -dm_coord_space 0 -dm_plex_transform_type refine_tobox -dm_plex_check_all
116: test:
117: suffix: ref_tri_tobox
118: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_refine 2
120: test:
121: suffix: box_tri_tobox
122: requires: triangle parmetis
123: nsize: {{1 3 5}}
124: args: -dm_plex_box_faces 3,3 -dm_refine 2 -petscpartitioner_type parmetis
126: test:
127: suffix: ref_tet_tobox
128: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron -dm_refine 2
130: test:
131: suffix: box_tet_tobox
132: requires: ctetgen parmetis
133: nsize: {{1 3 5}}
134: args: -dm_plex_dim 3 -dm_plex_box_faces 3,3,3 -dm_refine 2 -petscpartitioner_type parmetis
136: test:
137: suffix: ref_trip_tobox
138: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism -dm_refine 2
140: test:
141: suffix: ref_ttrip_tobox
142: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism -dm_refine 2
144: test:
145: suffix: ref_tquadp_tobox
146: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism -dm_refine 2
148: testset:
149: args: -dm_coord_space 0 -label_mesh -post_label_dm_extrude 2 -post_label_dm_plex_check_all -dm_view ::ascii_info_detail
151: test:
152: suffix: extrude_quad
153: args: -dm_plex_simplex 0
155: TEST*/