#include <gmsh_io.h>
Inherits MeshInput< MeshBase >, and MeshOutput< MeshBase >.
GmshIO (MeshBase &mesh)
GmshIO (const MeshBase &mesh)
virtual void read (const std::string &name)
virtual void write (const std::string &name)
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
bool & binary ()
virtual void write_equation_systems (const std::string &, const EquationSystems &)
MeshBase & mesh ()
void skip_comment_lines (std::istream &in, const char comment_start)
const MeshBase & mesh () const
virtual void read_mesh (std::istream &in)
virtual void write_mesh (std::ostream &out)
void write_post (const std::string &, const std::vector< Number > *=NULL, const std::vector< std::string > *=NULL)
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:
Martin Lüthi, 2005 (support for reading meshes and writing results)
Definition at line 47 of file gmsh_io.h.
Definition at line 142 of file gmsh_io.h.
:
MeshInput<MeshBase> (mesh),
MeshOutput<MeshBase> (mesh),
_binary (false)
{}
Definition at line 134 of file gmsh_io.h.
:
MeshOutput<MeshBase> (mesh),
_binary (false)
{
}
Definition at line 149 of file gmsh_io.h.
References _binary.
Referenced by write_post().
{
return _binary;
}
Referenced by GMVIO::_read_materials(), GMVIO::_read_nodes(), GMVIO::_read_one_cell(), GMVIO::add_cell_centered_data(), GMVIO::copy_nodal_solution(), UNVIO::element_in(), TetGenIO::element_in(), UNVIO::element_out(), UNVIO::node_in(), TetGenIO::node_in(), UNVIO::node_out(), XdrIO::read(), VTKIO::read(), GMVIO::read(), ExodusII_IO::read(), LegacyXdrIO::read_ascii(), LegacyXdrIO::read_binary(), UCDIO::read_implementation(), LegacyXdrIO::read_mesh(), read_mesh(), XdrIO::read_serialized_bcs(), XdrIO::read_serialized_connectivity(), XdrIO::read_serialized_nodes(), OFFIO::read_stream(), MatlabIO::read_stream(), VTKIO::solution_to_vtk(), XdrIO::write(), VTKIO::write(), TetGenIO::write(), GMVIO::write_ascii_new_impl(), GMVIO::write_ascii_old_impl(), GMVIO::write_binary(), GMVIO::write_discontinuous_gmv(), UNVIO::write_implementation(), UCDIO::write_implementation(), LegacyXdrIO::write_mesh(), write_mesh(), write_post(), XdrIO::write_serialized_bcs(), XdrIO::write_serialized_connectivity(), XdrIO::write_serialized_nodes(), and LegacyXdrIO::write_soln().
Referenced by PostscriptIO::write(), FroIO::write(), MEDITIO::write_ascii(), EnsightIO::write_geometry_ascii(), EnsightIO::write_scalar_ascii(), GnuPlotIO::write_solution(), DivaIO::write_stream(), and EnsightIO::write_vector_ascii().
Implements MeshInput< MeshBase >.
Definition at line 312 of file gmsh_io.C.
References read_mesh().
Referenced by UnstructuredMesh::read().
{
std::ifstream in (name.c_str());
this->read_mesh (in);
}
Read the element block
If the element dimension is smaller than the mesh dimension, this is a boundary element and will be added to mesh.boundary_info.
Because the elements might not yet exist, the sides are put on hold until the elements are created, and inserted once reading elements is finished
add the boundary element nodes to the set of nodes
If the element yet another dimension, just read in the nodes and throw them away
If any lower dimensional elements have been found in the file, try to add them to the mesh.boundary_info as sides and nodes with the respecitve id's (called 'physical' in Gmsh).
Definition at line 319 of file gmsh_io.C.
References MeshBase::active_elements_begin(), MeshBase::active_elements_end(), MeshBase::add_elem(), MeshBase::add_point(), MeshBase::boundary_info, Elem::build(), Elem::build_side(), MeshBase::clear(), MeshInput< MeshBase >::mesh(), MeshBase::mesh_dimension(), Elem::n_nodes(), Elem::n_sides(), Elem::neighbor(), MeshBase::node_ptr(), libMesh::processor_id(), MeshBase::reserve_elem(), MeshBase::reserve_nodes(), DofObject::set_id(), Elem::set_node(), Elem::subdomain_id(), and Elem::type().
Referenced by read().
{
// This is a serial-only process for now;
// the Mesh should be read on processor 0 and
// broadcast later
libmesh_assert(libMesh::processor_id() == 0);
libmesh_assert(in.good());
// initialize the map with element types
init_eletypes();
// clear any data in the mesh
MeshBase& mesh = MeshInput<MeshBase>::mesh();
mesh.clear();
const unsigned int dim = mesh.mesh_dimension();
// some variables
const int bufLen = 256;
char buf[bufLen+1];
int format=0, size=0;
Real version = 1.0;
// map to hold the node numbers for translation
// note the the nodes can be non-consecutive
std::map<unsigned int, unsigned int> nodetrans;
{
while (!in.eof()) {
in >> buf;
if (!std::strncmp(buf,'$MeshFormat',11))
{
in >> version >> format >> size;
if(version != 2.0){
std::cerr << 'Error: Wrong msh file version ' << version << ';
libmesh_error();
}
if(format){
std::cerr << 'Error: Unknown data format for mesh;
libmesh_error();
}
}
// read the node block
else if (!std::strncmp(buf,'$NOD',4) ||
!std::strncmp(buf,'$NOE',4) ||
!std::strncmp(buf,'$Nodes',6)
)
{
unsigned int numNodes = 0;
in >> numNodes;
mesh.reserve_nodes (numNodes);
// read in the nodal coordinates and form points.
Real x, y, z;
unsigned int id;
// add the nodal coordinates to the mesh
for (unsigned int i=0; i<numNodes; ++i)
{
in >> id >> x >> y >> z;
mesh.add_point (Point(x, y, z), i);
nodetrans[id] = i;
}
// read the $ENDNOD delimiter
in >> buf;
}
else if (!std::strncmp(buf,'$ELM',4) ||
!std::strncmp(buf,'$Elements',9)
)
{
unsigned int numElem = 0;
std::vector< boundaryElementInfo > boundary_elem;
// read how many elements are there, and reserve space in the mesh
in >> numElem;
mesh.reserve_elem (numElem);
// read the elements
unsigned int elem_id_counter = 0;
for (unsigned int iel=0; iel<numElem; ++iel)
{
unsigned int id, type, physical, elementary,
/* partition = 1,*/ nnodes, ntags;
// note - partition was assigned but never used - BSK
if(version <= 1.0)
{
in >> id >> type >> physical >> elementary >> nnodes;
}
else
{
in >> id >> type >> ntags;
elementary = physical = /* partition = */ 1;
for(unsigned int j = 0; j < ntags; j++)
{
int tag;
in >> tag;
if(j == 0)
physical = tag;
else if(j == 1)
elementary = tag;
// else if(j == 2)
// partition = tag;
// ignore any other tags for now
}
}
// consult the import element table which element to build
const elementDefinition& eletype = eletypes_imp[type];
nnodes = eletype.nnodes;
// only elements that match the mesh dimension are added
// if the element dimension is one less than dim, the nodes and
// sides are added to the mesh.boundary_info
if (eletype.dim == dim)
{
// add the elements to the mesh
Elem* elem = Elem::build(eletype.type).release();
elem->set_id(elem_id_counter);
mesh.add_elem(elem);
// different to iel, lower dimensional elems aren't added
elem_id_counter++;
// check number of nodes. We cannot do that for version 2.0
if (version <= 1.0)
{
if (elem->n_nodes() != nnodes)
{
std::cerr << 'Number of nodes for element ' << id
<< ' of type ' << eletypes_imp[type].type
<< ' (Gmsh type ' << type
<< ') does not match Libmesh definition. '
<< 'I expected ' << elem->n_nodes()
<< ' nodes, but got ' << nnodes << ';
libmesh_error();
}
}
// add node pointers to the elements
int nod = 0;
// if there is a node translation table, use it
if (eletype.nodes.size() > 0)
for (unsigned int i=0; i<nnodes; i++)
{
in >> nod;
elem->set_node(eletype.nodes[i]) = mesh.node_ptr(nodetrans[nod]);
}
else
{
for (unsigned int i=0; i<nnodes; i++)
{
in >> nod;
elem->set_node(i) = mesh.node_ptr(nodetrans[nod]);
}
}
// Finally, set the subdomain ID to physical
elem->subdomain_id() = physical;
} // if element.dim == dim
// if this is a boundary
else if (eletype.dim == dim-1)
{
boundaryElementInfo binfo;
std::set<unsigned int>::iterator iter = binfo.nodes.begin();
int nod = 0;
for (unsigned int i=0; i<nnodes; i++)
{
in >> nod;
mesh.boundary_info->add_node(nodetrans[nod], physical);
binfo.nodes.insert(iter, nodetrans[nod]);
}
binfo.id = physical;
boundary_elem.push_back(binfo);
}
else
{
int nod = 0;
for (unsigned int i=0; i<nnodes; i++)
in >> nod;
}
}//element loop
// read the $ENDELM delimiter
in >> buf;
if (boundary_elem.size() > 0)
{
// create a index of the boundary nodes to easily locate which
// element might have that boundary
std::map<unsigned int, std::vector<unsigned int> > node_index;
for (unsigned int i=0; i<boundary_elem.size(); i++)
{
boundaryElementInfo binfo = boundary_elem[i];
std::set<unsigned int>::iterator iter = binfo.nodes.begin();
for (;iter!= binfo.nodes.end(); iter++)
node_index[*iter].push_back(i);
}
MeshBase::const_element_iterator it = mesh.active_elements_begin();
const MeshBase::const_element_iterator end = mesh.active_elements_end();
// iterate over all elements and see which boundary element has
// the same set of nodes as on of the boundary elements previously read
for ( ; it != end; ++it)
{
const Elem* elem = *it;
for (unsigned int s=0; s<elem->n_sides(); s++)
if (elem->neighbor(s) == NULL)
{
AutoPtr<Elem> side (elem->build_side(s));
std::set<unsigned int> side_nodes;
std::set<unsigned int>::iterator iter = side_nodes.begin();
// make a set with all nodes from this side
// this allows for easy comparison
for (unsigned int ns=0; ns<side->n_nodes(); ns++)
side_nodes.insert(iter, side->node(ns));
// See whether one of the side node occurs in the list
// of tagged nodes. If we would loop over all side
// nodes, we would just get multiple hits, so taking
// node 0 is enough to do the job
unsigned int sn = side->node(0);
if (node_index.count(sn) > 0)
{
// Loop over all tagged ('physical') 'sides' which
// contain the node sn (typically just 1 to
// three). For each of these the set of nodes is
// compared to the current element's side nodes
for (unsigned int n=0; n<node_index[sn].size(); n++)
{
unsigned int bidx = node_index[sn][n];
if (boundary_elem[bidx].nodes == side_nodes)
mesh.boundary_info->add_side(elem, s, boundary_elem[bidx].id);
}
}
} // if elem->neighbor(s) == NULL
} // element loop
} // if boundary_elem.size() > 0
} // if $ELM
} // while !in.eof()
}
}
Referenced by TetGenIO::read(), and UCDIO::read_implementation().
Implements MeshOutput< MeshBase >.
Definition at line 594 of file gmsh_io.C.
References libMesh::processor_id(), and write_mesh().
Referenced by UnstructuredMesh::write().
{
if (libMesh::processor_id() == 0)
{
// Open the output file stream
std::ofstream out (name.c_str());
// Make sure it opened correctly
if (!out.good())
libmesh_file_error(name.c_str());
this->write_mesh (out);
}
}
Definition at line 620 of file gmsh_io.C.
References MeshBase::active_elements_begin(), MeshBase::active_elements_end(), DofObject::id(), MeshInput< MeshBase >::mesh(), MeshBase::n_active_elem(), Elem::n_nodes(), MeshBase::n_nodes(), Elem::node(), MeshBase::node(), DofObject::processor_id(), and Elem::type().
Referenced by write().
{
// Be sure that the stream is valid.
libmesh_assert (out.good());
// initialize the map with element types
init_eletypes();
// Get a const reference to the mesh
const MeshBase& mesh = MeshOutput<MeshBase>::mesh();
// Note: we are using version 2.0 of the gmsh output format.
{
// Write the file header.
out << '$MeshFormat;
out << '2.0 0 ' << sizeof(Real) << ';
out << '$EndMeshFormat;
}
{
// write the nodes in (n x y z) format
out << '$Nodes;
out << mesh.n_nodes() << ';
for (unsigned int v=0; v<mesh.n_nodes(); v++)
out << mesh.node(v).id()+1 << ' '
<< mesh.node(v)(0) << ' '
<< mesh.node(v)(1) << ' '
<< mesh.node(v)(2) << ';
out << '$EndNodes;
}
{
// write the connectivity
out << '$Elements;
out << mesh.n_active_elem() << ';
MeshBase::const_element_iterator it = mesh.active_elements_begin();
const MeshBase::const_element_iterator end = mesh.active_elements_end();
// loop over the elements
for ( ; it != end; ++it)
{
const Elem* elem = *it;
// Make sure we have a valid entry for
// the current element type.
libmesh_assert (eletypes_exp.count(elem->type()));
// consult the export element table
const elementDefinition& eletype = eletypes_exp[elem->type()];
// The element mapper better not require any more nodes
// than are present in the current element!
libmesh_assert (eletype.nodes.size() <= elem->n_nodes());
// elements ids are 1 based in Gmsh
out << elem->id()+1 << ' ';
// element type
out << eletype.exptype;
// write the number of tags and
// tag1 (physical entity), and tag2 (geometric entity)
out << ' 3 1 1 ';
// write the partition the element belongs to
out << elem->processor_id()+1 << ' ';
// if there is a node translation table, use it
if (eletype.nodes.size() > 0)
for (unsigned int i=0; i < elem->n_nodes(); i++)
out << elem->node(eletype.nodes[i])+1 << ' '; // gmsh is 1-based
// otherwise keep the same node order
else
for (unsigned int i=0; i < elem->n_nodes(); i++)
out << elem->node(i)+1 << ' '; // gmsh is 1-based
out << ';
} // element loop
out << '$EndElements;
}
}
Reimplemented from MeshOutput< MeshBase >.
Definition at line 610 of file gmsh_io.C.
References libMesh::processor_id(), and write_post().
{
//this->_binary = true;
if (libMesh::processor_id() == 0)
this->write_post (fname, &soln, &names);
}
Definition at line 705 of file gmsh_io.C.
References MeshBase::active_elements_begin(), MeshBase::active_elements_end(), binary(), libMeshEnums::EDGE2, libMeshEnums::EDGE3, libMeshEnums::EDGE4, libMeshEnums::HEX20, libMeshEnums::HEX27, libMeshEnums::HEX8, libmesh_real(), MeshInput< MeshBase >::mesh(), MeshBase::n_nodes(), Elem::n_vertices(), Elem::node(), Elem::point(), libMeshEnums::PRISM15, libMeshEnums::PRISM18, libMeshEnums::PRISM6, libMesh::processor_id(), libMeshEnums::PYRAMID5, libMeshEnums::QUAD4, libMeshEnums::QUAD8, libMeshEnums::QUAD9, libMeshEnums::TET10, libMeshEnums::TET4, libMeshEnums::TRI3, and libMeshEnums::TRI6.
Referenced by write_nodal_data().
{
// Should only do this on processor 0!
libmesh_assert (libMesh::processor_id() == 0);
// Create an output stream
std::ofstream out(fname.c_str());
// Make sure it opened correctly
if (!out.good())
libmesh_file_error(fname.c_str());
// initialize the map with element types
init_eletypes();
// create a character buffer
char buf[80];
// Get a constant reference to the mesh.
const MeshBase& mesh = MeshOutput<MeshBase>::mesh();
// write the data
if ((solution_names != NULL) && (v != NULL))
{
const unsigned int n_vars = solution_names->size();
if (!(v->size() == mesh.n_nodes()*n_vars))
std::cerr << 'ERROR: v->size()=' << v->size()
<< ', mesh.n_nodes()=' << mesh.n_nodes()
<< ', n_vars=' << n_vars
<< ', mesh.n_nodes()*n_vars=' << mesh.n_nodes()*n_vars
<< ';
libmesh_assert (v->size() == mesh.n_nodes()*n_vars);
// write the header
out << '$PostFormat;
if (this->binary())
out << '1.2 1 ' << sizeof(double) << ';
else
out << '1.2 0 ' << sizeof(double) << ';
out << '$EndPostFormat;
// Loop over the elements to see how much of each type there are
unsigned int n_points=0, n_lines=0, n_triangles=0, n_quadrangles=0,
n_tetrahedra=0, n_hexahedra=0, n_prisms=0, n_pyramids=0;
unsigned int n_scalar=0, n_vector=0, n_tensor=0;
unsigned int nb_text2d=0, nb_text2d_chars=0, nb_text3d=0, nb_text3d_chars=0;
{
MeshBase::const_element_iterator it = mesh.active_elements_begin();
const MeshBase::const_element_iterator end = mesh.active_elements_end();
for ( ; it != end; ++it)
{
const ElemType elemtype = (*it)->type();
switch (elemtype)
{
case EDGE2:
case EDGE3:
case EDGE4:
{
n_lines += 1;
break;
}
case TRI3:
case TRI6:
{
n_triangles += 1;
break;
}
case QUAD4:
case QUAD8:
case QUAD9:
{
n_quadrangles += 1;
break;
}
case TET4:
case TET10:
{
n_tetrahedra += 1;
break;
}
case HEX8:
case HEX20:
case HEX27:
{
n_hexahedra += 1;
break;
}
case PRISM6:
case PRISM15:
case PRISM18:
{
n_prisms += 1;
break;
}
case PYRAMID5:
{
n_pyramids += 1;
break;
}
default:
{
std::cerr << 'ERROR: Not existant element type '
<< (*it)->type() << std::endl;
libmesh_error();
}
}
}
}
// create a view for each variable
for (unsigned int ivar=0; ivar < n_vars; ivar++)
{
std::string varname = (*solution_names)[ivar];
// at the moment, we just write out scalar quantities
// later this should be made configurable through
// options to the writer class
n_scalar = 1;
// write the variable as a view, and the number of time steps
out << '$View << varname << ' ' << 1 << ';
// write how many of each geometry type are written
out << n_points * n_scalar << ' '
<< n_points * n_vector << ' '
<< n_points * n_tensor << ' '
<< n_lines * n_scalar << ' '
<< n_lines * n_vector << ' '
<< n_lines * n_tensor << ' '
<< n_triangles * n_scalar << ' '
<< n_triangles * n_vector << ' '
<< n_triangles * n_tensor << ' '
<< n_quadrangles * n_scalar << ' '
<< n_quadrangles * n_vector << ' '
<< n_quadrangles * n_tensor << ' '
<< n_tetrahedra * n_scalar << ' '
<< n_tetrahedra * n_vector << ' '
<< n_tetrahedra * n_tensor << ' '
<< n_hexahedra * n_scalar << ' '
<< n_hexahedra * n_vector << ' '
<< n_hexahedra * n_tensor << ' '
<< n_prisms * n_scalar << ' '
<< n_prisms * n_vector << ' '
<< n_prisms * n_tensor << ' '
<< n_pyramids * n_scalar << ' '
<< n_pyramids * n_vector << ' '
<< n_pyramids * n_tensor << ' '
<< nb_text2d << ' '
<< nb_text2d_chars << ' '
<< nb_text3d << ' '
<< nb_text3d_chars << ';
// if binary, write a marker to identify the endianness of the file
if (this->binary())
{
const int one = 1;
std::memcpy(buf, &one, sizeof(int));
out.write(buf, sizeof(int));
}
// the time steps (there is just 1 at the moment)
if (this->binary())
{
double one = 1;
std::memcpy(buf, &one, sizeof(double));
out.write(buf, sizeof(double));
}
else
out << '1;
// Loop over the elements and write out the data
MeshBase::const_element_iterator it = mesh.active_elements_begin();
const MeshBase::const_element_iterator end = mesh.active_elements_end();
for ( ; it != end; ++it)
{
const Elem* elem = *it;
// this is quite crappy, but I did not invent that file format!
for (unsigned int d=0; d<3; d++) // loop over the dimensions
{
for (unsigned int n=0; n < elem->n_vertices(); n++) // loop over vertices
{
const Point vertex = elem->point(n);
if (this->binary())
{
double tmp = vertex(d);
std::memcpy(buf, &tmp, sizeof(double));
out.write(reinterpret_cast<char *>(buf), sizeof(double));
}
else
out << vertex(d) << ' ';
}
if (!this->binary())
out << ';
}
// now finally write out the data
for (unsigned int i=0; i < elem->n_vertices(); i++) // loop over vertices
if (this->binary())
{
#ifdef LIBMESH_USE_COMPLEX_NUMBERS
std::cout << 'WARNING: Gmsh::write_post does not fully support '
<< 'complex numbers. Will only write the real part of '
<< 'variable ' << varname << std::endl;
#endif
double tmp = libmesh_real((*v)[elem->node(i)*n_vars + ivar]);
std::memcpy(buf, &tmp, sizeof(double));
out.write(reinterpret_cast<char *>(buf), sizeof(double));
}
else
{
#ifdef LIBMESH_USE_COMPLEX_NUMBERS
std::cout << 'WARNING: Gmsh::write_post does not fully support '
<< 'complex numbers. Will only write the real part of '
<< 'variable ' << varname << std::endl;
#endif
out << libmesh_real((*v)[elem->node(i)*n_vars + ivar]) << ';
}
}
if (this->binary())
out << ';
out << '$EndView;
} // end variable loop (writing the views)
}
}
Definition at line 123 of file gmsh_io.h.
Referenced by binary().
Generated automatically by Doxygen for libMesh from the source code.