libMesh::AbaqusIO Class Reference

#include <abaqus_io.h>

Inheritance diagram for libMesh::AbaqusIO:

Public Member Functions

 AbaqusIO (MeshBase &mesh)
 
virtual ~AbaqusIO ()
 
virtual void read (const std::string &name)
 

Public Attributes

bool build_sidesets_from_nodesets
 

Protected Member Functions

MeshBasemesh ()
 
void skip_comment_lines (std::istream &in, const char comment_start)
 

Protected Attributes

std::vector< bool > elems_of_dimension
 

Private Types

typedef std::map< std::string,
std::vector< dof_id_type > > 
container_t
 
typedef std::map< std::string,
std::vector< std::pair
< dof_id_type, unsigned > > > 
sideset_container_t
 

Private Member Functions

void read_nodes ()
 
void read_elements (std::string upper)
 
std::string parse_label (std::string line, std::string label_name)
 
void read_ids (std::string set_name, container_t &container)
 
void assign_subdomain_ids ()
 
void read_sideset (std::string sideset_name, sideset_container_t &container)
 
void assign_boundary_node_ids ()
 
void assign_sideset_ids ()
 
void process_and_discard_comments ()
 

Private Attributes

container_t _nodeset_ids
 
container_t _elemset_ids
 
sideset_container_t _sideset_ids
 
std::ifstream _in
 
std::set< ElemType > _elem_types
 
std::map< dof_id_type,
dof_id_type
_abaqus_to_libmesh_elem_mapping
 
std::map< dof_id_type,
dof_id_type
_abaqus_to_libmesh_node_mapping
 
bool _already_seen_part
 

Detailed Description

The AbaqusIO class is a preliminary implementation for reading Abaqus mesh files in ASCII format.

Author
John W. Peterson, 2011.

Definition at line 38 of file abaqus_io.h.

Member Typedef Documentation

typedef std::map<std::string, std::vector<dof_id_type> > libMesh::AbaqusIO::container_t
private

The type of data structure used to store Node and Elemset IDs.

Definition at line 68 of file abaqus_io.h.

typedef std::map<std::string, std::vector<std::pair<dof_id_type, unsigned> > > libMesh::AbaqusIO::sideset_container_t
private

Type of the data structure for storing the (elem ID, side) pairs defining sidesets. These come from the *Surface sections of the input file.

Definition at line 75 of file abaqus_io.h.

Constructor & Destructor Documentation

libMesh::AbaqusIO::AbaqusIO ( MeshBase mesh)
explicit

Constructor. Takes a writeable reference to a mesh object.

Definition at line 170 of file abaqus_io.C.

170  :
171  MeshInput<MeshBase> (mesh_in),
173  _already_seen_part(false)
174  {
175  }
libMesh::AbaqusIO::~AbaqusIO ( )
virtual

Destructor.

Definition at line 180 of file abaqus_io.C.

181  {
182  }

Member Function Documentation

void libMesh::AbaqusIO::assign_boundary_node_ids ( )
private

This function assigns boundary IDs to node sets based on the alphabetical order in which the sets are labelled in the Abaqus file. We choose the alphabetical ordering simply because Abaqus does not provide a numerical one within the file.

Definition at line 854 of file abaqus_io.C.

References _abaqus_to_libmesh_node_mapping, _nodeset_ids, libMesh::MeshBase::boundary_info, libMesh::err, libMesh::MeshInput< MT >::mesh(), libMesh::MeshBase::node_ptr(), and libMesh::out.

Referenced by read().

855  {
856  // Get a reference to the mesh we are reading
857  MeshBase& the_mesh = MeshInput<MeshBase>::mesh();
858 
859  // Iterate over the container of nodesets
860  container_t::iterator it=_nodeset_ids.begin();
861  for (unsigned current_id=0; it != _nodeset_ids.end(); ++it, ++current_id)
862  {
863  libMesh::out << "Assigning node boundary ID " << current_id << " to nodeset '"
864  << (*it).first
865  << "'." << std::endl;
866 
867  // Get a reference to the current vector of nodeset ID values
868  std::vector<dof_id_type>& nodeset_ids = (*it).second;
869 
870  for (std::size_t i=0; i<nodeset_ids.size(); ++i)
871  {
872  // Map the Abaqus global node ID to the libmesh node ID
873  dof_id_type libmesh_global_node_id = _abaqus_to_libmesh_node_mapping[nodeset_ids[i]];
874 
875  // Get node pointer from the mesh
876  Node* node = the_mesh.node_ptr(libmesh_global_node_id);
877 
878  if (node == NULL)
879  {
880  libMesh::err << "Error! Mesh returned NULL node pointer!" << std::endl;
881  libmesh_error();
882  }
883 
884  // Add this node with the current_id (which is determined by the
885  // alphabetical ordering of the map) to the BoundaryInfo object
886  the_mesh.boundary_info->add_node(node, current_id);
887  }
888  }
889 
890  } // assign_boundary_node_ids()
void libMesh::AbaqusIO::assign_sideset_ids ( )
private

Called at the end of the read() function, assigns any sideset IDs found when reading the file to the BoundaryInfo object.

Definition at line 895 of file abaqus_io.C.

References _abaqus_to_libmesh_elem_mapping, _sideset_ids, libMesh::MeshBase::boundary_info, libMesh::MeshBase::elem(), libMesh::Utility::enum_to_string(), libMesh::err, libMesh::MeshInput< MT >::mesh(), libMesh::out, and libMesh::Elem::type().

Referenced by read().

896  {
897  // Get a reference to the mesh we are reading
898  MeshBase& the_mesh = MeshInput<MeshBase>::mesh();
899 
900  // initialize the eletypes map (eletypes is a file-global variable)
901  init_eletypes();
902 
903  // Iterate over the container of sidesets
904  sideset_container_t::iterator it=_sideset_ids.begin();
905  for (unsigned current_id=0; it != _sideset_ids.end(); ++it, ++current_id)
906  {
907  libMesh::out << "Assigning sideset ID " << current_id << " to sideset '"
908  << (*it).first
909  << "'." << std::endl;
910 
911  // Get a reference to the current vector of nodeset ID values
912  std::vector<std::pair<dof_id_type,unsigned> >& sideset_ids = (*it).second;
913 
914  for (std::size_t i=0; i<sideset_ids.size(); ++i)
915  {
916  // sideset_ids is a vector of pairs (elem id, side id). Pull them out
917  // now to make the code below more readable.
918  dof_id_type abaqus_elem_id = sideset_ids[i].first;
919  unsigned abaqus_side_number = sideset_ids[i].second;
920 
921  // Map the Abaqus element ID to LibMesh numbering
922  dof_id_type libmesh_elem_id = _abaqus_to_libmesh_elem_mapping[ abaqus_elem_id ];
923 
924  // Get pointer to that element
925  Elem* elem = the_mesh.elem(libmesh_elem_id);
926 
927  // Check that the pointer returned from the Mesh is non-NULL
928  if (elem == NULL)
929  {
930  libMesh::err << "Mesh returned NULL pointer for Elem " << libmesh_elem_id << std::endl;
931  libmesh_error();
932  }
933 
934  // Grab a reference to the element definition for this element type
935  const ElementDefinition& eledef = eletypes[elem->type()];
936 
937  // If the element definition was not found, the call above would have
938  // created one with an uninitialized struct. Check for that here...
939  if (eledef.abaqus_zero_based_side_id_to_libmesh_side_id.size() == 0)
940  {
941  libMesh::err << "No Abaqus->LibMesh mapping information for ElemType " << Utility::enum_to_string(elem->type()) << "!" << std::endl;
942  libmesh_error();
943  }
944 
945  // Add this node with the current_id (which is determined by the
946  // alphabetical ordering of the map). Side numbers in Abaqus are 1-based,
947  // so we subtract 1 here before passing the abaqus side number to the
948  // mapping array
949  the_mesh.boundary_info->add_side(elem,
950  eledef.abaqus_zero_based_side_id_to_libmesh_side_id[abaqus_side_number-1],
951  current_id);
952  }
953  }
954  } // assign_sideset_ids()
void libMesh::AbaqusIO::assign_subdomain_ids ( )
private

This function is called after all the elements have been read and assigns element subdomain IDs.

The IDs are simply chosen in the order in which the elset labels are stored in the map (roughly alphabetically). To make this easy on people who are planning to use Exodus output, we'll assign different geometric elements to different (but related) subdomains, i.e. assuming there are E elemsets:

Elemset 0, Geometric Type 0: ID 0 Elemset 0, Geometric Type 1: ID 0+E ...

Elemset 0, Geometric Type N: ID 0+N*E

Elemset 1, Geometric Type 0: ID 1 Elemset 1, Geometric Type 1: ID 1+E ... Elemset 1, Geometric Type N: ID 1+N*E etc.

Definition at line 799 of file abaqus_io.C.

References _abaqus_to_libmesh_elem_mapping, _elem_types, _elemset_ids, libMesh::MeshBase::elem(), libMesh::err, libMesh::MeshInput< MT >::mesh(), libMesh::Elem::subdomain_id(), and libMesh::Elem::type().

Referenced by read().

800  {
801  // Get a reference to the mesh we are reading
802  MeshBase& the_mesh = MeshInput<MeshBase>::mesh();
803 
804  // The number of elemsets we've found while reading
805  std::size_t n_elemsets = _elemset_ids.size();
806 
807  // Fill in a temporary map with (ElemType, index) pairs based on the _elem_types set. This
808  // will allow us to easily look up this index in the loop below.
809  std::map<ElemType, unsigned> elem_types_map;
810  {
811  unsigned ctr=0;
812  for (std::set<ElemType>::iterator it=_elem_types.begin(); it!=_elem_types.end(); ++it)
813  elem_types_map[*it] = ctr++;
814  }
815 
816  // Loop over each Elemset and assign subdomain IDs to Mesh elements
817  {
818  // The elemset_id counter assigns a logical numbering to the _elemset_ids keys
819  container_t::iterator it=_elemset_ids.begin();
820  for (unsigned elemset_id=0; it != _elemset_ids.end(); ++it, ++elemset_id)
821  {
822  // Grab a reference to the vector of IDs
823  std::vector<dof_id_type>& id_vector = (*it).second;
824 
825  // Loop over this vector
826  for (std::size_t i=0; i<id_vector.size(); ++i)
827  {
828  // Map the id_vector[i]'th element ID (Abaqus numbering) to LibMesh numbering
829  dof_id_type libmesh_elem_id = _abaqus_to_libmesh_elem_mapping[ id_vector[i] ];
830 
831  // Get pointer to that element
832  Elem* elem = the_mesh.elem(libmesh_elem_id);
833 
834  if (elem == NULL)
835  {
836  libMesh::err << "Mesh returned NULL pointer for Elem " << libmesh_elem_id << std::endl;
837  libmesh_error();
838  }
839 
840  // Compute the proper subdomain ID, based on the formula in the
841  // documentation for this function.
842  subdomain_id_type computed_id = elemset_id + (elem_types_map[elem->type()] * n_elemsets);
843 
844  // Assign this ID to the element in question
845  elem->subdomain_id() = computed_id;
846  }
847  }
848  }
849  } // assign_subdomain_ids()
MeshBase & libMesh::MeshInput< MeshBase >::mesh ( )
protectedinherited

Returns the object as a writeable reference.

Referenced by libMesh::GMVIO::_read_one_cell(), libMesh::VTKIO::cells_to_vtk(), libMesh::TetGenIO::element_in(), libMesh::UNVIO::element_in(), libMesh::UNVIO::element_out(), libMesh::TetGenIO::node_in(), libMesh::UNVIO::node_in(), libMesh::UNVIO::node_out(), libMesh::VTKIO::nodes_to_vtk(), libMesh::ExodusII_IO::read(), libMesh::Nemesis_IO::read(), libMesh::GMVIO::read(), libMesh::CheckpointIO::read(), libMesh::XdrIO::read(), libMesh::VTKIO::read(), libMesh::LegacyXdrIO::read_ascii(), libMesh::CheckpointIO::read_bcs(), libMesh::CheckpointIO::read_connectivity(), libMesh::UCDIO::read_implementation(), libMesh::GmshIO::read_mesh(), libMesh::CheckpointIO::read_nodes(), libMesh::CheckpointIO::read_nodesets(), libMesh::XdrIO::read_serialized_bcs(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::OFFIO::read_stream(), libMesh::MatlabIO::read_stream(), libMesh::CheckpointIO::read_subdomain_names(), libMesh::TetGenIO::write(), libMesh::ExodusII_IO::write(), libMesh::Nemesis_IO::write(), libMesh::CheckpointIO::write(), libMesh::XdrIO::write(), libMesh::GMVIO::write_ascii_new_impl(), libMesh::GMVIO::write_ascii_old_impl(), libMesh::CheckpointIO::write_bcs(), libMesh::GMVIO::write_binary(), libMesh::CheckpointIO::write_connectivity(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::UCDIO::write_implementation(), libMesh::UNVIO::write_implementation(), libMesh::GmshIO::write_mesh(), libMesh::UCDIO::write_nodal_data(), libMesh::VTKIO::write_nodal_data(), libMesh::Nemesis_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_common(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::CheckpointIO::write_nodes(), libMesh::CheckpointIO::write_nodesets(), libMesh::XdrIO::write_parallel(), libMesh::GmshIO::write_post(), libMesh::XdrIO::write_serialized_bcs(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), and libMesh::CheckpointIO::write_subdomain_names().

std::string libMesh::AbaqusIO::parse_label ( std::string  line,
std::string  label_name 
)
private

This function parses a label of the form foo=bar from a comma-delimited line of the form ..., foo=bar, ... The input to the function in this case would be foo, the output would be bar

Definition at line 679 of file abaqus_io.C.

References end.

Referenced by read(), and read_elements().

680  {
681  // Do all string comparisons in upper-case
682  std::string upper_line(line), upper_label_name(label_name);
683  std::transform(upper_line.begin(), upper_line.end(), upper_line.begin(), ::toupper);
684  std::transform(upper_label_name.begin(), upper_label_name.end(), upper_label_name.begin(), ::toupper);
685 
686  // Get index of start of "label="
687  size_t label_index = upper_line.find(upper_label_name + "=");
688 
689  if (label_index != std::string::npos)
690  {
691  // Location of the first comma following "label="
692  size_t comma_index = upper_line.find(",", label_index);
693 
694  // Construct iterators from which to build the sub-string.
695  // Note the +1 is to skip past the "=" which follows the label name
696  std::string::iterator
697  beg = line.begin() + label_name.size() + 1 + label_index,
698  end = (comma_index == std::string::npos) ? line.end() : line.begin()+comma_index;
699 
700  return std::string(beg, end);
701  }
702 
703  // The label index was not found, return the empty string
704  return std::string("");
705  } // parse_label()
void libMesh::AbaqusIO::process_and_discard_comments ( )
private

Any of the various sections can start with some number of lines of comments, which start with "**". This function discards any lines of comments that it finds from the stream, leaving trailing data intact.

Definition at line 958 of file abaqus_io.C.

References _in.

Referenced by read().

959  {
960  std::string dummy;
961  while (true)
962  {
963  // We assume we are at the beginning of a line that may be
964  // comments or may be data. We need to only discard the line if
965  // it begins with **, but we must avoid calling std::getline()
966  // since there's no way to put that back.
967  if (_in.peek() == '*')
968  {
969  // The first character was a star, so actually read it from the stream.
970  _in.get();
971 
972  // Peek at the next character...
973  if (_in.peek() == '*')
974  {
975  // OK, second character was star also, by definition this
976  // line must be a comment! Read the rest of the line and discard!
977  std::getline(_in, dummy);
978 
979  // Debugging:
980  // libMesh::out << "Read comment line: " << dummy << std::endl;
981  }
982  else
983  {
984  // The second character was _not_ a star, so put back the first star
985  // we pulled out so that the line can be parsed correctly by somebody
986  // else!
987  _in.unget();
988 
989  // Finally, break out of the while loop, we are done parsing comments
990  break;
991  }
992  }
993  else
994  {
995  // First character was not *, so this line must be data! Break out of the
996  // while loop!
997  break;
998  }
999  }
1000 
1001  } // process_and_discard_comments()
void libMesh::AbaqusIO::read ( const std::string &  name)
virtual

This method implements reading a mesh from a specified file.

Implements libMesh::MeshInput< MeshBase >.

Definition at line 187 of file abaqus_io.C.

References _already_seen_part, _elemset_ids, _in, _nodeset_ids, _sideset_ids, assign_boundary_node_ids(), assign_sideset_ids(), assign_subdomain_ids(), libMesh::MeshBase::boundary_info, build_sidesets_from_nodesets, libMesh::MeshBase::clear(), libMesh::err, libMesh::libmesh_assert(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshBase::n_elem(), libMesh::MeshBase::n_nodes(), libMesh::out, parse_label(), process_and_discard_comments(), read_elements(), read_ids(), read_nodes(), and read_sideset().

Referenced by libMesh::UnstructuredMesh::read().

188  {
189  // Get a reference to the mesh we are reading
190  MeshBase& the_mesh = MeshInput<MeshBase>::mesh();
191 
192  // Clear any existing mesh data
193  the_mesh.clear();
194 
195  // Open stream for reading
196  _in.open(fname.c_str());
197  libmesh_assert(_in.good());
198 
199  // Read file line-by-line... this is based on a set of different
200  // test input files. I have not looked at the full input file
201  // specs for Abaqus.
202  std::string s;
203  while (true)
204  {
205  // Try to read something. This may set EOF!
206  std::getline(_in, s);
207 
208  if (_in)
209  {
210  // Process s...
211  //
212  // There are many sections in Abaqus files, we read some
213  // but others are just ignored... Some sections may occur
214  // more than once. For example for a hybrid grid, you
215  // will have multiple *Element sections...
216 
217  // Some Abaqus files use all upper-case for section names,
218  // so we will just convert s to uppercase
219  std::string upper(s);
220  std::transform(upper.begin(), upper.end(), upper.begin(), ::toupper);
221 
222  // 0.) Look for the "*Part" section
223  if (upper.find("*PART") == 0)
224  {
225  // libMesh::out << "Found parts section!" << std::endl;
226 
227  if (_already_seen_part)
228  {
229  libMesh::err << "We currently don't support reading Abaqus files with multiple PART sections" << std::endl;
230  libmesh_error();
231  }
232 
233  _already_seen_part = true;
234  }
235 
236  // 1.) Look for the "*Nodes" section
237  if (upper.find("*NODE") == 0)
238  {
239  // Process any lines of comments that may be present
241 
242  // Read a block of nodes
243  this->read_nodes();
244  }
245 
246 
247 
248  // 2.) Look for the "*Element" section
249  else if (upper.find("*ELEMENT,") == 0)
250  {
251  // Process any lines of comments that may be present
253 
254  // Read a block of elements
255  this->read_elements(upper);
256  }
257 
258 
259 
260  // 3.) Look for a Nodeset section
261  else if (upper.find("*NSET") == 0)
262  {
263  std::string nset_name = this->parse_label(s, "nset");
264 
265  // I haven't seen an unnamed elset yet, but let's detect it
266  // just in case...
267  if (nset_name == "")
268  {
269  libMesh::err << "Unnamed nset encountered!" << std::endl;
270  libmesh_error();
271  }
272 
273  // Process any lines of comments that may be present
275 
276  // Read the IDs, storing them in _nodeset_ids
277  this->read_ids(nset_name, _nodeset_ids);
278  } // *Nodeset
279 
280 
281 
282  // 4.) Look for an Elset section
283  else if (upper.find("*ELSET") == 0)
284  {
285  std::string elset_name = this->parse_label(s, "elset");
286 
287  // I haven't seen an unnamed elset yet, but let's detect it
288  // just in case...
289  if (elset_name == "")
290  {
291  libMesh::err << "Unnamed elset encountered!" << std::endl;
292  libmesh_error();
293  }
294 
295  // Debugging
296  // libMesh::out << "Processing ELSET: " << elset_name << std::endl;
297 
298  // Process any lines of comments that may be present
300 
301  // Read the IDs, storing them in _elemset_ids
302  this->read_ids(elset_name, _elemset_ids);
303  } // *Elset
304 
305 
306 
307  // 5.) Look for a Surface section. Need to be a little
308  // careful, since there are also "surface interaction"
309  // sections we don't want to read here.
310  else if (upper.find("*SURFACE,") == 0)
311  {
312  // libMesh::out << "Found SURFACE section: " << s << std::endl;
313 
314  // Get the name from the Name=Foo label. This will be the map key.
315  std::string sideset_name = this->parse_label(s, "name");
316 
317  // Print name of section we just found
318  // libMesh::out << "Found surface section named: " << sideset_name << std::endl;
319 
320  // Process any lines of comments that may be present
322 
323  // Read the sideset IDs
324  this->read_sideset(sideset_name, _sideset_ids);
325 
326  // Debugging: print status of most recently read sideset
327  // libMesh::out << "Read " << _sideset_ids[sideset_name].size() << " sides in " << sideset_name << std::endl;
328  }
329 
330  continue;
331  } // if (_in)
332 
333  // If !file, check to see if EOF was set. If so, break out
334  // of while loop.
335  if (_in.eof())
336  break;
337 
338  // If !in and !in.eof(), stream is in a bad state!
339  libMesh::err << "Stream is bad!\n";
340  libMesh::err << "Perhaps the file: " << fname << " does not exist?" << std::endl;
341  libmesh_error();
342  } // while
343 
344 
345  //
346  // We've read everything we can from the file at this point. Now
347  // do some more processing.
348  //
349  libMesh::out << "Mesh contains "
350  << the_mesh.n_elem()
351  << " elements, and "
352  << the_mesh.n_nodes()
353  << " nodes." << std::endl;
354 
355  // TODO: Remove these or write a function to do it?
356 // {
357 // container_t::iterator it=_nodeset_ids.begin();
358 // for (; it != _nodeset_ids.end(); ++it)
359 // {
360 // libMesh::out << "Node set '" << (*it).first << "' contains " << (*it).second.size() << " ID(s)." << std::endl;
361 // }
362 // }
363 //
364 // {
365 // container_t::iterator it=_elemset_ids.begin();
366 // for (; it != _elemset_ids.end(); ++it)
367 // {
368 // libMesh::out << "Elem set '" << (*it).first << "' contains " << (*it).second.size() << " ID(s)." << std::endl;
369 // }
370 // }
371 
372 
373  // Set element IDs based on the element sets.
374  this->assign_subdomain_ids();
375 
376  // Assign nodeset values to the BoundaryInfo object
377  this->assign_boundary_node_ids();
378 
379  // Assign sideset values in the BoundaryInfo object
380  this->assign_sideset_ids();
381 
382  // Abaqus files only contain nodesets by default. To be useful in
383  // applying most types of BCs in libmesh, we will definitely need
384  // sidesets. So we can call the new BoundaryInfo function which
385  // generates sidesets from nodesets.
387  the_mesh.boundary_info->build_side_list_from_node_list();
388 
389  } // read()
void libMesh::AbaqusIO::read_elements ( std::string  upper)
private

This function parses a block of elements in the Abaqus file. You must pass it an upper-cased version of the string declaring this section, which is typically something like: *ELEMENT, TYPE=CPS3 so that it can determine the type of elements to read.

Definition at line 480 of file abaqus_io.C.

References _abaqus_to_libmesh_elem_mapping, _abaqus_to_libmesh_node_mapping, _elem_types, _elemset_ids, _in, libMesh::MeshBase::add_elem(), libMesh::Elem::build(), libMesh::err, libMeshEnums::HEX20, libMeshEnums::HEX8, libMesh::DofObject::id(), libMeshEnums::INVALID_ELEM, libMesh::MeshInput< MT >::mesh(), libMesh::MeshBase::node_ptr(), parse_label(), libMeshEnums::PRISM15, libMeshEnums::PRISM6, libMeshEnums::QUAD4, libMesh::MeshBase::set_mesh_dimension(), libMesh::Elem::set_node(), libMeshEnums::TET10, libMeshEnums::TET4, and libMeshEnums::TRI3.

Referenced by read().

481  {
482  // Some *Element sections also specify an Elset name on the same line.
483  // Look for one here.
484  std::string elset_name = this->parse_label(upper, "elset");
485 
486  // Get a reference to the mesh we are reading
487  MeshBase& the_mesh = MeshInput<MeshBase>::mesh();
488 
489  // initialize the eletypes map (eletypes is a file-global variable)
490  init_eletypes();
491 
492  ElemType elem_type = INVALID_ELEM;
493  unsigned n_nodes_per_elem = 0;
494 
495  // Within s, we should have "type=XXXX"
496  if (upper.find("CPE4") != std::string::npos ||
497  upper.find("CPS4") != std::string::npos)
498  {
499  elem_type = QUAD4;
500  n_nodes_per_elem = 4;
501  the_mesh.set_mesh_dimension(2);
502  }
503  else if (upper.find("CPS3") != std::string::npos)
504  {
505  elem_type = TRI3;
506  n_nodes_per_elem = 3;
507  the_mesh.set_mesh_dimension(2);
508  }
509  else if (upper.find("C3D8") != std::string::npos)
510  {
511  elem_type = HEX8;
512  n_nodes_per_elem = 8;
513  the_mesh.set_mesh_dimension(3);
514  }
515  else if (upper.find("C3D4") != std::string::npos)
516  {
517  elem_type = TET4;
518  n_nodes_per_elem = 4;
519  the_mesh.set_mesh_dimension(3);
520  }
521  else if (upper.find("C3D20") != std::string::npos)
522  {
523  elem_type = HEX20;
524  n_nodes_per_elem = 20;
525  the_mesh.set_mesh_dimension(3);
526  }
527  else if (upper.find("C3D6") != std::string::npos)
528  {
529  elem_type = PRISM6;
530  n_nodes_per_elem = 6;
531  the_mesh.set_mesh_dimension(3);
532  }
533  else if (upper.find("C3D15") != std::string::npos)
534  {
535  elem_type = PRISM15;
536  n_nodes_per_elem = 15;
537  the_mesh.set_mesh_dimension(3);
538  }
539  else if (upper.find("C3D10") != std::string::npos)
540  {
541  elem_type = TET10;
542  n_nodes_per_elem = 10;
543  the_mesh.set_mesh_dimension(3);
544  }
545  else
546  {
547  libMesh::err << "Unrecognized element type: " << upper << std::endl;
548  libmesh_error();
549  }
550 
551  // Insert the elem type we detected into the set of all elem types for this mesh
552  _elem_types.insert(elem_type);
553 
554  // For reading in line endings
555  std::string dummy;
556 
557  // Grab a reference to the element definition for this element type
558  const ElementDefinition& eledef = eletypes[elem_type];
559 
560  // If the element definition was not found, the call above would have
561  // created one with an uninitialized struct. Check for that here...
562  if (eledef.abaqus_zero_based_node_id_to_libmesh_node_id.size() == 0)
563  {
564  // libMesh::err << "No Abaqus->LibMesh mapping information for ElemType "
565  // << Utility::enum_to_string(elem_type)
566  // << "!"
567  // << std::endl;
568  libmesh_error();
569  }
570 
571  // We will read elements until the next line begins with *, since that will be the
572  // next section.
573  while (_in.peek() != '*' && _in.peek() != EOF)
574  {
575  // Read the element ID, it is the first number on each line. It is
576  // followed by a comma, so read that also. We will need this ID later
577  // when we try to assign subdomain IDs
578  dof_id_type abaqus_elem_id = 0;
579  char c;
580  _in >> abaqus_elem_id >> c;
581 
582  // Debugging:
583  // libMesh::out << "Reading data for element " << abaqus_elem_id << std::endl;
584 
585  // Add an element of the appropriate type to the Mesh.
586  Elem* elem = the_mesh.add_elem(Elem::build(elem_type).release());
587 
588  // Associate the ID returned from libmesh with the abaqus element ID
589  //_libmesh_to_abaqus_elem_mapping[elem->id()] = abaqus_elem_id;
590  _abaqus_to_libmesh_elem_mapping[abaqus_elem_id] = elem->id();
591 
592  // The count of the total number of IDs read for the current element.
593  unsigned id_count=0;
594 
595  // Continue reading line-by-line until we have read enough nodes for this element
596  while (id_count < n_nodes_per_elem)
597  {
598  // Read entire line (up to carriage return) of comma-separated values
599  std::string csv_line;
600  std::getline(_in, csv_line);
601 
602  // Create a stream object out of the current line
603  std::stringstream line_stream(csv_line);
604 
605  // Process the comma-separated values
606  std::string cell;
607  while (std::getline(line_stream, cell, ','))
608  {
609  // FIXME: factor out this strtol stuff into a utility function.
610  char* endptr;
611  long abaqus_global_node_id = std::strtol(cell.c_str(), &endptr, /*base=*/10);
612 
613  if (abaqus_global_node_id!=0 || cell.c_str() != endptr)
614  {
615  // Use the global node number mapping to determine the corresponding libmesh global node id
616  dof_id_type libmesh_global_node_id = _abaqus_to_libmesh_node_mapping[abaqus_global_node_id];
617 
618  // Grab the node pointer from the mesh for this ID
619  Node* node = the_mesh.node_ptr(libmesh_global_node_id);
620 
621  // Debugging:
622  // libMesh::out << "Assigning global node id: " << abaqus_global_node_id
623  // << "(Abaqus), " << node->id() << "(LibMesh)" << std::endl;
624 
625  // If node_ptr() returns NULL, it may mean we have not yet read the
626  // *Nodes section, though I assumed that always came before the *Elements section...
627  if (node == NULL)
628  {
629  libMesh::err << "Error! Mesh returned NULL Node pointer.\n";
630  libMesh::err << "Either no node exists with ID " << libmesh_global_node_id
631  << " or perhaps this input file has *Elements defined before *Nodes?" << std::endl;
632  libmesh_error();
633  }
634 
635  // Note: id_count is the zero-based abaqus (elem local) node index. We therefore map
636  // it to a libmesh elem local node index using the element definition map
637  unsigned libmesh_elem_local_node_id =
638  eledef.abaqus_zero_based_node_id_to_libmesh_node_id[id_count];
639 
640  // Set this node pointer within the element.
641  elem->set_node(libmesh_elem_local_node_id) = node;
642 
643  // Debugging:
644  // libMesh::out << "Setting elem " << elem->id()
645  // << ", local node " << libmesh_elem_local_node_id
646  // << " to global node " << node->id() << std::endl;
647 
648  // Increment the count of IDs read for this element
649  id_count++;
650  } // end if strtol success
651  } // end while getline(',')
652  } // end while (id_count)
653 
654  // Ensure that we read *exactly* as many nodes as we were expecting to, no more.
655  if (id_count != n_nodes_per_elem)
656  {
657  libMesh::err << "Error: Needed to read "
658  << n_nodes_per_elem
659  << " nodes, but read "
660  << id_count
661  << " instead!" << std::endl;
662  libmesh_error();
663  }
664 
665  // If we are recording Elset IDs, add this element to the correct set for later processing.
666  // Make sure to add it with the Abaqus ID, not the libmesh one!
667  if (elset_name != "")
668  {
669  // Debugging:
670  // libMesh::out << "Adding Elem " << abaqus_elem_id << " to Elmset " << elset_name << std::endl;
671  _elemset_ids[elset_name].push_back(abaqus_elem_id);
672  }
673  } // end while (peek)
674  } // read_elements()
void libMesh::AbaqusIO::read_ids ( std::string  set_name,
container_t container 
)
private

This function reads all the IDs for the current node or element set of the given name, storing them in the passed map using the name as key.

Definition at line 710 of file abaqus_io.C.

References _in.

Referenced by read().

711  {
712  // Debugging
713  // libMesh::out << "Reading ids for set: " << set_name << std::endl;
714 
715  // Grab a reference to a vector that will hold all the IDs
716  std::vector<dof_id_type>& id_storage = container[set_name];
717 
718  // Read until the start of another section is detected, or EOF is encountered
719  while (_in.peek() != '*' && _in.peek() != EOF)
720  {
721  // Read entire comma-separated line into a string
722  std::string csv_line;
723  std::getline(_in, csv_line);
724 
725  // On that line, use std::getline again to parse each
726  // comma-separated entry.
727  std::string cell;
728  std::stringstream line_stream(csv_line);
729  while (std::getline(line_stream, cell, ','))
730  {
731  // If no conversion can be performed by strtol, 0 is returned.
732  //
733  // If endptr is not NULL, strtol() stores the address of the
734  // first invalid character in *endptr. If there were no
735  // digits at all, however, strtol() stores the original
736  // value of str in *endptr.
737  char* endptr;
738 
739  // FIXME - this needs to be updated for 64-bit inputs
740  long id = std::strtol(cell.c_str(), &endptr, /*base=*/10);
741 
742  // Note that lists of comma-separated values in abaqus also
743  // *end* with a comma, so the last call to getline on a given
744  // line will get an empty string, which we must detect.
745  if (id!=0 || cell.c_str() != endptr)
746  {
747  // Debugging
748  // libMesh::out << "Read id: " << id << std::endl;
749 
750  // 'cell' is now a string with an integer id in it
751  id_storage.push_back( id );
752  }
753  }
754  }
755 
756  // Status message
757  // libMesh::out << "Read " << id_storage.size() << " ID(s) for the set " << set_name << std::endl;
758  } // read_ids()
void libMesh::AbaqusIO::read_nodes ( )
private

This function parses a block of nodes in the Abaqus file once such a block has been found.

Definition at line 397 of file abaqus_io.C.

References _abaqus_to_libmesh_node_mapping, _in, libMesh::MeshBase::add_point(), libMesh::MeshInput< MT >::mesh(), and libMesh::Real.

Referenced by read().

398  {
399  // Get a reference to the mesh we are reading
400  MeshBase& the_mesh = MeshInput<MeshBase>::mesh();
401 
402  // Debugging: print node count
403  // libMesh::out << "Before read_nodes(), mesh contains "
404  // << the_mesh.n_elem()
405  // << " elements, and "
406  // << the_mesh.n_nodes()
407  // << " nodes." << std::endl;
408 
409  // In the input file I have, Abaqus neither tells what
410  // the mesh dimension is nor how many nodes it has...
411 
412  // The node line format is:
413  // id, x, y, z
414  // and you do have to parse out the commas.
415  // The z-coordinate will only be present for 3D meshes
416 
417  // Temporary variables for parsing lines of text
418  dof_id_type abaqus_node_id=0;
419  Real x=0, y=0, z=0;
420  char c;
421  std::string dummy;
422 
423  // Defines the sequential node numbering used by libmesh
424  dof_id_type libmesh_node_id = 0;
425 
426  // We will read nodes until the next line begins with *, since that will be the
427  // next section.
428  // TODO: Is Abaqus guaranteed to start the line with '*' or can there be leading white space?
429  while (_in.peek() != '*' && _in.peek() != EOF)
430  {
431  // Re-Initialize variables to be read in from file
432  abaqus_node_id=0;
433  x = y = z = 0.;
434 
435  // Note: we assume *at least* 2D points here, should we worry about
436  // trying to read 1D Abaqus meshes?
437  _in >> abaqus_node_id >> c >> x >> c >> y;
438 
439  // Peek at the next character. If it is a comma, then there is another
440  // value to read!
441  if (_in.peek() == ',')
442  _in >> c >> z;
443 
444  // Debugging: Print what we just read in.
445  // libMesh::out << "node_id=" << node_id
446  // << ", x=" << x
447  // << ", y=" << y
448  // << ", z=" << z
449  // << std::endl;
450 
451  // Read (and discard) the rest of the line, including the newline.
452  // This is required so that our 'peek()' at the beginning of this
453  // loop doesn't read the newline character, for example.
454  std::getline(_in, dummy);
455 
456  // Set up the abaqus -> libmesh node mapping. This is usually just the
457  // "off-by-one" map.
458  _abaqus_to_libmesh_node_mapping[abaqus_node_id] = libmesh_node_id;
459 
460  // Add the point to the mesh using libmesh's numbering,
461  // and post-increment the libmesh node counter.
462  the_mesh.add_point(Point(x,y,z), libmesh_node_id++);
463  } // while
464 
465  // Debugging: print node count. Note: in serial mesh, this count may
466  // be off by one, since Abaqus uses one-based numbering, and libmesh
467  // just reports the length of its _nodes vector for the number of nodes.
468  // libMesh::out << "After read_nodes(), mesh contains "
469  // << the_mesh.n_elem()
470  // << " elements, and "
471  // << the_mesh.n_nodes()
472  // << " nodes." << std::endl;
473 
474  } // read_nodes()
void libMesh::AbaqusIO::read_sideset ( std::string  sideset_name,
sideset_container_t container 
)
private

This function reads a sideset from the input file. This is defined by a "*Surface" section in the file, and then a list of element ID and side IDs for the set.

Definition at line 763 of file abaqus_io.C.

References _in.

Referenced by read().

764  {
765  // Grab a reference to a vector that will hold all the IDs
766  std::vector<std::pair<dof_id_type, unsigned> >& id_storage = container[sideset_name];
767 
768  // Variables for storing values read in from file
769  dof_id_type elem_id=0;
770  unsigned side_id=0;
771  char c;
772  std::string dummy;
773 
774  // Read until the start of another section is detected, or EOF is encountered
775  while (_in.peek() != '*' && _in.peek() != EOF)
776  {
777  // The strings are of the form: "391, S2"
778 
779  // Read the element ID and the leading comma
780  _in >> elem_id >> c;
781 
782  // Read another character (the 'S') and finally the side ID
783  _in >> c >> side_id;
784 
785  // Debugging: print status
786  // libMesh::out << "Read elem_id=" << elem_id << ", side_id=" << side_id << std::endl;
787 
788  // Store this pair of data in the vector
789  id_storage.push_back( std::make_pair(elem_id, side_id) );
790 
791  // Extract remaining characters on line including newline
792  std::getline(_in, dummy);
793  } // while
794  }
void libMesh::MeshInput< MeshBase >::skip_comment_lines ( std::istream &  in,
const char  comment_start 
)
protectedinherited

Reads input from in, skipping all the lines that start with the character comment_start.

Referenced by libMesh::TetGenIO::read(), and libMesh::UCDIO::read_implementation().

Member Data Documentation

std::map<dof_id_type, dof_id_type> libMesh::AbaqusIO::_abaqus_to_libmesh_elem_mapping
private

Map from libmesh element number -> abaqus element number, and the converse.

Definition at line 185 of file abaqus_io.h.

Referenced by assign_sideset_ids(), assign_subdomain_ids(), and read_elements().

std::map<dof_id_type, dof_id_type> libMesh::AbaqusIO::_abaqus_to_libmesh_node_mapping
private

Map from abaqus node number -> sequential, 0-based libmesh node numbering. Note that in every Abaqus file I've ever seen the node numbers were 1-based, sequential, and all in order, so that this map is probably overkill. Nevertheless, it is the most general solution in case we come across a weird Abaqus file some day.

Definition at line 194 of file abaqus_io.h.

Referenced by assign_boundary_node_ids(), read_elements(), and read_nodes().

bool libMesh::AbaqusIO::_already_seen_part
private

This flag gets set to true after the first "*PART" section we see. If it is still true when we see a second PART section, we will print an error message... we don't currently handle input files with multiple parts.

Definition at line 202 of file abaqus_io.h.

Referenced by read().

std::set<ElemType> libMesh::AbaqusIO::_elem_types
private

A set of the different geometric element types detected when reading the mesh.

Definition at line 178 of file abaqus_io.h.

Referenced by assign_subdomain_ids(), and read_elements().

container_t libMesh::AbaqusIO::_elemset_ids
private

Definition at line 166 of file abaqus_io.h.

Referenced by assign_subdomain_ids(), read(), and read_elements().

std::ifstream libMesh::AbaqusIO::_in
private

Stream object used to interact with the file

Definition at line 172 of file abaqus_io.h.

Referenced by process_and_discard_comments(), read(), read_elements(), read_ids(), read_nodes(), and read_sideset().

container_t libMesh::AbaqusIO::_nodeset_ids
private

Abaqus writes nodesets and elemsets with labels. As we read them in, we'll use these maps to provide a natural ordering for them.

Definition at line 165 of file abaqus_io.h.

Referenced by assign_boundary_node_ids(), and read().

sideset_container_t libMesh::AbaqusIO::_sideset_ids
private

Definition at line 167 of file abaqus_io.h.

Referenced by assign_sideset_ids(), and read().

bool libMesh::AbaqusIO::build_sidesets_from_nodesets

Default false. Abaqus files have only nodesets in them by default. Set this flag to true if you want libmesh to automatically generate sidesets from Abaqus' nodesets.

Definition at line 62 of file abaqus_io.h.

Referenced by read().


The documentation for this class was generated from the following files:

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:58 UTC

Hosted By:
SourceForge.net Logo