Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
SphereDecomp.cpp
Go to the documentation of this file.
00001 #include "SphereDecomp.hpp"
00002 #include "moab/MeshTopoUtil.hpp"
00003 #include "moab/Range.hpp"
00004 #include "moab/CN.hpp"
00005 #include <cmath>
00006 #include <cassert>
00007 #include <iostream>
00008 
00009 #define RR \
00010     if( MB_SUCCESS != result ) return result
00011 
00012 const char* SUBDIV_VERTICES_TAG_NAME = "subdiv_vertices";
00013 
00014 using namespace moab;
00015 
00016 SphereDecomp::SphereDecomp( Interface* impl )
00017 {
00018     mbImpl = impl;
00019 }
00020 
00021 ErrorCode SphereDecomp::build_sphere_mesh( const char* sphere_radii_tag_name, EntityHandle* hex_set )
00022 {
00023     ErrorCode result = mbImpl->tag_get_handle( sphere_radii_tag_name, 1, MB_TYPE_DOUBLE, sphereRadiiTag );RR;
00024 
00025     // need to make sure all interior edges and faces are created
00026     Range all_verts;
00027     result = mbImpl->get_entities_by_type( 0, MBVERTEX, all_verts );RR;
00028     MeshTopoUtil mtu( mbImpl );
00029     result = mtu.construct_aentities( all_verts );RR;
00030 
00031     // create tag to hold vertices
00032     result = mbImpl->tag_get_handle( SUBDIV_VERTICES_TAG_NAME, 9, MB_TYPE_HANDLE, subdivVerticesTag,
00033                                      MB_TAG_DENSE | MB_TAG_EXCL );RR;
00034 
00035     // compute nodal positions for each dimension element
00036     result = compute_nodes( 1 );RR;
00037     result = compute_nodes( 2 );RR;
00038     result = compute_nodes( 3 );RR;
00039 
00040     // build hex elements
00041     std::vector< EntityHandle > sphere_hexes, interstic_hexes;
00042     result = build_hexes( sphere_hexes, interstic_hexes );RR;
00043 
00044     result = mbImpl->tag_delete( subdivVerticesTag );RR;
00045 
00046     if( NULL != hex_set )
00047     {
00048         if( 0 == *hex_set )
00049         {
00050             EntityHandle this_set;
00051             // make a new set
00052             result = mbImpl->create_meshset( MESHSET_SET, this_set );RR;
00053             *hex_set = this_set;
00054         }
00055 
00056         // save all the hexes to this set
00057         result = mbImpl->add_entities( *hex_set, &sphere_hexes[0], sphere_hexes.size() );RR;
00058         result = mbImpl->add_entities( *hex_set, &interstic_hexes[0], interstic_hexes.size() );RR;
00059     }
00060 
00061     return result;
00062 }
00063 
00064 ErrorCode SphereDecomp::compute_nodes( const int dim )
00065 {
00066     // get facets of that dimension
00067     Range these_ents;
00068     const EntityType the_types[4] = { MBVERTEX, MBEDGE, MBTRI, MBTET };
00069 
00070     ErrorCode result = mbImpl->get_entities_by_dimension( 0, dim, these_ents );RR;
00071     assert( mbImpl->type_from_handle( *these_ents.begin() ) == the_types[dim] &&
00072             mbImpl->type_from_handle( *these_ents.rbegin() ) == the_types[dim] );
00073 
00074     EntityHandle subdiv_vertices[9];
00075     MeshTopoUtil mtu( mbImpl );
00076     double avg_pos[3], vert_pos[12], new_vert_pos[12], new_new_vert_pos[3];
00077     double radii[4], unitv[3];
00078     int num_verts = CN::VerticesPerEntity( the_types[dim] );
00079 
00080     for( Range::iterator rit = these_ents.begin(); rit != these_ents.end(); ++rit )
00081     {
00082 
00083         // get vertices
00084         const EntityHandle* connect;
00085         int num_connect;
00086         result = mbImpl->get_connectivity( *rit, connect, num_connect );RR;
00087 
00088         // compute center
00089         result = mtu.get_average_position( connect, num_connect, avg_pos );RR;
00090 
00091         // create center vertex
00092         result = mbImpl->create_vertex( avg_pos, subdiv_vertices[num_verts] );RR;
00093 
00094         // get coords of other vertices
00095         result = mbImpl->get_coords( connect, num_connect, vert_pos );RR;
00096 
00097         // get radii associated with each vertex
00098         result = mbImpl->tag_get_data( sphereRadiiTag, connect, num_connect, radii );RR;
00099 
00100         // compute subdiv vertex position for each vertex
00101         for( int i = 0; i < num_verts; i++ )
00102         {
00103             for( int j = 0; j < 3; j++ )
00104                 unitv[j] = avg_pos[j] - vert_pos[3 * i + j];
00105             double vlength = sqrt( unitv[0] * unitv[0] + unitv[1] * unitv[1] + unitv[2] * unitv[2] );
00106             if( vlength < radii[i] )
00107             {
00108                 std::cout << "Radius too large at vertex " << i << std::endl;
00109                 result = MB_FAILURE;
00110                 continue;
00111             }
00112 
00113             for( int j = 0; j < 3; j++ )
00114                 unitv[j] /= vlength;
00115 
00116             for( int j = 0; j < 3; j++ )
00117                 new_vert_pos[3 * i + j] = vert_pos[3 * i + j] + radii[i] * unitv[j];
00118 
00119             // create vertex at this position
00120             ErrorCode tmp_result = mbImpl->create_vertex( &new_vert_pos[3 * i], subdiv_vertices[i] );
00121             if( MB_SUCCESS != tmp_result ) result = tmp_result;
00122         }
00123 
00124         if( MB_SUCCESS != result ) return result;
00125 
00126         // compute subdiv vertex positions for vertices inside spheres; just mid-pt between
00127         // previous subdiv vertex and corner vertex
00128         for( int i = 0; i < num_verts; i++ )
00129         {
00130             for( int j = 0; j < 3; j++ )
00131                 new_new_vert_pos[j] = .5 * ( vert_pos[3 * i + j] + new_vert_pos[3 * i + j] );
00132 
00133             result = mbImpl->create_vertex( new_new_vert_pos, subdiv_vertices[num_verts + 1 + i] );
00134         }
00135 
00136         // set the tag
00137         result = mbImpl->tag_set_data( subdivVerticesTag, &( *rit ), 1, subdiv_vertices );RR;
00138     }
00139 
00140     return result;
00141 }
00142 
00143 ErrorCode SphereDecomp::build_hexes( std::vector< EntityHandle >& sphere_hexes,
00144                                      std::vector< EntityHandle >& interstic_hexes )
00145 {
00146     // build hexes inside each tet element separately
00147     Range tets;
00148     ErrorCode result = mbImpl->get_entities_by_type( 0, MBTET, tets );RR;
00149 
00150     for( Range::iterator vit = tets.begin(); vit != tets.end(); ++vit )
00151     {
00152         result = subdivide_tet( *vit, sphere_hexes, interstic_hexes );RR;
00153     }
00154 
00155     return MB_SUCCESS;
00156 }
00157 
00158 ErrorCode SphereDecomp::subdivide_tet( EntityHandle tet,
00159                                        std::vector< EntityHandle >& sphere_hexes,
00160                                        std::vector< EntityHandle >& interstic_hexes )
00161 {
00162     // 99: (#subdiv_verts/entity=9) * (#edges=6 + #faces=4 + 1=tet)
00163     EntityHandle subdiv_verts[99];
00164 
00165     // get tet connectivity
00166     std::vector< EntityHandle > tet_conn;
00167     ErrorCode result = mbImpl->get_connectivity( &tet, 1, tet_conn );RR;
00168 
00169     for( int dim = 1; dim <= 3; dim++ )
00170     {
00171         // get entities of this dimension
00172         std::vector< EntityHandle > ents;
00173         if( dim != 3 )
00174         {
00175             result = mbImpl->get_adjacencies( &tet, 1, dim, false, ents );RR;
00176         }
00177         else
00178             ents.push_back( tet );
00179 
00180         // for each, get subdiv verts & put into vector
00181         for( std::vector< EntityHandle >::iterator vit = ents.begin(); vit != ents.end(); ++vit )
00182         {
00183             result = retrieve_subdiv_verts( tet, *vit, &tet_conn[0], dim, subdiv_verts );RR;
00184         }
00185     }
00186 
00187     // ok, subdiv_verts are in canonical order; now create the hexes, using pre-computed templates
00188 
00189     // Templates are specified in terms of the vertices making up each hex; vertices are specified
00190     // by specifying the facet index and type they resolve, and the index of that vertex in that
00191     // facet's subdivision vertices list.
00192 
00193     // Each facet is subdivided into:
00194     // - a mid vertex
00195     // - one vertex for each corner vertex on the facet (located on a line between the mid vertex
00196     // and
00197     //   the corresponding corner vertex, a distance equal to the sphere radius away from the corner
00198     //   vertex)
00199     // - one vertex midway between each corner vertex and the corresponding "sphere surface" vertex
00200     // For edges, tris and tets this gives 5, 7 and 9 subdivision vertices, respectively.
00201     // Subdivision vertices appear in the list in the order: sphere surface vertices, mid vertex,
00202     // sphere interior vertices.  In each of those sub lists, vertices are listed in the canonical
00203     // order of the corresponding corner vertices for that facet.
00204 
00205     // Subdivision vertices for facetes are indexed by listing the facet type they resolve (EDGE,
00206     // FACE, TET), the index of that facet (integer = 0..5, 0..3, 0 for edges, tris, tet, resp), and
00207     // subdivision index (AINDEX..EINDEX for edges, AINDEX..GINDEX for tris, AINDEX..IINDEX for
00208     // tets).
00209 
00210     // Subdivision vertices for all facets of a tet are stored in one subdivision vertex vector, in
00211     // order of increasing facet dimension and index (index varies fastest).  The ESV, FSV, and TSV
00212     // macros are used to compute the indices into that vector for various parameters.  The CV macro
00213     // is used to index into the tet connectivity vector.
00214 
00215     // Subdivision templates for splitting the tet into 28 hexes were derived by hand, and are
00216     // listed below (using the indexing scheme described above).
00217 
00218 #define EDGE        0
00219 #define FACE        1
00220 #define TET         2
00221 #define AINDEX      0
00222 #define BINDEX      1
00223 #define CINDEX      2
00224 #define DINDEX      3
00225 #define EINDEX      4
00226 #define FINDEX      5
00227 #define GINDEX      6
00228 #define HINDEX      7
00229 #define IINDEX      8
00230 #define V0INDEX     0
00231 #define V1INDEX     1
00232 #define V2INDEX     2
00233 #define V3INDEX     3
00234 #define CV( a )     tet_conn[a]
00235 #define ESV( a, b ) subdiv_verts[(a)*9 + ( b )]
00236 #define FSV( a, b ) subdiv_verts[54 + (a)*9 + ( b )]
00237 #define TSV( a, b ) subdiv_verts[90 + (a)*9 + ( b )]
00238 
00239     EntityHandle this_connect[8], this_hex;
00240 
00241     // first, interstices hexes, three per vertex/spherical surface
00242     // V0:
00243     int i             = 0;
00244     this_connect[i++] = ESV( 0, AINDEX );
00245     this_connect[i++] = ESV( 0, CINDEX );
00246     this_connect[i++] = FSV( 3, DINDEX );
00247     this_connect[i++] = FSV( 3, AINDEX );
00248     this_connect[i++] = FSV( 0, AINDEX );
00249     this_connect[i++] = FSV( 0, DINDEX );
00250     this_connect[i++] = TSV( 0, EINDEX );
00251     this_connect[i++] = TSV( 0, AINDEX );
00252     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00253     interstic_hexes.push_back( this_hex );
00254 
00255     i                 = 0;
00256     this_connect[i++] = FSV( 0, AINDEX );
00257     this_connect[i++] = FSV( 0, DINDEX );
00258     this_connect[i++] = TSV( 0, EINDEX );
00259     this_connect[i++] = TSV( 0, AINDEX );
00260     this_connect[i++] = ESV( 3, AINDEX );
00261     this_connect[i++] = ESV( 3, CINDEX );
00262     this_connect[i++] = FSV( 2, DINDEX );
00263     this_connect[i++] = FSV( 2, AINDEX );
00264     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00265     interstic_hexes.push_back( this_hex );
00266 
00267     i                 = 0;
00268     this_connect[i++] = FSV( 3, AINDEX );
00269     this_connect[i++] = FSV( 3, DINDEX );
00270     this_connect[i++] = ESV( 2, CINDEX );
00271     this_connect[i++] = ESV( 2, BINDEX );
00272     this_connect[i++] = TSV( 0, AINDEX );
00273     this_connect[i++] = TSV( 0, EINDEX );
00274     this_connect[i++] = FSV( 2, DINDEX );
00275     this_connect[i++] = FSV( 2, AINDEX );
00276     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00277     interstic_hexes.push_back( this_hex );
00278 
00279     // V1:
00280     i                 = 0;
00281     this_connect[i++] = ESV( 0, CINDEX );
00282     this_connect[i++] = ESV( 0, BINDEX );
00283     this_connect[i++] = FSV( 3, CINDEX );
00284     this_connect[i++] = FSV( 3, DINDEX );
00285     this_connect[i++] = FSV( 0, DINDEX );
00286     this_connect[i++] = FSV( 0, BINDEX );
00287     this_connect[i++] = TSV( 0, BINDEX );
00288     this_connect[i++] = TSV( 0, EINDEX );
00289     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00290     interstic_hexes.push_back( this_hex );
00291 
00292     i                 = 0;
00293     this_connect[i++] = FSV( 0, DINDEX );
00294     this_connect[i++] = FSV( 0, BINDEX );
00295     this_connect[i++] = TSV( 0, BINDEX );
00296     this_connect[i++] = TSV( 0, EINDEX );
00297     this_connect[i++] = ESV( 4, CINDEX );
00298     this_connect[i++] = ESV( 4, AINDEX );
00299     this_connect[i++] = FSV( 1, AINDEX );
00300     this_connect[i++] = FSV( 1, DINDEX );
00301     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00302     interstic_hexes.push_back( this_hex );
00303 
00304     i                 = 0;
00305     this_connect[i++] = FSV( 1, DINDEX );
00306     this_connect[i++] = FSV( 1, AINDEX );
00307     this_connect[i++] = TSV( 0, BINDEX );
00308     this_connect[i++] = TSV( 0, EINDEX );
00309     this_connect[i++] = ESV( 1, CINDEX );
00310     this_connect[i++] = ESV( 1, AINDEX );
00311     this_connect[i++] = FSV( 3, CINDEX );
00312     this_connect[i++] = FSV( 3, DINDEX );
00313     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00314     interstic_hexes.push_back( this_hex );
00315 
00316     // V2:
00317     i                 = 0;
00318     this_connect[i++] = FSV( 3, DINDEX );
00319     this_connect[i++] = ESV( 1, CINDEX );
00320     this_connect[i++] = ESV( 1, BINDEX );
00321     this_connect[i++] = FSV( 3, BINDEX );
00322     this_connect[i++] = TSV( 0, EINDEX );
00323     this_connect[i++] = FSV( 1, DINDEX );
00324     this_connect[i++] = FSV( 1, BINDEX );
00325     this_connect[i++] = TSV( 0, CINDEX );
00326     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00327     interstic_hexes.push_back( this_hex );
00328 
00329     i                 = 0;
00330     this_connect[i++] = TSV( 0, EINDEX );
00331     this_connect[i++] = FSV( 1, DINDEX );
00332     this_connect[i++] = FSV( 1, BINDEX );
00333     this_connect[i++] = TSV( 0, CINDEX );
00334     this_connect[i++] = FSV( 2, DINDEX );
00335     this_connect[i++] = ESV( 5, CINDEX );
00336     this_connect[i++] = ESV( 5, AINDEX );
00337     this_connect[i++] = FSV( 2, CINDEX );
00338     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00339     interstic_hexes.push_back( this_hex );
00340 
00341     i                 = 0;
00342     this_connect[i++] = TSV( 0, CINDEX );
00343     this_connect[i++] = FSV( 2, CINDEX );
00344     this_connect[i++] = ESV( 2, AINDEX );
00345     this_connect[i++] = FSV( 3, BINDEX );
00346     this_connect[i++] = TSV( 0, EINDEX );
00347     this_connect[i++] = FSV( 2, DINDEX );
00348     this_connect[i++] = ESV( 2, CINDEX );
00349     this_connect[i++] = FSV( 3, DINDEX );
00350     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00351     interstic_hexes.push_back( this_hex );
00352 
00353     // V3:
00354     i                 = 0;
00355     this_connect[i++] = TSV( 0, EINDEX );
00356     this_connect[i++] = FSV( 1, DINDEX );
00357     this_connect[i++] = ESV( 5, CINDEX );
00358     this_connect[i++] = FSV( 2, DINDEX );
00359     this_connect[i++] = TSV( 0, DINDEX );
00360     this_connect[i++] = FSV( 1, CINDEX );
00361     this_connect[i++] = ESV( 5, BINDEX );
00362     this_connect[i++] = FSV( 2, BINDEX );
00363     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00364     interstic_hexes.push_back( this_hex );
00365 
00366     i                 = 0;
00367     this_connect[i++] = FSV( 0, DINDEX );
00368     this_connect[i++] = ESV( 4, CINDEX );
00369     this_connect[i++] = FSV( 1, DINDEX );
00370     this_connect[i++] = TSV( 0, EINDEX );
00371     this_connect[i++] = FSV( 0, CINDEX );
00372     this_connect[i++] = ESV( 4, BINDEX );
00373     this_connect[i++] = FSV( 1, CINDEX );
00374     this_connect[i++] = TSV( 0, DINDEX );
00375     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00376     interstic_hexes.push_back( this_hex );
00377 
00378     i                 = 0;
00379     this_connect[i++] = ESV( 3, CINDEX );
00380     this_connect[i++] = FSV( 0, DINDEX );
00381     this_connect[i++] = TSV( 0, EINDEX );
00382     this_connect[i++] = FSV( 2, DINDEX );
00383     this_connect[i++] = ESV( 3, BINDEX );
00384     this_connect[i++] = FSV( 0, CINDEX );
00385     this_connect[i++] = TSV( 0, DINDEX );
00386     this_connect[i++] = FSV( 2, BINDEX );
00387     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00388     interstic_hexes.push_back( this_hex );
00389 
00390     // now, the sphere interiors, four hexes per vertex sphere
00391 
00392     // V0:
00393     i                 = 0;
00394     this_connect[i++] = CV( V0INDEX );
00395     this_connect[i++] = ESV( 0, DINDEX );
00396     this_connect[i++] = FSV( 3, EINDEX );
00397     this_connect[i++] = ESV( 2, EINDEX );
00398     this_connect[i++] = ESV( 3, DINDEX );
00399     this_connect[i++] = FSV( 0, EINDEX );
00400     this_connect[i++] = TSV( 0, FINDEX );
00401     this_connect[i++] = FSV( 2, EINDEX );
00402     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00403     sphere_hexes.push_back( this_hex );
00404 
00405     i                 = 0;
00406     this_connect[i++] = ESV( 0, DINDEX );
00407     this_connect[i++] = ESV( 0, AINDEX );
00408     this_connect[i++] = FSV( 3, AINDEX );
00409     this_connect[i++] = FSV( 3, EINDEX );
00410     this_connect[i++] = FSV( 0, EINDEX );
00411     this_connect[i++] = FSV( 0, AINDEX );
00412     this_connect[i++] = TSV( 0, AINDEX );
00413     this_connect[i++] = TSV( 0, FINDEX );
00414     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00415     sphere_hexes.push_back( this_hex );
00416 
00417     i                 = 0;
00418     this_connect[i++] = FSV( 3, EINDEX );
00419     this_connect[i++] = FSV( 3, AINDEX );
00420     this_connect[i++] = ESV( 2, BINDEX );
00421     this_connect[i++] = ESV( 2, EINDEX );
00422     this_connect[i++] = TSV( 0, FINDEX );
00423     this_connect[i++] = TSV( 0, AINDEX );
00424     this_connect[i++] = FSV( 2, AINDEX );
00425     this_connect[i++] = FSV( 2, EINDEX );
00426     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00427     sphere_hexes.push_back( this_hex );
00428 
00429     i                 = 0;
00430     this_connect[i++] = TSV( 0, FINDEX );
00431     this_connect[i++] = TSV( 0, AINDEX );
00432     this_connect[i++] = FSV( 2, AINDEX );
00433     this_connect[i++] = FSV( 2, EINDEX );
00434     this_connect[i++] = FSV( 0, EINDEX );
00435     this_connect[i++] = FSV( 0, AINDEX );
00436     this_connect[i++] = ESV( 3, AINDEX );
00437     this_connect[i++] = ESV( 3, DINDEX );
00438     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00439     sphere_hexes.push_back( this_hex );
00440 
00441     // V1:
00442     i                 = 0;
00443     this_connect[i++] = CV( V1INDEX );
00444     this_connect[i++] = ESV( 1, DINDEX );
00445     this_connect[i++] = FSV( 3, GINDEX );
00446     this_connect[i++] = ESV( 0, EINDEX );
00447     this_connect[i++] = ESV( 4, DINDEX );
00448     this_connect[i++] = FSV( 1, EINDEX );
00449     this_connect[i++] = TSV( 0, GINDEX );
00450     this_connect[i++] = FSV( 0, FINDEX );
00451     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00452     sphere_hexes.push_back( this_hex );
00453 
00454     i                 = 0;
00455     this_connect[i++] = FSV( 3, GINDEX );
00456     this_connect[i++] = ESV( 1, DINDEX );
00457     this_connect[i++] = ESV( 1, AINDEX );
00458     this_connect[i++] = FSV( 3, CINDEX );
00459     this_connect[i++] = TSV( 0, GINDEX );
00460     this_connect[i++] = FSV( 1, EINDEX );
00461     this_connect[i++] = FSV( 1, AINDEX );
00462     this_connect[i++] = TSV( 0, BINDEX );
00463     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00464     sphere_hexes.push_back( this_hex );
00465 
00466     i                 = 0;
00467     this_connect[i++] = TSV( 0, GINDEX );
00468     this_connect[i++] = FSV( 1, EINDEX );
00469     this_connect[i++] = FSV( 1, AINDEX );
00470     this_connect[i++] = TSV( 0, BINDEX );
00471     this_connect[i++] = FSV( 0, FINDEX );
00472     this_connect[i++] = ESV( 4, DINDEX );
00473     this_connect[i++] = ESV( 4, AINDEX );
00474     this_connect[i++] = FSV( 0, BINDEX );
00475     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00476     sphere_hexes.push_back( this_hex );
00477 
00478     i                 = 0;
00479     this_connect[i++] = ESV( 0, BINDEX );
00480     this_connect[i++] = ESV( 0, EINDEX );
00481     this_connect[i++] = FSV( 3, GINDEX );
00482     this_connect[i++] = FSV( 3, CINDEX );
00483     this_connect[i++] = FSV( 0, BINDEX );
00484     this_connect[i++] = FSV( 0, FINDEX );
00485     this_connect[i++] = TSV( 0, GINDEX );
00486     this_connect[i++] = TSV( 0, BINDEX );
00487     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00488     sphere_hexes.push_back( this_hex );
00489 
00490     // V2:
00491     i                 = 0;
00492     this_connect[i++] = ESV( 1, BINDEX );
00493     this_connect[i++] = ESV( 1, EINDEX );
00494     this_connect[i++] = FSV( 3, FINDEX );
00495     this_connect[i++] = FSV( 3, BINDEX );
00496     this_connect[i++] = FSV( 1, BINDEX );
00497     this_connect[i++] = FSV( 1, FINDEX );
00498     this_connect[i++] = TSV( 0, HINDEX );
00499     this_connect[i++] = TSV( 0, CINDEX );
00500     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00501     sphere_hexes.push_back( this_hex );
00502 
00503     i                 = 0;
00504     this_connect[i++] = FSV( 3, FINDEX );
00505     this_connect[i++] = ESV( 1, EINDEX );
00506     this_connect[i++] = CV( V2INDEX );
00507     this_connect[i++] = ESV( 2, DINDEX );
00508     this_connect[i++] = TSV( 0, HINDEX );
00509     this_connect[i++] = FSV( 1, FINDEX );
00510     this_connect[i++] = ESV( 5, DINDEX );
00511     this_connect[i++] = FSV( 2, GINDEX );
00512     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00513     sphere_hexes.push_back( this_hex );
00514 
00515     i                 = 0;
00516     this_connect[i++] = TSV( 0, HINDEX );
00517     this_connect[i++] = FSV( 1, FINDEX );
00518     this_connect[i++] = ESV( 5, DINDEX );
00519     this_connect[i++] = FSV( 2, GINDEX );
00520     this_connect[i++] = TSV( 0, CINDEX );
00521     this_connect[i++] = FSV( 1, BINDEX );
00522     this_connect[i++] = ESV( 5, AINDEX );
00523     this_connect[i++] = FSV( 2, CINDEX );
00524     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00525     sphere_hexes.push_back( this_hex );
00526 
00527     i                 = 0;
00528     this_connect[i++] = FSV( 3, BINDEX );
00529     this_connect[i++] = FSV( 3, FINDEX );
00530     this_connect[i++] = ESV( 2, DINDEX );
00531     this_connect[i++] = ESV( 2, AINDEX );
00532     this_connect[i++] = TSV( 0, CINDEX );
00533     this_connect[i++] = TSV( 0, HINDEX );
00534     this_connect[i++] = FSV( 2, GINDEX );
00535     this_connect[i++] = FSV( 2, CINDEX );
00536     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00537     sphere_hexes.push_back( this_hex );
00538 
00539     // V3:
00540     i                 = 0;
00541     this_connect[i++] = FSV( 0, CINDEX );
00542     this_connect[i++] = ESV( 4, BINDEX );
00543     this_connect[i++] = FSV( 1, CINDEX );
00544     this_connect[i++] = TSV( 0, DINDEX );
00545     this_connect[i++] = FSV( 0, GINDEX );
00546     this_connect[i++] = ESV( 4, EINDEX );
00547     this_connect[i++] = FSV( 1, GINDEX );
00548     this_connect[i++] = TSV( 0, IINDEX );
00549     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00550     sphere_hexes.push_back( this_hex );
00551 
00552     i                 = 0;
00553     this_connect[i++] = ESV( 3, BINDEX );
00554     this_connect[i++] = FSV( 0, CINDEX );
00555     this_connect[i++] = TSV( 0, DINDEX );
00556     this_connect[i++] = FSV( 2, BINDEX );
00557     this_connect[i++] = ESV( 3, EINDEX );
00558     this_connect[i++] = FSV( 0, GINDEX );
00559     this_connect[i++] = TSV( 0, IINDEX );
00560     this_connect[i++] = FSV( 2, FINDEX );
00561     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00562     sphere_hexes.push_back( this_hex );
00563 
00564     i                 = 0;
00565     this_connect[i++] = TSV( 0, DINDEX );
00566     this_connect[i++] = FSV( 1, CINDEX );
00567     this_connect[i++] = ESV( 5, BINDEX );
00568     this_connect[i++] = FSV( 2, BINDEX );
00569     this_connect[i++] = TSV( 0, IINDEX );
00570     this_connect[i++] = FSV( 1, GINDEX );
00571     this_connect[i++] = ESV( 5, EINDEX );
00572     this_connect[i++] = FSV( 2, FINDEX );
00573     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00574     sphere_hexes.push_back( this_hex );
00575 
00576     i                 = 0;
00577     this_connect[i++] = FSV( 0, GINDEX );
00578     this_connect[i++] = ESV( 4, EINDEX );
00579     this_connect[i++] = FSV( 1, GINDEX );
00580     this_connect[i++] = TSV( 0, IINDEX );
00581     this_connect[i++] = ESV( 3, EINDEX );
00582     this_connect[i++] = CV( V3INDEX );
00583     this_connect[i++] = ESV( 5, EINDEX );
00584     this_connect[i++] = FSV( 2, FINDEX );
00585     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00586     sphere_hexes.push_back( this_hex );
00587 
00588     return result;
00589 }
00590 
00591 ErrorCode SphereDecomp::retrieve_subdiv_verts( EntityHandle tet,
00592                                                EntityHandle this_ent,
00593                                                const EntityHandle* tet_conn,
00594                                                const int dim,
00595                                                EntityHandle* subdiv_verts )
00596 {
00597     // get the subdiv verts for this entity
00598     ErrorCode result;
00599 
00600     // if it's a tet, just put them on the end & return
00601     if( tet == this_ent )
00602     {
00603         result = mbImpl->tag_get_data( subdivVerticesTag, &this_ent, 1, &subdiv_verts[90] );
00604         return MB_SUCCESS;
00605     }
00606 
00607     // if it's a sub-entity, need to find index, relative orientation, and offset
00608     // get connectivity of sub-entity
00609     std::vector< EntityHandle > this_conn;
00610     result = mbImpl->get_connectivity( &this_ent, 1, this_conn );RR;
00611 
00612     // get relative orientation
00613     std::vector< int > conn_tet_indices( this_conn.size() );
00614     for( size_t i = 0; i < this_conn.size(); ++i )
00615         conn_tet_indices[i] = std::find( tet_conn, tet_conn + 4, this_conn[i] ) - tet_conn;
00616     int sense, side_no, offset;
00617     int success = CN::SideNumber( MBTET, &conn_tet_indices[0], this_conn.size(), dim, side_no, sense, offset );
00618     if( -1 == success ) return MB_FAILURE;
00619 
00620     // start of this entity's subdiv_verts; edges go first, then preceding sides, then this one;
00621     // this assumes 6 edges/tet
00622     EntityHandle* subdiv_start = &subdiv_verts[( ( dim - 1 ) * 6 + side_no ) * 9];
00623 
00624     // get subdiv_verts and put them into proper place
00625     result = mbImpl->tag_get_data( subdivVerticesTag, &this_ent, 1, subdiv_start );
00626 
00627     // could probably do this more elegantly, but isn't worth it
00628 #define SWITCH( a, b )                        \
00629     {                                         \
00630         EntityHandle tmp_handle = a;          \
00631         ( a )                   = b;          \
00632         ( b )                   = tmp_handle; \
00633     }
00634     switch( dim )
00635     {
00636         case 1:
00637             if( offset != 0 || sense == -1 )
00638             {
00639                 SWITCH( subdiv_start[0], subdiv_start[1] );
00640                 SWITCH( subdiv_start[3], subdiv_start[4] );
00641             }
00642             break;
00643         case 2:
00644             // rotate first
00645             if( 0 != offset )
00646             {
00647                 std::rotate( subdiv_start, subdiv_start + offset, subdiv_start + 3 );
00648                 std::rotate( subdiv_start + 4, subdiv_start + 4 + offset, subdiv_start + 7 );
00649             }
00650             // now flip, if necessary
00651             if( -1 == sense )
00652             {
00653                 SWITCH( subdiv_start[1], subdiv_start[2] );
00654                 SWITCH( subdiv_start[5], subdiv_start[6] );
00655             }
00656             break;
00657         default:
00658             return MB_FAILURE;
00659     }
00660 
00661     // ok, we're done
00662     return MB_SUCCESS;
00663 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines