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