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