xref: /petsc/src/dm/impls/plex/tutorials/ex17.c (revision 3f02e49b19195914bf17f317a25cb39636853415)
1 static const char help[] = "Test of CAD functionality";
2 
3 #include <petscdmplexegads.h>
4 
5 static PetscErrorCode ComputeVolume(DM dm)
6 {
7   DMLabel     bodyLabel, faceLabel, edgeLabel;
8   double      surface = 0., volume = 0., vol;
9   PetscInt    dim, pStart, pEnd, pid;
10   const char *name;
11 
12   PetscFunctionBeginUser;
13   PetscCall(DMGetDimension(dm, &dim));
14   if (dim < 2) PetscFunctionReturn(0);
15   PetscCall(DMGetCoordinatesLocalSetUp(dm));
16   PetscCall(DMGetLabel(dm, "EGADS Body ID", &bodyLabel));
17   PetscCall(DMGetLabel(dm, "EGADS Face ID", &faceLabel));
18   PetscCall(DMGetLabel(dm, "EGADS Edge ID", &edgeLabel));
19 
20   PetscCall(DMPlexGetHeightStratum(dm, 0, &pStart, &pEnd));
21   for (PetscInt p = pStart; p < pEnd; ++p) {
22     PetscCall(DMLabelGetValue(dim == 2 ? faceLabel : bodyLabel, p, &pid));
23     if (pid >= 0) {
24       PetscCall(DMPlexComputeCellGeometryFVM(dm, p, &vol, NULL, NULL));
25       volume += vol;
26     }
27   }
28   PetscCall(DMPlexGetHeightStratum(dm, 1, &pStart, &pEnd));
29   for (PetscInt p = pStart; p < pEnd; ++p) {
30     PetscCall(DMLabelGetValue(dim == 2 ? edgeLabel : faceLabel, p, &pid));
31     if (pid >= 0) {
32       PetscCall(DMPlexComputeCellGeometryFVM(dm, p, &vol, NULL, NULL));
33       surface += vol;
34     }
35   }
36 
37   PetscCall(PetscObjectGetName((PetscObject)dm, &name));
38   PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "DM %s: Surface Area = %.6e Volume = %.6e\n", name ? name : "", surface, volume));
39   PetscFunctionReturn(0);
40 }
41 
42 int main(int argc, char *argv[])
43 {
44   DM dm;
45 
46   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
47   PetscCall(DMCreate(PETSC_COMM_WORLD, &dm));
48   PetscCall(DMSetType(dm, DMPLEX));
49   PetscCall(DMSetFromOptions(dm));
50   PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
51 
52   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)dm, "ref1_"));
53   PetscCall(DMSetFromOptions(dm));
54   PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
55 
56   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)dm, "inf1_"));
57   PetscCall(DMPlexInflateToGeomModel(dm, PETSC_TRUE));
58   PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
59 
60   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)dm, "ref2_"));
61   PetscCall(DMSetFromOptions(dm));
62   PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
63 
64   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)dm, "inf2_"));
65   PetscCall(DMPlexInflateToGeomModel(dm, PETSC_TRUE));
66   PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
67 
68   PetscCall(ComputeVolume(dm));
69   PetscCall(DMDestroy(&dm));
70 
71   PetscCall(PetscFinalize());
72   return 0;
73 }
74 
75 /*TEST
76 
77   build:
78     requires: egads tetgen
79 
80   testset:
81     requires: datafilespath
82     args: -bd_dm_distribute 0 -bd_dm_plex_name "CAD Surface" -bd_dm_plex_check_all -dm_plex_geom_print_model \
83           -bd_dm_view -bd_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID" \
84           -bd_dm_generator tetgen -dm_plex_name "CAD Mesh" -dm_view hdf5:sphere_volume.h5 \
85           -ref1_dm_refine 1 -ref1_dm_view hdf5:sphere_volume.h5 \
86           -inf1_dm_view hdf5:sphere_volume_inflate1.h5 \
87           -ref2_dm_refine 1 -ref2_dm_view hdf5:sphere_volume.h5 \
88           -inf2_dm_view hdf5:sphere_volume_inflate2.h5
89 
90     test:
91       suffix: sphere_egadslite
92       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/sphere_example.egadslite
93 
94     test:
95       suffix: sphere_egadslite_tess
96       TODO: broken
97       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/sphere_example.egadslite \
98             -dm_plex_geom_tess_model 1
99 
100     test:
101       suffix: sphere_egads
102       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/sphere_example.egads
103 
104     test:
105       suffix: sphere_egads_tess
106       TODO: broken
107       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/sphere_example.egads \
108             -dm_plex_geom_tess_model 1
109 
110     test:
111       suffix: cylinder_egads
112       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/cylinder_example.egads
113 
114     test:
115       suffix: cylinder_igs
116       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/cylinder_example.igs
117 
118     test:
119       suffix: cylinder_igs_tess
120       TODO: broken
121       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/cylinder_example.igs \
122             -dm_plex_geom_tess_model 1 -ref1_dm_refine 0 -ref2_dm_refine 0 -bd_dm_plex_view_labels
123 
124     test:
125       suffix: nozzle_stp
126       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/nozzle_example.stp
127 
128     test:
129       suffix: nozzle_stp_tess
130       TODO: broken
131       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/nozzle_example.stp \
132             -dm_plex_geom_tess_model 1 -ref1_dm_refine 0 -ref2_dm_refine 0 -bd_dm_plex_view_labels
133 
134     test:
135       suffix: nozzle_igs
136       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/nozzle_example.igs
137 
138     test:
139       suffix: nozzle_igs_tess
140       TODO: broken
141       args: -dm_plex_boundary_filename ${DATAFILESPATH}/meshes/cad/nozzle_example.igs \
142             -dm_plex_geom_tess_model 1 -ref1_dm_refine 0 -ref2_dm_refine 0 -bd_dm_plex_view_labels
143 
144 TEST*/
145