MOAB: Mesh Oriented datABase  (version 5.4.1)
ReadNASTRAN.cpp
Go to the documentation of this file.
00001 /**
00002  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
00003  * storing and accessing finite element mesh data.
00004  *
00005  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
00006  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
00007  * retains certain rights in this software.
00008  *
00009  * This library is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU Lesser General Public
00011  * License as published by the Free Software Foundation; either
00012  * version 2.1 of the License, or (at your option) any later version.
00013  *
00014  */
00015 
00016 #include "ReadNASTRAN.hpp"
00017 
00018 #include <iostream>
00019 #include <sstream>
00020 #include <fstream>
00021 #include <vector>
00022 #include <cstdlib>
00023 #include <cassert>
00024 #include <cmath>
00025 
00026 #include "moab/Interface.hpp"
00027 #include "moab/ReadUtilIface.hpp"
00028 #include "Internals.hpp"  // For MB_START_ID
00029 #include "moab/Range.hpp"
00030 #include "moab/FileOptions.hpp"
00031 #include "FileTokenizer.hpp"
00032 #include "MBTagConventions.hpp"
00033 #include "moab/CN.hpp"
00034 
00035 namespace moab
00036 {
00037 
00038 ReaderIface* ReadNASTRAN::factory( Interface* iface )
00039 {
00040     return new ReadNASTRAN( iface );
00041 }
00042 
00043 // Constructor
00044 ReadNASTRAN::ReadNASTRAN( Interface* impl ) : MBI( impl )
00045 {
00046     assert( NULL != impl );
00047     MBI->query_interface( readMeshIface );
00048     assert( NULL != readMeshIface );
00049 }
00050 
00051 // Destructor
00052 ReadNASTRAN::~ReadNASTRAN()
00053 {
00054     if( readMeshIface )
00055     {
00056         MBI->release_interface( readMeshIface );
00057         readMeshIface = 0;
00058     }
00059 }
00060 
00061 ErrorCode ReadNASTRAN::read_tag_values( const char* /*file_name*/,
00062                                         const char* /*tag_name*/,
00063                                         const FileOptions& /*opts*/,
00064                                         std::vector< int >& /*tag_values_out*/,
00065                                         const SubsetList* /*subset_list*/ )
00066 {
00067     return MB_NOT_IMPLEMENTED;
00068 }
00069 
00070 // Load the file as called by the Interface function
00071 ErrorCode ReadNASTRAN::load_file( const char* filename,
00072                                   const EntityHandle* /* file_set */,
00073                                   const FileOptions& /* opts */,
00074                                   const ReaderIface::SubsetList* subset_list,
00075                                   const Tag* file_id_tag )
00076 {
00077     // At this time there is no support for reading a subset of the file
00078     if( subset_list )
00079     {
00080         MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "Reading subset of files not supported for NASTRAN" );
00081     }
00082 
00083     nodeIdMap.clear();
00084     elemIdMap.clear();
00085 
00086     bool debug = false;
00087     if( debug ) std::cout << "begin ReadNASTRAN::load_file" << std::endl;
00088     ErrorCode result;
00089 
00090     // Count the entities of each type in the file. This is used to allocate the node array.
00091     int entity_count[MBMAXTYPE];
00092     for( int i = 0; i < MBMAXTYPE; i++ )
00093         entity_count[i] = 0;
00094 
00095     /* Determine the line_format of the first line. Assume that the entire file
00096        has the same format. */
00097     std::string line;
00098     std::ifstream file( filename );
00099     if( !getline( file, line ) ) return MB_FILE_DOES_NOT_EXIST;
00100     line_format format;
00101     result = determine_line_format( line, format );
00102     if( MB_SUCCESS != result ) return result;
00103 
00104     /* Count the number of each entity in the file. This allows one to allocate
00105        a sequential array of vertex handles. */
00106     while( !file.eof() )
00107     {
00108         // Cut the line into fields as determined by the line format.
00109         // Use a vector to allow for an unknown number of tokens (continue lines).
00110         // Continue lines are not implemented.
00111         std::vector< std::string > tokens;
00112         tokens.reserve( 10 );  // assume 10 fields to avoid extra vector resizing
00113         result = tokenize_line( line, format, tokens );
00114         if( MB_SUCCESS != result ) return result;
00115 
00116         // Process the tokens of the line. The first token describes the entity type.
00117         EntityType type;
00118         result = determine_entity_type( ( tokens.empty() ) ? "" : tokens.front(), type );
00119         if( MB_SUCCESS != result ) return result;
00120         entity_count[type]++;
00121         getline( file, line );
00122     }
00123 
00124     if( debug )
00125     {
00126         for( int i = 0; i < MBMAXTYPE; i++ )
00127         {
00128             std::cout << "entity_count[" << i << "]=" << entity_count[i] << std::endl;
00129         }
00130     }
00131 
00132     // Keep list of material sets
00133     std::vector< Range > materials;
00134 
00135     // Now that the number of vertices is known, create the vertices.
00136     EntityHandle start_vert = 0;
00137     std::vector< double* > coord_arrays( 3 );
00138     result = readMeshIface->get_node_coords( 3, entity_count[0], MB_START_ID, start_vert, coord_arrays );
00139     if( MB_SUCCESS != result ) return result;
00140     if( 0 == start_vert ) return MB_FAILURE;  // check for NULL
00141     int id, vert_index = 0;
00142     if( debug ) std::cout << "allocated coord arrays" << std::endl;
00143 
00144     // Read the file again to create the entities.
00145     file.clear();     // Clear eof state from object
00146     file.seekg( 0 );  // Rewind file
00147     while( !file.eof() )
00148     {
00149         getline( file, line );
00150 
00151         // Cut the line into fields as determined by the line format.
00152         // Use a vector to allow for an unknown number of tokens (continue lines).
00153         // Continue lines are not implemented.
00154         std::vector< std::string > tokens;
00155         tokens.reserve( 10 );  // assume 10 fields to avoid extra vector resizing
00156         result = tokenize_line( line, format, tokens );
00157         if( MB_SUCCESS != result ) return result;
00158 
00159         // Process the tokens of the line. The first token describes the entity type.
00160         EntityType type;
00161         result = determine_entity_type( tokens.front(), type );
00162         if( MB_SUCCESS != result ) return result;
00163 
00164         // Create the entity.
00165         if( MBVERTEX == type )
00166         {
00167             double* coords[3] = { coord_arrays[0] + vert_index, coord_arrays[1] + vert_index,
00168                                   coord_arrays[2] + vert_index };
00169             result            = read_node( tokens, debug, coords, id );
00170             if( MB_SUCCESS != result ) return result;
00171             if( !nodeIdMap.insert( id, start_vert + vert_index, 1 ).second ) return MB_FAILURE;  // Duplicate IDs!
00172             ++vert_index;
00173         }
00174         else
00175         {
00176             result = read_element( tokens, materials, type, debug );
00177             if( MB_SUCCESS != result ) return result;
00178         }
00179     }
00180 
00181     result = create_materials( materials );
00182     if( MB_SUCCESS != result ) return result;
00183 
00184     result = assign_ids( file_id_tag );
00185     if( MB_SUCCESS != result ) return result;
00186 
00187     file.close();
00188     nodeIdMap.clear();
00189     elemIdMap.clear();
00190     return MB_SUCCESS;
00191 }
00192 
00193 /* Determine the type of NASTRAN line: small field, large field, or free field.
00194    small field: each line has 10 fields of 8 characters
00195    large field: 1x8, 4x16, 1x8. Field 1 must have an asterisk following the character string
00196    free field: each line entry must be separated by a comma
00197    Implementation tries to avoid more searches than necessary. */
00198 ErrorCode ReadNASTRAN::determine_line_format( const std::string& line, line_format& format )
00199 {
00200     std::string::size_type found_asterisk = line.find( "*" );
00201     if( std::string::npos != found_asterisk )
00202     {
00203         format = LARGE_FIELD;
00204         return MB_SUCCESS;
00205     }
00206     else
00207     {
00208         std::string::size_type found_comma = line.find( "," );
00209         if( std::string::npos != found_comma )
00210         {
00211             format = FREE_FIELD;
00212             return MB_SUCCESS;
00213         }
00214         else
00215         {
00216             format = SMALL_FIELD;
00217             return MB_SUCCESS;
00218         }
00219     }
00220 }
00221 
00222 /* Tokenize the line. Continue-lines have not been implemented. */
00223 ErrorCode ReadNASTRAN::tokenize_line( const std::string& line,
00224                                       const line_format format,
00225                                       std::vector< std::string >& tokens )
00226 {
00227     size_t line_size = line.size();
00228     switch( format )
00229     {
00230         case SMALL_FIELD: {
00231             // Expect 10 fields of 8 characters.
00232             // The sample file does not have all 10 fields in each line
00233             const int field_length = 8;
00234             unsigned int n_tokens  = line_size / field_length;
00235             for( unsigned int i = 0; i < n_tokens; i++ )
00236             {
00237                 tokens.push_back( line.substr( i * field_length, field_length ) );
00238             }
00239             break;
00240         }
00241         case LARGE_FIELD:
00242             return MB_NOT_IMPLEMENTED;
00243         case FREE_FIELD:
00244             return MB_NOT_IMPLEMENTED;
00245         default:
00246             return MB_FAILURE;
00247     }
00248 
00249     return MB_SUCCESS;
00250 }
00251 
00252 ErrorCode ReadNASTRAN::determine_entity_type( const std::string& first_token, EntityType& type )
00253 {
00254     if( 0 == first_token.compare( "GRID    " ) )
00255         type = MBVERTEX;
00256     else if( 0 == first_token.compare( "CTETRA  " ) )
00257         type = MBTET;
00258     else if( 0 == first_token.compare( "CPENTA  " ) )
00259         type = MBPRISM;
00260     else if( 0 == first_token.compare( "CHEXA   " ) )
00261         type = MBHEX;
00262     else
00263         return MB_NOT_IMPLEMENTED;
00264 
00265     return MB_SUCCESS;
00266 }
00267 
00268 /* Some help from Jason:
00269    Nastran floats must contain a decimal point, may contain
00270    a leading '-' and may contain an exponent. The 'E' is optional
00271    when specifying an exponent.  A '-' or '+' at any location other
00272    than the first position indicates an exponent.  For a positive
00273    exponent, either a '+' or an 'E' must be specified.  For a
00274    negative exponent, the 'E' is option and the '-' is always specified.
00275    Examples for the real value 7.0 from mcs2006 quick reference guide:
00276    7.0  .7E1  0.7+1  .70+1  7.E+0  70.-1
00277 
00278    From the test file created in SC/Tetra:
00279    GRID           1       03.9804546.9052-15.6008-1
00280    has the coordinates: ( 3.980454, 6.9052e-1, 5.6008e-1 )
00281    GRID      200005       04.004752-3.985-15.4955-1
00282    has the coordinates: ( 4.004752, -3.985e-1, 5.4955e-1 ) */
00283 ErrorCode ReadNASTRAN::get_real( const std::string& token, double& real )
00284 {
00285     std::string significand = token;
00286     std::string exponent    = "0";
00287 
00288     // Cut off the first digit because a "-" could be here indicating a negative
00289     // number. Instead we are looking for a negative exponent.
00290     std::string back_token = token.substr( 1 );
00291 
00292     // A minus that is not the first digit is always a negative exponent
00293     std::string::size_type found_minus = back_token.find( "-" );
00294     if( std::string::npos != found_minus )
00295     {
00296         // separate the significand from the exponent at the "-"
00297         exponent    = token.substr( found_minus + 1 );
00298         significand = token.substr( 0, found_minus + 1 );
00299 
00300         // If the significand has an "E", remove it
00301         if( std::string::npos != significand.find( "E" ) )
00302             // Assume the "E" is at the end of the significand.
00303             significand = significand.substr( 1, significand.size() - 2 );
00304 
00305         // If a minus does not exist past the 1st digit, but an "E" or "+" does, then
00306         // it is a positive exponent. First look for an "E",
00307     }
00308     else
00309     {
00310         std::string::size_type found_E = token.find( "E" );
00311         if( std::string::npos != found_E )
00312         {
00313             significand = token.substr( 0, found_E - 1 );
00314             exponent    = token.substr( found_E + 1 );
00315             // If there is a "+" on the exponent, cut it off
00316             std::size_t found_plus = exponent.find( "+" );
00317             if( std::string::npos != found_plus )
00318             {
00319                 exponent = exponent.substr( found_plus + 1 );
00320             }
00321         }
00322         else
00323         {
00324             // If there is a "+" on the exponent, cut it off
00325             std::size_t found_plus = token.find( "+" );
00326             if( std::string::npos != found_plus )
00327             {
00328                 significand = token.substr( 0, found_plus - 1 );
00329                 exponent    = token.substr( found_plus + 1 );
00330             }
00331         }
00332     }
00333 
00334     // Now assemble the real number
00335     double signi = atof( significand.c_str() );
00336     double expon = atof( exponent.c_str() );
00337 
00338     if( HUGE_VAL == signi || HUGE_VAL == expon ) return MB_FAILURE;
00339 
00340     real = signi * pow( 10, expon );
00341 
00342     return MB_SUCCESS;
00343 }
00344 
00345 /* It has been determined that this line is a vertex. Read the rest of
00346    the line and create the vertex. */
00347 ErrorCode ReadNASTRAN::read_node( const std::vector< std::string >& tokens,
00348                                   const bool debug,
00349                                   double* coords[3],
00350                                   int& id )
00351 {
00352     // Read the node's id (unique)
00353     ErrorCode result;
00354     id = atoi( tokens[1].c_str() );
00355 
00356     // Read the node's coordinate system number
00357     // "0" or blank refers to the basic coordinate system.
00358     int coord_system = atoi( tokens[2].c_str() );
00359     if( 0 != coord_system )
00360     {
00361         std::cerr << "ReadNASTRAN: alternative coordinate systems not implemented" << std::endl;
00362         return MB_NOT_IMPLEMENTED;
00363     }
00364 
00365     // Read the coordinates
00366     for( unsigned int i = 0; i < 3; i++ )
00367     {
00368         result = get_real( tokens[i + 3], *coords[i] );
00369         if( MB_SUCCESS != result ) return result;
00370         if( debug ) std::cout << "read_node: coords[" << i << "]=" << coords[i] << std::endl;
00371     }
00372 
00373     return MB_SUCCESS;
00374 }
00375 
00376 /* The type of element has already been identified. Read the rest of the
00377    line and create the element. Assume that all of the nodes have already
00378    been created. */
00379 ErrorCode ReadNASTRAN::read_element( const std::vector< std::string >& tokens,
00380                                      std::vector< Range >& materials,
00381                                      const EntityType element_type,
00382                                      const bool /*debug*/ )
00383 {
00384     // Read the element's id (unique) and material set
00385     ErrorCode result;
00386     int id       = atoi( tokens[1].c_str() );
00387     int material = atoi( tokens[2].c_str() );
00388 
00389     // Resize materials list if necessary. This code is somewhat complicated
00390     // so as to avoid copying of Ranges
00391     if( material >= (int)materials.size() )
00392     {
00393         if( (int)materials.capacity() < material )
00394             materials.resize( material + 1 );
00395         else
00396         {
00397             std::vector< Range > new_mat( material + 1 );
00398             for( size_t i = 0; i < materials.size(); ++i )
00399                 new_mat[i].swap( materials[i] );
00400             materials.swap( new_mat );
00401         }
00402     }
00403 
00404     // The size of the connectivity array depends on the element type
00405     int n_conn = CN::VerticesPerEntity( element_type );
00406     EntityHandle conn_verts[27];
00407     assert( n_conn <= (int)( sizeof( conn_verts ) / sizeof( EntityHandle ) ) );
00408 
00409     // Read the connected node ids from the file
00410     for( int i = 0; i < n_conn; i++ )
00411     {
00412         int n         = atoi( tokens[3 + i].c_str() );
00413         conn_verts[i] = nodeIdMap.find( n );
00414         if( !conn_verts[i] )  // invalid vertex id
00415             return MB_FAILURE;
00416     }
00417 
00418     // Create the element and set the global id from the NASTRAN file
00419     EntityHandle element;
00420     result = MBI->create_element( element_type, conn_verts, n_conn, element );
00421     if( MB_SUCCESS != result ) return result;
00422     elemIdMap.insert( id, element, 1 );
00423 
00424     materials[material].insert( element );
00425     return MB_SUCCESS;
00426 }
00427 
00428 ErrorCode ReadNASTRAN::create_materials( const std::vector< Range >& materials )
00429 {
00430     ErrorCode result;
00431     Tag material_tag;
00432     int negone = -1;
00433     result = MBI->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, material_tag, MB_TAG_SPARSE | MB_TAG_CREAT,
00434                                   &negone );
00435     if( MB_SUCCESS != result ) return result;
00436 
00437     for( size_t i = 0; i < materials.size(); ++i )
00438     {
00439         if( materials[i].empty() ) continue;
00440 
00441         // Merge with existing or create new?  Original code effectively
00442         // created new by only merging with existing in current file set,
00443         // so do the same here. - j.kraftcheck
00444 
00445         EntityHandle handle;
00446         result = MBI->create_meshset( MESHSET_SET, handle );
00447         if( MB_SUCCESS != result ) return result;
00448 
00449         result = MBI->add_entities( handle, materials[i] );
00450         if( MB_SUCCESS != result ) return result;
00451 
00452         int id = i;
00453         result = MBI->tag_set_data( material_tag, &handle, 1, &id );
00454         if( MB_SUCCESS != result ) return result;
00455     }
00456 
00457     return MB_SUCCESS;
00458 }
00459 
00460 ErrorCode ReadNASTRAN::assign_ids( const Tag* file_id_tag )
00461 {
00462     // Create tag
00463     ErrorCode result;
00464     Tag id_tag = MBI->globalId_tag();
00465 
00466     RangeMap< int, EntityHandle >::iterator i;
00467     for( int t = 0; t < 2; ++t )
00468     {
00469         RangeMap< int, EntityHandle >& fileIdMap = t ? elemIdMap : nodeIdMap;
00470         for( i = fileIdMap.begin(); i != fileIdMap.end(); ++i )
00471         {
00472             Range range( i->value, i->value + i->count - 1 );
00473 
00474             result = readMeshIface->assign_ids( id_tag, range, i->begin );
00475             if( MB_SUCCESS != result ) return result;
00476 
00477             if( file_id_tag && *file_id_tag != id_tag )
00478             {
00479                 result = readMeshIface->assign_ids( *file_id_tag, range, i->begin );
00480                 if( MB_SUCCESS != result ) return result;
00481             }
00482         }
00483     }
00484 
00485     return MB_SUCCESS;
00486 }
00487 
00488 }  // namespace moab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines