libMesh
Classes | Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Static Private Member Functions | Private Attributes | Static Private Attributes | List of all members
libMesh::VTKIO Class Reference

This class implements reading and writing meshes in the VTK format. More...

#include <vtk_io.h>

Inheritance diagram for libMesh::VTKIO:
[legend]

Classes

struct  ElementMaps
 Helper object that holds a map from VTK to libMesh element types and vice-versa. More...
 

Public Member Functions

 VTKIO (MeshBase &mesh)
 Constructor. More...
 
 VTKIO (const MeshBase &mesh)
 Constructor. More...
 
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) libmesh_override
 This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided. More...
 
virtual void read (const std::string &) libmesh_override
 This method implements reading a mesh from a specified file in VTK format. More...
 
virtual void write (const std::string &) libmesh_override
 Output the mesh without solutions to a .pvtu file. More...
 
void set_compression (bool b)
 Setter for compression flag. More...
 
vtkUnstructuredGrid * get_vtk_grid ()
 Get a pointer to the VTK unstructured grid data structure. More...
 
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
 This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object. More...
 
virtual void write_nodal_data (const std::string &, const NumericVector< Number > &, const std::vector< std::string > &)
 This method should be overridden by "parallel" output formats for writing nodal data. More...
 
unsigned intascii_precision ()
 Return/set the precision to use when writing ASCII files. More...
 

Protected Member Functions

MeshBasemesh ()
 
void set_n_partitions (unsigned int n_parts)
 Sets the number of partitions in the mesh. More...
 
void skip_comment_lines (std::istream &in, const char comment_start)
 Reads input from in, skipping all the lines that start with the character comment_start. More...
 
const MeshBasemesh () const
 

Protected Attributes

std::vector< bool > elems_of_dimension
 A vector of bools describing what dimension elements have been encountered when reading a mesh. More...
 
const bool _is_parallel_format
 Flag specifying whether this format is parallel-capable. More...
 
const bool _serial_only_needed_on_proc_0
 Flag specifying whether this format can be written by only serializing the mesh to processor zero. More...
 

Private Member Functions

void nodes_to_vtk ()
 write the nodes from the mesh into a vtkUnstructuredGrid More...
 
void cells_to_vtk ()
 write the cells from the mesh into a vtkUnstructuredGrid More...
 

Static Private Member Functions

static ElementMaps build_element_maps ()
 Static function used to construct the _element_maps struct. More...
 

Private Attributes

vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
 Write the system vectors to vtk. More...
 
bool _compress
 Flag to indicate whether the output should be compressed. More...
 
std::map< dof_id_type, dof_id_type_local_node_map
 maps global node id to node id of partition More...
 

Static Private Attributes

static ElementMaps _element_maps = VTKIO::build_element_maps()
 ElementMaps object that is built statically and used by all instances of this class. More...
 

Detailed Description

This class implements reading and writing meshes in the VTK format.

Format description: cf. VTK home page.

This class will not have any functionality unless VTK is detected during configure and hence LIBMESH_HAVE_VTK is defined.

Author
Wout Ruijter
John W. Peterson
Date
2007

Definition at line 59 of file vtk_io.h.

Constructor & Destructor Documentation

libMesh::VTKIO::VTKIO ( MeshBase mesh)
explicit

Constructor.

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

Definition at line 74 of file vtk_io.C.

References _compress, _element_maps, build_element_maps(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::Quality::name(), write(), and write_nodal_data().

74  :
75  MeshInput<MeshBase> (mesh, /*is_parallel_format=*/true),
76  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
77 #ifdef LIBMESH_HAVE_VTK
78  ,_compress(false)
79 #endif
80 {
81 }
82 
83 
84 
85 // Constructor for writing
86 VTKIO::VTKIO (const MeshBase & mesh) :
87  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
88 #ifdef LIBMESH_HAVE_VTK
89  ,_compress(false)
90 #endif
91 {
92 }
93 
94 
95 
96 // Output the mesh without solutions to a .pvtu file
97 void VTKIO::write (const std::string & name)
98 {
99  std::vector<Number> soln;
100  std::vector<std::string> names;
101  this->write_nodal_data(name, soln, names);
102 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
bool _compress
Flag to indicate whether the output should be compressed.
Definition: vtk_io.h:156
MeshOutput(const bool is_parallel_format=false, const bool serial_only_needed_on_proc_0=false)
Default constructor.
VTKIO(MeshBase &mesh)
Constructor.
Definition: vtk_io.C:74
libMesh::VTKIO::VTKIO ( const MeshBase mesh)
explicit

Constructor.

Takes a read-only reference to a mesh object. This is the constructor required to write a mesh.

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

VTKIO::ElementMaps libMesh::VTKIO::build_element_maps ( )
staticprivate

Static function used to construct the _element_maps struct.

Definition at line 114 of file vtk_io.C.

References libMesh::VTKIO::ElementMaps::associate(), libMesh::EDGE2, libMesh::EDGE3, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI3SUBDIVISION, libMesh::TRI6, and libMesh::VTKIO::ElementMaps::writing_map.

Referenced by VTKIO().

115 {
116  // Object to be filled up
117  ElementMaps em;
118 
119  em.associate(EDGE2, VTK_LINE);
120  em.associate(EDGE3, VTK_QUADRATIC_EDGE);
121  em.associate(TRI3, VTK_TRIANGLE);
122  em.associate(TRI6, VTK_QUADRATIC_TRIANGLE);
123  em.associate(QUAD4, VTK_QUAD);
124  em.associate(QUAD8, VTK_QUADRATIC_QUAD);
125  em.associate(TET4, VTK_TETRA);
126  em.associate(TET10, VTK_QUADRATIC_TETRA);
127  em.associate(HEX8, VTK_HEXAHEDRON);
128  em.associate(HEX20, VTK_QUADRATIC_HEXAHEDRON);
129  em.associate(HEX27, VTK_TRIQUADRATIC_HEXAHEDRON);
130  em.associate(PRISM6, VTK_WEDGE);
131  em.associate(PRISM15, VTK_QUADRATIC_WEDGE);
132  em.associate(PRISM18, VTK_BIQUADRATIC_QUADRATIC_WEDGE);
133  em.associate(PYRAMID5, VTK_PYRAMID);
134 
135  // VTK_BIQUADRATIC_QUAD has been around since VTK 5.0
136 #if VTK_MAJOR_VERSION > 5 || (VTK_MAJOR_VERSION == 5 && VTK_MINOR_VERSION > 0)
137  em.associate(QUAD9, VTK_BIQUADRATIC_QUAD);
138 #endif
139 
140  // TRI3SUBDIVISION is for writing only
141  em.writing_map[TRI3SUBDIVISION] = VTK_TRIANGLE;
142 
143  return em;
144 }
void libMesh::VTKIO::cells_to_vtk ( )
private

write the cells from the mesh into a vtkUnstructuredGrid

Definition at line 404 of file vtk_io.C.

References _element_maps, _local_node_map, _vtk_grid, libMesh::MeshBase::active_local_element_ptr_range(), libMesh::VTKIO::ElementMaps::find(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::n_active_local_elem(), libMesh::MeshBase::point(), read(), libMesh::VTK, and write_nodal_data().

Referenced by write_nodal_data().

405 {
406  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
407 
408  vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
409  vtkSmartPointer<vtkIdList> pts = vtkSmartPointer<vtkIdList>::New();
410 
411  std::vector<int> types(mesh.n_active_local_elem());
412 
413  vtkSmartPointer<vtkIntArray> elem_id = vtkSmartPointer<vtkIntArray>::New();
414  elem_id->SetName("libmesh_elem_id");
415  elem_id->SetNumberOfComponents(1);
416 
417  vtkSmartPointer<vtkIntArray> subdomain_id = vtkSmartPointer<vtkIntArray>::New();
418  subdomain_id->SetName("subdomain_id");
419  subdomain_id->SetNumberOfComponents(1);
420 
421  vtkSmartPointer<vtkIntArray> elem_proc_id = vtkSmartPointer<vtkIntArray>::New();
422  elem_proc_id->SetName("processor_id");
423  elem_proc_id->SetNumberOfComponents(1);
424 
425  unsigned active_element_counter = 0;
426  for (const auto & elem : mesh.active_local_element_ptr_range())
427  {
428  pts->SetNumberOfIds(elem->n_nodes());
429 
430  // get the connectivity for this element
431  std::vector<dof_id_type> conn;
432  elem->connectivity(0, VTK, conn);
433 
434  for (std::size_t i=0; i<conn.size(); ++i)
435  {
436  // If the node ID is not found in the _local_node_map, we'll
437  // add it to the _vtk_grid. NOTE[JWP]: none of the examples
438  // I have actually enters this section of code...
439  if (_local_node_map.find(conn[i]) == _local_node_map.end())
440  {
441  dof_id_type global_node_id = elem->node_id(i);
442 
443  const Point & the_node = mesh.point(global_node_id);
444 
445  // InsertNextPoint accepts either a double or float array of length 3.
446  double pt[3] = {0., 0., 0.};
447  for (unsigned int d=0; d<LIBMESH_DIM; ++d)
448  pt[d] = the_node(d);
449 
450  // Insert the point into the _vtk_grid
451  vtkIdType local = _vtk_grid->GetPoints()->InsertNextPoint(pt);
452 
453  // Update the _local_node_map with the ID returned by VTK
454  _local_node_map[global_node_id] =
455  cast_int<dof_id_type>(local);
456  }
457 
458  // Otherwise, the node ID was found in the _local_node_map, so
459  // insert it into the vtkIdList.
460  pts->InsertId(i, _local_node_map[conn[i]]);
461  }
462 
463  vtkIdType vtkcellid = cells->InsertNextCell(pts);
464  types[active_element_counter] = cast_int<int>(_element_maps.find(elem->type()));
465 
466  elem_id->InsertTuple1(vtkcellid, elem->id());
467  subdomain_id->InsertTuple1(vtkcellid, elem->subdomain_id());
468  elem_proc_id->InsertTuple1(vtkcellid, elem->processor_id());
469  ++active_element_counter;
470  } // end loop over active elements
471 
472  _vtk_grid->SetCells(&types[0], cells);
473  _vtk_grid->GetCellData()->AddArray(elem_id);
474  _vtk_grid->GetCellData()->AddArray(subdomain_id);
475  _vtk_grid->GetCellData()->AddArray(elem_proc_id);
476 }
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:161
vtkIdType find(ElemType libmesh_type)
Definition: vtk_io.h:178
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:151
const MT & mesh() const
Definition: mesh_output.h:216
static ElementMaps _element_maps
ElementMaps object that is built statically and used by all instances of this class.
Definition: vtk_io.h:207
uint8_t dof_id_type
Definition: id_types.h:64
vtkUnstructuredGrid * libMesh::VTKIO::get_vtk_grid ( )

Get a pointer to the VTK unstructured grid data structure.

Definition at line 347 of file vtk_io.C.

References _vtk_grid.

348 {
349  return _vtk_grid;
350 }
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:151
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::XdrIO::read_header(), libMesh::UCDIO::read_implementation(), libMesh::UNVIO::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(), 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().

const MeshBase & libMesh::MeshOutput< MeshBase >::mesh ( ) const
protectedinherited
void libMesh::VTKIO::nodes_to_vtk ( )
private

write the nodes from the mesh into a vtkUnstructuredGrid

Definition at line 361 of file vtk_io.C.

References _local_node_map, _vtk_grid, libMesh::DofObject::id(), libMesh::MeshBase::local_node_ptr_range(), libMesh::MeshOutput< MT >::mesh(), and libMesh::MeshBase::n_local_nodes().

Referenced by write_nodal_data().

362 {
363  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
364 
365  // containers for points and coordinates of points
366  vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
367  vtkSmartPointer<vtkDoubleArray> pcoords = vtkSmartPointer<vtkDoubleArray>::New();
368  // if this grid is to be used in VTK then the dimension of the points should be 3
369  pcoords->SetNumberOfComponents(LIBMESH_DIM);
370  pcoords->Allocate(3*mesh.n_local_nodes());
371  points->SetNumberOfPoints(mesh.n_local_nodes()); // it seems that it needs this to prevent a segfault
372 
373  unsigned int local_node_counter = 0;
374 
375  for (const auto & node_ptr : mesh.local_node_ptr_range())
376  {
377  const Node & node = *node_ptr;
378 
379  double pnt[3] = {0, 0, 0};
380  for (unsigned int i=0; i<LIBMESH_DIM; ++i)
381  pnt[i] = node(i);
382 
383  // Fill mapping between global and local node numbers
384  _local_node_map[node.id()] = local_node_counter;
385 
386  // add point
387 #if VTK_VERSION_LESS_THAN(7,1,0)
388  pcoords->InsertNextTupleValue(pnt);
389 #else
390  pcoords->InsertNextTuple(pnt);
391 #endif
392  ++local_node_counter;
393  }
394 
395  // add coordinates to points
396  points->SetData(pcoords);
397 
398  // add points to grid
399  _vtk_grid->SetPoints(points);
400 }
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:161
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:151
const MT & mesh() const
Definition: mesh_output.h:216
void libMesh::VTKIO::read ( const std::string &  name)
virtual

This method implements reading a mesh from a specified file in VTK format.

FIXME: This is known to write nonsense on AMR meshes and it strips the imaginary parts of complex Numbers.

As with other Mesh IO classes, this interface is still still "available" when !LIBMESH_HAVE_VTK, however, it will throw a runtime error.

This function is not currently used by anything, so it is commented out, and may eventually be removed entirely.

Implements libMesh::MeshInput< MeshBase >.

Definition at line 148 of file vtk_io.C.

References _element_maps, _vtk_grid, libMesh::MeshBase::add_elem(), libMesh::MeshBase::add_point(), libMesh::Elem::build(), libMesh::MeshBase::clear(), libMesh::Elem::connectivity(), libMesh::Elem::dim(), libMesh::MeshInput< MeshBase >::elems_of_dimension, libMesh::VTKIO::ElementMaps::find(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::Elem::n_nodes(), libMesh::MeshBase::node_ptr(), libMesh::processor_id(), libMesh::DofObject::set_id(), libMesh::MeshBase::set_mesh_dimension(), libMesh::Elem::set_node(), and libMesh::VTK.

Referenced by cells_to_vtk(), and libMesh::NameBasedIO::read().

149 {
150  // This is a serial-only process for now;
151  // the Mesh should be read on processor 0 and
152  // broadcast later
153  libmesh_assert_equal_to (MeshOutput<MeshBase>::mesh().processor_id(), 0);
154 
155  // Keep track of what kinds of elements this file contains
156  elems_of_dimension.clear();
157  elems_of_dimension.resize(4, false);
158 
159  // Use a typedef, because these names are just crazy
160  typedef vtkSmartPointer<vtkXMLUnstructuredGridReader> MyReader;
161  MyReader reader = MyReader::New();
162 
163  // Pass the filename along to the reader
164  reader->SetFileName(name.c_str());
165 
166  // Force reading
167  reader->Update();
168 
169  // read in the grid
170  _vtk_grid = reader->GetOutput();
171 
172  // Get a reference to the mesh
173  MeshBase & mesh = MeshInput<MeshBase>::mesh();
174 
175  // Clear out any pre-existing data from the Mesh
176  mesh.clear();
177 
178  // Get the number of points from the _vtk_grid object
179  const unsigned int vtk_num_points = static_cast<unsigned int>(_vtk_grid->GetNumberOfPoints());
180 
181  // always numbered nicely so we can loop like this
182  for (unsigned int i=0; i<vtk_num_points; ++i)
183  {
184  // add to the id map
185  // and add the actual point
186  double pnt[3];
187  _vtk_grid->GetPoint(static_cast<vtkIdType>(i), pnt);
188  Point xyz(pnt[0], pnt[1], pnt[2]);
189  mesh.add_point(xyz, i);
190  }
191 
192  // Get the number of cells from the _vtk_grid object
193  const unsigned int vtk_num_cells = static_cast<unsigned int>(_vtk_grid->GetNumberOfCells());
194 
195  vtkSmartPointer<vtkGenericCell> cell = vtkSmartPointer<vtkGenericCell>::New();
196  for (unsigned int i=0; i<vtk_num_cells; ++i)
197  {
198  _vtk_grid->GetCell(i, cell);
199 
200  // Get the libMesh element type corresponding to this VTK element type.
201  ElemType libmesh_elem_type = _element_maps.find(cell->GetCellType());
202  Elem * elem = Elem::build(libmesh_elem_type).release();
203 
204  // get the straightforward numbering from the VTK cells
205  for (unsigned int j=0; j<elem->n_nodes(); ++j)
206  elem->set_node(j) =
207  mesh.node_ptr(cast_int<dof_id_type>(cell->GetPointId(j)));
208 
209  // then get the connectivity
210  std::vector<dof_id_type> conn;
211  elem->connectivity(0, VTK, conn);
212 
213  // then reshuffle the nodes according to the connectivity, this
214  // two-time-assign would evade the definition of the vtk_mapping
215  for (std::size_t j=0; j<conn.size(); ++j)
216  elem->set_node(j) = mesh.node_ptr(conn[j]);
217 
218  elem->set_id(i);
219 
220  elems_of_dimension[elem->dim()] = true;
221 
222  mesh.add_elem(elem);
223  } // end loop over VTK cells
224 
225  // Set the mesh dimension to the largest encountered for an element
226  for (unsigned char i=0; i!=4; ++i)
227  if (elems_of_dimension[i])
228  mesh.set_mesh_dimension(i);
229 
230 #if LIBMESH_DIM < 3
231  if (mesh.mesh_dimension() > LIBMESH_DIM)
232  libmesh_error_msg("Cannot open dimension " \
233  << mesh.mesh_dimension() \
234  << " mesh file when configured without " \
235  << mesh.mesh_dimension() \
236  << "D support.");
237 #endif // LIBMESH_DIM < 3
238 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:238
std::vector< bool > elems_of_dimension
A vector of bools describing what dimension elements have been encountered when reading a mesh...
Definition: mesh_input.h:97
vtkIdType find(ElemType libmesh_type)
Definition: vtk_io.h:178
ElemType
Defines an enum for geometric element types.
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:151
const MeshBase & mesh() const
static ElementMaps _element_maps
ElementMaps object that is built statically and used by all instances of this class.
Definition: vtk_io.h:207
processor_id_type processor_id()
Definition: libmesh_base.h:96
void libMesh::VTKIO::set_compression ( bool  b)

Setter for compression flag.

Definition at line 354 of file vtk_io.C.

References _compress.

355 {
356  this->_compress = b;
357 }
bool _compress
Flag to indicate whether the output should be compressed.
Definition: vtk_io.h:156
void libMesh::MeshInput< MeshBase >::set_n_partitions ( unsigned int  n_parts)
protectedinherited

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

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

virtual void libMesh::VTKIO::write ( const std::string &  )
virtual

Output the mesh without solutions to a .pvtu file.

As with other Mesh IO classes, this interface is still still "available" when !LIBMESH_HAVE_VTK, however, it will throw a runtime error.

Implements libMesh::MeshOutput< MeshBase >.

Referenced by main(), VTKIO(), and libMesh::NameBasedIO::write().

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::VTKIO::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.

As with other Mesh IO classes, this interface is still still "available" when !LIBMESH_HAVE_VTK, however, it will throw a runtime error.

Reimplemented from libMesh::MeshOutput< MeshBase >.

Definition at line 242 of file vtk_io.C.

References _compress, _local_node_map, _vtk_grid, cells_to_vtk(), data, libMesh::err, libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::n_nodes(), libMesh::n_processors(), nodes_to_vtk(), and libMesh::processor_id().

Referenced by cells_to_vtk(), VTKIO(), and libMesh::NameBasedIO::write_nodal_data().

245 {
246  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
247 
248  // Warn that the .pvtu file extension should be used. Paraview
249  // recognizes this, and it works in both serial and parallel. Only
250  // warn about this once.
251  if (fname.substr(fname.rfind("."), fname.size()) != ".pvtu")
252  libmesh_do_once(libMesh::err << "The .pvtu extension should be used when writing VTK files in libMesh.");
253 
254  // If there are variable names being written, the solution vector
255  // should not be empty, it should have been broadcast to all
256  // processors by the MeshOutput base class, since VTK is a parallel
257  // format. Verify this before going further.
258  if (!names.empty() && soln.empty())
259  libmesh_error_msg("Empty soln vector in VTKIO::write_nodal_data().");
260 
261  // we only use Unstructured grids
262  _vtk_grid = vtkSmartPointer<vtkUnstructuredGrid>::New();
263  vtkSmartPointer<vtkXMLPUnstructuredGridWriter> writer = vtkSmartPointer<vtkXMLPUnstructuredGridWriter>::New();
264 
265  // add nodes to the grid and update _local_node_map
266  _local_node_map.clear();
267  this->nodes_to_vtk();
268 
269  // add cells to the grid
270  this->cells_to_vtk();
271 
272  // add nodal solutions to the grid, if solutions are given
273  if (names.size() > 0)
274  {
275  std::size_t num_vars = names.size();
276  dof_id_type num_nodes = mesh.n_nodes();
277 
278  for (std::size_t variable=0; variable<num_vars; ++variable)
279  {
280  vtkSmartPointer<vtkDoubleArray> data = vtkSmartPointer<vtkDoubleArray>::New();
281  data->SetName(names[variable].c_str());
282 
283  // number of local and ghost nodes
284  data->SetNumberOfValues(_local_node_map.size());
285 
286  // loop over all nodes and get the solution for the current
287  // variable, if the node is in the current partition
288  for (dof_id_type k=0; k<num_nodes; ++k)
289  {
290  std::map<dof_id_type, dof_id_type>::iterator local_node_it = _local_node_map.find(k);
291  if (local_node_it == _local_node_map.end())
292  continue; // not a local node
293 
294 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
295  libmesh_do_once (libMesh::err << "Only writing the real part for complex numbers!\n"
296  << "if you need this support contact " << LIBMESH_PACKAGE_BUGREPORT
297  << std::endl);
298  data->SetValue(local_node_it->second, soln[k*num_vars + variable].real());
299 #else
300  data->SetValue(local_node_it->second, soln[k*num_vars + variable]);
301 #endif
302  }
303  _vtk_grid->GetPointData()->AddArray(data);
304  }
305  }
306 
307  // Tell the writer how many partitions exist and on which processor
308  // we are currently
309  writer->SetNumberOfPieces(MeshOutput<MeshBase>::mesh().n_processors());
310  writer->SetStartPiece(MeshOutput<MeshBase>::mesh().processor_id());
311  writer->SetEndPiece(MeshOutput<MeshBase>::mesh().processor_id());
312 
313  // partitions overlap by one node
314  // FIXME: According to this document
315  // http://paraview.org/Wiki/images/5/51/SC07_tut107_ParaView_Handouts.pdf
316  // the ghosts are cells rather than nodes.
317  writer->SetGhostLevel(1);
318 
319  // VTK 6 replaces SetInput() with SetInputData(). See
320  // http://www.vtk.org/Wiki/VTK/VTK_6_Migration/Replacement_of_SetInput
321  // for the full explanation.
322 #if VTK_VERSION_LESS_THAN(6,0,0)
323  writer->SetInput(_vtk_grid);
324 #else
325  writer->SetInputData(_vtk_grid);
326 #endif
327 
328  writer->SetFileName(fname.c_str());
329  writer->SetDataModeToAscii();
330 
331  // compress the output, if desired (switches also to binary)
332  if (this->_compress)
333  {
334 #if !VTK_VERSION_LESS_THAN(5,6,0)
335  writer->SetCompressorTypeToZLib();
336 #else
337  libmesh_do_once(libMesh::err << "Compression not implemented with old VTK libs!" << std::endl;);
338 #endif
339  }
340 
341  writer->Write();
342 
343 }
OStreamProxy err
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:161
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:151
bool _compress
Flag to indicate whether the output should be compressed.
Definition: vtk_io.h:156
const MT & mesh() const
Definition: mesh_output.h:216
void nodes_to_vtk()
write the nodes from the mesh into a vtkUnstructuredGrid
Definition: vtk_io.C:361
void cells_to_vtk()
write the cells from the mesh into a vtkUnstructuredGrid
Definition: vtk_io.C:404
IterBase * data
Ideally this private member data should have protected access.
processor_id_type processor_id()
Definition: libmesh_base.h:96
uint8_t dof_id_type
Definition: id_types.h:64
processor_id_type n_processors()
Definition: libmesh_base.h:88
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

bool libMesh::VTKIO::_compress
private

Flag to indicate whether the output should be compressed.

Definition at line 156 of file vtk_io.h.

Referenced by set_compression(), VTKIO(), and write_nodal_data().

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

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

Definition at line 207 of file vtk_io.h.

Referenced by cells_to_vtk(), read(), and VTKIO().

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<dof_id_type, dof_id_type> libMesh::VTKIO::_local_node_map
private

maps global node id to node id of partition

Definition at line 161 of file vtk_io.h.

Referenced by cells_to_vtk(), nodes_to_vtk(), and write_nodal_data().

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.

vtkSmartPointer<vtkUnstructuredGrid> libMesh::VTKIO::_vtk_grid
private

Write the system vectors to vtk.

This function is not currently used by anything, so it is commented out, and may eventually be removed entirely. pointer to the VTK grid. the vtkSmartPointer will automatically initialize the value to null and keep track of reference counting.

Definition at line 151 of file vtk_io.h.

Referenced by cells_to_vtk(), get_vtk_grid(), nodes_to_vtk(), read(), and write_nodal_data().

std::vector<bool> libMesh::MeshInput< MeshBase >::elems_of_dimension
protectedinherited

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