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