Actual source code: ex37.c

  1: static const char help[] = "Test of CAD functionality";

  3: #include <petscdmplex.h>

  5: /* TODO
  6:   - Fix IGES
  7:   - Test tessellation using -dm_plex_egads_with_tess
  8: */

 10: typedef struct {
 11:   char      filename[PETSC_MAX_PATH_LEN];
 12:   PetscBool volumeMesh;
 13: } AppCtx;

 15: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 16: {
 18:   options->filename[0] = '\0';
 19:   options->volumeMesh  = PETSC_TRUE;

 21:   PetscOptionsBegin(comm, "", "EGADSPlex Problem Options", "EGADSLite");
 22:   PetscOptionsString("-filename", "The CAD file", "ex37.c", options->filename, options->filename, sizeof(options->filename), NULL);
 23:   PetscOptionsBool("-volume_mesh", "Create a volume mesh", "ex37.c", options->volumeMesh, &options->volumeMesh, NULL);
 24:   PetscOptionsEnd();
 25:   return 0;
 26: }

 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;

 37:   DMGetDimension(dm, &dim);
 38:   if (dim < 2) return 0;
 39:   DMGetLabel(dm, "EGADS Body ID", &bodyLabel);
 40:   DMGetLabel(dm, "EGADS Face ID", &faceLabel);
 41:   DMGetLabel(dm, "EGADS Edge ID", &edgeLabel);

 43:   DMPlexGetHeightStratum(dm, 0, &pStart, &pEnd);
 44:   for (p = pStart; p < pEnd; ++p) {
 45:     DMLabelGetValue(dim == 2 ? faceLabel : bodyLabel, p, &pid);
 46:     if (pid >= 0) {
 47:       DMPlexComputeCellGeometryFVM(dm, p, &vol, NULL, NULL);
 48:       volume += vol;
 49:     }
 50:   }
 51:   DMPlexGetHeightStratum(dm, 1, &pStart, &pEnd);
 52:   for (p = pStart; p < pEnd; ++p) {
 53:     DMLabelGetValue(dim == 2 ? edgeLabel : faceLabel, p, &pid);
 54:     if (pid >= 0) {
 55:       DMPlexComputeCellGeometryFVM(dm, p, &vol, NULL, NULL);
 56:       surface += vol;
 57:     }
 58:   }

 60:   PetscObjectGetName(obj, &name);
 61:   PetscPrintf(PetscObjectComm(obj), "DM %s: Surface Area = %.6e Volume = %.6e\n", name ? name : "", surface, volume);
 62:   return 0;
 63: }

 65: int main(int argc, char *argv[])
 66: {
 67:   DM     surface, dm;
 68:   AppCtx ctx;

 71:   PetscInitialize(&argc, &argv, NULL, help);
 72:   ProcessOptions(PETSC_COMM_WORLD, &ctx);
 73:   DMPlexCreateFromFile(PETSC_COMM_WORLD, ctx.filename, "ex37_plex", PETSC_TRUE, &surface);
 74:   PetscObjectSetName((PetscObject)surface, "CAD Surface");
 75:   PetscObjectSetOptionsPrefix((PetscObject)surface, "sur_");
 76:   DMSetFromOptions(surface);
 77:   DMViewFromOptions(surface, NULL, "-dm_view");
 78:   ComputeVolume(surface);

 80:   if (ctx.volumeMesh) {
 81:     DMPlexGenerate(surface, "tetgen", PETSC_TRUE, &dm);
 82:     PetscObjectSetName((PetscObject)dm, "CAD Mesh");
 83:     DMPlexSetRefinementUniform(dm, PETSC_TRUE);
 84:     DMViewFromOptions(dm, NULL, "-pre_dm_view");

 86:     DMPlexInflateToGeomModel(dm);
 87:     DMViewFromOptions(dm, NULL, "-inf_dm_view");

 89:     DMSetFromOptions(dm);
 90:     DMViewFromOptions(dm, NULL, "-dm_view");
 91:     ComputeVolume(dm);
 92:   }

 94:   DMDestroy(&dm);
 95:   DMDestroy(&surface);
 96:   PetscFinalize();
 97:   return 0;
 98: }

100: /*TEST

102:   build:
103:     requires: egads tetgen

105:   test:
106:     suffix: sphere_0
107:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egadslite -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_egads_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"

109:   test:
110:     suffix: sphere_egads
111:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egads -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_egads_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"

113:   test:
114:     suffix: sphere_iges
115:     requires: broken
116:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.igs -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_egads_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"

118:   test:
119:     suffix: sphere_step
120:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.stp -dm_refine 1 -sur_dm_view -dm_plex_check_all -dm_plex_egads_print_model -sur_dm_plex_view_labels "EGADS Body ID","EGADS Face ID","EGADS Edge ID"

122:   test:
123:     suffix: nozzle_0
124:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/nozzle.egadslite -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all

126:   test:
127:     suffix: nozzle_egads
128:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/nozzle.egads -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all

130:   test:
131:     suffix: nozzle_iges
132:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/nozzle.igs -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all

134:   test:
135:     suffix: nozzle_step
136:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/nozzle.stp -sur_dm_refine 1 -sur_dm_view -dm_plex_check_all

138: TEST*/