Actual source code: ex1.cxx
1: static char help[] = "Simple MOAB example\n\n";
3: #include <petscdmmoab.h>
4: #include "moab/ScdInterface.hpp"
6: typedef struct {
7: DM dm; /* DM implementation using the MOAB interface */
8: PetscLogEvent createMeshEvent;
9: /* Domain and mesh definition */
10: PetscInt dim;
11: char filename[PETSC_MAX_PATH_LEN];
12: char tagname[PETSC_MAX_PATH_LEN];
13: } AppCtx;
15: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
16: {
17: PetscBool flg;
19: PetscStrcpy(options->filename, "");
20: PetscStrcpy(options->tagname, "petsc_tag");
21: options->dim = -1;
23: PetscOptionsBegin(comm, "", "MOAB example options", "DMMOAB");
24: PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex1.cxx", options->dim, &options->dim, NULL, PETSC_DECIDE, 3);
25: PetscOptionsString("-filename", "The file containing the mesh", "ex1.cxx", options->filename, options->filename, sizeof(options->filename), NULL);
26: PetscOptionsString("-tagname", "The tag name from which to create a vector", "ex1.cxx", options->tagname, options->tagname, sizeof(options->tagname), &flg);
27: PetscOptionsEnd();
29: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
30: return 0;
31: }
33: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
34: {
35: moab::Interface *iface = NULL;
36: moab::Tag tag = NULL;
37: moab::Tag ltog_tag = NULL;
38: moab::Range range;
39: PetscInt tagsize;
40: moab::ErrorCode merr;
42: PetscLogEventBegin(user->createMeshEvent, 0, 0, 0, 0);
43: DMMoabCreateMoab(comm, iface, <og_tag, &range, dm);
44: std::cout << "Created DMMoab using DMMoabCreateMoab." << std::endl;
45: DMMoabGetInterface(*dm, &iface);
47: // load file and get entities of requested or max dimension
48: if (strlen(user->filename) > 0) {
49: merr = iface->load_file(user->filename);
50: MBERRNM(merr);
51: std::cout << "Read mesh from file " << user->filename << std::endl;
52: } else {
53: // make a simple structured mesh
54: moab::ScdInterface *scdi;
55: merr = iface->query_interface(scdi);
57: moab::ScdBox *box;
58: merr = scdi->construct_box(moab::HomCoord(0, 0, 0), moab::HomCoord(5, 5, 5), NULL, 0, box);
59: MBERRNM(merr);
60: user->dim = 3;
61: merr = iface->set_dimension(user->dim);
62: MBERRNM(merr);
63: std::cout << "Created structured 5x5x5 mesh." << std::endl;
64: }
65: if (-1 == user->dim) {
66: moab::Range tmp_range;
67: merr = iface->get_entities_by_handle(0, tmp_range);
68: MBERRNM(merr);
69: if (tmp_range.empty()) MBERRNM(moab::MB_FAILURE);
70: user->dim = iface->dimension_from_handle(*tmp_range.rbegin());
71: }
72: merr = iface->get_entities_by_dimension(0, user->dim, range);
73: MBERRNM(merr);
74: DMMoabSetLocalVertices(*dm, &range);
76: // get the requested tag and create if necessary
77: std::cout << "Creating tag with name: " << user->tagname << ";\n";
78: merr = iface->tag_get_handle(user->tagname, 1, moab::MB_TYPE_DOUBLE, tag, moab::MB_TAG_CREAT | moab::MB_TAG_DENSE);
79: MBERRNM(merr);
80: {
81: // initialize new tag with gids
82: std::vector<double> tag_vals(range.size());
83: merr = iface->tag_get_data(ltog_tag, range, tag_vals.data());
84: MBERRNM(merr); // read them into ints
85: double *dval = tag_vals.data();
86: int *ival = reinterpret_cast<int *>(dval); // "stretch" them into doubles, from the end
87: for (int i = tag_vals.size() - 1; i >= 0; i--) dval[i] = ival[i];
88: merr = iface->tag_set_data(tag, range, tag_vals.data());
89: MBERRNM(merr); // write them into doubles
90: }
91: merr = iface->tag_get_length(tag, tagsize);
92: MBERRNM(merr);
94: DMSetUp(*dm);
96: // create the dmmoab and initialize its data
97: PetscObjectSetName((PetscObject)*dm, "MOAB mesh");
98: PetscLogEventEnd(user->createMeshEvent, 0, 0, 0, 0);
99: user->dm = *dm;
100: return 0;
101: }
103: int main(int argc, char **argv)
104: {
105: AppCtx user; /* user-defined work context */
106: moab::ErrorCode merr;
107: Vec vec;
108: moab::Interface *mbImpl = NULL;
109: moab::Tag datatag = NULL;
112: PetscInitialize(&argc, &argv, NULL, help);
113: ProcessOptions(PETSC_COMM_WORLD, &user);
115: CreateMesh(PETSC_COMM_WORLD, &user, &user.dm); /* create the MOAB dm and the mesh */
117: DMMoabGetInterface(user.dm, &mbImpl);
118: merr = mbImpl->tag_get_handle(user.tagname, datatag);
119: MBERRNM(merr);
120: DMMoabCreateVector(user.dm, datatag, NULL, PETSC_TRUE, PETSC_FALSE, &vec); /* create a vec from user-input tag */
122: std::cout << "Created VecMoab from existing tag." << std::endl;
123: VecDestroy(&vec);
124: std::cout << "Destroyed VecMoab." << std::endl;
125: DMDestroy(&user.dm);
126: std::cout << "Destroyed DMMoab." << std::endl;
127: PetscFinalize();
128: return 0;
129: }
131: /*TEST
133: build:
134: requires: moab
136: test:
138: TEST*/