Loading [MathJax]/extensions/tex2jax.js
Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ReadNASTRAN.cpp
Go to the documentation of this file.
1 /** 2  * MOAB, a Mesh-Oriented datABase, is a software component for creating, 3  * storing and accessing finite element mesh data. 4  * 5  * Copyright 2004 Sandia Corporation. Under the terms of Contract 6  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government 7  * retains certain rights in this software. 8  * 9  * This library is free software; you can redistribute it and/or 10  * modify it under the terms of the GNU Lesser General Public 11  * License as published by the Free Software Foundation; either 12  * version 2.1 of the License, or (at your option) any later version. 13  * 14  */ 15  16 #include "ReadNASTRAN.hpp" 17  18 #include <iostream> 19 #include <sstream> 20 #include <fstream> 21 #include <vector> 22 #include <cstdlib> 23 #include <cassert> 24 #include <cmath> 25  26 #include "moab/Interface.hpp" 27 #include "moab/ReadUtilIface.hpp" 28 #include "Internals.hpp" // For MB_START_ID 29 #include "moab/Range.hpp" 30 #include "moab/FileOptions.hpp" 31 #include "FileTokenizer.hpp" 32 #include "MBTagConventions.hpp" 33 #include "moab/CN.hpp" 34  35 namespace moab 36 { 37  38 ReaderIface* ReadNASTRAN::factory( Interface* iface ) 39 { 40  return new ReadNASTRAN( iface ); 41 } 42  43 // Constructor 44 ReadNASTRAN::ReadNASTRAN( Interface* impl ) : MBI( impl ) 45 { 46  assert( NULL != impl ); 47  MBI->query_interface( readMeshIface ); 48  assert( NULL != readMeshIface ); 49 } 50  51 // Destructor 52 ReadNASTRAN::~ReadNASTRAN() 53 { 54  if( readMeshIface ) 55  { 56  MBI->release_interface( readMeshIface ); 57  readMeshIface = 0; 58  } 59 } 60  61 ErrorCode ReadNASTRAN::read_tag_values( const char* /*file_name*/, 62  const char* /*tag_name*/, 63  const FileOptions& /*opts*/, 64  std::vector< int >& /*tag_values_out*/, 65  const SubsetList* /*subset_list*/ ) 66 { 67  return MB_NOT_IMPLEMENTED; 68 } 69  70 // Load the file as called by the Interface function 71 ErrorCode ReadNASTRAN::load_file( const char* filename, 72  const EntityHandle* /* file_set */, 73  const FileOptions& /* opts */, 74  const ReaderIface::SubsetList* subset_list, 75  const Tag* file_id_tag ) 76 { 77  // At this time there is no support for reading a subset of the file 78  if( subset_list ) 79  { 80  MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "Reading subset of files not supported for NASTRAN" ); 81  } 82  83  nodeIdMap.clear(); 84  elemIdMap.clear(); 85  86  bool debug = false; 87  if( debug ) std::cout << "begin ReadNASTRAN::load_file" << std::endl; 88  ErrorCode result; 89  90  // Count the entities of each type in the file. This is used to allocate the node array. 91  int entity_count[MBMAXTYPE]; 92  for( int i = 0; i < MBMAXTYPE; i++ ) 93  entity_count[i] = 0; 94  95  /* Determine the line_format of the first line. Assume that the entire file 96  has the same format. */ 97  std::string line; 98  std::ifstream file( filename ); 99  if( !getline( file, line ) ) return MB_FILE_DOES_NOT_EXIST; 100  line_format format; 101  result = determine_line_format( line, format ); 102  if( MB_SUCCESS != result ) return result; 103  104  /* Count the number of each entity in the file. This allows one to allocate 105  a sequential array of vertex handles. */ 106  while( !file.eof() ) 107  { 108  // Cut the line into fields as determined by the line format. 109  // Use a vector to allow for an unknown number of tokens (continue lines). 110  // Continue lines are not implemented. 111  std::vector< std::string > tokens; 112  tokens.reserve( 10 ); // assume 10 fields to avoid extra vector resizing 113  result = tokenize_line( line, format, tokens ); 114  if( MB_SUCCESS != result ) return result; 115  116  // Process the tokens of the line. The first token describes the entity type. 117  EntityType type; 118  result = determine_entity_type( ( tokens.empty() ) ? "" : tokens.front(), type ); 119  if( MB_SUCCESS != result ) return result; 120  entity_count[type]++; 121  getline( file, line ); 122  } 123  124  if( debug ) 125  { 126  for( int i = 0; i < MBMAXTYPE; i++ ) 127  { 128  std::cout << "entity_count[" << i << "]=" << entity_count[i] << std::endl; 129  } 130  } 131  132  // Keep list of material sets 133  std::vector< Range > materials; 134  135  // Now that the number of vertices is known, create the vertices. 136  EntityHandle start_vert = 0; 137  std::vector< double* > coord_arrays( 3 ); 138  result = readMeshIface->get_node_coords( 3, entity_count[0], MB_START_ID, start_vert, coord_arrays ); 139  if( MB_SUCCESS != result ) return result; 140  if( 0 == start_vert ) return MB_FAILURE; // check for NULL 141  int id, vert_index = 0; 142  if( debug ) std::cout << "allocated coord arrays" << std::endl; 143  144  // Read the file again to create the entities. 145  file.clear(); // Clear eof state from object 146  file.seekg( 0 ); // Rewind file 147  while( !file.eof() ) 148  { 149  getline( file, line ); 150  151  // Cut the line into fields as determined by the line format. 152  // Use a vector to allow for an unknown number of tokens (continue lines). 153  // Continue lines are not implemented. 154  std::vector< std::string > tokens; 155  tokens.reserve( 10 ); // assume 10 fields to avoid extra vector resizing 156  result = tokenize_line( line, format, tokens ); 157  if( MB_SUCCESS != result ) return result; 158  159  // Process the tokens of the line. The first token describes the entity type. 160  EntityType type; 161  result = determine_entity_type( tokens.front(), type ); 162  if( MB_SUCCESS != result ) return result; 163  164  // Create the entity. 165  if( MBVERTEX == type ) 166  { 167  double* coords[3] = { coord_arrays[0] + vert_index, coord_arrays[1] + vert_index, 168  coord_arrays[2] + vert_index }; 169  result = read_node( tokens, debug, coords, id ); 170  if( MB_SUCCESS != result ) return result; 171  if( !nodeIdMap.insert( id, start_vert + vert_index, 1 ).second ) return MB_FAILURE; // Duplicate IDs! 172  ++vert_index; 173  } 174  else 175  { 176  result = read_element( tokens, materials, type, debug ); 177  if( MB_SUCCESS != result ) return result; 178  } 179  } 180  181  result = create_materials( materials ); 182  if( MB_SUCCESS != result ) return result; 183  184  result = assign_ids( file_id_tag ); 185  if( MB_SUCCESS != result ) return result; 186  187  file.close(); 188  nodeIdMap.clear(); 189  elemIdMap.clear(); 190  return MB_SUCCESS; 191 } 192  193 /* Determine the type of NASTRAN line: small field, large field, or free field. 194  small field: each line has 10 fields of 8 characters 195  large field: 1x8, 4x16, 1x8. Field 1 must have an asterisk following the character string 196  free field: each line entry must be separated by a comma 197  Implementation tries to avoid more searches than necessary. */ 198 ErrorCode ReadNASTRAN::determine_line_format( const std::string& line, line_format& format ) 199 { 200  std::string::size_type found_asterisk = line.find( "*" ); 201  if( std::string::npos != found_asterisk ) 202  { 203  format = LARGE_FIELD; 204  return MB_SUCCESS; 205  } 206  else 207  { 208  std::string::size_type found_comma = line.find( "," ); 209  if( std::string::npos != found_comma ) 210  { 211  format = FREE_FIELD; 212  return MB_SUCCESS; 213  } 214  else 215  { 216  format = SMALL_FIELD; 217  return MB_SUCCESS; 218  } 219  } 220 } 221  222 /* Tokenize the line. Continue-lines have not been implemented. */ 223 ErrorCode ReadNASTRAN::tokenize_line( const std::string& line, 224  const line_format format, 225  std::vector< std::string >& tokens ) 226 { 227  size_t line_size = line.size(); 228  switch( format ) 229  { 230  case SMALL_FIELD: { 231  // Expect 10 fields of 8 characters. 232  // The sample file does not have all 10 fields in each line 233  const int field_length = 8; 234  unsigned int n_tokens = line_size / field_length; 235  for( unsigned int i = 0; i < n_tokens; i++ ) 236  { 237  tokens.push_back( line.substr( i * field_length, field_length ) ); 238  } 239  break; 240  } 241  case LARGE_FIELD: 242  return MB_NOT_IMPLEMENTED; 243  case FREE_FIELD: 244  return MB_NOT_IMPLEMENTED; 245  default: 246  return MB_FAILURE; 247  } 248  249  return MB_SUCCESS; 250 } 251  252 ErrorCode ReadNASTRAN::determine_entity_type( const std::string& first_token, EntityType& type ) 253 { 254  if( 0 == first_token.compare( "GRID " ) ) 255  type = MBVERTEX; 256  else if( 0 == first_token.compare( "CTETRA " ) ) 257  type = MBTET; 258  else if( 0 == first_token.compare( "CPENTA " ) ) 259  type = MBPRISM; 260  else if( 0 == first_token.compare( "CHEXA " ) ) 261  type = MBHEX; 262  else 263  return MB_NOT_IMPLEMENTED; 264  265  return MB_SUCCESS; 266 } 267  268 /* Some help from Jason: 269  Nastran floats must contain a decimal point, may contain 270  a leading '-' and may contain an exponent. The 'E' is optional 271  when specifying an exponent. A '-' or '+' at any location other 272  than the first position indicates an exponent. For a positive 273  exponent, either a '+' or an 'E' must be specified. For a 274  negative exponent, the 'E' is option and the '-' is always specified. 275  Examples for the real value 7.0 from mcs2006 quick reference guide: 276  7.0 .7E1 0.7+1 .70+1 7.E+0 70.-1 277  278  From the test file created in SC/Tetra: 279  GRID 1 03.9804546.9052-15.6008-1 280  has the coordinates: ( 3.980454, 6.9052e-1, 5.6008e-1 ) 281  GRID 200005 04.004752-3.985-15.4955-1 282  has the coordinates: ( 4.004752, -3.985e-1, 5.4955e-1 ) */ 283 ErrorCode ReadNASTRAN::get_real( const std::string& token, double& real ) 284 { 285  std::string significand = token; 286  std::string exponent = "0"; 287  288  // Cut off the first digit because a "-" could be here indicating a negative 289  // number. Instead we are looking for a negative exponent. 290  std::string back_token = token.substr( 1 ); 291  292  // A minus that is not the first digit is always a negative exponent 293  std::string::size_type found_minus = back_token.find( "-" ); 294  if( std::string::npos != found_minus ) 295  { 296  // separate the significand from the exponent at the "-" 297  exponent = token.substr( found_minus + 1 ); 298  significand = token.substr( 0, found_minus + 1 ); 299  300  // If the significand has an "E", remove it 301  if( std::string::npos != significand.find( "E" ) ) 302  // Assume the "E" is at the end of the significand. 303  significand = significand.substr( 1, significand.size() - 2 ); 304  305  // If a minus does not exist past the 1st digit, but an "E" or "+" does, then 306  // it is a positive exponent. First look for an "E", 307  } 308  else 309  { 310  std::string::size_type found_E = token.find( "E" ); 311  if( std::string::npos != found_E ) 312  { 313  significand = token.substr( 0, found_E - 1 ); 314  exponent = token.substr( found_E + 1 ); 315  // If there is a "+" on the exponent, cut it off 316  std::size_t found_plus = exponent.find( "+" ); 317  if( std::string::npos != found_plus ) 318  { 319  exponent = exponent.substr( found_plus + 1 ); 320  } 321  } 322  else 323  { 324  // If there is a "+" on the exponent, cut it off 325  std::size_t found_plus = token.find( "+" ); 326  if( std::string::npos != found_plus ) 327  { 328  significand = token.substr( 0, found_plus - 1 ); 329  exponent = token.substr( found_plus + 1 ); 330  } 331  } 332  } 333  334  // Now assemble the real number 335  double signi = atof( significand.c_str() ); 336  double expon = atof( exponent.c_str() ); 337  338  if( HUGE_VAL == signi || HUGE_VAL == expon ) return MB_FAILURE; 339  340  real = signi * pow( 10, expon ); 341  342  return MB_SUCCESS; 343 } 344  345 /* It has been determined that this line is a vertex. Read the rest of 346  the line and create the vertex. */ 347 ErrorCode ReadNASTRAN::read_node( const std::vector< std::string >& tokens, 348  const bool debug, 349  double* coords[3], 350  int& id ) 351 { 352  // Read the node's id (unique) 353  ErrorCode result; 354  id = atoi( tokens[1].c_str() ); 355  356  // Read the node's coordinate system number 357  // "0" or blank refers to the basic coordinate system. 358  int coord_system = atoi( tokens[2].c_str() ); 359  if( 0 != coord_system ) 360  { 361  std::cerr << "ReadNASTRAN: alternative coordinate systems not implemented" << std::endl; 362  return MB_NOT_IMPLEMENTED; 363  } 364  365  // Read the coordinates 366  for( unsigned int i = 0; i < 3; i++ ) 367  { 368  result = get_real( tokens[i + 3], *coords[i] ); 369  if( MB_SUCCESS != result ) return result; 370  if( debug ) std::cout << "read_node: coords[" << i << "]=" << coords[i] << std::endl; 371  } 372  373  return MB_SUCCESS; 374 } 375  376 /* The type of element has already been identified. Read the rest of the 377  line and create the element. Assume that all of the nodes have already 378  been created. */ 379 ErrorCode ReadNASTRAN::read_element( const std::vector< std::string >& tokens, 380  std::vector< Range >& materials, 381  const EntityType element_type, 382  const bool /*debug*/ ) 383 { 384  // Read the element's id (unique) and material set 385  ErrorCode result; 386  int id = atoi( tokens[1].c_str() ); 387  int material = atoi( tokens[2].c_str() ); 388  389  // Resize materials list if necessary. This code is somewhat complicated 390  // so as to avoid copying of Ranges 391  if( material >= (int)materials.size() ) 392  { 393  if( (int)materials.capacity() < material ) 394  materials.resize( material + 1 ); 395  else 396  { 397  std::vector< Range > new_mat( material + 1 ); 398  for( size_t i = 0; i < materials.size(); ++i ) 399  new_mat[i].swap( materials[i] ); 400  materials.swap( new_mat ); 401  } 402  } 403  404  // The size of the connectivity array depends on the element type 405  int n_conn = CN::VerticesPerEntity( element_type ); 406  EntityHandle conn_verts[27]; 407  assert( n_conn <= (int)( sizeof( conn_verts ) / sizeof( EntityHandle ) ) ); 408  409  // Read the connected node ids from the file 410  for( int i = 0; i < n_conn; i++ ) 411  { 412  int n = atoi( tokens[3 + i].c_str() ); 413  conn_verts[i] = nodeIdMap.find( n ); 414  if( !conn_verts[i] ) // invalid vertex id 415  return MB_FAILURE; 416  } 417  418  // Create the element and set the global id from the NASTRAN file 419  EntityHandle element; 420  result = MBI->create_element( element_type, conn_verts, n_conn, element ); 421  if( MB_SUCCESS != result ) return result; 422  elemIdMap.insert( id, element, 1 ); 423  424  materials[material].insert( element ); 425  return MB_SUCCESS; 426 } 427  428 ErrorCode ReadNASTRAN::create_materials( const std::vector< Range >& materials ) 429 { 430  ErrorCode result; 431  Tag material_tag; 432  int negone = -1; 433  result = MBI->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, material_tag, MB_TAG_SPARSE | MB_TAG_CREAT, 434  &negone ); 435  if( MB_SUCCESS != result ) return result; 436  437  for( size_t i = 0; i < materials.size(); ++i ) 438  { 439  if( materials[i].empty() ) continue; 440  441  // Merge with existing or create new? Original code effectively 442  // created new by only merging with existing in current file set, 443  // so do the same here. - j.kraftcheck 444  445  EntityHandle handle; 446  result = MBI->create_meshset( MESHSET_SET, handle ); 447  if( MB_SUCCESS != result ) return result; 448  449  result = MBI->add_entities( handle, materials[i] ); 450  if( MB_SUCCESS != result ) return result; 451  452  int id = i; 453  result = MBI->tag_set_data( material_tag, &handle, 1, &id ); 454  if( MB_SUCCESS != result ) return result; 455  } 456  457  return MB_SUCCESS; 458 } 459  460 ErrorCode ReadNASTRAN::assign_ids( const Tag* file_id_tag ) 461 { 462  // Create tag 463  ErrorCode result; 464  Tag id_tag = MBI->globalId_tag(); 465  466  RangeMap< int, EntityHandle >::iterator i; 467  for( int t = 0; t < 2; ++t ) 468  { 469  RangeMap< int, EntityHandle >& fileIdMap = t ? elemIdMap : nodeIdMap; 470  for( i = fileIdMap.begin(); i != fileIdMap.end(); ++i ) 471  { 472  Range range( i->value, i->value + i->count - 1 ); 473  474  result = readMeshIface->assign_ids( id_tag, range, i->begin ); 475  if( MB_SUCCESS != result ) return result; 476  477  if( file_id_tag && *file_id_tag != id_tag ) 478  { 479  result = readMeshIface->assign_ids( *file_id_tag, range, i->begin ); 480  if( MB_SUCCESS != result ) return result; 481  } 482  } 483  } 484  485  return MB_SUCCESS; 486 } 487  488 } // namespace moab