moab::ReadNCDF Class Reference

Output Exodus File for VERDE. More...

#include <ReadNCDF.hpp>

Inheritance diagram for moab::ReadNCDF:

Collaboration diagram for moab::ReadNCDF:

Public Member Functions
void	tokenize (const std::string &str, std::vector< std::string > &tokens, const char *delimiters)
ErrorCode	load_file (const char *file_name, const EntityHandle *file_set, const FileOptions &opts, const SubsetList *subset_list=0, const Tag *file_id_tag=0)
	load an ExoII file More...
ErrorCode	read_tag_values (const char *file_name, const char *tag_name, const FileOptions &opts, std::vector< int > &tag_values_out, const SubsetList *subset_list=0)
	Read tag values from a file. More...
	ReadNCDF (Interface *impl=NULL)
	Constructor. More...
virtual	~ReadNCDF ()
	Destructor. More...
ErrorCode	update (const char *exodus_file_name, const FileOptions &opts, const int num_blocks, const int *blocks_to_load, const EntityHandle file_set)
Public Member Functions inherited from moab::ReaderIface
virtual	~ReaderIface ()

Static Public Member Functions
static ReaderIface *	factory (Interface *)

Private Member Functions
bool	dimension_exists (const char *attrib_name)
void	reset ()
ErrorCode	read_exodus_header ()
	read the header from the ExoII file More...
ErrorCode	read_nodes (const Tag *file_id_tag)
	read the nodes More...
ErrorCode	read_face_blocks_headers ()
ErrorCode	read_block_headers (const int *blocks_to_load, const int num_blocks)
	read block headers, containing info about element type, number, etc. More...
ErrorCode	read_polyhedra_faces ()
ErrorCode	read_elements (const Tag *file_id_tag)
	read the element blocks More...
ErrorCode	read_global_ids ()
	read in the global element ids More...
ErrorCode	read_nodesets ()
	read the nodesets into meshsets More...
ErrorCode	read_sidesets ()
	read the sidesets (does nothing for now) More...
int	exodus_file ()
	exodus file bound to this object More...
int	number_dimensions ()
	number of dimensions in this exo file More...
ErrorCode	get_type (char *exo_element_type, EntityType &elem_type)
	map a character exodusII element type to a TSTT type & topology More...
ErrorCode	get_type (EntityType &elem_type, std::string &exo_element_type)
ErrorCode	read_qa_records (EntityHandle file_set)
ErrorCode	read_qa_information (std::vector< std::string > &qa_record_list)
ErrorCode	read_qa_string (char *string, int record_number, int record_position)
ErrorCode	create_ss_elements (int *element_ids, int *side_list, int num_sides, int num_dist_factors, std::vector< EntityHandle > &entities_to_add, std::vector< EntityHandle > &reverse_entities, std::vector< double > &dist_factor_vector, int ss_seq_id)
ErrorCode	find_side_element_type (const int element_id, ExoIIElementType &type, ReadBlockData &block_data, int &df_index, int side_id)
ErrorCode	create_sideset_element (const std::vector< EntityHandle > &, EntityType, EntityHandle &, int &)
	creates an element with the given connectivity More...
int	get_number_nodes (EntityHandle handle)

Private Attributes
ReadUtilIface *	readMeshIface
Interface *	mdbImpl
	interface instance More...
int	ncFile
int	CPU_WORD_SIZE
int	IO_WORD_SIZE
int	vertexOffset
	int to offset vertex ids with More...
std::string	exodusFile
	file name More...
int	numberNodes_loading
	number of nodes in the current exoII file More...
int	numberElements_loading
	number of elements in the current exoII file More...
int	numberElementBlocks_loading
	number of blocks in the current exoII file More...
int	numberFaceBlocks_loading
	number of face blocks in the current exoII file (used for polyhedra) More...
int	numberNodeSets_loading
	number of nodesets in the current exoII file More...
int	numberSideSets_loading
	number of sidesets in the current exoII file More...
int	numberDimensions_loading
EntityHandle	mCurrentMeshHandle
	Meshset Handle for the mesh that is currently being read. More...
std::vector< char >	nodesInLoadedBlocks
std::vector< ReadBlockData >	blocksLoading
std::vector< EntityHandle >	polyfaces
std::vector< ReadFaceBlockData >	faceBlocksLoading
Tag	mMaterialSetTag
	Cached tags for reading. Note that all these tags are defined when the core is initialized. More...
Tag	mDirichletSetTag
Tag	mNeumannSetTag
Tag	mHasMidNodesTag
Tag	mDistFactorTag
Tag	mGlobalIdTag
Tag	mQaRecordTag
int	max_line_length
int	max_str_length
Range	initRange
	range of entities in initial mesh, before this read More...

Detailed Description

Output Exodus File for VERDE.

Definition at line 73 of file ReadNCDF.hpp.

Constructor & Destructor Documentation

ReadNCDF()

moab::ReadNCDF::ReadNCDF

(

Interface *

impl = NULL

)

Constructor.

Get and cache predefined tag handles

Definition at line 140 of file ReadNCDF.cpp.

                                    : mdbImpl( impl ), max_line_length( -1 ), max_str_length( -1 )
{
    assert( impl != NULL );
    reset();
 
    impl->query_interface( readMeshIface );
 
    // Initialize in case tag_get_handle fails below
    mMaterialSetTag  = 0;
    mDirichletSetTag = 0;
    mNeumannSetTag   = 0;
    mHasMidNodesTag  = 0;
    mDistFactorTag   = 0;
    mQaRecordTag     = 0;
    mGlobalIdTag     = 0;
 
    //! Get and cache predefined tag handles
    ErrorCode result;
    const int negone = -1;
    result           = impl->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mMaterialSetTag,
                                             MB_TAG_SPARSE | MB_TAG_CREAT, &negone );
    assert( MB_SUCCESS == result );
    result = impl->tag_get_handle( DIRICHLET_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mDirichletSetTag,
                                   MB_TAG_SPARSE | MB_TAG_CREAT, &negone );
    assert( MB_SUCCESS == result );
    result = impl->tag_get_handle( NEUMANN_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mNeumannSetTag,
                                   MB_TAG_SPARSE | MB_TAG_CREAT, &negone );
    assert( MB_SUCCESS == result );
    const int mids[] = { -1, -1, -1, -1 };
    result           = impl->tag_get_handle( HAS_MID_NODES_TAG_NAME, 4, MB_TYPE_INTEGER, mHasMidNodesTag,
                                             MB_TAG_SPARSE | MB_TAG_CREAT, mids );
    assert( MB_SUCCESS == result );
    result = impl->tag_get_handle( "distFactor", 0, MB_TYPE_DOUBLE, mDistFactorTag,
                                   MB_TAG_SPARSE | MB_TAG_VARLEN | MB_TAG_CREAT );
    assert( MB_SUCCESS == result );
    result = impl->tag_get_handle( "qaRecord", 0, MB_TYPE_OPAQUE, mQaRecordTag,
                                   MB_TAG_SPARSE | MB_TAG_VARLEN | MB_TAG_CREAT );
    assert( MB_SUCCESS == result );
    mGlobalIdTag = impl->globalId_tag();
 
    assert( MB_SUCCESS == result );
#ifdef NDEBUG
    if( MB_SUCCESS == result )
    {
    };  // Line to avoid compiler warning about unused variable
#endif
    ncFile = 0;
}

References DIRICHLET_SET_TAG_NAME, ErrorCode, moab::Interface::globalId_tag(), HAS_MID_NODES_TAG_NAME, MATERIAL_SET_TAG_NAME, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_SPARSE, MB_TAG_VARLEN, MB_TYPE_DOUBLE, MB_TYPE_INTEGER, MB_TYPE_OPAQUE, mDirichletSetTag, mDistFactorTag, mGlobalIdTag, mHasMidNodesTag, mMaterialSetTag, mNeumannSetTag, mQaRecordTag, ncFile, NEUMANN_SET_TAG_NAME, moab::Interface::query_interface(), readMeshIface, reset(), and moab::Interface::tag_get_handle().

Referenced by factory().

~ReadNCDF()

moab::ReadNCDF::~ReadNCDF ( )

virtual

Destructor.

Definition at line 209 of file ReadNCDF.cpp.

{
    mdbImpl->release_interface( readMeshIface );
}

References mdbImpl, readMeshIface, and moab::Interface::release_interface().

Member Function Documentation

create_sideset_element()

ErrorCode moab::ReadNCDF::create_sideset_element ( const std::vector< EntityHandle > &  connectivity, EntityType  type, EntityHandle &  handle, int &  sense  )

private

creates an element with the given connectivity

Definition at line 1597 of file ReadNCDF.cpp.

{
    // Get adjacent entities
    ErrorCode error = MB_SUCCESS;
    int to_dim      = CN::Dimension( type );
    // for matching purposes , consider only corner nodes
    // basically, assume everything is matching if the corner nodes are matching, and
    // the number of total nodes is the same
    int nb_corner_nodes = CN::VerticesPerEntity( type );
    std::vector< EntityHandle > adj_ent;
    mdbImpl->get_adjacencies( &( connectivity[0] ), 1, to_dim, false, adj_ent );
 
    // For each entity, see if we can find a match
    // If we find a match, return it
    bool match_found = false;
    std::vector< EntityHandle > match_conn;
    // By default, assume positive sense
    sense = 1;
    for( unsigned int i = 0; i < adj_ent.size() && match_found == false; i++ )
    {
        // Get the connectivity
        error = mdbImpl->get_connectivity( &( adj_ent[i] ), 1, match_conn );
        if( error != MB_SUCCESS ) continue;
 
        // Make sure they have the same number of vertices (higher order elements ?)
        if( match_conn.size() != connectivity.size() ) continue;
 
        // Find a matching node
        std::vector< EntityHandle >::iterator iter;
        std::vector< EntityHandle >::iterator end_corner = match_conn.begin() + nb_corner_nodes;
        iter                                             = std::find( match_conn.begin(), end_corner, connectivity[0] );
        if( iter == end_corner ) continue;
 
        // Rotate to match connectivity
        // rotate only corner nodes, ignore the rest from now on
        std::rotate( match_conn.begin(), iter, end_corner );
 
        bool they_match = true;
        for( int j = 1; j < nb_corner_nodes; j++ )
        {
            if( connectivity[j] != match_conn[j] )
            {
                they_match = false;
                break;
            }
        }
 
        // If we didn't get a match
        if( !they_match )
        {
            // Try the opposite sense
            they_match = true;
 
            int k = nb_corner_nodes - 1;
            for( int j = 1; j < nb_corner_nodes; )
            {
                if( connectivity[j] != match_conn[k] )
                {
                    they_match = false;
                    break;
                }
                ++j;
                --k;
            }
            // If they matched here, sense is reversed
            if( they_match ) sense = -1;
        }
        match_found = they_match;
        if( match_found ) handle = adj_ent[i];
    }
 
    // If we didn't find a match, create an element
    if( !match_found ) error = mdbImpl->create_element( type, &connectivity[0], connectivity.size(), handle );
 
    return error;
}

References moab::Interface::create_element(), moab::CN::Dimension(), moab::error(), ErrorCode, moab::Interface::get_adjacencies(), moab::Interface::get_connectivity(), MB_SUCCESS, mdbImpl, and moab::CN::VerticesPerEntity().

Referenced by create_ss_elements().

create_ss_elements()

ErrorCode moab::ReadNCDF::create_ss_elements ( int *  element_ids, int *  side_list, int  num_sides, int  num_dist_factors, std::vector< EntityHandle > &  entities_to_add, std::vector< EntityHandle > &  reverse_entities, std::vector< double > &  dist_factor_vector, int  ss_seq_id  )

private

Definition at line 1357 of file ReadNCDF.cpp.

