MOAB: Mesh Oriented datABase  (version 5.2.1)
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 <math.h>
00006 #include <assert.h>
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, std::vector< EntityHandle >& sphere_hexes,
00159                                        std::vector< EntityHandle >& interstic_hexes )
00160 {
00161     // 99: (#subdiv_verts/entity=9) * (#edges=6 + #faces=4 + 1=tet)
00162     EntityHandle subdiv_verts[99];
00163 
00164     // get tet connectivity
00165     std::vector< EntityHandle > tet_conn;
00166     ErrorCode result = mbImpl->get_connectivity( &tet, 1, tet_conn );RR;
00167 
00168     for( int dim = 1; dim <= 3; dim++ )
00169     {
00170         // get entities of this dimension
00171         std::vector< EntityHandle > ents;
00172         if( dim != 3 )
00173         {
00174             result = mbImpl->get_adjacencies( &tet, 1, dim, false, ents );RR;
00175         }
00176         else
00177             ents.push_back( tet );
00178 
00179         // for each, get subdiv verts & put into vector
00180         for( std::vector< EntityHandle >::iterator vit = ents.begin(); vit != ents.end(); ++vit )
00181         {
00182             result = retrieve_subdiv_verts( tet, *vit, &tet_conn[0], dim, subdiv_verts );RR;
00183         }
00184     }
00185 
00186     // ok, subdiv_verts are in canonical order; now create the hexes, using pre-computed templates
00187 
00188     // Templates are specified in terms of the vertices making up each hex; vertices are specified
00189     // by specifying the facet index and type they resolve, and the index of that vertex in that
00190     // facet's subdivision vertices list.
00191 
00192     // Each facet is subdivided into:
00193     // - a mid vertex
00194     // - one vertex for each corner vertex on the facet (located on a line between the mid vertex
00195     // and
00196     //   the corresponding corner vertex, a distance equal to the sphere radius away from the corner
00197     //   vertex)
00198     // - one vertex midway between each corner vertex and the corresponding "sphere surface" vertex
00199     // For edges, tris and tets this gives 5, 7 and 9 subdivision vertices, respectively.
00200     // Subdivision vertices appear in the list in the order: sphere surface vertices, mid vertex,
00201     // sphere interior vertices.  In each of those sub lists, vertices are listed in the canonical
00202     // order of the corresponding corner vertices for that facet.
00203 
00204     // Subdivision vertices for facetes are indexed by listing the facet type they resolve (EDGE,
00205     // FACE, TET), the index of that facet (integer = 0..5, 0..3, 0 for edges, tris, tet, resp), and
00206     // subdivision index (AINDEX..EINDEX for edges, AINDEX..GINDEX for tris, AINDEX..IINDEX for
00207     // tets).
00208 
00209     // Subdivision vertices for all facets of a tet are stored in one subdivision vertex vector, in
00210     // order of increasing facet dimension and index (index varies fastest).  The ESV, FSV, and TSV
00211     // macros are used to compute the indices into that vector for various parameters.  The CV macro
00212     // is used to index into the tet connectivity vector.
00213 
00214     // Subdivision templates for splitting the tet into 28 hexes were derived by hand, and are
00215     // listed below (using the indexing scheme described above).
00216 
00217 #define EDGE        0
00218 #define FACE        1
00219 #define TET         2
00220 #define AINDEX      0
00221 #define BINDEX      1
00222 #define CINDEX      2
00223 #define DINDEX      3
00224 #define EINDEX      4
00225 #define FINDEX      5
00226 #define GINDEX      6
00227 #define HINDEX      7
00228 #define IINDEX      8
00229 #define V0INDEX     0
00230 #define V1INDEX     1
00231 #define V2INDEX     2
00232 #define V3INDEX     3
00233 #define CV( a )     tet_conn[a]
00234 #define ESV( a, b ) subdiv_verts[a * 9 + b]
00235 #define FSV( a, b ) subdiv_verts[54 + a * 9 + b]
00236 #define TSV( a, b ) subdiv_verts[90 + a * 9 + b]
00237 
00238     EntityHandle this_connect[8], this_hex;
00239 
00240     // first, interstices hexes, three per vertex/spherical surface
00241     // V0:
00242     int i             = 0;
00243     this_connect[i++] = ESV( 0, AINDEX );
00244     this_connect[i++] = ESV( 0, CINDEX );
00245     this_connect[i++] = FSV( 3, DINDEX );
00246     this_connect[i++] = FSV( 3, AINDEX );
00247     this_connect[i++] = FSV( 0, AINDEX );
00248     this_connect[i++] = FSV( 0, DINDEX );
00249     this_connect[i++] = TSV( 0, EINDEX );
00250     this_connect[i++] = TSV( 0, AINDEX );
00251     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00252     interstic_hexes.push_back( this_hex );
00253 
00254     i                 = 0;
00255     this_connect[i++] = FSV( 0, AINDEX );
00256     this_connect[i++] = FSV( 0, DINDEX );
00257     this_connect[i++] = TSV( 0, EINDEX );
00258     this_connect[i++] = TSV( 0, AINDEX );
00259     this_connect[i++] = ESV( 3, AINDEX );
00260     this_connect[i++] = ESV( 3, CINDEX );
00261     this_connect[i++] = FSV( 2, DINDEX );
00262     this_connect[i++] = FSV( 2, AINDEX );
00263     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00264     interstic_hexes.push_back( this_hex );
00265 
00266     i                 = 0;
00267     this_connect[i++] = FSV( 3, AINDEX );
00268     this_connect[i++] = FSV( 3, DINDEX );
00269     this_connect[i++] = ESV( 2, CINDEX );
00270     this_connect[i++] = ESV( 2, BINDEX );
00271     this_connect[i++] = TSV( 0, AINDEX );
00272     this_connect[i++] = TSV( 0, EINDEX );
00273     this_connect[i++] = FSV( 2, DINDEX );
00274     this_connect[i++] = FSV( 2, AINDEX );
00275     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00276     interstic_hexes.push_back( this_hex );
00277 
00278     // V1:
00279     i                 = 0;
00280     this_connect[i++] = ESV( 0, CINDEX );
00281     this_connect[i++] = ESV( 0, BINDEX );
00282     this_connect[i++] = FSV( 3, CINDEX );
00283     this_connect[i++] = FSV( 3, DINDEX );
00284     this_connect[i++] = FSV( 0, DINDEX );
00285     this_connect[i++] = FSV( 0, BINDEX );
00286     this_connect[i++] = TSV( 0, BINDEX );
00287     this_connect[i++] = TSV( 0, EINDEX );
00288     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00289     interstic_hexes.push_back( this_hex );
00290 
00291     i                 = 0;
00292     this_connect[i++] = FSV( 0, DINDEX );
00293     this_connect[i++] = FSV( 0, BINDEX );
00294     this_connect[i++] = TSV( 0, BINDEX );
00295     this_connect[i++] = TSV( 0, EINDEX );
00296     this_connect[i++] = ESV( 4, CINDEX );
00297     this_connect[i++] = ESV( 4, AINDEX );
00298     this_connect[i++] = FSV( 1, AINDEX );
00299     this_connect[i++] = FSV( 1, DINDEX );
00300     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00301     interstic_hexes.push_back( this_hex );
00302 
00303     i                 = 0;
00304     this_connect[i++] = FSV( 1, DINDEX );
00305     this_connect[i++] = FSV( 1, AINDEX );
00306     this_connect[i++] = TSV( 0, BINDEX );
00307     this_connect[i++] = TSV( 0, EINDEX );
00308     this_connect[i++] = ESV( 1, CINDEX );
00309     this_connect[i++] = ESV( 1, AINDEX );
00310     this_connect[i++] = FSV( 3, CINDEX );
00311     this_connect[i++] = FSV( 3, DINDEX );
00312     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00313     interstic_hexes.push_back( this_hex );
00314 
00315     // V2:
00316     i                 = 0;
00317     this_connect[i++] = FSV( 3, DINDEX );
00318     this_connect[i++] = ESV( 1, CINDEX );
00319     this_connect[i++] = ESV( 1, BINDEX );
00320     this_connect[i++] = FSV( 3, BINDEX );
00321     this_connect[i++] = TSV( 0, EINDEX );
00322     this_connect[i++] = FSV( 1, DINDEX );
00323     this_connect[i++] = FSV( 1, BINDEX );
00324     this_connect[i++] = TSV( 0, CINDEX );
00325     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00326     interstic_hexes.push_back( this_hex );
00327 
00328     i                 = 0;
00329     this_connect[i++] = TSV( 0, EINDEX );
00330     this_connect[i++] = FSV( 1, DINDEX );
00331     this_connect[i++] = FSV( 1, BINDEX );
00332     this_connect[i++] = TSV( 0, CINDEX );
00333     this_connect[i++] = FSV( 2, DINDEX );
00334     this_connect[i++] = ESV( 5, CINDEX );
00335     this_connect[i++] = ESV( 5, AINDEX );
00336     this_connect[i++] = FSV( 2, CINDEX );
00337     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00338     interstic_hexes.push_back( this_hex );
00339 
00340     i                 = 0;
00341     this_connect[i++] = TSV( 0, CINDEX );
00342     this_connect[i++] = FSV( 2, CINDEX );
00343     this_connect[i++] = ESV( 2, AINDEX );
00344     this_connect[i++] = FSV( 3, BINDEX );
00345     this_connect[i++] = TSV( 0, EINDEX );
00346     this_connect[i++] = FSV( 2, DINDEX );
00347     this_connect[i++] = ESV( 2, CINDEX );
00348     this_connect[i++] = FSV( 3, DINDEX );
00349     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00350     interstic_hexes.push_back( this_hex );
00351 
00352     // V3:
00353     i                 = 0;
00354     this_connect[i++] = TSV( 0, EINDEX );
00355     this_connect[i++] = FSV( 1, DINDEX );
00356     this_connect[i++] = ESV( 5, CINDEX );
00357     this_connect[i++] = FSV( 2, DINDEX );
00358     this_connect[i++] = TSV( 0, DINDEX );
00359     this_connect[i++] = FSV( 1, CINDEX );
00360     this_connect[i++] = ESV( 5, BINDEX );
00361     this_connect[i++] = FSV( 2, BINDEX );
00362     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00363     interstic_hexes.push_back( this_hex );
00364 
00365     i                 = 0;
00366     this_connect[i++] = FSV( 0, DINDEX );
00367     this_connect[i++] = ESV( 4, CINDEX );
00368     this_connect[i++] = FSV( 1, DINDEX );
00369     this_connect[i++] = TSV( 0, EINDEX );
00370     this_connect[i++] = FSV( 0, CINDEX );
00371     this_connect[i++] = ESV( 4, BINDEX );
00372     this_connect[i++] = FSV( 1, CINDEX );
00373     this_connect[i++] = TSV( 0, DINDEX );
00374     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00375     interstic_hexes.push_back( this_hex );
00376 
00377     i                 = 0;
00378     this_connect[i++] = ESV( 3, CINDEX );
00379     this_connect[i++] = FSV( 0, DINDEX );
00380     this_connect[i++] = TSV( 0, EINDEX );
00381     this_connect[i++] = FSV( 2, DINDEX );
00382     this_connect[i++] = ESV( 3, BINDEX );
00383     this_connect[i++] = FSV( 0, CINDEX );
00384     this_connect[i++] = TSV( 0, DINDEX );
00385     this_connect[i++] = FSV( 2, BINDEX );
00386     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00387     interstic_hexes.push_back( this_hex );
00388 
00389     // now, the sphere interiors, four hexes per vertex sphere
00390 
00391     // V0:
00392     i                 = 0;
00393     this_connect[i++] = CV( V0INDEX );
00394     this_connect[i++] = ESV( 0, DINDEX );
00395     this_connect[i++] = FSV( 3, EINDEX );
00396     this_connect[i++] = ESV( 2, EINDEX );
00397     this_connect[i++] = ESV( 3, DINDEX );
00398     this_connect[i++] = FSV( 0, EINDEX );
00399     this_connect[i++] = TSV( 0, FINDEX );
00400     this_connect[i++] = FSV( 2, EINDEX );
00401     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00402     sphere_hexes.push_back( this_hex );
00403 
00404     i                 = 0;
00405     this_connect[i++] = ESV( 0, DINDEX );
00406     this_connect[i++] = ESV( 0, AINDEX );
00407     this_connect[i++] = FSV( 3, AINDEX );
00408     this_connect[i++] = FSV( 3, EINDEX );
00409     this_connect[i++] = FSV( 0, EINDEX );
00410     this_connect[i++] = FSV( 0, AINDEX );
00411     this_connect[i++] = TSV( 0, AINDEX );
00412     this_connect[i++] = TSV( 0, FINDEX );
00413     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00414     sphere_hexes.push_back( this_hex );
00415 
00416     i                 = 0;
00417     this_connect[i++] = FSV( 3, EINDEX );
00418     this_connect[i++] = FSV( 3, AINDEX );
00419     this_connect[i++] = ESV( 2, BINDEX );
00420     this_connect[i++] = ESV( 2, EINDEX );
00421     this_connect[i++] = TSV( 0, FINDEX );
00422     this_connect[i++] = TSV( 0, AINDEX );
00423     this_connect[i++] = FSV( 2, AINDEX );
00424     this_connect[i++] = FSV( 2, EINDEX );
00425     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00426     sphere_hexes.push_back( this_hex );
00427 
00428     i                 = 0;
00429     this_connect[i++] = TSV( 0, FINDEX );
00430     this_connect[i++] = TSV( 0, AINDEX );
00431     this_connect[i++] = FSV( 2, AINDEX );
00432     this_connect[i++] = FSV( 2, EINDEX );
00433     this_connect[i++] = FSV( 0, EINDEX );
00434     this_connect[i++] = FSV( 0, AINDEX );
00435     this_connect[i++] = ESV( 3, AINDEX );
00436     this_connect[i++] = ESV( 3, DINDEX );
00437     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00438     sphere_hexes.push_back( this_hex );
00439 
00440     // V1:
00441     i                 = 0;
00442     this_connect[i++] = CV( V1INDEX );
00443     this_connect[i++] = ESV( 1, DINDEX );
00444     this_connect[i++] = FSV( 3, GINDEX );
00445     this_connect[i++] = ESV( 0, EINDEX );
00446     this_connect[i++] = ESV( 4, DINDEX );
00447     this_connect[i++] = FSV( 1, EINDEX );
00448     this_connect[i++] = TSV( 0, GINDEX );
00449     this_connect[i++] = FSV( 0, FINDEX );
00450     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00451     sphere_hexes.push_back( this_hex );
00452 
00453     i                 = 0;
00454     this_connect[i++] = FSV( 3, GINDEX );
00455     this_connect[i++] = ESV( 1, DINDEX );
00456     this_connect[i++] = ESV( 1, AINDEX );
00457     this_connect[i++] = FSV( 3, CINDEX );
00458     this_connect[i++] = TSV( 0, GINDEX );
00459     this_connect[i++] = FSV( 1, EINDEX );
00460     this_connect[i++] = FSV( 1, AINDEX );
00461     this_connect[i++] = TSV( 0, BINDEX );
00462     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00463     sphere_hexes.push_back( this_hex );
00464 
00465     i                 = 0;
00466     this_connect[i++] = TSV( 0, GINDEX );
00467     this_connect[i++] = FSV( 1, EINDEX );
00468     this_connect[i++] = FSV( 1, AINDEX );
00469     this_connect[i++] = TSV( 0, BINDEX );
00470     this_connect[i++] = FSV( 0, FINDEX );
00471     this_connect[i++] = ESV( 4, DINDEX );
00472     this_connect[i++] = ESV( 4, AINDEX );
00473     this_connect[i++] = FSV( 0, BINDEX );
00474     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00475     sphere_hexes.push_back( this_hex );
00476 
00477     i                 = 0;
00478     this_connect[i++] = ESV( 0, BINDEX );
00479     this_connect[i++] = ESV( 0, EINDEX );
00480     this_connect[i++] = FSV( 3, GINDEX );
00481     this_connect[i++] = FSV( 3, CINDEX );
00482     this_connect[i++] = FSV( 0, BINDEX );
00483     this_connect[i++] = FSV( 0, FINDEX );
00484     this_connect[i++] = TSV( 0, GINDEX );
00485     this_connect[i++] = TSV( 0, BINDEX );
00486     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00487     sphere_hexes.push_back( this_hex );
00488 
00489     // V2:
00490     i                 = 0;
00491     this_connect[i++] = ESV( 1, BINDEX );
00492     this_connect[i++] = ESV( 1, EINDEX );
00493     this_connect[i++] = FSV( 3, FINDEX );
00494     this_connect[i++] = FSV( 3, BINDEX );
00495     this_connect[i++] = FSV( 1, BINDEX );
00496     this_connect[i++] = FSV( 1, FINDEX );
00497     this_connect[i++] = TSV( 0, HINDEX );
00498     this_connect[i++] = TSV( 0, CINDEX );
00499     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00500     sphere_hexes.push_back( this_hex );
00501 
00502     i                 = 0;
00503     this_connect[i++] = FSV( 3, FINDEX );
00504     this_connect[i++] = ESV( 1, EINDEX );
00505     this_connect[i++] = CV( V2INDEX );
00506     this_connect[i++] = ESV( 2, DINDEX );
00507     this_connect[i++] = TSV( 0, HINDEX );
00508     this_connect[i++] = FSV( 1, FINDEX );
00509     this_connect[i++] = ESV( 5, DINDEX );
00510     this_connect[i++] = FSV( 2, GINDEX );
00511     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00512     sphere_hexes.push_back( this_hex );
00513 
00514     i                 = 0;
00515     this_connect[i++] = TSV( 0, HINDEX );
00516     this_connect[i++] = FSV( 1, FINDEX );
00517     this_connect[i++] = ESV( 5, DINDEX );
00518     this_connect[i++] = FSV( 2, GINDEX );
00519     this_connect[i++] = TSV( 0, CINDEX );
00520     this_connect[i++] = FSV( 1, BINDEX );
00521     this_connect[i++] = ESV( 5, AINDEX );
00522     this_connect[i++] = FSV( 2, CINDEX );
00523     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00524     sphere_hexes.push_back( this_hex );
00525 
00526     i                 = 0;
00527     this_connect[i++] = FSV( 3, BINDEX );
00528     this_connect[i++] = FSV( 3, FINDEX );
00529     this_connect[i++] = ESV( 2, DINDEX );
00530     this_connect[i++] = ESV( 2, AINDEX );
00531     this_connect[i++] = TSV( 0, CINDEX );
00532     this_connect[i++] = TSV( 0, HINDEX );
00533     this_connect[i++] = FSV( 2, GINDEX );
00534     this_connect[i++] = FSV( 2, CINDEX );
00535     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00536     sphere_hexes.push_back( this_hex );
00537 
00538     // V3:
00539     i                 = 0;
00540     this_connect[i++] = FSV( 0, CINDEX );
00541     this_connect[i++] = ESV( 4, BINDEX );
00542     this_connect[i++] = FSV( 1, CINDEX );
00543     this_connect[i++] = TSV( 0, DINDEX );
00544     this_connect[i++] = FSV( 0, GINDEX );
00545     this_connect[i++] = ESV( 4, EINDEX );
00546     this_connect[i++] = FSV( 1, GINDEX );
00547     this_connect[i++] = TSV( 0, IINDEX );
00548     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00549     sphere_hexes.push_back( this_hex );
00550 
00551     i                 = 0;
00552     this_connect[i++] = ESV( 3, BINDEX );
00553     this_connect[i++] = FSV( 0, CINDEX );
00554     this_connect[i++] = TSV( 0, DINDEX );
00555     this_connect[i++] = FSV( 2, BINDEX );
00556     this_connect[i++] = ESV( 3, EINDEX );
00557     this_connect[i++] = FSV( 0, GINDEX );
00558     this_connect[i++] = TSV( 0, IINDEX );
00559     this_connect[i++] = FSV( 2, FINDEX );
00560     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00561     sphere_hexes.push_back( this_hex );
00562 
00563     i                 = 0;
00564     this_connect[i++] = TSV( 0, DINDEX );
00565     this_connect[i++] = FSV( 1, CINDEX );
00566     this_connect[i++] = ESV( 5, BINDEX );
00567     this_connect[i++] = FSV( 2, BINDEX );
00568     this_connect[i++] = TSV( 0, IINDEX );
00569     this_connect[i++] = FSV( 1, GINDEX );
00570     this_connect[i++] = ESV( 5, EINDEX );
00571     this_connect[i++] = FSV( 2, FINDEX );
00572     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00573     sphere_hexes.push_back( this_hex );
00574 
00575     i                 = 0;
00576     this_connect[i++] = FSV( 0, GINDEX );
00577     this_connect[i++] = ESV( 4, EINDEX );
00578     this_connect[i++] = FSV( 1, GINDEX );
00579     this_connect[i++] = TSV( 0, IINDEX );
00580     this_connect[i++] = ESV( 3, EINDEX );
00581     this_connect[i++] = CV( V3INDEX );
00582     this_connect[i++] = ESV( 5, EINDEX );
00583     this_connect[i++] = FSV( 2, FINDEX );
00584     result            = mbImpl->create_element( MBHEX, this_connect, 8, this_hex );RR;
00585     sphere_hexes.push_back( this_hex );
00586 
00587     return result;
00588 }
00589 
00590 ErrorCode SphereDecomp::retrieve_subdiv_verts( EntityHandle tet, EntityHandle this_ent, const EntityHandle* tet_conn,
00591                                                const int dim, EntityHandle* subdiv_verts )
00592 {
00593     // get the subdiv verts for this entity
00594     ErrorCode result;
00595 
00596     // if it's a tet, just put them on the end & return
00597     if( tet == this_ent )
00598     {
00599         result = mbImpl->tag_get_data( subdivVerticesTag, &this_ent, 1, &subdiv_verts[90] );
00600         return MB_SUCCESS;
00601     }
00602 
00603     // if it's a sub-entity, need to find index, relative orientation, and offset
00604     // get connectivity of sub-entity
00605     std::vector< EntityHandle > this_conn;
00606     result = mbImpl->get_connectivity( &this_ent, 1, this_conn );RR;
00607 
00608     // get relative orientation
00609     std::vector< int > conn_tet_indices( this_conn.size() );
00610     for( size_t i = 0; i < this_conn.size(); ++i )
00611         conn_tet_indices[i] = std::find( tet_conn, tet_conn + 4, this_conn[i] ) - tet_conn;
00612     int sense, side_no, offset;
00613     int success = CN::SideNumber( MBTET, &conn_tet_indices[0], this_conn.size(), dim, side_no, sense, offset );
00614     if( -1 == success ) return MB_FAILURE;
00615 
00616     // start of this entity's subdiv_verts; edges go first, then preceding sides, then this one;
00617     // this assumes 6 edges/tet
00618     EntityHandle* subdiv_start = &subdiv_verts[( ( dim - 1 ) * 6 + side_no ) * 9];
00619 
00620     // get subdiv_verts and put them into proper place
00621     result = mbImpl->tag_get_data( subdivVerticesTag, &this_ent, 1, subdiv_start );
00622 
00623     // could probably do this more elegantly, but isn't worth it
00624 #define SWITCH( a, b )                        \
00625     {                                         \
00626         EntityHandle tmp_handle = a;          \
00627         a                       = b;          \
00628         b                       = tmp_handle; \
00629     }
00630     switch( dim )
00631     {
00632         case 1:
00633             if( offset != 0 || sense == -1 )
00634             {
00635                 SWITCH( subdiv_start[0], subdiv_start[1] );
00636                 SWITCH( subdiv_start[3], subdiv_start[4] );
00637             }
00638             break;
00639         case 2:
00640             // rotate first
00641             if( 0 != offset )
00642             {
00643                 std::rotate( subdiv_start, subdiv_start + offset, subdiv_start + 3 );
00644                 std::rotate( subdiv_start + 4, subdiv_start + 4 + offset, subdiv_start + 7 );
00645             }
00646             // now flip, if necessary
00647             if( -1 == sense )
00648             {
00649                 SWITCH( subdiv_start[1], subdiv_start[2] );
00650                 SWITCH( subdiv_start[5], subdiv_start[6] );
00651             }
00652             break;
00653         default:
00654             return MB_FAILURE;
00655     }
00656 
00657     // ok, we're done
00658     return MB_SUCCESS;
00659 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines