MOAB: Mesh Oriented datABase  (version 5.4.1)
MsqMOAB.cpp
Go to the documentation of this file.
00001 /** \file MsqMOAB.cpp
00002  *  \brief
00003  *  \author Vijay Mahadevan
00004  */
00005 #include "MsqMOAB.hpp"
00006 #include "MsqError.hpp"
00007 #include "MeshInterface.hpp"
00008 #include "MsqDebug.hpp"
00009 #include "MsqVertex.hpp"
00010 #include <cassert>
00011 #include "MsqIBase.hpp"
00012 #include <algorithm>
00013 
00014 #ifdef IMESH_MAJOR_VERSION
00015 #define IMESH_VERSION_ATLEAST( MAJOR, MINOR ) 1000 * IMESH_MAJOR_VERSION + IMESH_MINOR_VERSION <= 1000 * MAJOR + MINOR
00016 #else
00017 #define IMESH_VERSION_ATLEAST( MAJOR, MINOR ) 0
00018 #endif
00019 
00020 namespace MBMesquite
00021 {
00022 
00023 /*************************************************************************
00024  *                          Mesh Definition
00025  ************************************************************************/
00026 
00027 MsqMOAB::MsqMOAB( moab::Core* mesh,
00028                   moab::EntityHandle meshset,
00029                   moab::EntityType type,
00030                   MsqError& err,
00031                   const moab::Tag* fixed_tag,
00032                   const moab::Tag* slaved_tag )
00033     : meshInstance( mesh ), inputSetType( moab::MBMAXTYPE ), inputSet( 0 ), byteTag( 0 ), createdByteTag( false ),
00034       geometricDimension( 0 )
00035 {
00036     init_active_mesh( mesh, err, fixed_tag, slaved_tag );MSQ_ERRRTN( err );
00037     set_active_set( meshset, type, err );MSQ_ERRRTN( err );
00038 }
00039 
00040 MsqMOAB::MsqMOAB( moab::Core* mesh,
00041                   moab::EntityType type,
00042                   MsqError& err,
00043                   const moab::Tag* fixed_tag,
00044                   const moab::Tag* slaved_tag )
00045     : meshInstance( mesh ), inputSetType( moab::MBMAXTYPE ), inputSet( 0 ), byteTag( 0 ), createdByteTag( false ),
00046       geometricDimension( 0 )
00047 {
00048     init_active_mesh( mesh, err, fixed_tag, slaved_tag );MSQ_ERRRTN( err );
00049 
00050     moab::EntityHandle root_set = 0 /*meshInstance->get_root_set()*/;
00051     set_active_set( root_set, type, err );MSQ_ERRRTN( err );
00052 }
00053 
00054 MsqMOAB::~MsqMOAB()
00055 {
00056     moab::ErrorCode ierr;
00057 
00058     if( createdByteTag )
00059     {
00060         ierr = meshInstance->tag_delete( byteTag );MB_CHK_ERR_RET( ierr );
00061     }
00062 }
00063 
00064 moab::Core* MsqMOAB::get_interface() const
00065 {
00066     return meshInstance;
00067 }
00068 
00069 moab::EntityHandle MsqMOAB::get_entity_set() const
00070 {
00071     return inputSet;
00072 }
00073 
00074 moab::DataType MsqMOAB::check_valid_flag_tag( moab::Tag tag, const char* /*which_flag*/, MsqError& err )
00075 {
00076     moab::ErrorCode ierr;
00077     int size;
00078     std::string name;
00079     moab::DataType type = moab::MB_MAX_DATA_TYPE;
00080 
00081     ierr = meshInstance->tag_get_data_type( tag, type );MB_CHK_ERR_RET_VAL( ierr, type );
00082     ierr = meshInstance->tag_get_name( tag, name );MB_CHK_ERR_RET_VAL( ierr, type );
00083     ierr = meshInstance->tag_get_length( tag, size );MB_CHK_ERR_RET_VAL( ierr, type );
00084 
00085     err.set_error( MBMesquite::MsqError::NO_ERROR );
00086     return type;
00087 }
00088 
00089 void MsqMOAB::init_active_mesh( moab::Core* /*mesh*/,
00090                                 MsqError& err,
00091                                 const moab::Tag* fixed_tag,
00092                                 const moab::Tag* slaved_tag )
00093 {
00094     moab::ErrorCode ierr;
00095 
00096     // Initialize topology map
00097     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00098 
00099     const size_t mapsize = sizeof( topologyMap ) / sizeof( MBMesquite::EntityTopology );
00100 
00101     if( mapsize < moab::MBMAXTYPE )
00102     {
00103         MSQ_SETERR( err )
00104         ( "MsqMOAB needs to be updated for new iMesh element topologies.", MsqError::INTERNAL_ERROR );
00105     }
00106 
00107     for( size_t i = 0; i <= moab::MBMAXTYPE; ++i )
00108     {
00109         topologyMap[i] = MBMesquite::MIXED;
00110     }
00111 
00112     topologyMap[moab::MBTRI]     = MBMesquite::TRIANGLE;
00113     topologyMap[moab::MBQUAD]    = MBMesquite::QUADRILATERAL;
00114     topologyMap[moab::MBTET]     = MBMesquite::TETRAHEDRON;
00115     topologyMap[moab::MBHEX]     = MBMesquite::HEXAHEDRON;
00116     topologyMap[moab::MBPRISM]   = MBMesquite::PRISM;
00117     topologyMap[moab::MBPYRAMID] = MBMesquite::PYRAMID;
00118 
00119     // Check that fixed tag is valid
00120     haveFixedTag = false;
00121 
00122     if( fixed_tag )
00123     {
00124         fixedTagType = check_valid_flag_tag( *fixed_tag, "fixed", err );MSQ_ERRRTN( err );
00125         haveFixedTag = true;
00126         fixedTag     = *fixed_tag;
00127     }
00128 
00129     // Check that slaved tag is valid
00130     haveSlavedTag = false;
00131 
00132     if( slaved_tag )
00133     {
00134         slavedTagType = check_valid_flag_tag( *slaved_tag, "slaved", err );MSQ_ERRRTN( err );
00135         haveSlavedTag = true;
00136         slavedTag     = *slaved_tag;
00137     }
00138 
00139     // Get/create tag for vertex byte
00140     ierr = meshInstance->tag_get_handle( VERTEX_BYTE_TAG_NAME, byteTag );
00141 
00142     if( moab::MB_SUCCESS != ierr )
00143     {
00144         int defval = 0;
00145         ierr       = meshInstance->tag_get_handle( VERTEX_BYTE_TAG_NAME, 1, moab::MB_TYPE_INTEGER, byteTag,
00146                                                    moab::MB_TAG_CREAT | moab::MB_TAG_DENSE, &defval,
00147                                                    &createdByteTag );  // MB_CHK_ERR_RET(ierr);
00148     }
00149     else
00150     {
00151         int size;
00152         moab::DataType type;
00153         ierr = meshInstance->tag_get_length( byteTag, size );MB_CHK_ERR_RET( ierr );
00154         ierr = meshInstance->tag_get_data_type( byteTag, type );MB_CHK_ERR_RET( ierr );
00155     }
00156 
00157     ierr = meshInstance->get_dimension( geometricDimension );MB_CHK_ERR_RET( ierr );
00158     err.set_error( MBMesquite::MsqError::NO_ERROR );
00159 }
00160 
00161 void MsqMOAB::set_fixed_tag( moab::Tag tag, MsqError& err )
00162 {
00163     moab::DataType t = check_valid_flag_tag( tag, "fixed", err );MSQ_ERRRTN( err );
00164     fixedTag     = tag;
00165     fixedTagType = t;
00166     haveFixedTag = true;
00167 }
00168 
00169 void MsqMOAB::clear_fixed_tag()
00170 {
00171     haveFixedTag = false;
00172 }
00173 
00174 const moab::Tag* MsqMOAB::get_fixed_tag() const
00175 {
00176     return haveFixedTag ? &fixedTag : 0;
00177 }
00178 
00179 void MsqMOAB::set_slaved_tag( moab::Tag tag, MsqError& err )
00180 {
00181     moab::DataType t = check_valid_flag_tag( tag, "slaved", err );MSQ_ERRRTN( err );
00182     slavedTag     = tag;
00183     slavedTagType = t;
00184     haveSlavedTag = true;
00185 }
00186 
00187 void MsqMOAB::clear_slaved_tag()
00188 {
00189     haveSlavedTag = false;
00190 }
00191 
00192 const moab::Tag* MsqMOAB::get_slaved_tag() const
00193 {
00194     return haveSlavedTag ? &slavedTag : 0;
00195 }
00196 
00197 void MsqMOAB::set_active_set( moab::EntityHandle elem_set, moab::EntityType type_in, MsqError& err )
00198 {
00199     inputSetType = type_in;
00200     inputSet     = elem_set;
00201 
00202     // clear vertex byte
00203     std::vector< VertexHandle > verts;
00204     get_all_vertices( verts, err );MSQ_ERRRTN( err );
00205 
00206     if( !verts.empty() )
00207     {
00208         std::vector< unsigned char > zeros( verts.size(), 0 );
00209         vertices_set_byte( arrptr( verts ), arrptr( zeros ), verts.size(), err );MSQ_CHKERR( err );
00210     }
00211 }
00212 
00213 // Returns whether this mesh lies in a 2D or 3D coordinate system.
00214 int MsqMOAB::get_geometric_dimension( MBMesquite::MsqError& /*err*/ )
00215 {
00216     return geometricDimension;
00217 }
00218 
00219 //************ Vertex Properties ********************
00220 
00221 void MsqMOAB::get_flag_data( moab::Tag tag,
00222                              bool have_tag,
00223                              moab::DataType type,
00224                              const VertexHandle vert_array[],
00225                              std::vector< bool >& flag_array,
00226                              size_t num_vtx,
00227                              MsqError& err )
00228 {
00229     if( !num_vtx )
00230     {
00231         return;
00232     }
00233 
00234     if( !have_tag )
00235     {
00236         flag_array.clear();
00237         flag_array.resize( num_vtx, false );
00238         return;
00239     }
00240 
00241     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00242     flag_array.resize( num_vtx );
00243 
00244     assert( sizeof( VertexHandle ) == sizeof( moab::EntityHandle ) );
00245     const moab::EntityHandle* arr = reinterpret_cast< const moab::EntityHandle* >( vert_array );
00246 
00247     moab::ErrorCode ierr;
00248     int alloc = num_vtx;
00249     assert( (size_t)alloc == num_vtx );  // size_t can hold larger values than int if 64-bit
00250 
00251     if( type == moab::MB_TYPE_INTEGER )
00252     {
00253         std::vector< int > values( num_vtx );
00254         int* ptr = arrptr( values );
00255         ierr     = meshInstance->tag_get_data( tag, arr, num_vtx, ptr );MB_CHK_ERR_RET( ierr );
00256 
00257         for( size_t i = 0; i < num_vtx; ++i )
00258         {
00259             flag_array[i] = !!values[i];
00260         }
00261     }
00262     else if( type == moab::MB_TYPE_OPAQUE )
00263     {
00264         std::vector< char > values( num_vtx );
00265         void* ptr = arrptr( values );
00266         ierr      = meshInstance->tag_get_data( tag, arr, num_vtx, ptr );MB_CHK_ERR_RET( ierr );
00267 
00268         for( size_t i = 0; i < num_vtx; ++i )
00269         {
00270             flag_array[i] = !!values[i];
00271         }
00272     }
00273     else
00274     {
00275         MSQ_SETERR( err )( "Invalid tag type for vertex flag data", MsqError::INVALID_STATE );
00276         return;
00277     }
00278 
00279     err.set_error( MBMesquite::MsqError::NO_ERROR );
00280 }
00281 
00282 // Returns true or false, indicating whether the vertex
00283 // is allowed to be repositioned.  True indicates that the vertex
00284 // is fixed and cannot be moved.  Note that this is a read-only
00285 // property; this flag can't be modified by users of the
00286 // MBMesquite::Mesh interface.
00287 void MsqMOAB::vertices_get_fixed_flag( const VertexHandle vert_array[],
00288                                        std::vector< bool >& bool_array,
00289                                        size_t num_vtx,
00290                                        MsqError& err )
00291 {
00292     get_flag_data( fixedTag, haveFixedTag, fixedTagType, vert_array, bool_array, num_vtx, err );
00293 }
00294 
00295 void MsqMOAB::vertices_get_slaved_flag( const VertexHandle vert_array[],
00296                                         std::vector< bool >& bool_array,
00297                                         size_t num_vtx,
00298                                         MsqError& err )
00299 {
00300     get_flag_data( slavedTag, haveSlavedTag, slavedTagType, vert_array, bool_array, num_vtx, err );
00301 }
00302 
00303 // Get vertex coordinates
00304 void MsqMOAB::vertices_get_coordinates( const MBMesquite::Mesh::VertexHandle vert_array[],
00305                                         MsqVertex* coordinates,
00306                                         size_t num_vtx,
00307                                         MsqError& err )
00308 {
00309     if( !num_vtx )
00310     {
00311         return;
00312     }
00313 
00314     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00315     std::vector< double > dbl_store( 3 * num_vtx );
00316     double* dbl_array = arrptr( dbl_store );
00317 
00318     moab::ErrorCode ierr;
00319     // int junk = 3*num_vtx;
00320     assert( sizeof( VertexHandle ) == sizeof( moab::EntityHandle ) );
00321     const moab::EntityHandle* arr = reinterpret_cast< const moab::EntityHandle* >( vert_array );
00322     ierr                          = meshInstance->get_coords( arr, num_vtx, dbl_array );MB_CHK_ERR_RET( ierr );
00323 
00324     if( geometricDimension == 2 )
00325     {
00326         double* iter = dbl_array;
00327 
00328         for( size_t i = 0; i < num_vtx; ++i )
00329         {
00330             coordinates[i].x( *iter );
00331             ++iter;
00332             coordinates[i].y( *iter );
00333             ++iter;
00334             coordinates[i].z( 0 );
00335         }
00336     }
00337     else
00338     {
00339         double* iter = dbl_array;
00340 
00341         for( size_t i = 0; i < num_vtx; ++i )
00342         {
00343             coordinates[i].set( iter );
00344             iter += 3;
00345         }
00346     }
00347 
00348     err.set_error( MBMesquite::MsqError::NO_ERROR );
00349 }
00350 
00351 void MsqMOAB::vertex_set_coordinates( MBMesquite::Mesh::VertexHandle vertex, const Vector3D& coords, MsqError& err )
00352 {
00353     moab::ErrorCode ierr;
00354     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00355     moab::EntityHandle* bh = reinterpret_cast< moab::EntityHandle* >( &vertex );
00356     double xyz[3]          = { coords[0], coords[1], coords[2] };
00357     ierr                   = meshInstance->set_coords( bh, 1, xyz );MB_CHK_ERR_RET( ierr );
00358     err.set_error( MBMesquite::MsqError::NO_ERROR );
00359 }
00360 
00361 // Each vertex has a byte-sized flag that can be used to store
00362 // flags.  This byte's value is neither set nor used by the mesh
00363 // implementation.  It is intended to be used by Mesquite algorithms.
00364 // Until a vertex's byte has been explicitly set, its value is 0.
00365 void MsqMOAB::vertex_set_byte( MBMesquite::Mesh::VertexHandle vertex, unsigned char byte, MsqError& err )
00366 {
00367     moab::ErrorCode ierr;
00368     int value = byte;
00369     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00370     moab::EntityHandle* bh = reinterpret_cast< moab::EntityHandle* >( &vertex );
00371     ierr                   = meshInstance->tag_set_data( byteTag, bh, 1, &value );MB_CHK_ERR_RET( ierr );
00372     err.set_error( MBMesquite::MsqError::NO_ERROR );
00373 }
00374 
00375 void MsqMOAB::vertices_set_byte( const VertexHandle* vert_array,
00376                                  const unsigned char* byte_array,
00377                                  size_t array_size,
00378                                  MsqError& err )
00379 {
00380     if( !array_size )
00381     {
00382         return;
00383     }
00384 
00385     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00386     std::vector< int > data( array_size );
00387     std::copy( byte_array, byte_array + array_size, data.begin() );
00388     moab::ErrorCode ierr;
00389     assert( sizeof( VertexHandle ) == sizeof( moab::EntityHandle ) );
00390     const moab::EntityHandle* arr = reinterpret_cast< const moab::EntityHandle* >( vert_array );
00391     ierr                          = meshInstance->tag_set_data( byteTag, arr, array_size, arrptr( data ) );MB_CHK_ERR_RET( ierr );
00392     err.set_error( MBMesquite::MsqError::NO_ERROR );
00393 }
00394 
00395 // Retrieve the byte value for the specified vertex or vertices.
00396 // The byte value is 0 if it has not yet been set via one of the
00397 // *_set_byte() functions.
00398 void MsqMOAB::vertex_get_byte( MBMesquite::Mesh::VertexHandle vertex, unsigned char* byte, MsqError& err )
00399 {
00400     moab::ErrorCode ierr;
00401     int value;
00402     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00403     moab::EntityHandle* bh = reinterpret_cast< moab::EntityHandle* >( &vertex );
00404     ierr                   = meshInstance->tag_get_data( byteTag, bh, 1, &value );MB_CHK_ERR_RET( ierr );
00405     *byte = value;
00406     err.set_error( MBMesquite::MsqError::NO_ERROR );
00407 }
00408 
00409 void MsqMOAB::vertices_get_byte( const VertexHandle* vert_array,
00410                                  unsigned char* byte_array,
00411                                  size_t array_size,
00412                                  MsqError& err )
00413 {
00414     if( !array_size )
00415     {
00416         return;
00417     }
00418 
00419     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00420     std::vector< int > data( array_size );
00421     moab::ErrorCode ierr;
00422     int* ptr = arrptr( data );
00423     assert( sizeof( VertexHandle ) == sizeof( moab::EntityHandle ) );
00424     const moab::EntityHandle* arr = reinterpret_cast< const moab::EntityHandle* >( vert_array );
00425     ierr                          = meshInstance->tag_get_data( byteTag, arr, array_size, ptr );MB_CHK_ERR_RET( ierr );
00426     std::copy( data.begin(), data.end(), byte_array );
00427     err.set_error( MBMesquite::MsqError::NO_ERROR );
00428 }
00429 
00430 //**************** Topology *****************
00431 
00432 void MsqMOAB::get_adjacent_entities( const moab::EntityHandle* source,
00433                                      size_t num_source,
00434                                      moab::EntityType target_type,
00435                                      std::vector< EntityHandle >& target,
00436                                      std::vector< size_t >& offsets,
00437                                      MsqError& err )
00438 {
00439     if( num_source == 0 )
00440     {
00441         target.clear();
00442         offsets.clear();
00443         offsets.reserve( 1 );
00444         offsets.push_back( 0 );
00445         return;
00446     }
00447 
00448     err.set_error( MBMesquite::MsqError::UNKNOWN_ERROR );
00449 
00450     moab::ErrorCode ierr;
00451     int num_adj = 0, num_offset = 0;
00452     unsigned iadjoff = 0;
00453 
00454     assert( sizeof( size_t ) >= sizeof( int ) );
00455     offsets.resize( num_source + 1 );
00456     int* ptr2;
00457     bool expand = false;
00458 
00459     if( sizeof( size_t ) > sizeof( int ) )
00460     {
00461         ptr2       = (int*)malloc( sizeof( int ) * ( num_source + 1 ) );
00462         expand     = true;
00463         num_offset = num_source + 1;
00464     }
00465     else
00466     {
00467         // sizeof(int) == sizeof(size_t)
00468         ptr2       = (int*)arrptr( offsets );
00469         num_offset = offsets.size();
00470     }
00471 
00472     // std::cout << "Target capacity: " << target.capacity() << " and num sources = " <<  num_source
00473     // << std::endl;
00474 
00475     assert( sizeof( moab::EntityHandle ) == sizeof( EntityHandle ) );
00476     bool have_adj = false;
00477 
00478     // If passed vector has allocated storage, try to use existing space
00479     if( target.capacity() >= num_source )
00480     {
00481         target.resize( target.capacity() );
00482         // target.clear();
00483         moab::EntityHandle* ptr = reinterpret_cast< moab::EntityHandle* >( arrptr( target ) );
00484 
00485         ptr2[0] = 0;
00486 
00487         for( unsigned is = 0; is < num_source; ++is )
00488         {
00489             moab::Range adjents;
00490 
00491             if( target_type == moab::MBVERTEX )
00492             {
00493                 ierr = meshInstance->get_adjacencies( &source[is], 1, 0, false, adjents,
00494                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00495             }
00496             else if( target_type == moab::MBEDGE )
00497             {
00498                 ierr = meshInstance->get_adjacencies( &source[is], 1, 1, true, adjents,
00499                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00500             }
00501             else if( target_type <= moab::MBPOLYGON )
00502             {
00503                 ierr = meshInstance->get_adjacencies( &source[is], 1, 2, true, adjents,
00504                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00505             }
00506             else if( target_type < moab::MBENTITYSET )
00507             {
00508                 ierr = meshInstance->get_adjacencies( &source[is], 1, 3, true, adjents,
00509                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00510             }
00511             else /* Either EntitySet or MBMaxType -- Failures */
00512             {
00513                 MSQ_SETERR( err )
00514                 ( process_itaps_error( 1 ),
00515                   MsqError::NOT_IMPLEMENTED );  // "Invalid Target entity type specified"
00516                 return;
00517             }
00518 
00519             ptr2[is + 1] = iadjoff + adjents.size();
00520 
00521             for( unsigned iadj = 0; iadj < adjents.size(); ++iadj, ++iadjoff )
00522             {
00523                 ptr[iadjoff] = adjents[iadj];
00524             }
00525 
00526             // target.push_back( static_cast<EntityHandle>(&adjents[iadj]) );
00527             // std::cout << "1. Source entity [ " << is << "]: n(adjacencies) = " << offsets[is+1]
00528             // << std::endl;
00529         }
00530 
00531         if( moab::MB_SUCCESS == ierr )
00532         {
00533             have_adj = true;
00534             num_adj  = iadjoff;
00535         }
00536     }
00537 
00538     // If implementation passed back a size, try that
00539     if( !have_adj && num_adj && (unsigned)num_adj > target.capacity() )
00540     {
00541         target.resize( num_adj );
00542         // int junk1 = target.capacity(), junk3 = offsets.size();
00543         moab::EntityHandle* ptr = reinterpret_cast< moab::EntityHandle* >( arrptr( target ) );
00544 
00545         ptr2[0] = 0;
00546 
00547         for( unsigned is = 0; is < num_source; ++is )
00548         {
00549             moab::Range adjents;
00550 
00551             if( target_type == moab::MBVERTEX )
00552             {
00553                 ierr = meshInstance->get_adjacencies( &source[is], 1, 0, false, adjents,
00554                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00555             }
00556             else if( target_type == moab::MBEDGE )
00557             {
00558                 ierr = meshInstance->get_adjacencies( &source[is], 1, 1, true, adjents,
00559                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00560             }
00561             else if( target_type <= moab::MBPOLYGON )
00562             {
00563                 ierr = meshInstance->get_adjacencies( &source[is], 1, 2, true, adjents,
00564                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00565             }
00566             else if( target_type < moab::MBENTITYSET )
00567             {
00568                 ierr = meshInstance->get_adjacencies( &source[is], 1, 3, true, adjents,
00569                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00570             }
00571             else /* Either EntitySet or MBMaxType -- Failures */
00572             {
00573                 MSQ_SETERR( err )
00574                 ( process_itaps_error( 1 ),
00575                   MsqError::NOT_IMPLEMENTED );  // "Invalid Target entity type specified"
00576                 return;
00577             }
00578 
00579             ptr2[is + 1] = iadjoff + adjents.size();
00580 
00581             for( unsigned iadj = 0; iadj < adjents.size(); ++iadj, ++iadjoff )
00582             {
00583                 ptr[iadjoff] = adjents[iadj];
00584             }
00585         }
00586 
00587         if( moab::MB_SUCCESS == ierr )
00588         {
00589             have_adj = true;
00590         }
00591     }
00592 
00593     // Try with empty sidl array, and copy into elements vector
00594     if( !have_adj )
00595     {
00596         // If implementation passed back a size, try that
00597         std::vector< moab::EntityHandle > mArray;
00598 
00599         ptr2[0] = iadjoff;
00600 
00601         for( unsigned is = 0; is < num_source; ++is )
00602         {
00603             moab::Range adjents;
00604 
00605             if( target_type == moab::MBVERTEX )
00606             {
00607                 ierr = meshInstance->get_adjacencies( &source[is], 1, 0, false, adjents,
00608                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00609             }
00610             else if( target_type == moab::MBEDGE )
00611             {
00612                 ierr = meshInstance->get_adjacencies( &source[is], 1, 1, true, adjents,
00613                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00614             }
00615             else if( target_type <= moab::MBPOLYGON )
00616             {
00617                 ierr = meshInstance->get_adjacencies( &source[is], 1, 2, true, adjents,
00618                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00619             }
00620             else if( target_type < moab::MBENTITYSET )
00621             {
00622                 ierr = meshInstance->get_adjacencies( &source[is], 1, 3, true, adjents,
00623                                                       moab::Interface::INTERSECT );  // MB_CHK_ERR_RET(ierr);
00624             }
00625             else /* Either EntitySet or MBMaxType -- Failures */
00626             {
00627                 MSQ_SETERR( err )
00628                 ( process_itaps_error( 1 ),
00629                   MsqError::NOT_IMPLEMENTED );  // "Invalid Target entity type specified"
00630                 return;
00631             }
00632 
00633             ptr2[is + 1] = iadjoff + adjents.size();
00634 
00635             for( unsigned iadj = 0; iadj < adjents.size(); ++iadj, ++iadjoff )
00636             {
00637                 mArray.push_back( adjents[iadj] );
00638             }
00639         }
00640 
00641         num_adj = mArray.size();
00642         target.resize( num_adj );
00643         std::copy( mArray.begin(), mArray.end(), reinterpret_cast< moab::EntityHandle* >( arrptr( target ) ) );
00644         mArray.clear();
00645     }
00646 
00647     if( expand )
00648     {
00649         for( size_t i = num_offset; i > 0; --i )
00650         {
00651             offsets[i - 1] = static_cast< size_t >( ptr2[i - 1] );
00652         }
00653 
00654         free( ptr2 );
00655     }
00656 
00657     // assert( (unsigned)num_offset == offsets.size() );
00658     err.set_error( MBMesquite::MsqError::NO_ERROR );
00659 }
00660 
00661 void MsqMOAB::vertices_get_attached_elements( const VertexHandle* vertices,
00662                                               size_t num_vertex,
00663                                               std::vector< ElementHandle >& elements,
00664                                               std::vector< size_t >& offsets,
00665                                               MsqError& err )
00666 {
00667     // moab::ErrorCode ierr;
00668     bool cont;
00669     assert( sizeof( EntityHandle ) == sizeof( moab::EntityHandle ) );
00670     const moab::EntityHandle* verts = reinterpret_cast< const moab::EntityHandle* >( vertices );
00671     get_adjacent_entities( verts, num_vertex, inputSetType, elements, offsets, err );MSQ_ERRRTN( err );
00672 
00673     moab::EntityHandle root_set = 0 /*meshInstance->get_root_set()*/;
00674 
00675     // Remove all elements not in inputSet
00676     if( root_set != inputSet )
00677     {
00678         std::vector< size_t >::iterator offset_iter = offsets.begin();
00679         size_t read_idx, write_idx;
00680         moab::EntityHandle* bh = reinterpret_cast< moab::EntityHandle* >( arrptr( elements ) );
00681 
00682         for( read_idx = write_idx = 0; read_idx < elements.size(); ++read_idx )
00683         {
00684             if( *offset_iter == read_idx )
00685             {
00686                 *offset_iter = write_idx;
00687                 ++offset_iter;
00688             }
00689 
00690             cont = meshInstance->contains_entities( inputSet, &bh[read_idx], 1 );
00691 
00692             if( cont )
00693             {
00694                 elements[write_idx++] = elements[read_idx];
00695             }
00696         }
00697 
00698         *offset_iter = write_idx;
00699         elements.resize( write_idx );
00700     }
00701 }
00702 
00703 //**************** Element Topology *****************
00704 
00705 /** Get connectivity
00706  *\param elements - Array of length num_elems containing elements
00707  *                  handles of elements for which connectivity is to
00708  *                  be queried.
00709  *\param vertices - Array of vertex handles in connectivity list.
00710  *\param offsets  - Indices into \ref vertex_handles, one per element
00711  */
00712 void MsqMOAB::elements_get_attached_vertices( const ElementHandle* elements,
00713                                               size_t num_elems,
00714                                               std::vector< VertexHandle >& vertices,
00715                                               std::vector< size_t >& offsets,
00716                                               MBMesquite::MsqError& err )
00717 {
00718     moab::ErrorCode ierr;
00719     assert( sizeof( moab::EntityHandle ) == sizeof( EntityHandle ) );
00720     const moab::EntityHandle* elems = reinterpret_cast< const moab::EntityHandle* >( elements );
00721     // get_adjacent_entities( elems, num_elems, moab::MBVERTEX, vertices, offsets, err );
00722     // MSQ_CHKERR(err);
00723     offsets.resize( num_elems + 1 );
00724     offsets[0] = 0;
00725     vertices.clear();
00726     std::vector< moab::EntityHandle > mbverts;
00727 
00728     for( unsigned ie = 0; ie < num_elems; ++ie )
00729     {
00730         const moab::EntityHandle* conn;
00731         int nconn;
00732         ierr = meshInstance->get_connectivity( elems[ie], conn, nconn, false );MB_CHK_ERR_RET( ierr );
00733         offsets[ie + 1] = offsets[ie] + nconn;
00734 
00735         for( int iconn = 0; iconn < nconn; ++iconn )
00736         {
00737             mbverts.push_back( conn[iconn] );
00738         }
00739     }
00740 
00741     vertices.resize( mbverts.size() );
00742     moab::EntityHandle* verts = reinterpret_cast< moab::EntityHandle* >( arrptr( vertices ) );
00743 
00744     for( size_t iverts = 0; iverts < mbverts.size(); ++iverts )
00745     {
00746         verts[iverts] = mbverts[iverts];
00747     }
00748 
00749     mbverts.clear();
00750     err.set_error( MBMesquite::MsqError::NO_ERROR );
00751 }
00752 
00753 void MsqMOAB::get_all_elements( std::vector< ElementHandle >& elements, MsqError& err )
00754 {
00755     moab::ErrorCode ierr;
00756 
00757     if( inputSetType == moab::MBMAXTYPE )
00758     {
00759         int num_vol, num_face;
00760 
00761         ierr = meshInstance->get_number_entities_by_dimension( inputSet, 2, num_face, false );MB_CHK_ERR_RET( ierr );
00762         ierr = meshInstance->get_number_entities_by_dimension( inputSet, 3, num_vol, false );MB_CHK_ERR_RET( ierr );
00763         elements.resize( num_face + num_vol );
00764 
00765         if( elements.empty() )
00766         {
00767             return;
00768         }
00769 
00770         moab::EntityHandle* ptr = reinterpret_cast< moab::EntityHandle* >( arrptr( elements ) );
00771 
00772         if( num_face )
00773         {
00774             std::vector< moab::EntityHandle > faces;
00775             ierr = meshInstance->get_entities_by_dimension( inputSet, 2, faces, false );MB_CHK_ERR_RET( ierr );
00776             assert( faces.size() - num_face == 0 );
00777             std::copy( faces.begin(), faces.end(), ptr );
00778         }
00779 
00780         if( num_vol )
00781         {
00782             ptr += num_face;
00783             std::vector< moab::EntityHandle > regions;
00784             ierr = meshInstance->get_entities_by_dimension( inputSet, 3, regions, false );MB_CHK_ERR_RET( ierr );
00785             assert( regions.size() - num_vol == 0 );
00786             std::copy( regions.begin(), regions.end(), ptr );
00787         }
00788     }
00789     else
00790     {
00791         int count;
00792         ierr = meshInstance->get_number_entities_by_type( inputSet, inputSetType, count, false );MB_CHK_ERR_RET( ierr );
00793 
00794         if( !count )
00795         {
00796             return;
00797         }
00798 
00799         elements.resize( count );
00800 
00801         moab::EntityHandle* ptr = reinterpret_cast< moab::EntityHandle* >( arrptr( elements ) );
00802         std::vector< moab::EntityHandle > entities;
00803         ierr = meshInstance->get_entities_by_type( inputSet, inputSetType, entities, false );MB_CHK_ERR_RET( ierr );
00804         assert( entities.size() - count == 0 );
00805         std::copy( entities.begin(), entities.end(), ptr );
00806     }
00807 
00808     err.set_error( MBMesquite::MsqError::NO_ERROR );
00809 }
00810 
00811 void MsqMOAB::get_all_vertices( std::vector< VertexHandle >& vertices, MsqError& err )
00812 {
00813     std::vector< ElementHandle > elems;
00814     get_all_elements( elems, err );MSQ_CHKERR( err );
00815 
00816     if( elems.empty() )
00817     {
00818         return;
00819     }
00820 
00821     std::vector< size_t > offsets;
00822     elements_get_attached_vertices( arrptr( elems ), elems.size(), vertices, offsets, err );MSQ_CHKERR( err );
00823 
00824     std::sort( vertices.begin(), vertices.end() );
00825     vertices.erase( std::unique( vertices.begin(), vertices.end() ), vertices.end() );
00826 }
00827 
00828 // Returns the topologies of the given entities.  The "entity_topologies"
00829 // array must be at least "num_elements" in size.
00830 void MsqMOAB::elements_get_topologies( const ElementHandle* element_handle_array,
00831                                        EntityTopology* element_topologies,
00832                                        size_t num_elements,
00833                                        MsqError& err )
00834 {
00835     if( !num_elements )
00836     {
00837         return;
00838     }
00839 
00840     assert( sizeof( ElementHandle ) == sizeof( moab::EntityHandle ) );
00841     const moab::EntityHandle* arr = reinterpret_cast< const moab::EntityHandle* >( element_handle_array );
00842 
00843     for( size_t i = 0; i < num_elements; ++i )
00844     {
00845         element_topologies[i] = topologyMap[meshInstance->type_from_handle( arr[i] )];
00846     }
00847 
00848     err.set_error( MBMesquite::MsqError::NO_ERROR );
00849 }
00850 
00851 //**************** Memory Management ****************
00852 // Tells the mesh that the client is finished with a given
00853 // entity handle.
00854 void MsqMOAB::release_entity_handles( const MBMesquite::Mesh::EntityHandle* /*handle_array*/,
00855                                       size_t /*num_handles*/,
00856                                       MsqError& /*err*/ )
00857 {
00858     // Do nothing
00859 }
00860 
00861 // Instead of deleting a Mesh when you think you are done,
00862 // call release().  In simple cases, the implementation could
00863 // just call the destructor.  More sophisticated implementations
00864 // may want to keep the Mesh object to live longer than Mesquite
00865 // is using it.
00866 void MsqMOAB::release() {}
00867 
00868 //**************** Tags ****************
00869 TagHandle MsqMOAB::tag_create( const std::string& name, TagType type, unsigned length, const void*, MsqError& err )
00870 {
00871     moab::DataType moab_type;
00872 
00873     switch( type )
00874     {
00875         case MBMesquite::Mesh::BYTE:
00876             moab_type = moab::MB_TYPE_OPAQUE;
00877             break;
00878 
00879         case MBMesquite::Mesh::INT:
00880             moab_type = moab::MB_TYPE_INTEGER;
00881             break;
00882 
00883         case MBMesquite::Mesh::DOUBLE:
00884             moab_type = moab::MB_TYPE_DOUBLE;
00885             break;
00886 
00887         case MBMesquite::Mesh::HANDLE:
00888             moab_type = moab::MB_TYPE_HANDLE;
00889             break;
00890 
00891         default:
00892             MSQ_SETERR( err )( "Invalid tag type", MsqError::INVALID_ARG );
00893             return 0;
00894     }
00895 
00896     moab::ErrorCode ierr;
00897     moab::Tag result = 0;
00898     ierr             = meshInstance->tag_get_handle( name.c_str(), length, moab_type, result,
00899                                                      moab::MB_TAG_CREAT | moab::MB_TAG_DENSE );MB_CHK_ERR_RET_VAL( ierr, result );
00900 
00901     err.set_error( MBMesquite::MsqError::NO_ERROR );
00902     return static_cast< TagHandle >( result );
00903 }
00904 
00905 void MsqMOAB::tag_delete( TagHandle handle, MsqError& err )
00906 {
00907     moab::ErrorCode ierr = meshInstance->tag_delete( static_cast< moab::Tag >( handle ) );MB_CHK_ERR_RET( ierr );
00908     err.set_error( MBMesquite::MsqError::NO_ERROR );
00909 }
00910 
00911 TagHandle MsqMOAB::tag_get( const std::string& name, MsqError& err )
00912 {
00913     moab::ErrorCode ierr;
00914     moab::Tag handle = 0;
00915     ierr             = meshInstance->tag_get_handle( name.c_str(), handle );MB_CHK_ERR_RET_VAL( ierr, handle );
00916 
00917     err.set_error( MBMesquite::MsqError::NO_ERROR );
00918     return static_cast< TagHandle >( handle );
00919 }
00920 
00921 void MsqMOAB::tag_properties( TagHandle handle,
00922                               std::string& name_out,
00923                               TagType& type_out,
00924                               unsigned& length_out,
00925                               MsqError& err )
00926 {
00927     std::string buffer;
00928     moab::ErrorCode ierr;
00929     moab::DataType itype;
00930     int size;
00931 
00932     moab::Tag th = static_cast< moab::Tag >( handle );
00933     ierr         = meshInstance->tag_get_name( th, buffer );MB_CHK_ERR_RET( ierr );
00934     ierr = meshInstance->tag_get_length( th, size );MB_CHK_ERR_RET( ierr );
00935     ierr = meshInstance->tag_get_data_type( th, itype );MB_CHK_ERR_RET( ierr );
00936 
00937     name_out   = buffer;
00938     length_out = size;
00939 
00940     switch( itype )
00941     {
00942         case moab::MB_TYPE_OPAQUE:
00943             type_out = MBMesquite::Mesh::BYTE;
00944             break;
00945 
00946         case moab::MB_TYPE_INTEGER:
00947             type_out = MBMesquite::Mesh::INT;
00948             break;
00949 
00950         case moab::MB_TYPE_DOUBLE:
00951             type_out = MBMesquite::Mesh::DOUBLE;
00952             break;
00953 
00954         case moab::MB_TYPE_HANDLE:
00955             type_out = MBMesquite::Mesh::HANDLE;
00956             break;
00957 
00958         default:
00959             MSQ_SETERR( err )( "Unsupported iMesh tag type", MsqError::NOT_IMPLEMENTED );
00960             return;
00961     }
00962 
00963     err.set_error( MBMesquite::MsqError::NO_ERROR );
00964 }
00965 
00966 void MsqMOAB::tag_set_element_data( TagHandle tag,
00967                                     size_t num_elems,
00968                                     const ElementHandle* array,
00969                                     const void* data,
00970                                     MsqError& err )
00971 {
00972     tag_set_data( tag, num_elems, array, data, err );
00973 }
00974 
00975 void MsqMOAB::tag_set_vertex_data( TagHandle tag,
00976                                    size_t num_elems,
00977                                    const VertexHandle* array,
00978                                    const void* data,
00979                                    MsqError& err )
00980 {
00981     tag_set_data( tag, num_elems, array, data, err );
00982 }
00983 
00984 void MsqMOAB::tag_set_data( TagHandle tag,
00985                             size_t num_elems,
00986                             const EntityHandle* array,
00987                             const void* data,
00988                             MsqError& err )
00989 {
00990     moab::ErrorCode ierr;
00991     int size;
00992     moab::Tag th = static_cast< moab::Tag >( tag );
00993     ierr         = meshInstance->tag_get_length( th, size );MB_CHK_ERR_RET( ierr );
00994 
00995     assert( sizeof( EntityHandle ) == sizeof( moab::EntityHandle ) );
00996     const moab::EntityHandle* arr = reinterpret_cast< const moab::EntityHandle* >( array );
00997     ierr                          = meshInstance->tag_set_data( th, arr, num_elems, data );MB_CHK_ERR_RET( ierr );
00998 
00999     err.set_error( MBMesquite::MsqError::NO_ERROR );
01000 }
01001 
01002 void MsqMOAB::tag_get_element_data( TagHandle tag,
01003                                     size_t num_elems,
01004                                     const ElementHandle* array,
01005                                     void* data,
01006                                     MsqError& err )
01007 {
01008     tag_get_data( tag, num_elems, array, data, err );
01009 }
01010 
01011 void MsqMOAB::tag_get_vertex_data( TagHandle tag,
01012                                    size_t num_verts,
01013                                    const VertexHandle* array,
01014                                    void* data,
01015                                    MsqError& err )
01016 {
01017     tag_get_data( tag, num_verts, array, data, err );
01018 }
01019 
01020 void MsqMOAB::tag_get_data( TagHandle tag, size_t num_elems, const EntityHandle* array, void* data, MsqError& err )
01021 {
01022     moab::ErrorCode ierr;
01023 
01024 #if IMESH_VERSION_ATLEAST( 1, 1 )
01025     void* ptr = data;
01026 #else
01027     char* ptr = static_cast< char* >( data );
01028 #endif
01029 
01030     assert( sizeof( EntityHandle ) == sizeof( moab::EntityHandle ) );
01031     moab::Tag th                  = static_cast< moab::Tag >( tag );
01032     const moab::EntityHandle* arr = reinterpret_cast< const moab::EntityHandle* >( array );
01033     ierr                          = meshInstance->tag_get_data( th, arr, num_elems, ptr );MB_CHK_ERR_RET( ierr );
01034 
01035     err.set_error( MBMesquite::MsqError::NO_ERROR );
01036 }
01037 
01038 }  // namespace MBMesquite
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines