MOAB: Mesh Oriented datABase  (version 5.3.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     {
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
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines