libMesh::UNVIO Class Reference

#include <unv_io.h>

Inheritance diagram for libMesh::UNVIO:

Public Member Functions

 UNVIO (MeshBase &mesh)
 
 UNVIO (const MeshBase &mesh)
 
virtual ~UNVIO ()
 
virtual void read (const std::string &) libmesh_override
 
virtual void write (const std::string &) libmesh_override
 
bool & verbose ()
 
void read_dataset (std::string file_name)
 
const std::vector< Number > * get_data (Node *node) const
 
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
 
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
 
virtual void write_nodal_data (const std::string &, const NumericVector< Number > &, const std::vector< std::string > &)
 
unsigned int & ascii_precision ()
 

Protected Member Functions

MeshBasemesh ()
 
void set_n_partitions (unsigned int n_parts)
 
void skip_comment_lines (std::istream &in, const char comment_start)
 
const MeshBasemesh () const
 

Protected Attributes

std::vector< bool > elems_of_dimension
 
const bool _is_parallel_format
 
const bool _serial_only_needed_on_proc_0
 

Private Member Functions

void read_implementation (std::istream &in_stream)
 
void write_implementation (std::ostream &out_stream)
 
void nodes_in (std::istream &in_file)
 
void elements_in (std::istream &in_file)
 
void groups_in (std::istream &in_file)
 
void nodes_out (std::ostream &out_file)
 
void elements_out (std::ostream &out_file)
 
unsigned char max_elem_dimension_seen ()
 
bool need_D_to_e (std::string &number)
 

Private Attributes

bool _verbose
 
std::map< dof_id_type, Node * > _unv_node_id_to_libmesh_node_ptr
 
std::map< unsigned, unsigned > _unv_elem_id_to_libmesh_elem_id
 
std::map< Node *, std::vector< Number > > _node_data
 

Static Private Attributes

static const std::string _nodes_dataset_label = "2411"
 
static const std::string _elements_dataset_label = "2412"
 
static const std::string _groups_dataset_label = "2467"
 

Detailed Description

The UNVIO class implements the Ideas UNV universal file format. This class enables both reading and writing UNV files.

Author history

Author
Tammo Kaschner
Daniel Dreyer
Benjamin S. Kirk
John W. Peterson
Date
2003, 2004, 2014

Definition at line 52 of file unv_io.h.

Constructor & Destructor Documentation

libMesh::UNVIO::UNVIO ( MeshBase mesh)

Constructor. Takes a writable reference to a mesh object. This is the constructor required to read a mesh.

Definition at line 64 of file unv_io.C.

64  :
65  MeshInput<MeshBase> (mesh),
66  MeshOutput<MeshBase>(mesh),
67  _verbose (false)
68 {
69 }
bool _verbose
Definition: unv_io.h:182
libMesh::UNVIO::UNVIO ( const MeshBase mesh)

Constructor. Takes a reference to a constant mesh object. This constructor will only allow us to write the mesh.

Definition at line 73 of file unv_io.C.

73  :
74  MeshOutput<MeshBase> (mesh),
75  _verbose (false)
76 {
77 }
bool _verbose
Definition: unv_io.h:182
libMesh::UNVIO::~UNVIO ( )
virtual

Destructor.

Definition at line 81 of file unv_io.C.

82 {
83 }

Member Function Documentation

unsigned int& libMesh::MeshOutput< MeshBase >::ascii_precision ( )
inherited

Return/set the precision to use when writing ASCII files.

By default we use numeric_limits<Real>::digits10 + 2, which should be enough to write out to ASCII and get the exact same Real back when reading in.

Referenced by libMesh::TecplotIO::write_ascii(), libMesh::GMVIO::write_ascii_new_impl(), and libMesh::GMVIO::write_ascii_old_impl().

void libMesh::UNVIO::elements_in ( std::istream &  in_file)
private

Method reads elements and stores them in std::vector<Elem *> _elements in the same order as they come in. Within UNVIO, element labels are ignored.

Definition at line 599 of file unv_io.C.

References _unv_elem_id_to_libmesh_elem_id, _unv_node_id_to_libmesh_node_ptr, libMesh::MeshBase::add_elem(), libMesh::Elem::dim(), libMesh::MeshInput< MeshBase >::elems_of_dimension, libmesh_nullptr, libMesh::MeshInput< MT >::mesh(), libMesh::MeshInput< MeshBase >::mesh(), n_nodes, libMesh::out, libMesh::DofObject::set_id(), libMesh::Elem::set_node(), and verbose().

Referenced by read_implementation().

600 {
601  LOG_SCOPE("elements_in()","UNVIO");
602 
603  if (this->verbose())
604  libMesh::out << " Reading elements" << std::endl;
605 
606  MeshBase & mesh = MeshInput<MeshBase>::mesh();
607 
608  // node label, we use an int here so we can read in a -1
609  int element_label;
610 
611  unsigned
612  n_nodes, // number of nodes on element
613  fe_descriptor_id, // FE descriptor id
614  phys_prop_tab_num, // physical property table number (not supported yet)
615  mat_prop_tab_num, // material property table number (not supported yet)
616  color; // color (not supported yet)
617 
618  // vector that temporarily holds the node labels defining element
619  std::vector<unsigned int> node_labels (21);
620 
621  // vector that assigns element nodes to their correct position
622  // for example:
623  // 44:plane stress | QUAD4
624  // linear quadrilateral |
625  // position in UNV-file | position in libmesh
626  // assign_elem_node[1] = 0
627  // assign_elem_node[2] = 3
628  // assign_elem_node[3] = 2
629  // assign_elem_node[4] = 1
630  //
631  // UNV is 1-based, we leave the 0th element of the vectors unused in order
632  // to prevent confusion, this way we can store elements with up to 20 nodes
633  unsigned int assign_elem_nodes[21];
634 
635  // Read elements until there are none left
636  unsigned ctr = 0;
637  while (true)
638  {
639  // read element label, break out when we read -1
640  in_file >> element_label;
641 
642  if (element_label == -1)
643  break;
644 
645  in_file >> fe_descriptor_id // read FE descriptor id
646  >> phys_prop_tab_num // (not supported yet)
647  >> mat_prop_tab_num // (not supported yet)
648  >> color // (not supported yet)
649  >> n_nodes; // read number of nodes on element
650 
651  // For "beam" type elements, the next three numbers are:
652  // .) beam orientation node number
653  // .) beam fore-end cross section number
654  // .) beam aft-end cross section number
655  // which we discard in libmesh. The "beam" type elements:
656  // 11 Rod
657  // 21 Linear beam
658  // 22 Tapered beam
659  // 23 Curved beam
660  // 24 Parabolic beam
661  // all have fe_descriptor_id < 25.
662  // http://www.sdrl.uc.edu/universal-file-formats-for-modal-analysis-testing-1/file-format-storehouse/unv_2412.htm
663  if (fe_descriptor_id < 25)
664  {
665  unsigned dummy;
666  in_file >> dummy >> dummy >> dummy;
667  }
668 
669  // read node labels (1-based)
670  for (unsigned int j=1; j<=n_nodes; j++)
671  in_file >> node_labels[j];
672 
673  // element pointer, to be allocated
674  Elem * elem = libmesh_nullptr;
675 
676  switch (fe_descriptor_id)
677  {
678  case 11: // Rod
679  {
680  elem = new Edge2;
681 
682  assign_elem_nodes[1]=0;
683  assign_elem_nodes[2]=1;
684  break;
685  }
686 
687  case 41: // Plane Stress Linear Triangle
688  case 91: // Thin Shell Linear Triangle
689  {
690  elem = new Tri3; // create new element
691 
692  assign_elem_nodes[1]=0;
693  assign_elem_nodes[2]=2;
694  assign_elem_nodes[3]=1;
695  break;
696  }
697 
698  case 42: // Plane Stress Quadratic Triangle
699  case 92: // Thin Shell Quadratic Triangle
700  {
701  elem = new Tri6; // create new element
702 
703  assign_elem_nodes[1]=0;
704  assign_elem_nodes[2]=5;
705  assign_elem_nodes[3]=2;
706  assign_elem_nodes[4]=4;
707  assign_elem_nodes[5]=1;
708  assign_elem_nodes[6]=3;
709  break;
710  }
711 
712  case 43: // Plane Stress Cubic Triangle
713  libmesh_error_msg("ERROR: UNV-element type 43: Plane Stress Cubic Triangle not supported.");
714 
715  case 44: // Plane Stress Linear Quadrilateral
716  case 94: // Thin Shell Linear Quadrilateral
717  {
718  elem = new Quad4; // create new element
719 
720  assign_elem_nodes[1]=0;
721  assign_elem_nodes[2]=3;
722  assign_elem_nodes[3]=2;
723  assign_elem_nodes[4]=1;
724  break;
725  }
726 
727  case 45: // Plane Stress Quadratic Quadrilateral
728  case 95: // Thin Shell Quadratic Quadrilateral
729  {
730  elem = new Quad8; // create new element
731 
732  assign_elem_nodes[1]=0;
733  assign_elem_nodes[2]=7;
734  assign_elem_nodes[3]=3;
735  assign_elem_nodes[4]=6;
736  assign_elem_nodes[5]=2;
737  assign_elem_nodes[6]=5;
738  assign_elem_nodes[7]=1;
739  assign_elem_nodes[8]=4;
740  break;
741  }
742 
743  case 300: // Thin Shell Quadratic Quadrilateral (nine nodes)
744  {
745  elem = new Quad9; // create new element
746 
747  assign_elem_nodes[1]=0;
748  assign_elem_nodes[2]=7;
749  assign_elem_nodes[3]=3;
750  assign_elem_nodes[4]=6;
751  assign_elem_nodes[5]=2;
752  assign_elem_nodes[6]=5;
753  assign_elem_nodes[7]=1;
754  assign_elem_nodes[8]=4;
755  assign_elem_nodes[9]=8;
756  break;
757  }
758 
759  case 46: // Plane Stress Cubic Quadrilateral
760  libmesh_error_msg("ERROR: UNV-element type 46: Plane Stress Cubic Quadrilateral not supported.");
761 
762  case 111: // Solid Linear Tetrahedron
763  {
764  elem = new Tet4; // create new element
765 
766  assign_elem_nodes[1]=0;
767  assign_elem_nodes[2]=1;
768  assign_elem_nodes[3]=2;
769  assign_elem_nodes[4]=3;
770  break;
771  }
772 
773  case 112: // Solid Linear Prism
774  {
775  elem = new Prism6; // create new element
776 
777  assign_elem_nodes[1]=0;
778  assign_elem_nodes[2]=1;
779  assign_elem_nodes[3]=2;
780  assign_elem_nodes[4]=3;
781  assign_elem_nodes[5]=4;
782  assign_elem_nodes[6]=5;
783  break;
784  }
785 
786  case 115: // Solid Linear Brick
787  {
788  elem = new Hex8; // create new element
789 
790  assign_elem_nodes[1]=0;
791  assign_elem_nodes[2]=4;
792  assign_elem_nodes[3]=5;
793  assign_elem_nodes[4]=1;
794  assign_elem_nodes[5]=3;
795  assign_elem_nodes[6]=7;
796  assign_elem_nodes[7]=6;
797  assign_elem_nodes[8]=2;
798  break;
799  }
800 
801  case 116: // Solid Quadratic Brick
802  {
803  elem = new Hex20; // create new element
804 
805  assign_elem_nodes[1]=0;
806  assign_elem_nodes[2]=12;
807  assign_elem_nodes[3]=4;
808  assign_elem_nodes[4]=16;
809  assign_elem_nodes[5]=5;
810  assign_elem_nodes[6]=13;
811  assign_elem_nodes[7]=1;
812  assign_elem_nodes[8]=8;
813 
814  assign_elem_nodes[9]=11;
815  assign_elem_nodes[10]=19;
816  assign_elem_nodes[11]=17;
817  assign_elem_nodes[12]=9;
818 
819  assign_elem_nodes[13]=3;
820  assign_elem_nodes[14]=15;
821  assign_elem_nodes[15]=7;
822  assign_elem_nodes[16]=18;
823  assign_elem_nodes[17]=6;
824  assign_elem_nodes[18]=14;
825  assign_elem_nodes[19]=2;
826  assign_elem_nodes[20]=10;
827  break;
828  }
829 
830  case 117: // Solid Cubic Brick
831  libmesh_error_msg("Error: UNV-element type 117: Solid Cubic Brick not supported.");
832 
833  case 118: // Solid Quadratic Tetrahedron
834  {
835  elem = new Tet10; // create new element
836 
837  assign_elem_nodes[1]=0;
838  assign_elem_nodes[2]=4;
839  assign_elem_nodes[3]=1;
840  assign_elem_nodes[4]=5;
841  assign_elem_nodes[5]=2;
842  assign_elem_nodes[6]=6;
843  assign_elem_nodes[7]=7;
844  assign_elem_nodes[8]=8;
845  assign_elem_nodes[9]=9;
846  assign_elem_nodes[10]=3;
847  break;
848  }
849 
850  default: // Unrecognized element type
851  libmesh_error_msg("ERROR: UNV-element type " << fe_descriptor_id << " not supported.");
852  }
853 
854  // nodes are being stored in element
855  for (dof_id_type j=1; j<=n_nodes; j++)
856  {
857  // Map the UNV node ID to the libmesh node ID
858  std::map<dof_id_type, Node *>::iterator it =
859  _unv_node_id_to_libmesh_node_ptr.find(node_labels[j]);
860 
861  if (it != _unv_node_id_to_libmesh_node_ptr.end())
862  elem->set_node(assign_elem_nodes[j]) = it->second;
863  else
864  libmesh_error_msg("ERROR: UNV node " << node_labels[j] << " not found!");
865  }
866 
867  elems_of_dimension[elem->dim()] = true;
868 
869  // Set the element's ID
870  elem->set_id(ctr);
871 
872  // Maintain a map from the libmesh (0-based) numbering to the
873  // UNV numbering.
874  _unv_elem_id_to_libmesh_elem_id[element_label] = ctr;
875 
876  // Add the element to the Mesh
877  mesh.add_elem(elem);
878 
879  // Increment the counter for the next iteration
880  ctr++;
881  } // end while (true)
882 }
bool & verbose()
Definition: unv_io.C:87
std::map< dof_id_type, Node * > _unv_node_id_to_libmesh_node_ptr
Definition: unv_io.h:188
std::vector< bool > elems_of_dimension
Definition: mesh_input.h:97
const class libmesh_nullptr_t libmesh_nullptr
const dof_id_type n_nodes
Definition: tecplot_io.C:67
std::map< unsigned, unsigned > _unv_elem_id_to_libmesh_elem_id
Definition: unv_io.h:208
OStreamProxy out(std::cout)
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::UNVIO::elements_out ( std::ostream &  out_file)
private

Outputs the element data to the file out_file. Do not use this directly, but through the proper write method.

Definition at line 956 of file unv_io.C.

References _elements_dataset_label, libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), end, libMesh::HEX20, libMesh::HEX8, libMesh::DofObject::id(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::Elem::n_nodes(), libMesh::Elem::node_ptr(), libMesh::out, libMesh::PRISM6, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI6, libMesh::Elem::type(), and verbose().

Referenced by write_implementation().

957 {
958  if (this->verbose())
959  libMesh::out << " Writing elements" << std::endl;
960 
961  // Write beginning of dataset
962  out_file << " -1\n"
963  << " "
965  << "\n";
966 
967  unsigned
968  fe_descriptor_id = 0, // FE descriptor id
969  phys_prop_tab_dummy = 2, // physical property (not supported yet)
970  mat_prop_tab_dummy = 1, // material property (not supported yet)
971  color_dummy = 7; // color (not supported yet)
972 
973  // vector that assigns element nodes to their correct position
974  // currently only elements with up to 20 nodes
975  //
976  // Example:
977  // QUAD4 | 44:plane stress
978  // | linear quad
979  // position in libMesh | UNV numbering
980  // (note: 0-based) | (note: 1-based)
981  //
982  // assign_elem_node[0] = 0
983  // assign_elem_node[1] = 3
984  // assign_elem_node[2] = 2
985  // assign_elem_node[3] = 1
986  unsigned int assign_elem_nodes[20];
987 
988  unsigned int n_elem_written=0;
989 
990  // A reference to the parent class's mesh
991  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
992 
993  MeshBase::const_element_iterator it = mesh.elements_begin();
994  const MeshBase::const_element_iterator end = mesh.elements_end();
995 
996  for (; it != end; ++it)
997  {
998  const Elem * elem = *it;
999 
1000  elem->n_nodes();
1001 
1002  switch (elem->type())
1003  {
1004 
1005  case TRI3:
1006  {
1007  fe_descriptor_id = 41; // Plane Stress Linear Triangle
1008  assign_elem_nodes[0] = 0;
1009  assign_elem_nodes[1] = 2;
1010  assign_elem_nodes[2] = 1;
1011  break;
1012  }
1013 
1014  case TRI6:
1015  {
1016  fe_descriptor_id = 42; // Plane Stress Quadratic Triangle
1017  assign_elem_nodes[0] = 0;
1018  assign_elem_nodes[1] = 5;
1019  assign_elem_nodes[2] = 2;
1020  assign_elem_nodes[3] = 4;
1021  assign_elem_nodes[4] = 1;
1022  assign_elem_nodes[5] = 3;
1023  break;
1024  }
1025 
1026  case QUAD4:
1027  {
1028  fe_descriptor_id = 44; // Plane Stress Linear Quadrilateral
1029  assign_elem_nodes[0] = 0;
1030  assign_elem_nodes[1] = 3;
1031  assign_elem_nodes[2] = 2;
1032  assign_elem_nodes[3] = 1;
1033  break;
1034  }
1035 
1036  case QUAD8:
1037  {
1038  fe_descriptor_id = 45; // Plane Stress Quadratic Quadrilateral
1039  assign_elem_nodes[0] = 0;
1040  assign_elem_nodes[1] = 7;
1041  assign_elem_nodes[2] = 3;
1042  assign_elem_nodes[3] = 6;
1043  assign_elem_nodes[4] = 2;
1044  assign_elem_nodes[5] = 5;
1045  assign_elem_nodes[6] = 1;
1046  assign_elem_nodes[7] = 4;
1047  break;
1048  }
1049 
1050  case QUAD9:
1051  {
1052  fe_descriptor_id = 300; // Plane Stress Quadratic Quadrilateral
1053  assign_elem_nodes[0] = 0;
1054  assign_elem_nodes[1] = 7;
1055  assign_elem_nodes[2] = 3;
1056  assign_elem_nodes[3] = 6;
1057  assign_elem_nodes[4] = 2;
1058  assign_elem_nodes[5] = 5;
1059  assign_elem_nodes[6] = 1;
1060  assign_elem_nodes[7] = 4;
1061  assign_elem_nodes[8] = 8;
1062  break;
1063  }
1064 
1065  case TET4:
1066  {
1067  fe_descriptor_id = 111; // Solid Linear Tetrahedron
1068  assign_elem_nodes[0] = 0;
1069  assign_elem_nodes[1] = 1;
1070  assign_elem_nodes[2] = 2;
1071  assign_elem_nodes[3] = 3;
1072  break;
1073  }
1074 
1075  case PRISM6:
1076  {
1077  fe_descriptor_id = 112; // Solid Linear Prism
1078  assign_elem_nodes[0] = 0;
1079  assign_elem_nodes[1] = 1;
1080  assign_elem_nodes[2] = 2;
1081  assign_elem_nodes[3] = 3;
1082  assign_elem_nodes[4] = 4;
1083  assign_elem_nodes[5] = 5;
1084  break;
1085  }
1086 
1087  case HEX8:
1088  {
1089  fe_descriptor_id = 115; // Solid Linear Brick
1090  assign_elem_nodes[0] = 0;
1091  assign_elem_nodes[1] = 4;
1092  assign_elem_nodes[2] = 5;
1093  assign_elem_nodes[3] = 1;
1094  assign_elem_nodes[4] = 3;
1095  assign_elem_nodes[5] = 7;
1096  assign_elem_nodes[6] = 6;
1097  assign_elem_nodes[7] = 2;
1098  break;
1099  }
1100 
1101  case HEX20:
1102  {
1103  fe_descriptor_id = 116; // Solid Quadratic Brick
1104  assign_elem_nodes[ 0] = 0;
1105  assign_elem_nodes[ 1] = 12;
1106  assign_elem_nodes[ 2] = 4;
1107  assign_elem_nodes[ 3] = 16;
1108  assign_elem_nodes[ 4] = 5;
1109  assign_elem_nodes[ 5] = 13;
1110  assign_elem_nodes[ 6] = 1;
1111  assign_elem_nodes[ 7] = 8;
1112 
1113  assign_elem_nodes[ 8] = 11;
1114  assign_elem_nodes[ 9] = 19;
1115  assign_elem_nodes[10] = 17;
1116  assign_elem_nodes[11] = 9;
1117 
1118  assign_elem_nodes[12] = 3;
1119  assign_elem_nodes[13] = 15;
1120  assign_elem_nodes[14] = 7;
1121  assign_elem_nodes[15] = 18;
1122  assign_elem_nodes[16] = 6;
1123  assign_elem_nodes[17] = 14;
1124  assign_elem_nodes[18] = 2;
1125  assign_elem_nodes[19] = 10;
1126 
1127 
1128  break;
1129  }
1130 
1131  case TET10:
1132  {
1133  fe_descriptor_id = 118; // Solid Quadratic Tetrahedron
1134  assign_elem_nodes[0] = 0;
1135  assign_elem_nodes[1] = 4;
1136  assign_elem_nodes[2] = 1;
1137  assign_elem_nodes[3] = 5;
1138  assign_elem_nodes[4] = 2;
1139  assign_elem_nodes[5] = 6;
1140  assign_elem_nodes[6] = 7;
1141  assign_elem_nodes[7] = 8;
1142  assign_elem_nodes[8] = 9;
1143  assign_elem_nodes[9] = 3;
1144  break;
1145  }
1146 
1147  default:
1148  libmesh_error_msg("ERROR: Element type = " << elem->type() << " not supported in " << "UNVIO!");
1149  }
1150 
1151  dof_id_type elem_id = elem->id();
1152 
1153  out_file << std::setw(10) << elem_id // element ID
1154  << std::setw(10) << fe_descriptor_id // type of element
1155  << std::setw(10) << phys_prop_tab_dummy // not supported
1156  << std::setw(10) << mat_prop_tab_dummy // not supported
1157  << std::setw(10) << color_dummy // not supported
1158  << std::setw(10) << elem->n_nodes() // No. of nodes per element
1159  << '\n';
1160 
1161  for (unsigned int j=0; j<elem->n_nodes(); j++)
1162  {
1163  // assign_elem_nodes[j]-th node: i.e., j loops over the
1164  // libMesh numbering, and assign_elem_nodes[j] over the
1165  // UNV numbering.
1166  const Node * node_in_unv_order = elem->node_ptr(assign_elem_nodes[j]);
1167 
1168  // new record after 8 id entries
1169  if (j==8 || j==16)
1170  out_file << '\n';
1171 
1172  // Write foreign label for this node
1173  dof_id_type node_id = node_in_unv_order->id();
1174 
1175  out_file << std::setw(10) << node_id;
1176  }
1177 
1178  out_file << '\n';
1179 
1180  n_elem_written++;
1181  }
1182 
1183  if (this->verbose())
1184  libMesh::out << " Finished writing " << n_elem_written << " elements" << std::endl;
1185 
1186  // Write end of dataset
1187  out_file << " -1\n";
1188 }
bool & verbose()
Definition: unv_io.C:87
static const std::string _elements_dataset_label
Definition: unv_io.h:198
IterBase * end
const MT & mesh() const
Definition: mesh_output.h:216
OStreamProxy out(std::cout)
uint8_t dof_id_type
Definition: id_types.h:64
const std::vector< Number > * libMesh::UNVIO::get_data ( Node node) const
Returns
A pointer the values associated with the node node, as read in by the read_dataset() method.

If no values exist for the node in question, a libmesh_nullptr is returned instead. It is up to the user to check the return value before using it.

Definition at line 1354 of file unv_io.C.

References _node_data, and libmesh_nullptr.

1355 {
1356  std::map<Node *, std::vector<Number> >::const_iterator
1357  it = _node_data.find(node);
1358 
1359  if (it == _node_data.end())
1360  return libmesh_nullptr;
1361  else
1362  return &(it->second);
1363 }
std::map< Node *, std::vector< Number > > _node_data
Definition: unv_io.h:214
const class libmesh_nullptr_t libmesh_nullptr
void libMesh::UNVIO::groups_in ( std::istream &  in_file)
private

Reads the "groups" section of the file. The format of the groups section is described here: http://www.sdrl.uc.edu/universal-file-formats-for-modal-analysis-testing-1/file-format-storehouse/unv_2467.htm

Definition at line 413 of file unv_io.C.

References _unv_elem_id_to_libmesh_elem_id, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::BoundaryInfo::add_side(), libMesh::Elem::build_side_ptr(), libMesh::Elem::dim(), libMesh::MeshBase::elem_ptr(), end, libMesh::err, libMesh::MeshBase::get_boundary_info(), libMesh::Elem::key(), max_elem_dimension_seen(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::Elem::n_sides(), side, libMesh::BoundaryInfo::sideset_name(), libMesh::Elem::subdomain_id(), and libMesh::MeshBase::subdomain_name().

Referenced by read_implementation().

414 {
415  // Grab reference to the Mesh, so we can add boundary info data to it
416  MeshBase & mesh = MeshInput<MeshBase>::mesh();
417 
418  // Record the max and min element dimension seen while reading the file.
419  unsigned char max_dim = this->max_elem_dimension_seen();
420 
421  // Container which stores lower-dimensional elements (based on
422  // Elem::key()) for later assignment of boundary conditions. We
423  // could use e.g. an unordered_set with Elems as keys for this, but
424  // this turns out to be much slower on the search side, since we
425  // have to build an entire side in order to search, rather than just
426  // calling elem->key(side) to compute a key.
427  typedef LIBMESH_BEST_UNORDERED_MULTIMAP<dof_id_type, Elem *> map_type;
428  map_type provide_bcs;
429 
430  // Read groups until there aren't any more to read...
431  while (true)
432  {
433  // If we read a -1, it means there is nothing else to read in this section.
434  int group_number;
435  in_file >> group_number;
436 
437  if (group_number == -1)
438  break;
439 
440  // The first record consists of 8 entries:
441  // Field 1 -- group number (that we just read)
442  // Field 2 -- active constraint set no. for group
443  // Field 3 -- active restraint set no. for group
444  // Field 4 -- active load set no. for group
445  // Field 5 -- active dof set no. for group
446  // Field 6 -- active temperature set no. for group
447  // Field 7 -- active contact set no. for group
448  // Field 8 -- number of entities in group
449  // Only the first and last of these are relevant to us...
450  unsigned num_entities;
451  std::string group_name;
452  {
453  unsigned dummy;
454  in_file >> dummy >> dummy >> dummy >> dummy >> dummy >> dummy
455  >> num_entities;
456 
457  // The second record has 1 field, the group name
458  in_file >> group_name;
459  }
460 
461  // The dimension of the elements in the group will determine
462  // whether this is a sideset group or a subdomain group.
463  bool
464  is_subdomain_group = false,
465  is_sideset_group = false;
466 
467  // Each subsequent line has 8 entries, there are two entity type
468  // codes and two tags per line that we need to read. In all my
469  // examples, the entity type code was always "8", which stands for
470  // "finite elements" but it's possible that we could eventually
471  // support other entity type codes...
472  // Field 1 -- entity type code
473  // Field 2 -- entity tag
474  // Field 3 -- entity node leaf id.
475  // Field 4 -- entity component/ ham id.
476  // Field 5 -- entity type code
477  // Field 6 -- entity tag
478  // Field 7 -- entity node leaf id.
479  // Field 8 -- entity component/ ham id.
480  {
481  unsigned entity_type_code, entity_tag, dummy;
482  for (unsigned entity=0; entity<num_entities; ++entity)
483  {
484  in_file >> entity_type_code >> entity_tag >> dummy >> dummy;
485 
486  if (entity_type_code != 8)
487  libMesh::err << "Warning, unrecognized entity type code = "
488  << entity_type_code
489  << " in group "
490  << group_name
491  << std::endl;
492 
493  // Try to find this ID in the map from UNV element ids to libmesh ids.
494  std::map<unsigned, unsigned>::iterator it =
495  _unv_elem_id_to_libmesh_elem_id.find(entity_tag);
496 
497  if (it != _unv_elem_id_to_libmesh_elem_id.end())
498  {
499  unsigned libmesh_elem_id = it->second;
500 
501  // Attempt to get a pointer to the elem listed in the group
502  Elem * group_elem = mesh.elem_ptr(libmesh_elem_id);
503 
504  // dim < max_dim means the Elem defines a boundary condition
505  if (group_elem->dim() < max_dim)
506  {
507  is_sideset_group = true;
508 
509  // We can only handle elements that are *exactly*
510  // one dimension lower than the max element
511  // dimension. Not sure if "edge" BCs in 3D
512  // actually make sense/are required...
513  if (group_elem->dim()+1 != max_dim)
514  libmesh_error_msg("ERROR: Expected boundary element of dimension " \
515  << max_dim-1 << " but got " << group_elem->dim());
516 
517  // Set the current group number as the lower-dimensional element's subdomain ID.
518  // We will use this later to set a boundary ID.
519  group_elem->subdomain_id() =
520  cast_int<subdomain_id_type>(group_number);
521 
522  // Store the lower-dimensional element in the provide_bcs container.
523  provide_bcs.insert(std::make_pair(group_elem->key(), group_elem));
524  }
525 
526  // dim == max_dim means this group defines a subdomain ID
527  else if (group_elem->dim() == max_dim)
528  {
529  is_subdomain_group = true;
530  group_elem->subdomain_id() =
531  cast_int<subdomain_id_type>(group_number);
532  }
533 
534  else
535  libmesh_error_msg("ERROR: Found an elem with dim=" \
536  << group_elem->dim() << " > " << "max_dim=" << +max_dim);
537  }
538  else
539  libMesh::err << "WARNING: UNV Element " << entity_tag << " not found while parsing groups." << std::endl;
540  } // end for (entity)
541  } // end scope
542 
543  // Associate this group_number with the group_name in the BoundaryInfo object.
544  if (is_sideset_group)
545  mesh.get_boundary_info().sideset_name
546  (cast_int<boundary_id_type>(group_number)) = group_name;
547 
548  if (is_subdomain_group)
549  mesh.subdomain_name
550  (cast_int<subdomain_id_type>(group_number)) = group_name;
551 
552  } // end while (true)
553 
554  // Loop over elements and try to assign boundary information
555  {
556  MeshBase::element_iterator it = mesh.active_elements_begin();
557  const MeshBase::element_iterator end = mesh.active_elements_end();
558  for ( ; it != end; ++it)
559  {
560  Elem * elem = *it;
561 
562  if (elem->dim() == max_dim)
563  {
564  // This is a max-dimension element that may require BCs.
565  // For each of its sides, including internal sides, we'll
566  // see if a lower-dimensional element provides boundary
567  // information for it. Note that we have not yet called
568  // find_neighbors(), so we can't use elem->neighbor(sn) in
569  // this algorithm...
570  for (unsigned short sn=0; sn<elem->n_sides(); sn++)
571  {
572  // Look for this key in the provide_bcs map
573  std::pair<map_type::const_iterator,
574  map_type::const_iterator>
575  range = provide_bcs.equal_range (elem->key(sn));
576 
577  // Add boundary information for each side in the range.
578  for (map_type::const_iterator iter = range.first;
579  iter != range.second; ++iter)
580  {
581  // Build a side to confirm the hash mapped to the correct side.
582  UniquePtr<Elem> side (elem->build_side_ptr(sn));
583 
584  // Get a pointer to the lower-dimensional element
585  Elem * lower_dim_elem = iter->second;
586 
587  // This was a hash, so it might not be perfect. Let's verify...
588  if (*lower_dim_elem == *side)
589  mesh.get_boundary_info().add_side(elem, sn, lower_dim_elem->subdomain_id());
590  }
591  }
592  }
593  }
594  }
595 }
unsigned short int side
Definition: xdr_io.C:49
IterBase * end
unsigned char max_elem_dimension_seen()
Definition: unv_io.C:380
OStreamProxy err(std::cerr)
std::map< unsigned, unsigned > _unv_elem_id_to_libmesh_elem_id
Definition: unv_io.h:208
unsigned char libMesh::UNVIO::max_elem_dimension_seen ( )
private
Returns
The maximum geometric element dimension encountered while reading the Mesh. Only valid after the elements have been read in and the elems_of_dimension array has been populated.

Definition at line 380 of file unv_io.C.

References libMesh::MeshInput< MeshBase >::elems_of_dimension.

Referenced by groups_in(), and read_implementation().

381 {
382  unsigned char max_dim = 0;
383 
384  unsigned char elem_dimensions_size = cast_int<unsigned char>
385  (elems_of_dimension.size());
386  // The elems_of_dimension array is 1-based in the UNV reader
387  for (unsigned char i=1; i<elem_dimensions_size; ++i)
388  if (elems_of_dimension[i])
389  max_dim = i;
390 
391  return max_dim;
392 }
std::vector< bool > elems_of_dimension
Definition: mesh_input.h:97
MeshBase & libMesh::MeshInput< MeshBase >::mesh ( )
protectedinherited
Returns
The object as a writable reference.

Referenced by libMesh::GMVIO::_read_one_cell(), libMesh::TetGenIO::element_in(), elements_in(), elements_out(), groups_in(), libMesh::TetGenIO::node_in(), nodes_in(), nodes_out(), libMesh::GMVIO::read(), libMesh::ExodusII_IO::read(), libMesh::CheckpointIO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_bcs(), libMesh::CheckpointIO::read_connectivity(), libMesh::CheckpointIO::read_header(), libMesh::UCDIO::read_implementation(), read_implementation(), libMesh::GmshIO::read_mesh(), libMesh::CheckpointIO::read_nodes(), libMesh::CheckpointIO::read_nodesets(), libMesh::CheckpointIO::read_remote_elem(), libMesh::XdrIO::read_serialized_bcs_helper(), 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::UCDIO::UCDIO(), libMesh::VTKIO::VTKIO(), libMesh::TetGenIO::write(), libMesh::ExodusII_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::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::UCDIO::write_implementation(), libMesh::GmshIO::write_mesh(), libMesh::UCDIO::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_nodesets(), libMesh::XdrIO::write_parallel(), libMesh::GmshIO::write_post(), libMesh::XdrIO::write_serialized_bcs_helper(), 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().

bool libMesh::UNVIO::need_D_to_e ( std::string &  number)
private

Replaces "1.1111D+00" with "1.1111e+00" if necessary. This function only needs to be called once per stream, one can assume that if one number needs rewriting, they all do.

Returns
true if the replacement occurs, false otherwise.

Definition at line 396 of file unv_io.C.

Referenced by read_dataset().

397 {
398  // find "D" in string, start looking at 6th element since dont expect a "D" earlier.
399  std::string::size_type position = number.find("D", 6);
400 
401  if (position != std::string::npos)
402  {
403  // replace "D" in string
404  number.replace(position, 1, "e");
405  return true;
406  }
407  else
408  return false;
409 }
void libMesh::UNVIO::nodes_in ( std::istream &  in_file)
private

Read nodes from file.

Definition at line 307 of file unv_io.C.

References _unv_node_id_to_libmesh_node_ptr, libMesh::MeshBase::add_point(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::out, and verbose().

Referenced by read_implementation().

308 {
309  LOG_SCOPE("nodes_in()","UNVIO");
310 
311  if (this->verbose())
312  libMesh::out << " Reading nodes" << std::endl;
313 
314  MeshBase & mesh = MeshInput<MeshBase>::mesh();
315 
316  // node label, we use an int here so we can read in a -1
317  int node_label;
318 
319  // always 3 coordinates in the UNV file, no matter
320  // what LIBMESH_DIM is.
321  Point xyz;
322 
323  // We'll just read the floating point values as strings even when
324  // there are no "D" characters in the file. This will make parsing
325  // the file a bit slower but the logic to do so will all be
326  // simplified and in one place...
327  std::string line;
328  std::istringstream coords_stream;
329 
330  // Continue reading nodes until there are none left
331  unsigned ctr = 0;
332  while (true)
333  {
334  // Read the node label
335  in_file >> node_label;
336 
337  // Break out of the while loop when we hit -1
338  if (node_label == -1)
339  break;
340 
341  // Read and discard the the rest of the node data on this line
342  // which we do not currently use:
343  // .) exp_coord_sys_num
344  // .) disp_coord_sys_num
345  // .) color
346  std::getline(in_file, line);
347 
348  // read the floating-point data
349  std::getline(in_file, line);
350 
351  // Replace any "D" characters used for exponents with "E"
352  size_t last_pos = 0;
353  while (true)
354  {
355  last_pos = line.find("D", last_pos);
356 
357  if (last_pos != std::string::npos)
358  line.replace(last_pos, 1, "E");
359  else
360  break;
361  }
362 
363  // Construct a stream object from the current line and then
364  // stream in the xyz values.
365  coords_stream.str(line);
366  coords_stream.clear(); // clear iostate bits! Important!
367  coords_stream >> xyz(0) >> xyz(1) >> xyz(2);
368 
369  // Add node to the Mesh, bump the counter.
370  Node * added_node = mesh.add_point(xyz, ctr++);
371 
372  // Maintain the mapping between UNV node ids and libmesh Node
373  // pointers.
374  _unv_node_id_to_libmesh_node_ptr[node_label] = added_node;
375  }
376 }
bool & verbose()
Definition: unv_io.C:87
std::map< dof_id_type, Node * > _unv_node_id_to_libmesh_node_ptr
Definition: unv_io.h:188
OStreamProxy out(std::cout)
void libMesh::UNVIO::nodes_out ( std::ostream &  out_file)
private

Outputs nodes to the file out_file. Do not use this directly, but through the proper write method.

Definition at line 889 of file unv_io.C.

References _nodes_dataset_label, end, libMesh::DofObject::id(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), libMesh::out, libMesh::Real, verbose(), and libMesh::x.

Referenced by write_implementation().

890 {
891  if (this->verbose())
892  libMesh::out << " Writing " << MeshOutput<MeshBase>::mesh().n_nodes() << " nodes" << std::endl;
893 
894  // Write beginning of dataset
895  out_file << " -1\n"
896  << " "
898  << '\n';
899 
900  unsigned int
901  exp_coord_sys_dummy = 0, // export coordinate sys. (not supported yet)
902  disp_coord_sys_dummy = 0, // displacement coordinate sys. (not supp. yet)
903  color_dummy = 0; // color(not supported yet)
904 
905  // A reference to the parent class's mesh
906  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
907 
908  // Use scientific notation with capital E and 16 digits for printing out the coordinates
909  out_file << std::scientific << std::setprecision(16) << std::uppercase;
910 
911  MeshBase::const_node_iterator nd = mesh.nodes_begin();
912  const MeshBase::const_node_iterator end = mesh.nodes_end();
913 
914  for (; nd != end; ++nd)
915  {
916  const Node * current_node = *nd;
917 
918  dof_id_type node_id = current_node->id();
919 
920  out_file << std::setw(10) << node_id
921  << std::setw(10) << exp_coord_sys_dummy
922  << std::setw(10) << disp_coord_sys_dummy
923  << std::setw(10) << color_dummy
924  << '\n';
925 
926  // The coordinates - always write out three coords regardless of LIBMESH_DIM
927  Real x = (*current_node)(0);
928 
929 #if LIBMESH_DIM > 1
930  Real y = (*current_node)(1);
931 #else
932  Real y = 0.;
933 #endif
934 
935 #if LIBMESH_DIM > 2
936  Real z = (*current_node)(2);
937 #else
938  Real z = 0.;
939 #endif
940 
941  out_file << std::setw(25) << x
942  << std::setw(25) << y
943  << std::setw(25) << z
944  << '\n';
945  }
946 
947  // Write end of dataset
948  out_file << " -1\n";
949 }
bool & verbose()
Definition: unv_io.C:87
IterBase * end
const MT & mesh() const
Definition: mesh_output.h:216
PetscErrorCode Vec x
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string _nodes_dataset_label
Definition: unv_io.h:193
OStreamProxy out(std::cout)
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::UNVIO::read ( const std::string &  file_name)
virtual

This method implements reading a mesh from a specified file.

Implements libMesh::MeshInput< MeshBase >.

Definition at line 94 of file unv_io.C.

References read_implementation().

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

95 {
96  if (file_name.rfind(".gz") < file_name.size())
97  {
98 #ifdef LIBMESH_HAVE_GZSTREAM
99 
100  igzstream in_stream (file_name.c_str());
101  this->read_implementation (in_stream);
102 
103 #else
104 
105  libmesh_error_msg("ERROR: You must have the zlib.h header files and libraries to read and write compressed streams.");
106 
107 #endif
108  return;
109  }
110 
111  else
112  {
113  std::ifstream in_stream (file_name.c_str());
114  this->read_implementation (in_stream);
115  return;
116  }
117 }
void read_implementation(std::istream &in_stream)
Definition: unv_io.C:120
void libMesh::UNVIO::read_dataset ( std::string  file_name)

Read a UNV data file containing a dataset of type "2414". For more info, see http://tinyurl.com/htcf6zm

Definition at line 1192 of file unv_io.C.

References _node_data, _unv_node_id_to_libmesh_node_ptr, need_D_to_e(), and libMesh::Real.

1193 {
1194  std::ifstream in_stream(file_name.c_str());
1195 
1196  if (!in_stream.good())
1197  libmesh_error_msg("Error opening UNV data file.");
1198 
1199  std::string olds, news, dummy;
1200 
1201  while (true)
1202  {
1203  in_stream >> olds >> news;
1204 
1205  // A "-1" followed by a number means the beginning of a dataset.
1206  while (((olds != "-1") || (news == "-1")) && !in_stream.eof())
1207  {
1208  olds = news;
1209  in_stream >> news;
1210  }
1211 
1212  if (in_stream.eof())
1213  break;
1214 
1215  // We only support reading datasets in the "2414" format.
1216  if (news == "2414")
1217  {
1218  // Ignore the rest of the current line and the next two
1219  // lines that contain analysis dataset label and name.
1220  for (unsigned int i=0; i<3; i++)
1221  std::getline(in_stream, dummy);
1222 
1223  // Read the dataset location, where
1224  // 1: Data at nodes
1225  // 2: Data on elements
1226  // Currently only data on nodes is supported.
1227  unsigned int dataset_location;
1228  in_stream >> dataset_location;
1229 
1230  // Currently only nodal datasets are supported.
1231  if (dataset_location != 1)
1232  libmesh_error_msg("ERROR: Currently only Data at nodes is supported.");
1233 
1234  // Ignore the rest of this line and the next five records.
1235  for (unsigned int i=0; i<6; i++)
1236  std::getline(in_stream, dummy);
1237 
1238  // These data are all of no interest to us...
1239  unsigned int
1240  model_type,
1241  analysis_type,
1242  data_characteristic,
1243  result_type;
1244 
1245  // The type of data (complex, real, float, double etc, see
1246  // below).
1247  unsigned int data_type;
1248 
1249  // The number of floating-point values per entity.
1250  unsigned int num_vals_per_node;
1251 
1252  in_stream >> model_type // not used here
1253  >> analysis_type // not used here
1254  >> data_characteristic // not used here
1255  >> result_type // not used here
1256  >> data_type
1257  >> num_vals_per_node;
1258 
1259  // Ignore the rest of record 9, and records 10-13.
1260  for (unsigned int i=0; i<5; i++)
1261  std::getline(in_stream, dummy);
1262 
1263  // Now get the foreign (aka UNV node) node number and
1264  // the respective nodal data.
1265  int f_n_id;
1266  std::vector<Number> values;
1267 
1268  while (true)
1269  {
1270  in_stream >> f_n_id;
1271 
1272  // If -1 then we have reached the end of the dataset.
1273  if (f_n_id == -1)
1274  break;
1275 
1276  // Resize the values vector (usually data in three
1277  // principle directions, i.e. num_vals_per_node = 3).
1278  values.resize(num_vals_per_node);
1279 
1280  // Read the values for the respective node.
1281  for (unsigned int data_cnt=0; data_cnt<num_vals_per_node; data_cnt++)
1282  {
1283  // Check what data type we are reading.
1284  // 2,4: Real
1285  // 5,6: Complex
1286  // other data types are not supported yet.
1287  // As again, these floats may also be written
1288  // using a 'D' instead of an 'e'.
1289  if (data_type == 2 || data_type == 4)
1290  {
1291  std::string buf;
1292  in_stream >> buf;
1293  this->need_D_to_e(buf);
1294 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1295  values[data_cnt] = Complex(std::atof(buf.c_str()), 0.);
1296 #else
1297  values[data_cnt] = std::atof(buf.c_str());
1298 #endif
1299  }
1300 
1301  else if (data_type == 5 || data_type == 6)
1302  {
1303 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1304  Real re_val, im_val;
1305 
1306  std::string buf;
1307  in_stream >> buf;
1308 
1309  if (this->need_D_to_e(buf))
1310  {
1311  re_val = std::atof(buf.c_str());
1312  in_stream >> buf;
1313  this->need_D_to_e(buf);
1314  im_val = std::atof(buf.c_str());
1315  }
1316  else
1317  {
1318  re_val = std::atof(buf.c_str());
1319  in_stream >> im_val;
1320  }
1321 
1322  values[data_cnt] = Complex(re_val,im_val);
1323 #else
1324 
1325  libmesh_error_msg("ERROR: Complex data only supported when libMesh is configured with --enable-complex!");
1326 #endif
1327  }
1328 
1329  else
1330  libmesh_error_msg("ERROR: Data type not supported.");
1331 
1332  } // end loop data_cnt
1333 
1334  // Get a pointer to the Node associated with the UNV node id.
1335  std::map<dof_id_type, Node *>::const_iterator it =
1336  _unv_node_id_to_libmesh_node_ptr.find(f_n_id);
1337 
1338  if (it == _unv_node_id_to_libmesh_node_ptr.end())
1339  libmesh_error_msg("UNV node id " << f_n_id << " was not found.");
1340 
1341  // Store the nodal values in our map which uses the
1342  // libMesh Node* as the key. We use operator[] here
1343  // because we want to create an empty vector if the
1344  // entry does not already exist.
1345  _node_data[it->second] = values;
1346  } // end while (true)
1347  } // end if (news == "2414")
1348  } // end while (true)
1349 }
std::map< dof_id_type, Node * > _unv_node_id_to_libmesh_node_ptr
Definition: unv_io.h:188
std::map< Node *, std::vector< Number > > _node_data
Definition: unv_io.h:214
MPI_Datatype data_type
Definition: parallel.h:162
std::complex< Real > Complex
bool need_D_to_e(std::string &number)
Definition: unv_io.C:396
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void libMesh::UNVIO::read_implementation ( std::istream &  in_stream)
private

The actual implementation of the read function. The public read interface simply decides which type of stream to pass the implementation.

Definition at line 120 of file unv_io.C.

References _elements_dataset_label, _groups_dataset_label, _nodes_dataset_label, libMesh::MeshBase::delete_elem(), libMesh::Elem::dim(), libMesh::MeshBase::elements_begin(), libMesh::MeshBase::elements_end(), elements_in(), libMesh::MeshInput< MeshBase >::elems_of_dimension, groups_in(), max_elem_dimension_seen(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshInput< MeshBase >::mesh(), nodes_in(), libMesh::out, and verbose().

Referenced by read().

121 {
122  // Keep track of what kinds of elements this file contains
123  elems_of_dimension.clear();
124  elems_of_dimension.resize(4, false);
125 
126  {
127  if (!in_stream.good())
128  libmesh_error_msg("ERROR: Input file not good.");
129 
130  // Flags to be set when certain sections are encountered
131  bool
132  found_node = false,
133  found_elem = false,
134  found_group = false;
135 
136  // strings for reading the file line by line
137  std::string
138  old_line,
139  current_line;
140 
141  while (true)
142  {
143  // Save off the old_line. This will provide extra reliability
144  // for detecting the beginnings of the different sections.
145  old_line = current_line;
146 
147  // Try to read something. This may set EOF!
148  std::getline(in_stream, current_line);
149 
150  // If the stream is still "valid", parse the line
151  if (in_stream)
152  {
153  // UNV files always have some amount of leading
154  // whitespace, let's not rely on exactly how much... This
155  // command deletes it.
156  current_line.erase(std::remove_if(current_line.begin(), current_line.end(), isspace),
157  current_line.end());
158 
159  // Parse the nodes section
160  if (current_line == _nodes_dataset_label &&
161  old_line == "-1")
162  {
163  found_node = true;
164  this->nodes_in(in_stream);
165  }
166 
167  // Parse the elements section
168  else if (current_line == _elements_dataset_label &&
169  old_line == "-1")
170  {
171  // The current implementation requires the nodes to
172  // have been read before reaching the elements
173  // section.
174  if (!found_node)
175  libmesh_error_msg("ERROR: The Nodes section must come before the Elements section of the UNV file!");
176 
177  found_elem = true;
178  this->elements_in(in_stream);
179  }
180 
181  // Parse the groups section
182  else if (current_line == _groups_dataset_label &&
183  old_line == "-1")
184  {
185  // The current implementation requires the nodes and
186  // elements to have already been read before reaching
187  // the groups section.
188  if (!found_node || !found_elem)
189  libmesh_error_msg("ERROR: The Nodes and Elements sections must come before the Groups section of the UNV file!");
190 
191  found_group = true;
192  this->groups_in(in_stream);
193  }
194 
195  // We can stop reading once we've found the nodes, elements,
196  // and group sections.
197  if (found_node && found_elem && found_group)
198  break;
199 
200  // If we made it here, we're not done yet, so keep reading
201  continue;
202  }
203 
204  // if (!in_stream) check to see if EOF was set. If so, break out of while loop.
205  if (in_stream.eof())
206  break;
207 
208  // If !in_stream and !in_stream.eof(), stream is in a bad state!
209  libmesh_error_msg("Stream is bad! Perhaps the file does not exist?");
210  } // end while (true)
211 
212  // By now we better have found the datasets for nodes and elements,
213  // otherwise the unv file is bad!
214  if (!found_node)
215  libmesh_error_msg("ERROR: Could not find nodes!");
216 
217  if (!found_elem)
218  libmesh_error_msg("ERROR: Could not find elements!");
219  }
220 
221 
222 
223  {
224  // Set the mesh dimension to the largest element dimension encountered
225  MeshInput<MeshBase>::mesh().set_mesh_dimension(max_elem_dimension_seen());
226 
227 #if LIBMESH_DIM < 3
228  if (MeshInput<MeshBase>::mesh().mesh_dimension() > LIBMESH_DIM)
229  libmesh_error_msg("Cannot open dimension " \
230  << MeshInput<MeshBase>::mesh().mesh_dimension() \
231  << " mesh file when configured without " \
232  << MeshInput<MeshBase>::mesh().mesh_dimension() \
233  << "D support." );
234 #endif
235 
236  // Delete any lower-dimensional elements that might have been
237  // added to the mesh strictly for setting BCs.
238  {
239  // Grab reference to the Mesh
240  MeshBase & mesh = MeshInput<MeshBase>::mesh();
241 
242  unsigned char max_dim = this->max_elem_dimension_seen();
243 
244  MeshBase::const_element_iterator el = mesh.elements_begin();
245  const MeshBase::const_element_iterator end_el = mesh.elements_end();
246 
247  for (; el != end_el; ++el)
248  {
249  Elem * elem = *el;
250 
251  if (elem->dim() < max_dim)
252  mesh.delete_elem(elem);
253  }
254  }
255 
256  if (this->verbose())
257  libMesh::out << " Finished." << std::endl << std::endl;
258  }
259 }
bool & verbose()
Definition: unv_io.C:87
static const std::string _elements_dataset_label
Definition: unv_io.h:198
std::vector< bool > elems_of_dimension
Definition: mesh_input.h:97
void nodes_in(std::istream &in_file)
Definition: unv_io.C:307
unsigned char max_elem_dimension_seen()
Definition: unv_io.C:380
static const std::string _groups_dataset_label
Definition: unv_io.h:203
static const std::string _nodes_dataset_label
Definition: unv_io.h:193
OStreamProxy out(std::cout)
void elements_in(std::istream &in_file)
Definition: unv_io.C:599
void groups_in(std::istream &in_file)
Definition: unv_io.C:413
void libMesh::MeshInput< MeshBase >::set_n_partitions ( unsigned int  n_parts)
inlineprotectedinherited

Sets the number of partitions in the mesh. Typically this gets done by the partitioner, but some parallel file formats begin "pre-partitioned".

Definition at line 91 of file mesh_input.h.

References libMesh::MeshInput< MT >::mesh().

Referenced by libMesh::Nemesis_IO::read(), and libMesh::XdrIO::read().

91 { this->mesh().set_n_partitions() = n_parts; }
unsigned int & set_n_partitions()
Definition: mesh_base.h:1302
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().

bool & libMesh::UNVIO::verbose ( )

Set the flag indicating if we should be verbose.

Definition at line 87 of file unv_io.C.

References _verbose.

Referenced by elements_in(), elements_out(), nodes_in(), nodes_out(), and read_implementation().

88 {
89  return _verbose;
90 }
bool _verbose
Definition: unv_io.h:182
void libMesh::UNVIO::write ( const std::string &  file_name)
virtual

This method implements writing a mesh to a specified file.

Implements libMesh::MeshOutput< MeshBase >.

Definition at line 265 of file unv_io.C.

References write_implementation().

Referenced by libMesh::NameBasedIO::write().

266 {
267  if (file_name.rfind(".gz") < file_name.size())
268  {
269 #ifdef LIBMESH_HAVE_GZSTREAM
270 
271  ogzstream out_stream(file_name.c_str());
272  this->write_implementation (out_stream);
273 
274 #else
275 
276  libmesh_error_msg("ERROR: You must have the zlib.h header files and libraries to read and write compressed streams.");
277 
278 #endif
279 
280  return;
281  }
282 
283  else
284  {
285  std::ofstream out_stream (file_name.c_str());
286  this->write_implementation (out_stream);
287  return;
288  }
289 }
void write_implementation(std::ostream &out_stream)
Definition: unv_io.C:294
virtual void libMesh::MeshOutput< MeshBase >::write_equation_systems ( const std::string &  ,
const EquationSystems ,
const std::set< std::string > *  system_names = libmesh_nullptr 
)
virtualinherited

This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object.

Reimplemented in libMesh::NameBasedIO.

Referenced by libMesh::Nemesis_IO::write_timestep(), and libMesh::ExodusII_IO::write_timestep().

void libMesh::UNVIO::write_implementation ( std::ostream &  out_stream)
private

The actual implementation of the write function. The public write interface simply decides which type of stream to pass the implementation.

Definition at line 294 of file unv_io.C.

References elements_out(), and nodes_out().

Referenced by write().

295 {
296  if (!out_file.good())
297  libmesh_error_msg("ERROR: Output file not good.");
298 
299  // write the nodes, then the elements
300  this->nodes_out (out_file);
301  this->elements_out (out_file);
302 }
void nodes_out(std::ostream &out_file)
Definition: unv_io.C:889
void elements_out(std::ostream &out_file)
Definition: unv_io.C:956
virtual void libMesh::MeshOutput< MeshBase >::write_nodal_data ( const std::string &  ,
const std::vector< Number > &  ,
const std::vector< std::string > &   
)
inlinevirtualinherited

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided.

Reimplemented in libMesh::ExodusII_IO, libMesh::NameBasedIO, libMesh::GmshIO, libMesh::Nemesis_IO, libMesh::VTKIO, libMesh::UCDIO, libMesh::GMVIO, libMesh::MEDITIO, libMesh::GnuPlotIO, and libMesh::TecplotIO.

Definition at line 96 of file mesh_output.h.

References libMesh::MeshOutput< MT >::ascii_precision(), libMesh::MeshOutput< MT >::mesh(), and libMesh::MeshOutput< MT >::write_nodal_data().

99  { libmesh_not_implemented(); }
virtual void libMesh::MeshOutput< MeshBase >::write_nodal_data ( const std::string &  ,
const NumericVector< Number > &  ,
const std::vector< std::string > &   
)
virtualinherited

This method should be overridden by "parallel" output formats for writing nodal data. Instead of getting a localized copy of the nodal solution vector, it is passed a NumericVector of type=PARALLEL which is in node-major order i.e. (u0,v0,w0, u1,v1,w1, u2,v2,w2, u3,v3,w3, ...) and contains n_nodes*n_vars total entries. Then, it is up to the individual I/O class to extract the required solution values from this vector and write them in parallel.

If not implemented, localizes the parallel vector into a std::vector and calls the other version of this function.

Reimplemented in libMesh::Nemesis_IO.

Member Data Documentation

const std::string libMesh::UNVIO::_elements_dataset_label = "2412"
staticprivate

label for the element dataset

Definition at line 198 of file unv_io.h.

Referenced by elements_out(), and read_implementation().

const std::string libMesh::UNVIO::_groups_dataset_label = "2467"
staticprivate

label for the groups dataset

Definition at line 203 of file unv_io.h.

Referenced by read_implementation().

const bool libMesh::MeshOutput< MeshBase >::_is_parallel_format
protectedinherited

Flag specifying whether this format is parallel-capable. If this is false (default) I/O is only permitted when the mesh has been serialized.

Definition at line 141 of file mesh_output.h.

Referenced by libMesh::FroIO::write(), libMesh::PostscriptIO::write(), and libMesh::EnsightIO::write().

std::map<Node *, std::vector<Number> > libMesh::UNVIO::_node_data
private

Map from libMesh Node* to data at that node, as read in by the read_dataset() function.

Definition at line 214 of file unv_io.h.

Referenced by get_data(), and read_dataset().

const std::string libMesh::UNVIO::_nodes_dataset_label = "2411"
staticprivate

label for the node dataset

Definition at line 193 of file unv_io.h.

Referenced by nodes_out(), and read_implementation().

const bool libMesh::MeshOutput< MeshBase >::_serial_only_needed_on_proc_0
protectedinherited

Flag specifying whether this format can be written by only serializing the mesh to processor zero

If this is false (default) the mesh will be serialized to all processors

Definition at line 150 of file mesh_output.h.

std::map<unsigned, unsigned> libMesh::UNVIO::_unv_elem_id_to_libmesh_elem_id
private

Map UNV element IDs to libmesh element IDs.

Definition at line 208 of file unv_io.h.

Referenced by elements_in(), and groups_in().

std::map<dof_id_type, Node *> libMesh::UNVIO::_unv_node_id_to_libmesh_node_ptr
private

Maps UNV node IDs to libMesh Node*s. Used when reading. Even if the libMesh Mesh is renumbered, this map should continue to be valid.

Definition at line 188 of file unv_io.h.

Referenced by elements_in(), nodes_in(), and read_dataset().

bool libMesh::UNVIO::_verbose
private

should be be verbose?

Definition at line 182 of file unv_io.h.

Referenced by verbose().


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