MOAB: Mesh Oriented datABase  (version 5.3.0)
cub_file_test.cpp
Go to the documentation of this file.
00001 #include "TestUtil.hpp"
00002 #include "moab/Core.hpp"
00003 #include "MBTagConventions.hpp"
00004 #include "moab/CN.hpp"
00005 #include "moab/Range.hpp"
00006 #include "moab/GeomTopoTool.hpp"
00007 #include <cmath>
00008 #include <algorithm>
00009 
00010 using namespace moab;
00011 
00012 /**\brief Input test file: test.cub
00013  * Cubit 10.2 file.
00014  * File contains:
00015  * Two merged 10x10x10 bricks sharing a single surface (surface 6).
00016  * - Each brick is meshed with 8 5x5x5 hexes.
00017  * - Each surface is meshed with 4 5x5 quads.
00018  * - Each curve is meshed with 2 5-unit edges.
00019  * A single block containing both bricks.
00020  * Two side sets
00021  * - sideset 1: surfaces 1 and 7
00022  * - sideset 2: surfaces 5 and 11
00023  * Two node sets:
00024  * - nodeset 1: surfaces 2 and 8
00025  * - nodeset 2: surfaces 3 and 9
00026  *
00027  * Surfaces:
00028  *                    2          8
00029  *                   /          /
00030  *          o----------o----------o
00031  *         /.      /  /.      /  /|
00032  *        / .  (5)/  / . (11)/  / |
00033  *       /  .    L  /  .    L  /  |
00034  *      o----------o----------o   |
00035  *      |   .      |   .      |(12)
00036  *   4--|-> o . . .|. .o. . . | . o
00037  *      |  . (1)   |  . (9)   |  /
00038  *      | .   ^    | .   ^    | /
00039  *      |.    |    |.    |    |/
00040  *      o----------o----------o
00041  *            |          |
00042  *            3          9
00043  *
00044  * Curves:
00045  *
00046  *     o----8-----o----20----o
00047  *    /.         /.         /|
00048  *  12 .       11 .       24 |
00049  *  /  7       /  5       /  17
00050  * o----2-----o----14----o   |
00051  * |   .      |   .      |   |
00052  * |   o . .6.|. .o. . 18| . o
00053  * 3  .       1  .      13  /
00054  * | 9        | 10       | 22
00055  * |.         |.         |/
00056  * o----4-----o----16----o
00057  *
00058  * Vertices:
00059  *
00060  *     8----------5----------13
00061  *    /.         /.         /|
00062  *   / .        / .        / |
00063  *  /  .       /  .       /  |
00064  * 3----------2----------10  |
00065  * |   .      |   .      |   |
00066  * |   7 . . .|. .6. . . | . 14
00067  * |  .       |  .       |  /
00068  * | .        | .        | /
00069  * |.         |.         |/
00070  * 4----------1----------9
00071  */
00072 
00073 /* Input test file: ho_test.cub
00074  *
00075  * File is expected to contain at least one block for every
00076  * supported higher-order element type.  The coordinates of
00077  * every higher-order node are expected to be the mean of the
00078  * adjacent corner vertices of the element.
00079  */
00080 static const std::string input_file_1 = std::string( TestDir + "/io/test.cub" );
00081 static const std::string ho_file      = std::string( TestDir + "/io/ho_test.cub" );
00082 static const std::string cubit12_file = std::string( TestDir + "/io/cubtest12.cub" );
00083 static const std::string cubit14_file = std::string( TestDir + "/io/cubtest14.cub" );
00084 
00085 void read_file( Interface& moab, const std::string& input_file );
00086 
00087 // Check that adjacent lower-order entities have
00088 // higher-order nodes consitent with input entity.
00089 void check_adj_ho_nodes( Interface& moab, EntityHandle entity );
00090 
00091 // Check that element has expected higher-order nodes
00092 // and that each higher-order node is at the center
00093 // of the sub-entity it is on.
00094 void check_ho_element( Interface& moab, EntityHandle entity, int mid_nodes[4] );
00095 
00096 // Validate elements of specified type.
00097 // Looks for a block containing the specified entity type
00098 // and with the specified mid-node flags set in its
00099 // HAS_MID_NODES_TAG.
00100 void test_ho_elements( EntityType type, int num_nodes );
00101 
00102 void test_vertices();
00103 
00104 void test_edges();
00105 
00106 void test_quads();
00107 
00108 void test_hexes();
00109 
00110 void test_geometric_topology();
00111 
00112 void test_geometric_sets();
00113 
00114 void test_blocks();
00115 
00116 void test_side_sets();
00117 
00118 void test_node_sets();
00119 
00120 void test_tri6()
00121 {
00122     test_ho_elements( MBTRI, 6 );
00123 }
00124 void test_tri7()
00125 {
00126     test_ho_elements( MBTRI, 7 );
00127 }
00128 
00129 void test_quad5()
00130 {
00131     test_ho_elements( MBQUAD, 5 );
00132 }
00133 void test_quad8()
00134 {
00135     test_ho_elements( MBQUAD, 8 );
00136 }
00137 void test_quad9()
00138 {
00139     test_ho_elements( MBQUAD, 9 );
00140 }
00141 
00142 void test_tet8()
00143 {
00144     test_ho_elements( MBTET, 8 );
00145 }
00146 void test_tet10()
00147 {
00148     test_ho_elements( MBTET, 10 );
00149 }
00150 void test_tet14()
00151 {
00152     test_ho_elements( MBTET, 14 );
00153 }
00154 
00155 void test_hex9()
00156 {
00157     test_ho_elements( MBHEX, 9 );
00158 }
00159 void test_hex20()
00160 {
00161     test_ho_elements( MBHEX, 20 );
00162 }
00163 void test_hex27()
00164 {
00165     test_ho_elements( MBHEX, 27 );
00166 }
00167 
00168 void test_multiple_files();
00169 
00170 void test_cubit12();
00171 void test_cubit14();
00172 
00173 int main()
00174 {
00175     int result = 0;
00176 
00177     result += RUN_TEST( test_vertices );
00178     result += RUN_TEST( test_edges );
00179     result += RUN_TEST( test_quads );
00180     result += RUN_TEST( test_hexes );
00181     result += RUN_TEST( test_geometric_topology );
00182     result += RUN_TEST( test_geometric_sets );
00183     result += RUN_TEST( test_blocks );
00184     result += RUN_TEST( test_side_sets );
00185     result += RUN_TEST( test_node_sets );
00186     result += RUN_TEST( test_tri6 );
00187     result += RUN_TEST( test_tri7 );
00188     result += RUN_TEST( test_quad5 );
00189     result += RUN_TEST( test_quad8 );
00190     result += RUN_TEST( test_quad9 );
00191     result += RUN_TEST( test_tet8 );
00192     result += RUN_TEST( test_tet10 );
00193     result += RUN_TEST( test_tet14 );
00194     result += RUN_TEST( test_hex9 );
00195     result += RUN_TEST( test_hex20 );
00196     result += RUN_TEST( test_hex27 );
00197     result += RUN_TEST( test_multiple_files );
00198     result += RUN_TEST( test_cubit12 );
00199     result += RUN_TEST( test_cubit14 );
00200     return result;
00201 }
00202 
00203 void read_file( Interface& moab, const std::string& input_file )
00204 {
00205     ErrorCode rval = moab.load_file( input_file.c_str() );CHECK_ERR( rval );
00206 }
00207 
00208 void test_vertices()
00209 {
00210     /* Node coordinates, in order of node ID, beginning with 1. */
00211     const double node_coords[] = {
00212         5,  -5, 5,                                                   // 1
00213         5,  5,  5,  5,  0,  5,  -5, 5,  5,  0,  5,  5,               // 5
00214         -5, -5, 5,  -5, 0,  5,  0,  -5, 5,  0,  0,  5,  5,  5,  -5,  // 10
00215         5,  -5, -5, 5,  0,  -5, -5, -5, -5, 0,  -5, -5, -5, 5,  -5,  // 15
00216         -5, 0,  -5, 0,  5,  -5, 0,  0,  -5, -5, -5, 0,  5,  -5, 0,   // 20
00217         0,  -5, 0,  -5, 5,  0,  -5, 0,  0,  5,  5,  0,  0,  5,  0,   // 25
00218         5,  0,  0,  0,  0,  0,  15, -5, 5,  15, 5,  5,  15, 0,  5,   // 30
00219         10, 5,  5,  10, -5, 5,  10, 0,  5,  15, 5,  -5, 15, -5, -5,  // 35
00220         15, 0,  -5, 10, -5, -5, 10, 5,  -5, 10, 0,  -5, 15, -5, 0,   // 40
00221         10, -5, 0,  15, 5,  0,  10, 5,  0,  15, 0,  0,  10, 0,  0    // 45
00222     };
00223 
00224     ErrorCode rval;
00225     Core mb_impl;
00226     Interface& mb = mb_impl;
00227     read_file( mb, input_file_1 );
00228 
00229     // get vertex handles and check correct number of vertices
00230     const size_t num_nodes = sizeof( node_coords ) / ( 3 * sizeof( double ) );
00231     Range verts;
00232     rval = mb.get_entities_by_type( 0, MBVERTEX, verts );CHECK_ERR( rval );
00233     CHECK_EQUAL( num_nodes, (size_t)verts.size() );
00234 
00235     // check global ids (should be 1 to 45 for vertices.)
00236     Tag gid_tag;
00237     rval = mb.tag_get_handle( "GLOBAL_ID", 1, MB_TYPE_INTEGER, gid_tag );CHECK_ERR( rval );
00238     std::vector< int > ids( num_nodes );
00239     rval = mb.tag_get_data( gid_tag, verts, &ids[0] );CHECK_ERR( rval );
00240     std::vector< int > sorted( ids );
00241     std::sort( sorted.begin(), sorted.end() );
00242     for( size_t i = 0; i < num_nodes; ++i )
00243         CHECK_EQUAL( (int)( i + 1 ), sorted[i] );
00244 
00245     // check coordinates of each vertex
00246     std::vector< double > coords( 3 * num_nodes );
00247     rval = mb.get_coords( verts, &coords[0] );CHECK_ERR( rval );
00248     for( size_t i = 0; i < num_nodes; ++i )
00249     {
00250         const double* exp = node_coords + 3 * ( ids[i] - 1 );
00251         const double* act = &coords[3 * i];
00252         CHECK_REAL_EQUAL( exp[0], act[0], 1e-8 );
00253         CHECK_REAL_EQUAL( exp[1], act[1], 1e-8 );
00254         CHECK_REAL_EQUAL( exp[2], act[2], 1e-8 );
00255     }
00256 }
00257 
00258 void test_element( const std::string& filename, EntityType type, int num_elem, int node_per_elem, const int* conn_list )
00259 {
00260     ErrorCode rval;
00261     Core mb_impl;
00262     Interface& mb = mb_impl;
00263     read_file( mb, filename );
00264 
00265     Range elems;
00266     rval = mb.get_entities_by_type( 0, type, elems );CHECK_ERR( rval );
00267     CHECK_EQUAL( num_elem, (int)elems.size() );
00268 
00269     // get global ids
00270     Tag gid_tag;
00271     rval = mb.tag_get_handle( "GLOBAL_ID", 1, MB_TYPE_INTEGER, gid_tag );CHECK_ERR( rval );
00272     std::vector< int > ids( num_elem );
00273     rval = mb.tag_get_data( gid_tag, elems, &ids[0] );CHECK_ERR( rval );
00274 
00275     // check that global ids are consecutive, beginning with 1
00276     std::vector< int > sorted( ids );
00277     std::sort( sorted.begin(), sorted.end() );
00278     for( int i = 0; i < num_elem; ++i )
00279         CHECK_EQUAL( i + 1, sorted[i] );
00280 
00281     // check connectivity of each element
00282     std::vector< int > conn_ids( node_per_elem );
00283     std::vector< EntityHandle > conn_h;
00284     Range::iterator j = elems.begin();
00285     for( int i = 0; i < num_elem; ++i, ++j )
00286     {
00287         conn_h.clear();
00288         rval = mb.get_connectivity( &*j, 1, conn_h );CHECK_ERR( rval );
00289         CHECK_EQUAL( node_per_elem, (int)conn_h.size() );
00290         rval = mb.tag_get_data( gid_tag, &conn_h[0], node_per_elem, &conn_ids[0] );CHECK_ERR( rval );
00291         const int* exp = conn_list + node_per_elem * ( ids[i] - 1 );
00292         for( int k = 0; k < node_per_elem; ++k )
00293             CHECK_EQUAL( exp[k], conn_ids[k] );
00294     }
00295 }
00296 
00297 void test_edges()
00298 {
00299     const int edge_conn[] = { 1,  3,                                                                           // 1
00300                               3,  2,  2,  5,  5,  4,  4,  7,  7,  6,  6,  8,  8,  1,  3,  9,  9,  8,           // 10
00301                               5,  9,  9,  7,  10, 12, 12, 11, 11, 14, 14, 13, 13, 16, 16, 15, 15, 17, 17, 10,  // 20
00302                               12, 18, 18, 17, 14, 18, 18, 16, 6,  19, 19, 13, 1,  20, 20, 11, 19, 21, 21, 8,   // 30
00303                               14, 21, 21, 20, 4,  22, 22, 15, 22, 23, 23, 7,  16, 23, 23, 19, 2,  24, 24, 10,  // 40
00304                               24, 25, 25, 5,  17, 25, 25, 22, 20, 26, 26, 3,  12, 26, 26, 24, 28, 30, 30, 29,  // 50
00305                               29, 31, 31, 2,  1,  32, 32, 28, 30, 33, 33, 32, 31, 33, 33, 3,  34, 36, 36, 35,  // 60
00306                               35, 37, 37, 11, 10, 38, 38, 34, 36, 39, 39, 38, 37, 39, 39, 12, 28, 40, 40, 35,  // 70
00307                               20, 41, 41, 32, 37, 41, 41, 40, 29, 42, 42, 34, 42, 43, 43, 31, 38, 43, 43, 24,  // 80
00308                               40, 44, 44, 30, 36, 44, 44, 42 };
00309     test_element( input_file_1, MBEDGE, 84, 2, edge_conn );
00310 }
00311 
00312 void test_quads()
00313 {
00314     const int quad_conn[] = {
00315         1,  3,  9,  8,  // 1
00316         3,  2,  5,  9,  8,  9,  7,  6,  9,  5,  4,  7,  10, 12, 18, 17, 12, 11, 14, 18,
00317         17, 18, 16, 15, 18, 14, 13, 16, 6,  19, 21, 8,  19, 13, 14, 21,  // 10
00318         8,  21, 20, 1,  21, 14, 11, 20, 4,  22, 23, 7,  22, 15, 16, 23, 7,  23, 19, 6,
00319         23, 16, 13, 19, 2,  24, 25, 5,  24, 10, 17, 25, 5,  25, 22, 4,  25, 17, 15, 22,  // 20
00320         1,  20, 26, 3,  20, 11, 12, 26, 3,  26, 24, 2,  26, 12, 10, 24, 28, 30, 33, 32,
00321         30, 29, 31, 33, 32, 33, 3,  1,  33, 31, 2,  3,  34, 36, 39, 38, 36, 35, 37, 39,  // 30
00322         38, 39, 12, 10, 39, 37, 11, 12, 1,  20, 41, 32, 20, 11, 37, 41, 32, 41, 40, 28,
00323         41, 37, 35, 40, 29, 42, 43, 31, 42, 34, 38, 43, 31, 43, 24, 2,  43, 38, 10, 24,  // 40
00324         28, 40, 44, 30, 40, 35, 36, 44, 30, 44, 42, 29, 44, 36, 34, 42,
00325     };
00326     test_element( input_file_1, MBQUAD, 44, 4, quad_conn );
00327 }
00328 
00329 void test_hexes()
00330 {
00331     const int hex_conn[] = { 6,  19, 23, 7,  8,  21, 27, 9,  19, 13, 16, 23, 21, 14, 18, 27, 7,  23, 22, 4,  9,  27,
00332                              25, 5,  23, 16, 15, 22, 27, 18, 17, 25, 8,  21, 27, 9,  1,  20, 26, 3,  21, 14, 18, 27,
00333                              20, 11, 12, 26, 9,  27, 25, 5,  3,  26, 24, 2,  27, 18, 17, 25, 26, 12, 10, 24, 1,  20,
00334                              26, 3,  32, 41, 45, 33, 20, 11, 12, 26, 41, 37, 39, 45, 3,  26, 24, 2,  33, 45, 43, 31,
00335                              26, 12, 10, 24, 45, 39, 38, 43, 32, 41, 45, 33, 28, 40, 44, 30, 41, 37, 39, 45, 40, 35,
00336                              36, 44, 33, 45, 43, 31, 30, 44, 42, 29, 45, 39, 38, 43, 44, 36, 34, 42 };
00337     test_element( input_file_1, MBHEX, 16, 8, hex_conn );
00338 }
00339 
00340 template < int L >
00341 std::vector< int > find_parents( const int parent_conn[][L], int num_parent, int id )
00342 {
00343     std::vector< int > results;
00344     for( int i = 0; i < num_parent; ++i )
00345     {
00346         for( int j = 0; j < L; ++j )
00347         {
00348             if( parent_conn[i][j] == id ) results.push_back( i + 1 );
00349         }
00350     }
00351     return results;
00352 }
00353 
00354 int check_geometric_set( Interface& moab, int dim, int id, const int* children, int num_children,
00355                          std::vector< int > parents )
00356 {
00357     ErrorCode rval;
00358     Tag gid_tag, dim_tag;
00359 
00360     rval = moab.tag_get_handle( "GLOBAL_ID", 1, MB_TYPE_INTEGER, gid_tag );CHECK_ERR( rval );
00361     rval = moab.tag_get_handle( "GEOM_DIMENSION", 1, MB_TYPE_INTEGER, dim_tag );CHECK_ERR( rval );
00362     void* tag_vals[] = { &dim, &id };
00363     Tag tags[]       = { dim_tag, gid_tag };
00364     Range ents;
00365     rval = moab.get_entities_by_type_and_tag( 0, MBENTITYSET, tags, tag_vals, 2, ents );CHECK_ERR( rval );
00366     CHECK_EQUAL( 1u, (unsigned)ents.size() );
00367 
00368     const EntityHandle geom = ents.front();
00369     std::vector< int > exp_rel, act_rel;
00370     std::vector< EntityHandle > rel;
00371 
00372     if( num_children )
00373     {
00374         exp_rel.resize( num_children );
00375         std::copy( children, children + num_children, exp_rel.begin() );
00376         std::sort( exp_rel.begin(), exp_rel.end() );
00377         rel.clear();
00378         rval = moab.get_child_meshsets( geom, rel );CHECK_ERR( rval );
00379         CHECK_EQUAL( num_children, (int)rel.size() );
00380         act_rel.resize( rel.size() );
00381         rval = moab.tag_get_data( gid_tag, &rel[0], rel.size(), &act_rel[0] );CHECK_ERR( rval );
00382         std::sort( act_rel.begin(), act_rel.end() );
00383         CHECK( exp_rel == act_rel );
00384     }
00385 
00386     if( !parents.empty() )
00387     {
00388         exp_rel = parents;
00389         std::sort( exp_rel.begin(), exp_rel.end() );
00390         rel.clear();
00391         rval = moab.get_parent_meshsets( geom, rel );CHECK_ERR( rval );
00392         CHECK_EQUAL( parents.size(), rel.size() );
00393         act_rel.resize( rel.size() );
00394         rval = moab.tag_get_data( gid_tag, &rel[0], rel.size(), &act_rel[0] );CHECK_ERR( rval );
00395         std::sort( act_rel.begin(), act_rel.end() );
00396         CHECK( exp_rel == act_rel );
00397     }
00398 
00399     return 0;
00400 }
00401 
00402 void test_geometric_topology()
00403 {
00404     Core mb_impl;
00405     Interface& mb = mb_impl;
00406     read_file( mb, input_file_1 );
00407     // expected geometric vertices, specified by global ID
00408     const int vertex_ids[] = { 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 13, 14 };
00409     // List of global IDs of surfacs in geometric volumes, indexed by ID-1
00410     const int volume_surfs[2][6] = { { 1, 2, 3, 4, 5, 6 }, { 7, 8, 9, 6, 11, 12 } };
00411     // List of global IDs of curves in geometric surfaces, indexed by ID-1
00412     // Curve IDs of zero indicates that corresponding surface doesn't exist.
00413     const int surf_curves[12][4] = { { 1, 2, 3, 4 },     { 5, 6, 7, 8 },    { 9, 6, 10, 4 },   { 11, 7, 9, 3 },
00414                                      { 12, 8, 11, 2 },   { 10, 5, 12, 1 },  { 13, 14, 1, 16 }, { 17, 18, 5, 20 },
00415                                      { 10, 18, 22, 16 }, { 0, 0, 0, 0 },  // no surf 10
00416                                      { 24, 20, 12, 14 }, { 22, 17, 24, 13 } };
00417     // List of global IDs of vertices in geometric curves, indexed by ID-1
00418     // Vertex IDs of zero indicates that corresponding curve doesn't exist.
00419     const int curve_verts[24][2] = { { 1, 2 },  { 2, 3 },   { 3, 4 },  { 4, 1 }, { 5, 6 },  // 5
00420                                      { 6, 7 },  { 7, 8 },   { 8, 5 },  { 4, 7 }, { 1, 6 },  // 10
00421                                      { 3, 8 },  { 2, 5 },   { 9, 10 },                      // 13
00422                                      { 10, 2 }, { 0, 0 },                                   // no curve 15
00423                                      { 1, 9 },  { 13, 14 }, { 14, 6 }, { 0, 0 },            // no curve 19
00424                                      { 5, 13 }, { 0, 0 },                                   // no curve 21
00425                                      { 9, 14 }, { 0, 0 },                                   // no curve 23
00426                                      { 10, 13 } };
00427 
00428     // check all vertices
00429     for( unsigned i = 0; i < ( sizeof( vertex_ids ) / sizeof( vertex_ids[0] ) ); ++i )
00430         check_geometric_set( mb, 0, vertex_ids[i], 0, 0, find_parents< 2 >( curve_verts, 24, vertex_ids[i] ) );
00431 
00432     // check all curves
00433     for( int i = 1; i <= 24; ++i )
00434         if( curve_verts[i - 1][0] )
00435             check_geometric_set( mb, 1, i, curve_verts[i - 1], 2, find_parents< 4 >( surf_curves, 12, i ) );
00436 
00437     // check all surfs
00438     for( int i = 1; i <= 12; ++i )
00439         if( surf_curves[i - 1][0] )
00440             check_geometric_set( mb, 2, i, surf_curves[i - 1], 4, find_parents< 6 >( volume_surfs, 2, i ) );
00441 
00442     // check all volumes
00443     std::vector< int > empty;
00444     for( int i = 1; i <= 2; ++i )
00445         check_geometric_set( mb, 3, i, volume_surfs[i - 1], 6, empty );
00446 }
00447 
00448 void test_geometric_sets()
00449 {
00450     ErrorCode rval;
00451     Core mb_impl;
00452     Interface& mb = mb_impl;
00453     read_file( mb, input_file_1 );
00454     Tag gid_tag, dim_tag;
00455     rval = mb.tag_get_handle( "GLOBAL_ID", 1, MB_TYPE_INTEGER, gid_tag );CHECK_ERR( rval );
00456     rval = mb.tag_get_handle( "GEOM_DIMENSION", 1, MB_TYPE_INTEGER, dim_tag );CHECK_ERR( rval );
00457 
00458     // verify mesh entity counts
00459     Range verts, curves, surfs, vols;
00460     int dim = 0;
00461     // Cppcheck warning (false positive): variable dim is assigned a value that is never used
00462     // Cppcheck warning (false positive): variable dim is reassigned a value before the old one has
00463     // been used
00464     const void* vals[] = { &dim };
00465     rval               = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &dim_tag, vals, 1, verts );CHECK_ERR( rval );
00466     dim  = 1;
00467     rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &dim_tag, vals, 1, curves );CHECK_ERR( rval );
00468     dim  = 2;
00469     rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &dim_tag, vals, 1, surfs );CHECK_ERR( rval );
00470     dim  = 3;
00471     rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &dim_tag, vals, 1, vols );CHECK_ERR( rval );
00472 
00473     CHECK_EQUAL( 12u, (unsigned)verts.size() );
00474     CHECK_EQUAL( 20u, (unsigned)curves.size() );
00475     CHECK_EQUAL( 11u, (unsigned)surfs.size() );
00476     CHECK_EQUAL( 2u, (unsigned)vols.size() );
00477 
00478     // check that each vertex has a single node, and that the
00479     // node is also contained in any parent curve
00480     Range ents;
00481     Range::iterator i;
00482     for( i = verts.begin(); i != verts.end(); ++i )
00483     {
00484         ents.clear();
00485         rval = mb.get_entities_by_handle( *i, ents );CHECK_ERR( rval );
00486         CHECK_EQUAL( 1u, (unsigned)ents.size() );
00487         CHECK( ents.all_of_type( MBVERTEX ) );
00488     }
00489 
00490     // check that each curve has one node and two edges
00491     for( i = curves.begin(); i != curves.end(); ++i )
00492     {
00493         ents.clear();
00494         rval = mb.get_entities_by_handle( *i, ents );CHECK_ERR( rval );
00495         CHECK_EQUAL( 1u, (unsigned)ents.num_of_type( MBVERTEX ) );
00496         CHECK_EQUAL( 2u, (unsigned)ents.num_of_type( MBEDGE ) );
00497         CHECK_EQUAL( 3u, (unsigned)ents.size() );
00498     }
00499 
00500     // check that each surface has 1 node, 4 edges, 4 quads
00501     for( i = surfs.begin(); i != surfs.end(); ++i )
00502     {
00503         ents.clear();
00504         rval = mb.get_entities_by_handle( *i, ents );CHECK_ERR( rval );
00505         CHECK_EQUAL( 1u, (unsigned)ents.num_of_type( MBVERTEX ) );
00506         CHECK_EQUAL( 4u, (unsigned)ents.num_of_type( MBEDGE ) );
00507         CHECK_EQUAL( 4u, (unsigned)ents.num_of_type( MBQUAD ) );
00508         CHECK_EQUAL( 9u, (unsigned)ents.size() );
00509     }
00510 
00511     // check that each volume has 1 node and 8 hexes.
00512     for( i = vols.begin(); i != vols.end(); ++i )
00513     {
00514         ents.clear();
00515         rval = mb.get_entities_by_handle( *i, ents );CHECK_ERR( rval );
00516         CHECK_EQUAL( 1u, (unsigned)ents.num_of_type( MBVERTEX ) );
00517         CHECK_EQUAL( 8u, (unsigned)ents.num_of_type( MBHEX ) );
00518         CHECK_EQUAL( 9u, (unsigned)ents.size() );
00519     }
00520 
00521     // Check that for each geometric entity, any contained vertices
00522     // are adjacent to some entity in one of its parents.
00523     Range parents, geom, nodes, tmp;
00524     for( int d = 0; d < 3; ++d )
00525     {
00526         const void* vals1[] = { &d };
00527         rval                = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &dim_tag, vals1, 1, geom );CHECK_ERR( rval );
00528 
00529         for( i = geom.begin(); i != geom.end(); ++i )
00530         {
00531             nodes.clear();
00532             ents.clear();
00533             parents.clear();
00534             rval = mb.get_entities_by_type( *i, MBVERTEX, nodes );CHECK_ERR( rval );
00535             rval = mb.get_parent_meshsets( *i, parents );CHECK_ERR( rval );
00536             for( Range::iterator j = parents.begin(); j != parents.end(); ++j )
00537             {
00538                 tmp.clear();
00539                 rval = mb.get_entities_by_dimension( *j, d + 1, tmp );CHECK_ERR( rval );
00540                 ents.merge( tmp );
00541             }
00542             tmp.clear();
00543             rval = mb.get_adjacencies( ents, 0, false, tmp, Interface::UNION );CHECK_ERR( rval );
00544             nodes = subtract( nodes, tmp );
00545             CHECK( nodes.empty() );
00546         }
00547     }
00548 }
00549 
00550 // expect one block containing entire mesh, with id == 1
00551 void test_blocks()
00552 {
00553     ErrorCode rval;
00554     Core mb_impl;
00555     Interface& mb = mb_impl;
00556     read_file( mb, input_file_1 );
00557     Tag mat_tag;
00558     rval = mb.tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mat_tag );CHECK_ERR( rval );
00559 
00560     Range blocks;
00561     rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &mat_tag, 0, 1, blocks );CHECK_ERR( rval );
00562     CHECK_EQUAL( 1u, (unsigned)blocks.size() );
00563     EntityHandle block = blocks.front();
00564     int id;
00565     rval = mb.tag_get_data( mat_tag, &block, 1, &id );CHECK_ERR( rval );
00566     CHECK_EQUAL( 1, id );
00567 
00568     Range block_hexes, mesh_hexes;
00569     rval = mb.get_entities_by_dimension( 0, 3, mesh_hexes );CHECK_ERR( rval );
00570     rval = mb.get_entities_by_dimension( block, 3, block_hexes, true );CHECK_ERR( rval );
00571     CHECK( mesh_hexes == block_hexes );
00572 }
00573 
00574 // Common code for test_side_sets and test_node sets
00575 //\param tag_name  NEUMANN_SET_TAG_NAME or DIRICHLET_SET_TAG_NAME
00576 //\param count     Number of expected sets
00577 //\param ids       Expected IDs of sets
00578 //\param set_surfs One list for each id in "ids" containing the
00579 //                 ids of the geometric surfaces expected to be
00580 //                 contained in the boundary condition set.
00581 void test_bc_sets( const char* tag_name, unsigned count, const int* ids, const std::vector< int > set_surfs[] )
00582 {
00583     ErrorCode rval;
00584     Core mb_impl;
00585     Interface& mb = mb_impl;
00586     read_file( mb, input_file_1 );
00587     Tag ss_tag, gid_tag, dim_tag;
00588     rval = mb.tag_get_handle( tag_name, 1, MB_TYPE_INTEGER, ss_tag );CHECK_ERR( rval );
00589     rval = mb.tag_get_handle( "GLOBAL_ID", 1, MB_TYPE_INTEGER, gid_tag );CHECK_ERR( rval );
00590     rval = mb.tag_get_handle( "GEOM_DIMENSION", 1, MB_TYPE_INTEGER, dim_tag );CHECK_ERR( rval );
00591 
00592     // check number of sidesets and IDs
00593     Range sidesets;
00594     rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &ss_tag, 0, 1, sidesets );CHECK_ERR( rval );
00595     CHECK_EQUAL( count, (unsigned)sidesets.size() );
00596     std::vector< EntityHandle > handles( count, 0 );
00597     for( Range::iterator i = sidesets.begin(); i != sidesets.end(); ++i )
00598     {
00599         int id;
00600         rval = mb.tag_get_data( ss_tag, &*i, 1, &id );CHECK_ERR( rval );
00601         unsigned idx;
00602         for( idx = 0; idx < count; ++idx )
00603             if( ids[idx] == id ) break;
00604         CHECK( idx != count );
00605         CHECK( handles[idx] == 0 );
00606         handles[idx] = *i;
00607     }
00608 
00609     // get surface faces
00610     std::vector< Range > exp( count );
00611     Range surfs, tmp;
00612     Tag tags[] = { dim_tag, gid_tag };
00613     for( unsigned i = 0; i < count; ++i )
00614     {
00615         exp[i].clear();
00616         surfs.clear();
00617         const int two = 2;
00618         for( unsigned j = 0; j < set_surfs[i].size(); ++j )
00619         {
00620             const void* vals[] = { &two, &set_surfs[i][j] };
00621             surfs.clear();
00622             rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, tags, vals, 2, surfs );CHECK_ERR( rval );
00623             CHECK_EQUAL( 1u, (unsigned)surfs.size() );
00624             tmp.clear();
00625             rval = mb.get_entities_by_dimension( surfs.front(), 2, tmp, true );CHECK_ERR( rval );
00626             exp[i].merge( tmp );
00627         }
00628     }
00629 
00630     // check each bc set
00631     Range act;
00632     for( unsigned i = 0; i < count; ++i )
00633     {
00634         act.clear();
00635         rval = mb.get_entities_by_dimension( handles[i], 2, act, true );CHECK_ERR( rval );
00636         CHECK( exp[i] == act );
00637     }
00638 }
00639 
00640 // expect two sidesets containing two geometric surfaces each:
00641 // sideset 1 : surfaces 1 and 7
00642 // sideset 2 : surfaces 5 and 11
00643 void test_side_sets()
00644 {
00645     int ids[] = { 1, 2 };
00646     std::vector< int > surfs[2];
00647     surfs[0].push_back( 1 );
00648     surfs[0].push_back( 7 );
00649     surfs[1].push_back( 5 );
00650     surfs[1].push_back( 11 );
00651     test_bc_sets( NEUMANN_SET_TAG_NAME, 2, ids, surfs );
00652 }
00653 
00654 void test_node_sets()
00655 {
00656     int ids[] = { 1, 2 };
00657     std::vector< int > surfs[2];
00658     surfs[0].push_back( 2 );
00659     surfs[0].push_back( 8 );
00660     surfs[1].push_back( 3 );
00661     surfs[1].push_back( 9 );
00662     test_bc_sets( DIRICHLET_SET_TAG_NAME, 2, ids, surfs );
00663 }
00664 
00665 static EntityHandle find_side( Interface& moab, EntityHandle entity, int side_dim, int side_num )
00666 {
00667     ErrorCode rval;
00668 
00669     std::vector< EntityHandle > adj;
00670     rval = moab.get_adjacencies( &entity, 1, side_dim, false, adj );CHECK_ERR( rval );
00671 
00672     int sub_ent_indices[4];
00673     CN::SubEntityVertexIndices( TYPE_FROM_HANDLE( entity ), side_dim, side_num, sub_ent_indices );
00674     EntityType subtype  = CN::SubEntityType( TYPE_FROM_HANDLE( entity ), side_dim, side_num );
00675     int sub_ent_corners = CN::VerticesPerEntity( subtype );
00676 
00677     const EntityHandle* conn;
00678     int conn_len;
00679     rval = moab.get_connectivity( entity, conn, conn_len );CHECK_ERR( rval );
00680 
00681     for( size_t i = 0; i < adj.size(); ++i )
00682     {
00683         if( TYPE_FROM_HANDLE( adj[i] ) != subtype ) continue;
00684 
00685         const EntityHandle* sub_conn;
00686         int sub_len;
00687         rval = moab.get_connectivity( adj[i], sub_conn, sub_len );CHECK_ERR( rval );
00688 
00689         int n = std::find( sub_conn, sub_conn + sub_len, conn[sub_ent_indices[0]] ) - sub_conn;
00690         if( n == sub_len )  // no vertex in common
00691             continue;
00692 
00693         // check forward direction
00694         int j;
00695         for( j = 1; j < sub_ent_corners; ++j )
00696             if( conn[sub_ent_indices[j]] != sub_conn[( j + n ) % sub_ent_corners] ) break;
00697         if( j == sub_ent_corners ) return adj[i];
00698 
00699         // check reverse direction
00700         for( j = 1; j < sub_ent_corners; ++j )
00701             if( conn[sub_ent_indices[j]] != sub_conn[( n + sub_ent_corners - j ) % sub_ent_corners] ) break;
00702         if( j == sub_ent_corners ) return adj[i];
00703     }
00704 
00705     // no match
00706     return 0;
00707 }
00708 
00709 void check_adj_ho_nodes( Interface& moab, EntityHandle entity )
00710 {
00711     EntityType type = TYPE_FROM_HANDLE( entity );
00712     const EntityHandle* conn;
00713     int conn_len;
00714     ErrorCode rval = moab.get_connectivity( entity, conn, conn_len );CHECK_ERR( rval );
00715 
00716     int ho[4];
00717     CN::HasMidNodes( type, conn_len, ho );
00718     for( int dim = CN::Dimension( type ) - 1; dim > 0; --dim )
00719     {
00720         if( !ho[dim] ) continue;
00721 
00722         for( int j = 0; j < CN::NumSubEntities( type, dim ); ++j )
00723         {
00724             EntityHandle side = find_side( moab, entity, dim, j );
00725             if( !side ) continue;
00726 
00727             const EntityHandle* side_conn;
00728             int side_len;
00729             rval = moab.get_connectivity( side, side_conn, side_len );CHECK_ERR( rval );
00730 
00731             int this_idx = CN::HONodeIndex( type, conn_len, dim, j );
00732             int side_idx = CN::HONodeIndex( TYPE_FROM_HANDLE( side ), side_len, dim, 0 );
00733             CHECK_EQUAL( side_conn[side_idx], conn[this_idx] );
00734         }
00735     }
00736 }
00737 
00738 // Check that element has expected higher-order nodes
00739 // and that each higher-order node is at the center
00740 // of the sub-entity it is on.
00741 void check_ho_element( Interface& moab, EntityHandle entity, int mid_nodes[4] )
00742 {
00743     // get element info
00744     const EntityType type = TYPE_FROM_HANDLE( entity );
00745     const EntityHandle* conn;
00746     int conn_len;
00747     ErrorCode rval = moab.get_connectivity( entity, conn, conn_len );CHECK_ERR( rval );
00748     std::vector< double > coords( 3 * conn_len );
00749     rval = moab.get_coords( conn, conn_len, &coords[0] );CHECK_ERR( rval );
00750 
00751     // calculate and verify expected number of mid nodes
00752     int num_nodes = CN::VerticesPerEntity( type );
00753     for( int d = 1; d <= CN::Dimension( type ); ++d )
00754         if( mid_nodes[d] ) num_nodes += CN::NumSubEntities( type, d );
00755     CHECK_EQUAL( num_nodes, conn_len );
00756 
00757     // verify that each higher-order node is at the center
00758     // of its respective sub-entity.
00759     for( int i = CN::VerticesPerEntity( type ); i < num_nodes; ++i )
00760     {
00761         // get sub-entity owning ho-node
00762         int sub_dim, sub_num;
00763         CN::HONodeParent( type, num_nodes, i, sub_dim, sub_num );
00764         // get corner vertex indices
00765         int sub_conn[8], num_sub;
00766         if( sub_dim < CN::Dimension( type ) )
00767         {
00768             CN::SubEntityVertexIndices( type, sub_dim, sub_num, sub_conn );
00769             EntityType sub_type = CN::SubEntityType( type, sub_dim, sub_num );
00770             num_sub             = CN::VerticesPerEntity( sub_type );
00771         }
00772         else
00773         {
00774             num_sub = CN::VerticesPerEntity( type );
00775             for( int j = 0; j < num_sub; ++j )
00776                 sub_conn[j] = j;
00777         }
00778         // calculate mean of corner vertices
00779         double mean[3] = { 0, 0, 0 };
00780         for( int j = 0; j < num_sub; ++j )
00781         {
00782             int co = 3 * sub_conn[j];
00783             mean[0] += coords[co];
00784             mean[1] += coords[co + 1];
00785             mean[2] += coords[co + 2];
00786         }
00787         mean[0] /= num_sub;
00788         mean[1] /= num_sub;
00789         mean[2] /= num_sub;
00790         // verify that higher-order node is at expected location
00791         CHECK_REAL_EQUAL( mean[0], coords[3 * i], 1e-6 );
00792         CHECK_REAL_EQUAL( mean[1], coords[3 * i + 1], 1e-6 );
00793         CHECK_REAL_EQUAL( mean[2], coords[3 * i + 2], 1e-6 );
00794     }
00795 }
00796 
00797 // Validate elements of specified type.
00798 // Looks for a block containing the specified entity type
00799 // and with the specified mid-node flags set in its
00800 // HAS_MID_NODES_TAG.
00801 void test_ho_elements( EntityType type, int num_nodes )
00802 {
00803     Core mb_impl;
00804     Interface& mb = mb_impl;
00805     read_file( mb, ho_file );
00806 
00807     ErrorCode rval;
00808     Tag ho_tag, block_tag;
00809     rval = mb.tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, block_tag );CHECK_ERR( rval );
00810     rval = mb.tag_get_handle( HAS_MID_NODES_TAG_NAME, 4, MB_TYPE_INTEGER, ho_tag );CHECK_ERR( rval );
00811 
00812     // get material sets with expected higher-order nodes
00813     Range blocks;
00814     int ho_flags[4];
00815     CN::HasMidNodes( type, num_nodes, ho_flags );
00816     Tag tags[2]   = { ho_tag, block_tag };
00817     void* vals[2] = { ho_flags, NULL };
00818     rval          = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, tags, vals, 2, blocks );CHECK_ERR( rval );
00819 
00820     Range::iterator i;
00821     Range entities;
00822     for( i = blocks.begin(); i != blocks.end(); ++i )
00823     {
00824         rval = mb.get_entities_by_type( *i, type, entities, true );CHECK_ERR( rval );
00825     }
00826     // test file should contain all types of HO entities
00827     CHECK( !entities.empty() );
00828     // test ho node positions -- expect to be at center of adj corners
00829     for( i = entities.begin(); i != entities.end(); ++i )
00830         check_ho_element( mb, *i, ho_flags );
00831     // test ho node handles consistent with adjacent entities
00832     for( i = entities.begin(); i != entities.end(); ++i )
00833         check_adj_ho_nodes( mb, *i );
00834 }
00835 
00836 void test_multiple_files()
00837 {
00838     // load two surface meshes, one at z=+5 at z=-5.
00839     ErrorCode rval;
00840     Core mb_impl;
00841     Interface& mb = mb_impl;
00842     Range file1_ents, file2_ents;
00843     read_file( mb, input_file_1 );
00844     mb.get_entities_by_handle( 0, file1_ents );
00845     read_file( mb, input_file_1 );
00846     mb.get_entities_by_handle( 0, file2_ents );
00847     file2_ents = subtract( file2_ents, file1_ents );
00848     EntityHandle file1, file2;
00849     mb.create_meshset( MESHSET_SET, file1 );
00850     mb.create_meshset( MESHSET_SET, file2 );
00851     mb.add_entities( file1, file1_ents );
00852     mb.add_entities( file2, file2_ents );
00853 
00854     // first check that we get the same number of verts from
00855     // each file and that they are distinct vertices
00856     Range file1_verts, file2_verts;
00857     rval = mb.get_entities_by_type( file1, MBVERTEX, file1_verts );CHECK_ERR( rval );
00858     CHECK( !file1_verts.empty() );
00859     rval = mb.get_entities_by_type( file2, MBVERTEX, file2_verts );CHECK_ERR( rval );
00860     CHECK( !file2_verts.empty() );
00861     CHECK_EQUAL( file1_verts.size(), file2_verts.size() );
00862     CHECK( intersect( file1_verts, file2_verts ).empty() );
00863 
00864     // now check that we get the same number of elements from
00865     // each file and that they are distinct
00866     Range file1_elems, file2_elems;
00867     rval = mb.get_entities_by_dimension( file1, 3, file1_elems );CHECK_ERR( rval );
00868     CHECK( !file1_elems.empty() );
00869     rval = mb.get_entities_by_dimension( file2, 3, file2_elems );CHECK_ERR( rval );
00870     CHECK( !file2_elems.empty() );
00871     CHECK_EQUAL( file1_elems.size(), file2_elems.size() );
00872     CHECK( intersect( file1_elems, file2_elems ).empty() );
00873 
00874     // now check that the connectivity for each element is
00875     // defined using the appropriate vertex instances
00876     Range file1_elem_verts, file2_elem_verts;
00877     rval = mb.get_adjacencies( file1_elems, 0, false, file1_elem_verts, Interface::UNION );CHECK_ERR( rval );
00878     rval = mb.get_adjacencies( file2_elems, 0, false, file2_elem_verts, Interface::UNION );CHECK_ERR( rval );
00879     CHECK_EQUAL( file1_elem_verts.size(), file2_elem_verts.size() );
00880     CHECK( intersect( file1_elem_verts, file1_verts ) == file1_elem_verts );
00881     CHECK( intersect( file2_elem_verts, file2_verts ) == file2_elem_verts );
00882 }
00883 
00884 void test_cubit12()
00885 {
00886     Core mb_impl;
00887     Interface& mb = mb_impl;
00888     read_file( mb, cubit12_file );
00889 }
00890 
00891 void test_cubit14()
00892 {
00893     Core mb_impl;
00894     Interface& mb = mb_impl;
00895     read_file( mb, cubit14_file );
00896     // check the global id for some geometry sets
00897     GeomTopoTool gtt( &mb_impl );
00898     Range ranges[5];
00899     ErrorCode rval = gtt.find_geomsets( ranges );CHECK_ERR( rval );
00900     EntityHandle set0 = ranges[0][0];  // does it have a global id > 0?
00901     Tag gid_tag;
00902     rval = mb.tag_get_handle( "GLOBAL_ID", 1, MB_TYPE_INTEGER, gid_tag );CHECK_ERR( rval );
00903 
00904     int val;
00905     rval = mb.tag_get_data( gid_tag, &set0, 1, &val );
00906     CHECK( val != 0 );
00907 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines