Actual source code: ex1.c
petsc-3.4.2 2013-07-02
1: static char help[] = "Run C version of TetGen to construct and refine a mesh\n\n";
3: #include <petscdmplex.h>
5: #if defined(PETSC_HAVE_CGNS)
6: #include <cgnslib.h>
7: #endif
9: typedef struct {
10: DM dm; /* REQUIRED in order to use SNES evaluation functions */
11: PetscInt debug; /* The debugging level */
12: PetscLogEvent createMeshEvent;
13: /* Domain and mesh definition */
14: PetscInt dim; /* The topological mesh dimension */
15: PetscBool interpolate; /* Generate intermediate mesh elements */
16: PetscBool refinementUniform; /* Uniformly refine the mesh */
17: PetscReal refinementLimit; /* The largest allowable cell volume */
18: PetscBool cellSimplex; /* Use simplices or hexes */
19: char filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */
20: } AppCtx;
24: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
25: {
29: options->debug = 0;
30: options->dim = 2;
31: options->interpolate = PETSC_FALSE;
32: options->refinementUniform = PETSC_FALSE;
33: options->refinementLimit = 0.0;
34: options->cellSimplex = PETSC_TRUE;
35: options->filename[0] = '\0';
37: PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");
38: PetscOptionsInt("-debug", "The debugging level", "ex1.c", options->debug, &options->debug, NULL);
39: PetscOptionsInt("-dim", "The topological mesh dimension", "ex1.c", options->dim, &options->dim, NULL);
40: PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex1.c", options->interpolate, &options->interpolate, NULL);
41: PetscOptionsBool("-refinement_uniform", "Uniformly refine the mesh", "ex1.c", options->refinementUniform, &options->refinementUniform, NULL);
42: PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex1.c", options->refinementLimit, &options->refinementLimit, NULL);
43: PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex1.c", options->cellSimplex, &options->cellSimplex, NULL);
44: PetscOptionsString("-filename", "The mesh file", "ex7.c", options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);
45: PetscOptionsEnd();
47: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
48: return(0);
49: };
53: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
54: {
55: PetscInt dim = user->dim;
56: PetscBool interpolate = user->interpolate;
57: PetscBool refinementUniform = user->refinementUniform;
58: PetscReal refinementLimit = user->refinementLimit;
59: PetscBool cellSimplex = user->cellSimplex;
60: const char *filename = user->filename;
61: const char *partitioner = "chaco";
62: size_t len;
63: PetscMPIInt rank;
67: PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
68: MPI_Comm_rank(comm, &rank);
69: PetscStrlen(filename, &len);
70: if (len) {
71: #if defined(PETSC_HAVE_CGNS)
72: int cgid = -1;
74: if (!rank) {
75: cg_open(filename, CG_MODE_READ, &cgid);
76: if (cgid <= 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_LIB, "cg_open(\"%s\",...) did not return a valid file ID", filename);
77: }
78: DMPlexCreateCGNS(comm, cgid, interpolate, dm);
79: if (!rank) {cg_close(cgid);}
80: #else
81: SETERRQ(comm, PETSC_ERR_SUP, "Loading meshes requires CGNS support. Reconfigure using --with-cgns-dir");
82: #endif
83: } else if (cellSimplex) {
84: DMPlexCreateBoxMesh(comm, dim, interpolate, dm);
85: } else {
86: const PetscInt cells[3] = {2, 2, 2};
88: DMPlexCreateHexBoxMesh(comm, dim, cells, dm);
89: }
90: {
91: DM refinedMesh = NULL;
92: DM distributedMesh = NULL;
94: /* Refine mesh using a volume constraint */
95: DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
96: DMPlexSetRefinementLimit(*dm, refinementLimit);
97: DMRefine(*dm, comm, &refinedMesh);
98: if (refinedMesh) {
99: DMDestroy(dm);
100: *dm = refinedMesh;
101: }
102: /* Distribute mesh over processes */
103: DMPlexDistribute(*dm, partitioner, 0, &distributedMesh);
104: if (distributedMesh) {
105: DMDestroy(dm);
106: *dm = distributedMesh;
107: }
108: if (refinementUniform) {
109: DMPlexSetRefinementUniform(*dm, refinementUniform);
110: DMRefine(*dm, comm, &refinedMesh);
111: if (refinedMesh) {
112: DMDestroy(dm);
113: *dm = refinedMesh;
114: }
115: }
116: }
117: PetscObjectSetName((PetscObject) *dm, "Simplical Mesh");
118: DMSetFromOptions(*dm);
119: PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
120: user->dm = *dm;
121: return(0);
122: }
126: int main(int argc, char **argv)
127: {
128: AppCtx user; /* user-defined work context */
131: PetscInitialize(&argc, &argv, NULL, help);
132: ProcessOptions(PETSC_COMM_WORLD, &user);
133: CreateMesh(PETSC_COMM_WORLD, &user, &user.dm);
134: DMDestroy(&user.dm);
135: PetscFinalize();
136: return 0;
137: }