Actual source code: ex2.c
1: static char help[] = "PETSc Annual Meeting 2025: Meshing Tutorial.\n\n\n";
3: #include <petsc.h>
5: /* Use FORCE=1 to run the PyVista tests */
7: static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm)
8: {
9: PetscFunctionBeginUser;
10: PetscCall(DMCreate(comm, dm));
11: PetscCall(DMSetType(*dm, DMPLEX));
12: PetscCall(DMSetFromOptions(*dm));
13: PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
14: PetscFunctionReturn(PETSC_SUCCESS);
15: }
17: int main(int argc, char **argv)
18: {
19: DM dm;
21: PetscFunctionBeginUser;
22: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
23: PetscCall(CreateMesh(PETSC_COMM_WORLD, &dm));
24: PetscCall(DMDestroy(&dm));
25: PetscCall(PetscFinalize());
26: return PETSC_SUCCESS;
27: }
29: /*TEST
31: # Draw a square with X, use with -draw_pause -1
32: test:
33: suffix: 0
34: requires: triangle x
35: args: -dm_view draw -draw_pause 0
36: output_file: output/empty.out
38: # Draw a square with PyVista
39: test:
40: suffix: 1
41: requires: triangle pyvista
42: args: -dm_view pyvista
43: output_file: output/empty.out
45: # Refine the square
46: test:
47: suffix: 2
48: requires: triangle pyvista
49: args: -dm_view pyvista -dm_refine 1
50: output_file: output/empty.out
52: # Refine the square three times
53: test:
54: suffix: 3
55: requires: triangle pyvista
56: args: -dm_view pyvista -dm_refine 3
57: output_file: output/empty.out
59: # Refine the cube three times
60: test:
61: suffix: 4
62: requires: ctetgen pyvista
63: args: -dm_view pyvista -dm_plex_dim 3 -dm_refine 3
64: output_file: output/empty.out
66: # Draw a sphere with PyVista (all we get is an icosahedron)
67: test:
68: suffix: 5
69: requires: pyvista
70: args: -dm_view pyvista -dm_plex_shape sphere
71: output_file: output/empty.out
73: # Draw the 2-sphere with PyVista with quadratic coordinates
74: test:
75: suffix: 5b
76: requires: pyvista
77: args: -dm_view pyvista -dm_plex_shape sphere \
78: -dm_coord_petscspace_degree 2
79: output_file: output/empty.out
81: # Draw the 2-sphere with PyVista using quads with quadratic coordinates
82: test:
83: suffix: 5c
84: requires: pyvista
85: args: -dm_view pyvista -dm_plex_shape sphere -dm_plex_simplex 0 \
86: -dm_coord_petscspace_degree 2
87: output_file: output/empty.out
89: # Refine the sphere three times
90: test:
91: suffix: 6
92: requires: pyvista
93: args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 3
94: output_file: output/empty.out
96: # Draw the 1-sphere with PyVista, which is a circle
97: test:
98: suffix: 7a
99: requires: pyvista
100: args: -dm_view pyvista -view_pyvista_line_width 5 -dm_plex_dim 1 -dm_plex_shape sphere
101: output_file: output/empty.out
103: # Draw the 1-sphere with PyVista with cubic coordinates
104: test:
105: suffix: 7b
106: requires: pyvista
107: args: -dm_view pyvista -view_pyvista_line_width 5 -dm_plex_dim 1 -dm_plex_shape sphere \
108: -dm_coord_petscspace_degree 3
109: output_file: output/empty.out
111: # Show the 3-sphere
112: test:
113: suffix: 7
114: args: -dm_view -dm_plex_shape sphere -dm_plex_dim 3
116: # Extrude the sphere
117: test:
118: suffix: 8
119: requires: pyvista
120: args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 1 \
121: -dm_plex_transform_type extrude -dm_plex_transform_extrude_layers 3 \
122: -dm_plex_transform_extrude_use_tensor 0 -dm_plex_transform_extrude_thickness 0.3
123: output_file: output/empty.out
125: # Extrude the sphere with cutaway
126: test:
127: suffix: 9
128: requires: pyvista
129: args: -dm_view pyvista -view_pyvista_clip 1.,0.,0. -dm_plex_shape sphere -dm_refine 1 \
130: -dm_plex_transform_type extrude -dm_plex_transform_extrude_layers 3 \
131: -dm_plex_transform_extrude_use_tensor 0 -dm_plex_transform_extrude_thickness 0.3
132: output_file: output/empty.out
134: # Extrude the refined sphere
135: test:
136: suffix: 10
137: requires: pyvista
138: args: -dm_view pyvista -view_pyvista_clip 1.,0.,0. -dm_plex_shape sphere -dm_refine 3 \
139: -dm_plex_option_phases ext_ \
140: -ext_dm_refine 1 -ext_dm_plex_transform_type extrude -ext_dm_plex_transform_extrude_layers 3 \
141: -ext_dm_plex_transform_extrude_use_tensor 0 -ext_dm_plex_transform_extrude_thickness 0.5
142: output_file: output/empty.out
144: # Extrude the refined sphere
145: test:
146: suffix: 11
147: requires: pyvista
148: args: -dm_view pyvista -view_pyvista_clip 1.,0.,0. -dm_plex_shape sphere -dm_refine 3 \
149: -dm_plex_option_phases ext_,ref_ \
150: -ext_dm_refine 1 -ext_dm_plex_transform_type extrude -ext_dm_plex_transform_extrude_layers 3 \
151: -ext_dm_plex_transform_extrude_use_tensor 0 -ext_dm_plex_transform_extrude_thickness 0.5 \
152: -ref_dm_refine 1
153: output_file: output/empty.out
155: # Refine the Schwartz surface
156: test:
157: suffix: 12
158: requires: pyvista
159: args: -dm_view pyvista -dm_plex_shape schwarz_p -dm_plex_tps_extent 3,2 -dm_plex_tps_refine 2
160: output_file: output/empty.out
162: # Refine and extrude the Schwartz surface with given thickness
163: test:
164: suffix: 13
165: requires: pyvista
166: args: -dm_view pyvista -dm_plex_shape schwarz_p -dm_plex_tps_extent 3,2,2 -dm_plex_tps_refine 3 \
167: -dm_plex_tps_thickness 0.4 -dm_plex_tps_layers 3
168: output_file: output/empty.out
170: # Filter the square
171: test:
172: suffix: 14
173: requires: triangle pyvista
174: args: -dm_view pyvista -dm_refine 1 \
175: -dm_plex_transform_type transform_filter -dm_plex_transform_active filter_cells -dm_plex_label_filter_cells 0,1,2
176: output_file: output/empty.out
178: # Filter a cap on the sphere
179: test:
180: suffix: 15
181: requires: pyvista
182: args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 4 \
183: -dm_plex_option_phases filt_ \
184: -filt_dm_refine 1 -filt_dm_plex_transform_type transform_filter \
185: -filt_dm_plex_transform_active filter_cells -dm_plex_label_filter_cells 0
186: output_file: output/empty.out
188: # Filter the sphere minus a cap
189: test:
190: suffix: 16
191: requires: pyvista
192: args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 4 \
193: -dm_plex_option_phases filt_ \
194: -filt_dm_refine 1 -filt_dm_plex_transform_type transform_filter \
195: -filt_dm_plex_transform_active filter_cells \
196: -dm_plex_label_filter_cells 0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18
197: output_file: output/empty.out
199: # Convert sphere to quads
200: test:
201: suffix: 17
202: requires: pyvista
203: args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 3 \
204: -dm_plex_option_phases conv_ -conv_dm_refine 1 -conv_dm_plex_transform_type refine_tobox
205: output_file: output/empty.out
207: # Load and refine the nozzle
208: test:
209: suffix: 18
210: requires: pyvista egads datafilespath
211: args: -dm_view pyvista -dm_plex_filename /Users/knepley/PETSc4/petsc/petsc-dev/share/petsc/datafiles/meshes/nozzle.stp -dm_refine 3
212: output_file: output/empty.out
214: # Load and refine the nozzle, and convert to quads
215: test:
216: suffix: 19
217: requires: pyvista egads datafilespath
218: args: -dm_view pyvista -dm_plex_filename /Users/knepley/PETSc4/petsc/petsc-dev/share/petsc/datafiles/meshes/nozzle.stp -dm_refine 3 \
219: -dm_plex_option_phases conv_ -conv_dm_refine 1 -conv_dm_plex_transform_type refine_tobox
220: output_file: output/empty.out
222: # Extrude a periodic mesh, here the circle
223: test:
224: suffix: 20
225: args: -dm_plex_dim 1 -dm_plex_box_faces 8 -dm_plex_box_bd periodic \
226: -dm_refine 1 -dm_plex_transform_type extrude -dm_plex_transform_extrude_layers 2 \
227: -dm_view
229: # Construct and refine a torus
230: test:
231: suffix: 21
232: args: -dm_plex_option_phases phase1_,phase2_,phase3_ \
233: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 \
234: -dm_refine_pre 0 \
235: -phase1_dm_refine 0 \
236: -phase2_dm_plex_transform_type extrude -phase2_dm_extrude 16 \
237: -phase2_dm_plex_transform_extrude_thickness 1 \
238: -phase2_dm_plex_transform_extrude_use_tensor 0 \
239: -phase2_dm_plex_transform_extrude_periodic \
240: -phase3_dm_coord_remap -phase3_dm_coord_map torus \
241: -dm_view
243: # Construct a torus with a custom poloidial slice
244: # Looks better with -phase2_dm_refine 2 -phase3_dm_extrude 32
245: test:
246: suffix: 22
247: requires: triangle
248: args: -dm_plex_option_phases phase1_,phase2_,phase3_,phase4_ \
249: -dm_plex_shape diiid -dm_refine 0 -bd_dm_plex_generate_triangle_opts zpQq30e \
250: -phase1_dm_coord_remap -phase1_dm_coord_map rotate -phase1_dm_coord_map_params 0,0,0,0,1.5707 \
251: -phase2_dm_refine 0 \
252: -phase3_dm_plex_transform_type extrude -phase3_dm_extrude 16 \
253: -phase3_dm_plex_transform_extrude_thickness 1 \
254: -phase3_dm_plex_transform_extrude_use_tensor 0 \
255: -phase3_dm_plex_transform_extrude_periodic \
256: -phase4_dm_coord_remap -phase4_dm_coord_map torus \
257: -dm_view
259: TEST*/