Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
mbex2.cpp
Go to the documentation of this file.
1 /// \file mbex2.cpp
2 ///
3 /// \author Milad Fatenejad
4 ///
5 /// \brief beginner tutorial, example 2: Demonstrates loading a mesh
6 /// from a file, finding coordinate locations, connectivity
7 /// information...
8 ///
9 /// In this example, we read in the VTK file (mbex1.vtk) generated
10 /// in example 1, pick one of the hexahedrons, find out the vertexes
11 /// that define it (connectivity), and find the coordinates of those
12 /// vertexes. Finally, we will rotate the mesh a little and write it
13 /// out to a new file.
14 
15 // The moab/Core.hpp header file is needed for all MOAB work...
16 #include "moab/Core.hpp"
17 
18 #include <iostream>
19 #include <iomanip>
20 #include <cmath>
21 
22 // ****************
23 // * *
24 // * main *
25 // * *
26 // ****************
27 int main()
28 {
29  moab::ErrorCode rval;
30  moab::Core mbcore;
31  moab::Interface& mbint = mbcore;
32 
33  // First, lets load the file from the previous example. Note that
34  // you must compile and run the previous example to get mbex1.vtk!
35  rval = mbint.load_file( "mbex1.vtk" );MB_CHK_SET_ERR( rval, "load_file failed" );
36 
37  // *********************
38  // * Introspection *
39  // *********************
40 
41  // We can now access everything about the mesh through mbint. Lets
42  // pick an element and find some information about it. In this case,
43  // we know that all of the elements in the mesh are hexahedrons, so
44  // we will ask MOAB for a range containing the handle for all hexes:
45  moab::Range hex_range;
46  rval = mbint.get_entities_by_type( 0, moab::MBHEX, hex_range );MB_CHK_SET_ERR( rval, "get_entities_by_type(HEX) failed" );
47 
48  std::cout << "Loaded a mesh containing: " << hex_range << std::endl;
49 
50  // Let's analyze one of the hexes (in this case, I picked hex three):
51  moab::EntityHandle handle = hex_range[3];
52 
53  // Find out the eight vertexes that define this particular hex:
54  moab::Range connectivity;
55  rval = mbint.get_connectivity( &handle, 1, connectivity );MB_CHK_SET_ERR( rval, "get_connectivity failed" );
56 
57  // Connectivity should now contain the handles for the eight
58  // vertexes of interest. Lets print the handle and coordinate of
59  // each vertex. Note how we iterate through the range just like with
60  // STL containers...
61  double coord[3];
63 
64  std::cout << std::setw( 6 ) << "Handle" << std::setw( 10 ) << "X" << std::setw( 10 ) << "Y" << std::setw( 10 )
65  << "Z" << std::endl;
66 
67  for( iter = connectivity.begin(); iter != connectivity.end(); ++iter )
68  {
69  rval = mbint.get_coords( &( *iter ), 1, coord );MB_CHK_SET_ERR( rval, "get_coords" );
70 
71  // Print the entity handle followed by the x, y, z coordinate of
72  // the vertex:
73  std::cout << std::setw( 6 ) << *iter << std::setw( 10 ) << coord[0] << std::setw( 10 ) << coord[1]
74  << std::setw( 10 ) << coord[2] << std::endl;
75  }
76 
77  // ***********************
78  // * Rotate the Mesh *
79  // ***********************
80 
81  // Now, let's rotate the mesh about the z-axis by pi/4 radians. Do
82  // to this, we will iterate through all of the vertexes
83 
84  // Start by getting a range containing all of the vertex handles:
85  moab::Range vertex_range;
86  rval = mbint.get_entities_by_type( 0, moab::MBVERTEX, vertex_range );MB_CHK_SET_ERR( rval, "get_entities_by_type(VERTEX) failed" );
87 
88  // Get the coordinates of all of the vertexes:
89  std::vector< double > vertex_coords( 3 * vertex_range.size() );
90  rval = mbint.get_coords( vertex_range, vertex_coords.data() );MB_CHK_SET_ERR( rval, "get_coords" );
91 
92  unsigned count = 0;
93  const double PI = 3.14159265359;
94  const double ANGLE = PI / 4;
95  for( iter = vertex_range.begin(); iter != vertex_range.end(); ++iter )
96  {
97 
98  // Save the old coordinates:
99  double x = vertex_coords[count + 0];
100  double y = vertex_coords[count + 1];
101  // double z = vertex_coords[count+2];
102 
103  // Apply the rotation:
104  vertex_coords[count + 0] = x * std::cos( ANGLE ) - y * std::sin( ANGLE );
105  vertex_coords[count + 1] = x * std::sin( ANGLE ) + y * std::cos( ANGLE );
106 
107  count += 3;
108  }
109 
110  // Now, the vertex_coords vector contains all of the updated
111  // coordinates. Let's push these back into MOAB:
112  mbint.set_coords( vertex_range, vertex_coords.data() );MB_CHK_SET_ERR( rval, "set_coords" );
113 
114  // **************************
115  // * Write Mesh to File *
116  // **************************
117 
118  rval = mbint.write_file( "mbex2.vtk" );MB_CHK_SET_ERR( rval, "write_file(mbex2.vtk) failed" );
119 
120  // Now you can open your favorite visualization tool (VisIt or ParaView)
121  // and look at the original (mbex1.vtk) and rotated (mbex2.vtk) meshes.
122 
123  return 0;
124 }