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 > &) 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 &) override
 This method implements reading a mesh from a specified file in VTK format. More...
 
virtual void write (const std::string &) 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...
 
bool is_parallel_format () const
 Returns true iff this mesh file format and input class are parallelized, so that all processors can read their share of the data at once. More...
 
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=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_discontinuous_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=nullptr)
 This method implements writing a mesh with discontinuous 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 may be overridden by "parallel" output formats for writing nodal data. More...
 
virtual void write_nodal_data (const std::string &, const EquationSystems &, const std::set< std::string > *)
 This method should be overridden by "parallel" output formats for writing nodal data. More...
 
virtual void write_nodal_data_discontinuous (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
 This method implements writing a mesh with discontinuous data to a specified file where the nodal data and variables names are provided. 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
 
virtual bool get_add_sides ()
 

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 and update the local_node_map. More...
 
void cells_to_vtk ()
 write the cells from the mesh into a vtkUnstructuredGrid More...
 
void node_values_to_vtk (const std::string &name, const std::vector< Real > &local_values)
 write the nodal values of soln to a vtkUnstructuredGrid More...
 
void get_local_node_values (std::vector< Number > &local_values, std::size_t variable, const std::vector< Number > &soln, const std::vector< std::string > &names)
 Extract the values of soln that correspond to the nodes. More...
 

Static Private Member Functions

static std::map< ElemMappingType, ElementMapsbuild_element_maps ()
 Static function used to construct _element_maps. 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 std::map< ElemMappingType, ElementMaps_element_maps = VTKIO::build_element_maps()
 ElementMaps objects that are 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 60 of file vtk_io.h.

Constructor & Destructor Documentation

◆ VTKIO() [1/2]

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 82 of file vtk_io.C.

82  :
83  MeshInput<MeshBase> (mesh, /*is_parallel_format=*/true),
84  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
85 #ifdef LIBMESH_HAVE_VTK
86  ,_compress(false)
87 #endif
88 {
89 }
90 
91 
92 
93 // Constructor for writing
94 VTKIO::VTKIO (const MeshBase & mesh) :
95  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
96 #ifdef LIBMESH_HAVE_VTK
97  ,_compress(false)
98 #endif
99 {
100 }
101 
102 
103 
104 // Output the mesh without solutions to a .pvtu file
105 void VTKIO::write (const std::string & name)
106 {
107  std::vector<Number> soln;
108  std::vector<std::string> names;
109  this->write_nodal_data(name, soln, names);
110 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
bool _compress
Flag to indicate whether the output should be compressed.
Definition: vtk_io.h:174
MeshOutput(const bool is_parallel_format=false, const bool serial_only_needed_on_proc_0=false)
Default constructor.
Definition: mesh_output.h:219
template class LIBMESH_EXPORT MeshOutput< MeshBase >
Definition: mesh_output.C:180
VTKIO(MeshBase &mesh)
Constructor.
Definition: vtk_io.C:82

◆ VTKIO() [2/2]

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

◆ ascii_precision()

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

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

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

Definition at line 269 of file mesh_output.h.

Referenced by libMesh::UNVIO::nodes_out(), libMesh::FroIO::write(), libMesh::STLIO::write(), libMesh::MEDITIO::write_ascii(), libMesh::TecplotIO::write_ascii(), libMesh::GMVIO::write_ascii_new_impl(), and libMesh::GMVIO::write_ascii_old_impl().

270 {
271  return _ascii_precision;
272 }
unsigned int _ascii_precision
Precision to use when writing ASCII files.
Definition: mesh_output.h:207

◆ build_element_maps()

std::map< ElemMappingType, VTKIO::ElementMaps > libMesh::VTKIO::build_element_maps ( )
staticprivate

Static function used to construct _element_maps.

Definition at line 122 of file vtk_io.C.

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

123 {
124  // Object to be filled up
125  std::map<ElemMappingType, VTKIO::ElementMaps> all_maps;
126  ElementMaps em; // Lagrange element maps
127 
128  em.associate(EDGE2, VTK_LINE);
129  em.associate(EDGE3, VTK_QUADRATIC_EDGE);
130  em.associate(TRI3, VTK_TRIANGLE);
131  em.associate(TRI6, VTK_QUADRATIC_TRIANGLE);
132  em.associate(QUAD4, VTK_QUAD);
133  em.associate(QUAD8, VTK_QUADRATIC_QUAD);
134  em.associate(TET4, VTK_TETRA);
135  em.associate(TET10, VTK_QUADRATIC_TETRA);
136  em.associate(HEX8, VTK_HEXAHEDRON);
137  em.associate(HEX20, VTK_QUADRATIC_HEXAHEDRON);
138  em.associate(HEX27, VTK_TRIQUADRATIC_HEXAHEDRON);
139  em.associate(PRISM6, VTK_WEDGE);
140  em.associate(PRISM15, VTK_QUADRATIC_WEDGE);
141  em.associate(PRISM18, VTK_BIQUADRATIC_QUADRATIC_WEDGE);
142  em.associate(PYRAMID5, VTK_PYRAMID);
143 
144  // VTK_BIQUADRATIC_QUAD has been around since VTK 5.0
145 #if VTK_MAJOR_VERSION > 5 || (VTK_MAJOR_VERSION == 5 && VTK_MINOR_VERSION > 0)
146  em.associate(QUAD9, VTK_BIQUADRATIC_QUAD);
147 #endif
148 
149  // TRI3SUBDIVISION is for writing only
150  em.writing_map[TRI3SUBDIVISION] = VTK_TRIANGLE;
151  all_maps[ElemMappingType::LAGRANGE_MAP] = em;
152 
153 
154  // VTK_BEZIER_* types were introduced in VTK 9.0
155 #if VTK_VERSION_LESS_THAN(9,0,0)
156  // Revert back to previous behavior when using an older version of VTK
158 #else
159  ElementMaps bem; // Rational Bernstein element maps
160  bem.associate(EDGE2, VTK_LINE);
161  bem.associate(EDGE3, VTK_BEZIER_CURVE);
162  bem.associate(TRI3, VTK_TRIANGLE);
163  bem.associate(TRI6, VTK_BEZIER_TRIANGLE);
164  bem.associate(QUAD4, VTK_QUAD);
165  bem.associate(QUAD8, VTK_QUADRATIC_QUAD);
166  bem.associate(QUAD9, VTK_BEZIER_QUADRILATERAL);
167  bem.associate(TET4, VTK_TETRA);
168  bem.associate(TET10, VTK_QUADRATIC_TETRA);
169  bem.associate(HEX8, VTK_HEXAHEDRON);
170  bem.associate(HEX20, VTK_QUADRATIC_HEXAHEDRON);
171  bem.associate(HEX27, VTK_BEZIER_HEXAHEDRON);
172  bem.associate(PRISM6, VTK_WEDGE);
173  bem.associate(PRISM15, VTK_QUADRATIC_WEDGE);
174  bem.associate(PRISM18, VTK_BEZIER_WEDGE);
175  bem.associate(PYRAMID5, VTK_PYRAMID);
176  bem.writing_map[TRI3SUBDIVISION] = VTK_TRIANGLE;
178 #endif
179 
180  return all_maps;
181 }

◆ cells_to_vtk()

void libMesh::VTKIO::cells_to_vtk ( )
private

write the cells from the mesh into a vtkUnstructuredGrid

Definition at line 564 of file vtk_io.C.

References _element_maps, _local_node_map, _vtk_grid, libMesh::MeshBase::default_mapping_type(), libMesh::libmesh_assert(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::n_active_local_elem(), libMesh::NODEELEM, libMesh::MeshBase::point(), and libMesh::VTK.

Referenced by write_nodal_data().

565 {
566  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
567 
568  auto& element_map = libmesh_map_find(_element_maps, mesh.default_mapping_type());
569  vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
570  vtkSmartPointer<vtkIdList> pts = vtkSmartPointer<vtkIdList>::New();
571 
572  std::vector<int> types(mesh.n_active_local_elem());
573 
574  // We already created this but we need to add more if we have any
575  // ghost nodes
576  vtkAbstractArray * abstract_node_id =
577  _vtk_grid->GetPointData()->GetAbstractArray("libmesh_node_id");
578  vtkIntArray * node_id = vtkIntArray::SafeDownCast(abstract_node_id);
579  libmesh_assert(node_id);
580 
581  vtkSmartPointer<vtkIntArray> elem_id = vtkSmartPointer<vtkIntArray>::New();
582  elem_id->SetName("libmesh_elem_id");
583  elem_id->SetNumberOfComponents(1);
584 
585  vtkSmartPointer<vtkIntArray> subdomain_id = vtkSmartPointer<vtkIntArray>::New();
586  subdomain_id->SetName("subdomain_id");
587  subdomain_id->SetNumberOfComponents(1);
588 
589  vtkSmartPointer<vtkIntArray> elem_proc_id = vtkSmartPointer<vtkIntArray>::New();
590  elem_proc_id->SetName("processor_id");
591  elem_proc_id->SetNumberOfComponents(1);
592 
593  unsigned active_element_counter = 0;
594  for (const auto & elem : mesh.active_local_element_ptr_range())
595  {
596  // When using rational bernstein these hold the weights
597  if ( elem->type() == NODEELEM )
598  continue;
599 
600  pts->SetNumberOfIds(elem->n_nodes());
601 
602  // get the connectivity for this element
603  std::vector<dof_id_type> conn;
604  elem->connectivity(0, VTK, conn);
605 
606  for (unsigned int i=0,
607  n_conn = cast_int<unsigned int>(conn.size());
608  i != n_conn; ++i)
609  {
610  // If the node ID is not found in the _local_node_map, we'll
611  // add it to the _vtk_grid. NOTE[JWP]: none of the examples
612  // I have actually enters this section of code...
613  if (!_local_node_map.count(conn[i]))
614  {
615  dof_id_type global_node_id = elem->node_id(i);
616 
617  const Point & the_node = mesh.point(global_node_id);
618 
619  // InsertNextPoint accepts either a double or float array of length 3.
620  double pt[3] = {0., 0., 0.};
621  for (unsigned int d=0; d<LIBMESH_DIM; ++d)
622  pt[d] = double(the_node(d));
623 
624  // Insert the point into the _vtk_grid
625  vtkIdType local = _vtk_grid->GetPoints()->InsertNextPoint(pt);
626 
627  // Update the _local_node_map with the ID returned by VTK
628  _local_node_map[global_node_id] =
629  cast_int<dof_id_type>(local);
630 
631  node_id->InsertTuple1(local, global_node_id);
632  }
633 
634  // Otherwise, the node ID was found in the _local_node_map, so
635  // insert it into the vtkIdList.
636  pts->InsertId(i, _local_node_map[conn[i]]);
637  }
638 
639  vtkIdType vtkcellid = cells->InsertNextCell(pts);
640  types[active_element_counter] = cast_int<int>(element_map.find(elem->type()));
641 
642  elem_id->InsertTuple1(vtkcellid, elem->id());
643  subdomain_id->InsertTuple1(vtkcellid, elem->subdomain_id());
644  elem_proc_id->InsertTuple1(vtkcellid, elem->processor_id());
645  ++active_element_counter;
646  } // end loop over active elements
647 
648  _vtk_grid->SetCells(types.data(), cells);
649  _vtk_grid->GetCellData()->AddArray(elem_id);
650  _vtk_grid->GetCellData()->AddArray(subdomain_id);
651  _vtk_grid->GetCellData()->AddArray(elem_proc_id);
652 }
const MT & mesh() const
Definition: mesh_output.h:259
static std::map< ElemMappingType, ElementMaps > _element_maps
ElementMaps objects that are built statically and used by all instances of this class.
Definition: vtk_io.h:215
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:179
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:169
ElemMappingType default_mapping_type() const
Returns the default master space to physical space mapping basis functions to be used on newly added ...
Definition: mesh_base.h:812
dof_id_type n_active_local_elem() const
Definition: mesh_base.h:565
libmesh_assert(ctx)
virtual const Point & point(const dof_id_type i) const =0
uint8_t dof_id_type
Definition: id_types.h:67

◆ get_add_sides()

virtual bool libMesh::MeshOutput< MeshBase >::get_add_sides ( )
inlineprotectedvirtualinherited
Returns
Whether or not added sides are expected to be output, to plot SIDE_DISCONTINUOUS data. Subclasses should override this if they are capable of plotting such data.

Reimplemented in libMesh::ExodusII_IO.

Definition at line 176 of file mesh_output.h.

176 { return false; }

◆ get_local_node_values()

void libMesh::VTKIO::get_local_node_values ( std::vector< Number > &  local_values,
std::size_t  variable,
const std::vector< Number > &  soln,
const std::vector< std::string > &  names 
)
private

Extract the values of soln that correspond to the nodes.

This method overwrites all values in local_values

Definition at line 673 of file vtk_io.C.

References _local_node_map, libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), and libMesh::MeshBase::n_nodes().

Referenced by write_nodal_data().

677 {
678  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
679  std::size_t num_vars = names.size();
680  dof_id_type num_nodes = mesh.n_nodes();
681 
682  local_values.clear();
683  local_values.resize(_local_node_map.size(), 0.0);
684 
685  // loop over all nodes and get the solution for the current
686  // variable, if the node is in the current partition
687  for (dof_id_type k=0; k<num_nodes; ++k)
688  if (const auto local_node_it = _local_node_map.find(k);
689  local_node_it != _local_node_map.end())
690  local_values[local_node_it->second] = soln[k*num_vars + variable];
691 }
const MT & mesh() const
Definition: mesh_output.h:259
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:179
virtual dof_id_type n_nodes() const =0
uint8_t dof_id_type
Definition: id_types.h:67

◆ get_vtk_grid()

vtkUnstructuredGrid * libMesh::VTKIO::get_vtk_grid ( )

Get a pointer to the VTK unstructured grid data structure.

Definition at line 470 of file vtk_io.C.

References _vtk_grid.

471 {
472  return _vtk_grid;
473 }
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:169

◆ is_parallel_format()

bool libMesh::MeshInput< MeshBase >::is_parallel_format ( ) const
inlineinherited

Returns true iff this mesh file format and input class are parallelized, so that all processors can read their share of the data at once.

Definition at line 87 of file mesh_input.h.

References libMesh::MeshInput< MT >::_is_parallel_format.

87 { return this->_is_parallel_format; }
const bool _is_parallel_format
Flag specifying whether this format is parallel-capable.
Definition: mesh_input.h:130

◆ mesh() [1/2]

MeshBase & libMesh::MeshInput< MeshBase >::mesh ( )
inlineprotectedinherited
Returns
The object as a writable reference.

Definition at line 178 of file mesh_input.h.

Referenced by libMesh::GMVIO::_read_one_cell(), cells_to_vtk(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::Nemesis_IO::copy_elemental_solution(), libMesh::ExodusII_IO::copy_nodal_solution(), libMesh::TetGenIO::element_in(), libMesh::UNVIO::elements_in(), libMesh::UNVIO::elements_out(), get_local_node_values(), libMesh::ExodusII_IO::get_sideset_data_indices(), libMesh::UNVIO::groups_in(), libMesh::TetGenIO::node_in(), libMesh::UNVIO::nodes_in(), libMesh::UNVIO::nodes_out(), nodes_to_vtk(), libMesh::Nemesis_IO::prepare_to_write_nodal_data(), libMesh::GMVIO::read(), libMesh::STLIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::ExodusII_IO::read(), libMesh::CheckpointIO::read(), read(), libMesh::STLIO::read_ascii(), libMesh::CheckpointIO::read_bcs(), libMesh::STLIO::read_binary(), libMesh::CheckpointIO::read_connectivity(), libMesh::ExodusII_IO::read_header(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::UCDIO::read_implementation(), libMesh::UNVIO::read_implementation(), libMesh::GmshIO::read_mesh(), libMesh::DynaIO::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::ExodusII_IO::read_sideset_data(), libMesh::OFFIO::read_stream(), libMesh::MatlabIO::read_stream(), libMesh::CheckpointIO::read_subdomain_names(), libMesh::STLIO::write(), libMesh::TetGenIO::write(), libMesh::Nemesis_IO::write(), libMesh::XdrIO::write(), libMesh::CheckpointIO::write(), libMesh::ExodusII_IO::write(), libMesh::GMVIO::write_ascii_new_impl(), libMesh::GMVIO::write_ascii_old_impl(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::Nemesis_IO::write_element_data(), libMesh::ExodusII_IO::write_element_data(), libMesh::ExodusII_IO::write_elemsets(), libMesh::UCDIO::write_header(), libMesh::UCDIO::write_implementation(), libMesh::UCDIO::write_interior_elems(), libMesh::GmshIO::write_mesh(), libMesh::UCDIO::write_nodal_data(), write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_common(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::UCDIO::write_nodes(), 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(), libMesh::ExodusII_IO::write_sideset_data(), libMesh::UCDIO::write_soln(), and libMesh::CheckpointIO::write_subdomain_names().

179 {
180  libmesh_error_msg_if(_obj == nullptr, "ERROR: _obj should not be nullptr!");
181  return *_obj;
182 }
MeshBase * _obj
A pointer to a non-const object object.
Definition: mesh_input.h:123

◆ mesh() [2/2]

const MeshBase & libMesh::MeshOutput< MeshBase >::mesh ( ) const
inlineprotectedinherited

◆ node_values_to_vtk()

void libMesh::VTKIO::node_values_to_vtk ( const std::string &  name,
const std::vector< Real > &  local_values 
)
private

write the nodal values of soln to a vtkUnstructuredGrid

Definition at line 654 of file vtk_io.C.

References _local_node_map, _vtk_grid, libMesh::index_range(), and libMesh::Quality::name().

Referenced by write_nodal_data().

656 {
657  vtkSmartPointer<vtkDoubleArray> data = vtkSmartPointer<vtkDoubleArray>::New();
658  data->SetName(name.c_str());
659 
660  libmesh_assert_equal_to(_local_node_map.size(), local_values.size());
661 
662  // number of local and ghost nodes
663  data->SetNumberOfValues(_local_node_map.size());
664 
665  // copy values into vtk
666  for (auto i : index_range(local_values)) {
667  data->SetValue(i, double(local_values[i]));
668  }
669 
670  _vtk_grid->GetPointData()->AddArray(data);
671 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:179
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:169
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:117

◆ nodes_to_vtk()

void libMesh::VTKIO::nodes_to_vtk ( )
private

write the nodes from the mesh into a vtkUnstructuredGrid and update the local_node_map.

Definition at line 484 of file vtk_io.C.

References _local_node_map, _vtk_grid, libMesh::MeshBase::default_mapping_data(), libMesh::MeshBase::default_mapping_type(), libMesh::DofObject::get_extra_datum(), libMesh::DofObject::id(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::n_local_nodes(), libMesh::RATIONAL_BERNSTEIN_MAP, libMesh::Real, and libMesh::MeshTools::weight().

Referenced by write_nodal_data().

485 {
486  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
487 
488  // containers for points and coordinates of points
489  vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
490  vtkSmartPointer<vtkDoubleArray> pcoords = vtkSmartPointer<vtkDoubleArray>::New();
491  // if this grid is to be used in VTK then the dimension of the points should be 3
492  pcoords->SetNumberOfComponents(LIBMESH_DIM);
493  pcoords->Allocate(3*mesh.n_local_nodes());
494  points->SetNumberOfPoints(mesh.n_local_nodes()); // it seems that it needs this to prevent a segfault
495 
496  // SetRationalWeights() was introduced in VTK 9.0
497 #if !VTK_VERSION_LESS_THAN(9,0,0)
498  bool have_weights = false;
499  int weight_index = 0;
500  vtkSmartPointer<vtkDoubleArray> rational_weights;
501 
503  {
504  rational_weights = vtkSmartPointer<vtkDoubleArray>::New();
505  rational_weights->SetName("RationalWeights");
506  rational_weights->SetNumberOfComponents(1);
507  weight_index = static_cast<int>(mesh.default_mapping_data());
508  have_weights = true;
509  }
510 #endif
511 
512  vtkSmartPointer<vtkIntArray> node_id = vtkSmartPointer<vtkIntArray>::New();
513  node_id->SetName("libmesh_node_id");
514  node_id->SetNumberOfComponents(1);
515 
516  unsigned int local_node_counter = 0;
517 
518  for (const auto & node_ptr : mesh.local_node_ptr_range())
519  {
520  const Node & node = *node_ptr;
521 
522  double pnt[3] = {0, 0, 0};
523  for (unsigned int i=0; i<LIBMESH_DIM; ++i)
524  pnt[i] = double(node(i));
525 
526  // Fill mapping between global and local node numbers
527  _local_node_map[node.id()] = local_node_counter;
528 
529  // add point
530 #if VTK_VERSION_LESS_THAN(7,1,0)
531  pcoords->InsertNextTupleValue(pnt);
532 #else
533  pcoords->InsertNextTuple(pnt);
534 #endif
535 #if !VTK_VERSION_LESS_THAN(9,0,0)
536  if (have_weights)
537  {
538  Real weight = node.get_extra_datum<Real>(weight_index);
539  rational_weights->InsertTuple1(local_node_counter, double(weight));
540  }
541 #endif
542 
543  node_id->InsertTuple1(local_node_counter, node.id());
544 
545  ++local_node_counter;
546  }
547 
548  // add coordinates to points
549  points->SetData(pcoords);
550 
551  // add points to grid
552  _vtk_grid->SetPoints(points);
553  _vtk_grid->GetPointData()->AddArray(node_id);
554 
555 #if !VTK_VERSION_LESS_THAN(9,0,0)
556  if (have_weights)
557  _vtk_grid->GetPointData()->SetRationalWeights(rational_weights);
558 
559 #endif
560 }
const MT & mesh() const
Definition: mesh_output.h:259
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:179
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:169
dof_id_type n_local_nodes() const
Definition: mesh_base.h:442
ElemMappingType default_mapping_type() const
Returns the default master space to physical space mapping basis functions to be used on newly added ...
Definition: mesh_base.h:812
dof_id_type weight(const MeshBase &mesh, const processor_id_type pid)
Definition: mesh_tools.C:437
unsigned char default_mapping_data() const
Returns any default data value used by the master space to physical space mapping.
Definition: mesh_base.h:830
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ read()

void libMesh::VTKIO::read ( const std::string &  name)
overridevirtual

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 185 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::MeshBase::default_mapping_type(), libMesh::MeshInput< MeshBase >::elems_of_dimension, libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::Quality::name(), libMesh::MeshBase::node_ptr(), libMesh::MeshBase::query_elem_ptr(), libMesh::MeshBase::query_node_ptr(), libMesh::MeshBase::set_mesh_dimension(), and libMesh::VTK.

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

186 {
187  // This is a serial-only process for now;
188  // the Mesh should be read on processor 0 and
189  // broadcast later
190  libmesh_assert_equal_to (MeshOutput<MeshBase>::mesh().processor_id(), 0);
191 
192  // Keep track of what kinds of elements this file contains
193  elems_of_dimension.clear();
194  elems_of_dimension.resize(4, false);
195 
196  // Use a typedef, because these names are just crazy
197  typedef vtkSmartPointer<vtkXMLPUnstructuredGridReader> MyReader;
198  MyReader reader = MyReader::New();
199 
200  // Pass the filename along to the reader
201  reader->SetFileName(name.c_str());
202 
203  // Force reading
204  reader->Update();
205 
206  // read in the grid
207  _vtk_grid = reader->GetOutput();
208 
209  // Get a reference to the mesh
210  MeshBase & mesh = MeshInput<MeshBase>::mesh();
211 
212  // Clear out any pre-existing data from the Mesh
213  mesh.clear();
214 
215  // Try to preserve any libMesh ids and subdomain ids we find in the
216  // file. This will be null if there are none, e.g. if a non-libMesh
217  // program wrote this file.
218 
219  vtkAbstractArray * abstract_elem_id =
220  _vtk_grid->GetCellData()->GetAbstractArray("libmesh_elem_id");
221  vtkAbstractArray * abstract_node_id =
222  _vtk_grid->GetPointData()->GetAbstractArray("libmesh_node_id");
223  vtkAbstractArray * abstract_subdomain_id =
224  _vtk_grid->GetCellData()->GetAbstractArray("subdomain_id");
225 
226  // Get ids as integers. This will be null if they are another data
227  // type, e.g. if a non-libMesh program used the names we thought
228  // were unique for different data.
229  vtkIntArray * elem_id = vtkIntArray::SafeDownCast(abstract_elem_id);
230  vtkIntArray * node_id = vtkIntArray::SafeDownCast(abstract_node_id);
231  vtkIntArray * subdomain_id = vtkIntArray::SafeDownCast(abstract_subdomain_id);
232 
233  if (abstract_elem_id && !elem_id)
234  libmesh_warning("Found non-integral libmesh_elem_id array; forced to ignore it.\n"
235  "This is technically valid but probably broken.");
236 
237  if (abstract_node_id && !node_id)
238  libmesh_warning("Found non-integral libmesh_node_id array; forced to ignore it.\n"
239  "This is technically valid but probably broken.");
240 
241  if (abstract_subdomain_id && !subdomain_id)
242  libmesh_warning("Found non-integral subdomain_id array; forced to ignore it.\n"
243  "This is technically valid but probably broken.");
244 
245  // Get the number of points from the _vtk_grid object
246  const unsigned int vtk_num_points = static_cast<unsigned int>(_vtk_grid->GetNumberOfPoints());
247 
248  // Map from VTK indexing to libMesh id if necessary
249  std::vector<dof_id_type> vtk_node_to_libmesh;
250  if (node_id)
251  vtk_node_to_libmesh.resize(vtk_num_points);
252 
253  // always numbered nicely so we can loop like this
254  for (unsigned int i=0; i<vtk_num_points; ++i)
255  {
256  // add to the id map
257  // and add the actual point
258  double pnt[3];
259  _vtk_grid->GetPoint(static_cast<vtkIdType>(i), pnt);
260  Point xyz(pnt[0], pnt[1], pnt[2]);
261 
262  if (node_id)
263  {
264  auto id = node_id->GetValue(i);
265 
266  // It would nice to distinguish between "duplicate nodes
267  // because one was ghosted in a parallel file segment" and
268  // "duplicate nodes because there was a bug", but I'm not
269  // sure how to do that with vtkXMLPUnstructuredGridReader
270  if (!mesh.query_node_ptr(id))
271  mesh.add_point(xyz, id);
272  vtk_node_to_libmesh[i] = id;
273  }
274  else
275  mesh.add_point(xyz, i);
276  }
277 
278  // Get the number of cells from the _vtk_grid object
279  const unsigned int vtk_num_cells = static_cast<unsigned int>(_vtk_grid->GetNumberOfCells());
280 
281  auto& element_map = libmesh_map_find(_element_maps, mesh.default_mapping_type());
282 
283  vtkSmartPointer<vtkGenericCell> cell = vtkSmartPointer<vtkGenericCell>::New();
284  for (unsigned int i=0; i<vtk_num_cells; ++i)
285  {
286  _vtk_grid->GetCell(i, cell);
287 
288  // Get the libMesh element type corresponding to this VTK element type.
289  ElemType libmesh_elem_type = element_map.find(cell->GetCellType());
290  auto elem = Elem::build(libmesh_elem_type);
291 
292  // get the straightforward numbering from the VTK cells
293  for (auto j : elem->node_index_range())
294  {
295  const auto vtk_point_id = cell->GetPointId(j);
296  const dof_id_type libmesh_node_id = node_id ?
297  vtk_node_to_libmesh[vtk_point_id] : vtk_point_id;
298 
299  elem->set_node(j, mesh.node_ptr(libmesh_node_id));
300  }
301 
302  // then get the connectivity
303  std::vector<dof_id_type> conn;
304  elem->connectivity(0, VTK, conn);
305 
306  // then reshuffle the nodes according to the connectivity, this
307  // two-time-assign would evade the definition of the vtk_mapping
308  for (unsigned int j=0,
309  n_conn = cast_int<unsigned int>(conn.size());
310  j != n_conn; ++j)
311  elem->set_node(j, mesh.node_ptr(conn[j]));
312 
313  if (elem_id)
314  {
315  auto id = elem_id->GetValue(i);
316  libmesh_error_msg_if
317  (mesh.query_elem_ptr(id), "Duplicate element id " << id <<
318  " found in libmesh_elem_ids");
319  elem->set_id(id);
320  }
321  else
322  elem->set_id(i);
323 
324  if (subdomain_id)
325  {
326  auto sbdid = subdomain_id->GetValue(i);
327  elem->subdomain_id() = sbdid;
328  }
329 
330  elems_of_dimension[elem->dim()] = true;
331 
332  mesh.add_elem(std::move(elem));
333  } // end loop over VTK cells
334 
335  // Set the mesh dimension to the largest encountered for an element
336  for (unsigned char i=0; i!=4; ++i)
337  if (elems_of_dimension[i])
339 
340 #if LIBMESH_DIM < 3
341  libmesh_error_msg_if(mesh.mesh_dimension() > LIBMESH_DIM,
342  "Cannot open dimension "
343  << mesh.mesh_dimension()
344  << " mesh file when configured without "
345  << mesh.mesh_dimension()
346  << "D support.");
347 #endif // LIBMESH_DIM < 3
348 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
const MeshBase & mesh() const
Definition: mesh_output.h:259
static std::map< ElemMappingType, ElementMaps > _element_maps
ElementMaps objects that are built statically and used by all instances of this class.
Definition: vtk_io.h:215
ElemType
Defines an enum for geometric element types.
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:107
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:169
virtual Node * add_point(const Point &p, const dof_id_type id=DofObject::invalid_id, const processor_id_type proc_id=DofObject::invalid_processor_id)=0
Add a new Node at Point p to the end of the vertex array, with processor_id procid.
ElemMappingType default_mapping_type() const
Returns the default master space to physical space mapping basis functions to be used on newly added ...
Definition: mesh_base.h:812
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
static std::unique_ptr< Elem > build(const ElemType type, Elem *p=nullptr)
Definition: elem.C:444
virtual const Node * query_node_ptr(const dof_id_type i) const =0
void set_mesh_dimension(unsigned char d)
Resets the logical dimension of the mesh.
Definition: mesh_base.h:275
virtual void clear()
Deletes all the element and node data that is currently stored.
Definition: mesh_base.C:920
virtual const Elem * query_elem_ptr(const dof_id_type i) const =0
unsigned int mesh_dimension() const
Definition: mesh_base.C:372
virtual const Node * node_ptr(const dof_id_type i) const =0
uint8_t dof_id_type
Definition: id_types.h:67

◆ set_compression()

void libMesh::VTKIO::set_compression ( bool  b)

Setter for compression flag.

Definition at line 477 of file vtk_io.C.

References _compress.

478 {
479  this->_compress = b;
480 }
bool _compress
Flag to indicate whether the output should be compressed.
Definition: vtk_io.h:174

◆ set_n_partitions()

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 101 of file mesh_input.h.

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

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

101 { this->mesh().set_n_partitions() = n_parts; }
unsigned int & set_n_partitions()
Definition: mesh_base.h:1865

◆ skip_comment_lines()

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.

Definition at line 187 of file mesh_input.h.

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

189 {
190  char c, line[256];
191 
192  while (in.get(c), c==comment_start)
193  in.getline (line, 255);
194 
195  // put back first character of
196  // first non-comment line
197  in.putback (c);
198 }

◆ write()

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

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(), MeshInputTest::testExodusFileMappings(), MeshInputTest::testVTKPreserveElemIds(), MeshInputTest::testVTKPreserveSubdomainIds(), and libMesh::NameBasedIO::write().

◆ write_discontinuous_equation_systems()

void libMesh::MeshOutput< MeshBase >::write_discontinuous_equation_systems ( const std::string &  fname,
const EquationSystems es,
const std::set< std::string > *  system_names = nullptr 
)
virtualinherited

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

Definition at line 89 of file mesh_output.C.

References libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_variable_names(), libMesh::EquationSystems::get_mesh(), libMesh::libmesh_assert(), and libMesh::out.

Referenced by libMesh::ExodusII_IO::write_timestep_discontinuous().

92 {
93  LOG_SCOPE("write_discontinuous_equation_systems()", "MeshOutput");
94 
95  // We may need to gather and/or renumber a DistributedMesh to output
96  // it, making that const qualifier in our constructor a dirty lie
97  MT & my_mesh = const_cast<MT &>(*_obj);
98 
99  // If we're asked to write data that's associated with a different
100  // mesh, output files full of garbage are the result.
101  libmesh_assert_equal_to(&es.get_mesh(), _obj);
102 
103  // A non-renumbered mesh may not have a contiguous numbering, and
104  // that needs to be fixed before we can build a solution vector.
105  if (my_mesh.max_elem_id() != my_mesh.n_elem() ||
106  my_mesh.max_node_id() != my_mesh.n_nodes())
107  {
108  // If we were allowed to renumber then we should have already
109  // been properly renumbered...
110  libmesh_assert(!my_mesh.allow_renumbering());
111 
112  libmesh_do_once(libMesh::out <<
113  "Warning: This MeshOutput subclass only supports meshes which are contiguously renumbered!"
114  << std::endl;);
115 
116  my_mesh.allow_renumbering(true);
117 
118  my_mesh.renumber_nodes_and_elements();
119 
120  // Not sure what good going back to false will do here, the
121  // renumbering horses have already left the barn...
122  my_mesh.allow_renumbering(false);
123  }
124 
125  MeshSerializer serialize(const_cast<MT &>(*_obj), !_is_parallel_format, _serial_only_needed_on_proc_0);
126 
127  // Build the list of variable names that will be written.
128  std::vector<std::string> names;
129  es.build_variable_names (names, nullptr, system_names);
130 
131  if (!_is_parallel_format)
132  {
133  // Build the nodal solution values & get the variable
134  // names from the EquationSystems object
135  std::vector<Number> soln;
136  es.build_discontinuous_solution_vector (soln, system_names,
137  nullptr, false, /* defaults */
138  this->get_add_sides());
139 
140  this->write_nodal_data_discontinuous (fname, soln, names);
141  }
142  else // _is_parallel_format
143  {
144  libmesh_not_implemented();
145  }
146 }
virtual void write_nodal_data_discontinuous(const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
This method implements writing a mesh with discontinuous data to a specified file where the nodal dat...
Definition: mesh_output.h:118
const MeshBase *const _obj
A pointer to a constant object.
Definition: mesh_output.h:202
const bool _is_parallel_format
Flag specifying whether this format is parallel-capable.
Definition: mesh_output.h:184
libmesh_assert(ctx)
OStreamProxy out
const bool _serial_only_needed_on_proc_0
Flag specifying whether this format can be written by only serializing the mesh to processor zero...
Definition: mesh_output.h:193

◆ write_equation_systems()

void libMesh::MeshOutput< MeshBase >::write_equation_systems ( const std::string &  fname,
const EquationSystems es,
const std::set< std::string > *  system_names = 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::ExodusII_IO, and libMesh::NameBasedIO.

Definition at line 31 of file mesh_output.C.

References libMesh::EquationSystems::build_solution_vector(), libMesh::EquationSystems::build_variable_names(), libMesh::EquationSystems::get_mesh(), libMesh::libmesh_assert(), and libMesh::out.

Referenced by libMesh::Nemesis_IO::write_timestep().

34 {
35  LOG_SCOPE("write_equation_systems()", "MeshOutput");
36 
37  // We may need to gather and/or renumber a DistributedMesh to output
38  // it, making that const qualifier in our constructor a dirty lie
39  MT & my_mesh = const_cast<MT &>(*_obj);
40 
41  // If we're asked to write data that's associated with a different
42  // mesh, output files full of garbage are the result.
43  libmesh_assert_equal_to(&es.get_mesh(), _obj);
44 
45  // A non-parallel format, non-renumbered mesh may not have a contiguous
46  // numbering, and that needs to be fixed before we can build a solution vector.
47  if (!_is_parallel_format &&
48  (my_mesh.max_elem_id() != my_mesh.n_elem() ||
49  my_mesh.max_node_id() != my_mesh.n_nodes()))
50  {
51  // If we were allowed to renumber then we should have already
52  // been properly renumbered...
53  libmesh_assert(!my_mesh.allow_renumbering());
54 
55  libmesh_do_once(libMesh::out <<
56  "Warning: This MeshOutput subclass only supports meshes which are contiguously renumbered!"
57  << std::endl;);
58 
59  my_mesh.allow_renumbering(true);
60 
61  my_mesh.renumber_nodes_and_elements();
62 
63  // Not sure what good going back to false will do here, the
64  // renumbering horses have already left the barn...
65  my_mesh.allow_renumbering(false);
66  }
67 
69  {
70  MeshSerializer serialize(const_cast<MT &>(*_obj), !_is_parallel_format, _serial_only_needed_on_proc_0);
71 
72  // Build the list of variable names that will be written.
73  std::vector<std::string> names;
74  es.build_variable_names (names, nullptr, system_names);
75 
76  // Build the nodal solution values & get the variable
77  // names from the EquationSystems object
78  std::vector<Number> soln;
79  es.build_solution_vector (soln, system_names,
80  this->get_add_sides());
81 
82  this->write_nodal_data (fname, soln, names);
83  }
84  else // _is_parallel_format
85  this->write_nodal_data (fname, es, system_names);
86 }
virtual void write_nodal_data(const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: mesh_output.h:109
const MeshBase *const _obj
A pointer to a constant object.
Definition: mesh_output.h:202
const bool _is_parallel_format
Flag specifying whether this format is parallel-capable.
Definition: mesh_output.h:184
libmesh_assert(ctx)
OStreamProxy out
const bool _serial_only_needed_on_proc_0
Flag specifying whether this format can be written by only serializing the mesh to processor zero...
Definition: mesh_output.h:193

◆ write_nodal_data() [1/3]

void libMesh::VTKIO::write_nodal_data ( const std::string &  fname,
const std::vector< Number > &  soln,
const std::vector< std::string > &  names 
)
overridevirtual

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 352 of file vtk_io.C.

References _compress, _local_node_map, _vtk_grid, cells_to_vtk(), libMesh::ParallelObject::comm(), libMesh::err, TIMPI::Communicator::get(), get_local_node_values(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::ParallelObject::n_processors(), node_values_to_vtk(), nodes_to_vtk(), and libMesh::ParallelObject::processor_id().

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

355 {
356  // Warn that the .pvtu file extension should be used. Paraview
357  // recognizes this, and it works in both serial and parallel. Only
358  // warn about this once.
359  if (fname.substr(fname.rfind("."), fname.size()) != ".pvtu")
360  libmesh_do_once(libMesh::err << "The .pvtu extension should be used when writing VTK files in libMesh.");
361 
362  // If there are variable names being written, the solution vector
363  // should not be empty, it should have been broadcast to all
364  // processors by the MeshOutput base class, since VTK is a parallel
365  // format. Verify this before going further.
366  libmesh_error_msg_if(!names.empty() && soln.empty(),
367  "Empty soln vector in VTKIO::write_nodal_data().");
368 
369  // Get a reference to the mesh
370  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
371 
372  // we only use Unstructured grids
373  _vtk_grid = vtkSmartPointer<vtkUnstructuredGrid>::New();
374  vtkSmartPointer<vtkXMLPUnstructuredGridWriter> writer = vtkSmartPointer<vtkXMLPUnstructuredGridWriter>::New();
375 #ifdef LIBMESH_HAVE_MPI
376  // Set VTK to the same communicator as libMesh
377  vtkSmartPointer<vtkMPICommunicator> vtk_comm = vtkSmartPointer<vtkMPICommunicator>::New();
378  MPI_Comm mpi_comm = mesh.comm().get();
379  vtkMPICommunicatorOpaqueComm vtk_opaque_comm(&mpi_comm);
380  vtk_comm->InitializeExternal(&vtk_opaque_comm);
381 
382  vtkSmartPointer<vtkMPIController> vtk_mpi_ctrl = vtkSmartPointer<vtkMPIController>::New();
383  vtk_mpi_ctrl->SetCommunicator(vtk_comm);
384 
385  writer->SetController(vtk_mpi_ctrl);
386 #endif
387 
388  // add nodes to the grid and update _local_node_map
389  _local_node_map.clear();
390  this->nodes_to_vtk();
391 
392  // add cells to the grid
393  this->cells_to_vtk();
394 
395  // add nodal solutions to the grid, if solutions are given
396  if (names.size() > 0)
397  {
398  std::size_t num_vars = names.size();
399  std::vector<Number> local_values;
400 
401 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
402  std::vector<Real> local_real_values;
403 #endif
404 
405  for (std::size_t variable=0; variable<num_vars; ++variable)
406  {
407  get_local_node_values(local_values, variable, soln, names);
408 
409 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
410  // write real part
411  local_real_values.resize(local_values.size());
412  std::transform(local_values.begin(), local_values.end(),
413  local_real_values.begin(),
414  [](Number x) { return x.real(); });
415  node_values_to_vtk(names[variable] + "_real", local_real_values);
416 
417  // write imaginary part
418  local_real_values.resize(local_values.size());
419  std::transform(local_values.begin(), local_values.end(),
420  local_real_values.begin(),
421  [](Number x) { return x.imag(); });
422  node_values_to_vtk(names[variable] + "_imag", local_real_values);
423 #else
424  node_values_to_vtk(names[variable], local_values);
425 #endif
426  }
427  }
428 
429  // Tell the writer how many partitions exist and on which processor
430  // we are currently
431  writer->SetNumberOfPieces(mesh.n_processors());
432  writer->SetStartPiece(mesh.processor_id());
433  writer->SetEndPiece(mesh.processor_id());
434 
435  // partitions overlap by one node
436  // FIXME: According to this document
437  // http://paraview.org/Wiki/images/5/51/SC07_tut107_ParaView_Handouts.pdf
438  // the ghosts are cells rather than nodes.
439  writer->SetGhostLevel(1);
440 
441  // VTK 6 replaces SetInput() with SetInputData(). See
442  // http://www.vtk.org/Wiki/VTK/VTK_6_Migration/Replacement_of_SetInput
443  // for the full explanation.
444 #if VTK_VERSION_LESS_THAN(6,0,0)
445  writer->SetInput(_vtk_grid);
446 #else
447  writer->SetInputData(_vtk_grid);
448 #endif
449 
450  writer->SetFileName(fname.c_str());
451  writer->SetDataModeToAscii();
452 
453  // compress the output, if desired (switches also to binary)
454  if (this->_compress)
455  {
456 #if !VTK_VERSION_LESS_THAN(5,6,0)
457  writer->SetCompressorTypeToZLib();
458  writer->SetDataModeToBinary();
459 #else
460  libmesh_do_once(libMesh::err << "Compression not implemented with old VTK libs!" << std::endl;);
461 #endif
462  }
463 
464  writer->Write();
465 
466 }
OStreamProxy err
const MT & mesh() const
Definition: mesh_output.h:259
void node_values_to_vtk(const std::string &name, const std::vector< Real > &local_values)
write the nodal values of soln to a vtkUnstructuredGrid
Definition: vtk_io.C:654
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:179
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:169
bool _compress
Flag to indicate whether the output should be compressed.
Definition: vtk_io.h:174
const Parallel::Communicator & comm() const
void nodes_to_vtk()
write the nodes from the mesh into a vtkUnstructuredGrid and update the local_node_map.
Definition: vtk_io.C:484
processor_id_type n_processors() const
communicator & get()
void get_local_node_values(std::vector< Number > &local_values, std::size_t variable, const std::vector< Number > &soln, const std::vector< std::string > &names)
Extract the values of soln that correspond to the nodes.
Definition: vtk_io.C:673
void cells_to_vtk()
write the cells from the mesh into a vtkUnstructuredGrid
Definition: vtk_io.C:564
processor_id_type processor_id() const

◆ write_nodal_data() [2/3]

void libMesh::MeshOutput< MeshBase >::write_nodal_data ( const std::string &  fname,
const NumericVector< Number > &  parallel_soln,
const std::vector< std::string > &  names 
)
virtualinherited

This method may 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.

Definition at line 149 of file mesh_output.C.

References libMesh::NumericVector< T >::localize().

152 {
153  // This is the fallback implementation for parallel I/O formats that
154  // do not yet implement proper writing in parallel, and instead rely
155  // on the full solution vector being available on all processors.
156  std::vector<Number> soln;
157  parallel_soln.localize(soln);
158  this->write_nodal_data(fname, soln, names);
159 }
virtual void write_nodal_data(const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: mesh_output.h:109
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.

◆ write_nodal_data() [3/3]

void libMesh::MeshOutput< MeshBase >::write_nodal_data ( const std::string &  fname,
const EquationSystems es,
const std::set< std::string > *  system_names 
)
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 directly uses EquationSystems current_local_solution vectors to look up nodal values.

If not implemented, reorders the solutions into a nodal-only NumericVector and calls the above version of this function.

Reimplemented in libMesh::Nemesis_IO.

Definition at line 162 of file mesh_output.C.

References libMesh::EquationSystems::build_parallel_solution_vector(), and libMesh::EquationSystems::build_variable_names().

165 {
166  std::vector<std::string> names;
167  es.build_variable_names (names, nullptr, system_names);
168 
169  std::unique_ptr<NumericVector<Number>> parallel_soln =
170  es.build_parallel_solution_vector(system_names);
171 
172  this->write_nodal_data (fname, *parallel_soln, names);
173 }
virtual void write_nodal_data(const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: mesh_output.h:109

◆ write_nodal_data_discontinuous()

virtual void libMesh::MeshOutput< MeshBase >::write_nodal_data_discontinuous ( const std::string &  ,
const std::vector< Number > &  ,
const std::vector< std::string > &   
)
inlinevirtualinherited

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

Reimplemented in libMesh::ExodusII_IO.

Definition at line 118 of file mesh_output.h.

121  { libmesh_not_implemented(); }

Member Data Documentation

◆ _compress

bool libMesh::VTKIO::_compress
private

Flag to indicate whether the output should be compressed.

Definition at line 174 of file vtk_io.h.

Referenced by set_compression(), and write_nodal_data().

◆ _element_maps

std::map< ElemMappingType, VTKIO::ElementMaps > libMesh::VTKIO::_element_maps = VTKIO::build_element_maps()
staticprivate

ElementMaps objects that are built statically and used by all instances of this class.

Definition at line 215 of file vtk_io.h.

Referenced by cells_to_vtk(), and read().

◆ _is_parallel_format

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 184 of file mesh_output.h.

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

◆ _local_node_map

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 179 of file vtk_io.h.

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

◆ _serial_only_needed_on_proc_0

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 193 of file mesh_output.h.

◆ _vtk_grid

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 169 of file vtk_io.h.

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

◆ elems_of_dimension

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

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