|
MOAB: Mesh Oriented datABase
(version 5.2.1)
|
|
MOAB: Mesh Oriented datABase
(version 5.2.1)
|
00001 /** 00002 * MOAB, a Mesh-Oriented datABase, is a software component for creating, 00003 * storing and accessing finite element mesh data. 00004 * 00005 * Copyright 2004 Sandia Corporation. Under the terms of Contract 00006 * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government 00007 * retains certain rights in this software. 00008 * 00009 * This library is free software; you can redistribute it and/or 00010 * modify it under the terms of the GNU Lesser General Public 00011 * License as published by the Free Software Foundation; either 00012 * version 2.1 of the License, or (at your option) any later version. 00013 * 00014 */ 00015 00016 #include "moab/HalfFacetRep.hpp" 00017 #include "Internals.hpp" 00018 #include <iostream> 00019 #include <assert.h> 00020 #include <vector> 00021 #include <map> 00022 #include "MBTagConventions.hpp" 00023 #include "moab/ScdInterface.hpp" 00024 #ifdef MOAB_HAVE_MPI 00025 #include "moab/ParallelComm.hpp" 00026 #endif 00027 00028 namespace moab 00029 { 00030 00031 inline EntityHandle CREATE_HALFFACET( const unsigned lid, const EntityID id ) 00032 { 00033 assert( id <= MB_END_ID && lid < 6 ); 00034 return ( ( (HFacet)lid ) << MB_ID_WIDTH ) | id; 00035 } 00036 inline EntityID FID_FROM_HALFFACET( HFacet handle ) 00037 { 00038 return ( handle & MB_ID_MASK ); 00039 } 00040 inline int LID_FROM_HALFFACET( HFacet handle ) 00041 { 00042 return static_cast< int >( handle >> MB_ID_WIDTH ); 00043 } 00044 00045 HalfFacetRep::HalfFacetRep( Core* impl, ParallelComm* comm, moab::EntityHandle rset, bool filter_ghosts ) 00046 : thismeshtype( CURVE ), mb( impl ), pcomm( comm ), _rset( rset ), _filterghost( filter_ghosts ) 00047 { 00048 assert( NULL != impl ); 00049 mInitAHFmaps = false; 00050 chk_mixed = false; 00051 is_mixed = false; 00052 } 00053 00054 HalfFacetRep::~HalfFacetRep() {} 00055 00056 MESHTYPE HalfFacetRep::get_mesh_type( int nverts, int nedges, int nfaces, int ncells ) 00057 { 00058 MESHTYPE mesh_type = CURVE; 00059 00060 if( nverts && nedges && ( !nfaces ) && ( !ncells ) ) 00061 mesh_type = CURVE; 00062 else if( nverts && !nedges && nfaces && !ncells ) 00063 mesh_type = SURFACE; 00064 else if( nverts && nedges && nfaces && !ncells ) 00065 mesh_type = SURFACE_MIXED; 00066 else if( nverts && !nedges && !nfaces && ncells ) 00067 mesh_type = VOLUME; 00068 else if( nverts && nedges && !nfaces && ncells ) 00069 mesh_type = VOLUME_MIXED_1; 00070 else if( nverts && !nedges && nfaces && ncells ) 00071 mesh_type = VOLUME_MIXED_2; 00072 else if( nverts && nedges && nfaces && ncells ) 00073 mesh_type = VOLUME_MIXED; 00074 00075 return mesh_type; 00076 } 00077 00078 const HalfFacetRep::adj_matrix HalfFacetRep::adjMatrix[7] = { 00079 // Stores the adjacency matrix for each mesh type. 00080 // CURVE 00081 { { { 0, 1, 0, 0 }, { 1, 1, 0, 0 }, { 0, 0, 0, 0 }, { 0, 0, 0, 0 } } }, 00082 00083 // SURFACE 00084 { { { 0, 0, 1, 0 }, { 0, 0, 0, 0 }, { 1, 0, 1, 0 }, { 0, 0, 0, 0 } } }, 00085 00086 // SURFACE_MIXED 00087 { { { 0, 1, 1, 0 }, { 1, 1, 1, 0 }, { 1, 1, 1, 0 }, { 0, 0, 0, 0 } } }, 00088 00089 // VOLUME 00090 { { { 0, 0, 0, 1 }, { 0, 0, 0, 0 }, { 0, 0, 0, 0 }, { 1, 0, 0, 1 } } }, 00091 00092 // VOLUME_MIXED_1 00093 { { { 0, 1, 0, 1 }, { 1, 1, 0, 1 }, { 0, 0, 0, 0 }, { 1, 1, 0, 1 } } }, 00094 00095 // VOLUME_MIXED_2 00096 { { { 0, 0, 1, 1 }, { 0, 0, 0, 0 }, { 1, 0, 1, 1 }, { 1, 0, 1, 1 } } }, 00097 00098 // VOLUME_MIXED 00099 { { { 0, 1, 1, 1 }, { 1, 1, 1, 1 }, { 1, 1, 1, 1 }, { 1, 1, 1, 1 } } } 00100 }; 00101 00102 int HalfFacetRep::get_index_for_meshtype( MESHTYPE mesh_type ) 00103 { 00104 int index = 0; 00105 switch( mesh_type ) 00106 { 00107 case CURVE: 00108 index = 0; 00109 break; 00110 case SURFACE: 00111 index = 1; 00112 break; 00113 case SURFACE_MIXED: 00114 index = 2; 00115 break; 00116 case VOLUME: 00117 index = 3; 00118 break; 00119 case VOLUME_MIXED_1: 00120 index = 4; 00121 break; 00122 case VOLUME_MIXED_2: 00123 index = 5; 00124 break; 00125 case VOLUME_MIXED: 00126 index = 6; 00127 break; 00128 } 00129 return index; 00130 } 00131 00132 bool HalfFacetRep::check_mixed_entity_type() 00133 { 00134 if( !chk_mixed ) 00135 { 00136 chk_mixed = true; 00137 00138 ErrorCode error; 00139 Range felems, celems; 00140 00141 error = mb->get_entities_by_dimension( this->_rset, 2, felems );MB_CHK_ERR( error ); 00142 00143 if( felems.size() ) 00144 { 00145 Range tris, quad, poly; 00146 tris = felems.subset_by_type( MBTRI ); 00147 quad = felems.subset_by_type( MBQUAD ); 00148 poly = felems.subset_by_type( MBPOLYGON ); 00149 if( ( tris.size() && quad.size() ) || ( tris.size() && poly.size() ) || ( quad.size() && poly.size() ) ) 00150 is_mixed = true; 00151 if( poly.size() ) is_mixed = true; 00152 00153 if( is_mixed ) return is_mixed; 00154 } 00155 00156 error = mb->get_entities_by_dimension( this->_rset, 3, celems );MB_CHK_ERR( error ); 00157 if( celems.size() ) 00158 { 00159 Range tet, pyr, prism, hex, polyhed; 00160 tet = celems.subset_by_type( MBTET ); 00161 pyr = celems.subset_by_type( MBPYRAMID ); 00162 prism = celems.subset_by_type( MBPRISM ); 00163 hex = celems.subset_by_type( MBHEX ); 00164 polyhed = celems.subset_by_type( MBPOLYHEDRON ); 00165 if( ( tet.size() && pyr.size() ) || ( tet.size() && prism.size() ) || ( tet.size() && hex.size() ) || 00166 ( tet.size() && polyhed.size() ) || ( pyr.size() && prism.size() ) || ( pyr.size() && hex.size() ) || 00167 ( pyr.size() && polyhed.size() ) || ( prism.size() && hex.size() ) || 00168 ( prism.size() && polyhed.size() ) || ( hex.size() && polyhed.size() ) ) 00169 is_mixed = true; 00170 00171 if( polyhed.size() ) is_mixed = true; 00172 } 00173 00174 ScdInterface* scdi = NULL; 00175 error = mb->query_interface( scdi );MB_CHK_ERR( error ); 00176 if( scdi ) 00177 { 00178 Range boxes; 00179 error = scdi->find_boxes( boxes );MB_CHK_ERR( error ); 00180 00181 if( !boxes.empty() ) is_mixed = true; 00182 } 00183 } 00184 return is_mixed; 00185 } 00186 00187 /******************************************************* 00188 * initialize * 00189 ******************************************************/ 00190 00191 ErrorCode HalfFacetRep::initialize() 00192 { 00193 ErrorCode error; 00194 00195 if( !mInitAHFmaps ) 00196 { 00197 mInitAHFmaps = true; 00198 #ifdef MOAB_HAVE_MPI 00199 if( pcomm && _filterghost ) 00200 { 00201 moab::Range _averts, _aedgs, _afacs, _acels; 00202 error = mb->get_entities_by_dimension( this->_rset, 0, _averts, true );MB_CHK_ERR( error ); 00203 error = mb->get_entities_by_dimension( this->_rset, 1, _aedgs, true );MB_CHK_ERR( error ); 00204 error = mb->get_entities_by_dimension( this->_rset, 2, _afacs, true );MB_CHK_ERR( error ); 00205 error = mb->get_entities_by_dimension( this->_rset, 3, _acels, true );MB_CHK_ERR( error ); 00206 00207 // filter based on parallel status 00208 error = pcomm->filter_pstatus( _averts, PSTATUS_GHOST, PSTATUS_NOT, -1, &_verts );MB_CHK_ERR( error ); 00209 error = pcomm->filter_pstatus( _aedgs, PSTATUS_GHOST, PSTATUS_NOT, -1, &_edges );MB_CHK_ERR( error ); 00210 error = pcomm->filter_pstatus( _afacs, PSTATUS_GHOST, PSTATUS_NOT, -1, &_faces );MB_CHK_ERR( error ); 00211 error = pcomm->filter_pstatus( _acels, PSTATUS_GHOST, PSTATUS_NOT, -1, &_cells );MB_CHK_ERR( error ); 00212 } 00213 else 00214 { 00215 error = mb->get_entities_by_dimension( this->_rset, 0, _verts, true );MB_CHK_ERR( error ); 00216 error = mb->get_entities_by_dimension( this->_rset, 1, _edges, true );MB_CHK_ERR( error ); 00217 error = mb->get_entities_by_dimension( this->_rset, 2, _faces, true );MB_CHK_ERR( error ); 00218 error = mb->get_entities_by_dimension( this->_rset, 3, _cells, true );MB_CHK_ERR( error ); 00219 } 00220 #else 00221 error = mb->get_entities_by_dimension( this->_rset, 0, _verts, true );MB_CHK_ERR( error ); 00222 error = mb->get_entities_by_dimension( this->_rset, 1, _edges, true );MB_CHK_ERR( error ); 00223 error = mb->get_entities_by_dimension( this->_rset, 2, _faces, true );MB_CHK_ERR( error ); 00224 error = mb->get_entities_by_dimension( this->_rset, 3, _cells, true );MB_CHK_ERR( error ); 00225 00226 #endif 00227 00228 int nverts = _verts.size(); 00229 int nedges = _edges.size(); 00230 int nfaces = _faces.size(); 00231 int ncells = _cells.size(); 00232 00233 MESHTYPE mesh_type = get_mesh_type( nverts, nedges, nfaces, ncells ); 00234 thismeshtype = mesh_type; 00235 00236 // Initialize mesh type specific maps 00237 if( thismeshtype == CURVE ) 00238 { 00239 error = init_curve();MB_CHK_ERR( error ); 00240 } 00241 else if( thismeshtype == SURFACE ) 00242 { 00243 error = init_surface();MB_CHK_ERR( error ); 00244 } 00245 else if( thismeshtype == SURFACE_MIXED ) 00246 { 00247 error = init_curve();MB_CHK_ERR( error ); 00248 error = init_surface();MB_CHK_ERR( error ); 00249 } 00250 else if( thismeshtype == VOLUME ) 00251 { 00252 error = init_volume();MB_CHK_ERR( error ); 00253 } 00254 else if( thismeshtype == VOLUME_MIXED_1 ) 00255 { 00256 error = init_curve();MB_CHK_ERR( error ); 00257 error = init_volume();MB_CHK_ERR( error ); 00258 } 00259 else if( thismeshtype == VOLUME_MIXED_2 ) 00260 { 00261 error = init_surface();MB_CHK_ERR( error ); 00262 error = init_volume();MB_CHK_ERR( error ); 00263 } 00264 else if( thismeshtype == VOLUME_MIXED ) 00265 { 00266 error = init_curve();MB_CHK_ERR( error ); 00267 error = init_surface();MB_CHK_ERR( error ); 00268 error = init_volume();MB_CHK_ERR( error ); 00269 } 00270 } 00271 return MB_SUCCESS; 00272 } 00273 00274 ErrorCode HalfFacetRep::deinitialize() 00275 { 00276 return MB_SUCCESS; 00277 } 00278 00279 ErrorCode HalfFacetRep::init_curve() 00280 { 00281 ErrorCode error; 00282 00283 int nv = ID_FROM_HANDLE( *( _verts.end() - 1 ) ); 00284 int ne = ID_FROM_HANDLE( *( _edges.end() - 1 ) ); 00285 00286 v2hv.resize( nv, 0 ); 00287 sibhvs.resize( ne * 2, 0 ); 00288 00289 error = determine_sibling_halfverts( _verts, _edges );MB_CHK_ERR( error ); 00290 error = determine_incident_halfverts( _edges );MB_CHK_ERR( error ); 00291 00292 return MB_SUCCESS; 00293 } 00294 00295 ErrorCode HalfFacetRep::init_surface() 00296 { 00297 ErrorCode error; 00298 EntityType ftype = mb->type_from_handle( *_faces.begin() ); 00299 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 00300 00301 int nv = ID_FROM_HANDLE( *( _verts.end() - 1 ) ); 00302 int nf = ID_FROM_HANDLE( *( _faces.end() - 1 ) ); 00303 00304 v2he.resize( nv, 0 ); 00305 sibhes.resize( nf * nepf, 0 ); 00306 00307 // Construct ahf maps 00308 error = determine_sibling_halfedges( _faces );MB_CHK_ERR( error ); 00309 error = determine_incident_halfedges( _faces );MB_CHK_ERR( error ); 00310 00311 // Initialize queues for storing face and local id's during local search 00312 for( int i = 0; i < MAXSIZE; i++ ) 00313 { 00314 queue_fid[i] = 0; 00315 queue_lid[i] = 0; 00316 trackfaces[i] = 0; 00317 } 00318 00319 return MB_SUCCESS; 00320 } 00321 00322 ErrorCode HalfFacetRep::init_volume() 00323 { 00324 ErrorCode error; 00325 00326 // Initialize std::map between cell-types and their index in lConnMap3D 00327 cell_index[MBTET] = 0; 00328 cell_index[MBPYRAMID] = 1; 00329 cell_index[MBPRISM] = 2; 00330 cell_index[MBHEX] = 3; 00331 00332 int index = get_index_in_lmap( *_cells.begin() ); 00333 int nfpc = lConnMap3D[index].num_faces_in_cell; 00334 int nv = ID_FROM_HANDLE( *( _verts.end() - 1 ) ); 00335 int nc = ID_FROM_HANDLE( *( _cells.end() - 1 ) ); 00336 ; 00337 00338 v2hf.resize( nv, 0 ); 00339 sibhfs.resize( nc * nfpc, 0 ); 00340 00341 // Construct the maps 00342 error = determine_sibling_halffaces( _cells );MB_CHK_ERR( error ); 00343 error = determine_incident_halffaces( _cells );MB_CHK_ERR( error ); 00344 00345 // Initialize queues for storing face and local id's during local search 00346 for( int i = 0; i < MAXSIZE; i++ ) 00347 { 00348 Stkcells[i] = 0; 00349 cellq[i] = 0; 00350 trackcells[i] = 0; 00351 } 00352 00353 return MB_SUCCESS; 00354 } 00355 00356 ////////////////////////////////////////////////// 00357 ErrorCode HalfFacetRep::print_tags( int dim ) 00358 { 00359 if( dim == 1 ) 00360 { 00361 EntityHandle start_edge = *_edges.begin(); 00362 std::cout << "start_edge = " << start_edge << std::endl; 00363 std::cout << "<SIBHVS_EID,SIBHVS_LVID>" << std::endl; 00364 00365 for( Range::iterator i = _edges.begin(); i != _edges.end(); ++i ) 00366 { 00367 EntityHandle eid[2]; 00368 int lvid[2]; 00369 int eidx = ID_FROM_HANDLE( *i ) - 1; 00370 HFacet hf1 = sibhvs[2 * eidx]; 00371 HFacet hf2 = sibhvs[2 * eidx + 1]; 00372 eid[0] = fid_from_halfacet( hf1, MBEDGE ); 00373 eid[1] = fid_from_halfacet( hf2, MBEDGE ); 00374 lvid[0] = lid_from_halffacet( hf1 ); 00375 lvid[1] = lid_from_halffacet( hf2 ); 00376 std::cout << "Entity = " << *i << " :: <" << eid[0] << "," << lvid[0] << ">" 00377 << " " 00378 << "<" << eid[1] << "," << lvid[1] << ">" << std::endl; 00379 } 00380 00381 std::cout << "<V2HV_EID, V2HV_LVID>" << std::endl; 00382 00383 for( Range::iterator i = _verts.begin(); i != _verts.end(); ++i ) 00384 { 00385 int vidx = ID_FROM_HANDLE( *i ) - 1; 00386 HFacet hf = v2hv[vidx]; 00387 EntityHandle eid = fid_from_halfacet( hf, MBEDGE ); 00388 int lvid = lid_from_halffacet( hf ); 00389 std::cout << "Vertex = " << *i << " :: <" << eid << "," << lvid << ">" << std::endl; 00390 } 00391 } 00392 else if( dim == 2 ) 00393 { 00394 EntityType ftype = mb->type_from_handle( *_faces.begin() ); 00395 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 00396 EntityHandle start_face = *_faces.begin(); 00397 std::cout << "start_face = " << start_face << std::endl; 00398 std::cout << "<SIBHES_FID,SIBHES_LEID>" << std::endl; 00399 00400 for( Range::iterator i = _faces.begin(); i != _faces.end(); ++i ) 00401 { 00402 int fidx = ID_FROM_HANDLE( *i ) - 1; 00403 std::cout << "Entity = " << *i; 00404 for( int j = 0; j < nepf; j++ ) 00405 { 00406 HFacet hf = sibhes[nepf * fidx + j]; 00407 EntityHandle sib = fid_from_halfacet( hf, ftype ); 00408 int lid = lid_from_halffacet( hf ); 00409 std::cout << " :: <" << sib << "," << lid << ">" 00410 << " "; 00411 } 00412 std::cout << std::endl; 00413 } 00414 00415 std::cout << "<V2HE_FID, V2HE_LEID>" << std::endl; 00416 00417 for( Range::iterator i = _verts.begin(); i != _verts.end(); ++i ) 00418 { 00419 int vidx = ID_FROM_HANDLE( *i ) - 1; 00420 HFacet hf = v2he[vidx]; 00421 EntityHandle fid = fid_from_halfacet( hf, ftype ); 00422 int lid = lid_from_halffacet( hf ); 00423 std::cout << "Vertex = " << *i << " :: <" << fid << "," << lid << ">" << std::endl; 00424 } 00425 } 00426 else if( dim == 3 ) 00427 { 00428 EntityType ctype = mb->type_from_handle( *_cells.begin() ); 00429 00430 int index = get_index_in_lmap( *_cells.begin() ); 00431 int nfpc = lConnMap3D[index].num_faces_in_cell; 00432 EntityHandle start_cell = *_cells.begin(); 00433 std::cout << "start_cell = " << start_cell << std::endl; 00434 std::cout << "<SIBHES_CID,SIBHES_LFID>" << std::endl; 00435 00436 for( Range::iterator i = _cells.begin(); i != _cells.end(); ++i ) 00437 { 00438 int cidx = ID_FROM_HANDLE( *i ) - 1; 00439 std::cout << "Entity = " << *i; 00440 for( int j = 0; j < nfpc; j++ ) 00441 { 00442 HFacet hf = sibhfs[nfpc * cidx + j]; 00443 EntityHandle sib = fid_from_halfacet( hf, ctype ); 00444 int lid = lid_from_halffacet( hf ); 00445 std::cout << " :: <" << sib << "," << lid << ">" 00446 << " "; 00447 } 00448 std::cout << std::endl; 00449 } 00450 00451 std::cout << "<V2HF_CID, V2HF_LFID>" << std::endl; 00452 EntityHandle cid; 00453 int lid; 00454 00455 for( Range::iterator i = _verts.begin(); i != _verts.end(); ++i ) 00456 { 00457 int vidx = ID_FROM_HANDLE( *i ) - 1; 00458 HFacet hf = v2hf[vidx]; 00459 00460 if( hf == 0 && ( v2hfs.find( *i ) != v2hfs.end() ) ) 00461 { 00462 std::pair< std::multimap< EntityHandle, HFacet >::iterator, 00463 std::multimap< EntityHandle, HFacet >::iterator > 00464 it_hfs; 00465 it_hfs = v2hfs.equal_range( *i ); 00466 00467 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hfs.first; it != it_hfs.second; ++it ) 00468 { 00469 cid = fid_from_halfacet( it->second, ctype ); 00470 lid = lid_from_halffacet( hf ); 00471 00472 std::cout << "Vertex = " << *i << " :: <" << cid << "," << lid << ">" << std::endl; 00473 } 00474 } 00475 else 00476 { 00477 cid = fid_from_halfacet( hf, ctype ); 00478 lid = lid_from_halffacet( hf ); 00479 std::cout << "Vertex = " << *i << " :: <" << cid << "," << lid << ">" << std::endl; 00480 } 00481 } 00482 } 00483 return MB_SUCCESS; 00484 } 00485 00486 /********************************************************** 00487 * User interface for adjacency functions * 00488 ********************************************************/ 00489 00490 ErrorCode HalfFacetRep::get_adjacencies( const EntityHandle source_entity, const unsigned int target_dimension, 00491 std::vector< EntityHandle >& target_entities ) 00492 { 00493 00494 ErrorCode error; 00495 00496 unsigned int source_dimension = mb->dimension_from_handle( source_entity ); 00497 assert( ( source_dimension <= target_dimension ) || ( source_dimension > target_dimension ) ); 00498 00499 if( mInitAHFmaps == false ) 00500 { 00501 error = initialize();MB_CHK_ERR( error ); 00502 } 00503 00504 int mindex = get_index_for_meshtype( thismeshtype ); 00505 int adj_possible = adjMatrix[mindex].val[source_dimension][target_dimension]; 00506 00507 if( adj_possible ) 00508 { 00509 if( source_dimension < target_dimension ) 00510 { 00511 error = get_up_adjacencies( source_entity, target_dimension, target_entities );MB_CHK_ERR( error ); 00512 } 00513 else if( source_dimension == target_dimension ) 00514 { 00515 error = get_neighbor_adjacencies( source_entity, target_entities );MB_CHK_ERR( error ); 00516 } 00517 else 00518 { 00519 error = get_down_adjacencies( source_entity, target_dimension, target_entities );MB_CHK_ERR( error ); 00520 } 00521 } 00522 else 00523 return MB_SUCCESS; 00524 00525 return MB_SUCCESS; 00526 } 00527 00528 ErrorCode HalfFacetRep::get_up_adjacencies( EntityHandle ent, int out_dim, std::vector< EntityHandle >& adjents, 00529 std::vector< int >* lids ) 00530 { 00531 ErrorCode error; 00532 int in_dim = mb->dimension_from_handle( ent ); 00533 assert( ( in_dim >= 0 && in_dim <= 2 ) && ( out_dim > in_dim ) ); 00534 00535 if( in_dim == 0 ) 00536 { 00537 if( out_dim == 1 ) 00538 { 00539 error = get_up_adjacencies_1d( ent, adjents, lids );MB_CHK_ERR( error ); 00540 } 00541 else if( out_dim == 2 ) 00542 { 00543 error = get_up_adjacencies_vert_2d( ent, adjents );MB_CHK_ERR( error ); 00544 } 00545 else if( out_dim == 3 ) 00546 { 00547 error = get_up_adjacencies_vert_3d( ent, adjents );MB_CHK_ERR( error ); 00548 } 00549 } 00550 00551 else if( ( in_dim == 1 ) && ( out_dim == 2 ) ) 00552 { 00553 error = get_up_adjacencies_2d( ent, adjents, lids );MB_CHK_ERR( error ); 00554 } 00555 else if( ( in_dim == 1 ) && ( out_dim == 3 ) ) 00556 { 00557 error = get_up_adjacencies_edg_3d( ent, adjents, lids );MB_CHK_ERR( error ); 00558 } 00559 else if( ( in_dim == 2 ) && ( out_dim == 3 ) ) 00560 { 00561 error = get_up_adjacencies_face_3d( ent, adjents, lids );MB_CHK_ERR( error ); 00562 } 00563 return MB_SUCCESS; 00564 } 00565 00566 ErrorCode HalfFacetRep::get_neighbor_adjacencies( EntityHandle ent, std::vector< EntityHandle >& adjents ) 00567 { 00568 ErrorCode error; 00569 int in_dim = mb->dimension_from_handle( ent ); 00570 assert( in_dim >= 1 && in_dim <= 3 ); 00571 00572 if( in_dim == 1 ) 00573 { 00574 error = get_neighbor_adjacencies_1d( ent, adjents );MB_CHK_ERR( error ); 00575 } 00576 00577 else if( in_dim == 2 ) 00578 { 00579 error = get_neighbor_adjacencies_2d( ent, adjents );MB_CHK_ERR( error ); 00580 } 00581 else if( in_dim == 3 ) 00582 { 00583 error = get_neighbor_adjacencies_3d( ent, adjents );MB_CHK_ERR( error ); 00584 } 00585 return MB_SUCCESS; 00586 } 00587 00588 ErrorCode HalfFacetRep::get_down_adjacencies( EntityHandle ent, int out_dim, std::vector< EntityHandle >& adjents ) 00589 { 00590 ErrorCode error; 00591 int in_dim = mb->dimension_from_handle( ent ); 00592 assert( ( in_dim >= 2 && in_dim <= 3 ) && ( out_dim < in_dim ) ); 00593 00594 if( ( in_dim == 2 ) && ( out_dim == 1 ) ) 00595 { 00596 error = get_down_adjacencies_2d( ent, adjents );MB_CHK_ERR( error ); 00597 } 00598 else if( ( in_dim == 3 ) && ( out_dim == 1 ) ) 00599 { 00600 error = get_down_adjacencies_edg_3d( ent, adjents );MB_CHK_ERR( error ); 00601 } 00602 else if( ( in_dim == 3 ) && ( out_dim == 2 ) ) 00603 { 00604 error = get_down_adjacencies_face_3d( ent, adjents );MB_CHK_ERR( error ); 00605 } 00606 return MB_SUCCESS; 00607 } 00608 00609 ErrorCode HalfFacetRep::count_subentities( Range& edges, Range& faces, Range& cells, int* nedges, int* nfaces ) 00610 { 00611 ErrorCode error; 00612 if( edges.size() && !faces.size() && !cells.size() ) 00613 { 00614 nedges[0] = edges.size(); 00615 nfaces[0] = 0; 00616 } 00617 else if( faces.size() && !cells.size() ) 00618 { 00619 nedges[0] = find_total_edges_2d( faces ); 00620 nfaces[0] = 0; 00621 } 00622 else if( cells.size() ) 00623 { 00624 error = find_total_edges_faces_3d( cells, nedges, nfaces );MB_CHK_ERR( error ); 00625 } 00626 return MB_SUCCESS; 00627 } 00628 00629 /******************************************************** 00630 * 1D: sibhvs, v2hv, incident and neighborhood queries * 00631 *********************************************************/ 00632 ErrorCode HalfFacetRep::determine_sibling_halfverts( Range& verts, Range& edges ) 00633 { 00634 ErrorCode error; 00635 00636 // Step 1: Create an index list storing the starting position for each vertex 00637 int nv = verts.size(); 00638 std::vector< int > is_index( nv + 1 ); 00639 for( int i = 0; i < nv + 1; i++ ) 00640 is_index[i] = 0; 00641 00642 for( Range::iterator eid = edges.begin(); eid != edges.end(); ++eid ) 00643 { 00644 const EntityHandle* conn; 00645 int num_conn = 0; 00646 error = mb->get_connectivity( *eid, conn, num_conn, true );MB_CHK_ERR( error ); 00647 00648 int index = verts.index( conn[0] ); 00649 is_index[index + 1] += 1; 00650 index = verts.index( conn[1] ); 00651 is_index[index + 1] += 1; 00652 } 00653 is_index[0] = 0; 00654 00655 for( int i = 0; i < nv; i++ ) 00656 is_index[i + 1] = is_index[i] + is_index[i + 1]; 00657 00658 // Step 2: Define two arrays v2hv_eid, v2hv_lvid storing every half-facet on a vertex 00659 std::vector< EntityHandle > v2hv_map_eid( 2 * edges.size() ); 00660 std::vector< int > v2hv_map_lvid( 2 * edges.size() ); 00661 00662 for( Range::iterator eid = edges.begin(); eid != edges.end(); ++eid ) 00663 { 00664 const EntityHandle* conn; 00665 int num_conn = 0; 00666 error = mb->get_connectivity( *eid, conn, num_conn, true );MB_CHK_ERR( error ); 00667 00668 for( int j = 0; j < 2; j++ ) 00669 { 00670 int v = verts.index( conn[j] ); 00671 v2hv_map_eid[is_index[v]] = *eid; 00672 v2hv_map_lvid[is_index[v]] = j; 00673 is_index[v] += 1; 00674 } 00675 } 00676 00677 for( int i = nv - 2; i >= 0; i-- ) 00678 is_index[i + 1] = is_index[i]; 00679 is_index[0] = 0; 00680 00681 // Step 3: Fill up sibling half-verts map 00682 for( Range::iterator vid = verts.begin(); vid != verts.end(); ++vid ) 00683 { 00684 int v = verts.index( *vid ); 00685 int last = is_index[v + 1] - 1; 00686 if( last > is_index[v] ) 00687 { 00688 EntityHandle prev_eid = v2hv_map_eid[last]; 00689 int prev_lvid = v2hv_map_lvid[last]; 00690 00691 for( int i = is_index[v]; i <= last; i++ ) 00692 { 00693 EntityHandle cur_eid = v2hv_map_eid[i]; 00694 int cur_lvid = v2hv_map_lvid[i]; 00695 00696 int pidx = ID_FROM_HANDLE( prev_eid ) - 1; 00697 sibhvs[2 * pidx + prev_lvid] = create_halffacet( cur_eid, cur_lvid ); 00698 00699 prev_eid = cur_eid; 00700 prev_lvid = cur_lvid; 00701 } 00702 } 00703 } 00704 00705 return MB_SUCCESS; 00706 } 00707 00708 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 00709 ErrorCode HalfFacetRep::determine_incident_halfverts( Range& edges ) 00710 { 00711 ErrorCode error; 00712 00713 for( Range::iterator e_it = edges.begin(); e_it != edges.end(); ++e_it ) 00714 { 00715 EntityHandle cur_eid = *e_it; 00716 const EntityHandle* conn; 00717 int num_conn = 0; 00718 error = mb->get_connectivity( *e_it, conn, num_conn, true );MB_CHK_ERR( error ); 00719 00720 for( int i = 0; i < 2; ++i ) 00721 { 00722 EntityHandle v = conn[i]; 00723 int vidx = ID_FROM_HANDLE( v ) - 1; 00724 00725 HFacet hf = v2hv[vidx]; 00726 EntityHandle eid = fid_from_halfacet( hf, MBEDGE ); 00727 if( eid == 0 ) { v2hv[vidx] = create_halffacet( cur_eid, i ); } 00728 } 00729 } 00730 00731 return MB_SUCCESS; 00732 } 00733 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 00734 ErrorCode HalfFacetRep::get_up_adjacencies_1d( EntityHandle vid, std::vector< EntityHandle >& adjents, 00735 std::vector< int >* lvids ) 00736 { 00737 adjents.clear(); 00738 adjents.reserve( 20 ); 00739 00740 if( lvids != NULL ) lvids->reserve( 20 ); 00741 00742 int vidx = ID_FROM_HANDLE( vid ) - 1; 00743 HFacet hf = v2hv[vidx]; 00744 00745 EntityHandle start_eid = fid_from_halfacet( hf, MBEDGE ); 00746 int start_lid = lid_from_halffacet( hf ); 00747 00748 EntityHandle eid = start_eid; 00749 int lid = start_lid; 00750 00751 if( eid != 0 ) 00752 { 00753 adjents.push_back( eid ); 00754 if( lvids != NULL ) lvids->push_back( lid ); 00755 00756 while( eid != 0 ) 00757 { 00758 int eidx = ID_FROM_HANDLE( eid ) - 1; 00759 HFacet shf = sibhvs[2 * eidx + lid]; 00760 eid = fid_from_halfacet( shf, MBEDGE ); 00761 lid = lid_from_halffacet( shf ); 00762 00763 if( ( !eid ) || ( eid == start_eid ) ) break; 00764 00765 adjents.push_back( eid ); 00766 if( lvids != NULL ) lvids->push_back( lid ); 00767 } 00768 } 00769 00770 return MB_SUCCESS; 00771 } 00772 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 00773 ErrorCode HalfFacetRep::get_neighbor_adjacencies_1d( EntityHandle eid, std::vector< EntityHandle >& adjents ) 00774 { 00775 adjents.clear(); 00776 adjents.reserve( 20 ); 00777 00778 EntityHandle sibhv_eid; 00779 int sibhv_lid; 00780 int eidx = ID_FROM_HANDLE( eid ) - 1; 00781 00782 for( int lid = 0; lid < 2; ++lid ) 00783 { 00784 HFacet shf = sibhvs[2 * eidx + lid]; 00785 sibhv_eid = fid_from_halfacet( shf, MBEDGE ); 00786 sibhv_lid = lid_from_halffacet( shf ); 00787 00788 if( sibhv_eid != 0 ) 00789 { 00790 adjents.push_back( sibhv_eid ); 00791 00792 eidx = ID_FROM_HANDLE( sibhv_eid ) - 1; 00793 HFacet nhf = sibhvs[2 * eidx + sibhv_lid]; 00794 EntityHandle hv_eid = fid_from_halfacet( nhf, MBEDGE ); 00795 int hv_lid = lid_from_halffacet( nhf ); 00796 00797 while( hv_eid != 0 ) 00798 { 00799 if( hv_eid != eid ) adjents.push_back( hv_eid ); 00800 00801 eidx = ID_FROM_HANDLE( hv_eid ) - 1; 00802 HFacet hf = sibhvs[2 * eidx + hv_lid]; 00803 EntityHandle edge = fid_from_halfacet( hf, MBEDGE ); 00804 if( edge == sibhv_eid ) break; 00805 00806 hv_eid = edge; 00807 hv_lid = lid_from_halffacet( hf ); 00808 } 00809 } 00810 } 00811 00812 return MB_SUCCESS; 00813 } 00814 00815 /******************************************************* 00816 * 2D: sibhes, v2he, incident and neighborhood queries * 00817 ********************************************************/ 00818 const HalfFacetRep::LocalMaps2D HalfFacetRep::lConnMap2D[2] = { 00819 // Triangle 00820 { 3, { 1, 2, 0 }, { 2, 0, 1 } }, 00821 // Quad 00822 { 4, { 1, 2, 3, 0 }, { 3, 0, 1, 2 } } 00823 }; 00824 00825 /////////////////////////////////////////////////////////////////////////////////////////////////////////////// 00826 ErrorCode HalfFacetRep::determine_sibling_halfedges( Range& faces ) 00827 { 00828 ErrorCode error; 00829 EntityHandle start_face = *faces.begin(); 00830 EntityType ftype = mb->type_from_handle( start_face ); 00831 int nfaces = faces.size(); 00832 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 00833 00834 // Step 1: Create an index list storing the starting position for each vertex 00835 int nv = _verts.size(); 00836 std::vector< int > is_index( nv + 1 ); 00837 for( int i = 0; i < nv + 1; i++ ) 00838 is_index[i] = 0; 00839 00840 int index; 00841 00842 for( Range::iterator fid = faces.begin(); fid != faces.end(); ++fid ) 00843 { 00844 const EntityHandle* conn; 00845 error = mb->get_connectivity( *fid, conn, nepf, true );MB_CHK_ERR( error ); 00846 00847 for( int i = 0; i < nepf; i++ ) 00848 { 00849 index = _verts.index( conn[i] ); 00850 is_index[index + 1] += 1; 00851 } 00852 } 00853 is_index[0] = 0; 00854 00855 for( int i = 0; i < nv; i++ ) 00856 is_index[i + 1] = is_index[i] + is_index[i + 1]; 00857 00858 // Step 2: Define two arrays v2hv_eid, v2hv_lvid storing every half-facet on a vertex 00859 std::vector< EntityHandle > v2nv( nepf * nfaces ); 00860 std::vector< EntityHandle > v2he_map_fid( nepf * nfaces ); 00861 std::vector< int > v2he_map_leid( nepf * nfaces ); 00862 00863 for( Range::iterator fid = faces.begin(); fid != faces.end(); ++fid ) 00864 { 00865 const EntityHandle* conn; 00866 error = mb->get_connectivity( *fid, conn, nepf, true );MB_CHK_ERR( error ); 00867 00868 for( int j = 0; j < nepf; j++ ) 00869 { 00870 int v = _verts.index( conn[j] ); 00871 int nidx = lConnMap2D[ftype - 2].next[j]; 00872 v2nv[is_index[v]] = conn[nidx]; 00873 v2he_map_fid[is_index[v]] = *fid; 00874 v2he_map_leid[is_index[v]] = j; 00875 is_index[v] += 1; 00876 } 00877 } 00878 00879 for( int i = nv - 2; i >= 0; i-- ) 00880 is_index[i + 1] = is_index[i]; 00881 is_index[0] = 0; 00882 00883 // Step 3: Fill up sibling half-verts map 00884 for( Range::iterator fid = faces.begin(); fid != faces.end(); ++fid ) 00885 { 00886 const EntityHandle* conn; 00887 error = mb->get_connectivity( *fid, conn, nepf, true );MB_CHK_ERR( error ); 00888 00889 int fidx = ID_FROM_HANDLE( *fid ) - 1; 00890 for( int k = 0; k < nepf; k++ ) 00891 { 00892 HFacet hf = sibhes[nepf * fidx + k]; 00893 EntityHandle sibfid = fid_from_halfacet( hf, ftype ); 00894 00895 if( sibfid != 0 ) continue; 00896 00897 int nidx = lConnMap2D[ftype - 2].next[k]; 00898 int v = _verts.index( conn[k] ); 00899 int vn = _verts.index( conn[nidx] ); 00900 00901 EntityHandle first_fid = *fid; 00902 int first_leid = k; 00903 00904 EntityHandle prev_fid = *fid; 00905 int prev_leid = k; 00906 00907 for( index = is_index[vn]; index <= is_index[vn + 1] - 1; index++ ) 00908 { 00909 if( v2nv[index] == conn[k] ) 00910 { 00911 EntityHandle cur_fid = v2he_map_fid[index]; 00912 int cur_leid = v2he_map_leid[index]; 00913 00914 int pidx = ID_FROM_HANDLE( prev_fid ) - 1; 00915 sibhes[nepf * pidx + prev_leid] = create_halffacet( cur_fid, cur_leid ); 00916 00917 prev_fid = cur_fid; 00918 prev_leid = cur_leid; 00919 } 00920 } 00921 00922 for( index = is_index[v]; index <= is_index[v + 1] - 1; index++ ) 00923 { 00924 if( ( v2nv[index] == conn[nidx] ) && ( v2he_map_fid[index] != *fid ) ) 00925 { 00926 00927 EntityHandle cur_fid = v2he_map_fid[index]; 00928 int cur_leid = v2he_map_leid[index]; 00929 00930 int pidx = ID_FROM_HANDLE( prev_fid ) - 1; 00931 sibhes[nepf * pidx + prev_leid] = create_halffacet( cur_fid, cur_leid ); 00932 00933 prev_fid = cur_fid; 00934 prev_leid = cur_leid; 00935 } 00936 } 00937 00938 if( prev_fid != first_fid ) 00939 { 00940 00941 int pidx = ID_FROM_HANDLE( prev_fid ) - 1; 00942 sibhes[nepf * pidx + prev_leid] = create_halffacet( first_fid, first_leid ); 00943 } 00944 } 00945 } 00946 00947 return MB_SUCCESS; 00948 } 00949 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 00950 ErrorCode HalfFacetRep::determine_incident_halfedges( Range& faces ) 00951 { 00952 ErrorCode error; 00953 EntityType ftype = mb->type_from_handle( *faces.begin() ); 00954 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 00955 00956 std::vector< char > markEdges( nepf * faces.size(), 0 ); 00957 00958 for( Range::iterator it = faces.begin(); it != faces.end(); ++it ) 00959 { 00960 EntityHandle fid = *it; 00961 const EntityHandle* conn; 00962 error = mb->get_connectivity( fid, conn, nepf, true );MB_CHK_ERR( error ); 00963 00964 for( int i = 0; i < nepf; ++i ) 00965 { 00966 EntityHandle v = conn[i]; 00967 int vidx = ID_FROM_HANDLE( v ) - 1; 00968 HFacet hf = v2he[vidx]; 00969 00970 if( hf == 0 && ( v2hes.empty() || ( v2hes.find( v ) == v2hes.end() ) ) ) 00971 { 00972 // This is the first time a half-facet is assigned to a vertex. 00973 HFacet nwhf = 0; 00974 error = mark_halfedges( v, fid, i, faces, markEdges, nwhf );MB_CHK_ERR( error ); 00975 00976 if( nwhf == 0 ) nwhf = create_halffacet( fid, i ); 00977 00978 v2he[vidx] = nwhf; 00979 } 00980 else if( hf != 0 && !markEdges[nepf * faces.index( fid ) + i] ) 00981 { 00982 // This is the first time a non-manifold vertex is encountered. Copy the existing he 00983 // in v2he[v] to the multimap. 00984 v2hes.insert( std::pair< EntityHandle, HFacet >( v, hf ) ); 00985 HFacet nwhf = 0; 00986 error = mark_halfedges( v, fid, i, faces, markEdges, nwhf );MB_CHK_ERR( error ); 00987 00988 if( nwhf == 0 ) nwhf = create_halffacet( fid, i ); 00989 00990 v2hes.insert( std::pair< EntityHandle, HFacet >( v, nwhf ) ); 00991 v2he[vidx] = 0; 00992 } 00993 else if( hf == 0 && ( !v2hes.empty() ) && ( v2hes.find( v ) != v2hes.end() ) && 00994 !markEdges[nepf * faces.index( fid ) + i] ) 00995 { 00996 // This is check if reached if the vertex is non-manifold and has encountered a 00997 // half-facet to a new component. 00998 HFacet nwhf = 0; 00999 error = mark_halfedges( v, fid, i, faces, markEdges, nwhf );MB_CHK_ERR( error ); 01000 01001 if( nwhf == 0 ) nwhf = create_halffacet( fid, i ); 01002 01003 v2hes.insert( std::pair< EntityHandle, HFacet >( v, nwhf ) ); 01004 } 01005 } 01006 } 01007 01008 // error = print_tags(2); 01009 01010 return MB_SUCCESS; 01011 } 01012 01013 /////////////////////////////////////////////////////////////////// 01014 ErrorCode HalfFacetRep::mark_halfedges( EntityHandle vid, EntityHandle he_fid, int he_lid, Range& faces, 01015 std::vector< char >& markHEdgs, HFacet& bnd_hf ) 01016 { 01017 ErrorCode error; 01018 EntityType ftype = mb->type_from_handle( he_fid ); 01019 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01020 01021 int qsize = 0, count = -1; 01022 int num_qvals = 0; 01023 01024 error = gather_halfedges( vid, he_fid, he_lid, &qsize, &count );MB_CHK_ERR( error ); 01025 01026 while( num_qvals < qsize ) 01027 { 01028 EntityHandle curfid = queue_fid[num_qvals]; 01029 int curlid = queue_lid[num_qvals]; 01030 num_qvals += 1; 01031 01032 int fidx = ID_FROM_HANDLE( curfid ) - 1; 01033 01034 const EntityHandle* conn; 01035 error = mb->get_connectivity( curfid, conn, nepf, true );MB_CHK_ERR( error ); 01036 01037 if( !markHEdgs[nepf * faces.index( curfid ) + curlid] && ( conn[curlid] == vid ) ) 01038 { 01039 markHEdgs[nepf * faces.index( curfid ) + curlid] = 1; 01040 HFacet hf = sibhes[nepf * fidx + curlid]; 01041 EntityHandle sibfid = fid_from_halfacet( hf, ftype ); 01042 if( sibfid == 0 ) bnd_hf = create_halffacet( curfid, curlid ); 01043 } 01044 01045 EntityHandle he2_fid = 0; 01046 int he2_lid = 0; 01047 error = another_halfedge( vid, curfid, curlid, &he2_fid, &he2_lid );MB_CHK_ERR( error ); 01048 01049 if( !markHEdgs[nepf * faces.index( curfid ) + he2_lid] && ( conn[he2_lid] == vid ) ) 01050 { 01051 markHEdgs[nepf * faces.index( curfid ) + he2_lid] = 1; 01052 HFacet hf = sibhes[nepf * fidx + he2_lid]; 01053 EntityHandle sibfid = fid_from_halfacet( hf, ftype ); 01054 if( sibfid == 0 ) bnd_hf = create_halffacet( he2_fid, he2_lid ); 01055 } 01056 01057 bool val = find_match_in_array( he2_fid, trackfaces, count ); 01058 01059 if( val ) continue; 01060 01061 count += 1; 01062 trackfaces[count] = he2_fid; 01063 01064 error = get_up_adjacencies_2d( he2_fid, he2_lid, &qsize, &count );MB_CHK_ERR( error ); 01065 } 01066 01067 // Change the visited faces to false, also empty the queue 01068 for( int i = 0; i <= qsize; i++ ) 01069 { 01070 queue_fid[i] = 0; 01071 queue_lid[i] = 0; 01072 } 01073 01074 for( int i = 0; i <= count; i++ ) 01075 trackfaces[i] = 0; 01076 return MB_SUCCESS; 01077 } 01078 01079 /////////////////////////////////////////////////////////////////// 01080 ErrorCode HalfFacetRep::get_up_adjacencies_vert_2d( EntityHandle vid, std::vector< EntityHandle >& adjents ) 01081 { 01082 ErrorCode error; 01083 EntityType ftype = mb->type_from_handle( *_faces.begin() ); 01084 01085 int vidx = ID_FROM_HANDLE( vid ) - 1; 01086 HFacet hf = v2he[vidx]; 01087 01088 std::vector< EntityHandle > start_fids; 01089 std::vector< int > start_lids; 01090 01091 if( hf == 0 && ( v2hes.find( vid ) != v2hes.end() ) ) 01092 { 01093 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 01094 it_hes; 01095 it_hes = v2hes.equal_range( vid ); 01096 01097 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 01098 { 01099 start_fids.push_back( fid_from_halfacet( it->second, ftype ) ); 01100 start_lids.push_back( lid_from_halffacet( it->second ) ); 01101 } 01102 } 01103 else if( hf != 0 ) 01104 { 01105 start_fids.push_back( fid_from_halfacet( hf, ftype ) ); 01106 start_lids.push_back( lid_from_halffacet( hf ) ); 01107 } 01108 01109 if( start_fids.empty() ) return MB_SUCCESS; 01110 01111 int qsize = 0, count = -1; 01112 int num_qvals = 0; 01113 01114 adjents.reserve( (int)start_fids.size() ); 01115 01116 for( int i = 0; i < (int)start_fids.size(); i++ ) 01117 { 01118 adjents.push_back( start_fids[i] ); 01119 error = gather_halfedges( vid, start_fids[i], start_lids[i], &qsize, &count );MB_CHK_ERR( error ); 01120 } 01121 01122 while( num_qvals < qsize ) 01123 { 01124 EntityHandle curfid = queue_fid[num_qvals]; 01125 int curlid = queue_lid[num_qvals]; 01126 num_qvals += 1; 01127 01128 EntityHandle he2_fid = 0; 01129 int he2_lid = 0; 01130 error = another_halfedge( vid, curfid, curlid, &he2_fid, &he2_lid );MB_CHK_ERR( error ); 01131 01132 bool val = find_match_in_array( he2_fid, trackfaces, count ); 01133 01134 if( val ) continue; 01135 01136 count += 1; 01137 trackfaces[count] = he2_fid; 01138 01139 error = get_up_adjacencies_2d( he2_fid, he2_lid, &qsize, &count );MB_CHK_ERR( error ); 01140 01141 adjents.push_back( he2_fid ); 01142 } 01143 01144 // Change the visited faces to false, also empty the queue 01145 for( int i = 0; i <= qsize; i++ ) 01146 { 01147 queue_fid[i] = 0; 01148 queue_lid[i] = 0; 01149 } 01150 01151 for( int i = 0; i <= count; i++ ) 01152 trackfaces[i] = 0; 01153 01154 return MB_SUCCESS; 01155 } 01156 01157 ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 01158 ErrorCode HalfFacetRep::get_up_adjacencies_2d( EntityHandle eid, std::vector< EntityHandle >& adjents, 01159 std::vector< int >* leids ) 01160 { 01161 01162 // Given an explicit edge eid, find the incident faces. 01163 ErrorCode error; 01164 EntityHandle he_fid = 0; 01165 int he_lid = 0; 01166 01167 // Step 1: Given an explicit edge, find a corresponding half-edge from the surface mesh. 01168 bool found = find_matching_halfedge( eid, &he_fid, &he_lid ); 01169 01170 // Step 2: If there is a corresponding half-edge, collect all sibling half-edges and store the 01171 // incident faces. 01172 if( found ) 01173 { 01174 error = get_up_adjacencies_2d( he_fid, he_lid, true, adjents, leids );MB_CHK_ERR( error ); 01175 } 01176 01177 return MB_SUCCESS; 01178 } 01179 01180 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 01181 ErrorCode HalfFacetRep::get_up_adjacencies_2d( EntityHandle fid, int leid, bool add_inent, 01182 std::vector< EntityHandle >& adj_ents, std::vector< int >* adj_leids, 01183 std::vector< int >* adj_orients ) 01184 { 01185 // Given an implicit half-edge <fid, leid>, find the incident half-edges. 01186 ErrorCode error; 01187 01188 EntityType ftype = mb->type_from_handle( fid ); 01189 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01190 01191 if( !fid ) return MB_FAILURE; 01192 adj_ents.reserve( 20 ); 01193 01194 bool local_id = false; 01195 bool orient = false; 01196 if( adj_leids != NULL ) 01197 { 01198 local_id = true; 01199 adj_leids->reserve( 20 ); 01200 } 01201 if( adj_orients != NULL ) 01202 { 01203 orient = true; 01204 adj_orients->reserve( 20 ); 01205 } 01206 01207 if( add_inent ) 01208 { 01209 adj_ents.push_back( fid ); 01210 if( local_id ) adj_leids->push_back( leid ); 01211 } 01212 01213 EntityHandle fedge[2] = { 0, 0 }; 01214 01215 if( orient ) 01216 { 01217 // get connectivity and match their directions 01218 const EntityHandle* fid_conn; 01219 error = mb->get_connectivity( fid, fid_conn, nepf, true );MB_CHK_ERR( error ); 01220 01221 int nidx = lConnMap2D[ftype - 2].next[leid]; 01222 fedge[0] = fid_conn[leid]; 01223 fedge[1] = fid_conn[nidx]; 01224 } 01225 01226 int fidx = ID_FROM_HANDLE( fid ) - 1; 01227 HFacet hf = sibhes[nepf * fidx + leid]; 01228 EntityHandle curfid = fid_from_halfacet( hf, ftype ); 01229 int curlid = lid_from_halffacet( hf ); 01230 01231 while( ( curfid != fid ) && ( curfid != 0 ) ) 01232 { // Should not go into the loop when no sibling exists 01233 adj_ents.push_back( curfid ); 01234 01235 if( local_id ) adj_leids->push_back( curlid ); 01236 01237 if( orient ) 01238 { 01239 // get connectivity and match their directions 01240 const EntityHandle* conn; 01241 error = mb->get_connectivity( curfid, conn, nepf, true );MB_CHK_ERR( error ); 01242 01243 int nidx = lConnMap2D[ftype - 2].next[curlid]; 01244 01245 if( ( fedge[0] == conn[curlid] ) && ( fedge[1] == conn[nidx] ) ) 01246 adj_orients->push_back( 1 ); 01247 else if( ( fedge[1] == conn[curlid] ) && ( fedge[0] == conn[nidx] ) ) 01248 adj_orients->push_back( 0 ); 01249 } 01250 01251 int cidx = ID_FROM_HANDLE( curfid ) - 1; 01252 hf = sibhes[nepf * cidx + curlid]; 01253 curfid = fid_from_halfacet( hf, ftype ); 01254 curlid = lid_from_halffacet( hf ); 01255 } 01256 01257 return MB_SUCCESS; 01258 } 01259 01260 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 01261 ErrorCode HalfFacetRep::get_up_adjacencies_2d( EntityHandle fid, int lid, int* qsize, int* count ) 01262 { 01263 EntityType ftype = mb->type_from_handle( fid ); 01264 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01265 01266 int fidx = ID_FROM_HANDLE( fid ) - 1; 01267 HFacet hf = sibhes[nepf * fidx + lid]; 01268 EntityHandle curfid = fid_from_halfacet( hf, ftype ); 01269 int curlid = lid_from_halffacet( hf ); 01270 01271 if( curfid == 0 ) 01272 { 01273 int index = 0; 01274 bool found_ent = find_match_in_array( fid, queue_fid, qsize[0] - 1, true, &index ); 01275 01276 if( ( !found_ent ) || ( ( found_ent ) && ( queue_lid[index] != lid ) ) ) 01277 { 01278 queue_fid[qsize[0]] = fid; 01279 queue_lid[qsize[0]] = lid; 01280 qsize[0] += 1; 01281 } 01282 } 01283 01284 while( ( curfid != fid ) && ( curfid != 0 ) ) 01285 { 01286 bool val = find_match_in_array( curfid, trackfaces, count[0] ); 01287 01288 if( !val ) 01289 { 01290 queue_fid[qsize[0]] = curfid; 01291 queue_lid[qsize[0]] = curlid; 01292 qsize[0] += 1; 01293 } 01294 01295 int cidx = ID_FROM_HANDLE( curfid ) - 1; 01296 hf = sibhes[nepf * cidx + curlid]; 01297 curfid = fid_from_halfacet( hf, ftype ); 01298 curlid = lid_from_halffacet( hf ); 01299 } 01300 01301 return MB_SUCCESS; 01302 } 01303 01304 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 01305 bool HalfFacetRep::find_matching_halfedge( EntityHandle eid, EntityHandle* hefid, int* helid ) 01306 { 01307 ErrorCode error; 01308 EntityType ftype = mb->type_from_handle( *_faces.begin() ); 01309 01310 const EntityHandle* conn; 01311 int num_conn = 0; 01312 error = mb->get_connectivity( eid, conn, num_conn, true );MB_CHK_ERR( error ); 01313 01314 EntityHandle vid = conn[0]; 01315 int vidx = ID_FROM_HANDLE( conn[0] ) - 1; 01316 HFacet hf = v2he[vidx]; 01317 01318 if( hf == 0 ) 01319 { 01320 vidx = ID_FROM_HANDLE( conn[1] ) - 1; 01321 hf = v2he[vidx]; 01322 01323 if( hf == 0 ) // The edge is either a dangling edge or attached to two non-manifold 01324 // vertices 01325 return MB_SUCCESS; 01326 01327 vid = conn[1]; 01328 } 01329 01330 EntityHandle fid = fid_from_halfacet( hf, ftype ); 01331 int lid = lid_from_halffacet( hf ); 01332 01333 bool found = false; 01334 int qsize = 0, count = -1; 01335 01336 error = gather_halfedges( vid, fid, lid, &qsize, &count );MB_CHK_ERR( error ); 01337 01338 found = collect_and_compare( vid, conn, &qsize, &count, hefid, helid );MB_CHK_ERR( error ); 01339 01340 // Change the visited faces to false 01341 for( int i = 0; i < qsize; i++ ) 01342 { 01343 queue_fid[i] = 0; 01344 queue_lid[i] = 0; 01345 } 01346 01347 for( int i = 0; i <= count; i++ ) 01348 trackfaces[i] = 0; 01349 01350 return found; 01351 } 01352 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 01353 ErrorCode HalfFacetRep::gather_halfedges( EntityHandle vid, EntityHandle he_fid, int he_lid, int* qsize, int* count ) 01354 { 01355 ErrorCode error; 01356 EntityHandle he2_fid = 0; 01357 int he2_lid = 0; 01358 01359 error = another_halfedge( vid, he_fid, he_lid, &he2_fid, &he2_lid );MB_CHK_ERR( error ); 01360 01361 queue_fid[*qsize] = he_fid; 01362 queue_lid[*qsize] = he_lid; 01363 *qsize += 1; 01364 01365 queue_fid[*qsize] = he2_fid; 01366 queue_lid[*qsize] = he2_lid; 01367 *qsize += 1; 01368 01369 *count += 1; 01370 trackfaces[*count] = he_fid; 01371 01372 error = get_up_adjacencies_2d( he_fid, he_lid, qsize, count );MB_CHK_ERR( error ); 01373 error = get_up_adjacencies_2d( he2_fid, he2_lid, qsize, count );MB_CHK_ERR( error ); 01374 01375 return MB_SUCCESS; 01376 } 01377 01378 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 01379 ErrorCode HalfFacetRep::another_halfedge( EntityHandle vid, EntityHandle he_fid, int he_lid, EntityHandle* he2_fid, 01380 int* he2_lid ) 01381 { 01382 ErrorCode error; 01383 EntityType ftype = mb->type_from_handle( he_fid ); 01384 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01385 01386 const EntityHandle* conn; 01387 error = mb->get_connectivity( he_fid, conn, nepf, true );MB_CHK_ERR( error ); 01388 01389 *he2_fid = he_fid; 01390 if( conn[he_lid] == vid ) 01391 *he2_lid = lConnMap2D[ftype - 2].prev[he_lid]; 01392 else 01393 *he2_lid = lConnMap2D[ftype - 2].next[he_lid]; 01394 01395 return MB_SUCCESS; 01396 } 01397 01398 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 01399 bool HalfFacetRep::collect_and_compare( const EntityHandle vid, const EntityHandle* edg_vert, int* qsize, int* count, 01400 EntityHandle* he_fid, int* he_lid ) 01401 { 01402 ErrorCode error; 01403 EntityType ftype = mb->type_from_handle( *_faces.begin() ); 01404 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01405 01406 bool found = false; 01407 int num_qvals = 0, counter = 0; 01408 01409 while( num_qvals < *qsize && counter < MAXSIZE ) 01410 { 01411 EntityHandle curfid = queue_fid[num_qvals]; 01412 int curlid = queue_lid[num_qvals]; 01413 num_qvals += 1; 01414 01415 const EntityHandle* conn; 01416 error = mb->get_connectivity( curfid, conn, nepf, true );MB_CHK_ERR( error ); 01417 01418 int id = lConnMap2D[ftype - 2].next[curlid]; 01419 if( ( ( conn[curlid] == edg_vert[0] ) && ( conn[id] == edg_vert[1] ) ) || 01420 ( ( conn[curlid] == edg_vert[1] ) && ( conn[id] == edg_vert[0] ) ) ) 01421 { 01422 *he_fid = curfid; 01423 *he_lid = curlid; 01424 found = true; 01425 break; 01426 } 01427 01428 bool val = find_match_in_array( curfid, trackfaces, count[0] ); 01429 01430 if( val ) continue; 01431 01432 count[0] += 1; 01433 trackfaces[*count] = curfid; 01434 01435 EntityHandle he2_fid; 01436 int he2_lid; 01437 error = another_halfedge( vid, curfid, curlid, &he2_fid, &he2_lid );MB_CHK_ERR( error ); 01438 error = get_up_adjacencies_2d( he2_fid, he2_lid, qsize, count );MB_CHK_ERR( error ); 01439 01440 counter += 1; 01441 } 01442 return found; 01443 } 01444 01445 /////////////////////////////////////////////////////////////////////////////////////////////////// 01446 ErrorCode HalfFacetRep::get_neighbor_adjacencies_2d( EntityHandle fid, std::vector< EntityHandle >& adjents ) 01447 { 01448 ErrorCode error; 01449 01450 if( fid != 0 ) 01451 { 01452 EntityType ftype = mb->type_from_handle( fid ); 01453 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01454 01455 for( int lid = 0; lid < nepf; ++lid ) 01456 { 01457 error = get_up_adjacencies_2d( fid, lid, false, adjents );MB_CHK_ERR( error ); 01458 } 01459 } 01460 01461 return MB_SUCCESS; 01462 } 01463 01464 ///////////////////////////////////////////////////////////////////////////////////////////////// 01465 ErrorCode HalfFacetRep::get_down_adjacencies_2d( EntityHandle fid, std::vector< EntityHandle >& adjents ) 01466 { 01467 // Returns explicit edges, if any, of the face 01468 ErrorCode error; 01469 adjents.reserve( 10 ); 01470 EntityType ftype = mb->type_from_handle( fid ); 01471 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01472 01473 const EntityHandle* conn; 01474 error = mb->get_connectivity( fid, conn, nepf, true );MB_CHK_ERR( error ); 01475 01476 std::vector< EntityHandle > temp; 01477 01478 // Loop over 2 vertices 01479 for( int i = 0; i < 2; i++ ) 01480 { 01481 // Choose the two adjacent vertices for a triangle, and two opposite vertices for a quad 01482 int l; 01483 if( ftype == MBTRI ) 01484 l = i; 01485 else 01486 l = 2 * i; 01487 01488 // Get the current, next and prev vertices 01489 int nidx = lConnMap2D[ftype - 2].next[l]; 01490 int pidx = lConnMap2D[ftype - 2].prev[l]; 01491 EntityHandle v = conn[l]; 01492 EntityHandle vnext = conn[nidx]; 01493 EntityHandle vprev = conn[pidx]; 01494 01495 // Get incident edges on v 01496 error = get_up_adjacencies_1d( v, temp );MB_CHK_ERR( error ); 01497 01498 // Loop over the incident edges and check if its end vertices match those in the face 01499 for( int k = 0; k < (int)temp.size(); k++ ) 01500 { 01501 const EntityHandle* econn; 01502 int num_conn = 0; 01503 error = mb->get_connectivity( temp[k], econn, num_conn, true );MB_CHK_ERR( error ); 01504 01505 if( ( econn[0] == v && econn[1] == vnext ) || ( econn[0] == v && econn[1] == vprev ) || 01506 ( econn[0] == vnext && econn[1] == v ) || ( econn[0] == vprev && econn[1] == v ) ) 01507 { 01508 bool found = false; 01509 for( int j = 0; j < (int)adjents.size(); j++ ) 01510 { 01511 if( adjents[j] == temp[k] ) 01512 { 01513 found = true; 01514 break; 01515 } 01516 } 01517 if( !found ) adjents.push_back( temp[k] ); 01518 } 01519 } 01520 } 01521 01522 return MB_SUCCESS; 01523 } 01524 01525 //////////////////////////////////////////////////////////////////////////////////////////////// 01526 int HalfFacetRep::find_total_edges_2d( Range& faces ) 01527 { 01528 ErrorCode error; 01529 EntityType ftype = mb->type_from_handle( *faces.begin() ); 01530 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01531 int nfaces = faces.size(); 01532 01533 int total_edges = nepf * nfaces; 01534 01535 std::vector< int > trackF( total_edges, 0 ); 01536 std::vector< EntityHandle > adj_fids; 01537 std::vector< int > adj_lids; 01538 01539 for( Range::iterator f = faces.begin(); f != faces.end(); ++f ) 01540 { 01541 for( int l = 0; l < nepf; l++ ) 01542 { 01543 01544 adj_fids.clear(); 01545 adj_lids.clear(); 01546 01547 int id = nepf * ( faces.index( *f ) ) + l; 01548 if( !trackF[id] ) 01549 { 01550 error = get_up_adjacencies_2d( *f, l, false, adj_fids, &adj_lids );MB_CHK_ERR( error ); 01551 01552 total_edges -= adj_fids.size(); 01553 01554 for( int i = 0; i < (int)adj_fids.size(); i++ ) 01555 trackF[nepf * ( faces.index( adj_fids[i] ) ) + adj_lids[i]] = 1; 01556 }; 01557 }; 01558 }; 01559 01560 return total_edges; 01561 } 01562 01563 ErrorCode HalfFacetRep::get_face_edges( EntityHandle fid, std::vector< EntityHandle >& edges ) 01564 { 01565 ErrorCode error; 01566 edges.clear(); 01567 01568 EntityType ftype = mb->type_from_handle( fid ); 01569 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 01570 01571 std::vector< EntityHandle > conn; 01572 error = mb->get_connectivity( &fid, 1, conn );MB_CHK_ERR( error ); 01573 01574 for( int i = 0; i < nepf; i++ ) 01575 { 01576 EntityHandle v0 = conn[i]; 01577 EntityHandle v1 = conn[lConnMap2D[ftype - 2].next[i]]; 01578 01579 std::vector< EntityHandle > e0, e1, ecom; 01580 error = get_up_adjacencies_1d( v0, e0 );MB_CHK_ERR( error ); 01581 error = get_up_adjacencies_1d( v1, e1 );MB_CHK_ERR( error ); 01582 01583 std::sort( e0.begin(), e0.end() ); 01584 std::sort( e1.begin(), e1.end() ); 01585 std::set_intersection( e0.begin(), e0.end(), e1.begin(), e1.end(), std::back_inserter( ecom ) ); 01586 assert( ecom.size() == 1 || ecom.size() == 0 ); 01587 if( ecom.size() == 0 ) 01588 edges.push_back( 0 ); 01589 else 01590 edges.push_back( ecom[0] ); 01591 } 01592 01593 return MB_SUCCESS; 01594 } 01595 01596 /******************************************************* 01597 * 3D: sibhfs, v2hf, incident and neighborhood queries * 01598 ********************************************************/ 01599 01600 int HalfFacetRep::get_index_in_lmap( EntityHandle cid ) 01601 { 01602 EntityType type = mb->type_from_handle( cid ); 01603 int index = cell_index.find( type )->second; 01604 return index; 01605 } 01606 01607 const HalfFacetRep::LocalMaps3D HalfFacetRep::lConnMap3D[4] = { 01608 // Tet 01609 { 4, 01610 6, 01611 4, 01612 { 3, 3, 3, 3 }, 01613 { { 0, 1, 3 }, { 1, 2, 3 }, { 2, 0, 3 }, { 0, 2, 1 } }, 01614 { 3, 3, 3, 3 }, 01615 { { 0, 2, 3 }, { 0, 1, 3 }, { 1, 2, 3 }, { 0, 1, 2 } }, 01616 { { 0, 1 }, { 1, 2 }, { 2, 0 }, { 0, 3 }, { 1, 3 }, { 2, 3 } }, 01617 { { 3, 0 }, { 3, 1 }, { 3, 2 }, { 0, 2 }, { 0, 1 }, { 1, 2 } }, 01618 { { 0, 4, 3 }, { 1, 5, 4 }, { 2, 3, 5 }, { 2, 1, 0 } }, 01619 { { -1, 0, 2, 3 }, { 0, -1, 1, 4 }, { 2, 1, -1, 5 }, { 3, 4, 5, -1 } }, 01620 { 3, 0, 1, 2 }, 01621 { 0, 1 }, 01622 { { 3, 1, 2, 3 }, { 2, 2, 3 }, { 1, 3 } } }, 01623 01624 // Pyramid: Note: In MOAB pyramid follows the CGNS convention. Look up src/MBCNArrays.hpp 01625 { 5, 01626 8, 01627 5, 01628 { 4, 3, 3, 3, 3 }, 01629 { { 0, 3, 2, 1 }, { 0, 1, 4 }, { 1, 2, 4 }, { 2, 3, 4 }, { 3, 0, 4 } }, 01630 { 3, 3, 3, 3, 4 }, 01631 { { 0, 1, 4 }, { 0, 1, 2 }, { 0, 2, 3 }, { 0, 3, 4 }, { 1, 2, 3, 4 } }, 01632 { { 0, 1 }, { 1, 2 }, { 2, 3 }, { 3, 0 }, { 0, 4 }, { 1, 4 }, { 2, 4 }, { 3, 4 } }, 01633 { { 0, 1 }, { 0, 2 }, { 0, 3 }, { 0, 4 }, { 1, 4 }, { 1, 2 }, { 2, 3 }, { 3, 4 } }, 01634 { { 3, 2, 1, 0 }, { 0, 5, 4 }, { 1, 6, 5 }, { 2, 7, 6 }, { 3, 4, 7 } }, 01635 { { -1, 0, -1, 3, 4 }, { 0, -1, 1, -1, 5 }, { -1, 1, -1, 2, 6 }, { 3, -1, 2, -1, 7 }, { 4, 5, 6, 7, -1 } }, 01636 { 3, 4, 2, 0 }, 01637 { 0, 4 }, 01638 { { 4, 0, 1, 2, 3 }, { 2, 1, 3 }, { 2, 1, 3 } } }, 01639 01640 // Prism 01641 { 6, 01642 9, 01643 5, 01644 { 4, 4, 4, 3, 3 }, 01645 { { 0, 1, 4, 3 }, { 1, 2, 5, 4 }, { 0, 3, 5, 2 }, { 0, 2, 1 }, { 3, 4, 5 } }, 01646 { 3, 3, 3, 3, 3, 3 }, 01647 { { 0, 2, 3 }, { 0, 1, 3 }, { 1, 2, 3 }, { 0, 2, 4 }, { 0, 1, 4 }, { 1, 4, 2 } }, 01648 { { 0, 1 }, { 1, 2 }, { 2, 0 }, { 0, 3 }, { 1, 4 }, { 2, 5 }, { 3, 4 }, { 4, 5 }, { 5, 3 } }, 01649 { { 0, 3 }, { 1, 3 }, { 2, 3 }, { 0, 2 }, { 0, 1 }, { 1, 2 }, { 0, 4 }, { 1, 4 }, { 2, 4 } }, 01650 { { 0, 4, 6, 3 }, { 1, 5, 7, 4 }, { 2, 3, 8, 5 }, { 2, 1, 0 }, { 6, 7, 8 } }, 01651 { { -1, 0, 2, 3, -1, -1 }, 01652 { 0, -1, 1, -1, 4, -1 }, 01653 { 2, 1, -1, -1, -1, 5 }, 01654 { 3, -1, -1, -1, 6, 8 }, 01655 { -1, 4, -1, 6, -1, 7 }, 01656 { -1, -1, 5, 8, 7, -1 } }, 01657 { 4, 0, 5, 4, 1 }, 01658 { 0, 5 }, 01659 { { 3, 1, 2, 3 }, { 3, 2, 4, 3 }, { 2, 3, 1 }, { 1, 2 } } }, 01660 01661 // Hex 01662 { 8, 01663 12, 01664 6, 01665 { 4, 4, 4, 4, 4, 4 }, 01666 { { 0, 1, 5, 4 }, { 1, 2, 6, 5 }, { 2, 3, 7, 6 }, { 3, 0, 4, 7 }, { 0, 3, 2, 1 }, { 4, 5, 6, 7 } }, 01667 { 3, 3, 3, 3, 3, 3, 3, 3 }, 01668 { { 0, 3, 4 }, { 0, 1, 4 }, { 1, 2, 4 }, { 2, 3, 4 }, { 0, 3, 5 }, { 0, 1, 5 }, { 1, 2, 5 }, { 2, 3, 5 } }, 01669 { { 0, 1 }, 01670 { 1, 2 }, 01671 { 2, 3 }, 01672 { 3, 0 }, 01673 { 0, 4 }, 01674 { 1, 5 }, 01675 { 2, 6 }, 01676 { 3, 7 }, 01677 { 4, 5 }, 01678 { 5, 6 }, 01679 { 6, 7 }, 01680 { 7, 4 } }, 01681 { { 0, 4 }, 01682 { 1, 4 }, 01683 { 2, 4 }, 01684 { 3, 4 }, 01685 { 0, 3 }, 01686 { 0, 1 }, 01687 { 1, 2 }, 01688 { 2, 3 }, 01689 { 0, 5 }, 01690 { 1, 5 }, 01691 { 2, 5 }, 01692 { 3, 5 } }, 01693 { { 0, 5, 8, 4 }, { 1, 6, 9, 5 }, { 2, 7, 10, 6 }, { 3, 4, 11, 7 }, { 3, 2, 1, 0 }, { 8, 9, 10, 11 } }, 01694 { { -1, 0, -1, 3, 4, -1, -1, -1 }, 01695 { 0, -1, 1, -1, -1, 5, -1, -1 }, 01696 { -1, 1, -1, 2, -1, -1, 6, -1 }, 01697 { 3, -1, 2, -1, -1, -1, -1, 7 }, 01698 { 4, -1, -1, -1, -1, 8, -1, 11 }, 01699 { -1, 5, -1, -1, 8, -1, 9, -1 }, 01700 { -1, -1, 6, -1, -1, 9, -1, 10 }, 01701 { -1, -1, -1, 7, 11, -1, 10, -1 } }, 01702 { 4, 0, 2, 5, 7 }, 01703 { 0, 6 }, 01704 { { 3, 1, 3, 4 }, { 3, 1, 3, 6 }, { 3, 1, 4, 6 }, { 3, 3, 6, 4 } } } 01705 }; 01706 01707 ////////////////////////////////////////////////////////////////////////////////// 01708 ErrorCode HalfFacetRep::determine_sibling_halffaces( Range& cells ) 01709 { 01710 ErrorCode error; 01711 EntityHandle start_cell = *cells.begin(); 01712 EntityType ctype = mb->type_from_handle( start_cell ); 01713 int index = get_index_in_lmap( start_cell ); 01714 int nvpc = lConnMap3D[index].num_verts_in_cell; 01715 int nfpc = lConnMap3D[index].num_faces_in_cell; 01716 01717 // Step 1: Create an index list storing the starting position for each vertex 01718 int nv = _verts.size(); 01719 std::vector< int > is_index( nv + 1 ); 01720 for( int i = 0; i < nv + 1; i++ ) 01721 is_index[i] = 0; 01722 01723 int vindex; 01724 01725 for( Range::iterator cid = cells.begin(); cid != cells.end(); ++cid ) 01726 { 01727 const EntityHandle* conn; 01728 error = mb->get_connectivity( *cid, conn, nvpc, true );MB_CHK_ERR( error ); 01729 01730 for( int i = 0; i < nfpc; ++i ) 01731 { 01732 int nvF = lConnMap3D[index].hf2v_num[i]; 01733 EntityHandle v = 0; 01734 for( int k = 0; k < nvF; k++ ) 01735 { 01736 int id = lConnMap3D[index].hf2v[i][k]; 01737 if( v <= conn[id] ) v = conn[id]; 01738 } 01739 vindex = _verts.index( v ); 01740 is_index[vindex + 1] += 1; 01741 } 01742 } 01743 is_index[0] = 0; 01744 01745 for( int i = 0; i < nv; i++ ) 01746 is_index[i + 1] = is_index[i] + is_index[i + 1]; 01747 01748 // Step 2: Define four arrays v2hv_eid, v2hv_lvid storing every half-facet on a vertex 01749 std::vector< EntityHandle > v2oe_v1( is_index[nv] ); 01750 std::vector< EntityHandle > v2oe_v2( is_index[nv] ); 01751 std::vector< EntityHandle > v2hf_map_cid( is_index[nv] ); 01752 std::vector< int > v2hf_map_lfid( is_index[nv] ); 01753 01754 for( Range::iterator cid = cells.begin(); cid != cells.end(); ++cid ) 01755 { 01756 const EntityHandle* conn; 01757 error = mb->get_connectivity( *cid, conn, nvpc, true );MB_CHK_ERR( error ); 01758 01759 for( int i = 0; i < nfpc; i++ ) 01760 { 01761 int nvF = lConnMap3D[index].hf2v_num[i]; 01762 std::vector< EntityHandle > vs( nvF ); 01763 EntityHandle vmax = 0; 01764 int lv = -1; 01765 for( int k = 0; k < nvF; k++ ) 01766 { 01767 int id = lConnMap3D[index].hf2v[i][k]; 01768 vs[k] = conn[id]; 01769 if( vmax <= conn[id] ) 01770 { 01771 vmax = conn[id]; 01772 lv = k; 01773 } 01774 } 01775 01776 int nidx = lConnMap2D[nvF - 3].next[lv]; 01777 int pidx = lConnMap2D[nvF - 3].prev[lv]; 01778 01779 int v = _verts.index( vmax ); 01780 v2oe_v1[is_index[v]] = vs[nidx]; 01781 v2oe_v2[is_index[v]] = vs[pidx]; 01782 v2hf_map_cid[is_index[v]] = *cid; 01783 v2hf_map_lfid[is_index[v]] = i; 01784 is_index[v] += 1; 01785 } 01786 } 01787 01788 for( int i = nv - 2; i >= 0; i-- ) 01789 is_index[i + 1] = is_index[i]; 01790 is_index[0] = 0; 01791 01792 // Step 3: Fill up sibling half-verts map 01793 for( Range::iterator cid = cells.begin(); cid != cells.end(); ++cid ) 01794 { 01795 const EntityHandle* conn; 01796 error = mb->get_connectivity( *cid, conn, nvpc, true );MB_CHK_ERR( error ); 01797 01798 int cidx = ID_FROM_HANDLE( *cid ) - 1; 01799 for( int i = 0; i < nfpc; i++ ) 01800 { 01801 HFacet hf = sibhfs[nfpc * cidx + i]; 01802 EntityHandle sibcid = fid_from_halfacet( hf, ctype ); 01803 01804 if( sibcid != 0 ) continue; 01805 01806 int nvF = lConnMap3D[index].hf2v_num[i]; 01807 std::vector< EntityHandle > vs( nvF ); 01808 EntityHandle vmax = 0; 01809 int lv = -1; 01810 for( int k = 0; k < nvF; k++ ) 01811 { 01812 int id = lConnMap3D[index].hf2v[i][k]; 01813 vs[k] = conn[id]; 01814 if( vmax <= conn[id] ) 01815 { 01816 vmax = conn[id]; 01817 lv = k; 01818 } 01819 } 01820 if( lv < 0 ) MB_SET_ERR( MB_FAILURE, "did not find local vertex " ); 01821 int nidx = lConnMap2D[nvF - 3].next[lv]; 01822 int pidx = lConnMap2D[nvF - 3].prev[lv]; 01823 01824 int v = _verts.index( vmax ); 01825 EntityHandle v1 = vs[pidx]; 01826 EntityHandle v2 = vs[nidx]; 01827 01828 for( int ind = is_index[v]; ind <= is_index[v + 1] - 1; ind++ ) 01829 { 01830 if( ( v2oe_v1[ind] == v1 ) && ( v2oe_v2[ind] == v2 ) ) 01831 { 01832 // Map to opposite hf 01833 EntityHandle cur_cid = v2hf_map_cid[ind]; 01834 int cur_lfid = v2hf_map_lfid[ind]; 01835 01836 sibhfs[nfpc * cidx + i] = create_halffacet( cur_cid, cur_lfid ); 01837 01838 int scidx = ID_FROM_HANDLE( cur_cid ) - 1; 01839 sibhfs[nfpc * scidx + cur_lfid] = create_halffacet( *cid, i ); 01840 } 01841 } 01842 } 01843 } 01844 01845 return MB_SUCCESS; 01846 } 01847 01848 ErrorCode HalfFacetRep::determine_incident_halffaces( Range& cells ) 01849 { 01850 ErrorCode error; 01851 int index = get_index_in_lmap( *cells.begin() ); 01852 int nvpc = lConnMap3D[index].num_verts_in_cell; 01853 01854 std::multimap< EntityHandle, EntityHandle > comps; 01855 HFacet nwhf; 01856 01857 for( Range::iterator cid = cells.begin(); cid != cells.end(); ++cid ) 01858 { 01859 EntityHandle cell = *cid; 01860 const EntityHandle* conn; 01861 error = mb->get_connectivity( *cid, conn, nvpc, true );MB_CHK_ERR( error ); 01862 01863 for( int i = 0; i < nvpc; ++i ) 01864 { 01865 EntityHandle v = conn[i]; 01866 int vidx = ID_FROM_HANDLE( v ) - 1; 01867 HFacet hf = v2hf[vidx]; 01868 01869 bool found = find_cell_in_component( v, cell, comps ); 01870 01871 if( hf == 0 && !found && ( v2hfs.empty() || ( v2hfs.find( v ) == v2hfs.end() ) ) ) 01872 { 01873 nwhf = 0; 01874 error = add_cells_of_single_component( v, cell, lConnMap3D[index].v2hf[i][0], comps, nwhf );MB_CHK_ERR( error ); 01875 01876 v2hf[vidx] = nwhf; 01877 } 01878 else if( hf != 0 && !found ) 01879 { 01880 nwhf = 0; 01881 error = add_cells_of_single_component( v, cell, lConnMap3D[index].v2hf[i][0], comps, nwhf );MB_CHK_ERR( error ); 01882 01883 v2hfs.insert( std::pair< EntityHandle, HFacet >( v, hf ) ); 01884 v2hfs.insert( std::pair< EntityHandle, HFacet >( v, nwhf ) ); 01885 v2hf[vidx] = 0; 01886 } 01887 else if( hf == 0 && !found && ( !v2hfs.empty() ) && ( v2hfs.find( v ) != v2hfs.end() ) ) 01888 { 01889 nwhf = 0; 01890 error = add_cells_of_single_component( v, cell, lConnMap3D[index].v2hf[i][0], comps, nwhf );MB_CHK_ERR( error ); 01891 v2hfs.insert( std::pair< EntityHandle, HFacet >( v, nwhf ) ); 01892 } 01893 } 01894 } 01895 01896 // error = print_tags(3); 01897 01898 return MB_SUCCESS; 01899 } 01900 01901 ErrorCode HalfFacetRep::determine_border_vertices( Range& cells, Tag isborder ) 01902 { 01903 ErrorCode error; 01904 EntityHandle start_cell = *cells.begin(); 01905 EntityType ctype = mb->type_from_handle( *cells.begin() ); 01906 int index = get_index_in_lmap( start_cell ); 01907 int nvpc = lConnMap3D[index].num_verts_in_cell; 01908 int nfpc = lConnMap3D[index].num_faces_in_cell; 01909 01910 int val = 1; 01911 01912 for( Range::iterator t = cells.begin(); t != cells.end(); ++t ) 01913 { 01914 01915 const EntityHandle* conn; 01916 error = mb->get_connectivity( *t, conn, nvpc, true );MB_CHK_ERR( error ); 01917 01918 int cidx = ID_FROM_HANDLE( *t ) - 1; 01919 for( int i = 0; i < nfpc; ++i ) 01920 { 01921 HFacet hf = sibhfs[nfpc * cidx + i]; 01922 EntityHandle sib_cid = fid_from_halfacet( hf, ctype ); 01923 01924 if( sib_cid == 0 ) 01925 { 01926 int nvF = lConnMap3D[index].hf2v_num[i]; 01927 01928 for( int j = 0; j < nvF; ++j ) 01929 { 01930 int ind = lConnMap3D[index].hf2v[i][j]; 01931 error = mb->tag_set_data( isborder, &conn[ind], 1, &val );MB_CHK_ERR( error ); 01932 } 01933 } 01934 } 01935 } 01936 01937 return MB_SUCCESS; 01938 } 01939 01940 ///////////////////////////////////////////////////////////////////////// 01941 ErrorCode HalfFacetRep::add_cells_of_single_component( EntityHandle vid, EntityHandle curcid, int curlid, 01942 std::multimap< EntityHandle, EntityHandle >& comps, HFacet& hf ) 01943 { 01944 ErrorCode error; 01945 EntityType ctype = mb->type_from_handle( curcid ); 01946 int index = get_index_in_lmap( curcid ); 01947 int nvpc = lConnMap3D[index].num_verts_in_cell; 01948 int nfpc = lConnMap3D[index].num_faces_in_cell; 01949 01950 int Stksize = 0, count = -1; 01951 Stkcells[0] = curcid; 01952 01953 hf = create_halffacet( curcid, curlid ); 01954 01955 EntityHandle cur_cid; 01956 while( Stksize >= 0 ) 01957 { 01958 cur_cid = Stkcells[Stksize]; 01959 Stksize -= 1; 01960 01961 bool found = find_match_in_array( cur_cid, trackcells, count ); 01962 if( !found ) 01963 { 01964 count += 1; 01965 trackcells[count] = cur_cid; 01966 01967 // Add the current cell 01968 comps.insert( std::pair< EntityHandle, EntityHandle >( vid, cur_cid ) ); 01969 } 01970 01971 // Connectivity of the cell 01972 const EntityHandle* conn; 01973 error = mb->get_connectivity( cur_cid, conn, nvpc, true );MB_CHK_ERR( error ); 01974 01975 // Local id of vid in the cell and the half-faces incident on it 01976 int lv = -1; 01977 for( int i = 0; i < nvpc; ++i ) 01978 { 01979 if( conn[i] == vid ) 01980 { 01981 lv = i; 01982 break; 01983 } 01984 } 01985 if( lv < 0 ) MB_SET_ERR( MB_FAILURE, "did not find local vertex " ); 01986 int nhf_thisv = lConnMap3D[index].v2hf_num[lv]; 01987 int cidx = ID_FROM_HANDLE( cur_cid ) - 1; 01988 01989 // Add new cells into the stack 01990 EntityHandle ngb; 01991 HFacet hf_ngb; 01992 for( int i = 0; i < nhf_thisv; ++i ) 01993 { 01994 int ind = lConnMap3D[index].v2hf[lv][i]; 01995 hf_ngb = sibhfs[nfpc * cidx + ind]; 01996 ngb = fid_from_halfacet( hf_ngb, ctype ); 01997 01998 if( ngb ) 01999 { 02000 bool found_ent = find_match_in_array( ngb, trackcells, count ); 02001 02002 if( !found_ent ) 02003 { 02004 Stksize += 1; 02005 Stkcells[Stksize] = ngb; 02006 } 02007 } 02008 else 02009 hf = create_halffacet( cur_cid, ind ); 02010 } 02011 } 02012 02013 // Change the visited faces to false 02014 for( int i = 0; i < Stksize; i++ ) 02015 Stkcells[i] = 0; 02016 02017 for( int i = 0; i <= count; i++ ) 02018 trackcells[i] = 0; 02019 02020 return MB_SUCCESS; 02021 } 02022 02023 bool HalfFacetRep::find_cell_in_component( EntityHandle vid, EntityHandle cell, 02024 std::multimap< EntityHandle, EntityHandle >& comps ) 02025 { 02026 bool found = false; 02027 02028 if( comps.empty() ) return found; 02029 02030 std::pair< std::multimap< EntityHandle, EntityHandle >::iterator, 02031 std::multimap< EntityHandle, EntityHandle >::iterator > 02032 rit; 02033 02034 rit = comps.equal_range( vid ); 02035 02036 for( std::multimap< EntityHandle, EntityHandle >::iterator it = rit.first; it != rit.second; ++it ) 02037 { 02038 if( it->second == cell ) 02039 { 02040 found = true; 02041 break; 02042 } 02043 } 02044 02045 return found; 02046 } 02047 02048 //////////////////////////////////////////////////////////////////////// 02049 ErrorCode HalfFacetRep::get_up_adjacencies_vert_3d( EntityHandle vid, std::vector< EntityHandle >& adjents ) 02050 { 02051 ErrorCode error; 02052 adjents.reserve( 20 ); 02053 EntityType ctype = mb->type_from_handle( *_cells.begin() ); 02054 02055 // Obtain a half-face/s incident on v 02056 int vidx = ID_FROM_HANDLE( vid ) - 1; 02057 HFacet hf = v2hf[vidx]; 02058 02059 std::vector< EntityHandle > start_cells; 02060 if( hf == 0 && ( v2hfs.find( vid ) != v2hfs.end() ) ) // Vertex is non-manifold 02061 { 02062 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02063 it_hes; 02064 it_hes = v2hfs.equal_range( vid ); 02065 02066 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02067 { 02068 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02069 } 02070 } 02071 else if( hf != 0 ) 02072 start_cells.push_back( fid_from_halfacet( hf, ctype ) ); 02073 02074 if( start_cells.empty() ) return MB_SUCCESS; 02075 02076 int index = get_index_in_lmap( start_cells[0] ); 02077 int nvpc = lConnMap3D[index].num_verts_in_cell; 02078 int nfpc = lConnMap3D[index].num_faces_in_cell; 02079 02080 for( int i = 0; i < (int)start_cells.size(); i++ ) 02081 cellq[i] = start_cells[i]; 02082 02083 int qsize = start_cells.size(); 02084 EntityHandle cur_cid; 02085 int num_qvals = 0; 02086 02087 while( num_qvals < qsize ) 02088 { 02089 cur_cid = cellq[num_qvals]; 02090 num_qvals += 1; 02091 02092 // Add the current cell to output adj vector 02093 adjents.push_back( cur_cid ); 02094 02095 // Connectivity of the cell 02096 const EntityHandle* conn; 02097 error = mb->get_connectivity( cur_cid, conn, nvpc, true );MB_CHK_ERR( error ); 02098 02099 // Local id of vid in the cell and the half-faces incident on it 02100 int lv = -1; 02101 for( int i = 0; i < nvpc; ++i ) 02102 { 02103 if( conn[i] == vid ) 02104 { 02105 lv = i; 02106 break; 02107 } 02108 }; 02109 if( lv < 0 ) MB_SET_ERR( MB_FAILURE, "did not find local vertex " ); 02110 // Number of local half-faces incident on the current vertex 02111 int nhf_thisv = lConnMap3D[index].v2hf_num[lv]; 02112 int cidx = ID_FROM_HANDLE( cur_cid ) - 1; 02113 02114 // Add new cells into the stack 02115 EntityHandle ngb; 02116 for( int i = 0; i < nhf_thisv; ++i ) 02117 { 02118 int ind = lConnMap3D[index].v2hf[lv][i]; 02119 hf = sibhfs[nfpc * cidx + ind]; 02120 ngb = fid_from_halfacet( hf, ctype ); 02121 02122 if( ngb ) 02123 { 02124 bool found_ent = find_match_in_array( ngb, cellq, qsize - 1 ); 02125 02126 if( !found_ent ) 02127 { 02128 cellq[qsize] = ngb; 02129 qsize += 1; 02130 } 02131 } 02132 } 02133 } 02134 02135 // Change the visited faces to false 02136 for( int i = 0; i < qsize; i++ ) 02137 cellq[i] = 0; 02138 02139 return MB_SUCCESS; 02140 } 02141 02142 ErrorCode HalfFacetRep::get_up_adjacencies_edg_3d( EntityHandle eid, std::vector< EntityHandle >& adjents, 02143 std::vector< int >* leids ) 02144 { 02145 ErrorCode error; 02146 EntityType ctype = mb->type_from_handle( *_cells.begin() ); 02147 EntityHandle start_cell = *_cells.begin(); 02148 int index = get_index_in_lmap( start_cell ); 02149 int nvpc = lConnMap3D[index].num_verts_in_cell; 02150 int nfpc = lConnMap3D[index].num_faces_in_cell; 02151 02152 adjents.reserve( 20 ); 02153 if( leids != NULL ) leids->reserve( 20 ); 02154 02155 // Find the edge vertices 02156 const EntityHandle* econn; 02157 int num_conn = 0; 02158 error = mb->get_connectivity( eid, econn, num_conn, true );MB_CHK_ERR( error ); 02159 02160 EntityHandle v_start = econn[0], v_end = econn[1]; 02161 int v1idx = ID_FROM_HANDLE( v_start ) - 1; 02162 int v2idx = ID_FROM_HANDLE( v_end ) - 1; 02163 02164 // Find an half-facets incident to each end vertex of the edge 02165 std::vector< EntityHandle > start_cells; 02166 HFacet hf1 = v2hf[v1idx]; 02167 HFacet hf2 = v2hf[v2idx]; 02168 02169 if( hf1 == 0 && !v2hfs.empty() ) 02170 { 02171 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02172 it_hes; 02173 it_hes = v2hfs.equal_range( v_start ); 02174 02175 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02176 { 02177 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02178 } 02179 } 02180 else if( hf1 != 0 ) 02181 { 02182 start_cells.push_back( fid_from_halfacet( hf1, ctype ) ); 02183 } 02184 02185 if( hf2 == 0 && !v2hfs.empty() ) 02186 { 02187 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02188 it_hes; 02189 it_hes = v2hfs.equal_range( v_end ); 02190 02191 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02192 { 02193 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02194 } 02195 } 02196 else if( hf2 != 0 ) 02197 { 02198 start_cells.push_back( fid_from_halfacet( hf2, ctype ) ); 02199 } 02200 02201 if( start_cells.empty() ) return MB_SUCCESS; 02202 02203 std::sort( start_cells.begin(), start_cells.end() ); 02204 std::vector< EntityHandle >::iterator last = std::unique( start_cells.begin(), start_cells.end() ); 02205 start_cells.erase( last, start_cells.end() ); 02206 02207 for( int i = 0; i < (int)start_cells.size(); i++ ) 02208 cellq[i] = start_cells[i]; 02209 02210 int qsize = start_cells.size(); 02211 int num_qvals = 0; 02212 02213 while( num_qvals < qsize ) 02214 { 02215 EntityHandle cell_id = cellq[num_qvals]; 02216 num_qvals += 1; 02217 02218 const EntityHandle* conn; 02219 error = mb->get_connectivity( cell_id, conn, nvpc, true );MB_CHK_ERR( error ); 02220 02221 int lv0 = -1, lv1 = -1, lv = -1; 02222 02223 // locate v_origin in poped out tet, check if v_end is in 02224 for( int i = 0; i < nvpc; i++ ) 02225 { 02226 if( v_start == conn[i] ) 02227 { 02228 lv0 = i; 02229 lv = lv0; 02230 } 02231 else if( v_end == conn[i] ) 02232 { 02233 lv1 = i; 02234 lv = lv1; 02235 } 02236 } 02237 02238 if( ( lv0 >= 0 ) && ( lv1 >= 0 ) ) 02239 { 02240 adjents.push_back( cell_id ); 02241 if( leids != NULL ) leids->push_back( lConnMap3D[index].lookup_leids[lv0][lv1] ); 02242 } 02243 02244 // push back new found unchecked incident tets of v_start 02245 int cidx = ID_FROM_HANDLE( cell_id ) - 1; 02246 if( lv < 0 ) MB_SET_ERR( MB_FAILURE, "did not find local vertex " ); 02247 int nhf_thisv = lConnMap3D[index].v2hf_num[lv]; 02248 02249 for( int i = 0; i < nhf_thisv; i++ ) 02250 { 02251 int ind = lConnMap3D[index].v2hf[lv][i]; 02252 HFacet hf = sibhfs[nfpc * cidx + ind]; 02253 EntityHandle ngb = fid_from_halfacet( hf, ctype ); 02254 02255 if( ngb ) 02256 { 02257 bool found_ent = find_match_in_array( ngb, &cellq[0], qsize - 1 ); 02258 02259 if( !found_ent ) 02260 { 02261 cellq[qsize] = ngb; 02262 qsize += 1; 02263 } 02264 } 02265 } 02266 } 02267 02268 for( int i = 0; i < qsize; i++ ) 02269 cellq[i] = 0; 02270 02271 return MB_SUCCESS; 02272 } 02273 02274 ErrorCode HalfFacetRep::get_up_adjacencies_edg_3d( EntityHandle cid, int leid, std::vector< EntityHandle >& adjents, 02275 std::vector< int >* leids, std::vector< int >* adj_orients ) 02276 { 02277 ErrorCode error; 02278 EntityType ctype = mb->type_from_handle( cid ); 02279 int index = get_index_in_lmap( cid ); 02280 int nvpc = lConnMap3D[index].num_verts_in_cell; 02281 int nfpc = lConnMap3D[index].num_faces_in_cell; 02282 adjents.clear(); 02283 adjents.reserve( 20 ); 02284 02285 if( leids != NULL ) 02286 { 02287 leids->clear(); 02288 leids->reserve( 20 ); 02289 } 02290 if( adj_orients != NULL ) 02291 { 02292 adj_orients->clear(); 02293 adj_orients->reserve( 20 ); 02294 } 02295 02296 const EntityHandle* econn; 02297 error = mb->get_connectivity( cid, econn, nvpc, true );MB_CHK_ERR( error ); 02298 02299 // Get the end vertices of the edge <cid,leid> 02300 int id = lConnMap3D[index].e2v[leid][0]; 02301 EntityHandle v_start = econn[id]; 02302 id = lConnMap3D[index].e2v[leid][1]; 02303 EntityHandle v_end = econn[id]; 02304 02305 int v1idx = ID_FROM_HANDLE( v_start ) - 1; 02306 int v2idx = ID_FROM_HANDLE( v_end ) - 1; 02307 02308 // Find an half-facets incident to each end vertex of the edge 02309 std::vector< EntityHandle > start_cells; 02310 HFacet hf1 = v2hf[v1idx]; 02311 HFacet hf2 = v2hf[v2idx]; 02312 02313 if( hf1 == 0 && !v2hfs.empty() ) 02314 { 02315 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02316 it_hes; 02317 it_hes = v2hfs.equal_range( v_start ); 02318 02319 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02320 { 02321 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02322 } 02323 } 02324 else if( hf1 != 0 ) 02325 { 02326 start_cells.push_back( fid_from_halfacet( hf1, ctype ) ); 02327 } 02328 02329 if( hf2 == 0 && !v2hfs.empty() ) 02330 { 02331 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02332 it_hes; 02333 it_hes = v2hfs.equal_range( v_end ); 02334 02335 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02336 { 02337 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02338 } 02339 } 02340 else if( hf2 != 0 ) 02341 { 02342 start_cells.push_back( fid_from_halfacet( hf2, ctype ) ); 02343 } 02344 02345 if( start_cells.empty() ) return MB_SUCCESS; 02346 02347 std::sort( start_cells.begin(), start_cells.end() ); 02348 std::vector< EntityHandle >::iterator last = std::unique( start_cells.begin(), start_cells.end() ); 02349 start_cells.erase( last, start_cells.end() ); 02350 02351 for( int i = 0; i < (int)start_cells.size(); i++ ) 02352 cellq[i] = start_cells[i]; 02353 02354 int qsize = start_cells.size(); 02355 int num_qvals = 0; 02356 02357 while( num_qvals < qsize ) 02358 { 02359 EntityHandle cell_id = cellq[num_qvals]; 02360 num_qvals += 1; 02361 02362 const EntityHandle* conn; 02363 error = mb->get_connectivity( cell_id, conn, nvpc, true );MB_CHK_ERR( error ); 02364 02365 int lv0 = -1, lv1 = -1, lv = -1; 02366 02367 // locate v_origin in poped out tet, check if v_end is in 02368 for( int i = 0; i < nvpc; i++ ) 02369 { 02370 if( v_start == conn[i] ) 02371 { 02372 lv0 = i; 02373 lv = lv0; 02374 } 02375 else if( v_end == conn[i] ) 02376 { 02377 lv1 = i; 02378 lv = lv1; 02379 } 02380 } 02381 02382 if( ( lv0 >= 0 ) && ( lv1 >= 0 ) ) 02383 { 02384 adjents.push_back( cell_id ); 02385 if( leids != NULL ) leids->push_back( lConnMap3D[index].lookup_leids[lv0][lv1] ); 02386 02387 if( adj_orients != NULL ) 02388 { 02389 int cur_leid = lConnMap3D[index].lookup_leids[lv0][lv1]; 02390 int id1 = lConnMap3D[index].e2v[cur_leid][0]; 02391 int id2 = lConnMap3D[index].e2v[cur_leid][1]; 02392 if( ( v_start == conn[id1] ) && ( v_end == conn[id2] ) ) 02393 adj_orients->push_back( 1 ); 02394 else if( ( v_start == conn[id2] ) && ( v_end == conn[id1] ) ) 02395 adj_orients->push_back( 0 ); 02396 } 02397 } 02398 if( lv < 0 ) MB_SET_ERR( MB_FAILURE, "did not find local vertex " ); 02399 // push back new found unchecked incident tets of v_start 02400 int cidx = ID_FROM_HANDLE( cell_id ) - 1; 02401 int nhf_thisv = lConnMap3D[index].v2hf_num[lv]; 02402 02403 for( int i = 0; i < nhf_thisv; i++ ) 02404 { 02405 int ind = lConnMap3D[index].v2hf[lv][i]; 02406 HFacet hf = sibhfs[nfpc * cidx + ind]; 02407 EntityHandle ngb = fid_from_halfacet( hf, ctype ); 02408 02409 if( ngb ) 02410 { 02411 bool found_ent = find_match_in_array( ngb, &cellq[0], qsize - 1 ); 02412 02413 if( !found_ent ) 02414 { 02415 cellq[qsize] = ngb; 02416 qsize += 1; 02417 } 02418 } 02419 } 02420 } 02421 02422 for( int i = 0; i < qsize; i++ ) 02423 cellq[i] = 0; 02424 02425 return MB_SUCCESS; 02426 } 02427 02428 ErrorCode HalfFacetRep::get_up_adjacencies_edg_3d_comp( EntityHandle cid, int leid, 02429 std::vector< EntityHandle >& adjents, std::vector< int >* leids, 02430 std::vector< int >* adj_orients ) 02431 { 02432 ErrorCode error; 02433 EntityType ctype = mb->type_from_handle( cid ); 02434 int index = get_index_in_lmap( cid ); 02435 int nvpc = lConnMap3D[index].num_verts_in_cell; 02436 int nfpc = lConnMap3D[index].num_faces_in_cell; 02437 adjents.clear(); 02438 adjents.reserve( 20 ); 02439 02440 if( leids != NULL ) 02441 { 02442 leids->clear(); 02443 leids->reserve( 20 ); 02444 } 02445 if( adj_orients != NULL ) 02446 { 02447 adj_orients->clear(); 02448 adj_orients->reserve( 20 ); 02449 } 02450 02451 const EntityHandle* econn; 02452 error = mb->get_connectivity( cid, econn, nvpc );MB_CHK_ERR( error ); 02453 02454 // Get the end vertices of the edge <cid,leid> 02455 int id = lConnMap3D[index].e2v[leid][0]; 02456 EntityHandle v_start = econn[id]; 02457 id = lConnMap3D[index].e2v[leid][1]; 02458 EntityHandle v_end = econn[id]; 02459 02460 int v1idx = ID_FROM_HANDLE( v_start ) - 1; 02461 int v2idx = ID_FROM_HANDLE( v_end ) - 1; 02462 02463 // Find an half-facets incident to each end vertex of the edge 02464 std::vector< EntityHandle > start_cells; 02465 HFacet hf1 = v2hf[v1idx]; 02466 HFacet hf2 = v2hf[v2idx]; 02467 02468 if( ( hf1 == 0 ) && ( v2hfs.find( v_start ) != v2hfs.end() ) && ( hf2 == 0 ) && 02469 ( v2hfs.find( v_end ) != v2hfs.end() ) ) 02470 { 02471 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02472 it_hes; 02473 it_hes = v2hfs.equal_range( v_start ); 02474 02475 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02476 { 02477 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02478 } 02479 } 02480 else 02481 return MB_SUCCESS; 02482 02483 if( start_cells.empty() ) return MB_SUCCESS; 02484 02485 // std::sort(start_cells.begin(), start_cells.end()); 02486 // std::vector<EntityHandle>::iterator last = std::unique(start_cells.begin(), 02487 // start_cells.end()); start_cells.erase(last, start_cells.end()); 02488 02489 for( int c = 0; c < (int)start_cells.size(); c++ ) 02490 { 02491 cellq[0] = start_cells[c]; 02492 02493 int qsize = 1; 02494 int num_qvals = 0; 02495 02496 while( num_qvals < qsize ) 02497 { 02498 EntityHandle cell_id = cellq[num_qvals]; 02499 num_qvals += 1; 02500 02501 const EntityHandle* conn; 02502 error = mb->get_connectivity( cell_id, conn, nvpc );MB_CHK_ERR( error ); 02503 02504 int lv0 = -1, lv1 = -1, lv = -1; 02505 02506 // locate v_origin in poped out tet, check if v_end is in 02507 for( int i = 0; i < nvpc; i++ ) 02508 { 02509 if( v_start == conn[i] ) 02510 { 02511 lv0 = i; 02512 lv = lv0; 02513 } 02514 else if( v_end == conn[i] ) 02515 { 02516 lv1 = i; 02517 lv = lv1; 02518 } 02519 } 02520 02521 if( ( lv0 >= 0 ) && ( lv1 >= 0 ) ) 02522 { 02523 adjents.push_back( cell_id ); 02524 if( leids != NULL ) leids->push_back( lConnMap3D[index].lookup_leids[lv0][lv1] ); 02525 02526 if( adj_orients != NULL ) 02527 { 02528 int cur_leid = lConnMap3D[index].lookup_leids[lv0][lv1]; 02529 int id1 = lConnMap3D[index].e2v[cur_leid][0]; 02530 int id2 = lConnMap3D[index].e2v[cur_leid][1]; 02531 if( ( v_start == conn[id1] ) && ( v_end == conn[id2] ) ) 02532 adj_orients->push_back( 1 ); 02533 else if( ( v_start == conn[id2] ) && ( v_end == conn[id1] ) ) 02534 adj_orients->push_back( 0 ); 02535 } 02536 02537 for( int i = 0; i < qsize; i++ ) 02538 cellq[i] = 0; 02539 02540 break; 02541 } 02542 02543 // push back new found unchecked incident tets of v_start 02544 int cidx = ID_FROM_HANDLE( cell_id ) - 1; 02545 int nhf_thisv = lConnMap3D[index].v2hf_num[lv]; 02546 02547 for( int i = 0; i < nhf_thisv; i++ ) 02548 { 02549 int ind = lConnMap3D[index].v2hf[lv][i]; 02550 HFacet hf = sibhfs[nfpc * cidx + ind]; 02551 EntityHandle ngb = fid_from_halfacet( hf, ctype ); 02552 02553 if( ngb ) 02554 { 02555 bool found_ent = find_match_in_array( ngb, &cellq[0], qsize - 1 ); 02556 02557 if( !found_ent ) 02558 { 02559 cellq[qsize] = ngb; 02560 qsize += 1; 02561 } 02562 } 02563 } 02564 } 02565 02566 for( int i = 0; i < qsize; i++ ) 02567 cellq[i] = 0; 02568 } 02569 02570 return MB_SUCCESS; 02571 } 02572 02573 ErrorCode HalfFacetRep::get_up_adjacencies_face_3d( EntityHandle fid, std::vector< EntityHandle >& adjents, 02574 std::vector< int >* lfids ) 02575 { 02576 ErrorCode error; 02577 02578 EntityHandle cid = 0; 02579 int lid = 0; 02580 bool found = find_matching_halfface( fid, &cid, &lid ); 02581 02582 if( found ) 02583 { 02584 error = get_up_adjacencies_face_3d( cid, lid, adjents, lfids );MB_CHK_ERR( error ); 02585 } 02586 02587 return MB_SUCCESS; 02588 } 02589 02590 ErrorCode HalfFacetRep::get_up_adjacencies_face_3d( EntityHandle cid, int lfid, std::vector< EntityHandle >& adjents, 02591 std::vector< int >* lfids ) 02592 { 02593 02594 EntityHandle start_cell = *_cells.begin(); 02595 EntityType ctype = mb->type_from_handle( start_cell ); 02596 int index = get_index_in_lmap( start_cell ); 02597 int nfpc = lConnMap3D[index].num_faces_in_cell; 02598 02599 adjents.reserve( 4 ); 02600 adjents.push_back( cid ); 02601 02602 if( lfids != NULL ) 02603 { 02604 lfids->reserve( 4 ); 02605 lfids->push_back( lfid ); 02606 } 02607 02608 int cidx = ID_FROM_HANDLE( cid ) - 1; 02609 HFacet hf = sibhfs[nfpc * cidx + lfid]; 02610 EntityHandle sibcid = fid_from_halfacet( hf, ctype ); 02611 int siblid = lid_from_halffacet( hf ); 02612 02613 if( sibcid != 0 ) 02614 { 02615 adjents.push_back( sibcid ); 02616 if( lfids != NULL ) lfids->push_back( siblid ); 02617 } 02618 02619 return MB_SUCCESS; 02620 } 02621 02622 bool HalfFacetRep::find_matching_implicit_edge_in_cell( EntityHandle eid, std::vector< EntityHandle >& cid, 02623 std::vector< int >& leid ) 02624 { 02625 ErrorCode error; 02626 EntityType ctype = mb->type_from_handle( *_cells.begin() ); 02627 EntityHandle start_cell = *_cells.begin(); 02628 int index = get_index_in_lmap( start_cell ); 02629 int nvpc = lConnMap3D[index].num_verts_in_cell; 02630 int nfpc = lConnMap3D[index].num_faces_in_cell; 02631 02632 // Find the edge vertices 02633 const EntityHandle* econn; 02634 int num_conn = 0; 02635 error = mb->get_connectivity( eid, econn, num_conn, true );MB_CHK_ERR( error ); 02636 02637 EntityHandle v_start = econn[0], v_end = econn[1]; 02638 int v1idx = ID_FROM_HANDLE( v_start ) - 1; 02639 int v2idx = ID_FROM_HANDLE( v_end ) - 1; 02640 02641 // Find an incident cell to v_start 02642 std::vector< EntityHandle > start_cells; 02643 HFacet hf1 = v2hf[v1idx]; 02644 HFacet hf2 = v2hf[v2idx]; 02645 02646 int ncomps1 = 0, ncomps2 = 0, ncomp; 02647 if( hf1 == 0 && !v2hfs.empty() ) 02648 { 02649 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02650 it_hes; 02651 it_hes = v2hfs.equal_range( v_start ); 02652 02653 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02654 { 02655 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02656 ncomps1 += 1; 02657 } 02658 } 02659 else if( hf1 != 0 ) 02660 { 02661 start_cells.push_back( fid_from_halfacet( hf1, ctype ) ); 02662 ncomps1 += 1; 02663 } 02664 02665 if( hf2 == 0 && !v2hfs.empty() ) 02666 { 02667 std::pair< std::multimap< EntityHandle, HFacet >::iterator, std::multimap< EntityHandle, HFacet >::iterator > 02668 it_hes; 02669 it_hes = v2hfs.equal_range( v_end ); 02670 02671 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 02672 { 02673 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02674 ncomps2 += 1; 02675 } 02676 } 02677 else if( hf2 != 0 ) 02678 { 02679 start_cells.push_back( fid_from_halfacet( hf2, ctype ) ); 02680 ncomps2 += 1; 02681 } 02682 02683 ncomp = std::min( ncomps1, ncomps2 ); 02684 02685 bool found = false; 02686 if( start_cells.empty() ) return found; 02687 02688 for( int i = 0; i < (int)start_cells.size(); i++ ) 02689 cellq[i] = start_cells[i]; 02690 02691 int qsize = start_cells.size(); 02692 int num_qvals = 0; 02693 02694 while( num_qvals < qsize ) 02695 { 02696 EntityHandle cell_id = cellq[num_qvals]; 02697 num_qvals += 1; 02698 02699 const EntityHandle* conn; 02700 error = mb->get_connectivity( cell_id, conn, nvpc, true );MB_CHK_ERR( error ); 02701 02702 int lv0 = -1, lv1 = -1, lv = -1; 02703 02704 // locate v_origin in poped out tet, check if v_end is in 02705 for( int i = 0; i < nvpc; i++ ) 02706 { 02707 if( v_start == conn[i] ) 02708 { 02709 lv0 = i; 02710 lv = lv0; 02711 } 02712 else if( v_end == conn[i] ) 02713 { 02714 lv1 = i; 02715 lv = lv1; 02716 } 02717 } 02718 02719 if( ( lv0 >= 0 ) && ( lv1 >= 0 ) ) 02720 { 02721 found = true; 02722 cid.push_back( cell_id ); 02723 leid.push_back( lConnMap3D[index].lookup_leids[lv0][lv1] ); 02724 02725 if( (int)cid.size() == ncomp ) break; 02726 } 02727 02728 // push back new found unchecked incident tets of v_start 02729 int cidx = ID_FROM_HANDLE( cell_id ) - 1; 02730 int nhf_thisv = lConnMap3D[index].v2hf_num[lv]; 02731 if( lv < 0 ) MB_SET_ERR( MB_FAILURE, "did not find local vertex " ); 02732 for( int i = 0; i < nhf_thisv; i++ ) 02733 { 02734 int ind = lConnMap3D[index].v2hf[lv][i]; 02735 HFacet hf = sibhfs[nfpc * cidx + ind]; 02736 EntityHandle ngb = fid_from_halfacet( hf, ctype ); 02737 02738 if( ngb ) 02739 { 02740 bool found_ent = find_match_in_array( ngb, &cellq[0], qsize - 1 ); 02741 02742 if( !found_ent ) 02743 { 02744 cellq[qsize] = ngb; 02745 qsize += 1; 02746 } 02747 } 02748 } 02749 } 02750 02751 for( int i = 0; i < qsize; i++ ) 02752 cellq[i] = 0; 02753 02754 return found; 02755 } 02756 02757 bool HalfFacetRep::find_matching_halfface( EntityHandle fid, EntityHandle* cid, int* lid ) 02758 { 02759 ErrorCode error; 02760 EntityHandle start_cell = *_cells.begin(); 02761 EntityType ctype = mb->type_from_handle( start_cell ); 02762 int index = get_index_in_lmap( start_cell ); 02763 int nvpc = lConnMap3D[index].num_verts_in_cell; 02764 int nfpc = lConnMap3D[index].num_faces_in_cell; 02765 EntityType ftype = mb->type_from_handle( fid ); 02766 int nvF = lConnMap2D[ftype - 2].num_verts_in_face; 02767 02768 const EntityHandle* fid_verts; 02769 error = mb->get_connectivity( fid, fid_verts, nvF, true );MB_CHK_ERR( error ); 02770 02771 std::vector< EntityHandle > start_cells; 02772 int vidx, locfv0 = -1; 02773 HFacet hf = 0; 02774 02775 for( int i = 0; i < nvF; i++ ) 02776 { 02777 vidx = ID_FROM_HANDLE( fid_verts[i] ) - 1; 02778 hf = v2hf[vidx]; 02779 if( hf != 0 ) 02780 { 02781 start_cells.push_back( fid_from_halfacet( hf, ctype ) ); 02782 locfv0 = i; 02783 break; 02784 } 02785 else if( hf == 0 && v2hfs.find( fid_verts[i] ) != v2hfs.end() ) 02786 { 02787 std::pair< std::multimap< EntityHandle, HFacet >::iterator, 02788 std::multimap< EntityHandle, HFacet >::iterator > 02789 it_hfs; 02790 it_hfs = v2hfs.equal_range( fid_verts[i] ); 02791 02792 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hfs.first; it != it_hfs.second; ++it ) 02793 { 02794 start_cells.push_back( fid_from_halfacet( it->second, ctype ) ); 02795 } 02796 locfv0 = i; 02797 break; 02798 } 02799 } 02800 02801 if( start_cells.empty() ) return false; 02802 02803 EntityHandle cur_cid; 02804 bool found = false; 02805 02806 int Stksize = 0, count = -1; 02807 02808 for( int i = 0; i < (int)start_cells.size(); i++ ) 02809 Stkcells[i] = start_cells[i]; 02810 02811 Stksize = start_cells.size() - 1; 02812 02813 while( Stksize >= 0 ) 02814 { 02815 cur_cid = Stkcells[Stksize]; 02816 Stksize -= 1; 02817 count += 1; 02818 trackcells[count] = cur_cid; 02819 02820 const EntityHandle* conn; 02821 error = mb->get_connectivity( cur_cid, conn, nvpc, true );MB_CHK_ERR( error ); 02822 02823 int lv[4] = { -1, -1, -1, -1 }; 02824 int cnt = 0; 02825 for( int i = 0; i < nvpc; i++ ) 02826 { 02827 for( int j = 0; j < nvF; j++ ) 02828 if( conn[i] == fid_verts[j] ) 02829 { 02830 lv[j] = i; 02831 cnt += 1; 02832 } 02833 } 02834 if( cnt == nvF ) // All face verts are part of the cell 02835 { 02836 found = true; 02837 int nhf_thisv = lConnMap3D[index].v2hf_num[lv[locfv0]]; 02838 int lfid = -1; 02839 for( int i = 0; i < nhf_thisv; ++i ) 02840 { 02841 lfid = lConnMap3D[index].v2hf[lv[locfv0]][i]; 02842 int lcnt = 0; 02843 for( int j = 0; j < nvF; ++j ) 02844 { 02845 for( int k = 0; k < nvF; ++k ) 02846 { 02847 if( lv[k] == lConnMap3D[index].hf2v[lfid][j] ) lcnt += 1; 02848 } 02849 } 02850 if( lcnt == nvF ) break; 02851 } 02852 cid[0] = cur_cid; 02853 lid[0] = lfid; 02854 02855 break; 02856 } 02857 else 02858 { 02859 // Add other cells that are incident on fid_verts[0] 02860 if( locfv0 < 0 || lv[locfv0] < 0 ) MB_SET_ERR( MB_FAILURE, "did not find local vertex " ); 02861 int nhf_thisv = lConnMap3D[index].v2hf_num[lv[locfv0]]; 02862 int cidx = ID_FROM_HANDLE( cur_cid ) - 1; 02863 02864 // Add new cells into the stack 02865 EntityHandle ngb; 02866 for( int i = 0; i < nhf_thisv; ++i ) 02867 { 02868 int ind = lConnMap3D[index].v2hf[lv[locfv0]][i]; 02869 hf = sibhfs[nfpc * cidx + ind]; 02870 ngb = fid_from_halfacet( hf, ctype ); 02871 02872 if( ngb ) 02873 { 02874 02875 bool found_ent = find_match_in_array( ngb, trackcells, count ); 02876 02877 if( !found_ent ) 02878 { 02879 Stksize += 1; 02880 Stkcells[Stksize] = ngb; 02881 } 02882 } 02883 } 02884 } 02885 } 02886 02887 // Change the visited faces to false 02888 for( int i = 0; i < Stksize; i++ ) 02889 Stkcells[i] = 0; 02890 02891 for( int i = 0; i <= count; i++ ) 02892 trackcells[i] = 0; 02893 02894 return found; 02895 } 02896 02897 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 02898 ErrorCode HalfFacetRep::get_neighbor_adjacencies_3d( EntityHandle cid, std::vector< EntityHandle >& adjents ) 02899 { 02900 adjents.reserve( 20 ); 02901 EntityType ctype = mb->type_from_handle( cid ); 02902 int index = get_index_in_lmap( cid ); 02903 int nfpc = lConnMap3D[index].num_faces_in_cell; 02904 int cidx = ID_FROM_HANDLE( cid ) - 1; 02905 02906 if( cid != 0 ) 02907 { 02908 for( int lfid = 0; lfid < nfpc; ++lfid ) 02909 { 02910 HFacet hf = sibhfs[nfpc * cidx + lfid]; 02911 EntityHandle sibcid = fid_from_halfacet( hf, ctype ); 02912 if( sibcid != 0 ) adjents.push_back( sibcid ); 02913 } 02914 } 02915 02916 return MB_SUCCESS; 02917 } 02918 02919 ///////////////////////////////////////////////////////////////////////////////////////////////// 02920 ErrorCode HalfFacetRep::get_down_adjacencies_edg_3d( EntityHandle cid, std::vector< EntityHandle >& adjents ) 02921 { 02922 // TODO: Try intersection without using std templates 02923 // Returns explicit edges, if any, of the face 02924 ErrorCode error; 02925 adjents.reserve( 20 ); 02926 int index = get_index_in_lmap( cid ); 02927 int nvpc = lConnMap3D[index].num_verts_in_cell; 02928 02929 const EntityHandle* conn; 02930 error = mb->get_connectivity( cid, conn, nvpc, true );MB_CHK_ERR( error ); 02931 02932 // Gather all the incident edges on each vertex of the face 02933 int ns = lConnMap3D[index].search_everts[0]; 02934 std::vector< EntityHandle > temp; 02935 for( int i = 0; i < ns; i++ ) 02936 { 02937 temp.clear(); 02938 int lv0 = lConnMap3D[index].search_everts[i + 1]; 02939 error = get_up_adjacencies_1d( conn[lv0], temp );MB_CHK_ERR( error ); 02940 02941 int nle = lConnMap3D[index].v2le[i][0]; 02942 int count = 0; 02943 for( int j = 0; j < (int)temp.size(); j++ ) 02944 { 02945 const EntityHandle* econn; 02946 int nvpe = 0; 02947 error = mb->get_connectivity( temp[j], econn, nvpe, true );MB_CHK_ERR( error ); 02948 02949 for( int k = 0; k < nle; k++ ) 02950 { 02951 int lv1 = lConnMap3D[index].v2le[i][k + 1]; 02952 if( ( ( econn[0] == conn[lv0] ) && ( econn[1] == conn[lv1] ) ) || 02953 ( ( econn[1] == conn[lv0] ) && ( econn[0] == conn[lv1] ) ) ) 02954 { 02955 adjents.push_back( temp[j] ); 02956 count += 1; 02957 } 02958 } 02959 if( count == nle ) break; 02960 } 02961 } 02962 return MB_SUCCESS; 02963 } 02964 02965 ErrorCode HalfFacetRep::get_down_adjacencies_face_3d( EntityHandle cid, std::vector< EntityHandle >& adjents ) 02966 { 02967 // Returns explicit face, if any of the cell 02968 ErrorCode error; 02969 adjents.reserve( 10 ); 02970 int index = get_index_in_lmap( cid ); 02971 int nvpc = lConnMap3D[index].num_verts_in_cell; 02972 int nfpc = lConnMap3D[index].num_faces_in_cell; 02973 02974 // Get the connectivity of the input cell 02975 const EntityHandle* conn; 02976 error = mb->get_connectivity( cid, conn, nvpc, true );MB_CHK_ERR( error ); 02977 02978 // Collect all the half-faces of the cell 02979 EntityHandle half_faces[6][4]; 02980 for( int i = 0; i < nfpc; i++ ) 02981 { 02982 int nvf = lConnMap3D[index].hf2v_num[i]; 02983 for( int j = 0; j < nvf; j++ ) 02984 { 02985 int ind = lConnMap3D[index].hf2v[i][j]; 02986 half_faces[i][j] = conn[ind]; 02987 } 02988 } 02989 02990 // Add two vertices for which the upward adjacencies are computed 02991 int search_verts[2]; 02992 search_verts[0] = lConnMap3D[index].search_fverts[0]; 02993 search_verts[1] = lConnMap3D[index].search_fverts[1]; 02994 02995 std::vector< EntityHandle > temp; 02996 temp.reserve( 20 ); 02997 for( int i = 0; i < 2; i++ ) 02998 { 02999 // Get the incident faces on the local vertex 03000 int lv = search_verts[i]; 03001 temp.clear(); 03002 error = get_up_adjacencies_vert_2d( conn[lv], temp );MB_CHK_ERR( error ); 03003 03004 if( temp.size() == 0 ) continue; 03005 03006 // Get the half-faces incident on the local vertex and match it with the obtained faces 03007 int nhfthisv = lConnMap3D[index].v2hf_num[lv]; 03008 for( int k = 0; k < (int)temp.size(); k++ ) 03009 { 03010 const EntityHandle* fid_verts; 03011 int fsize = 0; 03012 error = mb->get_connectivity( temp[k], fid_verts, fsize, true );MB_CHK_ERR( error ); 03013 03014 for( int j = 0; j < nhfthisv; j++ ) 03015 { 03016 // Collect all the vertices of this half-face 03017 int idx = lConnMap3D[index].v2hf[lv][j]; 03018 int nvF = lConnMap3D[index].hf2v_num[idx]; 03019 03020 if( fsize != nvF ) continue; 03021 03022 int direct, offset; 03023 bool they_match = CN::ConnectivityMatch( &half_faces[idx][0], &fid_verts[0], nvF, direct, offset ); 03024 03025 if( they_match ) 03026 { 03027 bool found = false; 03028 for( int p = 0; p < (int)adjents.size(); p++ ) 03029 { 03030 if( adjents[p] == temp[k] ) 03031 { 03032 found = true; 03033 break; 03034 } 03035 } 03036 if( !found ) adjents.push_back( temp[k] ); 03037 } 03038 } 03039 } 03040 } 03041 03042 return MB_SUCCESS; 03043 } 03044 //////////////////////////////////////////////////////////////////////////////// 03045 ErrorCode HalfFacetRep::find_total_edges_faces_3d( Range cells, int* nedges, int* nfaces ) 03046 { 03047 ErrorCode error; 03048 int index = get_index_in_lmap( *cells.begin() ); 03049 int nepc = lConnMap3D[index].num_edges_in_cell; 03050 int nfpc = lConnMap3D[index].num_faces_in_cell; 03051 int ncells = cells.size(); 03052 int total_edges = nepc * ncells; 03053 int total_faces = nfpc * ncells; 03054 03055 std::vector< int > trackE( total_edges, 0 ); 03056 std::vector< int > trackF( total_faces, 0 ); 03057 03058 std::vector< EntityHandle > inc_cids, sib_cids; 03059 std::vector< int > inc_leids, sib_lfids; 03060 03061 for( Range::iterator it = cells.begin(); it != cells.end(); ++it ) 03062 { 03063 // Count edges 03064 for( int i = 0; i < nepc; i++ ) 03065 { 03066 inc_cids.clear(); 03067 inc_leids.clear(); 03068 03069 int id = nepc * ( cells.index( *it ) ) + i; 03070 if( !trackE[id] ) 03071 { 03072 error = get_up_adjacencies_edg_3d( *it, i, inc_cids, &inc_leids );MB_CHK_ERR( error ); 03073 03074 total_edges -= inc_cids.size() - 1; 03075 for( int j = 0; j < (int)inc_cids.size(); j++ ) 03076 trackE[nepc * ( cells.index( inc_cids[j] ) ) + inc_leids[j]] = 1; 03077 } 03078 } 03079 03080 // Count faces 03081 for( int i = 0; i < nfpc; i++ ) 03082 { 03083 sib_cids.clear(); 03084 sib_lfids.clear(); 03085 03086 int id = nfpc * ( cells.index( *it ) ) + i; 03087 if( !trackF[id] ) 03088 { 03089 error = get_up_adjacencies_face_3d( *it, i, sib_cids, &sib_lfids );MB_CHK_ERR( error ); 03090 03091 if( sib_cids.size() == 1 ) continue; 03092 03093 total_faces -= sib_cids.size() - 1; 03094 trackF[nfpc * ( cells.index( sib_cids[1] ) ) + sib_lfids[1]] = 1; 03095 } 03096 } 03097 } 03098 03099 nedges[0] = total_edges; 03100 nfaces[0] = total_faces; 03101 03102 return MB_SUCCESS; 03103 } 03104 03105 bool HalfFacetRep::find_match_in_array( EntityHandle ent, EntityHandle* ent_list, int count, bool get_index, 03106 int* index ) 03107 { 03108 bool found = false; 03109 for( int i = 0; i <= count; i++ ) 03110 { 03111 if( ent == ent_list[i] ) 03112 { 03113 found = true; 03114 if( get_index ) *index = i; 03115 break; 03116 } 03117 } 03118 03119 return found; 03120 } 03121 03122 ErrorCode HalfFacetRep::get_half_facet_in_comp( EntityHandle cid, int leid, std::vector< EntityHandle >& ents, 03123 std::vector< int >& lids, std::vector< int >& lfids ) 03124 { 03125 ErrorCode error; 03126 ents.clear(); 03127 lids.clear(); 03128 EntityType ctype = mb->type_from_handle( cid ); 03129 int index = get_index_in_lmap( *_cells.begin() ); 03130 int nfpc = lConnMap3D[index].num_faces_in_cell; 03131 int nvpc = lConnMap3D[index].num_verts_in_cell; 03132 03133 // Get all incident cells 03134 std::vector< EntityHandle > adjents; 03135 std::vector< int > adjlids; 03136 error = get_up_adjacencies_edg_3d( cid, leid, adjents, &adjlids );MB_CHK_ERR( error ); 03137 03138 // Get the end vertices of the edge <cid,leid> 03139 const EntityHandle* econn; 03140 error = mb->get_connectivity( cid, econn, nvpc, true );MB_CHK_ERR( error ); 03141 int id = lConnMap3D[index].e2v[leid][0]; 03142 EntityHandle vstart = econn[id]; 03143 id = lConnMap3D[index].e2v[leid][1]; 03144 EntityHandle vend = econn[id]; 03145 03146 std::vector< int > mark( adjents.size(), 0 ); 03147 int count = 0; 03148 03149 for( int k = 0; k < (int)adjents.size(); k++ ) 03150 { 03151 03152 if( mark[k] != 0 ) continue; 03153 03154 count += 1; 03155 mark[k] = count; 03156 03157 // Loop over each half-face incident on this local edge 03158 for( int i = 0; i < 2; i++ ) 03159 { 03160 EntityHandle cur_cell = adjents[k]; 03161 int cur_leid = adjlids[k]; 03162 03163 int lface = i; 03164 03165 while( true ) 03166 { 03167 int cidx = ID_FROM_HANDLE( cur_cell ) - 1; 03168 int lfid = lConnMap3D[index].e2hf[cur_leid][lface]; 03169 03170 HFacet hf = sibhfs[nfpc * cidx + lfid]; 03171 cur_cell = fid_from_halfacet( hf, ctype ); 03172 lfid = lid_from_halffacet( hf ); 03173 03174 // Check if loop reached starting cell or a boundary 03175 if( ( cur_cell == adjents[k] ) || ( cur_cell == 0 ) ) break; 03176 03177 const EntityHandle* sib_conn; 03178 error = mb->get_connectivity( cur_cell, sib_conn, nvpc, true );MB_CHK_ERR( error ); 03179 03180 // Find the local edge id wrt to sibhf 03181 int nv_curF = lConnMap3D[index].hf2v_num[lfid]; 03182 int lv0 = -1, lv1 = -1, idx = -1; 03183 for( int j = 0; j < nv_curF; j++ ) 03184 { 03185 idx = lConnMap3D[index].hf2v[lfid][j]; 03186 if( vstart == sib_conn[idx] ) lv0 = idx; 03187 if( vend == sib_conn[idx] ) lv1 = idx; 03188 } 03189 03190 assert( ( lv0 >= 0 ) && ( lv1 >= 0 ) ); 03191 cur_leid = lConnMap3D[index].lookup_leids[lv0][lv1]; 03192 03193 int chk_lfid = lConnMap3D[index].e2hf[cur_leid][0]; 03194 03195 if( lfid == chk_lfid ) 03196 lface = 1; 03197 else 03198 lface = 0; 03199 03200 int ind = std::find( adjents.begin(), adjents.end(), cur_cell ) - adjents.begin(); 03201 mark[ind] = count; 03202 } 03203 03204 // Loop back 03205 if( cur_cell != 0 ) break; 03206 } 03207 } 03208 03209 // Loop over again to find cells on the boundary 03210 for( int c = 0; c < count; c++ ) 03211 { 03212 for( int i = 0; i < (int)adjents.size(); i++ ) 03213 { 03214 if( mark[i] == c + 1 ) 03215 { 03216 int cidx = ID_FROM_HANDLE( adjents[i] ) - 1; 03217 for( int j = 0; j < nfpc; j++ ) 03218 { 03219 HFacet hf = sibhfs[nfpc * cidx + j]; 03220 EntityHandle cell = fid_from_halfacet( hf, ctype ); 03221 if( cell == 0 ) 03222 { 03223 ents.push_back( adjents[i] ); 03224 lids.push_back( adjlids[i] ); 03225 lfids.push_back( j ); 03226 break; 03227 } 03228 } 03229 } 03230 } 03231 } 03232 03233 return MB_SUCCESS; 03234 } 03235 03236 /////////////////////////////////////////////////////////////////////////////////////////// 03237 ErrorCode HalfFacetRep::resize_hf_maps( EntityHandle start_vert, int nverts, EntityHandle start_edge, int nedges, 03238 EntityHandle start_face, int nfaces, EntityHandle start_cell, int ncells ) 03239 { 03240 int nwsz = 0, insz = 0; 03241 if( nedges ) 03242 { 03243 if( ID_FROM_HANDLE( ( *( _edges.end() - 1 ) + 1 ) ) != ID_FROM_HANDLE( start_edge ) ) 03244 nwsz = ( ID_FROM_HANDLE( start_edge ) - ID_FROM_HANDLE( *_edges.end() ) + nedges ) * 2; 03245 else 03246 nwsz = nedges * 2; 03247 insz = sibhvs.size(); 03248 sibhvs.resize( insz + nwsz, 0 ); 03249 03250 if( v2hv.empty() ) 03251 { 03252 if( ( !v2he.empty() ) ) 03253 insz = v2he.size(); 03254 else if( ( !v2hf.empty() ) ) 03255 insz = v2hf.size(); 03256 else 03257 MB_SET_ERR( MB_FAILURE, "Trying to resize ahf maps for a mesh with no edges, faces and cells" ); 03258 } 03259 else 03260 insz = v2hv.size(); 03261 03262 if( ID_FROM_HANDLE( *( _verts.end() - 1 ) + 1 ) != ID_FROM_HANDLE( start_vert ) ) 03263 nwsz = ID_FROM_HANDLE( start_vert ) - ID_FROM_HANDLE( *_verts.end() ) + nverts; 03264 else 03265 nwsz = nverts; 03266 v2hv.resize( insz + nwsz, 0 ); 03267 } 03268 03269 if( nfaces ) 03270 { 03271 EntityType ftype = mb->type_from_handle( *_faces.begin() ); 03272 int nepf = lConnMap2D[ftype - 2].num_verts_in_face; 03273 03274 if( ID_FROM_HANDLE( ( *( _faces.end() - 1 ) + 1 ) ) != ID_FROM_HANDLE( start_face ) ) 03275 nwsz = ( ID_FROM_HANDLE( start_face ) - ID_FROM_HANDLE( *_faces.end() ) + nfaces ) * nepf; 03276 else 03277 nwsz = nfaces * nepf; 03278 insz = sibhes.size(); 03279 sibhes.resize( insz + nwsz, 0 ); 03280 03281 if( ID_FROM_HANDLE( *( _verts.end() - 1 ) + 1 ) != ID_FROM_HANDLE( start_vert ) ) 03282 nwsz = ID_FROM_HANDLE( start_vert ) - ID_FROM_HANDLE( *_verts.end() ) + nverts; 03283 else 03284 nwsz = nverts; 03285 insz = v2he.size(); 03286 v2he.resize( insz + nwsz, 0 ); 03287 } 03288 03289 if( ncells ) 03290 { 03291 int index = get_index_in_lmap( *_cells.begin() ); 03292 int nfpc = lConnMap3D[index].num_faces_in_cell; 03293 03294 if( ID_FROM_HANDLE( ( *( _cells.end() - 1 ) + 1 ) ) != ID_FROM_HANDLE( start_cell ) ) 03295 nwsz = ( ID_FROM_HANDLE( start_cell ) - ID_FROM_HANDLE( *_cells.end() ) + ncells ) * nfpc; 03296 else 03297 nwsz = ncells * nfpc; 03298 insz = sibhfs.size(); 03299 sibhfs.resize( insz + nwsz, 0 ); 03300 03301 if( ID_FROM_HANDLE( *( _verts.end() - 1 ) + 1 ) != ID_FROM_HANDLE( start_vert ) ) 03302 nwsz = ID_FROM_HANDLE( start_vert ) - ID_FROM_HANDLE( *_verts.end() ) + nverts; 03303 else 03304 nwsz = nverts; 03305 insz = v2hf.size(); 03306 v2hf.resize( insz + nwsz, 0 ); 03307 } 03308 03309 return MB_SUCCESS; 03310 } 03311 03312 bool HalfFacetRep::check_nonmanifold_vertices( EntityType type, EntityHandle vid ) 03313 { 03314 bool status = false; 03315 if( type == MBTRI || type == MBQUAD ) 03316 { 03317 HFacet hf = v2he[ID_FROM_HANDLE( vid ) - 1]; 03318 if( hf == 0 && ( v2hes.find( vid ) != v2hes.end() ) ) status = true; 03319 } 03320 else if( type == MBTET || type == MBHEX ) 03321 { 03322 HFacet hf = v2hf[ID_FROM_HANDLE( vid ) - 1]; 03323 if( hf == 0 && ( v2hfs.find( vid ) != v2hfs.end() ) ) status = true; 03324 } 03325 else 03326 MB_SET_ERR( MB_FAILURE, "Requesting non-manifold vertex checks for either (1) 1D mesh or " 03327 "(2) not-implemented entity types" ); 03328 03329 return status; 03330 } 03331 03332 ErrorCode HalfFacetRep::get_sibling_map( EntityType type, EntityHandle ent, EntityHandle* sib_entids, int* sib_lids, 03333 int num_halffacets ) 03334 { 03335 03336 if( type == MBEDGE ) 03337 { 03338 if( num_halffacets != 2 ) MB_SET_ERR( MB_FAILURE, "Incorrect number of halfvertices." ); 03339 03340 int eidx = ID_FROM_HANDLE( ent ) - 1; 03341 for( int i = 0; i < 2; i++ ) 03342 { 03343 HFacet hf = sibhvs[2 * eidx + i]; 03344 sib_entids[i] = fid_from_halfacet( hf, MBEDGE ); 03345 sib_lids[i] = lid_from_halffacet( hf ); 03346 } 03347 } 03348 else if( type == MBTRI || type == MBQUAD ) 03349 { 03350 int nepf = lConnMap2D[type - 2].num_verts_in_face; 03351 03352 if( num_halffacets != nepf ) MB_SET_ERR( MB_FAILURE, "Incorrect number of halfedges." ); 03353 03354 int fidx = ID_FROM_HANDLE( ent ) - 1; 03355 for( int i = 0; i < nepf; i++ ) 03356 { 03357 HFacet hf = sibhes[nepf * fidx + i]; 03358 sib_entids[i] = fid_from_halfacet( hf, type ); 03359 sib_lids[i] = lid_from_halffacet( hf ); 03360 } 03361 } 03362 else 03363 { 03364 int idx = get_index_in_lmap( *_cells.begin() ); 03365 int nfpc = lConnMap3D[idx].num_faces_in_cell; 03366 03367 if( num_halffacets != nfpc ) MB_SET_ERR( MB_FAILURE, "Incorrect number of halffaces." ); 03368 03369 int cidx = ID_FROM_HANDLE( ent ) - 1; 03370 for( int i = 0; i < nfpc; i++ ) 03371 { 03372 HFacet hf = sibhfs[nfpc * cidx + i]; 03373 sib_entids[i] = fid_from_halfacet( hf, type ); 03374 sib_lids[i] = lid_from_halffacet( hf ); 03375 } 03376 } 03377 return MB_SUCCESS; 03378 } 03379 03380 ErrorCode HalfFacetRep::get_sibling_map( EntityType type, EntityHandle ent, int lid, EntityHandle& sib_entid, 03381 int& sib_lid ) 03382 { 03383 03384 if( type == MBEDGE ) 03385 { 03386 int eidx = ID_FROM_HANDLE( ent ) - 1; 03387 HFacet hf = sibhvs[2 * eidx + lid]; 03388 sib_entid = fid_from_halfacet( hf, MBEDGE ); 03389 sib_lid = lid_from_halffacet( hf ); 03390 } 03391 else if( type == MBTRI || type == MBQUAD ) 03392 { 03393 int nepf = lConnMap2D[type - 2].num_verts_in_face; 03394 int fidx = ID_FROM_HANDLE( ent ) - 1; 03395 HFacet hf = sibhes[nepf * fidx + lid]; 03396 sib_entid = fid_from_halfacet( hf, type ); 03397 sib_lid = lid_from_halffacet( hf ); 03398 } 03399 else 03400 { 03401 int idx = get_index_in_lmap( *_cells.begin() ); 03402 int nfpc = lConnMap3D[idx].num_faces_in_cell; 03403 int cidx = ID_FROM_HANDLE( ent ) - 1; 03404 HFacet hf = sibhfs[nfpc * cidx + lid]; 03405 sib_entid = fid_from_halfacet( hf, type ); 03406 sib_lid = lid_from_halffacet( hf ); 03407 } 03408 return MB_SUCCESS; 03409 } 03410 03411 ErrorCode HalfFacetRep::set_sibling_map( EntityType type, EntityHandle ent, EntityHandle* set_entids, int* set_lids, 03412 int num_halffacets ) 03413 { 03414 03415 if( type == MBEDGE ) 03416 { 03417 if( num_halffacets != 2 ) MB_SET_ERR( MB_FAILURE, "Incorrect number of halfvertices" ); 03418 03419 int eidx = ID_FROM_HANDLE( ent ) - 1; 03420 for( int i = 0; i < 2; i++ ) 03421 { 03422 sibhvs[2 * eidx + i] = create_halffacet( set_entids[i], set_lids[i] ); 03423 } 03424 } 03425 else if( type == MBTRI || type == MBQUAD ) 03426 { 03427 int nepf = lConnMap2D[type - 2].num_verts_in_face; 03428 if( num_halffacets != nepf ) MB_SET_ERR( MB_FAILURE, "Incorrect number of halfedges." ); 03429 03430 int fidx = ID_FROM_HANDLE( ent ) - 1; 03431 for( int i = 0; i < nepf; i++ ) 03432 { 03433 sibhes[nepf * fidx + i] = create_halffacet( set_entids[i], set_lids[i] ); 03434 } 03435 } 03436 else 03437 { 03438 int idx = get_index_in_lmap( *_cells.begin() ); 03439 int nfpc = lConnMap3D[idx].num_faces_in_cell; 03440 if( num_halffacets != nfpc ) MB_SET_ERR( MB_FAILURE, "Incorrect number of halffaces." ); 03441 03442 int cidx = ID_FROM_HANDLE( ent ) - 1; 03443 for( int i = 0; i < nfpc; i++ ) 03444 { 03445 sibhfs[nfpc * cidx + i] = create_halffacet( set_entids[i], set_lids[i] ); 03446 } 03447 } 03448 03449 return MB_SUCCESS; 03450 } 03451 03452 ErrorCode HalfFacetRep::set_sibling_map( EntityType type, EntityHandle ent, int lid, EntityHandle& set_entid, 03453 int& set_lid ) 03454 { 03455 03456 if( type == MBEDGE ) 03457 { 03458 int eidx = ID_FROM_HANDLE( ent ) - 1; 03459 sibhvs[2 * eidx + lid] = create_halffacet( set_entid, set_lid ); 03460 } 03461 else if( type == MBTRI || type == MBQUAD ) 03462 { 03463 int nepf = lConnMap2D[type - 2].num_verts_in_face; 03464 int fidx = ID_FROM_HANDLE( ent ) - 1; 03465 sibhes[nepf * fidx + lid] = create_halffacet( set_entid, set_lid ); 03466 } 03467 else 03468 { 03469 int idx = get_index_in_lmap( *_cells.begin() ); 03470 int nfpc = lConnMap3D[idx].num_faces_in_cell; 03471 int cidx = ID_FROM_HANDLE( ent ) - 1; 03472 sibhfs[nfpc * cidx + lid] = create_halffacet( set_entid, set_lid ); 03473 } 03474 03475 return MB_SUCCESS; 03476 } 03477 03478 ErrorCode HalfFacetRep::get_incident_map( EntityType type, EntityHandle vid, std::vector< EntityHandle >& inci_entid, 03479 std::vector< int >& inci_lid ) 03480 { 03481 inci_entid.clear(); 03482 inci_lid.clear(); 03483 03484 if( type == MBEDGE ) 03485 { 03486 HFacet hf = v2hv[ID_FROM_HANDLE( vid ) - 1]; 03487 inci_entid.push_back( fid_from_halfacet( hf, type ) ); 03488 inci_lid.push_back( lid_from_halffacet( hf ) ); 03489 } 03490 else if( type == MBTRI || type == MBQUAD ) 03491 { 03492 HFacet hf = v2he[ID_FROM_HANDLE( vid ) - 1]; 03493 if( hf == 0 && ( v2hes.find( vid ) != v2hes.end() ) ) 03494 { 03495 std::pair< std::multimap< EntityHandle, HFacet >::iterator, 03496 std::multimap< EntityHandle, HFacet >::iterator > 03497 it_hes; 03498 it_hes = v2hes.equal_range( vid ); 03499 03500 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 03501 { 03502 inci_entid.push_back( fid_from_halfacet( it->second, type ) ); 03503 inci_lid.push_back( lid_from_halffacet( it->second ) ); 03504 } 03505 } 03506 else if( hf != 0 ) 03507 { 03508 inci_entid.push_back( fid_from_halfacet( hf, type ) ); 03509 inci_lid.push_back( lid_from_halffacet( hf ) ); 03510 } 03511 else if( hf == 0 && ( v2hes.find( vid ) == v2hes.end() ) ) 03512 { 03513 inci_entid.push_back( fid_from_halfacet( hf, type ) ); 03514 inci_lid.push_back( lid_from_halffacet( hf ) ); 03515 } 03516 } 03517 else 03518 { 03519 HFacet hf = v2hf[ID_FROM_HANDLE( vid ) - 1]; 03520 if( hf == 0 && ( v2hfs.find( vid ) != v2hfs.end() ) ) 03521 { 03522 std::pair< std::multimap< EntityHandle, HFacet >::iterator, 03523 std::multimap< EntityHandle, HFacet >::iterator > 03524 it_hes; 03525 it_hes = v2hfs.equal_range( vid ); 03526 03527 for( std::multimap< EntityHandle, HFacet >::iterator it = it_hes.first; it != it_hes.second; ++it ) 03528 { 03529 inci_entid.push_back( fid_from_halfacet( it->second, type ) ); 03530 inci_lid.push_back( lid_from_halffacet( it->second ) ); 03531 } 03532 } 03533 else if( hf != 0 ) 03534 { 03535 inci_entid.push_back( fid_from_halfacet( hf, type ) ); 03536 inci_lid.push_back( lid_from_halffacet( hf ) ); 03537 } 03538 else if( hf == 0 && ( v2hfs.find( vid ) == v2hfs.end() ) ) 03539 { 03540 inci_entid.push_back( fid_from_halfacet( hf, type ) ); 03541 inci_lid.push_back( lid_from_halffacet( hf ) ); 03542 } 03543 } 03544 03545 return MB_SUCCESS; 03546 } 03547 03548 ErrorCode HalfFacetRep::set_incident_map( EntityType type, EntityHandle vid, std::vector< EntityHandle >& set_entid, 03549 std::vector< int >& set_lid ) 03550 { 03551 if( type == MBEDGE ) { v2hv[ID_FROM_HANDLE( vid ) - 1] = create_halffacet( set_entid[0], set_lid[0] ); } 03552 else if( type == MBTRI || type == MBQUAD ) 03553 { 03554 if( set_entid.size() == 1 ) 03555 v2he[ID_FROM_HANDLE( vid ) - 1] = create_halffacet( set_entid[0], set_lid[0] ); 03556 else 03557 { 03558 HFacet hf = 0; 03559 for( int i = 0; i < (int)set_entid.size(); i++ ) 03560 { 03561 hf = create_halffacet( set_entid[i], set_lid[i] ); 03562 v2hes.insert( std::pair< EntityHandle, HFacet >( vid, hf ) ); 03563 } 03564 } 03565 } 03566 else 03567 { 03568 if( set_entid.size() == 1 ) 03569 v2hf[ID_FROM_HANDLE( vid ) - 1] = create_halffacet( set_entid[0], set_lid[0] ); 03570 else 03571 { 03572 HFacet hf = v2hf[ID_FROM_HANDLE( vid ) - 1]; 03573 if( hf != 0 ) { v2hf[ID_FROM_HANDLE( vid ) - 1] = 0; } 03574 for( int i = 0; i < (int)set_entid.size(); i++ ) 03575 { 03576 hf = create_halffacet( set_entid[i], set_lid[i] ); 03577 v2hfs.insert( std::pair< EntityHandle, HFacet >( vid, hf ) ); 03578 } 03579 } 03580 } 03581 03582 return MB_SUCCESS; 03583 } 03584 03585 ErrorCode HalfFacetRep::get_entity_ranges( Range& verts, Range& edges, Range& faces, Range& cells ) 03586 { 03587 verts = _verts; 03588 edges = _edges; 03589 faces = _faces; 03590 cells = _cells; 03591 return MB_SUCCESS; 03592 } 03593 03594 ErrorCode HalfFacetRep::update_entity_ranges( EntityHandle fileset ) 03595 { 03596 ErrorCode error; 03597 03598 error = mb->get_entities_by_dimension( fileset, 0, _verts, true );MB_CHK_ERR( error ); 03599 03600 error = mb->get_entities_by_dimension( fileset, 1, _edges, true );MB_CHK_ERR( error ); 03601 03602 error = mb->get_entities_by_dimension( fileset, 2, _faces, true );MB_CHK_ERR( error ); 03603 03604 error = mb->get_entities_by_dimension( fileset, 3, _cells, true );MB_CHK_ERR( error ); 03605 03606 return MB_SUCCESS; 03607 } 03608 03609 void HalfFacetRep::get_memory_use( unsigned long long& entity_total, unsigned long long& memory_total ) 03610 { 03611 entity_total = memory_total = 0; 03612 // 1D 03613 if( !v2hv.empty() ) entity_total += v2hv.capacity() * sizeof( HFacet ) + sizeof( v2hv ); 03614 if( !sibhvs.empty() ) entity_total += sibhvs.capacity() * sizeof( HFacet ) + sizeof( sibhvs ); 03615 03616 // 2D 03617 if( !v2he.empty() ) entity_total += v2he.capacity() * sizeof( HFacet ) + sizeof( v2he ); 03618 if( !sibhes.empty() ) entity_total += sibhes.capacity() * sizeof( HFacet ) + sizeof( sibhes ); 03619 03620 // 3D 03621 if( !v2hf.empty() ) entity_total += v2hf.capacity() * sizeof( HFacet ) + sizeof( v2hf ); 03622 if( !sibhfs.empty() ) entity_total += sibhfs.capacity() * sizeof( HFacet ) + sizeof( sibhfs ); 03623 03624 memory_total = entity_total; 03625 } 03626 03627 HFacet HalfFacetRep::create_halffacet( EntityHandle handle, int lid ) 03628 { 03629 EntityID fid = ID_FROM_HANDLE( handle ); 03630 return CREATE_HALFFACET( lid, fid ); 03631 } 03632 03633 EntityHandle HalfFacetRep::fid_from_halfacet( const HFacet facet, EntityType type ) 03634 { 03635 EntityID id = FID_FROM_HALFFACET( facet ); 03636 EntityHandle handle = 0; 03637 if( id == 0 ) return handle; 03638 03639 ErrorCode error = mb->handle_from_id( type, id, handle );MB_CHK_ERR( error ); 03640 return handle; 03641 } 03642 03643 int HalfFacetRep::lid_from_halffacet( const HFacet facet ) 03644 { 03645 if( facet == 0 ) 03646 return 0; 03647 else 03648 return LID_FROM_HALFFACET( facet ); 03649 } 03650 03651 } // namespace moab
1.7.6.1 
