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

83  :
84  MeshInput<MeshBase> (mesh, /*is_parallel_format=*/true),
85  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
86 #ifdef LIBMESH_HAVE_VTK
87  ,_compress(false)
88 #endif
89 {
90 }
91 
92 
93 
94 // Constructor for writing
95 VTKIO::VTKIO (const MeshBase & mesh) :
96  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
97 #ifdef LIBMESH_HAVE_VTK
98  ,_compress(false)
99 #endif
100 {
101 }
102 
103 
104 
105 // Output the mesh without solutions to a .pvtu file
106 void VTKIO::write (const std::string & name)
107 {
108  std::vector<Number> soln;
109  std::vector<std::string> names;
110  this->write_nodal_data(name, soln, names);
111 }
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:83

◆ 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 123 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.

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

◆ cells_to_vtk()

void libMesh::VTKIO::cells_to_vtk ( )
private

write the cells from the mesh into a vtkUnstructuredGrid

Definition at line 565 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().

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

678 {
679  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
680  std::size_t num_vars = names.size();
681  dof_id_type num_nodes = mesh.n_nodes();
682 
683  local_values.clear();
684  local_values.resize(_local_node_map.size(), 0.0);
685 
686  // loop over all nodes and get the solution for the current
687  // variable, if the node is in the current partition
688  for (dof_id_type k=0; k<num_nodes; ++k)
689  if (const auto local_node_it = _local_node_map.find(k);
690  local_node_it != _local_node_map.end())
691  local_values[local_node_it->second] = soln[k*num_vars + variable];
692 }
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 471 of file vtk_io.C.

References _vtk_grid.

472 {
473  return _vtk_grid;
474 }
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 655 of file vtk_io.C.

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

Referenced by write_nodal_data().

657 {
658  vtkSmartPointer<vtkDoubleArray> data = vtkSmartPointer<vtkDoubleArray>::New();
659  data->SetName(name.c_str());
660 
661  libmesh_assert_equal_to(_local_node_map.size(), local_values.size());
662 
663  // number of local and ghost nodes
664  data->SetNumberOfValues(_local_node_map.size());
665 
666  // copy values into vtk
667  for (auto i : index_range(local_values)) {
668  data->SetValue(i, double(local_values[i]));
669  }
670 
671  _vtk_grid->GetPointData()->AddArray(data);
672 }
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 485 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().

486 {
487  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
488 
489  // containers for points and coordinates of points
490  vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
491  vtkSmartPointer<vtkDoubleArray> pcoords = vtkSmartPointer<vtkDoubleArray>::New();
492  // if this grid is to be used in VTK then the dimension of the points should be 3
493  pcoords->SetNumberOfComponents(LIBMESH_DIM);
494  pcoords->Allocate(3*mesh.n_local_nodes());
495  points->SetNumberOfPoints(mesh.n_local_nodes()); // it seems that it needs this to prevent a segfault
496 
497  // SetRationalWeights() was introduced in VTK 9.0
498 #if !VTK_VERSION_LESS_THAN(9,0,0)
499  bool have_weights = false;
500  int weight_index = 0;
501  vtkSmartPointer<vtkDoubleArray> rational_weights;
502 
504  {
505  rational_weights = vtkSmartPointer<vtkDoubleArray>::New();
506  rational_weights->SetName("RationalWeights");
507  rational_weights->SetNumberOfComponents(1);
508  weight_index = static_cast<int>(mesh.default_mapping_data());
509  have_weights = true;
510  }
511 #endif
512 
513  vtkSmartPointer<vtkIntArray> node_id = vtkSmartPointer<vtkIntArray>::New();
514  node_id->SetName("libmesh_node_id");
515  node_id->SetNumberOfComponents(1);
516 
517  unsigned int local_node_counter = 0;
518 
519  for (const auto & node_ptr : mesh.local_node_ptr_range())
520  {
521  const Node & node = *node_ptr;
522 
523  double pnt[3] = {0, 0, 0};
524  for (unsigned int i=0; i<LIBMESH_DIM; ++i)
525  pnt[i] = double(node(i));
526 
527  // Fill mapping between global and local node numbers
528  _local_node_map[node.id()] = local_node_counter;
529 
530  // add point
531 #if VTK_VERSION_LESS_THAN(7,1,0)
532  pcoords->InsertNextTupleValue(pnt);
533 #else
534  pcoords->InsertNextTuple(pnt);
535 #endif
536 #if !VTK_VERSION_LESS_THAN(9,0,0)
537  if (have_weights)
538  {
539  Real weight = node.get_extra_datum<Real>(weight_index);
540  rational_weights->InsertTuple1(local_node_counter, double(weight));
541  }
542 #endif
543 
544  node_id->InsertTuple1(local_node_counter, node.id());
545 
546  ++local_node_counter;
547  }
548 
549  // add coordinates to points
550  points->SetData(pcoords);
551 
552  // add points to grid
553  _vtk_grid->SetPoints(points);
554  _vtk_grid->GetPointData()->AddArray(node_id);
555 
556 #if !VTK_VERSION_LESS_THAN(9,0,0)
557  if (have_weights)
558  _vtk_grid->GetPointData()->SetRationalWeights(rational_weights);
559 
560 #endif
561 }
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 186 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().

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

References _compress.

479 {
480  this->_compress = b;
481 }
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:1871

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

References _compress, _local_node_map, _vtk_grid, cells_to_vtk(), libMesh::ParallelObject::comm(), libMesh::Utility::ends_with(), 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().

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