{
    // Determine entity type from element id
 
    // If there are dist. factors, create a vector to hold the array
    // and place this array as a tag onto the sideset meshset
 
    std::vector< double > temp_dist_factor_vector( num_dist_factors );
    std::vector< char > temp_string_storage( max_str_length + 1 );
    char* temp_string = &temp_string_storage[0];
    int temp_var;
    if( num_dist_factors )
    {
        INS_ID( temp_string, "dist_fact_ss%d", ss_seq_id, max_str_length + 1 );
        temp_var = -1;
        GET_1D_DBL_VAR( temp_string, temp_var, temp_dist_factor_vector );
        if( -1 == temp_var )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting dist fact variable" );
        }
    }
 
    EntityType subtype;
    int num_side_nodes, num_elem_nodes;
    const EntityHandle* nodes;
    std::vector< EntityHandle > connectivity;
    int side_node_idx[32];
 
    int df_index = 0;
    int sense    = 0;
    for( int i = 0; i < num_sides; i++ )
    {
        ExoIIElementType exoii_type;
        ReadBlockData block_data;
        block_data.elemType = EXOII_MAX_ELEM_TYPE;
 
        if( find_side_element_type( element_ids[i], exoii_type, block_data, df_index, side_list[i] ) != MB_SUCCESS )
            continue;  // Isn't being read in this time
 
        EntityType type = ExoIIUtil::ExoIIElementMBEntity[exoii_type];
 
        EntityHandle ent_handle = element_ids[i] - block_data.startExoId + block_data.startMBId;
 
        const int side_num = side_list[i] - 1;
        if( type == MBHEX )
        {
            // Get the nodes of the element
            if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
 
            CN::SubEntityNodeIndices( type, num_elem_nodes, 2, side_num, subtype, num_side_nodes, side_node_idx );
            if( num_side_nodes <= 0 ) return MB_FAILURE;
 
            connectivity.resize( num_side_nodes );
            for( int k = 0; k < num_side_nodes; ++k )
                connectivity[k] = nodes[side_node_idx[k]];
 
            if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) ) return MB_FAILURE;
            if( 1 == sense )
                entities_to_add.push_back( ent_handle );
            else
                reverse_entities.push_back( ent_handle );
 
            // Read in distribution factor array
            if( num_dist_factors )
            {
                for( int k = 0; k < 4; k++ )
                    dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
            }
        }
 
        // If it is a Tet
        else if( type == MBTET )
        {
            // Get the nodes of the element
            if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
 
            CN::SubEntityNodeIndices( type, num_elem_nodes, 2, side_num, subtype, num_side_nodes, side_node_idx );
            if( num_side_nodes <= 0 ) return MB_FAILURE;
 
            connectivity.resize( num_side_nodes );
            for( int k = 0; k < num_side_nodes; ++k )
                connectivity[k] = nodes[side_node_idx[k]];
 
            if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) ) return MB_FAILURE;
            if( 1 == sense )
                entities_to_add.push_back( ent_handle );
            else
                reverse_entities.push_back( ent_handle );
 
            // Read in distribution factor array
            if( num_dist_factors )
            {
                for( int k = 0; k < 3; k++ )
                    dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
            }
        }
        else if( type == MBQUAD && exoii_type >= EXOII_SHELL && exoii_type <= EXOII_SHELL9 )
        {
            // ent_handle = CREATE_HANDLE(MBQUAD, base_id, error);
 
            // Just use this quad
            if( side_list[i] == 1 )
            {
                if( 1 == sense )
                    entities_to_add.push_back( ent_handle );
                else
                    reverse_entities.push_back( ent_handle );
 
                if( num_dist_factors )
                {
                    for( int k = 0; k < 4; k++ )
                        dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
                }
 
                continue;
            }
            else if( side_list[i] == 2 )
            {
                reverse_entities.push_back( ent_handle );
 
                if( num_dist_factors )
                {
                    for( int k = 0; k < 4; k++ )
                        dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
                }
 
                continue;
            }
            else
            {
                // Get the nodes of the element
                if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
 
                CN::SubEntityNodeIndices( type, num_elem_nodes, 1, side_num - 2, subtype, num_side_nodes,
                                          side_node_idx );
                if( num_side_nodes <= 0 ) return MB_FAILURE;
 
                connectivity.resize( num_side_nodes );
                for( int k = 0; k < num_side_nodes; ++k )
                    connectivity[k] = nodes[side_node_idx[k]];
 
                if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) )
                    return MB_FAILURE;
                if( 1 == sense )
                    entities_to_add.push_back( ent_handle );
                else
                    reverse_entities.push_back( ent_handle );
 
                if( num_dist_factors )
                {
                    for( int k = 0; k < 2; k++ )
                        dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
                }
            }
        }
        // If it is a Quad
        else if( type == MBQUAD )
        {
            // Get the nodes of the element
            if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
 
            CN::SubEntityNodeIndices( type, num_elem_nodes, 1, side_num, subtype, num_side_nodes, side_node_idx );
            if( num_side_nodes <= 0 ) return MB_FAILURE;
 
            connectivity.resize( num_side_nodes );
            for( int k = 0; k < num_side_nodes; ++k )
                connectivity[k] = nodes[side_node_idx[k]];
 
            if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) ) return MB_FAILURE;
            if( 1 == sense )
                entities_to_add.push_back( ent_handle );
            else
                reverse_entities.push_back( ent_handle );
 
            // Read in distribution factor array
            if( num_dist_factors )
            {
                for( int k = 0; k < 2; k++ )
                    dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
            }
        }
        else if( type == MBTRI )
        {
            int side_offset = 0;
            if( number_dimensions() == 3 && side_list[i] <= 2 )
            {
                if( 1 == sense )
                    entities_to_add.push_back( ent_handle );
                else
                    reverse_entities.push_back( ent_handle );
                if( num_dist_factors )
                {
                    for( int k = 0; k < 3; k++ )
                        dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
                }
            }
            else
            {
                if( number_dimensions() == 3 )
                {
                    if( side_list[i] > 2 ) side_offset = 2;
                }
 
                // Get the nodes of the element
                if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
 
                CN::SubEntityNodeIndices( type, num_elem_nodes, 1, side_num - side_offset, subtype, num_side_nodes,
                                          side_node_idx );
                if( num_side_nodes <= 0 ) return MB_FAILURE;
 
                connectivity.resize( num_side_nodes );
                for( int k = 0; k < num_side_nodes; ++k )
                    connectivity[k] = nodes[side_node_idx[k]];
 
                if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) )
                    return MB_FAILURE;
                if( 1 == sense )
                    entities_to_add.push_back( ent_handle );
                else
                    reverse_entities.push_back( ent_handle );
 
                if( num_dist_factors )
                {
                    for( int k = 0; k < 2; k++ )
                        dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
                }
            }
        }
    }
 
    return MB_SUCCESS;
}

References create_sideset_element(), moab::ReadBlockData::elemType, moab::EXOII_MAX_ELEM_TYPE, moab::EXOII_SHELL, moab::EXOII_SHELL9, moab::ExoIIUtil::ExoIIElementMBEntity, find_side_element_type(), GET_1D_DBL_VAR, moab::Interface::get_connectivity(), INS_ID, max_str_length, MB_SET_ERR, MB_SUCCESS, MBHEX, MBQUAD, MBTET, MBTRI, mdbImpl, number_dimensions(), moab::ReadBlockData::startExoId, moab::ReadBlockData::startMBId, and moab::CN::SubEntityNodeIndices().

Referenced by read_sidesets().

dimension_exists()

bool moab::ReadNCDF::dimension_exists ( const char *  attrib_name )

private

exodus_file()

int moab::ReadNCDF::exodus_file ( )

private

exodus file bound to this object

factory()

ReaderIface * moab::ReadNCDF::factory ( Interface *  iface )

static

Definition at line 135 of file ReadNCDF.cpp.

{
    return new ReadNCDF( iface );
}

References iface, and ReadNCDF().

Referenced by moab::ReaderWriterSet::ReaderWriterSet().

find_side_element_type()

ErrorCode moab::ReadNCDF::find_side_element_type ( const int  element_id, ExoIIElementType &  type, ReadBlockData &  block_data, int &  df_index, int  side_id  )

private

Definition at line 1677 of file ReadNCDF.cpp.

{
    std::vector< ReadBlockData >::iterator iter, end_iter;
    iter      = blocksLoading.begin();
    end_iter  = blocksLoading.end();
    elem_type = EXOII_MAX_ELEM_TYPE;
 
    for( ; iter != end_iter; ++iter )
    {
        if( exodus_id >= ( *iter ).startExoId && exodus_id < ( *iter ).startExoId + ( *iter ).numElements )
        {
            elem_type = ( *iter ).elemType;
 
            // If we're not reading this block in
            if( !( *iter ).reading_in )
            {
                // Offset df_indes according to type
                if( elem_type >= EXOII_HEX && elem_type <= EXOII_HEX27 )
                    df_index += 4;
                else if( elem_type >= EXOII_TETRA && elem_type <= EXOII_TETRA14 )
                    df_index += 3;
                else if( elem_type >= EXOII_QUAD && elem_type <= EXOII_QUAD9 )
                    df_index += 2;
                else if( elem_type >= EXOII_SHELL && elem_type <= EXOII_SHELL9 )
                {
                    if( side_id == 1 || side_id == 2 )
                        df_index += 4;
                    else
                        df_index += 2;
                }
                else if( elem_type >= EXOII_TRI && elem_type <= EXOII_TRI7 )
                    df_index += 3;
 
                return MB_FAILURE;
            }
 
            block_data = *iter;
 
            return MB_SUCCESS;
        }
    }
    return MB_FAILURE;
}

References blocksLoading, moab::EXOII_HEX, moab::EXOII_HEX27, moab::EXOII_MAX_ELEM_TYPE, moab::EXOII_QUAD, moab::EXOII_QUAD9, moab::EXOII_SHELL, moab::EXOII_SHELL9, moab::EXOII_TETRA, moab::EXOII_TETRA14, moab::EXOII_TRI, moab::EXOII_TRI7, and MB_SUCCESS.

Referenced by create_ss_elements().

get_number_nodes()

int moab::ReadNCDF::get_number_nodes ( EntityHandle  handle )

private

get_type() [1/2]

ErrorCode moab::ReadNCDF::get_type ( char *  exo_element_type, EntityType &  elem_type  )

private

map a character exodusII element type to a TSTT type & topology

get_type() [2/2]

ErrorCode moab::ReadNCDF::get_type ( EntityType &  elem_type, std::string &  exo_element_type  )

private

load_file()

ErrorCode moab::ReadNCDF::load_file ( const char *  file_name, const EntityHandle *  file_set, const FileOptions &  opts, const SubsetList *  subset_list = 0, const Tag *  file_id_tag = 0  )

virtual

load an ExoII file

Implements moab::ReaderIface.

Definition at line 284 of file ReadNCDF.cpp.

{
    ErrorCode status;
    int fail;
 
    int num_blocks            = 0;
    const int* blocks_to_load = 0;
    if( subset_list )
    {
        if( subset_list->tag_list_length > 1 || !strcmp( subset_list->tag_list[0].tag_name, MATERIAL_SET_TAG_NAME ) )
        {
            MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "ExodusII reader supports subset read only by material ID" );
        }
        if( subset_list->num_parts )
        {
            MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "ExodusII reader does not support mesh partitioning" );
        }
        blocks_to_load = subset_list->tag_list[0].tag_values;
        num_blocks     = subset_list->tag_list[0].num_tag_values;
    }
 
    // This function directs the reading of an exoii file, but doesn't do any of
    // the actual work
 
    // See if opts has tdata.
    ErrorCode rval;
    std::string s;
    rval = opts.get_str_option( "tdata", s );
    if( MB_SUCCESS == rval && !s.empty() )
        return update( exodus_file_name, opts, num_blocks, blocks_to_load, *file_set );
 
    reset();
 
    // 0. Open the file.
 
    // open netcdf/exodus file
    fail = nc_open( exodus_file_name, 0, &ncFile );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FILE_DOES_NOT_EXIST, "ReadNCDF:: problem opening Netcdf/Exodus II file " << exodus_file_name );
    }
 
    // 1. Read the header
    status = read_exodus_header();
    if( MB_FAILURE == status ) return status;
 
    status = mdbImpl->get_entities_by_handle( 0, initRange );
    if( MB_FAILURE == status ) return status;
 
    // 2. Read the nodes unless they've already been read before
    status = read_nodes( file_id_tag );
    if( MB_FAILURE == status ) return status;
 
    // 3.
    // extra for polyhedra blocks
    if( numberFaceBlocks_loading > 0 )
    {
        status = read_face_blocks_headers();
        if( MB_FAILURE == status ) return status;
    }
 
    status = read_block_headers( blocks_to_load, num_blocks );
    if( MB_FAILURE == status ) return status;
 
    // 4. Read elements (might not read them, depending on active blocks)
    if( numberFaceBlocks_loading > 0 )
    {
        status = read_polyhedra_faces();
        if( MB_FAILURE == status ) return status;
    }
    status = read_elements( file_id_tag );
    if( MB_FAILURE == status ) return status;
 
    // 5. Read global ids
    status = read_global_ids();
    if( MB_FAILURE == status ) return status;
 
    // 6. Read nodesets
    status = read_nodesets();
    if( MB_FAILURE == status ) return status;
 
    // 7. Read sidesets
    status = read_sidesets();
    if( MB_FAILURE == status ) return status;
 
    // 8. Read qa records
    if( file_set )
    {
        status = read_qa_records( *file_set );
        if( MB_FAILURE == status ) return status;
    }
 
    // What about properties???
 
    // Close the file
    fail = nc_close( ncFile );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "Trouble closing file" );
    }
 
    ncFile = 0;
    return MB_SUCCESS;
}

References ErrorCode, moab::fail(), moab::Interface::get_entities_by_handle(), moab::FileOptions::get_str_option(), initRange, MATERIAL_SET_TAG_NAME, MB_FILE_DOES_NOT_EXIST, MB_SET_ERR, MB_SUCCESS, MB_UNSUPPORTED_OPERATION, mdbImpl, ncFile, moab::ReaderIface::SubsetList::num_parts, moab::ReaderIface::IDTag::num_tag_values, numberFaceBlocks_loading, read_block_headers(), read_elements(), read_exodus_header(), read_face_blocks_headers(), read_global_ids(), read_nodes(), read_nodesets(), read_polyhedra_faces(), read_qa_records(), read_sidesets(), reset(), moab::ReaderIface::SubsetList::tag_list, moab::ReaderIface::SubsetList::tag_list_length, moab::ReaderIface::IDTag::tag_name, moab::ReaderIface::IDTag::tag_values, and update().

number_dimensions()

int moab::ReadNCDF::number_dimensions ( )

inlineprivate

number of dimensions in this exo file

Definition at line 255 of file ReadNCDF.hpp.

{
    return numberDimensions_loading;
}

References numberDimensions_loading.

Referenced by create_ss_elements().

read_block_headers()

ErrorCode moab::ReadNCDF::read_block_headers ( const int *  blocks_to_load, const int  num_blocks  )

private

read block headers, containing info about element type, number, etc.

Definition at line 553 of file ReadNCDF.cpp.

{
    // Get the element block ids; keep this in a separate list,
    // which is not offset by blockIdOffset; this list used later for
    // reading block connectivity
 
    // Get the ids of all the blocks of this file we're reading in
    std::vector< int > block_ids( numberElementBlocks_loading );
    int nc_block_ids = -1;
    GET_1D_INT_VAR( "eb_prop1", nc_block_ids, block_ids );
    if( -1 == nc_block_ids )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting eb_prop1 variable" );
    }
 
    int exodus_id = 1;
 
    // If the active_block_id_list is NULL all blocks are active.
    if( NULL == blocks_to_load || 0 == num_blocks )
    {
        blocks_to_load = &block_ids[0];
    }
 
    std::vector< int > new_blocks( blocks_to_load, blocks_to_load + numberElementBlocks_loading );
 
    std::vector< int >::iterator iter, end_iter;
    iter     = block_ids.begin();
    end_iter = block_ids.end();
 
    // Read header information and initialize header-type block information
    int temp_dim;
    std::vector< char > temp_string_storage( max_str_length + 1 );
    char* temp_string = &temp_string_storage[0];
    int block_seq_id  = 1;
 
    for( ; iter != end_iter; ++iter, block_seq_id++ )
    {
        int num_elements;
 
        GET_DIMB( temp_dim, temp_string, "num_el_in_blk%d", block_seq_id, num_elements );
 
        // Don't read element type string for now, since it's an attrib
        // on the connectivity
        // Get the entity type corresponding to this ExoII element type
        // ExoIIElementType elem_type =
        // ExoIIUtil::static_element_name_to_type(element_type_string);
 
        // Tag each element block(mesh set) with enum for ElementType (SHELL, QUAD4, ....etc)
        ReadBlockData block_data;
        block_data.elemType    = EXOII_MAX_ELEM_TYPE;
        block_data.blockId     = *iter;
        block_data.startExoId  = exodus_id;
        block_data.numElements = num_elements;
 
        // If block is in 'blocks_to_load'----load it!
        if( std::find( new_blocks.begin(), new_blocks.end(), *iter ) != new_blocks.end() )
            block_data.reading_in = true;
        else
            block_data.reading_in = false;
 
        blocksLoading.push_back( block_data );
        exodus_id += num_elements;
    }
 
    return MB_SUCCESS;
}

References moab::ReadBlockData::blockId, blocksLoading, moab::ReadBlockData::elemType, moab::EXOII_MAX_ELEM_TYPE, GET_1D_INT_VAR, GET_DIMB, max_str_length, MB_SET_ERR, MB_SUCCESS, numberElementBlocks_loading, moab::ReadBlockData::numElements, moab::ReadBlockData::reading_in, and moab::ReadBlockData::startExoId.

Referenced by load_file(), and update().

read_elements()

ErrorCode moab::ReadNCDF::read_elements ( const Tag *  file_id_tag )

private

read the element blocks

Definition at line 729 of file ReadNCDF.cpp.

{
    // Read in elements
 
    int result = 0;
 
    // Initialize the nodeInLoadedBlocks vector
    if( nodesInLoadedBlocks.size() < (size_t)( numberNodes_loading + 1 ) )
        nodesInLoadedBlocks.resize( numberNodes_loading + 1 );  // this could be repeated?
    memset( &nodesInLoadedBlocks[0], 0, ( numberNodes_loading + 1 ) * sizeof( char ) );
 
    std::vector< ReadBlockData >::iterator this_it;
    this_it = blocksLoading.begin();
 
    int temp_dim;
    std::vector< char > temp_string_storage( max_str_length + 1 );
    char* temp_string = &temp_string_storage[0];
 
    int nc_var;
    int block_seq_id = 1;
    std::vector< int > dims;
    size_t start[2] = { 0, 0 }, count[2];
 
    for( ; this_it != blocksLoading.end(); ++this_it, block_seq_id++ )
    {
        // If this block isn't to be read in --- continue
        if( !( *this_it ).reading_in ) continue;
 
        // Get some information about this block
        int block_id       = ( *this_it ).blockId;
        EntityHandle* conn = 0;
 
        // Get the ncdf connect variable and the element type
        INS_ID( temp_string, "connect%d", block_seq_id, max_str_length + 1 );
        GET_VAR( temp_string, nc_var, dims );
        if( -1 == nc_var || 0 == nc_var )
        {  // try other var, for polyhedra blocks
            // it could be polyhedra block, defined by fbconn and NFACED attribute
            INS_ID( temp_string, "facconn%d", block_seq_id, max_str_length + 1 );
            GET_VAR( temp_string, nc_var, dims );
 
            if( -1 == nc_var || 0 == nc_var )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connec or faccon variable" );
            }
        }
        nc_type att_type;
        size_t att_len;
        int fail = nc_inq_att( ncFile, nc_var, "elem_type", &att_type, &att_len );
        if( NC_NOERR != fail )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type attribute" );
        }
        std::vector< char > dum_str( att_len + 1 );
        fail = nc_get_att_text( ncFile, nc_var, "elem_type", &dum_str[0] );
        if( NC_NOERR != fail )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type" );
        }
        ExoIIElementType elem_type = ExoIIUtil::static_element_name_to_type( &dum_str[0] );
        ( *this_it ).elemType      = elem_type;
 
        int verts_per_element    = ExoIIUtil::VerticesPerElement[( *this_it ).elemType];
        int number_nodes         = ( *this_it ).numElements * verts_per_element;
        const EntityType mb_type = ExoIIUtil::ExoIIElementMBEntity[( *this_it ).elemType];
 
        if( mb_type == MBPOLYGON )
        {
            // need to read another variable, num_nod_per_el, and
            int num_nodes_poly_conn;
            GET_DIMB( temp_dim, temp_string, "num_nod_per_el%d", block_seq_id, num_nodes_poly_conn );
            // read the connec1 array, which is of dimensions num_nodes_poly_conn (only one )
            std::vector< int > connec;
            connec.resize( num_nodes_poly_conn );
            count[0] = num_nodes_poly_conn;
 
            fail = nc_get_vara_int( ncFile, nc_var, start, count, &connec[0] );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon connectivity " );
            }
            // ebepecnt1
            std::vector< int > ebec;
            ebec.resize( this_it->numElements );
            // Get the ncdf connect variable and the element type
            INS_ID( temp_string, "ebepecnt%d", block_seq_id, max_str_length + 1 );
            GET_VAR( temp_string, nc_var, dims );
            count[0] = this_it->numElements;
            fail     = nc_get_vara_int( ncFile, nc_var, start, count, &ebec[0] );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon nodes per elem " );
            }
            EntityHandle ms_handle;
            if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, ms_handle ) != MB_SUCCESS )
                return MB_FAILURE;
            // create polygons one by one, and put them in the list
            // also need to read the number of nodes for each polygon, then create
            std::vector< EntityHandle > polyConn;
            int ix = 0;
            for( int i = 0; i < this_it->numElements; i++ )
            {
                polyConn.resize( ebec[i] );
                for( int j = 0; j < ebec[i]; j++ )
                {
                    nodesInLoadedBlocks[connec[ix]] = 1;
                    polyConn[j]                     = vertexOffset + connec[ix++];
                }
                EntityHandle newp;
                /*
                 *  ErrorCode create_element(const EntityType type,
                                             const EntityHandle *connectivity,
                                             const int num_vertices,
                                             EntityHandle &element_handle)
                 */
                if( mdbImpl->create_element( MBPOLYGON, &polyConn[0], ebec[i], newp ) != MB_SUCCESS ) return MB_FAILURE;
 
                // add the element in order
                this_it->polys.push_back( newp );
                if( mdbImpl->add_entities( ms_handle, &newp, 1 ) != MB_SUCCESS ) return MB_FAILURE;
            }
            // Set the block id with an offset
            if( mdbImpl->tag_set_data( mMaterialSetTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
            if( mdbImpl->tag_set_data( mGlobalIdTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
        }
        else if( mb_type == MBPOLYHEDRON )
        {
            // need to read another variable, num_fac_per_el
            int num_fac_per_el;
            GET_DIMB( temp_dim, temp_string, "num_fac_per_el%d", block_seq_id, num_fac_per_el );
            // read the fbconn array, which is of dimensions num_nod_per_fa (only one dimension)
            std::vector< int > facconn;
            facconn.resize( num_fac_per_el );
            count[0] = num_fac_per_el;
 
            fail = nc_get_vara_int( ncFile, nc_var, start, count, &facconn[0] );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon connectivity " );
            }
            // ebepecnt1
            std::vector< int > ebepecnt;
            ebepecnt.resize( this_it->numElements );
            // Get the ncdf connect variable and the element type
            INS_ID( temp_string, "ebepecnt%d", block_seq_id, max_str_length + 1 );
            GET_VAR( temp_string, nc_var, dims );
            count[0] = this_it->numElements;
            fail     = nc_get_vara_int( ncFile, nc_var, start, count, &ebepecnt[0] );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon nodes per elem " );
            }
            EntityHandle ms_handle;
            if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, ms_handle ) != MB_SUCCESS )
                return MB_FAILURE;
            // create polygons one by one, and put them in the list
            // also need to read the number of nodes for each polygon, then create
            std::vector< EntityHandle > polyConn;
            int ix = 0;
            for( int i = 0; i < this_it->numElements; i++ )
            {
                polyConn.resize( ebepecnt[i] );
                for( int j = 0; j < ebepecnt[i]; j++ )
                {
                    polyConn[j] = polyfaces[facconn[ix++] - 1];
                }
                EntityHandle newp;
                /*
                 *  ErrorCode create_element(const EntityType type,
                                             const EntityHandle *connectivity,
                                             const int num_vertices,
                                             EntityHandle &element_handle)
                 */
                if( mdbImpl->create_element( MBPOLYHEDRON, &polyConn[0], ebepecnt[i], newp ) != MB_SUCCESS )
                    return MB_FAILURE;
 
                // add the element in order
                this_it->polys.push_back( newp );
                if( mdbImpl->add_entities( ms_handle, &newp, 1 ) != MB_SUCCESS ) return MB_FAILURE;
            }
            // Set the block id with an offset
            if( mdbImpl->tag_set_data( mMaterialSetTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
            if( mdbImpl->tag_set_data( mGlobalIdTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
        }
        else  // this is regular
        {
            // Allocate an array to read in connectivity data
            readMeshIface->get_element_connect( this_it->numElements, verts_per_element, mb_type, this_it->startExoId,
                                                this_it->startMBId, conn );
 
            // Create a range for this sequence of elements
            EntityHandle start_range, end_range;
            start_range = ( *this_it ).startMBId;
            end_range   = start_range + ( *this_it ).numElements - 1;
 
            Range new_range( start_range, end_range );
            // Range<EntityHandle> new_range((*this_it).startMBId,
            //                              (*this_it).startMBId + (*this_it).numElements - 1);
 
            // Create a MBSet for this block and set the material tag
 
            EntityHandle ms_handle;
            if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, ms_handle ) != MB_SUCCESS )
                return MB_FAILURE;
 
            if( mdbImpl->add_entities( ms_handle, new_range ) != MB_SUCCESS ) return MB_FAILURE;
 
            int mid_nodes[4];
            CN::HasMidNodes( mb_type, verts_per_element, mid_nodes );
            if( mdbImpl->tag_set_data( mHasMidNodesTag, &ms_handle, 1, mid_nodes ) != MB_SUCCESS ) return MB_FAILURE;
 
            // Just a check because the following code won't work if this case fails
            assert( sizeof( EntityHandle ) >= sizeof( int ) );
 
            // tmp_ptr is of type int* and points at the same place as conn
            int* tmp_ptr = reinterpret_cast< int* >( conn );
 
            // Read the connectivity into that memory,  this will take only part of the array
            // 1/2 if sizeof(EntityHandle) == 64 bits.
            count[0] = this_it->numElements;
            count[1] = verts_per_element;
            fail     = nc_get_vara_int( ncFile, nc_var, start, count, tmp_ptr );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connectivity" );
            }
 
            // Convert from exodus indices to vertex handles.
            // Iterate backwards in case handles are larger than ints.
            for( int i = number_nodes - 1; i >= 0; --i )
            {
                if( (unsigned)tmp_ptr[i] >= nodesInLoadedBlocks.size() )
                {
                    MB_SET_ERR( MB_FAILURE, "Invalid node ID in block connectivity" );
                }
                nodesInLoadedBlocks[tmp_ptr[i]] = 1;
                conn[i]                         = static_cast< EntityHandle >( tmp_ptr[i] ) + vertexOffset;
            }
 
            // Adjust connectivity order if necessary
            const int* reorder = exodus_elem_order_map[mb_type][verts_per_element];
            if( reorder ) ReadUtilIface::reorder( reorder, conn, this_it->numElements, verts_per_element );
 
            readMeshIface->update_adjacencies( ( *this_it ).startMBId, ( *this_it ).numElements,
                                               ExoIIUtil::VerticesPerElement[( *this_it ).elemType], conn );
 
            if( result == -1 )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: error getting element connectivity for block " << block_id );
            }
 
            // Set the block id with an offset
            if( mdbImpl->tag_set_data( mMaterialSetTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
            if( mdbImpl->tag_set_data( mGlobalIdTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
 
            if( file_id_tag )
            {
                Range range;
                range.insert( this_it->startMBId, this_it->startMBId + this_it->numElements - 1 );
                readMeshIface->assign_ids( *file_id_tag, range, this_it->startExoId );
            }
        }  // end regular block (not polygon)
    }
 
    return MB_SUCCESS;
}

References moab::Interface::add_entities(), moab::ReadUtilIface::assign_ids(), blocksLoading, moab::Interface::create_element(), moab::Interface::create_meshset(), exodus_elem_order_map, moab::ExoIIUtil::ExoIIElementMBEntity, moab::fail(), GET_DIMB, moab::ReadUtilIface::get_element_connect(), GET_VAR, moab::CN::HasMidNodes(), INS_ID, moab::Range::insert(), max_str_length, MB_SET_ERR, MB_SUCCESS, MBPOLYGON, MBPOLYHEDRON, mdbImpl, MESHSET_SET, MESHSET_TRACK_OWNER, mGlobalIdTag, mHasMidNodesTag, mMaterialSetTag, ncFile, nodesInLoadedBlocks, numberNodes_loading, polyfaces, readMeshIface, moab::ReadUtilIface::reorder(), moab::ExoIIUtil::static_element_name_to_type(), moab::Interface::tag_set_data(), moab::ReadUtilIface::update_adjacencies(), vertexOffset, and moab::ExoIIUtil::VerticesPerElement.

Referenced by load_file().

read_exodus_header()

ErrorCode moab::ReadNCDF::read_exodus_header ( )

private

read the header from the ExoII file

Definition at line 393 of file ReadNCDF.cpp.

{
    CPU_WORD_SIZE = sizeof( double );  // With ExodusII version 2, all floats
    IO_WORD_SIZE  = sizeof( double );  // should be changed to doubles
 
    // NetCDF doesn't check its own limits on file read, so check
    // them here so it doesn't corrupt memory any more than absolutely
    // necessary.
    int ndims;
    int fail = nc_inq_ndims( ncFile, &ndims );
    if( NC_NOERR != fail || ndims > NC_MAX_DIMS )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF: File contains " << ndims << " dims but NetCDF library supports only "
                                                           << (int)NC_MAX_DIMS );
    }
    int nvars;
    fail = nc_inq_nvars( ncFile, &nvars );
    if( nvars > NC_MAX_VARS )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF: File contains " << nvars << " vars but NetCDF library supports only "
                                                           << (int)NC_MAX_VARS );
    }
 
    // Get the attributes
 
    // Get the word size, scalar value
    nc_type att_type;
    size_t att_len;
    fail = nc_inq_att( ncFile, NC_GLOBAL, "floating_point_word_size", &att_type, &att_len );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting floating_point_word_size attribute" );
    }
    if( att_type != NC_INT || att_len != 1 )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Word size didn't have type int or size 1" );
    }
    fail = nc_get_att_int( ncFile, NC_GLOBAL, "floating_point_word_size", &IO_WORD_SIZE );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Trouble getting word size" );
    }
 
    // Exodus version
    fail = nc_inq_att( ncFile, NC_GLOBAL, "version", &att_type, &att_len );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting version attribute" );
    }
    if( att_type != NC_FLOAT || att_len != 1 )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Version didn't have type float or size 1" );
    }
    // float version = temp_att->as_float(0);
 
    // Read in initial variables
    int temp_dim;
    GET_DIM( temp_dim, "num_dim", numberDimensions_loading );
    GET_DIM( temp_dim, "num_nodes", numberNodes_loading );
    GET_DIM( temp_dim, "num_elem", numberElements_loading );
    GET_DIM( temp_dim, "num_el_blk", numberElementBlocks_loading );
    GET_DIM( temp_dim, "num_fa_blk", numberFaceBlocks_loading );  // for polyhedra blocks
    GET_DIM( temp_dim, "num_elem", numberElements_loading );
    GET_DIM( temp_dim, "num_node_sets", numberNodeSets_loading );
    GET_DIM( temp_dim, "num_side_sets", numberSideSets_loading );
    GET_DIM( temp_dim, "len_string", max_str_length );
    GET_DIM( temp_dim, "len_line", max_line_length );
 
    // Title; why are we even bothering if we're not going to keep it???
    fail = nc_inq_att( ncFile, NC_GLOBAL, "title", &att_type, &att_len );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting title attribute" );
    }
    if( att_type != NC_CHAR )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: title didn't have type char" );
    }
    char* title = new char[att_len + 1];
    fail        = nc_get_att_text( ncFile, NC_GLOBAL, "title", title );
    if( NC_NOERR != fail )
    {
        delete[] title;
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: trouble getting title" );
    }
    delete[] title;
 
    return MB_SUCCESS;
}

