1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2020 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 #ifndef LIBMESH_VTK_IO_H
20 #define LIBMESH_VTK_IO_H
21 
22 // Local includes
23 #include "libmesh/libmesh_common.h"
24 #include "libmesh/mesh_input.h"
25 #include "libmesh/mesh_output.h"
26 #include "libmesh/utility.h"
27 
28 #ifdef LIBMESH_HAVE_VTK
29 // Ignore "deprecated...header" warning from strstream
30 #include "libmesh/ignore_warnings.h"
31 #include "vtkType.h"
32 #include "vtkSmartPointer.h"
33 #include "libmesh/restore_warnings.h"
34 #endif
35 
36 // C++ includes
37 #include <cstddef>
38 #include <map>
39 
40 // Forward declarations
41 class vtkUnstructuredGrid;
42 
43 namespace libMesh
44 {
45 
46 class MeshBase;
47 
48 /**
49  * This class implements reading and writing meshes in the VTK format.
50  * Format description:
51  * cf. <a href="http://www.vtk.org/">VTK home page</a>.
52  *
53  * This class will not have any functionality unless VTK is detected
54  * during configure and hence LIBMESH_HAVE_VTK is defined.
55  *
56  * \author Wout Ruijter
57  * \author John W. Peterson
58  * \date 2007
59  */
60 class VTKIO : public MeshInput<MeshBase>,
61               public MeshOutput<MeshBase>
62 {
63 public:
64   /**
65    * Constructor.  Takes a writable reference to a mesh object.
66    * This is the constructor required to read a mesh.
67    */
68   explicit
69   VTKIO (MeshBase & mesh);
70 
71   /**
72    * Constructor.  Takes a read-only reference to a mesh object.
73    * This is the constructor required to write a mesh.
74    */
75   explicit
76   VTKIO (const MeshBase & mesh);
77 
78   /**
79    * Bring in base class functionality for name resolution and to
80    * avoid warnings about hidden overloaded virtual functions.
81    */
82   using MeshOutput<MeshBase>::write_nodal_data;
83 
84   /**
85    * This method implements writing a mesh with nodal data to a
86    * specified file where the nodal data and variable names are provided.
87    *
88    * As with other Mesh IO classes, this interface is still still
89    * "available" when !LIBMESH_HAVE_VTK, however, it will throw a
90    * runtime error.
91    */
92   virtual void write_nodal_data (const std::string &,
93                                  const std::vector<Number> &,
94                                  const std::vector<std::string> &) override;
95 
96   /**
97    * This method implements reading a mesh from a specified file
98    * in VTK format.
99    *
100    * As with other Mesh IO classes, this interface is still still
101    * "available" when !LIBMESH_HAVE_VTK, however, it will throw a
102    * runtime error.
103    */
104   virtual void read (const std::string &) override;
105 
106   /**
107    * Output the mesh without solutions to a .pvtu file.
108    *
109    * As with other Mesh IO classes, this interface is still still
110    * "available" when !LIBMESH_HAVE_VTK, however, it will throw a
111    * runtime error.
112    */
113   virtual void write (const std::string &) override;
114 
115 #ifdef LIBMESH_HAVE_VTK
116 
117   /**
118    * Setter for compression flag
119    */
120   void set_compression(bool b);
121 
122   /**
123    * Get a pointer to the VTK unstructured grid data structure.
124    */
125   vtkUnstructuredGrid * get_vtk_grid();
126 
127 private:
128   /**
129    * write the nodes from the mesh into a vtkUnstructuredGrid and update the
130    * local_node_map.
131    */
132   void nodes_to_vtk();
133 
134   /**
135    * write the cells from the mesh into a vtkUnstructuredGrid
136    */
137   void cells_to_vtk();
138 
139   /**
140    * write the nodal values of soln to a vtkUnstructuredGrid
141    */
142   void node_values_to_vtk(const std::string & name,
143                           const std::vector<Real> & local_values);
144 
145   /**
146    * Extract the values of soln that correspond to the nodes
147    *
148    * This method overwrites all values in local_values
149    */
150   void get_local_node_values(std::vector<Number> & local_values,
151                              std::size_t variable,
152                              const std::vector<Number> & soln,
153                              const std::vector<std::string> & names);
154 
155   /**
156    * Write the system vectors to vtk
157    *
158    * This function is not currently used by anything, so it is commented
159    * out, and may eventually be removed entirely.
160    */
161   // void system_vectors_to_vtk(const EquationSystems & es,
162   //                            vtkUnstructuredGrid * & grid);
163 
164   /**
165    * pointer to the VTK grid. the vtkSmartPointer will automatically
166    * initialize the value to null and keep track of reference
167    * counting.
168    */
169   vtkSmartPointer<vtkUnstructuredGrid> _vtk_grid;
170 
171   /**
172    * Flag to indicate whether the output should be compressed
173    */
174   bool _compress;
175 
176   /**
177    * maps global node id to node id of partition
178    */
179   std::map<dof_id_type, dof_id_type> _local_node_map;
180 
181 
182   /**
183    * Helper object that holds a map from VTK to libMesh element types
184    * and vice-versa.
185    */
186   struct ElementMaps
187   {
188     // Associate libmesh_type with vtk_type (searchable in both directions).
associateElementMaps189     void associate(ElemType libmesh_type, vtkIdType vtk_type)
190     {
191       writing_map[libmesh_type] = vtk_type;
192       reading_map[vtk_type] = libmesh_type;
193     }
194 
195     // Find an entry in the writing map, or throw an error.
findElementMaps196     vtkIdType find(ElemType libmesh_type)
197     {
198       return libmesh_map_find(writing_map, libmesh_type);
199     }
200 
201     // Find an entry in the reading map, or throw an error.
findElementMaps202     ElemType find(vtkIdType vtk_type)
203     {
204       return libmesh_map_find(reading_map, vtk_type);
205     }
206 
207     std::map<ElemType, vtkIdType> writing_map;
208     std::map<vtkIdType, ElemType> reading_map;
209   };
210 
211   /**
212    * ElementMaps object that is built statically and used by
213    * all instances of this class.
214    */
215   static ElementMaps _element_maps;
216 
217   /**
218    * Static function used to construct the _element_maps struct.
219    */
220   static ElementMaps build_element_maps();
221 
222 #endif
223 };
224 
225 
226 
227 } // namespace libMesh
228 
229 
230 #endif // LIBMESH_VTK_IO_H
231