MOAB: Mesh Oriented datABase
(version 5.4.1)
|
00001 #include <cstdlib> // For exit() 00002 #include <vector> 00003 #include <map> 00004 #include <iostream> 00005 #include <string> 00006 #include <algorithm> 00007 00008 #include "moab/CN.hpp" 00009 #include "moab/Range.hpp" 00010 #include "moab/Interface.hpp" 00011 #include "MBTagConventions.hpp" 00012 #include "Internals.hpp" 00013 #include "moab/ReadUtilIface.hpp" 00014 #include "moab/FileOptions.hpp" 00015 #include "ReadCCMIO.hpp" 00016 #include "moab/MeshTopoUtil.hpp" 00017 00018 #include "ccmio.h" 00019 00020 /* 00021 * CCMIO file structure 00022 * 00023 * Root 00024 * State(kCCMIOState) 00025 * Processor* 00026 * Vertices 00027 * ->ReadVerticesx, ReadMap 00028 * Topology 00029 * Boundary faces*(kCCMIOBoundaryFaces) 00030 * ->ReadFaces, ReadFaceCells, ReadMap 00031 * Internal faces(kCCMIOInternalFaces) 00032 * Cells (kCCMIOCells) 00033 * ->ReadCells (mapID), ReadMap, ReadCells (cellTypes) 00034 * Solution 00035 * Phase 00036 * Field 00037 * FieldData 00038 * Problem(kCCMIOProblemDescription) 00039 * CellType* (kCCMIOCellType) 00040 * Index (GetEntityIndex), MaterialId(ReadOpti), MaterialType(ReadOptstr), 00041 * PorosityId(ReadOpti), SpinId(ReadOpti), GroupId(ReadOpti) 00042 * 00043 * MaterialType (CCMIOReadOptstr in readexample) 00044 * constants (see readexample) 00045 * lagrangian data (CCMIOReadLagrangianData) 00046 * vertices label (CCMIOEntityDescription) 00047 * restart info: char solver[], iteratoins, time, char timeUnits[], angle 00048 * (CCMIOReadRestartInfo, kCCMIORestartData), reference data? 00049 * phase: 00050 * field: char name[], dims, CCMIODataType datatype, char units[] 00051 * dims = kCCMIOScalar (CCMIOReadFieldDataf), 00052 * kCCMIOVector (CCMIOReadMultiDimensionalFieldData), 00053 * kCCMIOTensor 00054 * MonitoringSets: num, name (CellSet, VertexSet, BoundarySet, BlockSet, SplineSet, CoupleSet) 00055 * CCMIOGetProstarSet, CCMIOReadOpt1i, 00056 */ 00057 00058 enum DataType 00059 { 00060 kScalar, 00061 kVector, 00062 kVertex, 00063 kCell, 00064 kInternalFace, 00065 kBoundaryFace, 00066 kBoundaryData, 00067 kBoundaryFaceData, 00068 kCellType 00069 }; 00070 00071 namespace moab 00072 { 00073 00074 static int const kNValues = 10; // Number of values of each element to print 00075 static char const kDefaultState[] = "default"; 00076 static char const kUnitsName[] = "Units"; 00077 static int const kVertOffset = 2; 00078 static int const kCellInc = 4; 00079 00080 #define CHK_SET_CCMERR( ccm_err_code, ccm_err_msg ) \ 00081 { \ 00082 if( kCCMIONoErr != ( ccm_err_code ) && kCCMIONoFileErr != ( ccm_err_code ) && \ 00083 kCCMIONoNodeErr != ( ccm_err_code ) ) \ 00084 MB_SET_ERR( MB_FAILURE, ccm_err_msg ); \ 00085 } 00086 00087 ReaderIface* ReadCCMIO::factory( Interface* iface ) 00088 { 00089 return new ReadCCMIO( iface ); 00090 } 00091 00092 ReadCCMIO::ReadCCMIO( Interface* impl ) 00093 : mMaterialIdTag( 0 ), mMaterialTypeTag( 0 ), mRadiationTag( 0 ), mPorosityIdTag( 0 ), mSpinIdTag( 0 ), 00094 mGroupIdTag( 0 ), mColorIdxTag( 0 ), mProcessorIdTag( 0 ), mLightMaterialTag( 0 ), mFreeSurfaceMaterialTag( 0 ), 00095 mThicknessTag( 0 ), mProstarRegionNumberTag( 0 ), mBoundaryTypeTag( 0 ), mCreatingProgramTag( 0 ), mbImpl( impl ), 00096 hasSolution( false ) 00097 { 00098 assert( impl != NULL ); 00099 00100 impl->query_interface( readMeshIface ); 00101 00102 // Initialize in case tag_get_handle fails below 00103 mMaterialSetTag = 0; 00104 mDirichletSetTag = 0; 00105 mNeumannSetTag = 0; 00106 mHasMidNodesTag = 0; 00107 mGlobalIdTag = 0; 00108 mNameTag = 0; 00109 00110 //! Get and cache predefined tag handles 00111 const int negone = -1; 00112 ErrorCode result = impl->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mMaterialSetTag, 00113 MB_TAG_CREAT | MB_TAG_SPARSE, &negone );MB_CHK_SET_ERR_RET( result, "Failed to get MATERIAL_SET tag" ); 00114 00115 result = impl->tag_get_handle( DIRICHLET_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mDirichletSetTag, 00116 MB_TAG_CREAT | MB_TAG_SPARSE, &negone );MB_CHK_SET_ERR_RET( result, "Failed to get DIRICHLET_SET tag" ); 00117 00118 result = impl->tag_get_handle( NEUMANN_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mNeumannSetTag, 00119 MB_TAG_CREAT | MB_TAG_SPARSE, &negone );MB_CHK_SET_ERR_RET( result, "Failed to get NEUMANN_SET tag" ); 00120 00121 const int negonearr[] = { -1, -1, -1, -1 }; 00122 result = impl->tag_get_handle( HAS_MID_NODES_TAG_NAME, 4, MB_TYPE_INTEGER, mHasMidNodesTag, 00123 MB_TAG_CREAT | MB_TAG_SPARSE, negonearr );MB_CHK_SET_ERR_RET( result, "Failed to get HAS_MID_NODES tag" ); 00124 00125 mGlobalIdTag = impl->globalId_tag(); 00126 00127 result = 00128 impl->tag_get_handle( NAME_TAG_NAME, NAME_TAG_SIZE, MB_TYPE_OPAQUE, mNameTag, MB_TAG_CREAT | MB_TAG_SPARSE );MB_CHK_SET_ERR_RET( result, "Failed to get NAME tag" ); 00129 } 00130 00131 ReadCCMIO::~ReadCCMIO() 00132 { 00133 mbImpl->release_interface( readMeshIface ); 00134 } 00135 00136 ErrorCode ReadCCMIO::load_file( const char* file_name, 00137 const EntityHandle* file_set, 00138 const FileOptions& /* opts */, 00139 const ReaderIface::SubsetList* subset_list, 00140 const Tag* /* file_id_tag */ ) 00141 { 00142 CCMIOID rootID, problemID, stateID, processorID, verticesID, topologyID, solutionID; 00143 CCMIOError error = kCCMIONoErr; 00144 00145 if( subset_list ) 00146 { 00147 MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "Reading subset of files not supported for CCMOI data" ); 00148 } 00149 00150 CCMIOOpenFile( &error, file_name, kCCMIORead, &rootID ); 00151 CHK_SET_CCMERR( error, "Problem opening file" ); 00152 00153 // Get the file state 00154 ErrorCode rval = get_state( rootID, problemID, stateID );MB_CHK_SET_ERR( rval, "Failed to get state" ); 00155 00156 // Get processors 00157 std::vector< CCMIOSize_t > procs; 00158 bool has_solution = false; 00159 rval = get_processors( stateID, processorID, verticesID, topologyID, solutionID, procs, has_solution );MB_CHK_SET_ERR( rval, "Failed to get processors" ); 00160 00161 std::vector< CCMIOSize_t >::iterator vit; 00162 Range new_ents, *new_ents_ptr = NULL; 00163 if( file_set ) new_ents_ptr = &new_ents; 00164 00165 for( vit = procs.begin(); vit != procs.end(); ++vit ) 00166 { 00167 rval = read_processor( stateID, problemID, processorID, verticesID, topologyID, *vit, new_ents_ptr );MB_CHK_SET_ERR( rval, "Failed to read processors" ); 00168 } 00169 00170 // Load some meta-data 00171 rval = load_metadata( rootID, problemID, stateID, processorID, file_set );MB_CHK_SET_ERR( rval, "Failed to load some meta-data" ); 00172 00173 // Now, put all this into the file set, if there is one 00174 if( file_set ) 00175 { 00176 rval = mbImpl->add_entities( *file_set, new_ents );MB_CHK_SET_ERR( rval, "Failed to add new entities to file set" ); 00177 } 00178 00179 return rval; 00180 } 00181 00182 ErrorCode ReadCCMIO::get_state( CCMIOID rootID, CCMIOID& problemID, CCMIOID& stateID ) 00183 { 00184 CCMIOError error = kCCMIONoErr; 00185 00186 // First try default 00187 CCMIOGetState( &error, rootID, "default", &problemID, &stateID ); 00188 if( kCCMIONoErr != error ) 00189 { 00190 CCMIOSize_t i = CCMIOSIZEC( 0 ); 00191 CCMIOError tmp_error = kCCMIONoErr; 00192 CCMIONextEntity( &tmp_error, rootID, kCCMIOState, &i, &stateID ); 00193 if( kCCMIONoErr == tmp_error ) CCMIONextEntity( &error, rootID, kCCMIOProblemDescription, &i, &problemID ); 00194 } 00195 CHK_SET_CCMERR( error, "Couldn't find state" ); 00196 00197 return MB_SUCCESS; 00198 } 00199 00200 ErrorCode ReadCCMIO::load_metadata( CCMIOID rootID, 00201 CCMIOID problemID, 00202 CCMIOID /* stateID */, 00203 CCMIOID processorID, 00204 const EntityHandle* file_set ) 00205 { 00206 // Read the simulation title. 00207 CCMIOError error = kCCMIONoErr; 00208 ErrorCode rval = MB_SUCCESS; 00209 CCMIONode rootNode; 00210 if( kCCMIONoErr == CCMIOGetEntityNode( &error, rootID, &rootNode ) ) 00211 { 00212 char* name = NULL; 00213 CCMIOGetTitle( &error, rootNode, &name ); 00214 00215 if( NULL != name && strlen( name ) != 0 ) 00216 { 00217 // Make a tag for it and tag the read set 00218 Tag simname; 00219 rval = mbImpl->tag_get_handle( "Title", strlen( name ), MB_TYPE_OPAQUE, simname, 00220 MB_TAG_CREAT | MB_TAG_SPARSE );MB_CHK_SET_ERR( rval, "Simulation name tag not found or created" ); 00221 EntityHandle set = file_set ? *file_set : 0; 00222 rval = mbImpl->tag_set_data( simname, &set, 1, name );MB_CHK_SET_ERR( rval, "Problem setting simulation name tag" ); 00223 } 00224 if( name ) free( name ); 00225 } 00226 00227 // Creating program 00228 EntityHandle dumh = ( file_set ? *file_set : 0 ); 00229 rval = get_str_option( "CreatingProgram", dumh, mCreatingProgramTag, processorID );MB_CHK_SET_ERR( rval, "Trouble getting CreatingProgram tag" ); 00230 00231 rval = load_matset_data( problemID );MB_CHK_SET_ERR( rval, "Failure loading matset data" ); 00232 00233 rval = load_neuset_data( problemID );MB_CHK_SET_ERR( rval, "Failure loading neuset data" ); 00234 00235 return rval; 00236 } 00237 00238 ErrorCode ReadCCMIO::load_matset_data( CCMIOID problemID ) 00239 { 00240 // Make sure there are matsets 00241 if( newMatsets.empty() ) return MB_SUCCESS; 00242 00243 // ... walk through each cell type 00244 CCMIOSize_t i = CCMIOSIZEC( 0 ); 00245 CCMIOID next; 00246 CCMIOError error = kCCMIONoErr; 00247 00248 while( CCMIONextEntity( NULL, problemID, kCCMIOCellType, &i, &next ) == kCCMIONoErr ) 00249 { 00250 // Get index, corresponding set, and label with material set tag 00251 int mindex; 00252 CCMIOGetEntityIndex( &error, next, &mindex ); 00253 std::map< int, EntityHandle >::iterator mit = newMatsets.find( mindex ); 00254 if( mit == newMatsets.end() ) 00255 // No actual faces for this matset; continue to next 00256 continue; 00257 00258 EntityHandle dum_ent = mit->second; 00259 ErrorCode rval = mbImpl->tag_set_data( mMaterialSetTag, &dum_ent, 1, &mindex );MB_CHK_SET_ERR( rval, "Trouble setting material set tag" ); 00260 00261 // Set name 00262 CCMIOSize_t len; 00263 CCMIOEntityLabel( &error, next, &len, NULL ); 00264 std::vector< char > opt_string2( GETINT32( len ) + 1, '\0' ); 00265 CCMIOEntityLabel( &error, next, NULL, &opt_string2[0] ); 00266 if( opt_string2.size() >= NAME_TAG_SIZE ) 00267 opt_string2[NAME_TAG_SIZE - 1] = '\0'; 00268 else 00269 ( opt_string2.resize( NAME_TAG_SIZE, '\0' ) ); 00270 rval = mbImpl->tag_set_data( mNameTag, &dum_ent, 1, &opt_string2[0] );MB_CHK_SET_ERR( rval, "Trouble setting name tag for material set" ); 00271 00272 // Material id 00273 rval = get_int_option( "MaterialId", dum_ent, mMaterialIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting MaterialId tag" ); 00274 00275 rval = get_str_option( "MaterialType", dum_ent, mMaterialTypeTag, next );MB_CHK_SET_ERR( rval, "Trouble getting MaterialType tag" ); 00276 00277 rval = get_int_option( "Radiation", dum_ent, mRadiationTag, next );MB_CHK_SET_ERR( rval, "Trouble getting Radiation option" ); 00278 00279 rval = get_int_option( "PorosityId", dum_ent, mPorosityIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting PorosityId option" ); 00280 00281 rval = get_int_option( "SpinId", dum_ent, mSpinIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting SpinId option" ); 00282 00283 rval = get_int_option( "GroupId", dum_ent, mGroupIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting GroupId option" ); 00284 00285 rval = get_int_option( "ColorIdx", dum_ent, mColorIdxTag, next );MB_CHK_SET_ERR( rval, "Trouble getting ColorIdx option" ); 00286 00287 rval = get_int_option( "ProcessorId", dum_ent, mProcessorIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting ProcessorId option" ); 00288 00289 rval = get_int_option( "LightMaterial", dum_ent, mLightMaterialTag, next );MB_CHK_SET_ERR( rval, "Trouble getting LightMaterial option" ); 00290 00291 rval = get_int_option( "FreeSurfaceMaterial", dum_ent, mFreeSurfaceMaterialTag, next );MB_CHK_SET_ERR( rval, "Trouble getting FreeSurfaceMaterial option" ); 00292 00293 rval = get_dbl_option( "Thickness", dum_ent, mThicknessTag, next );MB_CHK_SET_ERR( rval, "Trouble getting Thickness option" ); 00294 } 00295 00296 return MB_SUCCESS; 00297 } 00298 00299 ErrorCode ReadCCMIO::get_int_option( const char* opt_str, EntityHandle seth, Tag& tag, CCMIOID node ) 00300 { 00301 int idum; 00302 ErrorCode rval; 00303 if( kCCMIONoErr == CCMIOReadOpti( NULL, node, opt_str, &idum ) ) 00304 { 00305 if( !tag ) 00306 { 00307 rval = mbImpl->tag_get_handle( opt_str, 1, MB_TYPE_INTEGER, tag, MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Failed to get tag handle" ); 00308 } 00309 00310 rval = mbImpl->tag_set_data( tag, &seth, 1, &idum );MB_CHK_SET_ERR( rval, "Failed to set tag data" ); 00311 } 00312 00313 return MB_SUCCESS; 00314 } 00315 00316 ErrorCode ReadCCMIO::get_dbl_option( const char* opt_str, EntityHandle seth, Tag& tag, CCMIOID node ) 00317 { 00318 float fdum; 00319 if( kCCMIONoErr == CCMIOReadOptf( NULL, node, opt_str, &fdum ) ) 00320 { 00321 ErrorCode rval; 00322 if( !tag ) 00323 { 00324 rval = mbImpl->tag_get_handle( opt_str, 1, MB_TYPE_DOUBLE, tag, MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Failed to get tag handle" ); 00325 } 00326 00327 double dum_dbl = fdum; 00328 rval = mbImpl->tag_set_data( tag, &seth, 1, &dum_dbl );MB_CHK_SET_ERR( rval, "Failed to set tag data" ); 00329 } 00330 00331 return MB_SUCCESS; 00332 } 00333 00334 ErrorCode ReadCCMIO::get_str_option( const char* opt_str, 00335 EntityHandle seth, 00336 Tag& tag, 00337 CCMIOID node, 00338 const char* other_tag_name ) 00339 { 00340 int len; 00341 CCMIOError error = kCCMIONoErr; 00342 std::vector< char > opt_string; 00343 if( kCCMIONoErr != CCMIOReadOptstr( NULL, node, opt_str, &len, NULL ) ) return MB_SUCCESS; 00344 00345 opt_string.resize( len ); 00346 CCMIOReadOptstr( &error, node, opt_str, &len, &opt_string[0] ); 00347 ErrorCode rval = MB_SUCCESS; 00348 if( !tag ) 00349 { 00350 rval = mbImpl->tag_get_handle( other_tag_name ? other_tag_name : opt_str, NAME_TAG_SIZE, MB_TYPE_OPAQUE, tag, 00351 MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Failed to get tag handle" ); 00352 } 00353 00354 if( opt_string.size() > NAME_TAG_SIZE ) 00355 opt_string[NAME_TAG_SIZE - 1] = '\0'; 00356 else 00357 ( opt_string.resize( NAME_TAG_SIZE, '\0' ) ); 00358 00359 rval = mbImpl->tag_set_data( tag, &seth, 1, &opt_string[0] );MB_CHK_SET_ERR( rval, "Failed to set tag data" ); 00360 00361 return MB_SUCCESS; 00362 } 00363 00364 ErrorCode ReadCCMIO::load_neuset_data( CCMIOID problemID ) 00365 { 00366 CCMIOSize_t i = CCMIOSIZEC( 0 ); 00367 CCMIOID next; 00368 00369 // Make sure there are matsets 00370 if( newNeusets.empty() ) return MB_SUCCESS; 00371 00372 while( CCMIONextEntity( NULL, problemID, kCCMIOBoundaryRegion, &i, &next ) == kCCMIONoErr ) 00373 { 00374 // Get index, corresponding set, and label with neumann set tag 00375 int mindex; 00376 CCMIOError error = kCCMIONoErr; 00377 CCMIOGetEntityIndex( &error, next, &mindex ); 00378 std::map< int, EntityHandle >::iterator mit = newNeusets.find( mindex ); 00379 if( mit == newNeusets.end() ) 00380 // No actual faces for this neuset; continue to next 00381 continue; 00382 00383 EntityHandle dum_ent = mit->second; 00384 ErrorCode rval = mbImpl->tag_set_data( mNeumannSetTag, &dum_ent, 1, &mindex );MB_CHK_SET_ERR( rval, "Trouble setting neumann set tag" ); 00385 00386 // Set name 00387 rval = get_str_option( "BoundaryName", dum_ent, mNameTag, next, NAME_TAG_NAME );MB_CHK_SET_ERR( rval, "Trouble creating BoundaryName tag" ); 00388 00389 // BoundaryType 00390 rval = get_str_option( "BoundaryType", dum_ent, mBoundaryTypeTag, next );MB_CHK_SET_ERR( rval, "Trouble creating BoundaryType tag" ); 00391 00392 // ProstarRegionNumber 00393 rval = get_int_option( "ProstarRegionNumber", dum_ent, mProstarRegionNumberTag, next );MB_CHK_SET_ERR( rval, "Trouble creating ProstarRegionNumber tag" ); 00394 } 00395 00396 return MB_SUCCESS; 00397 } 00398 00399 ErrorCode ReadCCMIO::read_processor( CCMIOID /* stateID */, 00400 CCMIOID problemID, 00401 CCMIOID processorID, 00402 CCMIOID verticesID, 00403 CCMIOID topologyID, 00404 CCMIOSize_t proc, 00405 Range* new_ents ) 00406 { 00407 ErrorCode rval; 00408 00409 // vert_map fields: s: none, i: gid, ul: vert handle, r: none 00410 // TupleList vert_map(0, 1, 1, 0, 0); 00411 TupleList vert_map; 00412 rval = read_vertices( proc, processorID, verticesID, topologyID, new_ents, vert_map );MB_CHK_SET_ERR( rval, "Failed to read vertices" ); 00413 00414 rval = read_cells( proc, problemID, verticesID, topologyID, vert_map, new_ents );MB_CHK_SET_ERR( rval, "Failed to read cells" ); 00415 00416 return rval; 00417 } 00418 00419 ErrorCode ReadCCMIO::read_cells( CCMIOSize_t /* proc */, 00420 CCMIOID problemID, 00421 CCMIOID /* verticesID */, 00422 CCMIOID topologyID, 00423 TupleList& vert_map, 00424 Range* new_ents ) 00425 { 00426 // Read the faces. 00427 // face_map fields: s:forward/reverse, i: cell id, ul: face handle, r: none 00428 ErrorCode rval; 00429 #ifdef TUPLE_LIST 00430 TupleList face_map( 1, 1, 1, 0, 0 ); 00431 #else 00432 TupleList face_map; 00433 SenseList sense_map; 00434 #endif 00435 rval = read_all_faces( topologyID, vert_map, face_map, 00436 #ifndef TUPLE_LIST 00437 sense_map, 00438 #endif 00439 new_ents );MB_CHK_SET_ERR( rval, "Failed to read all cells" ); 00440 00441 // Read the cell topology types, if any exist in the file 00442 std::map< int, int > cell_topo_types; 00443 rval = read_topology_types( topologyID, cell_topo_types );MB_CHK_SET_ERR( rval, "Problem reading cell topo types" ); 00444 00445 // Now construct the cells; sort the face map by cell ids first 00446 #ifdef TUPLE_LIST 00447 rval = face_map.sort( 1 );MB_CHK_SET_ERR( rval, "Couldn't sort face map by cell id" ); 00448 #endif 00449 std::vector< EntityHandle > new_cells; 00450 rval = construct_cells( face_map, 00451 #ifndef TUPLE_LIST 00452 sense_map, 00453 #endif 00454 vert_map, cell_topo_types, new_cells );MB_CHK_SET_ERR( rval, "Failed to construct cells" ); 00455 if( new_ents ) 00456 { 00457 Range::iterator rit = new_ents->end(); 00458 std::vector< EntityHandle >::reverse_iterator vit; 00459 for( vit = new_cells.rbegin(); vit != new_cells.rend(); ++vit ) 00460 rit = new_ents->insert( rit, *vit ); 00461 } 00462 00463 rval = read_gids_and_types( problemID, topologyID, new_cells );MB_CHK_SET_ERR( rval, "Failed to read gids and types" ); 00464 00465 return MB_SUCCESS; 00466 } 00467 00468 ErrorCode ReadCCMIO::read_topology_types( CCMIOID& topologyID, std::map< int, int >& cell_topo_types ) 00469 { 00470 CCMIOError error = kCCMIONoErr; 00471 CCMIOID cellID, mapID; 00472 CCMIOSize_t ncells; 00473 CCMIOGetEntity( &error, topologyID, kCCMIOCells, 0, &cellID ); 00474 CCMIOEntitySize( &error, cellID, &ncells, NULL ); 00475 int num_cells = GETINT32( ncells ); 00476 00477 // First, do a dummy read to see if we even have topo types in this mesh 00478 int dum_int; 00479 CCMIOReadOpt1i( &error, cellID, "CellTopologyType", &dum_int, CCMIOINDEXC( kCCMIOStart ), 00480 CCMIOINDEXC( kCCMIOStart ) + 1 ); 00481 if( kCCMIONoErr != error ) return MB_SUCCESS; 00482 00483 // OK, we have topo types; first get the map node 00484 std::vector< int > dum_ints( num_cells ); 00485 CCMIOReadCells( &error, cellID, &mapID, &dum_ints[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOStart ) + 1 ); 00486 CHK_SET_CCMERR( error, "Failed to get the map node" ); 00487 00488 // Now read the map 00489 CCMIOReadMap( &error, mapID, &dum_ints[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) ); 00490 CHK_SET_CCMERR( error, "Failed to get cell ids" ); 00491 int i; 00492 for( i = 0; i < num_cells; i++ ) 00493 cell_topo_types[dum_ints[i]] = 0; 00494 00495 // Now read the cell topo types for real, reusing cell_topo_types 00496 std::vector< int > topo_types( num_cells ); 00497 CCMIOReadOpt1i( &error, cellID, "CellTopologyType", &topo_types[0], CCMIOINDEXC( kCCMIOStart ), 00498 CCMIOINDEXC( kCCMIOEnd ) ); 00499 CHK_SET_CCMERR( error, "Failed to get cell topo types" ); 00500 for( i = 0; i < num_cells; i++ ) 00501 cell_topo_types[dum_ints[i]] = topo_types[i]; 00502 00503 return MB_SUCCESS; 00504 } 00505 00506 ErrorCode ReadCCMIO::read_gids_and_types( CCMIOID /* problemID */, 00507 CCMIOID topologyID, 00508 std::vector< EntityHandle >& cells ) 00509 { 00510 // Get the cells entity and number of cells 00511 CCMIOSize_t dum_cells; 00512 int num_cells; 00513 CCMIOError error = kCCMIONoErr; 00514 CCMIOID cellsID, mapID; 00515 CCMIOGetEntity( &error, topologyID, kCCMIOCells, 0, &cellsID ); 00516 CCMIOEntitySize( &error, cellsID, &dum_cells, NULL ); 00517 num_cells = GETINT32( dum_cells ); 00518 00519 // Check the number of cells against how many are in the cell array 00520 if( num_cells != (int)cells.size() ) MB_SET_ERR( MB_FAILURE, "Number of cells doesn't agree" ); 00521 00522 // Read the gid map and set global ids 00523 std::vector< int > cell_gids( num_cells ); 00524 CCMIOReadCells( &error, cellsID, &mapID, NULL, CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) ); 00525 CHK_SET_CCMERR( error, "Couldn't read cells" ); 00526 CCMIOReadMap( &error, mapID, &cell_gids[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) ); 00527 CHK_SET_CCMERR( error, "Couldn't read cell id map" ); 00528 00529 ErrorCode rval = mbImpl->tag_set_data( mGlobalIdTag, &cells[0], cells.size(), &cell_gids[0] );MB_CHK_SET_ERR( rval, "Couldn't set gids tag" ); 00530 00531 // Now read cell material types; reuse cell_gids 00532 CCMIOReadCells( &error, cellsID, NULL, &cell_gids[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) ); 00533 CHK_SET_CCMERR( error, "Trouble reading cell types" ); 00534 00535 // Create the matsets 00536 std::map< int, Range > matset_ents; 00537 for( int i = 0; i < num_cells; i++ ) 00538 matset_ents[cell_gids[i]].insert( cells[i] ); 00539 00540 for( std::map< int, Range >::iterator mit = matset_ents.begin(); mit != matset_ents.end(); ++mit ) 00541 { 00542 EntityHandle matset; 00543 rval = mbImpl->create_meshset( MESHSET_SET, matset );MB_CHK_SET_ERR( rval, "Couldn't create material set" ); 00544 newMatsets[mit->first] = matset; 00545 00546 rval = mbImpl->add_entities( matset, mit->second );MB_CHK_SET_ERR( rval, "Couldn't add entities to material set" ); 00547 } 00548 00549 return MB_SUCCESS; 00550 } 00551 00552 ErrorCode ReadCCMIO::construct_cells( TupleList& face_map, 00553 #ifndef TUPLE_LIST 00554 SenseList& sense_map, 00555 #endif 00556 TupleList& /* vert_map */, 00557 std::map< int, int >& cell_topo_types, 00558 std::vector< EntityHandle >& new_cells ) 00559 { 00560 std::vector< EntityHandle > facehs; 00561 std::vector< int > senses; 00562 EntityHandle cell; 00563 ErrorCode tmp_rval, rval = MB_SUCCESS; 00564 EntityType this_type = MBMAXTYPE; 00565 bool has_mid_nodes = false; 00566 #ifdef TUPLE_LIST 00567 unsigned int i = 0; 00568 while( i < face_map.n ) 00569 { 00570 // Pull out face handles bounding the same cell 00571 facehs.clear(); 00572 int this_id = face_map.get_int( i ); 00573 unsigned int inext = i; 00574 while( face_map.get_int( inext ) == this_id && inext <= face_map.n ) 00575 { 00576 inext++; 00577 EntityHandle face = face_map.get_ulong( inext ); 00578 facehs.push_back( face ); 00579 senses.push_back( face_map.get_short( inext ) ); 00580 } 00581 this_type = MBMAXTYPE; 00582 has_mid_nodes = false; 00583 #else 00584 std::map< int, std::vector< EntityHandle > >::iterator fmit; 00585 std::map< int, std::vector< int > >::iterator smit; 00586 std::map< int, int >::iterator typeit; 00587 for( fmit = face_map.begin(), smit = sense_map.begin(); fmit != face_map.end(); ++fmit, ++smit ) 00588 { 00589 // Pull out face handles bounding the same cell 00590 facehs.clear(); 00591 int this_id = ( *fmit ).first; 00592 facehs = ( *fmit ).second; 00593 senses.clear(); 00594 senses = ( *smit ).second; 00595 typeit = cell_topo_types.find( this_id ); 00596 if( typeit != cell_topo_types.end() ) 00597 { 00598 rval = ccmio_to_moab_type( typeit->second, this_type, has_mid_nodes ); 00599 } 00600 else 00601 { 00602 this_type = MBMAXTYPE; 00603 has_mid_nodes = false; 00604 } 00605 #endif 00606 tmp_rval = create_cell_from_faces( facehs, senses, this_type, has_mid_nodes, cell ); 00607 if( MB_SUCCESS != tmp_rval ) 00608 rval = tmp_rval; 00609 else 00610 { 00611 new_cells.push_back( cell ); 00612 // Tag cell with global id 00613 tmp_rval = mbImpl->tag_set_data( mGlobalIdTag, &cell, 1, &this_id ); 00614 if( MB_SUCCESS != tmp_rval ) rval = tmp_rval; 00615 } 00616 } 00617 00618 return rval; 00619 } 00620 00621 ErrorCode ReadCCMIO::ccmio_to_moab_type( int ccm_type, EntityType& moab_type, bool& has_mid_nodes ) 00622 { 00623 switch( ccm_type ) 00624 { 00625 case 1: 00626 moab_type = MBVERTEX; 00627 break; 00628 case 2: 00629 case 28: 00630 moab_type = MBEDGE; 00631 break; 00632 case 29: 00633 moab_type = MBMAXTYPE; 00634 break; 00635 case 3: 00636 case 4: 00637 moab_type = MBQUAD; 00638 break; 00639 case 11: 00640 case 21: 00641 moab_type = MBHEX; 00642 break; 00643 case 12: 00644 case 22: 00645 moab_type = MBPRISM; 00646 break; 00647 case 13: 00648 case 23: 00649 moab_type = MBTET; 00650 break; 00651 case 14: 00652 case 24: 00653 moab_type = MBPYRAMID; 00654 break; 00655 case 255: 00656 moab_type = MBPOLYHEDRON; 00657 break; 00658 default: 00659 moab_type = MBMAXTYPE; 00660 } 00661 00662 switch( ccm_type ) 00663 { 00664 case 28: 00665 case 4: 00666 case 21: 00667 case 22: 00668 case 23: 00669 case 24: 00670 has_mid_nodes = true; 00671 break; 00672 default: 00673 break; 00674 } 00675 00676 return MB_SUCCESS; 00677 } 00678 00679 ErrorCode ReadCCMIO::create_cell_from_faces( std::vector< EntityHandle >& facehs, 00680 std::vector< int >& senses, 00681 EntityType this_type, 00682 bool /* has_mid_nodes */, 00683 EntityHandle& cell ) 00684 { 00685 ErrorCode rval; 00686 00687 // Test up front to see if they're one type 00688 EntityType face_type = mbImpl->type_from_handle( facehs[0] ); 00689 bool same_type = true; 00690 for( std::vector< EntityHandle >::iterator vit = facehs.begin(); vit != facehs.end(); ++vit ) 00691 { 00692 if( face_type != mbImpl->type_from_handle( *vit ) ) 00693 { 00694 same_type = false; 00695 break; 00696 } 00697 } 00698 00699 std::vector< EntityHandle > verts; 00700 EntityType input_type = this_type; 00701 std::vector< EntityHandle > storage; 00702 MeshTopoUtil mtu( mbImpl ); 00703 00704 // Preset this to maxtype, so we get an affirmative choice in loop 00705 this_type = MBMAXTYPE; 00706 00707 if( ( MBTET == input_type || MBMAXTYPE == input_type ) && same_type && face_type == MBTRI && facehs.size() == 4 ) 00708 { 00709 // Try to get proper connectivity for tet 00710 00711 // Get connectivity of first face, and reverse it if sense is forward, since 00712 // base face always points into entity 00713 rval = mbImpl->get_connectivity( &facehs[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" ); 00714 if( senses[0] > 0 ) std::reverse( verts.begin(), verts.end() ); 00715 00716 // Get the 4th vertex through the next tri 00717 const EntityHandle* conn; 00718 int conn_size; 00719 rval = mbImpl->get_connectivity( facehs[1], conn, conn_size, true, &storage );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" ); 00720 int i = 0; 00721 while( std::find( verts.begin(), verts.end(), conn[i] ) != verts.end() && i < conn_size ) 00722 i++; 00723 00724 // If i is not at the end of the verts, found the apex; otherwise fall back to polyhedron 00725 if( conn_size != i ) 00726 { 00727 this_type = MBTET; 00728 verts.push_back( conn[i] ); 00729 } 00730 } 00731 else if( ( MBHEX == input_type || MBMAXTYPE == input_type ) && same_type && MBQUAD == face_type && 00732 facehs.size() == 6 ) 00733 { 00734 // Build hex from quads 00735 // Algorithm: 00736 // - verts = vertices from 1st quad 00737 // - Find quad q1 sharing verts[0] and verts[1] 00738 // - Find quad q2 sharing other 2 verts in q1 00739 // - Find v1 = opposite vert from verts[1] in q1 , v2 = opposite from verts[0] 00740 // - Get i = offset of v1 in verts2 of q2, rotate verts2 by i 00741 // - If verts2[(i + 1) % 4] != v2, flip verts2 by switching verts2[1] and verts2[3] 00742 // - append verts2 to verts 00743 00744 // Get the other vertices for this hex; need to find the quad with no common vertices 00745 Range tmp_faces, tmp_verts; 00746 // Get connectivity of first face, and reverse it if sense is forward, since 00747 // base face always points into entity 00748 rval = mbImpl->get_connectivity( &facehs[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" ); 00749 if( senses[0] > 0 ) std::reverse( verts.begin(), verts.end() ); 00750 00751 // Get q1, which shares 2 vertices with q0 00752 std::copy( facehs.begin(), facehs.end(), range_inserter( tmp_faces ) ); 00753 rval = mbImpl->get_adjacencies( &verts[0], 2, 2, false, tmp_faces ); 00754 if( MB_SUCCESS != rval || tmp_faces.size() != 2 ) MB_SET_ERR( MB_FAILURE, "Couldn't get adj face" ); 00755 tmp_faces.erase( facehs[0] ); 00756 EntityHandle q1 = *tmp_faces.begin(); 00757 // Get other 2 verts of q1 00758 rval = mbImpl->get_connectivity( &q1, 1, tmp_verts );MB_CHK_SET_ERR( rval, "Couldn't get adj verts" ); 00759 tmp_verts.erase( verts[0] ); 00760 tmp_verts.erase( verts[1] ); 00761 // Get q2 00762 std::copy( facehs.begin(), facehs.end(), range_inserter( tmp_faces ) ); 00763 rval = mbImpl->get_adjacencies( tmp_verts, 2, false, tmp_faces ); 00764 if( MB_SUCCESS != rval || tmp_faces.size() != 2 ) MB_SET_ERR( MB_FAILURE, "Couldn't get adj face" ); 00765 tmp_faces.erase( q1 ); 00766 EntityHandle q2 = *tmp_faces.begin(); 00767 // Get verts in q2 00768 rval = mbImpl->get_connectivity( &q2, 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get adj vertices" ); 00769 00770 // Get verts in q1 opposite from v[1] and v[0] in q0 00771 EntityHandle v0 = 0, v1 = 0; 00772 rval = mtu.opposite_entity( q1, verts[1], v0 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" ); 00773 rval = mtu.opposite_entity( q1, verts[0], v1 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" ); 00774 if( v0 && v1 ) 00775 { 00776 // Offset of v0 in q2, then rotate and flip 00777 unsigned int ioff = std::find( storage.begin(), storage.end(), v0 ) - storage.begin(); 00778 if( 4 == ioff ) MB_SET_ERR( MB_FAILURE, "Trouble finding offset" ); 00779 00780 if( storage[( ioff + 1 ) % 4] != v1 ) 00781 { 00782 std::reverse( storage.begin(), storage.end() ); 00783 ioff = std::find( storage.begin(), storage.end(), v0 ) - storage.begin(); 00784 } 00785 if( 0 != ioff ) std::rotate( storage.begin(), storage.begin() + ioff, storage.end() ); 00786 00787 // Copy into verts, and make hex 00788 std::copy( storage.begin(), storage.end(), std::back_inserter( verts ) ); 00789 this_type = MBHEX; 00790 } 00791 } 00792 00793 if( MBMAXTYPE == this_type && facehs.size() == 5 ) 00794 { 00795 // Some preliminaries 00796 std::vector< EntityHandle > tris, quads; 00797 for( unsigned int i = 0; i < 5; i++ ) 00798 { 00799 if( MBTRI == mbImpl->type_from_handle( facehs[i] ) ) 00800 tris.push_back( facehs[i] ); 00801 else if( MBQUAD == mbImpl->type_from_handle( facehs[i] ) ) 00802 quads.push_back( facehs[i] ); 00803 } 00804 00805 // Check for prisms 00806 if( 2 == tris.size() && 3 == quads.size() ) 00807 { 00808 // OK, we have the right number of tris and quads; try to find the proper verts 00809 00810 // Get connectivity of first tri, and reverse if necessary 00811 int index = std::find( facehs.begin(), facehs.end(), tris[0] ) - facehs.begin(); 00812 rval = mbImpl->get_connectivity( &tris[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" ); 00813 if( senses[index] > 0 ) std::reverse( verts.begin(), verts.end() ); 00814 00815 // Now align vertices of other tri, through a quad, similar to how we did hexes 00816 // Get q1, which shares 2 vertices with t0 00817 Range tmp_faces, tmp_verts; 00818 std::copy( facehs.begin(), facehs.end(), range_inserter( tmp_faces ) ); 00819 rval = mbImpl->get_adjacencies( &verts[0], 2, 2, false, tmp_faces ); 00820 if( MB_SUCCESS != rval || tmp_faces.size() != 2 ) MB_SET_ERR( MB_FAILURE, "Couldn't get adj face" ); 00821 tmp_faces.erase( tris[0] ); 00822 EntityHandle q1 = *tmp_faces.begin(); 00823 // Get verts in q1 00824 rval = mbImpl->get_connectivity( &q1, 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get adj vertices" ); 00825 00826 // Get verts in q1 opposite from v[1] and v[0] in q0 00827 EntityHandle v0 = 0, v1 = 0; 00828 rval = mtu.opposite_entity( q1, verts[1], v0 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" ); 00829 rval = mtu.opposite_entity( q1, verts[0], v1 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" ); 00830 if( v0 && v1 ) 00831 { 00832 // Offset of v0 in t2, then rotate and flip 00833 storage.clear(); 00834 rval = mbImpl->get_connectivity( &tris[1], 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" ); 00835 00836 index = std::find( facehs.begin(), facehs.end(), tris[1] ) - facehs.begin(); 00837 if( senses[index] < 0 ) std::reverse( storage.begin(), storage.end() ); 00838 unsigned int ioff = std::find( storage.begin(), storage.end(), v0 ) - storage.begin(); 00839 if( 3 == ioff ) MB_SET_ERR( MB_FAILURE, "Trouble finding offset" ); 00840 for( unsigned int i = 0; i < 3; i++ ) 00841 verts.push_back( storage[( ioff + i ) % 3] ); 00842 00843 this_type = MBPRISM; 00844 } 00845 } 00846 else if( tris.size() == 4 && quads.size() == 1 ) 00847 { 00848 // Check for pyramid 00849 // Get connectivity of first tri, and reverse if necessary 00850 int index = std::find( facehs.begin(), facehs.end(), quads[0] ) - facehs.begin(); 00851 rval = mbImpl->get_connectivity( &quads[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" ); 00852 if( senses[index] > 0 ) std::reverse( verts.begin(), verts.end() ); 00853 00854 // Get apex node 00855 rval = mbImpl->get_connectivity( &tris[0], 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" ); 00856 for( unsigned int i = 0; i < 3; i++ ) 00857 { 00858 if( std::find( verts.begin(), verts.end(), storage[i] ) == verts.end() ) 00859 { 00860 verts.push_back( storage[i] ); 00861 break; 00862 } 00863 } 00864 00865 if( 5 == verts.size() ) this_type = MBPYRAMID; 00866 } 00867 else 00868 { 00869 // Dummy else clause to stop in debugger 00870 this_type = MBMAXTYPE; 00871 } 00872 } 00873 00874 if( MBMAXTYPE != input_type && input_type != this_type && this_type != MBMAXTYPE ) 00875 std::cerr << "Warning: types disagree (cell_topo_type = " << CN::EntityTypeName( input_type ) 00876 << ", faces indicate type " << CN::EntityTypeName( this_type ) << std::endl; 00877 00878 if( MBMAXTYPE != input_type && this_type == MBMAXTYPE && input_type != MBPOLYHEDRON ) 00879 std::cerr << "Warning: couldn't find proper connectivity for specified topo_type = " 00880 << CN::EntityTypeName( input_type ) << std::endl; 00881 00882 // Now make the element; if we fell back to polyhedron, use faces, otherwise use verts 00883 if( MBPOLYHEDRON == this_type || MBMAXTYPE == this_type ) 00884 { 00885 rval = mbImpl->create_element( MBPOLYHEDRON, &facehs[0], facehs.size(), cell );MB_CHK_SET_ERR( rval, "create_element failed" ); 00886 } 00887 else 00888 { 00889 rval = mbImpl->create_element( this_type, &verts[0], verts.size(), cell );MB_CHK_SET_ERR( rval, "create_element failed" ); 00890 } 00891 00892 return MB_SUCCESS; 00893 } 00894 00895 ErrorCode ReadCCMIO::read_all_faces( CCMIOID topologyID, 00896 TupleList& vert_map, 00897 TupleList& face_map, 00898 #ifndef TUPLE_LIST 00899 SenseList& sense_map, 00900 #endif 00901 Range* new_faces ) 00902 { 00903 CCMIOSize_t index; 00904 CCMIOID faceID; 00905 ErrorCode rval; 00906 CCMIOError error = kCCMIONoErr; 00907 00908 // Get total # internal/bdy faces, size the face map accordingly 00909 #ifdef TUPLE_LIST 00910 index = CCMIOSIZEC( 0 ); 00911 int nbdy_faces = 0; 00912 CCMIOSize_t nf; 00913 error = kCCMIONoErr; 00914 while( kCCMIONoErr == CCMIONextEntity( NULL, topologyID, kCCMIOBoundaryFaces, &index, &faceID ) ) 00915 { 00916 CCMIOEntitySize( &error, faceID, &nf, NULL ); 00917 nbdy_faces += nf; 00918 } 00919 CCMIOGetEntity( &error, topologyID, kCCMIOInternalFaces, 0, &faceID ); 00920 CCMIOEntitySize( &error, faceID, &nf, NULL ); 00921 00922 int nint_faces = nf; 00923 face_map.resize( 2 * nint_faces + nbdy_faces ); 00924 #endif 00925 00926 // Get multiple blocks of bdy faces 00927 index = CCMIOSIZEC( 0 ); 00928 while( kCCMIONoErr == CCMIONextEntity( NULL, topologyID, kCCMIOBoundaryFaces, &index, &faceID ) ) 00929 { 00930 rval = read_faces( faceID, kCCMIOBoundaryFaces, vert_map, face_map, 00931 #ifndef TUPLE_LIST 00932 sense_map, 00933 #endif 00934 new_faces );MB_CHK_SET_ERR( rval, "Trouble reading boundary faces" ); 00935 } 00936 00937 // Now get internal faces 00938 CCMIOGetEntity( &error, topologyID, kCCMIOInternalFaces, 0, &faceID ); 00939 CHK_SET_CCMERR( error, "Couldn't get internal faces" ); 00940 00941 rval = read_faces( faceID, kCCMIOInternalFaces, vert_map, face_map, 00942 #ifndef TUPLE_LIST 00943 sense_map, 00944 #endif 00945 new_faces );MB_CHK_SET_ERR( rval, "Trouble reading internal faces" ); 00946 00947 return rval; 00948 } 00949 00950 ErrorCode ReadCCMIO::read_faces( CCMIOID faceID, 00951 CCMIOEntity bdy_or_int, 00952 TupleList& vert_map, 00953 TupleList& face_map, 00954 #ifndef TUPLE_LIST 00955 SenseList& sense_map, 00956 #endif 00957 Range* new_faces ) 00958 { 00959 if( kCCMIOInternalFaces != bdy_or_int && kCCMIOBoundaryFaces != bdy_or_int ) 00960 MB_SET_ERR( MB_FAILURE, "Face type isn't boundary or internal" ); 00961 00962 CCMIOSize_t dum_faces; 00963 CCMIOError error = kCCMIONoErr; 00964 CCMIOEntitySize( &error, faceID, &dum_faces, NULL ); 00965 int num_faces = GETINT32( dum_faces ); 00966 00967 // Get the size of the face connectivity array (not really a straight connect 00968 // array, has n, connect(n), ...) 00969 CCMIOSize_t farray_size = CCMIOSIZEC( 0 ); 00970 CCMIOReadFaces( &error, faceID, bdy_or_int, NULL, &farray_size, NULL, CCMIOINDEXC( kCCMIOStart ), 00971 CCMIOINDEXC( kCCMIOEnd ) ); 00972 CHK_SET_CCMERR( error, "Trouble reading face connectivity length" ); 00973 00974 // Allocate vectors for holding farray and cells for each face; use new for finer 00975 // control of de-allocation 00976 int num_sides = ( kCCMIOInternalFaces == bdy_or_int ? 2 : 1 ); 00977 int* farray = new int[GETINT32( farray_size )]; 00978 00979 // Read farray and make the faces 00980 CCMIOID mapID; 00981 CCMIOReadFaces( &error, faceID, bdy_or_int, &mapID, NULL, farray, CCMIOINDEXC( kCCMIOStart ), 00982 CCMIOINDEXC( kCCMIOEnd ) ); 00983 CHK_SET_CCMERR( error, "Trouble reading face connectivity" ); 00984 00985 std::vector< EntityHandle > face_handles; 00986 ErrorCode rval = make_faces( farray, vert_map, face_handles, num_faces );MB_CHK_SET_ERR( rval, "Failed to make the faces" ); 00987 00988 // Read face cells and make tuples 00989 int* face_cells; 00990 if( num_sides * num_faces < farray_size ) 00991 face_cells = new int[num_sides * num_faces]; 00992 else 00993 face_cells = farray; 00994 CCMIOReadFaceCells( &error, faceID, bdy_or_int, face_cells, CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) ); 00995 CHK_SET_CCMERR( error, "Trouble reading face cells" ); 00996 00997 int* tmp_ptr = face_cells; 00998 for( unsigned int i = 0; i < face_handles.size(); i++ ) 00999 { 01000 #ifdef TUPLE_LIST 01001 short forward = 1, reverse = -1; 01002 face_map.push_back( &forward, tmp_ptr++, &face_handles[i], NULL ); 01003 if( 2 == num_sides ) face_map.push_back( &reverse, tmp_ptr++, &face_handles[i], NULL ); 01004 #else 01005 face_map[*tmp_ptr].push_back( face_handles[i] ); 01006 sense_map[*tmp_ptr++].push_back( 1 ); 01007 if( 2 == num_sides ) 01008 { 01009 face_map[*tmp_ptr].push_back( face_handles[i] ); 01010 sense_map[*tmp_ptr++].push_back( -1 ); 01011 } 01012 #endif 01013 } 01014 01015 // Now read & set face gids, reuse face_cells 'cuz we know it's big enough 01016 CCMIOReadMap( &error, mapID, face_cells, CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) ); 01017 CHK_SET_CCMERR( error, "Trouble reading face gids" ); 01018 01019 rval = mbImpl->tag_set_data( mGlobalIdTag, &face_handles[0], face_handles.size(), face_cells );MB_CHK_SET_ERR( rval, "Couldn't set face global ids" ); 01020 01021 // Make a neumann set for these faces if they're all in a boundary face set 01022 if( kCCMIOBoundaryFaces == bdy_or_int ) 01023 { 01024 EntityHandle neuset; 01025 rval = mbImpl->create_meshset( MESHSET_SET, neuset );MB_CHK_SET_ERR( rval, "Failed to create neumann set" ); 01026 01027 // Don't trust entity index passed in 01028 int index; 01029 CCMIOGetEntityIndex( &error, faceID, &index ); 01030 newNeusets[index] = neuset; 01031 01032 rval = mbImpl->add_entities( neuset, &face_handles[0], face_handles.size() );MB_CHK_SET_ERR( rval, "Failed to add faces to neumann set" ); 01033 01034 // Now tag as neumann set; will add id later 01035 int dum_val = 0; 01036 rval = mbImpl->tag_set_data( mNeumannSetTag, &neuset, 1, &dum_val );MB_CHK_SET_ERR( rval, "Failed to tag neumann set" ); 01037 } 01038 01039 if( new_faces ) 01040 { 01041 std::sort( face_handles.begin(), face_handles.end() ); 01042 std::copy( face_handles.rbegin(), face_handles.rend(), range_inserter( *new_faces ) ); 01043 } 01044 01045 return MB_SUCCESS; 01046 } 01047 01048 ErrorCode ReadCCMIO::make_faces( int* farray, 01049 TupleList& vert_map, 01050 std::vector< EntityHandle >& new_faces, 01051 int num_faces ) 01052 { 01053 std::vector< EntityHandle > verts; 01054 ErrorCode tmp_rval = MB_SUCCESS, rval = MB_SUCCESS; 01055 01056 for( int i = 0; i < num_faces; i++ ) 01057 { 01058 int num_verts = *farray++; 01059 verts.resize( num_verts ); 01060 01061 // Fill in connectivity by looking up by gid in vert tuple_list 01062 for( int j = 0; j < num_verts; j++ ) 01063 { 01064 #ifdef TUPLE_LIST 01065 int tindex = vert_map.find( 1, farray[j] ); 01066 if( -1 == tindex ) 01067 { 01068 tmp_rval = MB_FAILURE; 01069 break; 01070 } 01071 verts[j] = vert_map.get_ulong( tindex, 0 ); 01072 #else 01073 verts[j] = ( vert_map[farray[j]] )[0]; 01074 #endif 01075 } 01076 farray += num_verts; 01077 01078 if( MB_SUCCESS == tmp_rval ) 01079 { 01080 // Make face 01081 EntityType ftype = ( 3 == num_verts ? MBTRI : ( 4 == num_verts ? MBQUAD : MBPOLYGON ) ); 01082 EntityHandle faceh; 01083 tmp_rval = mbImpl->create_element( ftype, &verts[0], num_verts, faceh ); 01084 if( faceh ) new_faces.push_back( faceh ); 01085 } 01086 01087 if( MB_SUCCESS != tmp_rval ) rval = tmp_rval; 01088 } 01089 01090 return rval; 01091 } 01092 01093 ErrorCode ReadCCMIO::read_vertices( CCMIOSize_t /* proc */, 01094 CCMIOID /* processorID */, 01095 CCMIOID verticesID, 01096 CCMIOID /* topologyID */, 01097 Range* verts, 01098 TupleList& vert_map ) 01099 { 01100 CCMIOError error = kCCMIONoErr; 01101 01102 // Pre-read the number of vertices, so we can pre-allocate & read directly in 01103 CCMIOSize_t nverts = CCMIOSIZEC( 0 ); 01104 CCMIOEntitySize( &error, verticesID, &nverts, NULL ); 01105 CHK_SET_CCMERR( error, "Couldn't get number of vertices" ); 01106 01107 // Get # dimensions 01108 CCMIOSize_t dims; 01109 float scale; 01110 CCMIOReadVerticesf( &error, verticesID, &dims, NULL, NULL, NULL, CCMIOINDEXC( 0 ), CCMIOINDEXC( 1 ) ); 01111 CHK_SET_CCMERR( error, "Couldn't get number of dimensions" ); 01112 01113 // Allocate vertex space 01114 EntityHandle node_handle = 0; 01115 std::vector< double* > arrays; 01116 readMeshIface->get_node_coords( 3, GETINT32( nverts ), MB_START_ID, node_handle, arrays ); 01117 01118 // Read vertex coords 01119 CCMIOID mapID; 01120 std::vector< double > tmp_coords( GETINT32( dims ) * GETINT32( nverts ) ); 01121 CCMIOReadVerticesd( &error, verticesID, &dims, &scale, &mapID, &tmp_coords[0], CCMIOINDEXC( 0 ), 01122 CCMIOINDEXC( 0 + nverts ) ); 01123 CHK_SET_CCMERR( error, "Trouble reading vertex coordinates" ); 01124 01125 // Copy interleaved coords into moab blocked coordinate space 01126 int i = 0, threei = 0; 01127 for( ; i < nverts; i++ ) 01128 { 01129 arrays[0][i] = tmp_coords[threei++]; 01130 arrays[1][i] = tmp_coords[threei++]; 01131 if( 3 == GETINT32( dims ) ) 01132 arrays[2][i] = tmp_coords[threei++]; 01133 else 01134 arrays[2][i] = 0.0; 01135 } 01136 01137 // Scale, if necessary 01138 if( 1.0 != scale ) 01139 { 01140 for( i = 0; i < nverts; i++ ) 01141 { 01142 arrays[0][i] *= scale; 01143 arrays[1][i] *= scale; 01144 if( 3 == GETINT32( dims ) ) arrays[2][i] *= scale; 01145 } 01146 } 01147 01148 // Read gids for vertices 01149 std::vector< int > gids( GETINT32( nverts ) ); 01150 CCMIOReadMap( &error, mapID, &gids[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) ); 01151 CHK_SET_CCMERR( error, "Trouble reading vertex global ids" ); 01152 01153 // Put new vertex handles into range, and set gids for them 01154 Range new_verts( node_handle, node_handle + nverts - 1 ); 01155 ErrorCode rval = mbImpl->tag_set_data( mGlobalIdTag, new_verts, &gids[0] );MB_CHK_SET_ERR( rval, "Couldn't set gids on vertices" ); 01156 01157 // Pack vert_map with global ids and handles for these vertices 01158 #ifdef TUPLE_LIST 01159 vert_map.resize( GETINT32( nverts ) ); 01160 for( i = 0; i < GETINT32( nverts ); i++ ) 01161 { 01162 vert_map.push_back( NULL, &gids[i], &node_handle, NULL ); 01163 #else 01164 for( i = 0; i < GETINT32( nverts ); i++ ) 01165 { 01166 ( vert_map[gids[i]] ).push_back( node_handle ); 01167 #endif 01168 node_handle += 1; 01169 } 01170 01171 if( verts ) verts->merge( new_verts ); 01172 01173 return MB_SUCCESS; 01174 } 01175 01176 ErrorCode ReadCCMIO::get_processors( CCMIOID stateID, 01177 CCMIOID& processorID, 01178 CCMIOID& verticesID, 01179 CCMIOID& topologyID, 01180 CCMIOID& solutionID, 01181 std::vector< CCMIOSize_t >& procs, 01182 bool& /* has_solution */ ) 01183 { 01184 CCMIOSize_t proc = CCMIOSIZEC( 0 ); 01185 CCMIOError error = kCCMIONoErr; 01186 01187 CCMIONextEntity( &error, stateID, kCCMIOProcessor, &proc, &processorID ); 01188 CHK_SET_CCMERR( error, "CCMIONextEntity() failed" ); 01189 if( CCMIOReadProcessor( NULL, processorID, &verticesID, &topologyID, NULL, &solutionID ) != kCCMIONoErr ) 01190 { 01191 // Maybe no solution; try again 01192 CCMIOReadProcessor( &error, processorID, &verticesID, &topologyID, NULL, NULL ); 01193 CHK_SET_CCMERR( error, "CCMIOReadProcessor() failed" ); 01194 hasSolution = false; 01195 } 01196 01197 procs.push_back( proc ); 01198 01199 return MB_SUCCESS; 01200 } 01201 01202 ErrorCode ReadCCMIO::read_tag_values( const char* /* file_name */, 01203 const char* /* tag_name */, 01204 const FileOptions& /* opts */, 01205 std::vector< int >& /* tag_values_out */, 01206 const SubsetList* /* subset_list */ ) 01207 { 01208 return MB_FAILURE; 01209 } 01210 01211 } // namespace moab