References CPU_WORD_SIZE, moab::fail(), GET_DIM, IO_WORD_SIZE, max_line_length, max_str_length, MB_SET_ERR, MB_SUCCESS, ncFile, numberDimensions_loading, numberElementBlocks_loading, numberElements_loading, numberFaceBlocks_loading, numberNodes_loading, numberNodeSets_loading, and numberSideSets_loading.

Referenced by load_file(), read_tag_values(), and update().

read_face_blocks_headers()

ErrorCode moab::ReadNCDF::read_face_blocks_headers ( )

private

Definition at line 620 of file ReadNCDF.cpp.

{
    // Get the ids of all the blocks of this file we're reading in
    std::vector< int > fblock_ids( numberFaceBlocks_loading );
    int nc_fblock_ids = -1;
    GET_1D_INT_VAR( "fa_prop1", nc_fblock_ids, fblock_ids );
    if( -1 == nc_fblock_ids )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting fa_prop1 variable" );
    }
 
    int temp_dim;
    std::vector< char > temp_string_storage( max_str_length + 1 );
    char* temp_string = &temp_string_storage[0];
    // int block_seq_id = 1;
    int exodus_id = 1;
 
    for( int i = 1; i <= numberFaceBlocks_loading; i++ )
    {
        int num_elements;
 
        GET_DIMB( temp_dim, temp_string, "num_fa_in_blk%d", i, num_elements );
 
        // Don't read element type string for now, since it's an attrib
        // on the connectivity
        // Get the entity type corresponding to this ExoII element type
        // ExoIIElementType elem_type =
        // ExoIIUtil::static_element_name_to_type(element_type_string);
 
        // Tag each element block(mesh set) with enum for ElementType (SHELL, QUAD4, ....etc)
        ReadFaceBlockData fblock_data;
        fblock_data.faceBlockId = i;
        fblock_data.startExoId  = exodus_id;
        fblock_data.numElements = num_elements;
 
        faceBlocksLoading.push_back( fblock_data );
        exodus_id += num_elements;
    }
    return MB_SUCCESS;
}

References moab::ReadFaceBlockData::faceBlockId, faceBlocksLoading, GET_1D_INT_VAR, GET_DIMB, max_str_length, MB_SET_ERR, MB_SUCCESS, numberFaceBlocks_loading, moab::ReadFaceBlockData::numElements, and moab::ReadFaceBlockData::startExoId.

Referenced by load_file().

read_global_ids()

ErrorCode moab::ReadNCDF::read_global_ids ( )

private

read in the global element ids

Definition at line 996 of file ReadNCDF.cpp.

{
    // Read in the map from the exodus file
    std::vector< int > ids( std::max( numberElements_loading, numberNodes_loading ) );
 
    int varid = -1;
    GET_1D_INT_VAR( "elem_map", varid, ids );
    if( -1 != varid )
    {
        std::vector< ReadBlockData >::iterator iter;
        int id_pos = 0;
        for( iter = blocksLoading.begin(); iter != blocksLoading.end(); ++iter )
        {
            if( iter->reading_in )
            {
                if( iter->elemType != EXOII_POLYGON && iter->elemType != EXOII_POLYHEDRON )
                {
                    if( iter->startMBId != 0 )
                    {
                        Range range( iter->startMBId, iter->startMBId + iter->numElements - 1 );
                        ErrorCode error = mdbImpl->tag_set_data( mGlobalIdTag, range, &ids[id_pos] );
                        if( error != MB_SUCCESS ) return error;
                        id_pos += iter->numElements;
                    }
                    else
                        return MB_FAILURE;
                }
                else  // polygons or polyhedrons; use id from elements
                {
                    for( std::vector< EntityHandle >::iterator eit = iter->polys.begin(); eit != iter->polys.end();
                         eit++, id_pos++ )
                    {
                        EntityHandle peh = *eit;
                        if( mdbImpl->tag_set_data( mGlobalIdTag, &peh, 1, &ids[id_pos] ) != MB_SUCCESS )
                            return MB_FAILURE;
                    }
                }
            }
        }
    }
 
    // Read in node map next
    varid = -1;
    GET_1D_INT_VAR( "node_num_map", varid, ids );
    if( -1 != varid )
    {
        Range range( MB_START_ID + vertexOffset, MB_START_ID + vertexOffset + numberNodes_loading - 1 );
        ErrorCode error = mdbImpl->tag_set_data( mGlobalIdTag, range, &ids[0] );
        if( MB_SUCCESS != error ) MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem setting node global ids" );
    }
 
    return MB_SUCCESS;
}

References blocksLoading, moab::error(), ErrorCode, moab::EXOII_POLYGON, moab::EXOII_POLYHEDRON, GET_1D_INT_VAR, MB_SET_ERR, MB_START_ID, MB_SUCCESS, mdbImpl, mGlobalIdTag, numberElements_loading, numberNodes_loading, moab::Interface::tag_set_data(), and vertexOffset.

Referenced by load_file().

read_nodes()

ErrorCode moab::ReadNCDF::read_nodes ( const Tag *  file_id_tag )

private

read the nodes

Definition at line 483 of file ReadNCDF.cpp.

{
    // Read the nodes into memory
 
    assert( 0 != ncFile );
 
    // Create a sequence to hold the node coordinates
    // Get the current number of entities and start at the next slot
 
    EntityHandle node_handle = 0;
    std::vector< double* > arrays;
    readMeshIface->get_node_coords( 3, numberNodes_loading, MB_START_ID, node_handle, arrays );
 
    vertexOffset = ID_FROM_HANDLE( node_handle ) - MB_START_ID;
 
    // Read in the coordinates
    int fail;
    int coord = 0;
    nc_inq_varid( ncFile, "coord", &coord );
 
    // Single var for all coords
    size_t start[2] = { 0, 0 }, count[2] = { 1, static_cast< size_t >( numberNodes_loading ) };
    if( coord )
    {
        for( int d = 0; d < numberDimensions_loading; ++d )
        {
            start[0] = d;
            fail     = nc_get_vara_double( ncFile, coord, start, count, arrays[d] );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting " << (char)( 'x' + d ) << " coord array" );
            }
        }
    }
    // Var for each coord
    else
    {
        char varname[] = "coord ";
        for( int d = 0; d < numberDimensions_loading; ++d )
        {
            varname[5] = 'x' + (char)d;
            fail       = nc_inq_varid( ncFile, varname, &coord );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting " << (char)( 'x' + d ) << " coord variable" );
            }
 
            fail = nc_get_vara_double( ncFile, coord, start, &count[1], arrays[d] );
            if( NC_NOERR != fail )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting " << (char)( 'x' + d ) << " coord array" );
            }
        }
    }
 
    // Zero out any coord values that are in database but not in file
    // (e.g. if MOAB has 3D coords but file is 2D then set Z coords to zero.)
    for( unsigned d = numberDimensions_loading; d < arrays.size(); ++d )
        std::fill( arrays[d], arrays[d] + numberNodes_loading, 0.0 );
 
    if( file_id_tag )
    {
        Range nodes;
        nodes.insert( node_handle, node_handle + numberNodes_loading - 1 );
        readMeshIface->assign_ids( *file_id_tag, nodes, vertexOffset );
    }
 
    return MB_SUCCESS;
}

References moab::ReadUtilIface::assign_ids(), moab::fail(), moab::ReadUtilIface::get_node_coords(), moab::ID_FROM_HANDLE(), moab::Range::insert(), MB_SET_ERR, MB_START_ID, MB_SUCCESS, ncFile, numberDimensions_loading, numberNodes_loading, readMeshIface, and vertexOffset.

Referenced by load_file().

read_nodesets()

ErrorCode moab::ReadNCDF::read_nodesets ( )

private

read the nodesets into meshsets

Definition at line 1050 of file ReadNCDF.cpp.

