1 /// \file mbex4.cpp
2 ///
3 /// \author Milad Fatenejad
4 ///
5 /// \brief beginner tutorial, example 4: Create a 2D structured mesh
6 /// and set some tag data on it
7 ///
8 /// In this example, we create a 2D structured mesh (actually a 3D
9 /// mesh made of quads) and then we will actually set data on the
10 /// mesh. Tags represent data that is attached to entities. In this
11 /// example, we will create two tags:
12 ///
13 /// -# A "temperature" tag that is a single double precision number
14 /// attached to each quad.
15 /// -# A "velocity" tag that is an array of 2 double precision numbers
16 /// attached to each vertex.
17 ///
18 /// We will write the mesh out to a file, and you can visualize the
19 /// data using your favorite tool.
20 ///
21 /// In this example, I am demonstrating these operations in the
22 /// clearest possible way - not using the most efficient method. MOAB
23 /// has been designed so that users do not need to sacrifice
24 /// performance - future examples will demonstrate how to
25 /// access/manipulate the mesh using the fastest methods possible.
26
27 // The moab/Core.hpp header file is needed for all MOAB work...
28 #include "moab/Core.hpp"
29
30 // The moab/ScdInterface.hpp contains code which defines the moab
31 // structured mesh interface
32 #include "moab/ScdInterface.hpp"
33
34 #include <iostream>
35 #include <cmath>
36
37 // ****************
38 // * *
39 // * main *
40 // * *
41 // ****************
42 int main()
43 {
44 moab::ErrorCode rval;
45 moab::Core mbint;
46
47 // ***********************
48 // * Create the Mesh *
49 // ***********************
50
51 // First, lets make the mesh. It will be a 100 by 100 uniform grid
52 // (there will be 100x100 quads, 101x101 vertexes) with dx = dy =
53 // 0.1. Unlike the previous example, we will first make the mesh,
54 // then set the coordinates one at a time.
55
56 const unsigned NI = 100;
57 const unsigned NJ = 100;
58
59 // moab::ScdInterface is the structured mesh interface class for
60 // MOAB.
61 moab::ScdInterface* scdint;
62
63 // Tell MOAB that our mesh is structured:
64 rval = mbint.query_interface( scdint );MB_CHK_SET_ERR( rval, "mbint.query_interface failed" );
65
66 // Create the mesh:
67 moab::ScdBox* scdbox = NULL;
68 rval = scdint->construct_box( moab::HomCoord( 0, 0, 0 ), moab::HomCoord( NI, NJ, 0 ), NULL, 0, scdbox );MB_CHK_SET_ERR( rval, "scdint->construct_box failed" );
69
70 // MOAB knows to make quads instead of hexes because the last start
71 // and end indexes are the same (0). Note that it is still a "3D"
72 // mesh because each vertex coordinate is still defined using three
73 // numbers although every element in the mesh is a quadrilateral.
74
75 // ******************************
76 // * Set Vertex Coordinates *
77 // ******************************
78
79 // The "NULL" and "0" arguments in the call to construct_box are
80 // where we could specify the vertex coordinates. Since we didn't
81 // give any coordinates, every vertex is given a position of 0,0,0
82 // by default. Now we will set the vertex coordinates...
83
84 const double DX = 0.1;
85 const double DY = 0.1;
86
87 for( unsigned i = 0; i < NI + 1; i++ )
88 for( unsigned j = 0; j < NJ + 1; j++ )
89 {
90 // First, get the entity handle:
91 moab::EntityHandle handle = scdbox->get_vertex( i, j );
92
93 // Compute the coordinate:
94 double coord[3] = { DX * i, DY * j, 0.0 };
95
96 // Change the coordinate of the vertex:
97 mbint.set_coords( &handle, 1, coord );
98 }
99
100 // *******************
101 // * Attach Tags *
102 // *******************
103
104 // The vertex coordinates have been defined, now let's attach some
105 // data to the mesh. In MOAB this is done using "tags". Tags are
106 // little bits of information that can be attached to any mesh
107 // entity. In our example, we want to create two tags. The
108 // "temperature" tag will be attached to each quad and will be 1
109 // double. The "velocity" tag will be attached to each vertex and
110 // will be an array of two doubles.
111
112 // zero and twozeros represent the initial tag
113
114 // Create the tags:
115 moab::Tag temp_tag;
116 double temp_default_value = 0.0;
117 rval = mbint.tag_get_handle( "temperature", 1, moab::MB_TYPE_DOUBLE, temp_tag,
118 moab::MB_TAG_DENSE | moab::MB_TAG_CREAT, &temp_default_value );MB_CHK_SET_ERR( rval, "mbint.tag_get_handle(temperature) failed" );
119
120 moab::Tag vel_tag;
121 double vel_default_value[2] = { 0.0, 0.0 };
122 rval = mbint.tag_get_handle( "velocity", 2, moab::MB_TYPE_DOUBLE, vel_tag, moab::MB_TAG_DENSE | moab::MB_TAG_CREAT,
123 vel_default_value );MB_CHK_SET_ERR( rval, "mbint.tag_get_handle(velocity) failed" );
124
125 // Note that when we created each tag, we specified two flags:
126 //
127 // The moab::MB_TAG_DENSE flag tells MOAB that this is a dense
128 // tag. Dense tags will get automatically assigned to entities which
129 // have continuous handles. For this example, this means that once
130 // we set a tag on one vertex, memory will be allocated for
131 // assigning the tag to all vertexes. The same is true of
132 // quads. Dense tags are much more efficient when assigning tags to
133 // lots of entities. If you only want to assign a tag to a few
134 // entities, it is more efficient to use sparse tags
135 // (moab::MB_TAG_SPARSE).
136 //
137 // The moab::MB_TAG_CREAT flag tells MOAB to create the tag if it
138 // doesn't already exist.
139
140 // The tags have now been created, now we have to attach them to
141 // entities and set their values. NOTE: I am going to do this in a
142 // manner which emphasizes clarity - this is not the most efficient
143 // approach - that will be saved for later tutorial examples.
144
145 // Loop through each quad and set the temperature:
146 for( unsigned i = 0; i < NI; i++ )
147 for( unsigned j = 0; j < NJ; j++ )
148 {
149 // Get the handle for this quad:
150 moab::EntityHandle handle = scdbox->get_element( i, j );
151
152 // Compute the temperature...
153 double xc = DX * ( i + 0.5 );
154 double yc = DY * ( j + 0.5 );
155 double r = std::sqrt( xc * xc + yc * yc );
156 double temperature = std::exp( -0.5 * r );
157
158 // Set the temperature on a single quad:
159 rval = mbint.tag_set_data( temp_tag, &handle, 1, &temperature );MB_CHK_SET_ERR( rval, "mbint.tag_set_data(temp_tag) failed" );
160 }
161
162 // Loop through each vertex and set the velocity:
163 for( unsigned i = 0; i < NI + 1; i++ )
164 for( unsigned j = 0; j < NJ + 1; j++ )
165 {
166 // Get the handle for this vertex:
167 moab::EntityHandle handle = scdbox->get_vertex( i, j );
168 double velocity[2] = { i, j };
169
170 // Set the velocity on a vertex:
171 rval = mbint.tag_set_data( vel_tag, &handle, 1, velocity );MB_CHK_SET_ERR( rval, "mbint.tag_set_data(vel_tag) failed" );
172 }
173
174 // ***************************
175 // * Write Mesh to Files *
176 // ***************************
177
178 // NOTE: Some visualization software (such as VisIt) may not
179 // interpret the velocity tag as a vector and you may not be able to
180 // plot it. But you should be able to plot the temperature on top of
181 // the mesh.
182
183 rval = mbint.write_file( "mbex4.vtk" );MB_CHK_SET_ERR( rval, "write_file(mbex4.vtk) failed" );
184
185 return 0;
186 }