libMesh::GmshIO Class Reference

#include <gmsh_io.h>

Inheritance diagram for libMesh::GmshIO:

Classes

struct  ElementDefinition
 
struct  ElementMaps
 

Public Member Functions

 GmshIO (MeshBase &mesh)
 
 GmshIO (const MeshBase &mesh)
 
virtual void read (const std::string &name) libmesh_override
 
virtual void write (const std::string &name) libmesh_override
 
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) libmesh_override
 
bool & binary ()
 
bool & write_lower_dimensional_elements ()
 
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 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_mesh (std::istream &in)
 
void write_mesh (std::ostream &out)
 
void write_post (const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
 

Static Private Member Functions

static ElementMaps build_element_maps ()
 

Private Attributes

bool _binary
 
bool _write_lower_dimensional_elements
 

Static Private Attributes

static ElementMaps _element_maps = GmshIO::build_element_maps()
 

Detailed Description

This class implements writing meshes in the Gmsh format. For a full description of the Gmsh format and to obtain the GMSH software see the Gmsh home page

Author
John W. Peterson
Date
2004, 2014
Author
Martin Luthi
Date
2005

Definition at line 50 of file gmsh_io.h.

Constructor & Destructor Documentation

libMesh::GmshIO::GmshIO ( MeshBase mesh)
explicit

Constructor. Takes a non-const Mesh reference which it will fill up with elements via the read() command.

Definition at line 120 of file gmsh_io.C.

120  :
121  MeshInput<MeshBase> (mesh),
122  MeshOutput<MeshBase> (mesh),
123  _binary (false),
125 {
126 }
bool _write_lower_dimensional_elements
Definition: gmsh_io.h:148
libMesh::GmshIO::GmshIO ( const MeshBase mesh)
explicit

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

Definition at line 111 of file gmsh_io.C.

111  :
112  MeshOutput<MeshBase>(mesh),
113  _binary(false),
115 {
116 }
bool _write_lower_dimensional_elements
Definition: gmsh_io.h:148

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().

bool & libMesh::GmshIO::binary ( )

Flag indicating whether or not to write a binary file. While binary files may end up being smaller than equivalent ASCII files, they will almost certainly take longer to write. The reason for this is that the ostream::write() function which is used to write "binary" data to streams, only takes a pointer to char as its first argument. This means if you want to write anything other than a buffer of chars, you first have to use a strange memcpy hack to get the data into the desired format. See the templated to_binary_stream() function below.

Definition at line 130 of file gmsh_io.C.

References _binary.

Referenced by write_post().

131 {
132  return _binary;
133 }
GmshIO::ElementMaps libMesh::GmshIO::build_element_maps ( )
staticprivate

A static function used to construct the _element_maps struct, statically.

Definition at line 42 of file gmsh_io.C.

References libMesh::GmshIO::ElementMaps::add_def(), libMesh::EDGE2, libMesh::EDGE3, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::GmshIO::ElementMaps::in, libMesh::GmshIO::ElementDefinition::nnodes, libMesh::NODEELEM, libMesh::GmshIO::ElementDefinition::nodes, libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::swap(), libMesh::TET10, libMesh::TET4, libMesh::TRI3, and libMesh::TRI6.

43 {
44  // Object to be filled up
45  ElementMaps em;
46 
47  // POINT (import only)
48  em.in.insert(std::make_pair(15, ElementDefinition(NODEELEM, 15, 0, 1)));
49 
50  // Add elements with trivial node mappings
51  em.add_def(ElementDefinition(EDGE2, 1, 1, 2));
52  em.add_def(ElementDefinition(EDGE3, 8, 1, 3));
53  em.add_def(ElementDefinition(TRI3, 2, 2, 3));
54  em.add_def(ElementDefinition(TRI6, 9, 2, 6));
55  em.add_def(ElementDefinition(QUAD4, 3, 2, 4));
56  em.add_def(ElementDefinition(QUAD8, 16, 2, 8));
57  em.add_def(ElementDefinition(QUAD9, 10, 2, 9));
58  em.add_def(ElementDefinition(HEX8, 5, 3, 8));
59  em.add_def(ElementDefinition(TET4, 4, 3, 4));
60  em.add_def(ElementDefinition(PRISM6, 6, 3, 6));
61  em.add_def(ElementDefinition(PYRAMID5, 7, 3, 5));
62 
63  // Add elements with non-trivial node mappings
64 
65  // HEX20
66  {
67  ElementDefinition eledef(HEX20, 17, 3, 20);
68  const unsigned int nodes[] = {0,1,2,3,4,5,6,7,8,11,12,9,13,10,14,15,16,19,17,18};
69  std::vector<unsigned int>(nodes, nodes+eledef.nnodes).swap(eledef.nodes); // swap trick
70  em.add_def(eledef);
71  }
72 
73  // HEX27
74  {
75  ElementDefinition eledef(HEX27, 12, 3, 27);
76  const unsigned int nodes[] = {0,1,2,3,4,5,6,7,8,11,12,9,13,10,14,
77  15,16,19,17,18,20,21,24,22,23,25,26};
78  std::vector<unsigned int>(nodes, nodes+eledef.nnodes).swap(eledef.nodes); // swap trick
79  em.add_def(eledef);
80  }
81 
82  // TET10
83  {
84  ElementDefinition eledef(TET10, 11, 3, 10);
85  const unsigned int nodes[] = {0,1,2,3,4,5,6,7,9,8};
86  std::vector<unsigned int>(nodes, nodes+eledef.nnodes).swap(eledef.nodes); // swap trick
87  em.add_def(eledef);
88  }
89 
90  // PRISM15
91  {
92  ElementDefinition eledef(PRISM15, 18, 3, 15);
93  const unsigned int nodes[] = {0,1,2,3,4,5,6,8,9,7,10,11,12,14,13};
94  std::vector<unsigned int>(nodes, nodes+eledef.nnodes).swap(eledef.nodes); // swap trick
95  em.add_def(eledef);
96  }
97 
98  // PRISM18
99  {
100  ElementDefinition eledef(PRISM18, 13, 3, 18);
101  const unsigned int nodes[] = {0,1,2,3,4,5,6,8,9,7,10,11,12,14,13,15,17,16};
102  std::vector<unsigned int>(nodes, nodes+eledef.nnodes).swap(eledef.nodes); // swap trick
103  em.add_def(eledef);
104  }
105 
106  return em;
107 }
X_input swap(X_system)
MeshBase & libMesh::MeshInput< MeshBase >::mesh ( )
protectedinherited
Returns
The object as a writable reference.

Referenced by libMesh::GMVIO::_read_one_cell(), libMesh::TetGenIO::element_in(), libMesh::UNVIO::elements_in(), libMesh::UNVIO::elements_out(), libMesh::UNVIO::groups_in(), libMesh::TetGenIO::node_in(), libMesh::UNVIO::nodes_in(), libMesh::UNVIO::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(), libMesh::UNVIO::read_implementation(), 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(), 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(), 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().

void libMesh::GmshIO::read ( const std::string &  name)
virtual

Reads in a mesh in the Gmsh *.msh format from the ASCII file given by name.

The user is responsible for calling Mesh::prepare_for_use() after reading the mesh and before using it.

Implements libMesh::MeshInput< MeshBase >.

Definition at line 144 of file gmsh_io.C.

References read_mesh().

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

145 {
146  std::ifstream in (name.c_str());
147  this->read_mesh (in);
148 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
void read_mesh(std::istream &in)
Definition: gmsh_io.C:152
void libMesh::GmshIO::read_mesh ( std::istream &  in)
private

Implementation of the read() function. This function is called by the public interface function and implements reading the file.

Definition at line 152 of file gmsh_io.C.

References _element_maps, libMesh::Elem::build(), libMesh::Elem::build_side_ptr(), libMesh::GmshIO::ElementDefinition::dim, libMesh::Elem::dim(), end, libMesh::err, libMesh::GmshIO::ElementMaps::in, libMesh::Elem::key(), libMesh::libmesh_assert(), std::max(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshInput< MeshBase >::mesh(), std::min(), libMesh::Elem::n_nodes(), libMesh::Elem::n_sides(), libMesh::GmshIO::ElementDefinition::nnodes, libMesh::Elem::node_id(), libMesh::GmshIO::ElementDefinition::nodes, libMesh::processor_id(), libMesh::Real, libMesh::DofObject::set_id(), libMesh::Elem::set_node(), side, libMesh::Elem::subdomain_id(), libMesh::GmshIO::ElementDefinition::type, and libMesh::x.

Referenced by read().

153 {
154  // This is a serial-only process for now;
155  // the Mesh should be read on processor 0 and
156  // broadcast later
157  libmesh_assert_equal_to (MeshOutput<MeshBase>::mesh().processor_id(), 0);
158 
159  libmesh_assert(in.good());
160 
161  // clear any data in the mesh
162  MeshBase & mesh = MeshInput<MeshBase>::mesh();
163  mesh.clear();
164 
165  // some variables
166  int format=0, size=0;
167  Real version = 1.0;
168 
169  // Keep track of lower-dimensional blocks which are not BCs, but
170  // actually blocks of lower-dimensional elements.
171  std::set<subdomain_id_type> lower_dimensional_blocks;
172 
173  // Mapping from physical id -> (physical dim, physical name) pairs.
174  // These can refer to either "sidesets" or "subdomains"; we need to
175  // wait until the Mesh has been read to know which is which. Note
176  // that we are using 'int' as the key here rather than
177  // subdomain_id_type or boundary_id_type, since at this point, it
178  // could be either.
179  typedef std::pair<unsigned, std::string> GmshPhysical;
180  std::map<int, GmshPhysical> gmsh_physicals;
181 
182  // map to hold the node numbers for translation
183  // note the the nodes can be non-consecutive
184  std::map<unsigned int, unsigned int> nodetrans;
185 
186  // For reading the file line by line
187  std::string s;
188 
189  while (true)
190  {
191  // Try to read something. This may set EOF!
192  std::getline(in, s);
193 
194  if (in)
195  {
196  // Process s...
197 
198  if (s.find("$MeshFormat") == static_cast<std::string::size_type>(0))
199  {
200  in >> version >> format >> size;
201  if ((version != 2.0) && (version != 2.1) && (version != 2.2))
202  {
203  // Some notes on gmsh mesh versions:
204  //
205  // Mesh version 2.0 goes back as far as I know. It's not explicitly
206  // mentioned here: http://www.geuz.org/gmsh/doc/VERSIONS.txt
207  //
208  // As of gmsh-2.4.0:
209  // bumped mesh version format to 2.1 (small change in the $PhysicalNames
210  // section, where the group dimension is now required);
211  // [Since we don't even parse the PhysicalNames section at the time
212  // of this writing, I don't think this change affects us.]
213  //
214  // Mesh version 2.2 tested by Manav Bhatia; no other
215  // libMesh code changes were required for support
216  libmesh_error_msg("Error: Unknown msh file version " << version);
217  }
218 
219  if (format)
220  libmesh_error_msg("Error: Unknown data format for mesh in Gmsh reader.");
221  }
222 
223  // Read and process the "PhysicalNames" section.
224  else if (s.find("$PhysicalNames") == static_cast<std::string::size_type>(0))
225  {
226  // The lines in the PhysicalNames section should look like the following:
227  // 2 1 "frac" lower_dimensional_block
228  // 2 3 "top"
229  // 2 4 "bottom"
230  // 3 2 "volume"
231 
232  // Read in the number of physical groups to expect in the file.
233  unsigned int num_physical_groups = 0;
234  in >> num_physical_groups;
235 
236  // Read rest of line including newline character.
237  std::getline(in, s);
238 
239  for (unsigned int i=0; i<num_physical_groups; ++i)
240  {
241  // Read an entire line of the PhysicalNames section.
242  std::getline(in, s);
243 
244  // Use an istringstream to extract the physical
245  // dimension, physical id, and physical name from
246  // this line.
247  std::istringstream s_stream(s);
248  unsigned phys_dim;
249  int phys_id;
250  std::string phys_name;
251  s_stream >> phys_dim >> phys_id >> phys_name;
252 
253  // Not sure if this is true for all Gmsh files, but
254  // my test file has quotes around the phys_name
255  // string. So let's erase any quotes now...
256  phys_name.erase(std::remove(phys_name.begin(), phys_name.end(), '"'), phys_name.end());
257 
258  // Record this ID for later assignment of subdomain/sideset names.
259  gmsh_physicals[phys_id] = std::make_pair(phys_dim, phys_name);
260 
261  // If 's' also contains the libmesh-specific string
262  // "lower_dimensional_block", add this block ID to
263  // the list of blocks which are not boundary
264  // conditions.
265  if (s.find("lower_dimensional_block") != std::string::npos)
266  {
267  lower_dimensional_blocks.insert(cast_int<subdomain_id_type>(phys_id));
268 
269  // The user has explicitly told us that this
270  // block is a subdomain, so set that association
271  // in the Mesh.
272  mesh.subdomain_name(cast_int<subdomain_id_type>(phys_id)) = phys_name;
273  }
274  }
275  }
276 
277  // read the node block
278  else if (s.find("$NOD") == static_cast<std::string::size_type>(0) ||
279  s.find("$NOE") == static_cast<std::string::size_type>(0) ||
280  s.find("$Nodes") == static_cast<std::string::size_type>(0))
281  {
282  unsigned int num_nodes = 0;
283  in >> num_nodes;
284  mesh.reserve_nodes (num_nodes);
285 
286  // read in the nodal coordinates and form points.
287  Real x, y, z;
288  unsigned int id;
289 
290  // add the nodal coordinates to the mesh
291  for (unsigned int i=0; i<num_nodes; ++i)
292  {
293  in >> id >> x >> y >> z;
294  mesh.add_point (Point(x, y, z), i);
295  nodetrans[id] = i;
296  }
297 
298  // read the $ENDNOD delimiter
299  std::getline(in, s);
300  }
301 
302 
303  // Read the element block
304  else if (s.find("$ELM") == static_cast<std::string::size_type>(0) ||
305  s.find("$Elements") == static_cast<std::string::size_type>(0))
306  {
307  // For reading the number of elements and the node ids from the stream
308  unsigned int
309  num_elem = 0,
310  node_id = 0;
311 
312  // read how many elements are there, and reserve space in the mesh
313  in >> num_elem;
314  mesh.reserve_elem (num_elem);
315 
316  // As of version 2.2, the format for each element line is:
317  // elm-number elm-type number-of-tags < tag > ... node-number-list
318  // From the Gmsh docs:
319  // * the first tag is the number of the
320  // physical entity to which the element belongs
321  // * the second is the number of the elementary geometrical
322  // entity to which the element belongs
323  // * the third is the number of mesh partitions to which the element
324  // belongs
325  // * The rest of the tags are the partition ids (negative
326  // partition ids indicate ghost cells). A zero tag is
327  // equivalent to no tag. Gmsh and most codes using the
328  // MSH 2 format require at least the first two tags
329  // (physical and elementary tags).
330 
331  // Keep track of all the element dimensions seen
332  std::vector<unsigned> elem_dimensions_seen(3);
333 
334  // read the elements
335  for (unsigned int iel=0; iel<num_elem; ++iel)
336  {
337  unsigned int
338  id, type,
339  physical=1, elementary=1,
340  nnodes=0, ntags;
341 
342  // Note: tag has to be an int because it could be negative,
343  // see above.
344  int tag;
345 
346  if (version <= 1.0)
347  in >> id >> type >> physical >> elementary >> nnodes;
348 
349  else
350  {
351  in >> id >> type >> ntags;
352 
353  if (ntags > 2)
354  libmesh_do_once(libMesh::err << "Warning, ntags=" << ntags << ", but we currently only support reading 2 flags." << std::endl;);
355 
356  for (unsigned int j = 0; j < ntags; j++)
357  {
358  in >> tag;
359  if (j == 0)
360  physical = tag;
361  else if (j == 1)
362  elementary = tag;
363  }
364  }
365 
366  // Consult the import element table to determine which element to build
367  std::map<unsigned int, GmshIO::ElementDefinition>::iterator eletypes_it = _element_maps.in.find(type);
368 
369  // Make sure we actually found something
370  if (eletypes_it == _element_maps.in.end())
371  libmesh_error_msg("Element type " << type << " not found!");
372 
373  // Get a reference to the ElementDefinition
374  const GmshIO::ElementDefinition & eletype = eletypes_it->second;
375 
376  // If we read nnodes, make sure it matches the number in eletype.nnodes
377  if (nnodes != 0 && nnodes != eletype.nnodes)
378  libmesh_error_msg("nnodes = " << nnodes << " and eletype.nnodes = " << eletype.nnodes << " do not match.");
379 
380  // Assign the value from the eletype object.
381  nnodes = eletype.nnodes;
382 
383  // Don't add 0-dimensional "point" elements to the
384  // Mesh. They should *always* be treated as boundary
385  // "nodeset" data.
386  if (eletype.dim > 0)
387  {
388  // Record this element dimension as being "seen".
389  // We will treat all elements with dimension <
390  // max(dimension) as specifying boundary conditions,
391  // but we won't know what max_elem_dimension_seen is
392  // until we read the entire file.
393  elem_dimensions_seen[eletype.dim-1] = 1;
394 
395  // Add the element to the mesh
396  {
397  Elem * elem = Elem::build(eletype.type).release();
398  elem->set_id(iel);
399  elem = mesh.add_elem(elem);
400 
401  // Make sure that the libmesh element we added has nnodes nodes.
402  if (elem->n_nodes() != nnodes)
403  libmesh_error_msg("Number of nodes for element " \
404  << id \
405  << " of type " << eletype.type \
406  << " (Gmsh type " << type \
407  << ") does not match Libmesh definition. " \
408  << "I expected " << elem->n_nodes() \
409  << " nodes, but got " << nnodes);
410 
411  // Add node pointers to the elements.
412  // If there is a node translation table, use it.
413  if (eletype.nodes.size() > 0)
414  for (unsigned int i=0; i<nnodes; i++)
415  {
416  in >> node_id;
417  elem->set_node(eletype.nodes[i]) = mesh.node_ptr(nodetrans[node_id]);
418  }
419  else
420  {
421  for (unsigned int i=0; i<nnodes; i++)
422  {
423  in >> node_id;
424  elem->set_node(i) = mesh.node_ptr(nodetrans[node_id]);
425  }
426  }
427 
428  // Finally, set the subdomain ID to physical. If this is a lower-dimension element, this ID will
429  // eventually go into the Mesh's BoundaryInfo object.
430  elem->subdomain_id() = static_cast<subdomain_id_type>(physical);
431  }
432  }
433 
434  // Handle 0-dimensional elements (points) by adding
435  // them to the BoundaryInfo object with
436  // boundary_id == physical.
437  else
438  {
439  // This seems like it should always be the same
440  // number as the 'id' we already read in on this
441  // line. At least it was in the example gmsh
442  // file I had...
443  in >> node_id;
444  mesh.get_boundary_info().add_node
445  (nodetrans[node_id],
446  static_cast<boundary_id_type>(physical));
447  }
448  } // element loop
449 
450  // read the $ENDELM delimiter
451  std::getline(in, s);
452 
453  // Record the max and min element dimension seen while reading the file.
454  unsigned char
455  max_elem_dimension_seen=1,
456  min_elem_dimension_seen=3;
457 
458  for (std::size_t i=0; i<elem_dimensions_seen.size(); ++i)
459  if (elem_dimensions_seen[i])
460  {
461  // Debugging
462  // libMesh::out << "Seen elements of dimension " << i+1 << std::endl;
463  max_elem_dimension_seen =
464  std::max(max_elem_dimension_seen, cast_int<unsigned char>(i+1));
465  min_elem_dimension_seen =
466  std::min(min_elem_dimension_seen, cast_int<unsigned char>(i+1));
467  }
468 
469  // Debugging:
470  // libMesh::out << "max_elem_dimension_seen=" << max_elem_dimension_seen << std::endl;
471  // libMesh::out << "min_elem_dimension_seen=" << min_elem_dimension_seen << std::endl;
472 
473  // If the difference between the max and min element dimension seen is larger than
474  // 1, (e.g. the file has 1D and 3D elements only) we don't handle this case.
475  if (max_elem_dimension_seen - min_elem_dimension_seen > 1)
476  libmesh_error_msg("Cannot handle meshes with dimension mismatch greater than 1.");
477 
478  // How many different element dimensions did we see while reading from file?
479  unsigned n_dims_seen = std::accumulate(elem_dimensions_seen.begin(),
480  elem_dimensions_seen.end(),
481  static_cast<unsigned>(0),
482  std::plus<unsigned>());
483 
484  // Have not yet tested a case where 1, 2, and 3D elements are all in the same Mesh,
485  // though it should theoretically be possible to handle.
486  if (n_dims_seen == 3)
487  libmesh_error_msg("Reading meshes with 1, 2, and 3D elements not currently supported.");
488 
489  // Set mesh_dimension based on the largest element dimension seen.
490  mesh.set_mesh_dimension(max_elem_dimension_seen);
491 
492  // Now that we know the maximum element dimension seen,
493  // we know whether the physical names are subdomain
494  // names or sideset names.
495  {
496  std::map<int, GmshPhysical>::iterator it = gmsh_physicals.begin();
497  for (; it != gmsh_physicals.end(); ++it)
498  {
499  // Extract data
500  int phys_id = it->first;
501  unsigned phys_dim = it->second.first;
502  std::string phys_name = it->second.second;
503 
504  // If the physical's dimension matches the largest
505  // dimension we've seen, it's a subdomain name.
506  if (phys_dim == max_elem_dimension_seen)
507  mesh.subdomain_name(cast_int<subdomain_id_type>(phys_id)) = phys_name;
508 
509  // Otherwise, if it's not a lower-dimensional
510  // block, it's a sideset name.
511  else if (phys_dim < max_elem_dimension_seen &&
512  !lower_dimensional_blocks.count(cast_int<boundary_id_type>(phys_id)))
513  mesh.get_boundary_info().sideset_name(cast_int<boundary_id_type>(phys_id)) = phys_name;
514  }
515  }
516 
517  if (n_dims_seen > 1)
518  {
519  // Store lower-dimensional elements in a map sorted
520  // by Elem::key(). We use a multimap for two reasons:
521  // 1.) The hash function is not guaranteed to be
522  // unique, so different lower-dimensional elements
523  // could theoretically hash to the same value,
524  // although this is pretty unlikely.
525  // 2.) The Gmsh file may contain multiple
526  // lower-dimensional elements for a single side in
527  // order to implement multiple boundary ids for a
528  // single side. These lower-dimensional elements
529  // will all hash to the same value, and we need to
530  // be able to store all of them.
531  typedef LIBMESH_BEST_UNORDERED_MULTIMAP<dof_id_type, Elem *> provide_container_t;
532  provide_container_t provide_bcs;
533 
534  // 1st loop over active elements - get info about lower-dimensional elements.
535  {
536  MeshBase::element_iterator it = mesh.active_elements_begin();
537  const MeshBase::element_iterator end = mesh.active_elements_end();
538  for ( ; it != end; ++it)
539  {
540  Elem * elem = *it;
541 
542  if (elem->dim() < max_elem_dimension_seen &&
543  !lower_dimensional_blocks.count(elem->subdomain_id()))
544  {
545  // To be consistent with the previous
546  // GmshIO behavior, add all the
547  // lower-dimensional elements' nodes to
548  // the Mesh's BoundaryInfo object with the
549  // lower-dimensional element's subdomain
550  // ID.
551  for (unsigned n=0; n<elem->n_nodes(); n++)
552  mesh.get_boundary_info().add_node(elem->node_id(n),
553  elem->subdomain_id());
554 
555  // Store this elem in a quickly-searchable
556  // container to use it to assign boundary
557  // conditions later.
558  provide_bcs.insert(std::make_pair(elem->key(), elem));
559  }
560  }
561  }
562 
563  // 2nd loop over active elements - use lower dimensional element data to set BCs for higher dimensional elements
564  {
565  MeshBase::element_iterator it = mesh.active_elements_begin();
566  const MeshBase::element_iterator end = mesh.active_elements_end();
567 
568  for ( ; it != end; ++it)
569  {
570  Elem * elem = *it;
571 
572  if (elem->dim() == max_elem_dimension_seen)
573  {
574  // This is a max-dimension element that
575  // may require BCs. For each of its
576  // sides, including internal sides, we'll
577  // see if one more more lower-dimensional elements
578  // provides boundary information for it.
579  // Note that we have not yet called
580  // find_neighbors(), so we can't use
581  // elem->neighbor(sn) in this algorithm...
582  for (unsigned short sn=0; sn<elem->n_sides(); sn++)
583  {
584  // Look for the current side in the provide_bcs multimap.
585  std::pair<provide_container_t::iterator,
586  provide_container_t::iterator>
587  rng = provide_bcs.equal_range(elem->key(sn));
588 
589  for (provide_container_t::iterator iter = rng.first;
590  iter != rng.second; ++iter)
591  {
592  // Construct the side for hash verification.
593  UniquePtr<Elem> side (elem->build_side_ptr(sn));
594 
595  // Construct the lower-dimensional element to compare to the side.
596  Elem * lower_dim_elem = iter->second;
597 
598  // This was a hash, so it might not be perfect. Let's verify...
599  if (*lower_dim_elem == *side)
600  {
601  // Add the lower-dimensional
602  // element's subdomain_id as a
603  // boundary_id for the
604  // higher-dimensional element.
605  boundary_id_type bid = cast_int<boundary_id_type>(lower_dim_elem->subdomain_id());
606  mesh.get_boundary_info().add_side(elem, sn, bid);
607  }
608  }
609  }
610  }
611  }
612  } // end 2nd loop over active elements
613 
614  // 3rd loop over active elements - Remove the lower-dimensional elements
615  {
616  MeshBase::element_iterator it = mesh.active_elements_begin();
617  const MeshBase::element_iterator end = mesh.active_elements_end();
618  for ( ; it != end; ++it)
619  {
620  Elem * elem = *it;
621 
622  if (elem->dim() < max_elem_dimension_seen &&
623  !lower_dimensional_blocks.count(elem->subdomain_id()))
624  mesh.delete_elem(elem);
625  }
626  } // end 3rd loop over active elements
627  } // end if (n_dims_seen > 1)
628  } // if $ELM
629 
630  continue;
631  } // if (in)
632 
633 
634  // If !in, check to see if EOF was set. If so, break out
635  // of while loop.
636  if (in.eof())
637  break;
638 
639  // If !in and !in.eof(), stream is in a bad state!
640  libmesh_error_msg("Stream is bad! Perhaps the file does not exist?");
641 
642  } // while true
643 }
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:234
TestClass subdomain_id_type
Definition: id_types.h:43
unsigned short int side
Definition: xdr_io.C:49
IterBase * end
const MeshBase & mesh() const
std::map< unsigned int, ElementDefinition > in
Definition: gmsh_io.h:186
long double max(long double a, double b)
libmesh_assert(j)
PetscErrorCode Vec x
int8_t boundary_id_type
Definition: id_types.h:51
OStreamProxy err(std::cerr)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static ElementMaps _element_maps
Definition: gmsh_io.h:193
processor_id_type processor_id()
Definition: libmesh_base.h:96
long double min(long double a, double b)
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().

void libMesh::GmshIO::write ( const std::string &  name)
virtual

This method implements writing a mesh to a specified file in the Gmsh *.msh format.

Implements libMesh::MeshOutput< MeshBase >.

Definition at line 647 of file gmsh_io.C.

References libMesh::processor_id(), and write_mesh().

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

648 {
650  {
651  // Open the output file stream
652  std::ofstream out_stream (name.c_str());
653 
654  // Make sure it opened correctly
655  if (!out_stream.good())
656  libmesh_file_error(name.c_str());
657 
658  this->write_mesh (out_stream);
659  }
660 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
const MeshBase & mesh() const
void write_mesh(std::ostream &out)
Definition: gmsh_io.C:676
processor_id_type processor_id()
Definition: libmesh_base.h:96
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().

bool & libMesh::GmshIO::write_lower_dimensional_elements ( )

Access to the flag which controls whether boundary elements are written to the Mesh file.

Definition at line 137 of file gmsh_io.C.

References _write_lower_dimensional_elements.

Referenced by write_mesh().

138 {
140 }
bool _write_lower_dimensional_elements
Definition: gmsh_io.h:148
void libMesh::GmshIO::write_mesh ( std::ostream &  out)
private

This method implements writing a mesh to a specified file. This will write an ASCII *.msh file.

Definition at line 676 of file gmsh_io.C.

References _element_maps, libMesh::Elem::build_side_ptr(), end, libMesh::GmshIO::ElementDefinition::gmsh_type, libMesh::DofObject::id(), libMesh::MeshTools::Generation::Private::idx(), libMesh::libmesh_assert(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::Elem::n_nodes(), libMesh::Elem::node_id(), libMesh::GmshIO::ElementDefinition::nodes, libMesh::GmshIO::ElementMaps::out, libMesh::DofObject::processor_id(), libMesh::Real, side, libMesh::Elem::subdomain_id(), libMesh::Elem::type(), and write_lower_dimensional_elements().

Referenced by write().

677 {
678  // Be sure that the stream is valid.
679  libmesh_assert (out_stream.good());
680 
681  // Get a const reference to the mesh
682  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
683 
684  // If requested, write out lower-dimensional elements for
685  // element-side-based boundary conditions.
686  unsigned int n_boundary_faces = 0;
688  n_boundary_faces = mesh.get_boundary_info().n_boundary_conds();
689 
690  // Note: we are using version 2.0 of the gmsh output format.
691 
692  // Write the file header.
693  out_stream << "$MeshFormat\n";
694  out_stream << "2.0 0 " << sizeof(Real) << '\n';
695  out_stream << "$EndMeshFormat\n";
696 
697  // write the nodes in (n x y z) format
698  out_stream << "$Nodes\n";
699  out_stream << mesh.n_nodes() << '\n';
700 
701  for (unsigned int v=0; v<mesh.n_nodes(); v++)
702  out_stream << mesh.node_ref(v).id()+1 << " "
703  << mesh.node_ref(v)(0) << " "
704  << mesh.node_ref(v)(1) << " "
705  << mesh.node_ref(v)(2) << '\n';
706  out_stream << "$EndNodes\n";
707 
708  {
709  // write the connectivity
710  out_stream << "$Elements\n";
711  out_stream << mesh.n_active_elem() + n_boundary_faces << '\n';
712 
713  MeshBase::const_element_iterator it = mesh.active_elements_begin();
714  const MeshBase::const_element_iterator end = mesh.active_elements_end();
715 
716  // loop over the elements
717  for ( ; it != end; ++it)
718  {
719  const Elem * elem = *it;
720 
721  // Make sure we have a valid entry for
722  // the current element type.
723  libmesh_assert (_element_maps.out.count(elem->type()));
724 
725  // consult the export element table
726  std::map<ElemType, ElementDefinition>::iterator def_it =
727  _element_maps.out.find(elem->type());
728 
729  // Assert that we found it
730  if (def_it == _element_maps.out.end())
731  libmesh_error_msg("Element type " << elem->type() << " not found in _element_maps.");
732 
733  // Get a reference to the ElementDefinition object
734  const ElementDefinition & eletype = def_it->second;
735 
736  // The element mapper better not require any more nodes
737  // than are present in the current element!
738  libmesh_assert_less_equal (eletype.nodes.size(), elem->n_nodes());
739 
740  // elements ids are 1 based in Gmsh
741  out_stream << elem->id()+1 << " ";
742  // element type
743  out_stream << eletype.gmsh_type;
744 
745  // write the number of tags (3) and their values:
746  // 1 (physical entity)
747  // 2 (geometric entity)
748  // 3 (partition entity)
749  out_stream << " 3 "
750  << static_cast<unsigned int>(elem->subdomain_id())
751  << " 0 "
752  << elem->processor_id()+1
753  << " ";
754 
755  // if there is a node translation table, use it
756  if (eletype.nodes.size() > 0)
757  for (unsigned int i=0; i < elem->n_nodes(); i++)
758  out_stream << elem->node_id(eletype.nodes[i])+1 << " "; // gmsh is 1-based
759  // otherwise keep the same node order
760  else
761  for (unsigned int i=0; i < elem->n_nodes(); i++)
762  out_stream << elem->node_id(i)+1 << " "; // gmsh is 1-based
763  out_stream << "\n";
764  } // element loop
765  }
766 
767  {
768  // A counter for writing surface elements to the Gmsh file
769  // sequentially. We start numbering them with a number strictly
770  // larger than the largest element ID in the mesh. Note: the
771  // MeshBase docs say "greater than or equal to" the maximum
772  // element id in the mesh, so technically we might need a +1 here,
773  // but all of the implementations return an ID strictly greater
774  // than the largest element ID in the Mesh.
775  unsigned int e_id = mesh.max_elem_id();
776 
777  // loop over the elements, writing out boundary faces
778  MeshBase::const_element_iterator it = mesh.active_elements_begin();
779  const MeshBase::const_element_iterator end = mesh.active_elements_end();
780 
781  if (n_boundary_faces)
782  {
783  // Construct the list of boundary sides
784  std::vector<dof_id_type> element_id_list;
785  std::vector<unsigned short int> side_list;
786  std::vector<boundary_id_type> bc_id_list;
787 
788  mesh.get_boundary_info().build_side_list(element_id_list, side_list, bc_id_list);
789 
790  // Loop over these lists, writing data to the file.
791  for (std::size_t idx=0; idx<element_id_list.size(); ++idx)
792  {
793  const Elem & elem = mesh.elem_ref(element_id_list[idx]);
794 
795  UniquePtr<const Elem> side = elem.build_side_ptr(side_list[idx]);
796 
797  // Map from libmesh elem type to gmsh elem type.
798  std::map<ElemType, ElementDefinition>::iterator def_it =
799  _element_maps.out.find(side->type());
800 
801  // If we didn't find it, that's an error
802  if (def_it == _element_maps.out.end())
803  libmesh_error_msg("Element type " << side->type() << " not found in _element_maps.");
804 
805  // consult the export element table
806  const GmshIO::ElementDefinition & eletype = def_it->second;
807 
808  // The element mapper better not require any more nodes
809  // than are present in the current element!
810  libmesh_assert_less_equal (eletype.nodes.size(), side->n_nodes());
811 
812  // elements ids are 1-based in Gmsh
813  out_stream << e_id+1 << " ";
814 
815  // element type
816  out_stream << eletype.gmsh_type;
817 
818  // write the number of tags:
819  // 1 (physical entity)
820  // 2 (geometric entity)
821  // 3 (partition entity)
822  out_stream << " 3 "
823  << bc_id_list[idx]
824  << " 0 "
825  << elem.processor_id()+1
826  << " ";
827 
828  // if there is a node translation table, use it
829  if (eletype.nodes.size() > 0)
830  for (unsigned int i=0; i < side->n_nodes(); i++)
831  out_stream << side->node_id(eletype.nodes[i])+1 << " "; // gmsh is 1-based
832 
833  // otherwise keep the same node order
834  else
835  for (unsigned int i=0; i < side->n_nodes(); i++)
836  out_stream << side->node_id(i)+1 << " "; // gmsh is 1-based
837 
838  // Go to the next line
839  out_stream << "\n";
840 
841  // increment this index too...
842  ++e_id;
843  }
844  }
845 
846  out_stream << "$EndElements\n";
847  }
848 }
bool & write_lower_dimensional_elements()
Definition: gmsh_io.C:137
unsigned short int side
Definition: xdr_io.C:49
IterBase * end
const MT & mesh() const
Definition: mesh_output.h:216
libmesh_assert(j)
std::map< ElemType, ElementDefinition > out
Definition: gmsh_io.h:185
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static ElementMaps _element_maps
Definition: gmsh_io.h:193
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int i, const unsigned int j)
void libMesh::GmshIO::write_nodal_data ( const std::string &  fname,
const std::vector< Number > &  soln,
const std::vector< std::string > &  names 
)
virtual

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

Reimplemented from libMesh::MeshOutput< MeshBase >.

Definition at line 664 of file gmsh_io.C.

References libMesh::processor_id(), and write_post().

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

667 {
668  LOG_SCOPE("write_nodal_data()", "GmshIO");
669 
671  this->write_post (fname, &soln, &names);
672 }
void write_post(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
Definition: gmsh_io.C:852
const MeshBase & mesh() const
processor_id_type processor_id()
Definition: libmesh_base.h:96
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.

void libMesh::GmshIO::write_post ( const std::string &  fname,
const std::vector< Number > *  v = libmesh_nullptr,
const std::vector< std::string > *  solution_names = libmesh_nullptr 
)
private

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are optionally provided. This will write an ASCII or binary *.pos file, depending on the binary flag.

Definition at line 852 of file gmsh_io.C.

References libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), binary(), libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, end, libMesh::err, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libmesh_nullptr, libMesh::libmesh_real(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::n_nodes(), n_vars, libMesh::Elem::n_vertices(), libMesh::Elem::node_id(), libMesh::out, libMesh::Elem::point(), libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::processor_id(), libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::TET10, libMesh::TET4, libMesh::TRI3, and libMesh::TRI6.

Referenced by write_nodal_data().

855 {
856 
857  // Should only do this on processor 0!
858  libmesh_assert_equal_to (MeshOutput<MeshBase>::mesh().processor_id(), 0);
859 
860  // Create an output stream
861  std::ofstream out_stream(fname.c_str());
862 
863  // Make sure it opened correctly
864  if (!out_stream.good())
865  libmesh_file_error(fname.c_str());
866 
867  // create a character buffer
868  char buf[80];
869 
870  // Get a constant reference to the mesh.
871  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
872 
873  // write the data
874  if ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr))
875  {
876  const unsigned int n_vars =
877  cast_int<unsigned int>(solution_names->size());
878 
879  if (!(v->size() == mesh.n_nodes()*n_vars))
880  libMesh::err << "ERROR: v->size()=" << v->size()
881  << ", mesh.n_nodes()=" << mesh.n_nodes()
882  << ", n_vars=" << n_vars
883  << ", mesh.n_nodes()*n_vars=" << mesh.n_nodes()*n_vars
884  << "\n";
885 
886  libmesh_assert_equal_to (v->size(), mesh.n_nodes()*n_vars);
887 
888  // write the header
889  out_stream << "$PostFormat\n";
890  if (this->binary())
891  out_stream << "1.2 1 " << sizeof(double) << "\n";
892  else
893  out_stream << "1.2 0 " << sizeof(double) << "\n";
894  out_stream << "$EndPostFormat\n";
895 
896  // Loop over the elements to see how much of each type there are
897  unsigned int n_points=0, n_lines=0, n_triangles=0, n_quadrangles=0,
898  n_tetrahedra=0, n_hexahedra=0, n_prisms=0, n_pyramids=0;
899  unsigned int n_scalar=0, n_vector=0, n_tensor=0;
900  unsigned int nb_text2d=0, nb_text2d_chars=0, nb_text3d=0, nb_text3d_chars=0;
901 
902  {
903  MeshBase::const_element_iterator it = mesh.active_elements_begin();
904  const MeshBase::const_element_iterator end = mesh.active_elements_end();
905 
906 
907  for ( ; it != end; ++it)
908  {
909  const ElemType elemtype = (*it)->type();
910 
911  switch (elemtype)
912  {
913  case EDGE2:
914  case EDGE3:
915  case EDGE4:
916  {
917  n_lines += 1;
918  break;
919  }
920  case TRI3:
921  case TRI6:
922  {
923  n_triangles += 1;
924  break;
925  }
926  case QUAD4:
927  case QUAD8:
928  case QUAD9:
929  {
930  n_quadrangles += 1;
931  break;
932  }
933  case TET4:
934  case TET10:
935  {
936  n_tetrahedra += 1;
937  break;
938  }
939  case HEX8:
940  case HEX20:
941  case HEX27:
942  {
943  n_hexahedra += 1;
944  break;
945  }
946  case PRISM6:
947  case PRISM15:
948  case PRISM18:
949  {
950  n_prisms += 1;
951  break;
952  }
953  case PYRAMID5:
954  {
955  n_pyramids += 1;
956  break;
957  }
958  default:
959  libmesh_error_msg("ERROR: Nonexistent element type " << (*it)->type());
960  }
961  }
962  }
963 
964  // create a view for each variable
965  for (unsigned int ivar=0; ivar < n_vars; ivar++)
966  {
967  std::string varname = (*solution_names)[ivar];
968 
969  // at the moment, we just write out scalar quantities
970  // later this should be made configurable through
971  // options to the writer class
972  n_scalar = 1;
973 
974  // write the variable as a view, and the number of time steps
975  out_stream << "$View\n" << varname << " " << 1 << "\n";
976 
977  // write how many of each geometry type are written
978  out_stream << n_points * n_scalar << " "
979  << n_points * n_vector << " "
980  << n_points * n_tensor << " "
981  << n_lines * n_scalar << " "
982  << n_lines * n_vector << " "
983  << n_lines * n_tensor << " "
984  << n_triangles * n_scalar << " "
985  << n_triangles * n_vector << " "
986  << n_triangles * n_tensor << " "
987  << n_quadrangles * n_scalar << " "
988  << n_quadrangles * n_vector << " "
989  << n_quadrangles * n_tensor << " "
990  << n_tetrahedra * n_scalar << " "
991  << n_tetrahedra * n_vector << " "
992  << n_tetrahedra * n_tensor << " "
993  << n_hexahedra * n_scalar << " "
994  << n_hexahedra * n_vector << " "
995  << n_hexahedra * n_tensor << " "
996  << n_prisms * n_scalar << " "
997  << n_prisms * n_vector << " "
998  << n_prisms * n_tensor << " "
999  << n_pyramids * n_scalar << " "
1000  << n_pyramids * n_vector << " "
1001  << n_pyramids * n_tensor << " "
1002  << nb_text2d << " "
1003  << nb_text2d_chars << " "
1004  << nb_text3d << " "
1005  << nb_text3d_chars << "\n";
1006 
1007  // if binary, write a marker to identify the endianness of the file
1008  if (this->binary())
1009  {
1010  const int one = 1;
1011  std::memcpy(buf, &one, sizeof(int));
1012  out_stream.write(buf, sizeof(int));
1013  }
1014 
1015  // the time steps (there is just 1 at the moment)
1016  if (this->binary())
1017  {
1018  double one = 1;
1019  std::memcpy(buf, &one, sizeof(double));
1020  out_stream.write(buf, sizeof(double));
1021  }
1022  else
1023  out_stream << "1\n";
1024 
1025  // Loop over the elements and write out the data
1026  MeshBase::const_element_iterator it = mesh.active_elements_begin();
1027  const MeshBase::const_element_iterator end = mesh.active_elements_end();
1028 
1029  for ( ; it != end; ++it)
1030  {
1031  const Elem * elem = *it;
1032 
1033  // this is quite crappy, but I did not invent that file format!
1034  for (unsigned int d=0; d<3; d++) // loop over the dimensions
1035  {
1036  for (unsigned int n=0; n < elem->n_vertices(); n++) // loop over vertices
1037  {
1038  const Point & vertex = elem->point(n);
1039  if (this->binary())
1040  {
1041  double tmp = vertex(d);
1042  std::memcpy(buf, &tmp, sizeof(double));
1043  out_stream.write(reinterpret_cast<char *>(buf), sizeof(double));
1044  }
1045  else
1046  out_stream << vertex(d) << " ";
1047  }
1048  if (!this->binary())
1049  out_stream << "\n";
1050  }
1051 
1052  // now finally write out the data
1053  for (unsigned int i=0; i < elem->n_vertices(); i++) // loop over vertices
1054  if (this->binary())
1055  {
1056 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1057  libMesh::out << "WARNING: Gmsh::write_post does not fully support "
1058  << "complex numbers. Will only write the real part of "
1059  << "variable " << varname << std::endl;
1060 #endif
1061  double tmp = libmesh_real((*v)[elem->node_id(i)*n_vars + ivar]);
1062  std::memcpy(buf, &tmp, sizeof(double));
1063  out_stream.write(reinterpret_cast<char *>(buf), sizeof(double));
1064  }
1065  else
1066  {
1067 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1068  libMesh::out << "WARNING: Gmsh::write_post does not fully support "
1069  << "complex numbers. Will only write the real part of "
1070  << "variable " << varname << std::endl;
1071 #endif
1072  out_stream << libmesh_real((*v)[elem->node_id(i)*n_vars + ivar]) << "\n";
1073  }
1074  }
1075  if (this->binary())
1076  out_stream << "\n";
1077  out_stream << "$EndView\n";
1078 
1079  } // end variable loop (writing the views)
1080  }
1081 }
T libmesh_real(T a)
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
const MeshBase & mesh() const
bool & binary()
Definition: gmsh_io.C:130
const unsigned int n_vars
Definition: tecplot_io.C:68
OStreamProxy err(std::cerr)
OStreamProxy out(std::cout)
processor_id_type processor_id()
Definition: libmesh_base.h:96

Member Data Documentation

bool libMesh::GmshIO::_binary
private

Flag to write binary data.

Definition at line 142 of file gmsh_io.h.

Referenced by binary().

GmshIO::ElementMaps libMesh::GmshIO::_element_maps = GmshIO::build_element_maps()
staticprivate

A static ElementMaps object that is built statically and used by all instances of this class.

Definition at line 193 of file gmsh_io.h.

Referenced by read_mesh(), and write_mesh().

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().

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.

bool libMesh::GmshIO::_write_lower_dimensional_elements
private

If true, lower-dimensional elements based on the boundary conditions get written to the output file.

Definition at line 148 of file gmsh_io.h.

Referenced by write_lower_dimensional_elements().


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