1 static char help[] = "PETSc Annual Meeting 2025: Meshing Tutorial.\n\n\n"; 2 3 #include <petsc.h> 4 5 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm) 6 { 7 PetscFunctionBeginUser; 8 PetscCall(DMCreate(comm, dm)); 9 PetscCall(DMSetType(*dm, DMPLEX)); 10 PetscCall(DMSetFromOptions(*dm)); 11 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 12 PetscFunctionReturn(PETSC_SUCCESS); 13 } 14 15 int main(int argc, char **argv) 16 { 17 DM dm; 18 19 PetscFunctionBeginUser; 20 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 21 PetscCall(CreateMesh(PETSC_COMM_WORLD, &dm)); 22 PetscCall(DMDestroy(&dm)); 23 PetscCall(PetscFinalize()); 24 return PETSC_SUCCESS; 25 } 26 27 /*TEST 28 29 # Draw a square with X, use with -draw_pause -1 30 test: 31 suffix: 0 32 requires: triangle x 33 args: -dm_view draw -draw_pause 0 34 35 # Draw a square with PyVista 36 test: 37 suffix: 1 38 requires: triangle pyvista 39 args: -dm_view pyvista 40 output_file: output/ex2_0.out 41 42 # Refine the square 43 test: 44 suffix: 2 45 requires: triangle pyvista 46 args: -dm_view pyvista -dm_refine 1 47 output_file: output/ex2_0.out 48 49 # Refine the square three times 50 test: 51 suffix: 3 52 requires: triangle pyvista 53 args: -dm_view pyvista -dm_refine 3 54 output_file: output/ex2_0.out 55 56 # Refine the cube three times 57 test: 58 suffix: 4 59 requires: ctetgen pyvista 60 args: -dm_view pyvista -dm_plex_dim 3 -dm_refine 3 61 output_file: output/ex2_0.out 62 63 # Draw a sphere with PyVista (all we get is an icosahedron) 64 test: 65 suffix: 5 66 requires: pyvista 67 args: -dm_view pyvista -dm_plex_shape sphere 68 output_file: output/ex2_0.out 69 70 # Refine the sphere three times 71 test: 72 suffix: 6 73 requires: pyvista 74 args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 3 75 output_file: output/ex2_0.out 76 77 # Show the 3-sphere 78 test: 79 suffix: 7 80 args: -dm_view -dm_plex_shape sphere -dm_plex_dim 3 81 82 # Extrude the sphere 83 test: 84 suffix: 8 85 requires: pyvista 86 args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 1 \ 87 -dm_plex_transform_type extrude -dm_plex_transform_extrude_layers 3 \ 88 -dm_plex_transform_extrude_use_tensor 0 -dm_plex_transform_extrude_thickness 0.3 89 output_file: output/ex2_0.out 90 91 # Extrude the sphere with cutaway 92 test: 93 suffix: 9 94 requires: pyvista 95 args: -dm_view pyvista -view_pyvista_clip 1.,0.,0. -dm_plex_shape sphere -dm_refine 1 \ 96 -dm_plex_transform_type extrude -dm_plex_transform_extrude_layers 3 \ 97 -dm_plex_transform_extrude_use_tensor 0 -dm_plex_transform_extrude_thickness 0.3 98 output_file: output/ex2_0.out 99 100 # Extrude the refined sphere 101 test: 102 suffix: 10 103 requires: pyvista 104 args: -dm_view pyvista -view_pyvista_clip 1.,0.,0. -dm_plex_shape sphere -dm_refine 3 \ 105 -dm_plex_option_phases ext_ \ 106 -ext_dm_refine 1 -ext_dm_plex_transform_type extrude -ext_dm_plex_transform_extrude_layers 3 \ 107 -ext_dm_plex_transform_extrude_use_tensor 0 -ext_dm_plex_transform_extrude_thickness 0.5 108 output_file: output/ex2_0.out 109 110 # Extrude the refined sphere 111 test: 112 suffix: 11 113 requires: pyvista 114 args: -dm_view pyvista -view_pyvista_clip 1.,0.,0. -dm_plex_shape sphere -dm_refine 3 \ 115 -dm_plex_option_phases ext_,ref_ \ 116 -ext_dm_refine 1 -ext_dm_plex_transform_type extrude -ext_dm_plex_transform_extrude_layers 3 \ 117 -ext_dm_plex_transform_extrude_use_tensor 0 -ext_dm_plex_transform_extrude_thickness 0.5 \ 118 -ref_dm_refine 1 119 output_file: output/ex2_0.out 120 121 # Refine the Schwartz surface 122 test: 123 suffix: 12 124 requires: pyvista 125 args: -dm_view pyvista -dm_plex_shape schwarz_p -dm_plex_tps_extent 3,2 -dm_plex_tps_refine 2 126 output_file: output/ex2_0.out 127 128 # Refine and extrude the Schwartz surface with given thickness 129 test: 130 suffix: 13 131 requires: pyvista 132 args: -dm_view pyvista -dm_plex_shape schwarz_p -dm_plex_tps_extent 3,2,2 -dm_plex_tps_refine 3 \ 133 -dm_plex_tps_thickness 0.4 -dm_plex_tps_layers 3 134 output_file: output/ex2_0.out 135 136 # Filter the square 137 test: 138 suffix: 14 139 requires: triangle pyvista 140 args: -dm_view pyvista -dm_refine 1 \ 141 -dm_plex_transform_type transform_filter -dm_plex_transform_active filter_cells -dm_plex_label_filter_cells 0,1,2 142 output_file: output/ex2_0.out 143 144 # Filter a cap on the sphere 145 test: 146 suffix: 15 147 requires: pyvista 148 args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 4 \ 149 -dm_plex_option_phases filt_ \ 150 -filt_dm_refine 1 -filt_dm_plex_transform_type transform_filter \ 151 -filt_dm_plex_transform_active filter_cells -dm_plex_label_filter_cells 0 152 output_file: output/ex2_0.out 153 154 # Filter the sphere minus a cap 155 test: 156 suffix: 16 157 requires: pyvista 158 args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 4 \ 159 -dm_plex_option_phases filt_ \ 160 -filt_dm_refine 1 -filt_dm_plex_transform_type transform_filter \ 161 -filt_dm_plex_transform_active filter_cells \ 162 -dm_plex_label_filter_cells 0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18 163 output_file: output/ex2_0.out 164 165 # Convert sphere to quads 166 test: 167 suffix: 17 168 requires: pyvista 169 args: -dm_view pyvista -dm_plex_shape sphere -dm_refine 3 \ 170 -dm_plex_option_phases conv_ -conv_dm_refine 1 -conv_dm_plex_transform_type refine_tobox 171 output_file: output/ex2_0.out 172 173 # Load and refine the nozzle 174 test: 175 suffix: 18 176 requires: pyvista egads datafilespath 177 args: -dm_view pyvista -dm_plex_filename /Users/knepley/PETSc4/petsc/petsc-dev/share/petsc/datafiles/meshes/nozzle.stp -dm_refine 3 178 output_file: output/ex2_0.out 179 180 # Load and refine the nozzle, and convert to quads 181 test: 182 suffix: 19 183 requires: pyvista egads datafilespath 184 args: -dm_view pyvista -dm_plex_filename /Users/knepley/PETSc4/petsc/petsc-dev/share/petsc/datafiles/meshes/nozzle.stp -dm_refine 3 \ 185 -dm_plex_option_phases conv_ -conv_dm_refine 1 -conv_dm_plex_transform_type refine_tobox 186 output_file: output/ex2_0.out 187 188 TEST*/ 189