libMesh
vtk_io.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 // C++ includes
20 #include <fstream>
21 
22 // Local includes
23 #include "libmesh/libmesh_config.h"
24 #include "libmesh/vtk_io.h"
25 #include "libmesh/mesh_base.h"
26 #include "libmesh/equation_systems.h"
27 #include "libmesh/numeric_vector.h"
28 #include "libmesh/system.h"
29 #include "libmesh/node.h"
30 #include "libmesh/elem.h"
31 
32 #ifdef LIBMESH_HAVE_VTK
33 
34 // I get a lot of "warning: extra ';' inside a class [-Wextra-semi]" from clang
35 // on VTK header files.
36 #include "libmesh/ignore_warnings.h"
37 
38 #include "vtkXMLUnstructuredGridReader.h"
39 #include "vtkXMLUnstructuredGridWriter.h"
40 #include "vtkXMLPUnstructuredGridWriter.h"
41 #include "vtkUnstructuredGrid.h"
42 #include "vtkIntArray.h"
43 #include "vtkCellArray.h"
44 #include "vtkCellData.h"
45 #include "vtkConfigure.h"
46 #include "vtkDoubleArray.h"
47 #include "vtkGenericCell.h"
48 #include "vtkPointData.h"
49 #include "vtkPoints.h"
50 #include "vtkSmartPointer.h"
51 
52 #include "libmesh/restore_warnings.h"
53 
54 // A convenient macro for comparing VTK versions. Returns 1 if the
55 // current VTK version is < major.minor.subminor and zero otherwise.
56 //
57 // It relies on the VTK version numbers detected during configure. Note that if
58 // LIBMESH_HAVE_VTK is not defined, none of the LIBMESH_DETECTED_VTK_VERSION_* variables will
59 // be defined either.
60 #define VTK_VERSION_LESS_THAN(major,minor,subminor) \
61  ((LIBMESH_DETECTED_VTK_VERSION_MAJOR < (major) || \
62  (LIBMESH_DETECTED_VTK_VERSION_MAJOR == (major) && (LIBMESH_DETECTED_VTK_VERSION_MINOR < (minor) || \
63  (LIBMESH_DETECTED_VTK_VERSION_MINOR == (minor) && \
64  LIBMESH_DETECTED_VTK_VERSION_SUBMINOR < (subminor))))) ? 1 : 0)
65 
66 #endif // LIBMESH_HAVE_VTK
67 
68 
69 
70 namespace libMesh
71 {
72 
73 // Constructor for reading
75  MeshInput<MeshBase> (mesh, /*is_parallel_format=*/true),
76  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
77 #ifdef LIBMESH_HAVE_VTK
78  ,_compress(false)
79 #endif
80 {
81 }
82 
83 
84 
85 // Constructor for writing
86 VTKIO::VTKIO (const MeshBase & mesh) :
87  MeshOutput<MeshBase>(mesh, /*is_parallel_format=*/true)
88 #ifdef LIBMESH_HAVE_VTK
89  ,_compress(false)
90 #endif
91 {
92 }
93 
94 
95 
96 // Output the mesh without solutions to a .pvtu file
97 void VTKIO::write (const std::string & name)
98 {
99  std::vector<Number> soln;
100  std::vector<std::string> names;
101  this->write_nodal_data(name, soln, names);
102 }
103 
104 
105 
106 // The rest of the file is wrapped in ifdef LIBMESH_HAVE_VTK except for
107 // a couple of "stub" functions at the bottom.
108 #ifdef LIBMESH_HAVE_VTK
109 
110 // Initialize the static _element_maps struct.
112 
113 // Static function which constructs the ElementMaps object.
115 {
116  // Object to be filled up
117  ElementMaps em;
118 
119  em.associate(EDGE2, VTK_LINE);
120  em.associate(EDGE3, VTK_QUADRATIC_EDGE);
121  em.associate(TRI3, VTK_TRIANGLE);
122  em.associate(TRI6, VTK_QUADRATIC_TRIANGLE);
123  em.associate(QUAD4, VTK_QUAD);
124  em.associate(QUAD8, VTK_QUADRATIC_QUAD);
125  em.associate(TET4, VTK_TETRA);
126  em.associate(TET10, VTK_QUADRATIC_TETRA);
127  em.associate(HEX8, VTK_HEXAHEDRON);
128  em.associate(HEX20, VTK_QUADRATIC_HEXAHEDRON);
129  em.associate(HEX27, VTK_TRIQUADRATIC_HEXAHEDRON);
130  em.associate(PRISM6, VTK_WEDGE);
131  em.associate(PRISM15, VTK_QUADRATIC_WEDGE);
132  em.associate(PRISM18, VTK_BIQUADRATIC_QUADRATIC_WEDGE);
133  em.associate(PYRAMID5, VTK_PYRAMID);
134 
135  // VTK_BIQUADRATIC_QUAD has been around since VTK 5.0
136 #if VTK_MAJOR_VERSION > 5 || (VTK_MAJOR_VERSION == 5 && VTK_MINOR_VERSION > 0)
137  em.associate(QUAD9, VTK_BIQUADRATIC_QUAD);
138 #endif
139 
140  // TRI3SUBDIVISION is for writing only
141  em.writing_map[TRI3SUBDIVISION] = VTK_TRIANGLE;
142 
143  return em;
144 }
145 
146 
147 
148 void VTKIO::read (const std::string & name)
149 {
150  // This is a serial-only process for now;
151  // the Mesh should be read on processor 0 and
152  // broadcast later
153  libmesh_assert_equal_to (MeshOutput<MeshBase>::mesh().processor_id(), 0);
154 
155  // Keep track of what kinds of elements this file contains
156  elems_of_dimension.clear();
157  elems_of_dimension.resize(4, false);
158 
159  // Use a typedef, because these names are just crazy
160  typedef vtkSmartPointer<vtkXMLUnstructuredGridReader> MyReader;
161  MyReader reader = MyReader::New();
162 
163  // Pass the filename along to the reader
164  reader->SetFileName(name.c_str());
165 
166  // Force reading
167  reader->Update();
168 
169  // read in the grid
170  _vtk_grid = reader->GetOutput();
171 
172  // Get a reference to the mesh
174 
175  // Clear out any pre-existing data from the Mesh
176  mesh.clear();
177 
178  // Get the number of points from the _vtk_grid object
179  const unsigned int vtk_num_points = static_cast<unsigned int>(_vtk_grid->GetNumberOfPoints());
180 
181  // always numbered nicely so we can loop like this
182  for (unsigned int i=0; i<vtk_num_points; ++i)
183  {
184  // add to the id map
185  // and add the actual point
186  double pnt[3];
187  _vtk_grid->GetPoint(static_cast<vtkIdType>(i), pnt);
188  Point xyz(pnt[0], pnt[1], pnt[2]);
189  mesh.add_point(xyz, i);
190  }
191 
192  // Get the number of cells from the _vtk_grid object
193  const unsigned int vtk_num_cells = static_cast<unsigned int>(_vtk_grid->GetNumberOfCells());
194 
195  vtkSmartPointer<vtkGenericCell> cell = vtkSmartPointer<vtkGenericCell>::New();
196  for (unsigned int i=0; i<vtk_num_cells; ++i)
197  {
198  _vtk_grid->GetCell(i, cell);
199 
200  // Get the libMesh element type corresponding to this VTK element type.
201  ElemType libmesh_elem_type = _element_maps.find(cell->GetCellType());
202  Elem * elem = Elem::build(libmesh_elem_type).release();
203 
204  // get the straightforward numbering from the VTK cells
205  for (unsigned int j=0; j<elem->n_nodes(); ++j)
206  elem->set_node(j) =
207  mesh.node_ptr(cast_int<dof_id_type>(cell->GetPointId(j)));
208 
209  // then get the connectivity
210  std::vector<dof_id_type> conn;
211  elem->connectivity(0, VTK, conn);
212 
213  // then reshuffle the nodes according to the connectivity, this
214  // two-time-assign would evade the definition of the vtk_mapping
215  for (std::size_t j=0; j<conn.size(); ++j)
216  elem->set_node(j) = mesh.node_ptr(conn[j]);
217 
218  elem->set_id(i);
219 
220  elems_of_dimension[elem->dim()] = true;
221 
222  mesh.add_elem(elem);
223  } // end loop over VTK cells
224 
225  // Set the mesh dimension to the largest encountered for an element
226  for (unsigned char i=0; i!=4; ++i)
227  if (elems_of_dimension[i])
228  mesh.set_mesh_dimension(i);
229 
230 #if LIBMESH_DIM < 3
231  if (mesh.mesh_dimension() > LIBMESH_DIM)
232  libmesh_error_msg("Cannot open dimension " \
233  << mesh.mesh_dimension() \
234  << " mesh file when configured without " \
235  << mesh.mesh_dimension() \
236  << "D support.");
237 #endif // LIBMESH_DIM < 3
238 }
239 
240 
241 
242 void VTKIO::write_nodal_data (const std::string & fname,
243  const std::vector<Number> & soln,
244  const std::vector<std::string> & names)
245 {
246  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
247 
248  // Warn that the .pvtu file extension should be used. Paraview
249  // recognizes this, and it works in both serial and parallel. Only
250  // warn about this once.
251  if (fname.substr(fname.rfind("."), fname.size()) != ".pvtu")
252  libmesh_do_once(libMesh::err << "The .pvtu extension should be used when writing VTK files in libMesh.");
253 
254  // If there are variable names being written, the solution vector
255  // should not be empty, it should have been broadcast to all
256  // processors by the MeshOutput base class, since VTK is a parallel
257  // format. Verify this before going further.
258  if (!names.empty() && soln.empty())
259  libmesh_error_msg("Empty soln vector in VTKIO::write_nodal_data().");
260 
261  // we only use Unstructured grids
262  _vtk_grid = vtkSmartPointer<vtkUnstructuredGrid>::New();
263  vtkSmartPointer<vtkXMLPUnstructuredGridWriter> writer = vtkSmartPointer<vtkXMLPUnstructuredGridWriter>::New();
264 
265  // add nodes to the grid and update _local_node_map
266  _local_node_map.clear();
267  this->nodes_to_vtk();
268 
269  // add cells to the grid
270  this->cells_to_vtk();
271 
272  // add nodal solutions to the grid, if solutions are given
273  if (names.size() > 0)
274  {
275  std::size_t num_vars = names.size();
276  dof_id_type num_nodes = mesh.n_nodes();
277 
278  for (std::size_t variable=0; variable<num_vars; ++variable)
279  {
280  vtkSmartPointer<vtkDoubleArray> data = vtkSmartPointer<vtkDoubleArray>::New();
281  data->SetName(names[variable].c_str());
282 
283  // number of local and ghost nodes
284  data->SetNumberOfValues(_local_node_map.size());
285 
286  // loop over all nodes and get the solution for the current
287  // variable, if the node is in the current partition
288  for (dof_id_type k=0; k<num_nodes; ++k)
289  {
290  std::map<dof_id_type, dof_id_type>::iterator local_node_it = _local_node_map.find(k);
291  if (local_node_it == _local_node_map.end())
292  continue; // not a local node
293 
294 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
295  libmesh_do_once (libMesh::err << "Only writing the real part for complex numbers!\n"
296  << "if you need this support contact " << LIBMESH_PACKAGE_BUGREPORT
297  << std::endl);
298  data->SetValue(local_node_it->second, soln[k*num_vars + variable].real());
299 #else
300  data->SetValue(local_node_it->second, soln[k*num_vars + variable]);
301 #endif
302  }
303  _vtk_grid->GetPointData()->AddArray(data);
304  }
305  }
306 
307  // Tell the writer how many partitions exist and on which processor
308  // we are currently
309  writer->SetNumberOfPieces(MeshOutput<MeshBase>::mesh().n_processors());
310  writer->SetStartPiece(MeshOutput<MeshBase>::mesh().processor_id());
311  writer->SetEndPiece(MeshOutput<MeshBase>::mesh().processor_id());
312 
313  // partitions overlap by one node
314  // FIXME: According to this document
315  // http://paraview.org/Wiki/images/5/51/SC07_tut107_ParaView_Handouts.pdf
316  // the ghosts are cells rather than nodes.
317  writer->SetGhostLevel(1);
318 
319  // VTK 6 replaces SetInput() with SetInputData(). See
320  // http://www.vtk.org/Wiki/VTK/VTK_6_Migration/Replacement_of_SetInput
321  // for the full explanation.
322 #if VTK_VERSION_LESS_THAN(6,0,0)
323  writer->SetInput(_vtk_grid);
324 #else
325  writer->SetInputData(_vtk_grid);
326 #endif
327 
328  writer->SetFileName(fname.c_str());
329  writer->SetDataModeToAscii();
330 
331  // compress the output, if desired (switches also to binary)
332  if (this->_compress)
333  {
334 #if !VTK_VERSION_LESS_THAN(5,6,0)
335  writer->SetCompressorTypeToZLib();
336 #else
337  libmesh_do_once(libMesh::err << "Compression not implemented with old VTK libs!" << std::endl;);
338 #endif
339  }
340 
341  writer->Write();
342 
343 }
344 
345 
346 
347 vtkUnstructuredGrid * VTKIO::get_vtk_grid()
348 {
349  return _vtk_grid;
350 }
351 
352 
353 
355 {
356  this->_compress = b;
357 }
358 
359 
360 
362 {
363  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
364 
365  // containers for points and coordinates of points
366  vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
367  vtkSmartPointer<vtkDoubleArray> pcoords = vtkSmartPointer<vtkDoubleArray>::New();
368  // if this grid is to be used in VTK then the dimension of the points should be 3
369  pcoords->SetNumberOfComponents(LIBMESH_DIM);
370  pcoords->Allocate(3*mesh.n_local_nodes());
371  points->SetNumberOfPoints(mesh.n_local_nodes()); // it seems that it needs this to prevent a segfault
372 
373  unsigned int local_node_counter = 0;
374 
375  for (const auto & node_ptr : mesh.local_node_ptr_range())
376  {
377  const Node & node = *node_ptr;
378 
379  double pnt[3] = {0, 0, 0};
380  for (unsigned int i=0; i<LIBMESH_DIM; ++i)
381  pnt[i] = node(i);
382 
383  // Fill mapping between global and local node numbers
384  _local_node_map[node.id()] = local_node_counter;
385 
386  // add point
387 #if VTK_VERSION_LESS_THAN(7,1,0)
388  pcoords->InsertNextTupleValue(pnt);
389 #else
390  pcoords->InsertNextTuple(pnt);
391 #endif
392  ++local_node_counter;
393  }
394 
395  // add coordinates to points
396  points->SetData(pcoords);
397 
398  // add points to grid
399  _vtk_grid->SetPoints(points);
400 }
401 
402 
403 
405 {
406  const MeshBase & mesh = MeshOutput<MeshBase>::mesh();
407 
408  vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
409  vtkSmartPointer<vtkIdList> pts = vtkSmartPointer<vtkIdList>::New();
410 
411  std::vector<int> types(mesh.n_active_local_elem());
412 
413  vtkSmartPointer<vtkIntArray> elem_id = vtkSmartPointer<vtkIntArray>::New();
414  elem_id->SetName("libmesh_elem_id");
415  elem_id->SetNumberOfComponents(1);
416 
417  vtkSmartPointer<vtkIntArray> subdomain_id = vtkSmartPointer<vtkIntArray>::New();
418  subdomain_id->SetName("subdomain_id");
419  subdomain_id->SetNumberOfComponents(1);
420 
421  vtkSmartPointer<vtkIntArray> elem_proc_id = vtkSmartPointer<vtkIntArray>::New();
422  elem_proc_id->SetName("processor_id");
423  elem_proc_id->SetNumberOfComponents(1);
424 
425  unsigned active_element_counter = 0;
426  for (const auto & elem : mesh.active_local_element_ptr_range())
427  {
428  pts->SetNumberOfIds(elem->n_nodes());
429 
430  // get the connectivity for this element
431  std::vector<dof_id_type> conn;
432  elem->connectivity(0, VTK, conn);
433 
434  for (std::size_t i=0; i<conn.size(); ++i)
435  {
436  // If the node ID is not found in the _local_node_map, we'll
437  // add it to the _vtk_grid. NOTE[JWP]: none of the examples
438  // I have actually enters this section of code...
439  if (_local_node_map.find(conn[i]) == _local_node_map.end())
440  {
441  dof_id_type global_node_id = elem->node_id(i);
442 
443  const Point & the_node = mesh.point(global_node_id);
444 
445  // InsertNextPoint accepts either a double or float array of length 3.
446  double pt[3] = {0., 0., 0.};
447  for (unsigned int d=0; d<LIBMESH_DIM; ++d)
448  pt[d] = the_node(d);
449 
450  // Insert the point into the _vtk_grid
451  vtkIdType local = _vtk_grid->GetPoints()->InsertNextPoint(pt);
452 
453  // Update the _local_node_map with the ID returned by VTK
454  _local_node_map[global_node_id] =
455  cast_int<dof_id_type>(local);
456  }
457 
458  // Otherwise, the node ID was found in the _local_node_map, so
459  // insert it into the vtkIdList.
460  pts->InsertId(i, _local_node_map[conn[i]]);
461  }
462 
463  vtkIdType vtkcellid = cells->InsertNextCell(pts);
464  types[active_element_counter] = cast_int<int>(_element_maps.find(elem->type()));
465 
466  elem_id->InsertTuple1(vtkcellid, elem->id());
467  subdomain_id->InsertTuple1(vtkcellid, elem->subdomain_id());
468  elem_proc_id->InsertTuple1(vtkcellid, elem->processor_id());
469  ++active_element_counter;
470  } // end loop over active elements
471 
472  _vtk_grid->SetCells(&types[0], cells);
473  _vtk_grid->GetCellData()->AddArray(elem_id);
474  _vtk_grid->GetCellData()->AddArray(subdomain_id);
475  _vtk_grid->GetCellData()->AddArray(elem_proc_id);
476 }
477 
478 
479 
487 // void VTKIO::system_vectors_to_vtk(const EquationSystems & es,
488 // vtkUnstructuredGrid *& grid)
489 // {
490 // if (MeshOutput<MeshBase>::mesh().processor_id() == 0)
491 // {
492 // std::map<std::string, std::vector<Number>> vecs;
493 // for (unsigned int i=0; i<es.n_systems(); ++i)
494 // {
495 // const System & sys = es.get_system(i);
496 // System::const_vectors_iterator v_end = sys.vectors_end();
497 // System::const_vectors_iterator it = sys.vectors_begin();
498 // for (; it!= v_end; ++it)
499 // {
500 // // for all vectors on this system
501 // std::vector<Number> values;
502 // // libMesh::out<<"it "<<it->first<<std::endl;
503 //
504 // it->second->localize_to_one(values, 0);
505 // // libMesh::out<<"finish localize"<<std::endl;
506 // vecs[it->first] = values;
507 // }
508 // }
509 //
510 // std::map<std::string, std::vector<Number>>::iterator it = vecs.begin();
511 //
512 // for (; it!=vecs.end(); ++it)
513 // {
514 // vtkSmartPointer<vtkDoubleArray> data = vtkSmartPointer<vtkDoubleArray>::New();
515 // data->SetName(it->first.c_str());
516 // libmesh_assert_equal_to (it->second.size(), es.get_mesh().n_nodes());
517 // data->SetNumberOfValues(it->second.size());
518 //
519 // for (std::size_t i=0; i<it->second.size(); ++i)
520 // {
521 // #ifdef LIBMESH_USE_COMPLEX_NUMBERS
522 // libmesh_do_once (libMesh::err << "Only writing the real part for complex numbers!\n"
523 // << "if you need this support contact " << LIBMESH_PACKAGE_BUGREPORT
524 // << std::endl);
525 // data->SetValue(i, it->second[i].real());
526 // #else
527 // data->SetValue(i, it->second[i]);
528 // #endif
529 //
530 // }
531 // grid->GetPointData()->AddArray(data);
532 // }
533 // }
534 // }
535 
536 
537 
538 #else // !LIBMESH_HAVE_VTK
539 
540 void VTKIO::read (const std::string & name)
541 {
542  libmesh_error_msg("Cannot read VTK file: " << name \
543  << "\nYou must have VTK installed and correctly configured to read VTK meshes.");
544 }
545 
546 
547 
548 void VTKIO::write_nodal_data (const std::string & fname,
549  const std::vector<Number> &,
550  const std::vector<std::string> &)
551 {
552  libmesh_error_msg("Cannot write VTK file: " << fname \
553  << "\nYou must have VTK installed and correctly configured to read VTK meshes.");
554 }
555 
556 
557 #endif // LIBMESH_HAVE_VTK
558 
559 
560 
561 } // namespace libMesh
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
void set_compression(bool b)
Setter for compression flag.
Definition: vtk_io.C:354
OStreamProxy err
virtual void connectivity(const unsigned int sc, const IOPackage iop, std::vector< dof_id_type > &conn) const =0
std::map< dof_id_type, dof_id_type > _local_node_map
maps global node id to node id of partition
Definition: vtk_io.h:161
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:1941
A Node is like a Point, but with more information.
Definition: node.h:52
virtual const Point & point(const dof_id_type i) const =0
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:238
virtual SimpleRange< node_iterator > local_node_ptr_range()=0
dof_id_type n_local_nodes() const
Definition: mesh_base.h:272
std::vector< bool > elems_of_dimension
A vector of bools describing what dimension elements have been encountered when reading a mesh...
Definition: mesh_input.h:97
vtkIdType find(ElemType libmesh_type)
Definition: vtk_io.h:178
ElemType
Defines an enum for geometric element types.
vtkSmartPointer< vtkUnstructuredGrid > _vtk_grid
Write the system vectors to vtk.
Definition: vtk_io.h:151
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
bool _compress
Flag to indicate whether the output should be compressed.
Definition: vtk_io.h:156
virtual const Node * node_ptr(const dof_id_type i) const =0
const MT & mesh() const
Definition: mesh_output.h:216
This class defines an abstract interface for Mesh output.
Definition: mesh_output.h:53
static ElementMaps _element_maps
ElementMaps object that is built statically and used by all instances of this class.
Definition: vtk_io.h:207
The libMesh namespace provides an interface to certain functionality in the library.
Helper object that holds a map from VTK to libMesh element types and vice-versa.
Definition: vtk_io.h:168
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.
This is the MeshBase class.
Definition: mesh_base.h:68
dof_id_type & set_id()
Definition: dof_object.h:641
virtual void read(const std::string &) libmesh_override
This method implements reading a mesh from a specified file in VTK format.
Definition: vtk_io.C:148
virtual unsigned int n_nodes() const =0
void nodes_to_vtk()
write the nodes from the mesh into a vtkUnstructuredGrid
Definition: vtk_io.C:361
virtual SimpleRange< element_iterator > active_local_element_ptr_range()=0
This class defines an abstract interface for Mesh input.
Definition: mesh_base.h:50
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
static ElementMaps build_element_maps()
Static function used to construct the _element_maps struct.
Definition: vtk_io.C:114
vtkUnstructuredGrid * get_vtk_grid()
Get a pointer to the VTK unstructured grid data structure.
Definition: vtk_io.C:347
void set_mesh_dimension(unsigned char d)
Resets the logical dimension of the mesh.
Definition: mesh_base.h:199
std::map< ElemType, vtkIdType > writing_map
Definition: vtk_io.h:199
virtual void clear()
Deletes all the data that are currently stored.
Definition: mesh_base.C:285
dof_id_type n_active_local_elem() const
Definition: mesh_base.h:389
virtual void write(const std::string &) libmesh_override
Output the mesh without solutions to a .pvtu file.
VTKIO(MeshBase &mesh)
Constructor.
Definition: vtk_io.C:74
void associate(ElemType libmesh_type, vtkIdType vtk_type)
Definition: vtk_io.h:171
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
void cells_to_vtk()
write the cells from the mesh into a vtkUnstructuredGrid
Definition: vtk_io.C:404
virtual unsigned int dim() const =0
IterBase * data
Ideally this private member data should have protected access.
dof_id_type id() const
Definition: dof_object.h:632
virtual void write_nodal_data(const std::string &, const std::vector< Number > &, const std::vector< std::string > &) libmesh_override
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: vtk_io.C:242
virtual dof_id_type n_nodes() const =0
processor_id_type processor_id()
Definition: libmesh_base.h:96
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
uint8_t dof_id_type
Definition: id_types.h:64
processor_id_type n_processors()
Definition: libmesh_base.h:88