1 static const char help[] = "Test for non-manifold interpolation";
2
3 #include <petscdmplex.h>
4
5 /*
6 Test 0:
7
8 3-------------7
9 /| /|
10 / | / |
11 / | / |
12 1-------------5 |
13 | | | |
14 | | | |
15 | | | |
16 | | | |
17 z 4---------|---8
18 ^ / | /
19 | y | /
20 |/ |/
21 2--->-x-------6-------------9
22
23 Test 1:
24
25 3-------------7
26 /| /|
27 / | / |
28 / | / |
29 1-------------5 |
30 | | | |
31 | | | |
32 | | | |
33 | | | |
34 z 4---------|---8
35 ^ / | / \
36 | y | / \
37 |/ |/ \
38 2--->-x-------6-------9
39 */
main(int argc,char ** argv)40 int main(int argc, char **argv)
41 {
42 DM dm, idm;
43 DMLabel ctLabel;
44 PetscInt testNum = 0;
45
46 PetscFunctionBeginUser;
47 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
48 PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Non-Manifold Options", "ex66");
49 PetscCall(PetscOptionsInt("-test_num", "Test number", "", testNum, &testNum, NULL));
50 PetscOptionsEnd();
51
52 PetscCall(DMCreate(PETSC_COMM_WORLD, &dm));
53 PetscCall(PetscObjectSetName((PetscObject)dm, "cubeline-fromdag"));
54 PetscCall(DMSetType(dm, DMPLEX));
55 PetscCall(DMSetDimension(dm, 3));
56
57 switch (testNum) {
58 case 0: {
59 // 9 vertices
60 // 1 edge
61 // 0 faces
62 // 1 volume
63 PetscInt num_points[4] = {9, 1, 0, 1};
64
65 // point 0 = hexahedron (defined by 8 vertices)
66 // points 1-9 = vertices
67 // point 10 = edged (defined by 2 vertices)
68 PetscInt cone_size[11] = {8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2};
69
70 // hexahedron defined by points
71 PetscInt cones[11] = {3, 4, 2, 1, 7, 5, 6, 8, 6, 9};
72 PetscInt cone_orientations[11] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
73 PetscScalar vertex_coords[3 * 9] = {0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 0, 1, 1, 1, 1, 1, 0, 2, 0, 0};
74
75 PetscCall(DMPlexCreateFromDAG(dm, 3, num_points, cone_size, cones, cone_orientations, vertex_coords));
76 } break;
77 case 1: {
78 // 9 vertices
79 // 0 edges
80 // 1 face
81 // 1 volume
82 PetscInt num_points[4] = {9, 0, 1, 1};
83
84 // point 0 = hexahedron (defined by 8 vertices)
85 // points 1-9 = vertices
86 // point 10 = triangle (defined by 3 vertices)
87 PetscInt cone_size[11] = {8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3};
88
89 // hexahedron defined by 8 points (point 0)
90 PetscInt cones[11] = {3, 4, 2, 1, 7, 5, 6, 8, 6, 9, 8};
91 PetscInt cone_orientations[11] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
92 PetscScalar vertex_coords[3 * 9] = {0, 0, 1, // point 1
93 0, 0, 0, // point 2
94 0, 1, 1, // point 3
95 0, 1, 0, // point 4
96 1, 0, 1, // point 5
97 1, 0, 0, // point 6
98 1, 1, 1, // point 7
99 1, 1, 0, // point 8
100 2, 0, 0}; // point 9
101
102 PetscCall(DMPlexCreateFromDAG(dm, 3, num_points, cone_size, cones, cone_orientations, vertex_coords));
103 } break;
104 default:
105 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Invalid test number %" PetscInt_FMT, testNum);
106 }
107 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
108
109 // TODO: make it work with a DM made from a msh file
110 // PetscCall(DMPlexCreateGmshFromFile(PETSC_COMM_WORLD, "cube-line.msh", PETSC_FALSE, &dm));
111 // PetscCall(PetscObjectSetName((PetscObject)dm, "msh"));
112
113 // Must set cell types
114 PetscCall(DMPlexGetCellTypeLabel(dm, &ctLabel));
115 switch (testNum) {
116 case 0:
117 PetscCall(DMLabelSetValue(ctLabel, 0, DM_POLYTOPE_HEXAHEDRON));
118 PetscCall(DMLabelSetValue(ctLabel, 1, DM_POLYTOPE_POINT));
119 PetscCall(DMLabelSetValue(ctLabel, 2, DM_POLYTOPE_POINT));
120 PetscCall(DMLabelSetValue(ctLabel, 3, DM_POLYTOPE_POINT));
121 PetscCall(DMLabelSetValue(ctLabel, 4, DM_POLYTOPE_POINT));
122 PetscCall(DMLabelSetValue(ctLabel, 5, DM_POLYTOPE_POINT));
123 PetscCall(DMLabelSetValue(ctLabel, 6, DM_POLYTOPE_POINT));
124 PetscCall(DMLabelSetValue(ctLabel, 7, DM_POLYTOPE_POINT));
125 PetscCall(DMLabelSetValue(ctLabel, 8, DM_POLYTOPE_POINT));
126 PetscCall(DMLabelSetValue(ctLabel, 9, DM_POLYTOPE_POINT));
127 PetscCall(DMLabelSetValue(ctLabel, 10, DM_POLYTOPE_SEGMENT));
128 break;
129 case 1:
130 PetscCall(DMLabelSetValue(ctLabel, 0, DM_POLYTOPE_HEXAHEDRON));
131 PetscCall(DMLabelSetValue(ctLabel, 1, DM_POLYTOPE_POINT));
132 PetscCall(DMLabelSetValue(ctLabel, 2, DM_POLYTOPE_POINT));
133 PetscCall(DMLabelSetValue(ctLabel, 3, DM_POLYTOPE_POINT));
134 PetscCall(DMLabelSetValue(ctLabel, 4, DM_POLYTOPE_POINT));
135 PetscCall(DMLabelSetValue(ctLabel, 5, DM_POLYTOPE_POINT));
136 PetscCall(DMLabelSetValue(ctLabel, 6, DM_POLYTOPE_POINT));
137 PetscCall(DMLabelSetValue(ctLabel, 7, DM_POLYTOPE_POINT));
138 PetscCall(DMLabelSetValue(ctLabel, 8, DM_POLYTOPE_POINT));
139 PetscCall(DMLabelSetValue(ctLabel, 9, DM_POLYTOPE_POINT));
140 PetscCall(DMLabelSetValue(ctLabel, 10, DM_POLYTOPE_TRIANGLE));
141 break;
142 default:
143 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Invalid test number %" PetscInt_FMT, testNum);
144 }
145
146 // interpolate (make sure to use -interp_dm_plex_stratify_celltype)
147 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)dm, "interp_"));
148 PetscCall(DMPlexInterpolate(dm, &idm));
149 PetscCall(DMDestroy(&dm));
150 dm = idm;
151
152 PetscCall(DMSetFromOptions(dm));
153 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
154
155 PetscCall(DMDestroy(&dm));
156 PetscCall(PetscFinalize());
157 return 0;
158 }
159
160 /*TEST
161
162 testset:
163 args: -interp_dm_plex_stratify_celltype -dm_view ::ascii_info_detail -interp_dm_view ::ascii_info_detail
164
165 test:
166 suffix: 0
167 test:
168 suffix: 1
169 args: -test_num 1
170
171 TEST*/
172