{
    // Read in the nodesets for the model
 
    if( 0 == numberNodeSets_loading ) return MB_SUCCESS;
    std::vector< int > id_array( numberNodeSets_loading );
 
    // Read in the nodeset ids
    int nc_var;
    GET_1D_INT_VAR( "ns_prop1", nc_var, id_array );
    if( -1 == nc_var )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting ns_prop1 variable" );
    }
 
    // Use a vector of ints to read node handles
    std::vector< int > node_handles;
 
    int i;
    std::vector< char > temp_string_storage( max_str_length + 1 );
    char* temp_string = &temp_string_storage[0];
    int temp_dim;
    for( i = 0; i < numberNodeSets_loading; i++ )
    {
        // Get nodeset parameters
        int number_nodes_in_set;
        int number_dist_factors_in_set;
 
        GET_DIMB( temp_dim, temp_string, "num_nod_ns%d", i + 1, number_nodes_in_set );
        GET_DIMB( temp_dim, temp_string, "num_df_ns%d", i + 1, number_dist_factors_in_set );
 
        // Need to new a vector to store dist. factors
        // This vector * gets stored as a tag on the sideset meshset
        std::vector< double > temp_dist_factor_vector( number_nodes_in_set );
        if( number_dist_factors_in_set != 0 )
        {
            INS_ID( temp_string, "dist_fact_ns%d", i + 1, max_str_length + 1 );
            GET_1D_DBL_VAR( temp_string, temp_dim, temp_dist_factor_vector );
            if( -1 == temp_dim )
            {
                MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting dist fact variable" );
            }
        }
 
        // Size new arrays and get ids and distribution factors
        if( node_handles.size() < (unsigned int)number_nodes_in_set )
        {
            node_handles.reserve( number_nodes_in_set );
            node_handles.resize( number_nodes_in_set );
        }
 
        INS_ID( temp_string, "node_ns%d", i + 1, max_str_length + 1 );
        int temp_var = -1;
        GET_1D_INT_VAR( temp_string, temp_var, node_handles );
        if( -1 == temp_var )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting nodeset node variable" );
        }
 
        // Maybe there is already a nodesets meshset here we can append to
        Range child_meshsets;
        if( mdbImpl->get_entities_by_handle( 0, child_meshsets ) != MB_SUCCESS ) return MB_FAILURE;
 
        child_meshsets = subtract( child_meshsets, initRange );
 
        Range::iterator iter, end_iter;
        iter     = child_meshsets.begin();
        end_iter = child_meshsets.end();
 
        EntityHandle ns_handle = 0;
        for( ; iter != end_iter; ++iter )
        {
            int nodeset_id;
            if( mdbImpl->tag_get_data( mDirichletSetTag, &( *iter ), 1, &nodeset_id ) != MB_SUCCESS ) continue;
 
            if( id_array[i] == nodeset_id )
            {
                // Found the meshset
                ns_handle = *iter;
                break;
            }
        }
 
        std::vector< EntityHandle > nodes_of_nodeset;
        if( ns_handle )
            if( mdbImpl->get_entities_by_handle( ns_handle, nodes_of_nodeset, true ) != MB_SUCCESS ) return MB_FAILURE;
 
        // Make these into entity handles
        // TODO: could we have read it into EntityHandle sized array in the first place?
        int j, temp;
        std::vector< EntityHandle > nodes;
        std::vector< double > dist_factor_vector;
        for( j = 0; j < number_nodes_in_set; j++ )
        {
            // See if this node is one we're currently reading in
            if( nodesInLoadedBlocks[node_handles[j]] == 1 )
            {
                // Make sure that it already isn't in a nodeset
                unsigned int node_id = CREATE_HANDLE( MBVERTEX, node_handles[j] + vertexOffset, temp );
                if( !ns_handle ||
                    std::find( nodes_of_nodeset.begin(), nodes_of_nodeset.end(), node_id ) == nodes_of_nodeset.end() )
                {
                    nodes.push_back( node_id );
 
                    if( number_dist_factors_in_set != 0 ) dist_factor_vector.push_back( temp_dist_factor_vector[j] );
                }
            }
        }
 
        // No nodes to add
        if( nodes.empty() ) continue;
 
        // If there was no meshset found --> create one
        if( ns_handle == 0 )
        {
            if( mdbImpl->create_meshset( MESHSET_ORDERED | MESHSET_TRACK_OWNER, ns_handle ) != MB_SUCCESS )
                return MB_FAILURE;
 
            // Set a tag signifying dirichlet bc
            // TODO: create this tag another way
 
            int nodeset_id = id_array[i];
            if( mdbImpl->tag_set_data( mDirichletSetTag, &ns_handle, 1, &nodeset_id ) != MB_SUCCESS ) return MB_FAILURE;
            if( mdbImpl->tag_set_data( mGlobalIdTag, &ns_handle, 1, &nodeset_id ) != MB_SUCCESS ) return MB_FAILURE;
 
            if( !dist_factor_vector.empty() )
            {
                int size         = dist_factor_vector.size();
                const void* data = &dist_factor_vector[0];
                if( mdbImpl->tag_set_by_ptr( mDistFactorTag, &ns_handle, 1, &data, &size ) != MB_SUCCESS )
                    return MB_FAILURE;
            }
        }
        else if( !dist_factor_vector.empty() )
        {
            // Append dist factors to vector
            const void* ptr = 0;
            int size        = 0;
            if( mdbImpl->tag_get_by_ptr( mDistFactorTag, &ns_handle, 1, &ptr, &size ) != MB_SUCCESS ) return MB_FAILURE;
 
            const double* data = reinterpret_cast< const double* >( ptr );
            dist_factor_vector.reserve( dist_factor_vector.size() + size );
            std::copy( data, data + size, std::back_inserter( dist_factor_vector ) );
            size = dist_factor_vector.size();
            ptr  = &dist_factor_vector[0];
            if( mdbImpl->tag_set_by_ptr( mDistFactorTag, &ns_handle, 1, &ptr, &size ) != MB_SUCCESS ) return MB_FAILURE;
        }
 
        // Add the nodes to the meshset
        if( mdbImpl->add_entities( ns_handle, &nodes[0], nodes.size() ) != MB_SUCCESS ) return MB_FAILURE;
    }
 
    return MB_SUCCESS;
}

References moab::Interface::add_entities(), moab::Range::begin(), moab::CREATE_HANDLE(), moab::Interface::create_meshset(), moab::Range::end(), GET_1D_DBL_VAR, GET_1D_INT_VAR, GET_DIMB, moab::Interface::get_entities_by_handle(), initRange, INS_ID, max_str_length, MB_SET_ERR, MB_SUCCESS, MBVERTEX, mdbImpl, mDirichletSetTag, mDistFactorTag, MESHSET_TRACK_OWNER, mGlobalIdTag, nodesInLoadedBlocks, numberNodeSets_loading, size, moab::subtract(), moab::Interface::tag_get_by_ptr(), moab::Interface::tag_get_data(), moab::Interface::tag_set_by_ptr(), moab::Interface::tag_set_data(), and vertexOffset.

Referenced by load_file().

read_polyhedra_faces()

ErrorCode moab::ReadNCDF::read_polyhedra_faces ( )

private

Definition at line 660 of file ReadNCDF.cpp.

{
    int temp_dim;
    std::vector< char > temp_string_storage( max_str_length + 1 );
    char* temp_string = &temp_string_storage[0];
    nodesInLoadedBlocks.resize( numberNodes_loading + 1 );
    size_t start[1] = { 0 }, count[1];  // one dim arrays only here!
    std::vector< ReadFaceBlockData >::iterator this_it;
    this_it = faceBlocksLoading.begin();
 
    int fblock_seq_id = 1;
    std::vector< int > dims;
    int nc_var, fail;
 
    for( ; this_it != faceBlocksLoading.end(); ++this_it, fblock_seq_id++ )
    {
        int num_fa_in_blk;
        GET_DIMB( temp_dim, temp_string, "num_fa_in_blk%d", fblock_seq_id, num_fa_in_blk );
        int num_nod_per_fa;
        GET_DIMB( temp_dim, temp_string, "num_nod_per_fa%d", fblock_seq_id, num_nod_per_fa );
        // Get the ncdf connect variable and the element type
        INS_ID( temp_string, "fbconn%d", fblock_seq_id, max_str_length + 1 );
        GET_VAR( temp_string, nc_var, dims );
        std::vector< int > fbconn;
        fbconn.resize( num_nod_per_fa );
        count[0] = num_nod_per_fa;
 
        fail = nc_get_vara_int( ncFile, nc_var, start, count, &fbconn[0] );
        if( NC_NOERR != fail )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting face polyhedra connectivity " );
        }
        std::vector< int > fbepecnt;
        fbepecnt.resize( num_fa_in_blk );
        INS_ID( temp_string, "fbepecnt%d", fblock_seq_id, max_str_length + 1 );
        GET_VAR( temp_string, nc_var, dims );
        count[0] = num_fa_in_blk;
        fail     = nc_get_vara_int( ncFile, nc_var, start, count, &fbepecnt[0] );
        if( NC_NOERR != fail )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting face polyhedra connectivity " );
        }
        // now we can create some polygons
        std::vector< EntityHandle > polyConn;
        int ix = 0;
        for( int i = 0; i < num_fa_in_blk; i++ )
        {
            polyConn.resize( fbepecnt[i] );
            for( int j = 0; j < fbepecnt[i]; j++ )
            {
                nodesInLoadedBlocks[fbconn[ix]] = 1;
                polyConn[j]                     = vertexOffset + fbconn[ix++];
            }
            EntityHandle newp;
            /*
             *  ErrorCode create_element(const EntityType type,
                                         const EntityHandle *connectivity,
                                         const int num_vertices,
                                         EntityHandle &element_handle)
             */
            if( mdbImpl->create_element( MBPOLYGON, &polyConn[0], fbepecnt[i], newp ) != MB_SUCCESS ) return MB_FAILURE;
 
            // add the element in order
            polyfaces.push_back( newp );
        }
    }
    return MB_SUCCESS;
}

References moab::Interface::create_element(), faceBlocksLoading, moab::fail(), GET_DIMB, GET_VAR, INS_ID, max_str_length, MB_SET_ERR, MB_SUCCESS, MBPOLYGON, mdbImpl, ncFile, nodesInLoadedBlocks, numberNodes_loading, polyfaces, and vertexOffset.

Referenced by load_file().

read_qa_information()

ErrorCode moab::ReadNCDF::read_qa_information ( std::vector< std::string > &  qa_record_list )

private

Definition at line 1751 of file ReadNCDF.cpp.

{
    // Inquire on the genesis file to find the number of qa records
 
    int number_records = 0;
 
    int temp_dim;
    GET_DIM( temp_dim, "num_qa_rec", number_records );
    std::vector< char > data( max_str_length + 1 );
 
    for( int i = 0; i < number_records; i++ )
    {
        for( int j = 0; j < 4; j++ )
        {
            data[max_str_length] = '\0';
            if( read_qa_string( &data[0], i, j ) != MB_SUCCESS ) return MB_FAILURE;
 
            qa_record_list.push_back( &data[0] );
        }
    }
 
    return MB_SUCCESS;
}

References GET_DIM, max_str_length, MB_SUCCESS, and read_qa_string().

Referenced by read_qa_records().

read_qa_records()

ErrorCode moab::ReadNCDF::read_qa_records ( EntityHandle  file_set )

private

Definition at line 1725 of file ReadNCDF.cpp.

{
    std::vector< std::string > qa_records;
    read_qa_information( qa_records );
 
    std::vector< char > tag_data;
    for( std::vector< std::string >::iterator i = qa_records.begin(); i != qa_records.end(); ++i )
    {
        std::copy( i->begin(), i->end(), std::back_inserter( tag_data ) );
        tag_data.push_back( '\0' );
    }
 
    // If there were qa_records...tag them to the mCurrentMeshHandle
    if( !tag_data.empty() )
    {
        const void* ptr = &tag_data[0];
        int size        = tag_data.size();
        if( mdbImpl->tag_set_by_ptr( mQaRecordTag, &file_set, 1, &ptr, &size ) != MB_SUCCESS )
        {
            return MB_FAILURE;
        }
    }
 
    return MB_SUCCESS;
}

References MB_SUCCESS, mdbImpl, mQaRecordTag, read_qa_information(), size, and moab::Interface::tag_set_by_ptr().

Referenced by load_file().

read_qa_string()

ErrorCode moab::ReadNCDF::read_qa_string ( char *  string, int  record_number, int  record_position  )

private

Definition at line 1775 of file ReadNCDF.cpp.

{
    std::vector< int > dims;
    int temp_var = -1;
    GET_VAR( "qa_records", temp_var, dims );
    if( -1 == temp_var )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting qa record variable" );
        return MB_FAILURE;
    }
    size_t count[3], start[3];
    start[0] = record_number;
    start[1] = record_position;
    start[2] = 0;
    count[0] = 1;
    count[1] = 1;
    count[2] = max_str_length;
    int fail = nc_get_vara_text( ncFile, temp_var, start, count, temp_string );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem setting current record number variable position" );
    }
    // Get the variable id in the exodus file
 
    return MB_SUCCESS;
}

References moab::fail(), GET_VAR, max_str_length, MB_SET_ERR, MB_SUCCESS, and ncFile.

Referenced by read_qa_information().

read_sidesets()

ErrorCode moab::ReadNCDF::read_sidesets ( )

private

read the sidesets (does nothing for now)

Definition at line 1205 of file ReadNCDF.cpp.

{
    // Uhhh if you read the same file (previously_read==true) then blow away the
    // sidesets pertaining to this file?  How do you do that? If you read a
    // new file make sure all the offsets are correct.
 
    // If not loading any sidesets -- exit
    if( 0 == numberSideSets_loading ) return MB_SUCCESS;
 
    // Read in the sideset ids
    std::vector< int > id_array( numberSideSets_loading );
    int temp_var = -1;
    GET_1D_INT_VAR( "ss_prop1", temp_var, id_array );
    if( -1 == temp_var )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting ss_prop1 variable" );
    }
 
    // Create a side set for each one
    int number_sides_in_set;
    int number_dist_factors_in_set;
 
    // Maybe there is already a sidesets meshset here we can append to
    Range child_meshsets;
    if( mdbImpl->get_entities_by_type( 0, MBENTITYSET, child_meshsets ) != MB_SUCCESS ) return MB_FAILURE;
 
    child_meshsets = subtract( child_meshsets, initRange );
 
    Range::iterator iter, end_iter;
 
    int i;
    std::vector< char > temp_string_storage( max_str_length + 1 );
    char* temp_string = &temp_string_storage[0];
    int temp_dim;
    for( i = 0; i < numberSideSets_loading; i++ )
    {
        // Get sideset parameters
        GET_DIMB( temp_dim, temp_string, "num_side_ss%d", i + 1, number_sides_in_set );
        GET_DIMB( temp_dim, temp_string, "num_df_ss%d", i + 1, number_dist_factors_in_set );
 
        // Size new arrays and get element and side lists
        std::vector< int > side_list( number_sides_in_set );
        std::vector< int > element_list( number_sides_in_set );
        INS_ID( temp_string, "side_ss%d", i + 1, max_str_length + 1 );
        temp_var = -1;
        GET_1D_INT_VAR( temp_string, temp_var, side_list );
        if( -1 == temp_var )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting sideset side variable" );
        }
 
        INS_ID( temp_string, "elem_ss%d", i + 1, max_str_length + 1 );
        temp_var = -1;
        GET_1D_INT_VAR( temp_string, temp_var, element_list );
        if( -1 == temp_var )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting sideset elem variable" );
        }
 
        std::vector< double > temp_dist_factor_vector;
        std::vector< EntityHandle > entities_to_add, reverse_entities;
        // Create the sideset entities
        if( create_ss_elements( &element_list[0], &side_list[0], number_sides_in_set, number_dist_factors_in_set,
                                entities_to_add, reverse_entities, temp_dist_factor_vector, i + 1 ) != MB_SUCCESS )
            return MB_FAILURE;
 
        // If there are elements to add
        if( !entities_to_add.empty() || !reverse_entities.empty() )
        {
            iter     = child_meshsets.begin();
            end_iter = child_meshsets.end();
 
            EntityHandle ss_handle = 0;
            for( ; iter != end_iter; ++iter )
            {
                int sideset_id;
                if( mdbImpl->tag_get_data( mNeumannSetTag, &( *iter ), 1, &sideset_id ) != MB_SUCCESS ) continue;
 
                if( id_array[i] == sideset_id )
                {
                    // Found the meshset
                    ss_handle = *iter;
                    break;
                }
            }
 
            // If we didn't find a sideset already
            if( ss_handle == 0 )
            {
                if( mdbImpl->create_meshset( MESHSET_ORDERED | MESHSET_TRACK_OWNER, ss_handle ) != MB_SUCCESS )
                    return MB_FAILURE;
 
                if( ss_handle == 0 ) return MB_FAILURE;
 
                int sideset_id = id_array[i];
                if( mdbImpl->tag_set_data( mNeumannSetTag, &ss_handle, 1, &sideset_id ) != MB_SUCCESS )
                    return MB_FAILURE;
                if( mdbImpl->tag_set_data( mGlobalIdTag, &ss_handle, 1, &sideset_id ) != MB_SUCCESS ) return MB_FAILURE;
 
                if( !reverse_entities.empty() )
                {
                    // Also make a reverse set to put in this set
                    EntityHandle reverse_set;
                    if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, reverse_set ) != MB_SUCCESS )
                        return MB_FAILURE;
 
                    // Add the reverse set to the sideset set and the entities to the reverse set
                    ErrorCode result = mdbImpl->add_entities( ss_handle, &reverse_set, 1 );
                    if( MB_SUCCESS != result ) return result;
 
                    result = mdbImpl->add_entities( reverse_set, &reverse_entities[0], reverse_entities.size() );
                    if( MB_SUCCESS != result ) return result;
 
                    // Set the reverse tag
                    Tag sense_tag;
                    int dum_sense = 0;
                    result        = mdbImpl->tag_get_handle( "NEUSET_SENSE", 1, MB_TYPE_INTEGER, sense_tag,
                                                             MB_TAG_SPARSE | MB_TAG_CREAT, &dum_sense );
                    if( result != MB_SUCCESS ) return result;
                    dum_sense = -1;
                    result    = mdbImpl->tag_set_data( sense_tag, &reverse_set, 1, &dum_sense );
                    if( result != MB_SUCCESS ) return result;
                }
            }
 
            if( mdbImpl->add_entities( ss_handle, ( entities_to_add.empty() ) ? NULL : &entities_to_add[0],
                                       entities_to_add.size() ) != MB_SUCCESS )
                return MB_FAILURE;
 
            // Distribution factor stuff
            if( number_dist_factors_in_set )
            {
                // If this sideset does not already has a distribution factor array...set one
                const void* ptr = 0;
                int size        = 0;
                if( MB_SUCCESS == mdbImpl->tag_get_by_ptr( mDistFactorTag, &ss_handle, 1, &ptr, &size ) )
                {
                    const double* data = reinterpret_cast< const double* >( ptr );
                    std::copy( data, data + size, std::back_inserter( temp_dist_factor_vector ) );
                }
 
                ptr  = &temp_dist_factor_vector[0];
                size = temp_dist_factor_vector.size();
                if( mdbImpl->tag_set_by_ptr( mDistFactorTag, &ss_handle, 1, &ptr, &size ) != MB_SUCCESS )
                    return MB_FAILURE;
            }
        }
    }
 
    return MB_SUCCESS;
}

