1 static const char help[] = "Test of CAD functionality";
2
3 #include <petscdmplex.h>
4
5 /* TODO
6 - Fix IGES
7 - Test tessellation using -dm_plex_egads_with_tess
8 */
9
10 typedef struct {
11 char filename[PETSC_MAX_PATH_LEN];
12 PetscBool volumeMesh;
13 } AppCtx;
14
ProcessOptions(MPI_Comm comm,AppCtx * options)15 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
16 {
17 PetscFunctionBeginUser;
18 options->filename[0] = '\0';
19 options->volumeMesh = PETSC_TRUE;
20
21 PetscOptionsBegin(comm, "", "EGADSPlex Problem Options", "EGADSlite");
22 PetscCall(PetscOptionsString("-filename", "The CAD file", "ex37.c", options->filename, options->filename, sizeof(options->filename), NULL));
23 PetscCall(PetscOptionsBool("-volume_mesh", "Create a volume mesh", "ex37.c", options->volumeMesh, &options->volumeMesh, NULL));
24 PetscOptionsEnd();
25 PetscFunctionReturn(PETSC_SUCCESS);
26 }
27
ComputeVolume(DM dm)28 static PetscErrorCode ComputeVolume(DM dm)
29 {
30 PetscObject obj = (PetscObject)dm;
31 DMLabel bodyLabel, faceLabel, edgeLabel;
32 double surface = 0., volume = 0., vol;
33 PetscInt dim, pStart, pEnd, p, pid;
34 const char *name;
35
36 PetscFunctionBeginUser;
37 PetscCall(DMGetDimension(dm, &dim));
38 if (dim < 2) PetscFunctionReturn(PETSC_SUCCESS);
39 PetscCall(DMGetLabel(dm, "EGADS Body ID", &bodyLabel));
40 PetscCall(DMGetLabel(dm, "EGADS Face ID", &faceLabel));
41 PetscCall(DMGetLabel(dm, "EGADS Edge ID", &edgeLabel));
42
43 PetscCall(DMPlexGetHeightStratum(dm, 0, &pStart, &pEnd));
44 for (p = pStart; p < pEnd; ++p) {
45 PetscCall(DMLabelGetValue(dim == 2 ? faceLabel : bodyLabel, p, &pid));
46 if (pid >= 0) {
47 PetscCall(DMPlexComputeCellGeometryFVM(dm, p, &vol, NULL, NULL));
48 volume += vol;
49 }
50 }
51 PetscCall(DMPlexGetHeightStratum(dm, 1, &pStart, &pEnd));
52 for (p = pStart; p < pEnd; ++p) {
53 PetscCall(DMLabelGetValue(dim == 2 ? edgeLabel : faceLabel, p, &pid));
54 if (pid >= 0) {
55 PetscCall(DMPlexComputeCellGeometryFVM(dm, p, &vol, NULL, NULL));
56 surface += vol;
57 }
58 }
59
60 PetscCall(PetscObjectGetName(obj, &name));
61 PetscCall(PetscPrintf(PetscObjectComm(obj), "DM %s: Surface Area = %.6e Volume = %.6e\n", name ? name : "", surface, volume));
62 PetscFunctionReturn(PETSC_SUCCESS);
63 }
64
main(int argc,char * argv[])65 int main(int argc, char *argv[])
66 {
67 DM surface, dm;
68 AppCtx ctx;
69
70 PetscFunctionBeginUser;
71 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
72 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &ctx));
73 PetscCall(DMPlexCreateFromFile(PETSC_COMM_WORLD, ctx.filename, "ex37_plex", PETSC_TRUE, &surface));
74 PetscCall(PetscObjectSetName((PetscObject)surface, "CAD Surface"));
75 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)surface, "sur_"));
76 PetscCall(DMSetFromOptions(surface));
77 PetscCall(DMViewFromOptions(surface, NULL, "-dm_view"));
78 PetscCall(ComputeVolume(surface));
79
80 if (ctx.volumeMesh) {
81 PetscCall(DMPlexGenerate(surface, "tetgen", PETSC_TRUE, &dm));
82 PetscCall(PetscObjectSetName((PetscObject)dm, "CAD Mesh"));
83 PetscCall(DMPlexSetRefinementUniform(dm, PETSC_TRUE));
84 PetscCall(DMViewFromOptions(dm, NULL, "-pre_dm_view"));
85
86 PetscCall(DMPlexInflateToGeomModel(dm, PETSC_TRUE));
87 PetscCall(DMViewFromOptions(dm, NULL, "-inf_dm_view"));
88
89 PetscCall(DMSetFromOptions(dm));
90 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
91 PetscCall(ComputeVolume(dm));
92 PetscCall(DMDestroy(&dm));
93 }
94 PetscCall(DMDestroy(&surface));
95 PetscCall(PetscFinalize());
96 return 0;
97 }
98
99 /*TEST
100
101 build:
102 requires: egads tetgen datafilespath
103
104 test:
105 suffix: sphere_0
106 args: -filename ${DATAFILESPATH}/meshes/cad/sphere_example.egadslite -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_geom_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"
107
108 test:
109 suffix: sphere_egads
110 args: -filename ${DATAFILESPATH}/meshes/cad/sphere_example.egads -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_geom_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"
111
112 test:
113 suffix: sphere_iges
114 TODO: broken
115 args: -filename ${DATAFILESPATH}/meshes/cad/sphere_example.igs -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_geom_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"
116
117 test:
118 suffix: sphere_step
119 args: -filename ${DATAFILESPATH}/meshes/cad/sphere_example.stp -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_geom_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"
120
121 test:
122 suffix: nozzle_0
123 args: -filename ${DATAFILESPATH}/meshes/cad/nozzle.egadslite -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all
124
125 test:
126 suffix: nozzle_egads
127 args: -filename ${DATAFILESPATH}/meshes/cad/nozzle.egads -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all
128
129 test:
130 suffix: nozzle_iges
131 args: -filename ${DATAFILESPATH}/meshes/cad/nozzle.igs -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all
132
133 test:
134 suffix: nozzle_step
135 args: -filename ${DATAFILESPATH}/meshes/cad/nozzle.stp -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all
136
137 TEST*/
138