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