moab
|
00001 /* \example structuredmesh structuredmesh.cpp 00002 * \brief Show creation and query of structured mesh through MOAB's structured mesh interface. 00003 * This is a serial example showing creation and query of a 3D structured mesh. A single 00004 * rectangular brick of hex elements is created, then referenced by its ijk parameterization. 00005 * 1D and 2D examples could be made simply by changing the dimension parameter and the number 00006 * of ijk parameters passed into the MOAB functions. 00007 * 00008 * This example: 00009 * 0. Instantiate MOAB and get the structured mesh interface 00010 * 1. Creates a IxJxK structured mesh, which includes I*J*K vertices and (I-1)*(J-1)*(K-1) hexes. 00011 * 2. Get the vertices and hexes from moab and check their numbers against I*J*K and (I-1)*(J-1)*(K-1), resp. 00012 * 3. Loop over elements in 3 nested loops over i, j, k; for each (i,j,k): 00013 * 3a. Get the element corresponding to (i,j,k) 00014 * 3b. Get the connectivity of the element 00015 * 3c. Get the coordinates of the vertices comprising that element 00016 * 4. Release the structured mesh interface and destroy the MOAB instance 00017 * 00018 * To run: ./structuredmesh <dim> <N> 00019 * (default values so can run w/ no user interaction) 00020 */ 00021 00022 #include "moab/Core.hpp" 00023 #include "moab/ScdInterface.hpp" 00024 #include <iostream> 00025 #include <vector> 00026 00027 using namespace moab; 00028 00029 int main(int argc, char **argv) 00030 { 00031 argv[0] = argv[argc - argc]; // To remove the warnings about unused parameters 00032 int I, J, K; 00033 // progoptions? 00034 std::cout << "Enter I, J, K... " << std::endl; 00035 std::cin >> I >> J >> K; 00036 00037 // 0. Instantiate MOAB and get the structured mesh interface 00038 Interface *mb = new Core(); 00039 ScdInterface *scdiface; 00040 ErrorCode rval = mb->query_interface(scdiface); // get a ScdInterface object through moab instance 00041 if (MB_SUCCESS != rval) return rval; 00042 00043 // 1. Creates a IxJxK structured mesh, which includes I*J*K vertices and (I-1)*(J-1)*(K-1) hexes. 00044 ScdBox *box; 00045 rval = scdiface->construct_box(HomCoord(0, 0, 0), HomCoord(I-1, J-1, K-1), // low, high box corners in parametric space 00046 NULL, 0, // NULL coords vector and 0 coords (don't specify coords for now) 00047 box); // box is the structured box object providing the parametric 00048 // structured mesh interface for this rectangle of elements 00049 if (MB_SUCCESS != rval) return rval; 00050 00051 // 2. Get the vertices and hexes from moab and check their numbers against I*J*K and (I-1)*(J-1)*(K-1), resp. 00052 Range verts, hexes; 00053 rval = mb->get_entities_by_dimension(0, 0, verts); // first '0' specifies "root set", or entire MOAB instance, second the entity dimension being requested 00054 if (MB_SUCCESS != rval) return rval; 00055 rval = mb->get_entities_by_dimension(0, 3, hexes); 00056 if (MB_SUCCESS != rval) return rval; 00057 00058 if ((I-1)*(J-1)*(K-1) == (int) hexes.size() && I*J*K == (int) verts.size()) 00059 std::cout << "Created " << hexes.size() << " hexes and " << verts.size() << " vertices." << std::endl; 00060 else 00061 std::cout << "Created the wrong number of vertices or hexes!" << std::endl; 00062 00063 // 3. Loop over elements in 3 nested loops over i, j, k; for each (i,j,k): 00064 std::vector<double> coords(8*3); 00065 std::vector<EntityHandle> connect; 00066 for (int k = 0; k < K-1; k++) { 00067 for (int j = 0; j < J-1; j++) { 00068 for (int i = 0; i < I-1; i++) { 00069 // 3a. Get the element corresponding to (i,j,k) 00070 EntityHandle ehandle = box->get_element(i, j, k); 00071 if (0 == ehandle) return MB_FAILURE; 00072 // 3b. Get the connectivity of the element 00073 rval = mb->get_connectivity(&ehandle, 1, connect); // get the connectivity, in canonical order 00074 if (MB_SUCCESS != rval) return rval; 00075 // 3c. Get the coordinates of the vertices comprising that element 00076 rval = mb->get_coords(connect.data(), connect.size(), coords.data()); // get the coordinates of those vertices 00077 if (MB_SUCCESS != rval) return rval; 00078 } 00079 } 00080 } 00081 00082 // 4. Release the structured mesh interface and destroy the MOAB instance 00083 mb->release_interface(scdiface); // tell MOAB we're done with the ScdInterface 00084 delete mb; 00085 00086 return 0; 00087 }