 |
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
|
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
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... | |
Output Exodus File for VERDE.
Definition at line 73 of file ReadNCDF.hpp.
| moab::ReadNCDF::ReadNCDF | ( | Interface * | impl = NULL | ) |
Constructor.
Get and cache predefined tag handles
Definition at line 140 of file ReadNCDF.cpp.
140 : mdbImpl( impl ), max_line_length( -1 ), max_str_length( -1 )
141 {
142 assert( impl != NULL );
143 reset();
144
145 impl->query_interface( readMeshIface );
146
147 // Initialize in case tag_get_handle fails below
148 mMaterialSetTag = 0;
149 mDirichletSetTag = 0;
150 mNeumannSetTag = 0;
151 mHasMidNodesTag = 0;
152 mDistFactorTag = 0;
153 mQaRecordTag = 0;
154 mGlobalIdTag = 0;
155
156 //! Get and cache predefined tag handles
157 ErrorCode result;
158 const int negone = -1;
159 result = impl->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mMaterialSetTag,
160 MB_TAG_SPARSE | MB_TAG_CREAT, &negone );
161 assert( MB_SUCCESS == result );
162 result = impl->tag_get_handle( DIRICHLET_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mDirichletSetTag,
163 MB_TAG_SPARSE | MB_TAG_CREAT, &negone );
164 assert( MB_SUCCESS == result );
165 result = impl->tag_get_handle( NEUMANN_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mNeumannSetTag,
166 MB_TAG_SPARSE | MB_TAG_CREAT, &negone );
167 assert( MB_SUCCESS == result );
168 const int mids[] = { -1, -1, -1, -1 };
169 result = impl->tag_get_handle( HAS_MID_NODES_TAG_NAME, 4, MB_TYPE_INTEGER, mHasMidNodesTag,
170 MB_TAG_SPARSE | MB_TAG_CREAT, mids );
171 assert( MB_SUCCESS == result );
172 result = impl->tag_get_handle( "distFactor", 0, MB_TYPE_DOUBLE, mDistFactorTag,
173 MB_TAG_SPARSE | MB_TAG_VARLEN | MB_TAG_CREAT );
174 assert( MB_SUCCESS == result );
175 result = impl->tag_get_handle( "qaRecord", 0, MB_TYPE_OPAQUE, mQaRecordTag,
176 MB_TAG_SPARSE | MB_TAG_VARLEN | MB_TAG_CREAT );
177 assert( MB_SUCCESS == result );
178 mGlobalIdTag = impl->globalId_tag();
179
180 assert( MB_SUCCESS == result );
181 #ifdef NDEBUG
182 if( MB_SUCCESS == result )
183 {
184 }; // Line to avoid compiler warning about unused variable
185 #endif
186 ncFile = 0;
187 }
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().
|
virtual |
Destructor.
Definition at line 209 of file ReadNCDF.cpp.
210 { 211 mdbImpl->release_interface( readMeshIface ); 212 }
References mdbImpl, readMeshIface, and moab::Interface::release_interface().
|
private |
creates an element with the given connectivity
Definition at line 1597 of file ReadNCDF.cpp.
1601 {
1602 // Get adjacent entities
1603 ErrorCode error = MB_SUCCESS;
1604 int to_dim = CN::Dimension( type );
1605 // for matching purposes , consider only corner nodes
1606 // basically, assume everything is matching if the corner nodes are matching, and
1607 // the number of total nodes is the same
1608 int nb_corner_nodes = CN::VerticesPerEntity( type );
1609 std::vector< EntityHandle > adj_ent;
1610 mdbImpl->get_adjacencies( &( connectivity[0] ), 1, to_dim, false, adj_ent );
1611
1612 // For each entity, see if we can find a match
1613 // If we find a match, return it
1614 bool match_found = false;
1615 std::vector< EntityHandle > match_conn;
1616 // By default, assume positive sense
1617 sense = 1;
1618 for( unsigned int i = 0; i < adj_ent.size() && match_found == false; i++ )
1619 {
1620 // Get the connectivity
1621 error = mdbImpl->get_connectivity( &( adj_ent[i] ), 1, match_conn );
1622 if( error != MB_SUCCESS ) continue;
1623
1624 // Make sure they have the same number of vertices (higher order elements ?)
1625 if( match_conn.size() != connectivity.size() ) continue;
1626
1627 // Find a matching node
1628 std::vector< EntityHandle >::iterator iter;
1629 std::vector< EntityHandle >::iterator end_corner = match_conn.begin() + nb_corner_nodes;
1630 iter = std::find( match_conn.begin(), end_corner, connectivity[0] );
1631 if( iter == end_corner ) continue;
1632
1633 // Rotate to match connectivity
1634 // rotate only corner nodes, ignore the rest from now on
1635 std::rotate( match_conn.begin(), iter, end_corner );
1636
1637 bool they_match = true;
1638 for( int j = 1; j < nb_corner_nodes; j++ )
1639 {
1640 if( connectivity[j] != match_conn[j] )
1641 {
1642 they_match = false;
1643 break;
1644 }
1645 }
1646
1647 // If we didn't get a match
1648 if( !they_match )
1649 {
1650 // Try the opposite sense
1651 they_match = true;
1652
1653 int k = nb_corner_nodes - 1;
1654 for( int j = 1; j < nb_corner_nodes; )
1655 {
1656 if( connectivity[j] != match_conn[k] )
1657 {
1658 they_match = false;
1659 break;
1660 }
1661 ++j;
1662 --k;
1663 }
1664 // If they matched here, sense is reversed
1665 if( they_match ) sense = -1;
1666 }
1667 match_found = they_match;
1668 if( match_found ) handle = adj_ent[i];
1669 }
1670
1671 // If we didn't find a match, create an element
1672 if( !match_found ) error = mdbImpl->create_element( type, &connectivity[0], connectivity.size(), handle );
1673
1674 return error;
1675 }
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().
|
private |
Definition at line 1357 of file ReadNCDF.cpp.
1365 {
1366 // Determine entity type from element id
1367
1368 // If there are dist. factors, create a vector to hold the array
1369 // and place this array as a tag onto the sideset meshset
1370
1371 std::vector< double > temp_dist_factor_vector( num_dist_factors );
1372 std::vector< char > temp_string_storage( max_str_length + 1 );
1373 char* temp_string = &temp_string_storage[0];
1374 int temp_var;
1375 if( num_dist_factors )
1376 {
1377 INS_ID( temp_string, "dist_fact_ss%d", ss_seq_id, max_str_length + 1 );
1378 temp_var = -1;
1379 GET_1D_DBL_VAR( temp_string, temp_var, temp_dist_factor_vector );
1380 if( -1 == temp_var )
1381 {
1382 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting dist fact variable" );
1383 }
1384 }
1385
1386 EntityType subtype;
1387 int num_side_nodes, num_elem_nodes;
1388 const EntityHandle* nodes;
1389 std::vector< EntityHandle > connectivity;
1390 int side_node_idx[32];
1391
1392 int df_index = 0;
1393 int sense = 0;
1394 for( int i = 0; i < num_sides; i++ )
1395 {
1396 ExoIIElementType exoii_type;
1397 ReadBlockData block_data;
1398 block_data.elemType = EXOII_MAX_ELEM_TYPE;
1399
1400 if( find_side_element_type( element_ids[i], exoii_type, block_data, df_index, side_list[i] ) != MB_SUCCESS )
1401 continue; // Isn't being read in this time
1402
1403 EntityType type = ExoIIUtil::ExoIIElementMBEntity[exoii_type];
1404
1405 EntityHandle ent_handle = element_ids[i] - block_data.startExoId + block_data.startMBId;
1406
1407 const int side_num = side_list[i] - 1;
1408 if( type == MBHEX )
1409 {
1410 // Get the nodes of the element
1411 if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
1412
1413 CN::SubEntityNodeIndices( type, num_elem_nodes, 2, side_num, subtype, num_side_nodes, side_node_idx );
1414 if( num_side_nodes <= 0 ) return MB_FAILURE;
1415
1416 connectivity.resize( num_side_nodes );
1417 for( int k = 0; k < num_side_nodes; ++k )
1418 connectivity[k] = nodes[side_node_idx[k]];
1419
1420 if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) ) return MB_FAILURE;
1421 if( 1 == sense )
1422 entities_to_add.push_back( ent_handle );
1423 else
1424 reverse_entities.push_back( ent_handle );
1425
1426 // Read in distribution factor array
1427 if( num_dist_factors )
1428 {
1429 for( int k = 0; k < 4; k++ )
1430 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1431 }
1432 }
1433
1434 // If it is a Tet
1435 else if( type == MBTET )
1436 {
1437 // Get the nodes of the element
1438 if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
1439
1440 CN::SubEntityNodeIndices( type, num_elem_nodes, 2, side_num, subtype, num_side_nodes, side_node_idx );
1441 if( num_side_nodes <= 0 ) return MB_FAILURE;
1442
1443 connectivity.resize( num_side_nodes );
1444 for( int k = 0; k < num_side_nodes; ++k )
1445 connectivity[k] = nodes[side_node_idx[k]];
1446
1447 if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) ) return MB_FAILURE;
1448 if( 1 == sense )
1449 entities_to_add.push_back( ent_handle );
1450 else
1451 reverse_entities.push_back( ent_handle );
1452
1453 // Read in distribution factor array
1454 if( num_dist_factors )
1455 {
1456 for( int k = 0; k < 3; k++ )
1457 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1458 }
1459 }
1460 else if( type == MBQUAD && exoii_type >= EXOII_SHELL && exoii_type <= EXOII_SHELL9 )
1461 {
1462 // ent_handle = CREATE_HANDLE(MBQUAD, base_id, error);
1463
1464 // Just use this quad
1465 if( side_list[i] == 1 )
1466 {
1467 if( 1 == sense )
1468 entities_to_add.push_back( ent_handle );
1469 else
1470 reverse_entities.push_back( ent_handle );
1471
1472 if( num_dist_factors )
1473 {
1474 for( int k = 0; k < 4; k++ )
1475 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1476 }
1477
1478 continue;
1479 }
1480 else if( side_list[i] == 2 )
1481 {
1482 reverse_entities.push_back( ent_handle );
1483
1484 if( num_dist_factors )
1485 {
1486 for( int k = 0; k < 4; k++ )
1487 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1488 }
1489
1490 continue;
1491 }
1492 else
1493 {
1494 // Get the nodes of the element
1495 if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
1496
1497 CN::SubEntityNodeIndices( type, num_elem_nodes, 1, side_num - 2, subtype, num_side_nodes,
1498 side_node_idx );
1499 if( num_side_nodes <= 0 ) return MB_FAILURE;
1500
1501 connectivity.resize( num_side_nodes );
1502 for( int k = 0; k < num_side_nodes; ++k )
1503 connectivity[k] = nodes[side_node_idx[k]];
1504
1505 if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) )
1506 return MB_FAILURE;
1507 if( 1 == sense )
1508 entities_to_add.push_back( ent_handle );
1509 else
1510 reverse_entities.push_back( ent_handle );
1511
1512 if( num_dist_factors )
1513 {
1514 for( int k = 0; k < 2; k++ )
1515 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1516 }
1517 }
1518 }
1519 // If it is a Quad
1520 else if( type == MBQUAD )
1521 {
1522 // Get the nodes of the element
1523 if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
1524
1525 CN::SubEntityNodeIndices( type, num_elem_nodes, 1, side_num, subtype, num_side_nodes, side_node_idx );
1526 if( num_side_nodes <= 0 ) return MB_FAILURE;
1527
1528 connectivity.resize( num_side_nodes );
1529 for( int k = 0; k < num_side_nodes; ++k )
1530 connectivity[k] = nodes[side_node_idx[k]];
1531
1532 if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) ) return MB_FAILURE;
1533 if( 1 == sense )
1534 entities_to_add.push_back( ent_handle );
1535 else
1536 reverse_entities.push_back( ent_handle );
1537
1538 // Read in distribution factor array
1539 if( num_dist_factors )
1540 {
1541 for( int k = 0; k < 2; k++ )
1542 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1543 }
1544 }
1545 else if( type == MBTRI )
1546 {
1547 int side_offset = 0;
1548 if( number_dimensions() == 3 && side_list[i] <= 2 )
1549 {
1550 if( 1 == sense )
1551 entities_to_add.push_back( ent_handle );
1552 else
1553 reverse_entities.push_back( ent_handle );
1554 if( num_dist_factors )
1555 {
1556 for( int k = 0; k < 3; k++ )
1557 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1558 }
1559 }
1560 else
1561 {
1562 if( number_dimensions() == 3 )
1563 {
1564 if( side_list[i] > 2 ) side_offset = 2;
1565 }
1566
1567 // Get the nodes of the element
1568 if( mdbImpl->get_connectivity( ent_handle, nodes, num_elem_nodes ) != MB_SUCCESS ) return MB_FAILURE;
1569
1570 CN::SubEntityNodeIndices( type, num_elem_nodes, 1, side_num - side_offset, subtype, num_side_nodes,
1571 side_node_idx );
1572 if( num_side_nodes <= 0 ) return MB_FAILURE;
1573
1574 connectivity.resize( num_side_nodes );
1575 for( int k = 0; k < num_side_nodes; ++k )
1576 connectivity[k] = nodes[side_node_idx[k]];
1577
1578 if( MB_SUCCESS != create_sideset_element( connectivity, subtype, ent_handle, sense ) )
1579 return MB_FAILURE;
1580 if( 1 == sense )
1581 entities_to_add.push_back( ent_handle );
1582 else
1583 reverse_entities.push_back( ent_handle );
1584
1585 if( num_dist_factors )
1586 {
1587 for( int k = 0; k < 2; k++ )
1588 dist_factor_vector.push_back( temp_dist_factor_vector[df_index++] );
1589 }
1590 }
1591 }
1592 }
1593
1594 return MB_SUCCESS;
1595 }
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().
|
private |
|
private |
exodus file bound to this object
|
static |
Definition at line 135 of file ReadNCDF.cpp.
136 {
137 return new ReadNCDF( iface );
138 }
References iface, and ReadNCDF().
Referenced by moab::ReaderWriterSet::ReaderWriterSet().
|
private |
Definition at line 1677 of file ReadNCDF.cpp.
1682 {
1683 std::vector< ReadBlockData >::iterator iter, end_iter;
1684 iter = blocksLoading.begin();
1685 end_iter = blocksLoading.end();
1686 elem_type = EXOII_MAX_ELEM_TYPE;
1687
1688 for( ; iter != end_iter; ++iter )
1689 {
1690 if( exodus_id >= ( *iter ).startExoId && exodus_id < ( *iter ).startExoId + ( *iter ).numElements )
1691 {
1692 elem_type = ( *iter ).elemType;
1693
1694 // If we're not reading this block in
1695 if( !( *iter ).reading_in )
1696 {
1697 // Offset df_indes according to type
1698 if( elem_type >= EXOII_HEX && elem_type <= EXOII_HEX27 )
1699 df_index += 4;
1700 else if( elem_type >= EXOII_TETRA && elem_type <= EXOII_TETRA14 )
1701 df_index += 3;
1702 else if( elem_type >= EXOII_QUAD && elem_type <= EXOII_QUAD9 )
1703 df_index += 2;
1704 else if( elem_type >= EXOII_SHELL && elem_type <= EXOII_SHELL9 )
1705 {
1706 if( side_id == 1 || side_id == 2 )
1707 df_index += 4;
1708 else
1709 df_index += 2;
1710 }
1711 else if( elem_type >= EXOII_TRI && elem_type <= EXOII_TRI7 )
1712 df_index += 3;
1713
1714 return MB_FAILURE;
1715 }
1716
1717 block_data = *iter;
1718
1719 return MB_SUCCESS;
1720 }
1721 }
1722 return MB_FAILURE;
1723 }
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().
|
private |
|
private |
map a character exodusII element type to a TSTT type & topology
|
private |
|
virtual |
load an ExoII file
Implements moab::ReaderIface.
Definition at line 284 of file ReadNCDF.cpp.
289 {
290 ErrorCode status;
291 int fail;
292
293 int num_blocks = 0;
294 const int* blocks_to_load = 0;
295 if( subset_list )
296 {
297 if( subset_list->tag_list_length > 1 || !strcmp( subset_list->tag_list[0].tag_name, MATERIAL_SET_TAG_NAME ) )
298 {
299 MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "ExodusII reader supports subset read only by material ID" );
300 }
301 if( subset_list->num_parts )
302 {
303 MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "ExodusII reader does not support mesh partitioning" );
304 }
305 blocks_to_load = subset_list->tag_list[0].tag_values;
306 num_blocks = subset_list->tag_list[0].num_tag_values;
307 }
308
309 // This function directs the reading of an exoii file, but doesn't do any of
310 // the actual work
311
312 // See if opts has tdata.
313 ErrorCode rval;
314 std::string s;
315 rval = opts.get_str_option( "tdata", s );
316 if( MB_SUCCESS == rval && !s.empty() )
317 return update( exodus_file_name, opts, num_blocks, blocks_to_load, *file_set );
318
319 reset();
320
321 // 0. Open the file.
322
323 // open netcdf/exodus file
324 fail = nc_open( exodus_file_name, 0, &ncFile );
325 if( NC_NOERR != fail )
326 {
327 MB_SET_ERR( MB_FILE_DOES_NOT_EXIST, "ReadNCDF:: problem opening Netcdf/Exodus II file " << exodus_file_name );
328 }
329
330 // 1. Read the header
331 status = read_exodus_header();
332 if( MB_FAILURE == status ) return status;
333
334 status = mdbImpl->get_entities_by_handle( 0, initRange );
335 if( MB_FAILURE == status ) return status;
336
337 // 2. Read the nodes unless they've already been read before
338 status = read_nodes( file_id_tag );
339 if( MB_FAILURE == status ) return status;
340
341 // 3.
342 // extra for polyhedra blocks
343 if( numberFaceBlocks_loading > 0 )
344 {
345 status = read_face_blocks_headers();
346 if( MB_FAILURE == status ) return status;
347 }
348
349 status = read_block_headers( blocks_to_load, num_blocks );
350 if( MB_FAILURE == status ) return status;
351
352 // 4. Read elements (might not read them, depending on active blocks)
353 if( numberFaceBlocks_loading > 0 )
354 {
355 status = read_polyhedra_faces();
356 if( MB_FAILURE == status ) return status;
357 }
358 status = read_elements( file_id_tag );
359 if( MB_FAILURE == status ) return status;
360
361 // 5. Read global ids
362 status = read_global_ids();
363 if( MB_FAILURE == status ) return status;
364
365 // 6. Read nodesets
366 status = read_nodesets();
367 if( MB_FAILURE == status ) return status;
368
369 // 7. Read sidesets
370 status = read_sidesets();
371 if( MB_FAILURE == status ) return status;
372
373 // 8. Read qa records
374 if( file_set )
375 {
376 status = read_qa_records( *file_set );
377 if( MB_FAILURE == status ) return status;
378 }
379
380 // What about properties???
381
382 // Close the file
383 fail = nc_close( ncFile );
384 if( NC_NOERR != fail )
385 {
386 MB_SET_ERR( MB_FAILURE, "Trouble closing file" );
387 }
388
389 ncFile = 0;
390 return MB_SUCCESS;
391 }
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().
|
inlineprivate |
number of dimensions in this exo file
Definition at line 255 of file ReadNCDF.hpp.
256 {
257 return numberDimensions_loading;
258 }
References numberDimensions_loading.
Referenced by create_ss_elements().
|
private |
read block headers, containing info about element type, number, etc.
Definition at line 553 of file ReadNCDF.cpp.
554 {
555 // Get the element block ids; keep this in a separate list,
556 // which is not offset by blockIdOffset; this list used later for
557 // reading block connectivity
558
559 // Get the ids of all the blocks of this file we're reading in
560 std::vector< int > block_ids( numberElementBlocks_loading );
561 int nc_block_ids = -1;
562 GET_1D_INT_VAR( "eb_prop1", nc_block_ids, block_ids );
563 if( -1 == nc_block_ids )
564 {
565 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting eb_prop1 variable" );
566 }
567
568 int exodus_id = 1;
569
570 // If the active_block_id_list is NULL all blocks are active.
571 if( NULL == blocks_to_load || 0 == num_blocks )
572 {
573 blocks_to_load = &block_ids[0];
574 }
575
576 std::vector< int > new_blocks( blocks_to_load, blocks_to_load + numberElementBlocks_loading );
577
578 std::vector< int >::iterator iter, end_iter;
579 iter = block_ids.begin();
580 end_iter = block_ids.end();
581
582 // Read header information and initialize header-type block information
583 int temp_dim;
584 std::vector< char > temp_string_storage( max_str_length + 1 );
585 char* temp_string = &temp_string_storage[0];
586 int block_seq_id = 1;
587
588 for( ; iter != end_iter; ++iter, block_seq_id++ )
589 {
590 int num_elements;
591
592 GET_DIMB( temp_dim, temp_string, "num_el_in_blk%d", block_seq_id, num_elements );
593
594 // Don't read element type string for now, since it's an attrib
595 // on the connectivity
596 // Get the entity type corresponding to this ExoII element type
597 // ExoIIElementType elem_type =
598 // ExoIIUtil::static_element_name_to_type(element_type_string);
599
600 // Tag each element block(mesh set) with enum for ElementType (SHELL, QUAD4, ....etc)
601 ReadBlockData block_data;
602 block_data.elemType = EXOII_MAX_ELEM_TYPE;
603 block_data.blockId = *iter;
604 block_data.startExoId = exodus_id;
605 block_data.numElements = num_elements;
606
607 // If block is in 'blocks_to_load'----load it!
608 if( std::find( new_blocks.begin(), new_blocks.end(), *iter ) != new_blocks.end() )
609 block_data.reading_in = true;
610 else
611 block_data.reading_in = false;
612
613 blocksLoading.push_back( block_data );
614 exodus_id += num_elements;
615 }
616
617 return MB_SUCCESS;
618 }
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 the element blocks
Definition at line 729 of file ReadNCDF.cpp.
730 {
731 // Read in elements
732
733 int result = 0;
734
735 // Initialize the nodeInLoadedBlocks vector
736 if( nodesInLoadedBlocks.size() < (size_t)( numberNodes_loading + 1 ) )
737 nodesInLoadedBlocks.resize( numberNodes_loading + 1 ); // this could be repeated?
738 memset( &nodesInLoadedBlocks[0], 0, ( numberNodes_loading + 1 ) * sizeof( char ) );
739
740 std::vector< ReadBlockData >::iterator this_it;
741 this_it = blocksLoading.begin();
742
743 int temp_dim;
744 std::vector< char > temp_string_storage( max_str_length + 1 );
745 char* temp_string = &temp_string_storage[0];
746
747 int nc_var;
748 int block_seq_id = 1;
749 std::vector< int > dims;
750 size_t start[2] = { 0, 0 }, count[2];
751
752 for( ; this_it != blocksLoading.end(); ++this_it, block_seq_id++ )
753 {
754 // If this block isn't to be read in --- continue
755 if( !( *this_it ).reading_in ) continue;
756
757 // Get some information about this block
758 int block_id = ( *this_it ).blockId;
759 EntityHandle* conn = 0;
760
761 // Get the ncdf connect variable and the element type
762 INS_ID( temp_string, "connect%d", block_seq_id, max_str_length + 1 );
763 GET_VAR( temp_string, nc_var, dims );
764 if( -1 == nc_var || 0 == nc_var )
765 { // try other var, for polyhedra blocks
766 // it could be polyhedra block, defined by fbconn and NFACED attribute
767 INS_ID( temp_string, "facconn%d", block_seq_id, max_str_length + 1 );
768 GET_VAR( temp_string, nc_var, dims );
769
770 if( -1 == nc_var || 0 == nc_var )
771 {
772 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connec or faccon variable" );
773 }
774 }
775 nc_type att_type;
776 size_t att_len;
777 int fail = nc_inq_att( ncFile, nc_var, "elem_type", &att_type, &att_len );
778 if( NC_NOERR != fail )
779 {
780 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type attribute" );
781 }
782 std::vector< char > dum_str( att_len + 1 );
783 fail = nc_get_att_text( ncFile, nc_var, "elem_type", &dum_str[0] );
784 if( NC_NOERR != fail )
785 {
786 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type" );
787 }
788 ExoIIElementType elem_type = ExoIIUtil::static_element_name_to_type( &dum_str[0] );
789 ( *this_it ).elemType = elem_type;
790
791 int verts_per_element = ExoIIUtil::VerticesPerElement[( *this_it ).elemType];
792 int number_nodes = ( *this_it ).numElements * verts_per_element;
793 const EntityType mb_type = ExoIIUtil::ExoIIElementMBEntity[( *this_it ).elemType];
794
795 if( mb_type == MBPOLYGON )
796 {
797 // need to read another variable, num_nod_per_el, and
798 int num_nodes_poly_conn;
799 GET_DIMB( temp_dim, temp_string, "num_nod_per_el%d", block_seq_id, num_nodes_poly_conn );
800 // read the connec1 array, which is of dimensions num_nodes_poly_conn (only one )
801 std::vector< int > connec;
802 connec.resize( num_nodes_poly_conn );
803 count[0] = num_nodes_poly_conn;
804
805 fail = nc_get_vara_int( ncFile, nc_var, start, count, &connec[0] );
806 if( NC_NOERR != fail )
807 {
808 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon connectivity " );
809 }
810 // ebepecnt1
811 std::vector< int > ebec;
812 ebec.resize( this_it->numElements );
813 // Get the ncdf connect variable and the element type
814 INS_ID( temp_string, "ebepecnt%d", block_seq_id, max_str_length + 1 );
815 GET_VAR( temp_string, nc_var, dims );
816 count[0] = this_it->numElements;
817 fail = nc_get_vara_int( ncFile, nc_var, start, count, &ebec[0] );
818 if( NC_NOERR != fail )
819 {
820 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon nodes per elem " );
821 }
822 EntityHandle ms_handle;
823 if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, ms_handle ) != MB_SUCCESS )
824 return MB_FAILURE;
825 // create polygons one by one, and put them in the list
826 // also need to read the number of nodes for each polygon, then create
827 std::vector< EntityHandle > polyConn;
828 int ix = 0;
829 for( int i = 0; i < this_it->numElements; i++ )
830 {
831 polyConn.resize( ebec[i] );
832 for( int j = 0; j < ebec[i]; j++ )
833 {
834 nodesInLoadedBlocks[connec[ix]] = 1;
835 polyConn[j] = vertexOffset + connec[ix++];
836 }
837 EntityHandle newp;
838 /*
839 * ErrorCode create_element(const EntityType type,
840 const EntityHandle *connectivity,
841 const int num_vertices,
842 EntityHandle &element_handle)
843 */
844 if( mdbImpl->create_element( MBPOLYGON, &polyConn[0], ebec[i], newp ) != MB_SUCCESS ) return MB_FAILURE;
845
846 // add the element in order
847 this_it->polys.push_back( newp );
848 if( mdbImpl->add_entities( ms_handle, &newp, 1 ) != MB_SUCCESS ) return MB_FAILURE;
849 }
850 // Set the block id with an offset
851 if( mdbImpl->tag_set_data( mMaterialSetTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
852 if( mdbImpl->tag_set_data( mGlobalIdTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
853 }
854 else if( mb_type == MBPOLYHEDRON )
855 {
856 // need to read another variable, num_fac_per_el
857 int num_fac_per_el;
858 GET_DIMB( temp_dim, temp_string, "num_fac_per_el%d", block_seq_id, num_fac_per_el );
859 // read the fbconn array, which is of dimensions num_nod_per_fa (only one dimension)
860 std::vector< int > facconn;
861 facconn.resize( num_fac_per_el );
862 count[0] = num_fac_per_el;
863
864 fail = nc_get_vara_int( ncFile, nc_var, start, count, &facconn[0] );
865 if( NC_NOERR != fail )
866 {
867 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon connectivity " );
868 }
869 // ebepecnt1
870 std::vector< int > ebepecnt;
871 ebepecnt.resize( this_it->numElements );
872 // Get the ncdf connect variable and the element type
873 INS_ID( temp_string, "ebepecnt%d", block_seq_id, max_str_length + 1 );
874 GET_VAR( temp_string, nc_var, dims );
875 count[0] = this_it->numElements;
876 fail = nc_get_vara_int( ncFile, nc_var, start, count, &ebepecnt[0] );
877 if( NC_NOERR != fail )
878 {
879 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting polygon nodes per elem " );
880 }
881 EntityHandle ms_handle;
882 if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, ms_handle ) != MB_SUCCESS )
883 return MB_FAILURE;
884 // create polygons one by one, and put them in the list
885 // also need to read the number of nodes for each polygon, then create
886 std::vector< EntityHandle > polyConn;
887 int ix = 0;
888 for( int i = 0; i < this_it->numElements; i++ )
889 {
890 polyConn.resize( ebepecnt[i] );
891 for( int j = 0; j < ebepecnt[i]; j++ )
892 {
893 polyConn[j] = polyfaces[facconn[ix++] - 1];
894 }
895 EntityHandle newp;
896 /*
897 * ErrorCode create_element(const EntityType type,
898 const EntityHandle *connectivity,
899 const int num_vertices,
900 EntityHandle &element_handle)
901 */
902 if( mdbImpl->create_element( MBPOLYHEDRON, &polyConn[0], ebepecnt[i], newp ) != MB_SUCCESS )
903 return MB_FAILURE;
904
905 // add the element in order
906 this_it->polys.push_back( newp );
907 if( mdbImpl->add_entities( ms_handle, &newp, 1 ) != MB_SUCCESS ) return MB_FAILURE;
908 }
909 // Set the block id with an offset
910 if( mdbImpl->tag_set_data( mMaterialSetTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
911 if( mdbImpl->tag_set_data( mGlobalIdTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
912 }
913 else // this is regular
914 {
915 // Allocate an array to read in connectivity data
916 readMeshIface->get_element_connect( this_it->numElements, verts_per_element, mb_type, this_it->startExoId,
917 this_it->startMBId, conn );
918
919 // Create a range for this sequence of elements
920 EntityHandle start_range, end_range;
921 start_range = ( *this_it ).startMBId;
922 end_range = start_range + ( *this_it ).numElements - 1;
923
924 Range new_range( start_range, end_range );
925 // Range<EntityHandle> new_range((*this_it).startMBId,
926 // (*this_it).startMBId + (*this_it).numElements - 1);
927
928 // Create a MBSet for this block and set the material tag
929
930 EntityHandle ms_handle;
931 if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, ms_handle ) != MB_SUCCESS )
932 return MB_FAILURE;
933
934 if( mdbImpl->add_entities( ms_handle, new_range ) != MB_SUCCESS ) return MB_FAILURE;
935
936 int mid_nodes[4];
937 CN::HasMidNodes( mb_type, verts_per_element, mid_nodes );
938 if( mdbImpl->tag_set_data( mHasMidNodesTag, &ms_handle, 1, mid_nodes ) != MB_SUCCESS ) return MB_FAILURE;
939
940 // Just a check because the following code won't work if this case fails
941 assert( sizeof( EntityHandle ) >= sizeof( int ) );
942
943 // tmp_ptr is of type int* and points at the same place as conn
944 int* tmp_ptr = reinterpret_cast< int* >( conn );
945
946 // Read the connectivity into that memory, this will take only part of the array
947 // 1/2 if sizeof(EntityHandle) == 64 bits.
948 count[0] = this_it->numElements;
949 count[1] = verts_per_element;
950 fail = nc_get_vara_int( ncFile, nc_var, start, count, tmp_ptr );
951 if( NC_NOERR != fail )
952 {
953 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connectivity" );
954 }
955
956 // Convert from exodus indices to vertex handles.
957 // Iterate backwards in case handles are larger than ints.
958 for( int i = number_nodes - 1; i >= 0; --i )
959 {
960 if( (unsigned)tmp_ptr[i] >= nodesInLoadedBlocks.size() )
961 {
962 MB_SET_ERR( MB_FAILURE, "Invalid node ID in block connectivity" );
963 }
964 nodesInLoadedBlocks[tmp_ptr[i]] = 1;
965 conn[i] = static_cast< EntityHandle >( tmp_ptr[i] ) + vertexOffset;
966 }
967
968 // Adjust connectivity order if necessary
969 const int* reorder = exodus_elem_order_map[mb_type][verts_per_element];
970 if( reorder ) ReadUtilIface::reorder( reorder, conn, this_it->numElements, verts_per_element );
971
972 readMeshIface->update_adjacencies( ( *this_it ).startMBId, ( *this_it ).numElements,
973 ExoIIUtil::VerticesPerElement[( *this_it ).elemType], conn );
974
975 if( result == -1 )
976 {
977 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: error getting element connectivity for block " << block_id );
978 }
979
980 // Set the block id with an offset
981 if( mdbImpl->tag_set_data( mMaterialSetTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
982 if( mdbImpl->tag_set_data( mGlobalIdTag, &ms_handle, 1, &block_id ) != MB_SUCCESS ) return MB_FAILURE;
983
984 if( file_id_tag )
985 {
986 Range range;
987 range.insert( this_it->startMBId, this_it->startMBId + this_it->numElements - 1 );
988 readMeshIface->assign_ids( *file_id_tag, range, this_it->startExoId );
989 }
990 } // end regular block (not polygon)
991 }
992
993 return MB_SUCCESS;
994 }
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().
|
private |
read the header from the ExoII file
Definition at line 393 of file ReadNCDF.cpp.
394 {
395 CPU_WORD_SIZE = sizeof( double ); // With ExodusII version 2, all floats
396 IO_WORD_SIZE = sizeof( double ); // should be changed to doubles
397
398 // NetCDF doesn't check its own limits on file read, so check
399 // them here so it doesn't corrupt memory any more than absolutely
400 // necessary.
401 int ndims;
402 int fail = nc_inq_ndims( ncFile, &ndims );
403 if( NC_NOERR != fail || ndims > NC_MAX_DIMS )
404 {
405 MB_SET_ERR( MB_FAILURE, "ReadNCDF: File contains " << ndims << " dims but NetCDF library supports only "
406 << (int)NC_MAX_DIMS );
407 }
408 int nvars;
409 fail = nc_inq_nvars( ncFile, &nvars );
410 if( nvars > NC_MAX_VARS )
411 {
412 MB_SET_ERR( MB_FAILURE, "ReadNCDF: File contains " << nvars << " vars but NetCDF library supports only "
413 << (int)NC_MAX_VARS );
414 }
415
416 // Get the attributes
417
418 // Get the word size, scalar value
419 nc_type att_type;
420 size_t att_len;
421 fail = nc_inq_att( ncFile, NC_GLOBAL, "floating_point_word_size", &att_type, &att_len );
422 if( NC_NOERR != fail )
423 {
424 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting floating_point_word_size attribute" );
425 }
426 if( att_type != NC_INT || att_len != 1 )
427 {
428 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Word size didn't have type int or size 1" );
429 }
430 fail = nc_get_att_int( ncFile, NC_GLOBAL, "floating_point_word_size", &IO_WORD_SIZE );
431 if( NC_NOERR != fail )
432 {
433 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Trouble getting word size" );
434 }
435
436 // Exodus version
437 fail = nc_inq_att( ncFile, NC_GLOBAL, "version", &att_type, &att_len );
438 if( NC_NOERR != fail )
439 {
440 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting version attribute" );
441 }
442 if( att_type != NC_FLOAT || att_len != 1 )
443 {
444 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Version didn't have type float or size 1" );
445 }
446 // float version = temp_att->as_float(0);
447
448 // Read in initial variables
449 int temp_dim;
450 GET_DIM( temp_dim, "num_dim", numberDimensions_loading );
451 GET_DIM( temp_dim, "num_nodes", numberNodes_loading );
452 GET_DIM( temp_dim, "num_elem", numberElements_loading );
453 GET_DIM( temp_dim, "num_el_blk", numberElementBlocks_loading );
454 GET_DIM( temp_dim, "num_fa_blk", numberFaceBlocks_loading ); // for polyhedra blocks
455 GET_DIM( temp_dim, "num_elem", numberElements_loading );
456 GET_DIM( temp_dim, "num_node_sets", numberNodeSets_loading );
457 GET_DIM( temp_dim, "num_side_sets", numberSideSets_loading );
458 GET_DIM( temp_dim, "len_string", max_str_length );
459 GET_DIM( temp_dim, "len_line", max_line_length );
460
461 // Title; why are we even bothering if we're not going to keep it???
462 fail = nc_inq_att( ncFile, NC_GLOBAL, "title", &att_type, &att_len );
463 if( NC_NOERR != fail )
464 {
465 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting title attribute" );
466 }
467 if( att_type != NC_CHAR )
468 {
469 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: title didn't have type char" );
470 }
471 char* title = new char[att_len + 1];
472 fail = nc_get_att_text( ncFile, NC_GLOBAL, "title", title );
473 if( NC_NOERR != fail )
474 {
475 delete[] title;
476 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: trouble getting title" );
477 }
478 delete[] title;
479
480 return MB_SUCCESS;
481 }
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().
|
private |
Definition at line 620 of file ReadNCDF.cpp.
621 {
622 // Get the ids of all the blocks of this file we're reading in
623 std::vector< int > fblock_ids( numberFaceBlocks_loading );
624 int nc_fblock_ids = -1;
625 GET_1D_INT_VAR( "fa_prop1", nc_fblock_ids, fblock_ids );
626 if( -1 == nc_fblock_ids )
627 {
628 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting fa_prop1 variable" );
629 }
630
631 int temp_dim;
632 std::vector< char > temp_string_storage( max_str_length + 1 );
633 char* temp_string = &temp_string_storage[0];
634 // int block_seq_id = 1;
635 int exodus_id = 1;
636
637 for( int i = 1; i <= numberFaceBlocks_loading; i++ )
638 {
639 int num_elements;
640
641 GET_DIMB( temp_dim, temp_string, "num_fa_in_blk%d", i, num_elements );
642
643 // Don't read element type string for now, since it's an attrib
644 // on the connectivity
645 // Get the entity type corresponding to this ExoII element type
646 // ExoIIElementType elem_type =
647 // ExoIIUtil::static_element_name_to_type(element_type_string);
648
649 // Tag each element block(mesh set) with enum for ElementType (SHELL, QUAD4, ....etc)
650 ReadFaceBlockData fblock_data;
651 fblock_data.faceBlockId = i;
652 fblock_data.startExoId = exodus_id;
653 fblock_data.numElements = num_elements;
654
655 faceBlocksLoading.push_back( fblock_data );
656 exodus_id += num_elements;
657 }
658 return MB_SUCCESS;
659 }
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().
|
private |
read in the global element ids
Definition at line 996 of file ReadNCDF.cpp.
997 {
998 // Read in the map from the exodus file
999 std::vector< int > ids( std::max( numberElements_loading, numberNodes_loading ) );
1000
1001 int varid = -1;
1002 GET_1D_INT_VAR( "elem_map", varid, ids );
1003 if( -1 != varid )
1004 {
1005 std::vector< ReadBlockData >::iterator iter;
1006 int id_pos = 0;
1007 for( iter = blocksLoading.begin(); iter != blocksLoading.end(); ++iter )
1008 {
1009 if( iter->reading_in )
1010 {
1011 if( iter->elemType != EXOII_POLYGON && iter->elemType != EXOII_POLYHEDRON )
1012 {
1013 if( iter->startMBId != 0 )
1014 {
1015 Range range( iter->startMBId, iter->startMBId + iter->numElements - 1 );
1016 ErrorCode error = mdbImpl->tag_set_data( mGlobalIdTag, range, &ids[id_pos] );
1017 if( error != MB_SUCCESS ) return error;
1018 id_pos += iter->numElements;
1019 }
1020 else
1021 return MB_FAILURE;
1022 }
1023 else // polygons or polyhedrons; use id from elements
1024 {
1025 for( std::vector< EntityHandle >::iterator eit = iter->polys.begin(); eit != iter->polys.end();
1026 eit++, id_pos++ )
1027 {
1028 EntityHandle peh = *eit;
1029 if( mdbImpl->tag_set_data( mGlobalIdTag, &peh, 1, &ids[id_pos] ) != MB_SUCCESS )
1030 return MB_FAILURE;
1031 }
1032 }
1033 }
1034 }
1035 }
1036
1037 // Read in node map next
1038 varid = -1;
1039 GET_1D_INT_VAR( "node_num_map", varid, ids );
1040 if( -1 != varid )
1041 {
1042 Range range( MB_START_ID + vertexOffset, MB_START_ID + vertexOffset + numberNodes_loading - 1 );
1043 ErrorCode error = mdbImpl->tag_set_data( mGlobalIdTag, range, &ids[0] );
1044 if( MB_SUCCESS != error ) MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem setting node global ids" );
1045 }
1046
1047 return MB_SUCCESS;
1048 }
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 the nodes
Definition at line 483 of file ReadNCDF.cpp.
484 {
485 // Read the nodes into memory
486
487 assert( 0 != ncFile );
488
489 // Create a sequence to hold the node coordinates
490 // Get the current number of entities and start at the next slot
491
492 EntityHandle node_handle = 0;
493 std::vector< double* > arrays;
494 readMeshIface->get_node_coords( 3, numberNodes_loading, MB_START_ID, node_handle, arrays );
495
496 vertexOffset = ID_FROM_HANDLE( node_handle ) - MB_START_ID;
497
498 // Read in the coordinates
499 int fail;
500 int coord = 0;
501 nc_inq_varid( ncFile, "coord", &coord );
502
503 // Single var for all coords
504 size_t start[2] = { 0, 0 }, count[2] = { 1, static_cast< size_t >( numberNodes_loading ) };
505 if( coord )
506 {
507 for( int d = 0; d < numberDimensions_loading; ++d )
508 {
509 start[0] = d;
510 fail = nc_get_vara_double( ncFile, coord, start, count, arrays[d] );
511 if( NC_NOERR != fail )
512 {
513 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting " << (char)( 'x' + d ) << " coord array" );
514 }
515 }
516 }
517 // Var for each coord
518 else
519 {
520 char varname[] = "coord ";
521 for( int d = 0; d < numberDimensions_loading; ++d )
522 {
523 varname[5] = 'x' + (char)d;
524 fail = nc_inq_varid( ncFile, varname, &coord );
525 if( NC_NOERR != fail )
526 {
527 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting " << (char)( 'x' + d ) << " coord variable" );
528 }
529
530 fail = nc_get_vara_double( ncFile, coord, start, &count[1], arrays[d] );
531 if( NC_NOERR != fail )
532 {
533 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting " << (char)( 'x' + d ) << " coord array" );
534 }
535 }
536 }
537
538 // Zero out any coord values that are in database but not in file
539 // (e.g. if MOAB has 3D coords but file is 2D then set Z coords to zero.)
540 for( unsigned d = numberDimensions_loading; d < arrays.size(); ++d )
541 std::fill( arrays[d], arrays[d] + numberNodes_loading, 0.0 );
542
543 if( file_id_tag )
544 {
545 Range nodes;
546 nodes.insert( node_handle, node_handle + numberNodes_loading - 1 );
547 readMeshIface->assign_ids( *file_id_tag, nodes, vertexOffset );
548 }
549
550 return MB_SUCCESS;
551 }
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().
|
private |
read the nodesets into meshsets
Definition at line 1050 of file ReadNCDF.cpp.
1051 {
1052 // Read in the nodesets for the model
1053
1054 if( 0 == numberNodeSets_loading ) return MB_SUCCESS;
1055 std::vector< int > id_array( numberNodeSets_loading );
1056
1057 // Read in the nodeset ids
1058 int nc_var;
1059 GET_1D_INT_VAR( "ns_prop1", nc_var, id_array );
1060 if( -1 == nc_var )
1061 {
1062 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting ns_prop1 variable" );
1063 }
1064
1065 // Use a vector of ints to read node handles
1066 std::vector< int > node_handles;
1067
1068 int i;
1069 std::vector< char > temp_string_storage( max_str_length + 1 );
1070 char* temp_string = &temp_string_storage[0];
1071 int temp_dim;
1072 for( i = 0; i < numberNodeSets_loading; i++ )
1073 {
1074 // Get nodeset parameters
1075 int number_nodes_in_set;
1076 int number_dist_factors_in_set;
1077
1078 GET_DIMB( temp_dim, temp_string, "num_nod_ns%d", i + 1, number_nodes_in_set );
1079 GET_DIMB( temp_dim, temp_string, "num_df_ns%d", i + 1, number_dist_factors_in_set );
1080
1081 // Need to new a vector to store dist. factors
1082 // This vector * gets stored as a tag on the sideset meshset
1083 std::vector< double > temp_dist_factor_vector( number_nodes_in_set );
1084 if( number_dist_factors_in_set != 0 )
1085 {
1086 INS_ID( temp_string, "dist_fact_ns%d", i + 1, max_str_length + 1 );
1087 GET_1D_DBL_VAR( temp_string, temp_dim, temp_dist_factor_vector );
1088 if( -1 == temp_dim )
1089 {
1090 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting dist fact variable" );
1091 }
1092 }
1093
1094 // Size new arrays and get ids and distribution factors
1095 if( node_handles.size() < (unsigned int)number_nodes_in_set )
1096 {
1097 node_handles.reserve( number_nodes_in_set );
1098 node_handles.resize( number_nodes_in_set );
1099 }
1100
1101 INS_ID( temp_string, "node_ns%d", i + 1, max_str_length + 1 );
1102 int temp_var = -1;
1103 GET_1D_INT_VAR( temp_string, temp_var, node_handles );
1104 if( -1 == temp_var )
1105 {
1106 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting nodeset node variable" );
1107 }
1108
1109 // Maybe there is already a nodesets meshset here we can append to
1110 Range child_meshsets;
1111 if( mdbImpl->get_entities_by_handle( 0, child_meshsets ) != MB_SUCCESS ) return MB_FAILURE;
1112
1113 child_meshsets = subtract( child_meshsets, initRange );
1114
1115 Range::iterator iter, end_iter;
1116 iter = child_meshsets.begin();
1117 end_iter = child_meshsets.end();
1118
1119 EntityHandle ns_handle = 0;
1120 for( ; iter != end_iter; ++iter )
1121 {
1122 int nodeset_id;
1123 if( mdbImpl->tag_get_data( mDirichletSetTag, &( *iter ), 1, &nodeset_id ) != MB_SUCCESS ) continue;
1124
1125 if( id_array[i] == nodeset_id )
1126 {
1127 // Found the meshset
1128 ns_handle = *iter;
1129 break;
1130 }
1131 }
1132
1133 std::vector< EntityHandle > nodes_of_nodeset;
1134 if( ns_handle )
1135 if( mdbImpl->get_entities_by_handle( ns_handle, nodes_of_nodeset, true ) != MB_SUCCESS ) return MB_FAILURE;
1136
1137 // Make these into entity handles
1138 // TODO: could we have read it into EntityHandle sized array in the first place?
1139 int j, temp;
1140 std::vector< EntityHandle > nodes;
1141 std::vector< double > dist_factor_vector;
1142 for( j = 0; j < number_nodes_in_set; j++ )
1143 {
1144 // See if this node is one we're currently reading in
1145 if( nodesInLoadedBlocks[node_handles[j]] == 1 )
1146 {
1147 // Make sure that it already isn't in a nodeset
1148 unsigned int node_id = CREATE_HANDLE( MBVERTEX, node_handles[j] + vertexOffset, temp );
1149 if( !ns_handle ||
1150 std::find( nodes_of_nodeset.begin(), nodes_of_nodeset.end(), node_id ) == nodes_of_nodeset.end() )
1151 {
1152 nodes.push_back( node_id );
1153
1154 if( number_dist_factors_in_set != 0 ) dist_factor_vector.push_back( temp_dist_factor_vector[j] );
1155 }
1156 }
1157 }
1158
1159 // No nodes to add
1160 if( nodes.empty() ) continue;
1161
1162 // If there was no meshset found --> create one
1163 if( ns_handle == 0 )
1164 {
1165 if( mdbImpl->create_meshset( MESHSET_ORDERED | MESHSET_TRACK_OWNER, ns_handle ) != MB_SUCCESS )
1166 return MB_FAILURE;
1167
1168 // Set a tag signifying dirichlet bc
1169 // TODO: create this tag another way
1170
1171 int nodeset_id = id_array[i];
1172 if( mdbImpl->tag_set_data( mDirichletSetTag, &ns_handle, 1, &nodeset_id ) != MB_SUCCESS ) return MB_FAILURE;
1173 if( mdbImpl->tag_set_data( mGlobalIdTag, &ns_handle, 1, &nodeset_id ) != MB_SUCCESS ) return MB_FAILURE;
1174
1175 if( !dist_factor_vector.empty() )
1176 {
1177 int size = dist_factor_vector.size();
1178 const void* data = &dist_factor_vector[0];
1179 if( mdbImpl->tag_set_by_ptr( mDistFactorTag, &ns_handle, 1, &data, &size ) != MB_SUCCESS )
1180 return MB_FAILURE;
1181 }
1182 }
1183 else if( !dist_factor_vector.empty() )
1184 {
1185 // Append dist factors to vector
1186 const void* ptr = 0;
1187 int size = 0;
1188 if( mdbImpl->tag_get_by_ptr( mDistFactorTag, &ns_handle, 1, &ptr, &size ) != MB_SUCCESS ) return MB_FAILURE;
1189
1190 const double* data = reinterpret_cast< const double* >( ptr );
1191 dist_factor_vector.reserve( dist_factor_vector.size() + size );
1192 std::copy( data, data + size, std::back_inserter( dist_factor_vector ) );
1193 size = dist_factor_vector.size();
1194 ptr = &dist_factor_vector[0];
1195 if( mdbImpl->tag_set_by_ptr( mDistFactorTag, &ns_handle, 1, &ptr, &size ) != MB_SUCCESS ) return MB_FAILURE;
1196 }
1197
1198 // Add the nodes to the meshset
1199 if( mdbImpl->add_entities( ns_handle, &nodes[0], nodes.size() ) != MB_SUCCESS ) return MB_FAILURE;
1200 }
1201
1202 return MB_SUCCESS;
1203 }
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().
|
private |
Definition at line 660 of file ReadNCDF.cpp.
661 {
662 int temp_dim;
663 std::vector< char > temp_string_storage( max_str_length + 1 );
664 char* temp_string = &temp_string_storage[0];
665 nodesInLoadedBlocks.resize( numberNodes_loading + 1 );
666 size_t start[1] = { 0 }, count[1]; // one dim arrays only here!
667 std::vector< ReadFaceBlockData >::iterator this_it;
668 this_it = faceBlocksLoading.begin();
669
670 int fblock_seq_id = 1;
671 std::vector< int > dims;
672 int nc_var, fail;
673
674 for( ; this_it != faceBlocksLoading.end(); ++this_it, fblock_seq_id++ )
675 {
676 int num_fa_in_blk;
677 GET_DIMB( temp_dim, temp_string, "num_fa_in_blk%d", fblock_seq_id, num_fa_in_blk );
678 int num_nod_per_fa;
679 GET_DIMB( temp_dim, temp_string, "num_nod_per_fa%d", fblock_seq_id, num_nod_per_fa );
680 // Get the ncdf connect variable and the element type
681 INS_ID( temp_string, "fbconn%d", fblock_seq_id, max_str_length + 1 );
682 GET_VAR( temp_string, nc_var, dims );
683 std::vector< int > fbconn;
684 fbconn.resize( num_nod_per_fa );
685 count[0] = num_nod_per_fa;
686
687 fail = nc_get_vara_int( ncFile, nc_var, start, count, &fbconn[0] );
688 if( NC_NOERR != fail )
689 {
690 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting face polyhedra connectivity " );
691 }
692 std::vector< int > fbepecnt;
693 fbepecnt.resize( num_fa_in_blk );
694 INS_ID( temp_string, "fbepecnt%d", fblock_seq_id, max_str_length + 1 );
695 GET_VAR( temp_string, nc_var, dims );
696 count[0] = num_fa_in_blk;
697 fail = nc_get_vara_int( ncFile, nc_var, start, count, &fbepecnt[0] );
698 if( NC_NOERR != fail )
699 {
700 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting face polyhedra connectivity " );
701 }
702 // now we can create some polygons
703 std::vector< EntityHandle > polyConn;
704 int ix = 0;
705 for( int i = 0; i < num_fa_in_blk; i++ )
706 {
707 polyConn.resize( fbepecnt[i] );
708 for( int j = 0; j < fbepecnt[i]; j++ )
709 {
710 nodesInLoadedBlocks[fbconn[ix]] = 1;
711 polyConn[j] = vertexOffset + fbconn[ix++];
712 }
713 EntityHandle newp;
714 /*
715 * ErrorCode create_element(const EntityType type,
716 const EntityHandle *connectivity,
717 const int num_vertices,
718 EntityHandle &element_handle)
719 */
720 if( mdbImpl->create_element( MBPOLYGON, &polyConn[0], fbepecnt[i], newp ) != MB_SUCCESS ) return MB_FAILURE;
721
722 // add the element in order
723 polyfaces.push_back( newp );
724 }
725 }
726 return MB_SUCCESS;
727 }
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().
|
private |
Definition at line 1751 of file ReadNCDF.cpp.
1752 {
1753 // Inquire on the genesis file to find the number of qa records
1754
1755 int number_records = 0;
1756
1757 int temp_dim;
1758 GET_DIM( temp_dim, "num_qa_rec", number_records );
1759 std::vector< char > data( max_str_length + 1 );
1760
1761 for( int i = 0; i < number_records; i++ )
1762 {
1763 for( int j = 0; j < 4; j++ )
1764 {
1765 data[max_str_length] = '\0';
1766 if( read_qa_string( &data[0], i, j ) != MB_SUCCESS ) return MB_FAILURE;
1767
1768 qa_record_list.push_back( &data[0] );
1769 }
1770 }
1771
1772 return MB_SUCCESS;
1773 }
References GET_DIM, max_str_length, MB_SUCCESS, and read_qa_string().
Referenced by read_qa_records().
|
private |
Definition at line 1725 of file ReadNCDF.cpp.
1726 {
1727 std::vector< std::string > qa_records;
1728 read_qa_information( qa_records );
1729
1730 std::vector< char > tag_data;
1731 for( std::vector< std::string >::iterator i = qa_records.begin(); i != qa_records.end(); ++i )
1732 {
1733 std::copy( i->begin(), i->end(), std::back_inserter( tag_data ) );
1734 tag_data.push_back( '\0' );
1735 }
1736
1737 // If there were qa_records...tag them to the mCurrentMeshHandle
1738 if( !tag_data.empty() )
1739 {
1740 const void* ptr = &tag_data[0];
1741 int size = tag_data.size();
1742 if( mdbImpl->tag_set_by_ptr( mQaRecordTag, &file_set, 1, &ptr, &size ) != MB_SUCCESS )
1743 {
1744 return MB_FAILURE;
1745 }
1746 }
1747
1748 return MB_SUCCESS;
1749 }
References MB_SUCCESS, mdbImpl, mQaRecordTag, read_qa_information(), size, and moab::Interface::tag_set_by_ptr().
Referenced by load_file().
|
private |
Definition at line 1775 of file ReadNCDF.cpp.
1776 {
1777 std::vector< int > dims;
1778 int temp_var = -1;
1779 GET_VAR( "qa_records", temp_var, dims );
1780 if( -1 == temp_var )
1781 {
1782 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting qa record variable" );
1783 return MB_FAILURE;
1784 }
1785 size_t count[3], start[3];
1786 start[0] = record_number;
1787 start[1] = record_position;
1788 start[2] = 0;
1789 count[0] = 1;
1790 count[1] = 1;
1791 count[2] = max_str_length;
1792 int fail = nc_get_vara_text( ncFile, temp_var, start, count, temp_string );
1793 if( NC_NOERR != fail )
1794 {
1795 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem setting current record number variable position" );
1796 }
1797 // Get the variable id in the exodus file
1798
1799 return MB_SUCCESS;
1800 }
References moab::fail(), GET_VAR, max_str_length, MB_SET_ERR, MB_SUCCESS, and ncFile.
Referenced by read_qa_information().
|
private |
read the sidesets (does nothing for now)
Definition at line 1205 of file ReadNCDF.cpp.
1206 {
1207 // Uhhh if you read the same file (previously_read==true) then blow away the
1208 // sidesets pertaining to this file? How do you do that? If you read a
1209 // new file make sure all the offsets are correct.
1210
1211 // If not loading any sidesets -- exit
1212 if( 0 == numberSideSets_loading ) return MB_SUCCESS;
1213
1214 // Read in the sideset ids
1215 std::vector< int > id_array( numberSideSets_loading );
1216 int temp_var = -1;
1217 GET_1D_INT_VAR( "ss_prop1", temp_var, id_array );
1218 if( -1 == temp_var )
1219 {
1220 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting ss_prop1 variable" );
1221 }
1222
1223 // Create a side set for each one
1224 int number_sides_in_set;
1225 int number_dist_factors_in_set;
1226
1227 // Maybe there is already a sidesets meshset here we can append to
1228 Range child_meshsets;
1229 if( mdbImpl->get_entities_by_type( 0, MBENTITYSET, child_meshsets ) != MB_SUCCESS ) return MB_FAILURE;
1230
1231 child_meshsets = subtract( child_meshsets, initRange );
1232
1233 Range::iterator iter, end_iter;
1234
1235 int i;
1236 std::vector< char > temp_string_storage( max_str_length + 1 );
1237 char* temp_string = &temp_string_storage[0];
1238 int temp_dim;
1239 for( i = 0; i < numberSideSets_loading; i++ )
1240 {
1241 // Get sideset parameters
1242 GET_DIMB( temp_dim, temp_string, "num_side_ss%d", i + 1, number_sides_in_set );
1243 GET_DIMB( temp_dim, temp_string, "num_df_ss%d", i + 1, number_dist_factors_in_set );
1244
1245 // Size new arrays and get element and side lists
1246 std::vector< int > side_list( number_sides_in_set );
1247 std::vector< int > element_list( number_sides_in_set );
1248 INS_ID( temp_string, "side_ss%d", i + 1, max_str_length + 1 );
1249 temp_var = -1;
1250 GET_1D_INT_VAR( temp_string, temp_var, side_list );
1251 if( -1 == temp_var )
1252 {
1253 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting sideset side variable" );
1254 }
1255
1256 INS_ID( temp_string, "elem_ss%d", i + 1, max_str_length + 1 );
1257 temp_var = -1;
1258 GET_1D_INT_VAR( temp_string, temp_var, element_list );
1259 if( -1 == temp_var )
1260 {
1261 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting sideset elem variable" );
1262 }
1263
1264 std::vector< double > temp_dist_factor_vector;
1265 std::vector< EntityHandle > entities_to_add, reverse_entities;
1266 // Create the sideset entities
1267 if( create_ss_elements( &element_list[0], &side_list[0], number_sides_in_set, number_dist_factors_in_set,
1268 entities_to_add, reverse_entities, temp_dist_factor_vector, i + 1 ) != MB_SUCCESS )
1269 return MB_FAILURE;
1270
1271 // If there are elements to add
1272 if( !entities_to_add.empty() || !reverse_entities.empty() )
1273 {
1274 iter = child_meshsets.begin();
1275 end_iter = child_meshsets.end();
1276
1277 EntityHandle ss_handle = 0;
1278 for( ; iter != end_iter; ++iter )
1279 {
1280 int sideset_id;
1281 if( mdbImpl->tag_get_data( mNeumannSetTag, &( *iter ), 1, &sideset_id ) != MB_SUCCESS ) continue;
1282
1283 if( id_array[i] == sideset_id )
1284 {
1285 // Found the meshset
1286 ss_handle = *iter;
1287 break;
1288 }
1289 }
1290
1291 // If we didn't find a sideset already
1292 if( ss_handle == 0 )
1293 {
1294 if( mdbImpl->create_meshset( MESHSET_ORDERED | MESHSET_TRACK_OWNER, ss_handle ) != MB_SUCCESS )
1295 return MB_FAILURE;
1296
1297 if( ss_handle == 0 ) return MB_FAILURE;
1298
1299 int sideset_id = id_array[i];
1300 if( mdbImpl->tag_set_data( mNeumannSetTag, &ss_handle, 1, &sideset_id ) != MB_SUCCESS )
1301 return MB_FAILURE;
1302 if( mdbImpl->tag_set_data( mGlobalIdTag, &ss_handle, 1, &sideset_id ) != MB_SUCCESS ) return MB_FAILURE;
1303
1304 if( !reverse_entities.empty() )
1305 {
1306 // Also make a reverse set to put in this set
1307 EntityHandle reverse_set;
1308 if( mdbImpl->create_meshset( MESHSET_SET | MESHSET_TRACK_OWNER, reverse_set ) != MB_SUCCESS )
1309 return MB_FAILURE;
1310
1311 // Add the reverse set to the sideset set and the entities to the reverse set
1312 ErrorCode result = mdbImpl->add_entities( ss_handle, &reverse_set, 1 );
1313 if( MB_SUCCESS != result ) return result;
1314
1315 result = mdbImpl->add_entities( reverse_set, &reverse_entities[0], reverse_entities.size() );
1316 if( MB_SUCCESS != result ) return result;
1317
1318 // Set the reverse tag
1319 Tag sense_tag;
1320 int dum_sense = 0;
1321 result = mdbImpl->tag_get_handle( "NEUSET_SENSE", 1, MB_TYPE_INTEGER, sense_tag,
1322 MB_TAG_SPARSE | MB_TAG_CREAT, &dum_sense );
1323 if( result != MB_SUCCESS ) return result;
1324 dum_sense = -1;
1325 result = mdbImpl->tag_set_data( sense_tag, &reverse_set, 1, &dum_sense );
1326 if( result != MB_SUCCESS ) return result;
1327 }
1328 }
1329
1330 if( mdbImpl->add_entities( ss_handle, ( entities_to_add.empty() ) ? NULL : &entities_to_add[0],
1331 entities_to_add.size() ) != MB_SUCCESS )
1332 return MB_FAILURE;
1333
1334 // Distribution factor stuff
1335 if( number_dist_factors_in_set )
1336 {
1337 // If this sideset does not already has a distribution factor array...set one
1338 const void* ptr = 0;
1339 int size = 0;
1340 if( MB_SUCCESS == mdbImpl->tag_get_by_ptr( mDistFactorTag, &ss_handle, 1, &ptr, &size ) )
1341 {
1342 const double* data = reinterpret_cast< const double* >( ptr );
1343 std::copy( data, data + size, std::back_inserter( temp_dist_factor_vector ) );
1344 }
1345
1346 ptr = &temp_dist_factor_vector[0];
1347 size = temp_dist_factor_vector.size();
1348 if( mdbImpl->tag_set_by_ptr( mDistFactorTag, &ss_handle, 1, &ptr, &size ) != MB_SUCCESS )
1349 return MB_FAILURE;
1350 }
1351 }
1352 }
1353
1354 return MB_SUCCESS;
1355 }
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().
|
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.
| 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.
219 {
220 if( subset_list )
221 {
222 MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "ExodusII reader supports subset read only by material ID" );
223 }
224
225 // Open netcdf/exodus file
226 int fail = nc_open( file_name, 0, &ncFile );
227 if( NC_NOWRITE != fail )
228 {
229 MB_SET_ERR( MB_FILE_DOES_NOT_EXIST, "ReadNCDF:: problem opening Netcdf/Exodus II file " << file_name );
230 }
231
232 // 1. Read the header
233 ErrorCode rval = read_exodus_header();
234 if( MB_FAILURE == rval ) return rval;
235
236 int count = 0;
237 const char* prop;
238 const char* blocks = "eb_prop1";
239 const char* nodesets = "ns_prop1";
240 const char* sidesets = "ss_prop1";
241
242 if( !strcmp( tag_name, MATERIAL_SET_TAG_NAME ) )
243 {
244 count = numberElementBlocks_loading;
245 prop = blocks;
246 }
247 else if( !strcmp( tag_name, DIRICHLET_SET_TAG_NAME ) )
248 {
249 count = numberNodeSets_loading;
250 prop = nodesets;
251 }
252 else if( !strcmp( tag_name, NEUMANN_SET_TAG_NAME ) )
253 {
254 count = numberSideSets_loading;
255 prop = sidesets;
256 }
257 else
258 {
259 ncFile = 0;
260 return MB_TAG_NOT_FOUND;
261 }
262
263 if( count )
264 {
265 int nc_var = -1;
266 GET_1D_INT_VAR( prop, nc_var, id_array );
267 if( !nc_var )
268 {
269 MB_SET_ERR( MB_FAILURE, "Problem getting prop variable" );
270 }
271 }
272
273 // Close the file
274 fail = nc_close( ncFile );
275 if( NC_NOERR != fail )
276 {
277 MB_SET_ERR( MB_FAILURE, "Trouble closing file" );
278 }
279
280 ncFile = 0;
281 return MB_SUCCESS;
282 }
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().
|
private |
Definition at line 189 of file ReadNCDF.cpp.
190 { 191 mCurrentMeshHandle = 0; 192 vertexOffset = 0; 193 194 numberNodes_loading = 0; 195 numberElements_loading = 0; 196 numberElementBlocks_loading = 0; 197 numberFaceBlocks_loading = 0; 198 numberNodeSets_loading = 0; 199 numberSideSets_loading = 0; 200 numberDimensions_loading = 0; 201 202 if( !blocksLoading.empty() ) blocksLoading.clear(); 203 204 if( !nodesInLoadedBlocks.empty() ) nodesInLoadedBlocks.clear(); 205 206 if( !faceBlocksLoading.empty() ) faceBlocksLoading.clear(); 207 }
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().
| void moab::ReadNCDF::tokenize | ( | const std::string & | str, |
| std::vector< std::string > & | tokens, | ||
| const char * | delimiters | ||
| ) |
Definition at line 2445 of file ReadNCDF.cpp.
2446 {
2447 std::string::size_type last = str.find_first_not_of( delimiters, 0 );
2448 std::string::size_type pos = str.find_first_of( delimiters, last );
2449 while( std::string::npos != pos && std::string::npos != last )
2450 {
2451 tokens.push_back( str.substr( last, pos - last ) );
2452 last = str.find_first_not_of( delimiters, pos );
2453 pos = str.find_first_of( delimiters, last );
2454 if( std::string::npos == pos ) pos = str.size();
2455 }
2456 }
Referenced by 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.
1807 {
1808 // Function : updating current database from new exodus_file.
1809 // Creator: Jane Hu
1810 // opts is currently designed as following
1811 // tdata = <var_name>[, time][,op][,destination]
1812 // where var_name show the tag name to be updated, this version just takes
1813 // coord.
1814 // time is the optional, and it gives time step of each of the mesh
1815 // info in exodus file. It start from 1.
1816 // op is the operation that is going to be performed on the var_name info.
1817 // currently support 'set'
1818 // destination shows where to store the updated info, currently assume it is
1819 // stored in the same database by replacing the old info if there's no input
1820 // for destination, or the destination data is given in exodus format and just
1821 // need to update the coordinates.
1822 // Assumptions:
1823 // 1. Assume the num_el_blk's in both DB and update exodus file are the same.
1824 // 2. Assume num_el_in_blk1...num_el_in_blk(num_el_blk) numbers are matching, may in
1825 // different order. example: num_el_in_blk11 = num_el_in_blk22 && num_el_in_blk12 =
1826 // num_el_in_blk21.
1827 // 3. In exodus file, get node_num_map
1828 // 4. loop through the node_num_map, use it to find the node in the cub file.
1829 // 5. Replace coord[0][n] with coordx[m] + vals_nod_var1(time_step, m) for all directions for
1830 // matching nodes.
1831 // Test: try to get hexes
1832
1833 // *******************************************************************
1834 // Move nodes to their deformed locations.
1835 // *******************************************************************
1836 ErrorCode rval;
1837 std::string s;
1838 rval = opts.get_str_option( "tdata", s );
1839 if( MB_SUCCESS != rval )
1840 {
1841 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem reading file options" );
1842 }
1843 std::vector< std::string > tokens;
1844 tokenize( s, tokens, "," );
1845
1846 // 1. Check for time step to find the match time
1847 int time_step = 1;
1848 if( tokens.size() > 1 && !tokens[1].empty() )
1849 {
1850 const char* time_s = tokens[1].c_str();
1851 char* endptr;
1852 long int pval = strtol( time_s, &endptr, 0 );
1853 std::string end = endptr;
1854 if( !end.empty() ) // Syntax error
1855 return MB_TYPE_OUT_OF_RANGE;
1856
1857 // Check for overflow (parsing long int, returning int)
1858 time_step = pval;
1859 if( pval != (long int)time_step ) return MB_TYPE_OUT_OF_RANGE;
1860 if( time_step <= 0 ) return MB_TYPE_OUT_OF_RANGE;
1861 }
1862
1863 // 2. Check for the operations, currently support set.
1864 const char* op;
1865 if( tokens.size() < 3 || ( tokens[2] != "set" && tokens[2] != "add" ) )
1866 {
1867 MB_SET_ERR( MB_TYPE_OUT_OF_RANGE, "ReadNCDF: invalid operation specified for update" );
1868 }
1869 op = tokens[2].c_str();
1870
1871 // Open netcdf/exodus file
1872 ncFile = 0;
1873 int fail = nc_open( exodus_file_name, 0, &ncFile );
1874 if( !ncFile )
1875 {
1876 MB_SET_ERR( MB_FILE_DOES_NOT_EXIST, "ReadNCDF:: problem opening Netcdf/Exodus II file " << exodus_file_name );
1877 }
1878
1879 rval = read_exodus_header();MB_CHK_ERR( rval );
1880
1881 // Check to make sure that the requested time step exists
1882 int ncdim = -1;
1883 int max_time_steps;
1884 GET_DIM( ncdim, "time_step", max_time_steps );
1885 if( -1 == ncdim )
1886 {
1887 std::cout << "ReadNCDF: could not get number of time steps" << std::endl;
1888 return MB_FAILURE;
1889 }
1890 std::cout << " Maximum time step=" << max_time_steps << std::endl;
1891 if( max_time_steps < time_step )
1892 {
1893 std::cout << "ReadNCDF: time step is greater than max_time_steps" << std::endl;
1894 return MB_FAILURE;
1895 }
1896
1897 // Get the time
1898 std::vector< double > times( max_time_steps );
1899 int nc_var = -1;
1900 GET_1D_DBL_VAR( "time_whole", nc_var, times );
1901 if( -1 == nc_var )
1902 {
1903 std::cout << "ReadNCDF: unable to get time variable" << std::endl;
1904 }
1905 else
1906 {
1907 std::cout << " Step " << time_step << " is at " << times[time_step - 1] << " seconds" << std::endl;
1908 }
1909
1910 // Read in the node_num_map.
1911 std::vector< int > ptr( numberNodes_loading );
1912
1913 int varid = -1;
1914 GET_1D_INT_VAR( "node_num_map", varid, ptr );
1915 if( -1 == varid )
1916 {
1917 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting node number map data" );
1918 }
1919
1920 // Read in the deformations.
1921 std::vector< std::vector< double > > deformed_arrays( 3 );
1922 std::vector< std::vector< double > > orig_coords( 3 );
1923 deformed_arrays[0].reserve( numberNodes_loading );
1924 deformed_arrays[1].reserve( numberNodes_loading );
1925 deformed_arrays[2].reserve( numberNodes_loading );
1926 orig_coords[0].reserve( numberNodes_loading );
1927 orig_coords[1].reserve( numberNodes_loading );
1928 orig_coords[2].reserve( numberNodes_loading );
1929 size_t start[2] = { static_cast< size_t >( time_step - 1 ), 0 },
1930 count[2] = { 1, static_cast< size_t >( numberNodes_loading ) };
1931 std::vector< int > dims;
1932 int coordx = -1, coordy = -1, coordz = -1;
1933 GET_VAR( "vals_nod_var1", coordx, dims );
1934 GET_VAR( "vals_nod_var2", coordy, dims );
1935 if( numberDimensions_loading == 3 ) GET_VAR( "vals_nod_var3", coordz, dims );
1936 if( -1 == coordx || -1 == coordy || ( numberDimensions_loading == 3 && -1 == coordz ) )
1937 {
1938 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting coords variable" );
1939 }
1940
1941 fail = nc_get_vara_double( ncFile, coordx, start, count, &deformed_arrays[0][0] );
1942 if( NC_NOERR != fail )
1943 {
1944 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting x deformation array" );
1945 }
1946
1947 fail = nc_get_vara_double( ncFile, coordy, start, count, &deformed_arrays[1][0] );
1948 if( NC_NOERR != fail )
1949 {
1950 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting y deformation array" );
1951 }
1952 if( numberDimensions_loading == 3 )
1953 {
1954 fail = nc_get_vara_double( ncFile, coordz, start, count, &deformed_arrays[2][0] );
1955 if( NC_NOERR != fail )
1956 {
1957 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting z deformation array" );
1958 }
1959 }
1960
1961 int coord1 = -1, coord2 = -1, coord3 = -1;
1962 GET_1D_DBL_VAR( "coordx", coord1, orig_coords[0] );
1963 if( -1 == coord1 )
1964 {
1965 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting x coord array" );
1966 }
1967 GET_1D_DBL_VAR( "coordy", coord2, orig_coords[1] );
1968 if( -1 == coord2 )
1969 {
1970 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting y coord array" );
1971 }
1972 if( numberDimensions_loading == 3 )
1973 {
1974 GET_1D_DBL_VAR( "coordz", coord3, orig_coords[2] );
1975 if( -1 == coord3 )
1976 {
1977 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting z coord array" );
1978 }
1979 }
1980
1981 // b. Deal with DB file : get node info. according to node_num_map.
1982 if( tokens[0] != "coord" && tokens[0] != "COORD" ) return MB_NOT_IMPLEMENTED;
1983
1984 if( strcmp( op, "set" ) != 0 && strcmp( op, " set" ) != 0 ) return MB_NOT_IMPLEMENTED;
1985
1986 // Two methods of matching nodes (id vs. proximity)
1987 const bool match_node_ids = true;
1988
1989 // Get nodes in cubit file
1990 Range cub_verts;
1991 rval = mdbImpl->get_entities_by_type( cub_file_set, MBVERTEX, cub_verts );MB_CHK_ERR( rval );
1992 std::cout << " cub_file_set contains " << cub_verts.size() << " nodes." << std::endl;
1993
1994 // Some accounting
1995 std::cout << " exodus file contains " << numberNodes_loading << " nodes." << std::endl;
1996 double max_magnitude = 0;
1997 double average_magnitude = 0;
1998 int found = 0;
1999 int lost = 0;
2000 std::map< int, EntityHandle > cub_verts_id_map;
2001 AdaptiveKDTree kdtree( mdbImpl );
2002 EntityHandle root;
2003
2004 // Should not use cub verts unless they have been matched. Place in a map
2005 // for fast handle_by_id lookup.
2006 std::map< int, EntityHandle > matched_cub_vert_id_map;
2007
2008 // Place cub verts in a map for searching by id
2009 if( match_node_ids )
2010 {
2011 std::vector< int > cub_ids( cub_verts.size() );
2012 rval = mdbImpl->tag_get_data( mGlobalIdTag, cub_verts, &cub_ids[0] );MB_CHK_ERR( rval );
2013 for( unsigned i = 0; i != cub_verts.size(); ++i )
2014 {
2015 cub_verts_id_map.insert( std::pair< int, EntityHandle >( cub_ids[i], cub_verts[i] ) );
2016 }
2017
2018 // Place cub verts in a kdtree for searching by proximity
2019 }
2020 else
2021 {
2022 FileOptions tree_opts( "MAX_PER_LEAF=1;SPLITS_PER_DIR=1;CANDIDATE_PLANE_SET=0" );
2023 rval = kdtree.build_tree( cub_verts, &root, &tree_opts );MB_CHK_ERR( rval );
2024 AdaptiveKDTreeIter tree_iter;
2025 rval = kdtree.get_tree_iterator( root, tree_iter );MB_CHK_ERR( rval );
2026 }
2027
2028 // For each exo vert, find the matching cub vert
2029 for( int i = 0; i < numberNodes_loading; ++i )
2030 {
2031 int exo_id = ptr[i];
2032 CartVect exo_coords( orig_coords[0][i], orig_coords[1][i], orig_coords[2][i] );
2033 EntityHandle cub_vert = 0;
2034 bool found_match = false;
2035
2036 // By id
2037 if( match_node_ids )
2038 {
2039 std::map< int, EntityHandle >::iterator i_iter;
2040 i_iter = cub_verts_id_map.find( exo_id );
2041 if( i_iter != cub_verts_id_map.end() )
2042 {
2043 found_match = true;
2044 cub_vert = i_iter->second;
2045 }
2046
2047 // By proximity
2048 }
2049 else
2050 {
2051 // The MAX_NODE_DIST is the farthest distance to search for a node.
2052 // For the 1/12th symmetry 85 pin model, the max node dist could not be less
2053 // than 1e-1 (March 26, 2010).
2054 const double MAX_NODE_DIST = 1e-1;
2055
2056 std::vector< EntityHandle > leaves;
2057 double min_dist = MAX_NODE_DIST;
2058 rval = kdtree.distance_search( exo_coords.array(), MAX_NODE_DIST, leaves );MB_CHK_ERR( rval );
2059 for( std::vector< EntityHandle >::const_iterator j = leaves.begin(); j != leaves.end(); ++j )
2060 {
2061 std::vector< EntityHandle > leaf_verts;
2062 rval = mdbImpl->get_entities_by_type( *j, MBVERTEX, leaf_verts );MB_CHK_ERR( rval );
2063 for( std::vector< EntityHandle >::const_iterator k = leaf_verts.begin(); k != leaf_verts.end(); ++k )
2064 {
2065 CartVect orig_cub_coords, difference;
2066 rval = mdbImpl->get_coords( &( *k ), 1, orig_cub_coords.array() );MB_CHK_ERR( rval );
2067 difference = orig_cub_coords - exo_coords;
2068 double dist = difference.length();
2069 if( dist < min_dist )
2070 {
2071 min_dist = dist;
2072 cub_vert = *k;
2073 }
2074 }
2075 }
2076 if( 0 != cub_vert ) found_match = true;
2077 }
2078
2079 // If a match is found, update it with the deformed coords from the exo file.
2080 if( found_match )
2081 {
2082 CartVect updated_exo_coords;
2083 matched_cub_vert_id_map.insert( std::pair< int, EntityHandle >( exo_id, cub_vert ) );
2084 updated_exo_coords[0] = orig_coords[0][i] + deformed_arrays[0][i];
2085 updated_exo_coords[1] = orig_coords[1][i] + deformed_arrays[1][i];
2086
2087 if( numberDimensions_loading == 3 ) updated_exo_coords[2] = orig_coords[2][i] + deformed_arrays[2][i];
2088 rval = mdbImpl->set_coords( &cub_vert, 1, updated_exo_coords.array() );MB_CHK_ERR( rval );
2089 ++found;
2090
2091 double magnitude =
2092 sqrt( deformed_arrays[0][i] * deformed_arrays[0][i] + deformed_arrays[1][i] * deformed_arrays[1][i] +
2093 deformed_arrays[2][i] * deformed_arrays[2][i] );
2094 if( magnitude > max_magnitude ) max_magnitude = magnitude;
2095 average_magnitude += magnitude;
2096 }
2097 else
2098 {
2099 ++lost;
2100 std::cout << "cannot match exo node " << exo_id << " " << exo_coords << std::endl;
2101 }
2102 }
2103
2104 // Exo verts that could not be matched have already been printed. Now print the
2105 // cub verts that could not be matched.
2106 if( matched_cub_vert_id_map.size() < cub_verts.size() )
2107 {
2108 Range unmatched_cub_verts = cub_verts;
2109 for( std::map< int, EntityHandle >::const_iterator i = matched_cub_vert_id_map.begin();
2110 i != matched_cub_vert_id_map.end(); ++i )
2111 {
2112 unmatched_cub_verts.erase( i->second );
2113 }
2114
2115 for( Range::const_iterator i = unmatched_cub_verts.begin(); i != unmatched_cub_verts.end(); ++i )
2116 {
2117 int cub_id;
2118 rval = mdbImpl->tag_get_data( mGlobalIdTag, &( *i ), 1, &cub_id );MB_CHK_ERR( rval );
2119
2120 CartVect cub_coords;
2121 rval = mdbImpl->get_coords( &( *i ), 1, cub_coords.array() );MB_CHK_ERR( rval );
2122 std::cout << "cannot match cub node " << cub_id << " " << cub_coords << std::endl;
2123 }
2124
2125 std::cout << " " << unmatched_cub_verts.size() << " nodes from the cub file could not be matched."
2126 << std::endl;
2127 }
2128
2129 // Summarize statistics
2130 std::cout << " " << found << " nodes from the exodus file were matched in the cub_file_set ";
2131 if( match_node_ids )
2132 {
2133 std::cout << "by id." << std::endl;
2134 }
2135 else
2136 {
2137 std::cout << "by proximity." << std::endl;
2138 }
2139
2140 // Fail if all of the nodes could not be matched.
2141 if( 0 != lost )
2142 {
2143 std::cout << "Error: " << lost << " nodes from the exodus file could not be matched." << std::endl;
2144 // return MB_FAILURE;
2145 }
2146 std::cout << " maximum node displacement magnitude: " << max_magnitude << " cm" << std::endl;
2147 std::cout << " average node displacement magnitude: " << average_magnitude / found << " cm" << std::endl;
2148
2149 // *******************************************************************
2150 // Remove dead elements from the MOAB instance.
2151 // *******************************************************************
2152
2153 // How many element variables are in the file?
2154 int n_elem_var;
2155 GET_DIM( ncdim, "num_elem_var", n_elem_var );
2156
2157 // Get element variable names
2158 varid = -1;
2159 int cstatus = nc_inq_varid( ncFile, "name_elem_var", &varid );
2160 std::vector< char > names_memory( n_elem_var * max_str_length );
2161 std::vector< char* > names( n_elem_var );
2162 for( int i = 0; i < n_elem_var; ++i )
2163 names[i] = &names_memory[i * max_str_length];
2164 if( cstatus != NC_NOERR || varid == -1 )
2165 {
2166 std::cout << "ReadNCDF: name_elem_var does not exist" << std::endl;
2167 return MB_FAILURE;
2168 }
2169 int status = nc_get_var_text( ncFile, varid, &names_memory[0] );
2170 if( NC_NOERR != status )
2171 {
2172 MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem getting element variable names" );
2173 }
2174
2175 // Is one of the element variable names "death_status"? If so, get its index
2176 // in the element variable array.
2177 int death_index;
2178 bool found_death_index = false;
2179 for( int i = 0; i < n_elem_var; ++i )
2180 {
2181 std::string temp( names[i] );
2182 std::string::size_type pos0 = temp.find( "death_status" );
2183 std::string::size_type pos1 = temp.find( "Death_Status" );
2184 std::string::size_type pos2 = temp.find( "DEATH_STATUS" );
2185 if( std::string::npos != pos0 || std::string::npos != pos1 || std::string::npos != pos2 )
2186 {
2187 found_death_index = true;
2188 death_index = i + 1; // NetCDF variables start with 1
2189 break;
2190 }
2191 }
2192 if( !found_death_index )
2193 {
2194 MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem getting index of death_status variable" );
2195 }
2196
2197 // The exodus header has already been read. This contains the number of element
2198 // blocks.
2199
2200 // Dead elements are listed by block. Read the block headers to determine how
2201 // many elements are in each block.
2202 rval = read_block_headers( blocks_to_load, num_blocks );
2203 if( MB_FAILURE == rval )
2204 {
2205 MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem reading block headers" );
2206 }
2207
2208 // Dead elements from the Exodus file can be located in the cub_file_set by id
2209 // or by connectivity. Currently, finding elements by id requires careful book
2210 // keeping when constructing the model in Cubit. To avoid this, one can match
2211 // elements by connectivity instead.
2212 const bool match_elems_by_connectivity = true;
2213
2214 // Get the element id map. The ids in the map are from the elements in the blocks.
2215 // elem_num_map(blk1 elem ids, blk2 elem ids, blk3 elem ids, ...)
2216 std::vector< int > elem_ids( numberNodes_loading );
2217 if( !match_elems_by_connectivity )
2218 {
2219 GET_1D_INT_VAR( "elem_num_map", varid, elem_ids );
2220 if( -1 == varid )
2221 {
2222 MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem getting element number map data" );
2223 }
2224 }
2225
2226 // For each block
2227 int first_elem_id_in_block = 0;
2228 int block_count = 1; // NetCDF variables start with 1
2229 int total_elems = 0;
2230 int total_dead_elems = 0;
2231 for( std::vector< ReadBlockData >::iterator i = blocksLoading.begin(); i != blocksLoading.end(); ++i )
2232 {
2233 // Get the ncdf connect variable
2234 std::string temp_string( "connect" );
2235 std::stringstream temp_ss;
2236 temp_ss << block_count;
2237 temp_string += temp_ss.str();
2238 temp_string += "\0";
2239 std::vector< int > ldims;
2240 GET_VAR( temp_string.c_str(), nc_var, ldims );
2241 if( !nc_var )
2242 {
2243 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connect variable" );
2244 }
2245 // The element type is an attribute of the connectivity variable
2246 nc_type att_type;
2247 size_t att_len;
2248 fail = nc_inq_att( ncFile, nc_var, "elem_type", &att_type, &att_len );
2249 if( NC_NOERR != fail )
2250 {
2251 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type attribute" );
2252 }
2253 std::vector< char > dum_str( att_len + 1 );
2254 fail = nc_get_att_text( ncFile, nc_var, "elem_type", &dum_str[0] );
2255 if( NC_NOERR != fail )
2256 {
2257 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting elem type" );
2258 }
2259 ExoIIElementType elem_type = ExoIIUtil::static_element_name_to_type( &dum_str[0] );
2260 ( *i ).elemType = elem_type;
2261 const EntityType mb_type = ExoIIUtil::ExoIIElementMBEntity[( *i ).elemType];
2262
2263 // Get the number of nodes per element
2264 unsigned int nodes_per_element = ExoIIUtil::VerticesPerElement[( *i ).elemType];
2265
2266 // Read the connectivity into that memory.
2267 // int exo_conn[i->numElements][nodes_per_element];
2268 int* exo_conn = new int[i->numElements * nodes_per_element];
2269 // NcBool status = temp_var->get(&exo_conn[0][0], i->numElements, nodes_per_element);
2270 size_t lstart[2] = { 0, 0 }, lcount[2];
2271 lcount[0] = i->numElements;
2272 lcount[1] = nodes_per_element;
2273 fail = nc_get_vara_int( ncFile, nc_var, lstart, lcount, exo_conn );
2274 if( NC_NOERR != fail )
2275 {
2276 delete[] exo_conn;
2277 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting connectivity" );
2278 }
2279
2280 // Get the death_status at the correct time step.
2281 std::vector< double > death_status( i->numElements ); // It seems wrong, but it uses doubles
2282 std::string array_name( "vals_elem_var" );
2283 temp_ss.str( "" ); // stringstream won't clear by temp.clear()
2284 temp_ss << death_index;
2285 array_name += temp_ss.str();
2286 array_name += "eb";
2287 temp_ss.str( "" ); // stringstream won't clear by temp.clear()
2288 temp_ss << block_count;
2289 array_name += temp_ss.str();
2290 array_name += "\0";
2291 GET_VAR( array_name.c_str(), nc_var, ldims );
2292 if( !nc_var )
2293 {
2294 MB_SET_ERR( MB_FAILURE, "ReadNCDF:: Problem getting death status variable" );
2295 }
2296 lstart[0] = time_step - 1;
2297 lstart[1] = 0;
2298 lcount[0] = 1;
2299 lcount[1] = i->numElements;
2300 status = nc_get_vara_double( ncFile, nc_var, lstart, lcount, &death_status[0] );
2301 if( NC_NOERR != status )
2302 {
2303 delete[] exo_conn;
2304 MB_SET_ERR( MB_FAILURE, "ReadNCDF: Problem setting time step for death_status" );
2305 }
2306
2307 // Look for dead elements. If there is too many dead elements and this starts
2308 // to take too long, I should put the elems in a kd-tree for more efficient
2309 // searching. Alternatively I could get the exo connectivity and match nodes.
2310 int dead_elem_counter = 0, missing_elem_counter = 0;
2311 for( int j = 0; j < i->numElements; ++j )
2312 {
2313 if( 1 != death_status[j] )
2314 {
2315 Range cub_elem, cub_nodes;
2316 if( match_elems_by_connectivity )
2317 {
2318 // Get exodus nodes for the element
2319 std::vector< int > elem_conn( nodes_per_element );
2320 for( unsigned int k = 0; k < nodes_per_element; ++k )
2321 {
2322 // elem_conn[k] = exo_conn[j][k];
2323 elem_conn[k] = exo_conn[j * nodes_per_element + k];
2324 }
2325 // Get the ids of the nodes (assume we are matching by id)
2326 // Remember that the exodus array locations start with 1 (not 0).
2327 std::vector< int > elem_conn_node_ids( nodes_per_element );
2328 for( unsigned int k = 0; k < nodes_per_element; ++k )
2329 {
2330 elem_conn_node_ids[k] = ptr[elem_conn[k] - 1];
2331 }
2332 // Get the cub_file_set nodes by id
2333 // The map is a log search and takes almost no time.
2334 // MOAB's linear tag search takes 5-10 minutes.
2335 for( unsigned int k = 0; k < nodes_per_element; ++k )
2336 {
2337 std::map< int, EntityHandle >::iterator k_iter;
2338 k_iter = matched_cub_vert_id_map.find( elem_conn_node_ids[k] );
2339
2340 if( k_iter == matched_cub_vert_id_map.end() )
2341 {
2342 std::cout << "ReadNCDF: Found no cub node with id=" << elem_conn_node_ids[k]
2343 << ", but expected to find only 1." << std::endl;
2344 break;
2345 }
2346 cub_nodes.insert( k_iter->second );
2347 }
2348
2349 if( nodes_per_element != cub_nodes.size() )
2350 {
2351 std::cout << "ReadNCDF: nodes_per_elemenet != cub_nodes.size()" << std::endl;
2352 delete[] exo_conn;
2353 return MB_INVALID_SIZE;
2354 }
2355
2356 // Get the cub_file_set element with the same nodes
2357 int to_dim = CN::Dimension( mb_type );
2358 rval = mdbImpl->get_adjacencies( cub_nodes, to_dim, false, cub_elem );
2359 if( MB_SUCCESS != rval )
2360 {
2361 std::cout << "ReadNCDF: could not get dead cub element" << std::endl;
2362 delete[] exo_conn;
2363 return rval;
2364 }
2365
2366 // Pronto/Presto renumbers elements, so matching cub and exo elements by
2367 // id is not possible at this time.
2368 }
2369 else
2370 {
2371 // Get dead element's id
2372 int elem_id = elem_ids[first_elem_id_in_block + j];
2373 void* id[] = { &elem_id };
2374 // Get the element by id
2375 rval = mdbImpl->get_entities_by_type_and_tag( cub_file_set, mb_type, &mGlobalIdTag, id, 1, cub_elem,
2376 Interface::INTERSECT );
2377 if( MB_SUCCESS != rval )
2378 {
2379 delete[] exo_conn;
2380 return rval;
2381 }
2382 }
2383
2384 if( 1 == cub_elem.size() )
2385 {
2386 // Delete the dead element from the cub file. It will be removed from sets
2387 // ONLY if they are tracking meshsets.
2388 rval = mdbImpl->remove_entities( cub_file_set, cub_elem );
2389 if( MB_SUCCESS != rval )
2390 {
2391 delete[] exo_conn;
2392 return rval;
2393 }
2394 rval = mdbImpl->delete_entities( cub_elem );
2395 if( MB_SUCCESS != rval )
2396 {
2397 delete[] exo_conn;
2398 return rval;
2399 }
2400 }
2401 else
2402 {
2403 std::cout << "ReadNCDF: Should have found 1 element with type=" << mb_type
2404 << " in cub_file_set, but instead found " << cub_elem.size() << std::endl;
2405 rval = mdbImpl->list_entities( cub_nodes );
2406 ++missing_elem_counter;
2407 delete[] exo_conn;
2408 return MB_FAILURE;
2409 }
2410 ++dead_elem_counter;
2411 }
2412 }
2413 // Print some statistics
2414 temp_ss.str( "" );
2415 temp_ss << i->blockId;
2416 total_dead_elems += dead_elem_counter;
2417 total_elems += i->numElements;
2418 std::cout << " Block " << temp_ss.str() << " has " << dead_elem_counter << "/" << i->numElements
2419 << " dead elements." << std::endl;
2420 if( 0 != missing_elem_counter )
2421 {
2422 std::cout << " " << missing_elem_counter
2423 << " dead elements in this block were not found in the cub_file_set." << std::endl;
2424 }
2425
2426 // Advance the pointers into element ids and block_count. Memory cleanup.
2427 first_elem_id_in_block += i->numElements;
2428 ++block_count;
2429 delete[] exo_conn;
2430 }
2431
2432 std::cout << " Total: " << total_dead_elems << "/" << total_elems << " dead elements." << std::endl;
2433
2434 // Close the file
2435 fail = nc_close( ncFile );
2436 if( NC_NOERR != fail )
2437 {
2438 MB_SET_ERR( MB_FAILURE, "Trouble closing file" );
2439 }
2440
2441 ncFile = 0;
2442 return MB_SUCCESS;
2443 }
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().
|
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().
|
private |
Definition at line 194 of file ReadNCDF.hpp.
Referenced by read_exodus_header().
|
private |
file name
Definition at line 201 of file ReadNCDF.hpp.
|
private |
Definition at line 236 of file ReadNCDF.hpp.
Referenced by read_face_blocks_headers(), read_polyhedra_faces(), and reset().
|
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().
|
private |
Definition at line 195 of file ReadNCDF.hpp.
Referenced by read_exodus_header().
|
private |
Definition at line 248 of file ReadNCDF.hpp.
Referenced by read_exodus_header().
|
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().
|
private |
Meshset Handle for the mesh that is currently being read.
Definition at line 224 of file ReadNCDF.hpp.
Referenced by reset().
|
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().
|
private |
Definition at line 241 of file ReadNCDF.hpp.
Referenced by read_nodesets(), and ReadNCDF().
|
private |
Definition at line 244 of file ReadNCDF.hpp.
Referenced by read_nodesets(), read_sidesets(), and ReadNCDF().
|
private |
Definition at line 245 of file ReadNCDF.hpp.
Referenced by read_elements(), read_global_ids(), read_nodesets(), read_sidesets(), ReadNCDF(), and update().
|
private |
Definition at line 243 of file ReadNCDF.hpp.
Referenced by read_elements(), and ReadNCDF().
|
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().
|
private |
Definition at line 242 of file ReadNCDF.hpp.
Referenced by read_sidesets(), and ReadNCDF().
|
private |
Definition at line 246 of file ReadNCDF.hpp.
Referenced by read_qa_records(), and ReadNCDF().
|
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().
|
private |
Definition at line 227 of file ReadNCDF.hpp.
Referenced by read_elements(), read_nodesets(), read_polyhedra_faces(), and reset().
|
private |
Definition at line 221 of file ReadNCDF.hpp.
Referenced by number_dimensions(), read_exodus_header(), read_nodes(), reset(), and update().
|
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().
|
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().
|
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().
|
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().
|
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().
|
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().
|
private |
Definition at line 233 of file ReadNCDF.hpp.
Referenced by read_elements(), and read_polyhedra_faces().
|
private |
Definition at line 108 of file ReadNCDF.hpp.
Referenced by read_elements(), read_nodes(), ReadNCDF(), and ~ReadNCDF().
|
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().
1.9.1.