References moab::Interface::add_entities(), moab::Range::begin(), moab::Interface::create_meshset(), create_ss_elements(), moab::Range::end(), ErrorCode, GET_1D_INT_VAR, GET_DIMB, moab::Interface::get_entities_by_type(), initRange, INS_ID, max_str_length, MB_SET_ERR, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_SPARSE, MB_TYPE_INTEGER, MBENTITYSET, mdbImpl, mDistFactorTag, MESHSET_SET, MESHSET_TRACK_OWNER, mGlobalIdTag, mNeumannSetTag, numberSideSets_loading, size, moab::subtract(), moab::Interface::tag_get_by_ptr(), moab::Interface::tag_get_data(), moab::Interface::tag_get_handle(), moab::Interface::tag_set_by_ptr(), and moab::Interface::tag_set_data().

Referenced by load_file().

read_tag_values()

ErrorCode moab::ReadNCDF::read_tag_values ( const char *  file_name, const char *  tag_name, const FileOptions &  opts, std::vector< int > &  tag_values_out, const SubsetList *  subset_list = 0  )

virtual

Read tag values from a file.

Read the list if all integer tag values from the file for a tag that is a single integer value per entity.

Parameters

file_name	The file to read.
tag_name	The tag for which to read values
tag_values_out	Output: The list of tag values.
subset_list	An array of tag name and value sets specifying the subset of the file to read. If multiple tags are specified, the sets that match all tags (intersection) should be read.
subset_list_length	The length of the 'subset_list' array.

Implements moab::ReaderIface.

Definition at line 214 of file ReadNCDF.cpp.

{
    if( subset_list )
    {
        MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "ExodusII reader supports subset read only by material ID" );
    }
 
    // Open netcdf/exodus file
    int fail = nc_open( file_name, 0, &ncFile );
    if( NC_NOWRITE != fail )
    {
        MB_SET_ERR( MB_FILE_DOES_NOT_EXIST, "ReadNCDF:: problem opening Netcdf/Exodus II file " << file_name );
    }
 
    // 1. Read the header
    ErrorCode rval = read_exodus_header();
    if( MB_FAILURE == rval ) return rval;
 
    int count = 0;
    const char* prop;
    const char* blocks   = "eb_prop1";
    const char* nodesets = "ns_prop1";
    const char* sidesets = "ss_prop1";
 
    if( !strcmp( tag_name, MATERIAL_SET_TAG_NAME ) )
    {
        count = numberElementBlocks_loading;
        prop  = blocks;
    }
    else if( !strcmp( tag_name, DIRICHLET_SET_TAG_NAME ) )
    {
        count = numberNodeSets_loading;
        prop  = nodesets;
    }
    else if( !strcmp( tag_name, NEUMANN_SET_TAG_NAME ) )
    {
        count = numberSideSets_loading;
        prop  = sidesets;
    }
    else
    {
        ncFile = 0;
        return MB_TAG_NOT_FOUND;
    }
 
    if( count )
    {
        int nc_var = -1;
        GET_1D_INT_VAR( prop, nc_var, id_array );
        if( !nc_var )
        {
            MB_SET_ERR( MB_FAILURE, "Problem getting prop variable" );
        }
    }
 
    // Close the file
    fail = nc_close( ncFile );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "Trouble closing file" );
    }
 
    ncFile = 0;
    return MB_SUCCESS;
}

References DIRICHLET_SET_TAG_NAME, ErrorCode, moab::fail(), GET_1D_INT_VAR, MATERIAL_SET_TAG_NAME, MB_FILE_DOES_NOT_EXIST, MB_SET_ERR, MB_SUCCESS, MB_TAG_NOT_FOUND, MB_UNSUPPORTED_OPERATION, ncFile, NEUMANN_SET_TAG_NAME, numberElementBlocks_loading, numberNodeSets_loading, numberSideSets_loading, and read_exodus_header().

reset()

void moab::ReadNCDF::reset ( )

private

Definition at line 189 of file ReadNCDF.cpp.

{
    mCurrentMeshHandle = 0;
    vertexOffset       = 0;
 
    numberNodes_loading         = 0;
    numberElements_loading      = 0;
    numberElementBlocks_loading = 0;
    numberFaceBlocks_loading    = 0;
    numberNodeSets_loading      = 0;
    numberSideSets_loading      = 0;
    numberDimensions_loading    = 0;
 
    if( !blocksLoading.empty() ) blocksLoading.clear();
 
    if( !nodesInLoadedBlocks.empty() ) nodesInLoadedBlocks.clear();
 
    if( !faceBlocksLoading.empty() ) faceBlocksLoading.clear();
}

References blocksLoading, faceBlocksLoading, mCurrentMeshHandle, nodesInLoadedBlocks, numberDimensions_loading, numberElementBlocks_loading, numberElements_loading, numberFaceBlocks_loading, numberNodes_loading, numberNodeSets_loading, numberSideSets_loading, and vertexOffset.

Referenced by load_file(), and ReadNCDF().

tokenize()

void moab::ReadNCDF::tokenize	(	const std::string &	str,
		std::vector< std::string > &	tokens,
		const char *	delimiters
	)

Definition at line 2445 of file ReadNCDF.cpp.

{
    std::string::size_type last = str.find_first_not_of( delimiters, 0 );
    std::string::size_type pos  = str.find_first_of( delimiters, last );
    while( std::string::npos != pos && std::string::npos != last )
    {
        tokens.push_back( str.substr( last, pos - last ) );
        last = str.find_first_not_of( delimiters, pos );
        pos  = str.find_first_of( delimiters, last );
        if( std::string::npos == pos ) pos = str.size();
    }
}

Referenced by update().

update()

ErrorCode moab::ReadNCDF::update	(	const char *	exodus_file_name,
		const FileOptions &	opts,
		const int	num_blocks,
		const int *	blocks_to_load,
		const EntityHandle	file_set
	)

Definition at line 1802 of file ReadNCDF.cpp.

{
    // Function : updating current database from new exodus_file.
    // Creator:   Jane Hu
    // opts is currently designed as following
    // tdata = <var_name>[, time][,op][,destination]
    // where var_name show the tag name to be updated, this version just takes
    // coord.
    // time is the optional, and it gives time step of each of the mesh
    // info in exodus file. It start from 1.
    // op is the operation that is going to be performed on the var_name info.
    // currently support 'set'
    // destination shows where to store the updated info, currently assume it is
    // stored in the same database by replacing the old info if there's no input
    // for destination, or the destination data is given in exodus format and just
    // need to update the coordinates.
    // Assumptions:
    // 1. Assume the num_el_blk's in both DB and update exodus file are the same.
    // 2. Assume num_el_in_blk1...num_el_in_blk(num_el_blk) numbers are matching, may in
    // different order. example: num_el_in_blk11 = num_el_in_blk22 && num_el_in_blk12 =
    // num_el_in_blk21.
    // 3. In exodus file, get node_num_map
    // 4. loop through the node_num_map, use it to find the node in the cub file.
    // 5. Replace coord[0][n] with coordx[m] + vals_nod_var1(time_step, m) for all directions for
    // matching nodes.
    //  Test: try to get hexes
 
    // *******************************************************************
    // Move nodes to their deformed locations.
    // *******************************************************************
    ErrorCode rval;
    std::string s;
    rval = opts.get_str_option( "tdata", s );
    if( MB_SUCCESS != rval )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem reading file options" );
    }
    std::vector< std::string > tokens;
    tokenize( s, tokens, "," );
 
    // 1. Check for time step to find the match time
    int time_step = 1;
    if( tokens.size() > 1 && !tokens[1].empty() )
    {
        const char* time_s = tokens[1].c_str();
        char* endptr;
        long int pval   = strtol( time_s, &endptr, 0 );
        std::string end = endptr;
        if( !end.empty() )  // Syntax error
            return MB_TYPE_OUT_OF_RANGE;
 
        // Check for overflow (parsing long int, returning int)
        time_step = pval;
        if( pval != (long int)time_step ) return MB_TYPE_OUT_OF_RANGE;
        if( time_step <= 0 ) return MB_TYPE_OUT_OF_RANGE;
    }
 
    // 2. Check for the operations, currently support set.
    const char* op;
    if( tokens.size() < 3 || ( tokens[2] != "set" && tokens[2] != "add" ) )
    {
        MB_SET_ERR( MB_TYPE_OUT_OF_RANGE, "ReadNCDF: invalid operation specified for update" );
    }
    op = tokens[2].c_str();
 
    // Open netcdf/exodus file
    ncFile   = 0;
    int fail = nc_open( exodus_file_name, 0, &ncFile );
    if( !ncFile )
    {
        MB_SET_ERR( MB_FILE_DOES_NOT_EXIST, "ReadNCDF:: problem opening Netcdf/Exodus II file " << exodus_file_name );
    }
 
    rval = read_exodus_header();MB_CHK_ERR( rval );
 
    // Check to make sure that the requested time step exists
    int ncdim = -1;
    int max_time_steps;
    GET_DIM( ncdim, "time_step", max_time_steps );
    if( -1 == ncdim )
    {
        std::cout << "ReadNCDF: could not get number of time steps" << std::endl;
        return MB_FAILURE;
    }
    std::cout << "  Maximum time step=" << max_time_steps << std::endl;
    if( max_time_steps < time_step )
    {
        std::cout << "ReadNCDF: time step is greater than max_time_steps" << std::endl;
        return MB_FAILURE;
    }
 
    // Get the time
    std::vector< double > times( max_time_steps );
    int nc_var = -1;
    GET_1D_DBL_VAR( "time_whole", nc_var, times );
    if( -1 == nc_var )
    {
        std::cout << "ReadNCDF: unable to get time variable" << std::endl;
    }
    else
    {
        std::cout << "  Step " << time_step << " is at " << times[time_step - 1] << " seconds" << std::endl;
    }
 
    // Read in the node_num_map.
    std::vector< int > ptr( numberNodes_loading );
 
    int varid = -1;
    GET_1D_INT_VAR( "node_num_map", varid, ptr );
    if( -1 == varid )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting node number map data" );
    }
 
    // Read in the deformations.
    std::vector< std::vector< double > > deformed_arrays( 3 );
    std::vector< std::vector< double > > orig_coords( 3 );
    deformed_arrays[0].reserve( numberNodes_loading );
    deformed_arrays[1].reserve( numberNodes_loading );
    deformed_arrays[2].reserve( numberNodes_loading );
    orig_coords[0].reserve( numberNodes_loading );
    orig_coords[1].reserve( numberNodes_loading );
    orig_coords[2].reserve( numberNodes_loading );
    size_t start[2] = { static_cast< size_t >( time_step - 1 ), 0 },
           count[2] = { 1, static_cast< size_t >( numberNodes_loading ) };
    std::vector< int > dims;
    int coordx = -1, coordy = -1, coordz = -1;
    GET_VAR( "vals_nod_var1", coordx, dims );
    GET_VAR( "vals_nod_var2", coordy, dims );
    if( numberDimensions_loading == 3 ) GET_VAR( "vals_nod_var3", coordz, dims );
    if( -1 == coordx || -1 == coordy || ( numberDimensions_loading == 3 && -1 == coordz ) )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting coords variable" );
    }
 
    fail = nc_get_vara_double( ncFile, coordx, start, count, &deformed_arrays[0][0] );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting x deformation array" );
    }
 
    fail = nc_get_vara_double( ncFile, coordy, start, count, &deformed_arrays[1][0] );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting y deformation array" );
    }
    if( numberDimensions_loading == 3 )
    {
        fail = nc_get_vara_double( ncFile, coordz, start, count, &deformed_arrays[2][0] );
        if( NC_NOERR != fail )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting z deformation array" );
        }
    }
 
    int coord1 = -1, coord2 = -1, coord3 = -1;
    GET_1D_DBL_VAR( "coordx", coord1, orig_coords[0] );
    if( -1 == coord1 )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting x coord array" );
    }
    GET_1D_DBL_VAR( "coordy", coord2, orig_coords[1] );
    if( -1 == coord2 )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting y coord array" );
    }
    if( numberDimensions_loading == 3 )
    {
        GET_1D_DBL_VAR( "coordz", coord3, orig_coords[2] );
        if( -1 == coord3 )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting z coord array" );
        }
    }
 
    // b. Deal with DB file : get node info. according to node_num_map.
    if( tokens[0] != "coord" && tokens[0] != "COORD" ) return MB_NOT_IMPLEMENTED;
 
    if( strcmp( op, "set" ) != 0 && strcmp( op, " set" ) != 0 ) return MB_NOT_IMPLEMENTED;
 
    // Two methods of matching nodes (id vs. proximity)
    const bool match_node_ids = true;
 
    // Get nodes in cubit file
    Range cub_verts;
    rval = mdbImpl->get_entities_by_type( cub_file_set, MBVERTEX, cub_verts );MB_CHK_ERR( rval );
    std::cout << "  cub_file_set contains " << cub_verts.size() << " nodes." << std::endl;
 
    // Some accounting
    std::cout << "  exodus file contains " << numberNodes_loading << " nodes." << std::endl;
    double max_magnitude     = 0;
    double average_magnitude = 0;
    int found                = 0;
    int lost                 = 0;
    std::map< int, EntityHandle > cub_verts_id_map;
    AdaptiveKDTree kdtree( mdbImpl );
    EntityHandle root;
 
    // Should not use cub verts unless they have been matched. Place in a map
    // for fast handle_by_id lookup.
    std::map< int, EntityHandle > matched_cub_vert_id_map;
 
    // Place cub verts in a map for searching by id
    if( match_node_ids )
    {
        std::vector< int > cub_ids( cub_verts.size() );
        rval = mdbImpl->tag_get_data( mGlobalIdTag, cub_verts, &cub_ids[0] );MB_CHK_ERR( rval );
        for( unsigned i = 0; i != cub_verts.size(); ++i )
        {
            cub_verts_id_map.insert( std::pair< int, EntityHandle >( cub_ids[i], cub_verts[i] ) );
        }
 
        // Place cub verts in a kdtree for searching by proximity
    }
    else
    {
        FileOptions tree_opts( "MAX_PER_LEAF=1;SPLITS_PER_DIR=1;CANDIDATE_PLANE_SET=0" );
        rval = kdtree.build_tree( cub_verts, &root, &tree_opts );MB_CHK_ERR( rval );
        AdaptiveKDTreeIter tree_iter;
        rval = kdtree.get_tree_iterator( root, tree_iter );MB_CHK_ERR( rval );
    }
 
    // For each exo vert, find the matching cub vert
    for( int i = 0; i < numberNodes_loading; ++i )
    {
        int exo_id = ptr[i];
        CartVect exo_coords( orig_coords[0][i], orig_coords[1][i], orig_coords[2][i] );
        EntityHandle cub_vert = 0;
        bool found_match      = false;
 
        // By id
        if( match_node_ids )
        {
            std::map< int, EntityHandle >::iterator i_iter;
            i_iter = cub_verts_id_map.find( exo_id );
            if( i_iter != cub_verts_id_map.end() )
            {
                found_match = true;
                cub_vert    = i_iter->second;
            }
 
            // By proximity
        }
        else
        {
            // The MAX_NODE_DIST is the farthest distance to  search for a node.
            // For the 1/12th symmetry 85 pin model, the max node dist could not be less
            // than 1e-1 (March 26, 2010).
            const double MAX_NODE_DIST = 1e-1;
 
            std::vector< EntityHandle > leaves;
            double min_dist = MAX_NODE_DIST;
            rval            = kdtree.distance_search( exo_coords.array(), MAX_NODE_DIST, leaves );MB_CHK_ERR( rval );
            for( std::vector< EntityHandle >::const_iterator j = leaves.begin(); j != leaves.end(); ++j )
            {
                std::vector< EntityHandle > leaf_verts;
                rval = mdbImpl->get_entities_by_type( *j, MBVERTEX, leaf_verts );MB_CHK_ERR( rval );
                for( std::vector< EntityHandle >::const_iterator k = leaf_verts.begin(); k != leaf_verts.end(); ++k )
                {
                    CartVect orig_cub_coords, difference;
                    rval = mdbImpl->get_coords( &( *k ), 1, orig_cub_coords.array() );MB_CHK_ERR( rval );
                    difference  = orig_cub_coords - exo_coords;
                    double dist = difference.length();
                    if( dist < min_dist )
                    {
                        min_dist = dist;
                        cub_vert = *k;
                    }
                }
            }
            if( 0 != cub_vert ) found_match = true;
        }
 
        // If a match is found, update it with the deformed coords from the exo file.
        if( found_match )
        {
            CartVect updated_exo_coords;
            matched_cub_vert_id_map.insert( std::pair< int, EntityHandle >( exo_id, cub_vert ) );
            updated_exo_coords[0] = orig_coords[0][i] + deformed_arrays[0][i];
            updated_exo_coords[1] = orig_coords[1][i] + deformed_arrays[1][i];
 
            if( numberDimensions_loading == 3 ) updated_exo_coords[2] = orig_coords[2][i] + deformed_arrays[2][i];
            rval = mdbImpl->set_coords( &cub_vert, 1, updated_exo_coords.array() );MB_CHK_ERR( rval );
            ++found;
 
            double magnitude =
                sqrt( deformed_arrays[0][i] * deformed_arrays[0][i] + deformed_arrays[1][i] * deformed_arrays[1][i] +
                      deformed_arrays[2][i] * deformed_arrays[2][i] );
            if( magnitude > max_magnitude ) max_magnitude = magnitude;
            average_magnitude += magnitude;
        }
        else
        {
            ++lost;
            std::cout << "cannot match exo node " << exo_id << " " << exo_coords << std::endl;
        }
    }
 
    // Exo verts that could not be matched have already been printed. Now print the
    // cub verts that could not be matched.
    if( matched_cub_vert_id_map.size() < cub_verts.size() )
    {
        Range unmatched_cub_verts = cub_verts;
        for( std::map< int, EntityHandle >::const_iterator i = matched_cub_vert_id_map.begin();
             i != matched_cub_vert_id_map.end(); ++i )
        {
            unmatched_cub_verts.erase( i->second );
        }
 
        for( Range::const_iterator i = unmatched_cub_verts.begin(); i != unmatched_cub_verts.end(); ++i )
        {
            int cub_id;
            rval = mdbImpl->tag_get_data( mGlobalIdTag, &( *i ), 1, &cub_id );MB_CHK_ERR( rval );
 
            CartVect cub_coords;
            rval = mdbImpl->get_coords( &( *i ), 1, cub_coords.array() );MB_CHK_ERR( rval );
            std::cout << "cannot match cub node " << cub_id << " " << cub_coords << std::endl;
        }
 
        std::cout << "  " << unmatched_cub_verts.size() << " nodes from the cub file could not be matched."
                  << std::endl;
    }
 
    // Summarize statistics
    std::cout << "  " << found << " nodes from the exodus file were matched in the cub_file_set ";
    if( match_node_ids )
    {
        std::cout << "by id." << std::endl;
    }
    else
    {
        std::cout << "by proximity." << std::endl;
    }
 
    // Fail if all of the nodes could not be matched.
    if( 0 != lost )
    {
        std::cout << "Error:  " << lost << " nodes from the exodus file could not be matched." << std::endl;
        // return MB_FAILURE;
    }
    std::cout << "  maximum node displacement magnitude: " << max_magnitude << " cm" << std::endl;
    std::cout << "  average node displacement magnitude: " << average_magnitude / found << " cm" << std::endl;
 
    // *******************************************************************
    // Remove dead elements from the MOAB instance.
    // *******************************************************************
 
    // How many element variables are in the file?
    int n_elem_var;
    GET_DIM( ncdim, "num_elem_var", n_elem_var );
 
    // Get element variable names
    varid       = -1;
    int cstatus = nc_inq_varid( ncFile, "name_elem_var", &varid );
    std::vector< char > names_memory( n_elem_var * max_str_length );
    std::vector< char* > names( n_elem_var );
    for( int i = 0; i < n_elem_var; ++i )
        names[i] = &names_memory[i * max_str_length];
    if( cstatus != NC_NOERR || varid == -1 )
    {
        std::cout << "ReadNCDF: name_elem_var does not exist" << std::endl;
        return MB_FAILURE;
    }
    int status = nc_get_var_text( ncFile, varid, &names_memory[0] );
    if( NC_NOERR != status )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem getting element variable names" );
    }
 
    // Is one of the element variable names "death_status"? If so, get its index
    // in the element variable array.
    int death_index;
    bool found_death_index = false;
    for( int i = 0; i < n_elem_var; ++i )
    {
        std::string temp( names[i] );
        std::string::size_type pos0 = temp.find( "death_status" );
        std::string::size_type pos1 = temp.find( "Death_Status" );
        std::string::size_type pos2 = temp.find( "DEATH_STATUS" );
        if( std::string::npos != pos0 || std::string::npos != pos1 || std::string::npos != pos2 )
        {
            found_death_index = true;
            death_index       = i + 1;  // NetCDF variables start with 1
            break;
        }
    }
    if( !found_death_index )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem getting index of death_status variable" );
    }
 
    // The exodus header has already been read. This contains the number of element
    // blocks.
 
    // Dead elements are listed by block. Read the block headers to determine how
    // many elements are in each block.
    rval = read_block_headers( blocks_to_load, num_blocks );
    if( MB_FAILURE == rval )
    {
        MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem reading block headers" );
    }
 
    // Dead elements from the Exodus file can be located in the cub_file_set by id
    // or by connectivity. Currently, finding elements by id requires careful book
    // keeping when constructing the model in Cubit. To avoid this, one can match
    // elements by connectivity instead.
    const bool match_elems_by_connectivity = true;
 
    // Get the element id map. The ids in the map are from the elements in the blocks.
    // elem_num_map(blk1 elem ids, blk2 elem ids, blk3 elem ids, ...)
    std::vector< int > elem_ids( numberNodes_loading );
    if( !match_elems_by_connectivity )
    {
        GET_1D_INT_VAR( "elem_num_map", varid, elem_ids );
        if( -1 == varid )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem getting element number map data" );
        }
    }
 
    // For each block
    int first_elem_id_in_block = 0;
    int block_count            = 1;  // NetCDF variables start with 1
    int total_elems            = 0;
    int total_dead_elems       = 0;
    for( std::vector< ReadBlockData >::iterator i = blocksLoading.begin(); i != blocksLoading.end(); ++i )
    {
        // Get the ncdf connect variable
        std::string temp_string( "connect" );
        std::stringstream temp_ss;
        temp_ss << block_count;
        temp_string += temp_ss.str();
        temp_string += "\0";
        std::vector< int > ldims;
        GET_VAR( temp_string.c_str(), nc_var, ldims );
        if( !nc_var )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connect variable" );
        }
        // The element type is an attribute of the connectivity variable
        nc_type att_type;
        size_t att_len;
        fail = nc_inq_att( ncFile, nc_var, "elem_type", &att_type, &att_len );
        if( NC_NOERR != fail )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type attribute" );
        }
        std::vector< char > dum_str( att_len + 1 );
        fail = nc_get_att_text( ncFile, nc_var, "elem_type", &dum_str[0] );
        if( NC_NOERR != fail )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type" );
        }
        ExoIIElementType elem_type = ExoIIUtil::static_element_name_to_type( &dum_str[0] );
        ( *i ).elemType            = elem_type;
        const EntityType mb_type   = ExoIIUtil::ExoIIElementMBEntity[( *i ).elemType];
 
        // Get the number of nodes per element
        unsigned int nodes_per_element = ExoIIUtil::VerticesPerElement[( *i ).elemType];
 
        // Read the connectivity into that memory.
        // int exo_conn[i->numElements][nodes_per_element];
        int* exo_conn = new int[i->numElements * nodes_per_element];
        // NcBool status = temp_var->get(&exo_conn[0][0], i->numElements, nodes_per_element);
        size_t lstart[2] = { 0, 0 }, lcount[2];
        lcount[0]        = i->numElements;
        lcount[1]        = nodes_per_element;
        fail             = nc_get_vara_int( ncFile, nc_var, lstart, lcount, exo_conn );
        if( NC_NOERR != fail )
        {
            delete[] exo_conn;
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connectivity" );
        }
 
        // Get the death_status at the correct time step.
        std::vector< double > death_status( i->numElements );  // It seems wrong, but it uses doubles
        std::string array_name( "vals_elem_var" );
        temp_ss.str( "" );  // stringstream won't clear by temp.clear()
        temp_ss << death_index;
        array_name += temp_ss.str();
        array_name += "eb";
        temp_ss.str( "" );  // stringstream won't clear by temp.clear()
        temp_ss << block_count;
        array_name += temp_ss.str();
        array_name += "\0";
        GET_VAR( array_name.c_str(), nc_var, ldims );
        if( !nc_var )
        {
            MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting death status variable" );
        }
        lstart[0] = time_step - 1;
        lstart[1] = 0;
        lcount[0] = 1;
        lcount[1] = i->numElements;
        status    = nc_get_vara_double( ncFile, nc_var, lstart, lcount, &death_status[0] );
        if( NC_NOERR != status )
        {
            delete[] exo_conn;
            MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem setting time step for death_status" );
        }
 
        // Look for dead elements. If there is too many dead elements and this starts
        // to take too long, I should put the elems in a kd-tree for more efficient
        // searching. Alternatively I could get the exo connectivity and match nodes.
        int dead_elem_counter = 0, missing_elem_counter = 0;
        for( int j = 0; j < i->numElements; ++j )
        {
            if( 1 != death_status[j] )
            {
                Range cub_elem, cub_nodes;
                if( match_elems_by_connectivity )
                {
                    // Get exodus nodes for the element
                    std::vector< int > elem_conn( nodes_per_element );
                    for( unsigned int k = 0; k < nodes_per_element; ++k )
                    {
                        // elem_conn[k] = exo_conn[j][k];
                        elem_conn[k] = exo_conn[j * nodes_per_element + k];
                    }
                    // Get the ids of the nodes (assume we are matching by id)
                    // Remember that the exodus array locations start with 1 (not 0).
                    std::vector< int > elem_conn_node_ids( nodes_per_element );
                    for( unsigned int k = 0; k < nodes_per_element; ++k )
                    {
                        elem_conn_node_ids[k] = ptr[elem_conn[k] - 1];
                    }
                    // Get the cub_file_set nodes by id
                    // The map is a log search and takes almost no time.
                    // MOAB's linear tag search takes 5-10 minutes.
                    for( unsigned int k = 0; k < nodes_per_element; ++k )
                    {
                        std::map< int, EntityHandle >::iterator k_iter;
                        k_iter = matched_cub_vert_id_map.find( elem_conn_node_ids[k] );
 
                        if( k_iter == matched_cub_vert_id_map.end() )
                        {
                            std::cout << "ReadNCDF: Found no cub node with id=" << elem_conn_node_ids[k]
                                      << ", but expected to find only 1." << std::endl;
                            break;
                        }
                        cub_nodes.insert( k_iter->second );
                    }
 
                    if( nodes_per_element != cub_nodes.size() )
                    {
                        std::cout << "ReadNCDF: nodes_per_elemenet != cub_nodes.size()" << std::endl;
                        delete[] exo_conn;
                        return MB_INVALID_SIZE;
                    }
 
                    // Get the cub_file_set element with the same nodes
                    int to_dim = CN::Dimension( mb_type );
                    rval       = mdbImpl->get_adjacencies( cub_nodes, to_dim, false, cub_elem );
                    if( MB_SUCCESS != rval )
                    {
                        std::cout << "ReadNCDF: could not get dead cub element" << std::endl;
                        delete[] exo_conn;
                        return rval;
                    }
 
                    // Pronto/Presto renumbers elements, so matching cub and exo elements by
                    // id is not possible at this time.
                }
                else
                {
                    // Get dead element's id
                    int elem_id = elem_ids[first_elem_id_in_block + j];
                    void* id[]  = { &elem_id };
                    // Get the element by id
                    rval = mdbImpl->get_entities_by_type_and_tag( cub_file_set, mb_type, &mGlobalIdTag, id, 1, cub_elem,
                                                                  Interface::INTERSECT );
                    if( MB_SUCCESS != rval )
                    {
                        delete[] exo_conn;
                        return rval;
                    }
                }
 
                if( 1 == cub_elem.size() )
                {
                    // Delete the dead element from the cub file. It will be removed from sets
                    // ONLY if they are tracking meshsets.
                    rval = mdbImpl->remove_entities( cub_file_set, cub_elem );
                    if( MB_SUCCESS != rval )
                    {
                        delete[] exo_conn;
                        return rval;
                    }
                    rval = mdbImpl->delete_entities( cub_elem );
                    if( MB_SUCCESS != rval )
                    {
                        delete[] exo_conn;
                        return rval;
                    }
                }
                else
                {
                    std::cout << "ReadNCDF: Should have found 1 element with  type=" << mb_type
                              << " in cub_file_set, but instead found " << cub_elem.size() << std::endl;
                    rval = mdbImpl->list_entities( cub_nodes );
                    ++missing_elem_counter;
                    delete[] exo_conn;
                    return MB_FAILURE;
                }
                ++dead_elem_counter;
            }
        }
        // Print some statistics
        temp_ss.str( "" );
        temp_ss << i->blockId;
        total_dead_elems += dead_elem_counter;
        total_elems += i->numElements;
        std::cout << "  Block " << temp_ss.str() << " has " << dead_elem_counter << "/" << i->numElements
                  << " dead elements." << std::endl;
        if( 0 != missing_elem_counter )
        {
            std::cout << "    " << missing_elem_counter
                      << " dead elements in this block were not found in the cub_file_set." << std::endl;
        }
 
        // Advance the pointers into element ids and block_count. Memory cleanup.
        first_elem_id_in_block += i->numElements;
        ++block_count;
        delete[] exo_conn;
    }
 
    std::cout << " Total: " << total_dead_elems << "/" << total_elems << " dead elements." << std::endl;
 
    // Close the file
    fail = nc_close( ncFile );
    if( NC_NOERR != fail )
    {
        MB_SET_ERR( MB_FAILURE, "Trouble closing file" );
    }
 
    ncFile = 0;
    return MB_SUCCESS;
}

