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 */1516#include"ReadNASTRAN.hpp"1718#include<iostream>19#include<sstream>20#include<fstream>21#include<vector>22#include<cstdlib>23#include<cassert>24#include<cmath>2526#include"moab/Interface.hpp"27#include"moab/ReadUtilIface.hpp"28#include"Internals.hpp"// For MB_START_ID29#include"moab/Range.hpp"30#include"moab/FileOptions.hpp"31#include"FileTokenizer.hpp"32#include"MBTagConventions.hpp"33#include"moab/CN.hpp"3435namespace moab
36 {
3738ReaderIface* ReadNASTRAN::factory( Interface* iface )
39 {
40returnnewReadNASTRAN( iface );
41 }
4243// Constructor44 ReadNASTRAN::ReadNASTRAN( Interface* impl ) : MBI( impl )
45 {
46assert( NULL != impl );
47 MBI->query_interface( readMeshIface );
48assert( NULL != readMeshIface );
49 }
5051// Destructor52 ReadNASTRAN::~ReadNASTRAN()
53 {
54if( readMeshIface )
55 {
56 MBI->release_interface( readMeshIface );
57 readMeshIface = 0;
58 }
59 }
6061ErrorCode ReadNASTRAN::read_tag_values( constchar* /*file_name*/,
62constchar* /*tag_name*/,
63const FileOptions& /*opts*/,
64 std::vector< int >& /*tag_values_out*/,
65const SubsetList* /*subset_list*/ )
66 {
67return MB_NOT_IMPLEMENTED;
68 }
6970// Load the file as called by the Interface function71ErrorCode ReadNASTRAN::load_file( constchar* filename,
72const EntityHandle* /* file_set */,
73const FileOptions& /* opts */,
74const ReaderIface::SubsetList* subset_list,
75const Tag* file_id_tag )
76 {
77// At this time there is no support for reading a subset of the file78if( subset_list )
79 {
80MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "Reading subset of files not supported for NASTRAN" );
81 }
8283 nodeIdMap.clear();
84 elemIdMap.clear();
8586bool debug = false;
87if( debug ) std::cout << "begin ReadNASTRAN::load_file" << std::endl;
88 ErrorCode result;
8990// Count the entities of each type in the file. This is used to allocate the node array.91int entity_count[MBMAXTYPE];
92for( int i = 0; i < MBMAXTYPE; i++ )
93 entity_count[i] = 0;
9495/* Determine the line_format of the first line. Assume that the entire file
96 has the same format. */97 std::string line;
98std::ifstream file( filename );
99if( !getline( file, line ) ) return MB_FILE_DOES_NOT_EXIST;
100 line_format format;
101 result = determine_line_format( line, format );
102if( MB_SUCCESS != result ) return result;
103104/* Count the number of each entity in the file. This allows one to allocate
105 a sequential array of vertex handles. */106while( !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 resizing113 result = tokenize_line( line, format, tokens );
114if( MB_SUCCESS != result ) return result;
115116// 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 );
119if( MB_SUCCESS != result ) return result;
120 entity_count[type]++;
121getline( file, line );
122 }
123124if( debug )
125 {
126for( int i = 0; i < MBMAXTYPE; i++ )
127 {
128 std::cout << "entity_count[" << i << "]=" << entity_count[i] << std::endl;
129 }
130 }
131132// Keep list of material sets133 std::vector< Range > materials;
134135// Now that the number of vertices is known, create the vertices.136 EntityHandle start_vert = 0;
137std::vector< double* > coord_arrays( 3 );
138 result = readMeshIface->get_node_coords( 3, entity_count[0], MB_START_ID, start_vert, coord_arrays );
139if( MB_SUCCESS != result ) return result;
140if( 0 == start_vert ) return MB_FAILURE; // check for NULL141int id, vert_index = 0;
142if( debug ) std::cout << "allocated coord arrays" << std::endl;
143144// Read the file again to create the entities.145 file.clear(); // Clear eof state from object146 file.seekg( 0 ); // Rewind file147while( !file.eof() )
148 {
149getline( file, line );
150151// 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 resizing156 result = tokenize_line( line, format, tokens );
157if( MB_SUCCESS != result ) return result;
158159// Process the tokens of the line. The first token describes the entity type.160 EntityType type;
161 result = determine_entity_type( tokens.front(), type );
162if( MB_SUCCESS != result ) return result;
163164// Create the entity.165if( MBVERTEX == type )
166 {
167double* 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 );
170if( MB_SUCCESS != result ) return result;
171if( !nodeIdMap.insert( id, start_vert + vert_index, 1 ).second ) return MB_FAILURE; // Duplicate IDs!172 ++vert_index;
173 }
174else175 {
176 result = read_element( tokens, materials, type, debug );
177if( MB_SUCCESS != result ) return result;
178 }
179 }
180181 result = create_materials( materials );
182if( MB_SUCCESS != result ) return result;
183184 result = assign_ids( file_id_tag );
185if( MB_SUCCESS != result ) return result;
186187 file.close();
188 nodeIdMap.clear();
189 elemIdMap.clear();
190return MB_SUCCESS;
191 }
192193/* 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. */198ErrorCode ReadNASTRAN::determine_line_format( const std::string& line, line_format& format )
199 {
200 std::string::size_type found_asterisk = line.find( "*" );
201if( std::string::npos != found_asterisk )
202 {
203 format = LARGE_FIELD;
204return MB_SUCCESS;
205 }
206else207 {
208 std::string::size_type found_comma = line.find( "," );
209if( std::string::npos != found_comma )
210 {
211 format = FREE_FIELD;
212return MB_SUCCESS;
213 }
214else215 {
216 format = SMALL_FIELD;
217return MB_SUCCESS;
218 }
219 }
220 }
221222/* Tokenize the line. Continue-lines have not been implemented. */223ErrorCode ReadNASTRAN::tokenize_line( const std::string& line,
224const line_format format,
225 std::vector< std::string >& tokens )
226 {
227size_t line_size = line.size();
228switch( format )
229 {
230case SMALL_FIELD: {
231// Expect 10 fields of 8 characters.232// The sample file does not have all 10 fields in each line233constint field_length = 8;
234unsignedint n_tokens = line_size / field_length;
235for( unsignedint i = 0; i < n_tokens; i++ )
236 {
237 tokens.push_back( line.substr( i * field_length, field_length ) );
238 }
239break;
240 }
241case LARGE_FIELD:
242return MB_NOT_IMPLEMENTED;
243case FREE_FIELD:
244return MB_NOT_IMPLEMENTED;
245default:
246return MB_FAILURE;
247 }
248249return MB_SUCCESS;
250 }
251252ErrorCode ReadNASTRAN::determine_entity_type( const std::string& first_token, EntityType& type )
253 {
254if( 0 == first_token.compare( "GRID " ) )
255 type = MBVERTEX;
256elseif( 0 == first_token.compare( "CTETRA " ) )
257 type = MBTET;
258elseif( 0 == first_token.compare( "CPENTA " ) )
259 type = MBPRISM;
260elseif( 0 == first_token.compare( "CHEXA " ) )
261 type = MBHEX;
262else263return MB_NOT_IMPLEMENTED;
264265return MB_SUCCESS;
266 }
267268/* 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 ) */283ErrorCode ReadNASTRAN::get_real( const std::string& token, double& real )
284 {
285 std::string significand = token;
286 std::string exponent = "0";
287288// Cut off the first digit because a "-" could be here indicating a negative289// number. Instead we are looking for a negative exponent.290 std::string back_token = token.substr( 1 );
291292// A minus that is not the first digit is always a negative exponent293 std::string::size_type found_minus = back_token.find( "-" );
294if( 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 );
299300// If the significand has an "E", remove it301if( 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 );
304305// If a minus does not exist past the 1st digit, but an "E" or "+" does, then306// it is a positive exponent. First look for an "E",307 }
308else309 {
310 std::string::size_type found_E = token.find( "E" );
311if( 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 off316 std::size_t found_plus = exponent.find( "+" );
317if( std::string::npos != found_plus )
318 {
319 exponent = exponent.substr( found_plus + 1 );
320 }
321 }
322else323 {
324// If there is a "+" on the exponent, cut it off325 std::size_t found_plus = token.find( "+" );
326if( std::string::npos != found_plus )
327 {
328 significand = token.substr( 0, found_plus - 1 );
329 exponent = token.substr( found_plus + 1 );
330 }
331 }
332 }
333334// Now assemble the real number335double signi = atof( significand.c_str() );
336double expon = atof( exponent.c_str() );
337338if( HUGE_VAL == signi || HUGE_VAL == expon ) return MB_FAILURE;
339340 real = signi * pow( 10, expon );
341342return MB_SUCCESS;
343 }
344345/* It has been determined that this line is a vertex. Read the rest of
346 the line and create the vertex. */347ErrorCode ReadNASTRAN::read_node( const std::vector< std::string >& tokens,
348constbool debug,
349double* coords[3],
350int& id )
351 {
352// Read the node's id (unique)353 ErrorCode result;
354 id = atoi( tokens[1].c_str() );
355356// Read the node's coordinate system number357// "0" or blank refers to the basic coordinate system.358int coord_system = atoi( tokens[2].c_str() );
359if( 0 != coord_system )
360 {
361 std::cerr << "ReadNASTRAN: alternative coordinate systems not implemented" << std::endl;
362return MB_NOT_IMPLEMENTED;
363 }
364365// Read the coordinates366for( unsignedint i = 0; i < 3; i++ )
367 {
368 result = get_real( tokens[i + 3], *coords[i] );
369if( MB_SUCCESS != result ) return result;
370if( debug ) std::cout << "read_node: coords[" << i << "]=" << coords[i] << std::endl;
371 }
372373return MB_SUCCESS;
374 }
375376/* 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. */379ErrorCode ReadNASTRAN::read_element( const std::vector< std::string >& tokens,
380 std::vector< Range >& materials,
381const EntityType element_type,
382constbool/*debug*/ )
383 {
384// Read the element's id (unique) and material set385 ErrorCode result;
386int id = atoi( tokens[1].c_str() );
387int material = atoi( tokens[2].c_str() );
388389// Resize materials list if necessary. This code is somewhat complicated390// so as to avoid copying of Ranges391if( material >= (int)materials.size() )
392 {
393if( (int)materials.capacity() < material )
394 materials.resize( material + 1 );
395else396 {
397 std::vector< Range > new_mat( material + 1 );
398for( size_t i = 0; i < materials.size(); ++i )
399 new_mat[i].swap( materials[i] );
400 materials.swap( new_mat );
401 }
402 }
403404// The size of the connectivity array depends on the element type405int n_conn = CN::VerticesPerEntity( element_type );
406 EntityHandle conn_verts[27];
407assert( n_conn <= (int)( sizeof( conn_verts ) / sizeof( EntityHandle ) ) );
408409// Read the connected node ids from the file410for( int i = 0; i < n_conn; i++ )
411 {
412int n = atoi( tokens[3 + i].c_str() );
413 conn_verts[i] = nodeIdMap.find( n );
414if( !conn_verts[i] ) // invalid vertex id415return MB_FAILURE;
416 }
417418// Create the element and set the global id from the NASTRAN file419 EntityHandle element;
420 result = MBI->create_element( element_type, conn_verts, n_conn, element );
421if( MB_SUCCESS != result ) return result;
422 elemIdMap.insert( id, element, 1 );
423424 materials[material].insert( element );
425return MB_SUCCESS;
426 }
427428ErrorCode ReadNASTRAN::create_materials( const std::vector< Range >& materials )
429 {
430 ErrorCode result;
431 Tag material_tag;
432int 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 );
435if( MB_SUCCESS != result ) return result;
436437for( size_t i = 0; i < materials.size(); ++i )
438 {
439if( materials[i].empty() ) continue;
440441// Merge with existing or create new? Original code effectively442// created new by only merging with existing in current file set,443// so do the same here. - j.kraftcheck444445 EntityHandle handle;
446 result = MBI->create_meshset( MESHSET_SET, handle );
447if( MB_SUCCESS != result ) return result;
448449 result = MBI->add_entities( handle, materials[i] );
450if( MB_SUCCESS != result ) return result;
451452int id = i;
453 result = MBI->tag_set_data( material_tag, &handle, 1, &id );
454if( MB_SUCCESS != result ) return result;
455 }
456457return MB_SUCCESS;
458 }
459460ErrorCode ReadNASTRAN::assign_ids( const Tag* file_id_tag )
461 {
462// Create tag463 ErrorCode result;
464 Tag id_tag = MBI->globalId_tag();
465466 RangeMap< int, EntityHandle >::iterator i;
467for( int t = 0; t < 2; ++t )
468 {
469 RangeMap< int, EntityHandle >& fileIdMap = t ? elemIdMap : nodeIdMap;
470for( i = fileIdMap.begin(); i != fileIdMap.end(); ++i )
471 {
472Range range( i->value, i->value + i->count - 1 );
473474 result = readMeshIface->assign_ids( id_tag, range, i->begin );
475if( MB_SUCCESS != result ) return result;
476477if( file_id_tag && *file_id_tag != id_tag )
478 {
479 result = readMeshIface->assign_ids( *file_id_tag, range, i->begin );
480if( MB_SUCCESS != result ) return result;
481 }
482 }
483 }
484485return MB_SUCCESS;
486 }
487488 } // namespace moab