References moab::CartVect::array(), moab::Range::begin(), blocksLoading, moab::AdaptiveKDTree::build_tree(), moab::Interface::delete_entities(), moab::CN::Dimension(), moab::AdaptiveKDTree::distance_search(), moab::Range::end(), moab::Range::erase(), ErrorCode, moab::ExoIIUtil::ExoIIElementMBEntity, moab::fail(), GET_1D_DBL_VAR, GET_1D_INT_VAR, moab::Interface::get_adjacencies(), moab::Interface::get_coords(), GET_DIM, moab::Interface::get_entities_by_type(), moab::Interface::get_entities_by_type_and_tag(), moab::FileOptions::get_str_option(), moab::AdaptiveKDTree::get_tree_iterator(), GET_VAR, moab::Range::insert(), moab::Interface::INTERSECT, moab::CartVect::length(), moab::Interface::list_entities(), max_str_length, MB_CHK_ERR, MB_FILE_DOES_NOT_EXIST, MB_INVALID_SIZE, MB_NOT_IMPLEMENTED, MB_SET_ERR, MB_SUCCESS, MB_TYPE_OUT_OF_RANGE, MBVERTEX, mdbImpl, mGlobalIdTag, ncFile, numberDimensions_loading, numberNodes_loading, read_block_headers(), read_exodus_header(), moab::Interface::remove_entities(), moab::Interface::set_coords(), moab::Range::size(), moab::ExoIIUtil::static_element_name_to_type(), moab::Interface::tag_get_data(), tokenize(), and moab::ExoIIUtil::VerticesPerElement.

Referenced by load_file().

Member Data Documentation

blocksLoading

std::vector< ReadBlockData > moab::ReadNCDF::blocksLoading

private

Definition at line 231 of file ReadNCDF.hpp.

Referenced by find_side_element_type(), read_block_headers(), read_elements(), read_global_ids(), reset(), and update().

CPU_WORD_SIZE

int moab::ReadNCDF::CPU_WORD_SIZE

private

Definition at line 194 of file ReadNCDF.hpp.

Referenced by read_exodus_header().

exodusFile

std::string moab::ReadNCDF::exodusFile

private

file name

Definition at line 201 of file ReadNCDF.hpp.

faceBlocksLoading

std::vector< ReadFaceBlockData > moab::ReadNCDF::faceBlocksLoading

private

Definition at line 236 of file ReadNCDF.hpp.

Referenced by read_face_blocks_headers(), read_polyhedra_faces(), and reset().

initRange

Range moab::ReadNCDF::initRange

private

range of entities in initial mesh, before this read

Definition at line 251 of file ReadNCDF.hpp.

Referenced by load_file(), read_nodesets(), and read_sidesets().

IO_WORD_SIZE

int moab::ReadNCDF::IO_WORD_SIZE

private

Definition at line 195 of file ReadNCDF.hpp.

Referenced by read_exodus_header().

max_line_length

int moab::ReadNCDF::max_line_length

private

Definition at line 248 of file ReadNCDF.hpp.

Referenced by read_exodus_header().

max_str_length

int moab::ReadNCDF::max_str_length

private

Definition at line 248 of file ReadNCDF.hpp.

Referenced by create_ss_elements(), read_block_headers(), read_elements(), read_exodus_header(), read_face_blocks_headers(), read_nodesets(), read_polyhedra_faces(), read_qa_information(), read_qa_string(), read_sidesets(), and update().

mCurrentMeshHandle

EntityHandle moab::ReadNCDF::mCurrentMeshHandle

private

Meshset Handle for the mesh that is currently being read.

Definition at line 224 of file ReadNCDF.hpp.

Referenced by reset().

mdbImpl

Interface* moab::ReadNCDF::mdbImpl

private

interface instance

Definition at line 190 of file ReadNCDF.hpp.

Referenced by create_sideset_element(), create_ss_elements(), load_file(), read_elements(), read_global_ids(), read_nodesets(), read_polyhedra_faces(), read_qa_records(), read_sidesets(), update(), and ~ReadNCDF().

mDirichletSetTag

Tag moab::ReadNCDF::mDirichletSetTag

private

Definition at line 241 of file ReadNCDF.hpp.

Referenced by read_nodesets(), and ReadNCDF().

mDistFactorTag

Tag moab::ReadNCDF::mDistFactorTag

private

Definition at line 244 of file ReadNCDF.hpp.

Referenced by read_nodesets(), read_sidesets(), and ReadNCDF().

mGlobalIdTag

Tag moab::ReadNCDF::mGlobalIdTag

private

Definition at line 245 of file ReadNCDF.hpp.

Referenced by read_elements(), read_global_ids(), read_nodesets(), read_sidesets(), ReadNCDF(), and update().

mHasMidNodesTag

Tag moab::ReadNCDF::mHasMidNodesTag

private

Definition at line 243 of file ReadNCDF.hpp.

Referenced by read_elements(), and ReadNCDF().

mMaterialSetTag

Tag moab::ReadNCDF::mMaterialSetTag

private

Cached tags for reading. Note that all these tags are defined when the core is initialized.

Definition at line 240 of file ReadNCDF.hpp.

Referenced by read_elements(), and ReadNCDF().

mNeumannSetTag

Tag moab::ReadNCDF::mNeumannSetTag

private

Definition at line 242 of file ReadNCDF.hpp.

Referenced by read_sidesets(), and ReadNCDF().

mQaRecordTag

Tag moab::ReadNCDF::mQaRecordTag

private

Definition at line 246 of file ReadNCDF.hpp.

Referenced by read_qa_records(), and ReadNCDF().

ncFile

int moab::ReadNCDF::ncFile

private

Definition at line 192 of file ReadNCDF.hpp.

Referenced by load_file(), read_elements(), read_exodus_header(), read_nodes(), read_polyhedra_faces(), read_qa_string(), read_tag_values(), ReadNCDF(), and update().

nodesInLoadedBlocks

std::vector< char > moab::ReadNCDF::nodesInLoadedBlocks

private

Definition at line 227 of file ReadNCDF.hpp.

Referenced by read_elements(), read_nodesets(), read_polyhedra_faces(), and reset().

numberDimensions_loading

int moab::ReadNCDF::numberDimensions_loading

private

Definition at line 221 of file ReadNCDF.hpp.

Referenced by number_dimensions(), read_exodus_header(), read_nodes(), reset(), and update().

numberElementBlocks_loading

int moab::ReadNCDF::numberElementBlocks_loading

private

number of blocks in the current exoII file

Definition at line 210 of file ReadNCDF.hpp.

Referenced by read_block_headers(), read_exodus_header(), read_tag_values(), and reset().

numberElements_loading

int moab::ReadNCDF::numberElements_loading

private

number of elements in the current exoII file

Definition at line 207 of file ReadNCDF.hpp.

Referenced by read_exodus_header(), read_global_ids(), and reset().

numberFaceBlocks_loading

int moab::ReadNCDF::numberFaceBlocks_loading

private

number of face blocks in the current exoII file (used for polyhedra)

Definition at line 213 of file ReadNCDF.hpp.

Referenced by load_file(), read_exodus_header(), read_face_blocks_headers(), and reset().

numberNodes_loading

int moab::ReadNCDF::numberNodes_loading

private

number of nodes in the current exoII file

Definition at line 204 of file ReadNCDF.hpp.

Referenced by read_elements(), read_exodus_header(), read_global_ids(), read_nodes(), read_polyhedra_faces(), reset(), and update().

numberNodeSets_loading

int moab::ReadNCDF::numberNodeSets_loading

private

number of nodesets in the current exoII file

Definition at line 216 of file ReadNCDF.hpp.

Referenced by read_exodus_header(), read_nodesets(), read_tag_values(), and reset().

numberSideSets_loading

int moab::ReadNCDF::numberSideSets_loading

private

number of sidesets in the current exoII file

Definition at line 219 of file ReadNCDF.hpp.

Referenced by read_exodus_header(), read_sidesets(), read_tag_values(), and reset().

polyfaces

std::vector< EntityHandle > moab::ReadNCDF::polyfaces

private

Definition at line 233 of file ReadNCDF.hpp.

Referenced by read_elements(), and read_polyhedra_faces().

readMeshIface

ReadUtilIface* moab::ReadNCDF::readMeshIface

private

Definition at line 108 of file ReadNCDF.hpp.

Referenced by read_elements(), read_nodes(), ReadNCDF(), and ~ReadNCDF().

vertexOffset

int moab::ReadNCDF::vertexOffset

private

int to offset vertex ids with

Definition at line 198 of file ReadNCDF.hpp.

Referenced by read_elements(), read_global_ids(), read_nodes(), read_nodesets(), read_polyhedra_faces(), and reset().

---

The documentation for this class was generated from the following files:

• ReadNCDF.hpp
• ReadNCDF.cpp