dune-vtk  0.2
types.hh
Go to the documentation of this file.
1 #pragma once
2 
3 #include <cstdint>
4 #include <map>
5 #include <string>
6 #include <vector>
7 
8 #include <dune/common/ftraits.hh>
9 #include <dune/common/typelist.hh>
10 #include <dune/common/version.hh>
11 #include <dune/geometry/type.hh>
12 #include <dune/grid/io/file/vtk/common.hh>
15 
16 namespace Dune
17 {
18  namespace Vtk
19  {
21  enum FormatTypes {
22  ASCII = 1<<0,
23  BINARY = 1<<1,
24  COMPRESSED = 1<<2,
26  };
27  std::string to_string (Vtk::FormatTypes);
28 
30  Vtk::FormatTypes formatTypeOf(Dune::VTK::OutputType);
31 
32 
34  enum class RangeTypes {
35  UNSPECIFIED, //< The output components are not restricted
36  AUTO, //< Detect the category automatically from number of components
37  SCALAR, //< Use exactly 1 component
38  VECTOR, //< Use exactly 3 components
39  TENSOR //< Use exactly 9 components
40  };
41  std::string to_string (Vtk::RangeTypes);
42 
43  // Map a dune-grid FieldInfo::Type to ValueTypes
44  Vtk::RangeTypes rangeTypeOf (Dune::VTK::FieldInfo::Type);
45 
46  // Map a number of components to a corresponding value type
47  Vtk::RangeTypes rangeTypeOf (int ncomps);
48 
49 
50  enum DataTypes {
51  UNKNOWN = 0,
56  FLOAT32 = 32,
57  FLOAT64 = 64
58  };
59  std::string to_string (Vtk::DataTypes);
60 
61  // Map a dune-grid Precision type to DataTypes
62  Vtk::DataTypes dataTypeOf (Dune::VTK::Precision);
63 
64  // Map a string to DataTypes
65  Vtk::DataTypes dataTypeOf (std::string);
66 
67  // Map the field_type of T to DataTypes
68  template <class T>
70  {
71  using F = typename FieldTraits<T>::field_type;
72  if constexpr (std::is_same_v<F, std::int8_t>) { return Vtk::DataTypes::INT8; }
73  if constexpr (std::is_same_v<F, std::uint8_t>) { return Vtk::DataTypes::UINT8; }
74  if constexpr (std::is_same_v<F, std::int16_t>) { return Vtk::DataTypes::INT16; }
75  if constexpr (std::is_same_v<F, std::uint16_t>) { return Vtk::DataTypes::UINT16; }
76  if constexpr (std::is_same_v<F, std::int32_t>) { return Vtk::DataTypes::INT32; }
77  if constexpr (std::is_same_v<F, std::uint32_t>) { return Vtk::DataTypes::UINT32; }
78  if constexpr (std::is_same_v<F, std::int64_t>) { return Vtk::DataTypes::INT64; }
79  if constexpr (std::is_same_v<F, std::uint64_t>) { return Vtk::DataTypes::UINT64; }
80  if constexpr (std::is_same_v<F, float>) { return Vtk::DataTypes::FLOAT32; }
81  if constexpr (std::is_same_v<F, double>) { return Vtk::DataTypes::FLOAT64; }
82  if constexpr (std::is_same_v<F, long double>) { return Vtk::DataTypes::FLOAT64; }
84  }
85 
86  template <class> struct NoConstraint : std::true_type {};
87 
89  template <template <class> class C = NoConstraint, class Caller>
90  void mapDataTypes (Vtk::DataTypes t, Caller caller)
91  {
92  switch (t) {
93  case INT8: if constexpr(C<std::int8_t>::value) caller(MetaType<std::int8_t>{}); break;
94  case UINT8: if constexpr(C<std::uint8_t>::value) caller(MetaType<std::uint8_t>{}); break;
95  case INT16: if constexpr(C<std::int16_t>::value) caller(MetaType<std::int16_t>{}); break;
96  case UINT16: if constexpr(C<std::uint16_t>::value) caller(MetaType<std::uint16_t>{}); break;
97  case INT32: if constexpr(C<std::int32_t>::value) caller(MetaType<std::int32_t>{}); break;
98  case UINT32: if constexpr(C<std::uint32_t>::value) caller(MetaType<std::uint32_t>{}); break;
99  case INT64: if constexpr(C<std::int64_t>::value) caller(MetaType<std::int64_t>{}); break;
100  case UINT64: if constexpr(C<std::uint64_t>::value) caller(MetaType<std::uint64_t>{}); break;
101  case FLOAT32: if constexpr(C<float>::value) caller(MetaType<float>{}); break;
102  case FLOAT64: if constexpr(C<double>::value) caller(MetaType<double>{}); break;
103  default:
104  VTK_ASSERT_MSG(false, "Unsupported type " + to_string(t));
105  break;
106  }
107  }
108 
110  template <template <class> class Constraint1 = NoConstraint,
111  template <class> class Constraint2 = NoConstraint,
112  class Caller>
113  void mapDataTypes (Vtk::DataTypes t1, Vtk::DataTypes t2, Caller caller)
114  {
115  mapDataTypes<Constraint1>(t1, [&](auto type1) {
116  mapDataTypes<Constraint2>(t2, [&](auto type2) {
117  caller(type1, type2);
118  });
119  });
120  }
121 
123  template <template <class> class Constraint1 = NoConstraint,
124  template <class> class Constraint2 = NoConstraint,
125  template <class> class Constraint3 = NoConstraint,
126  class Caller>
127  void mapDataTypes (Vtk::DataTypes t1, Vtk::DataTypes t2, Vtk::DataTypes t3, Caller caller)
128  {
129  mapDataTypes<Constraint1>(t1, [&](auto type1) {
130  mapDataTypes<Constraint2>(t2, [&](auto type2) {
131  mapDataTypes<Constraint3>(t3, [&](auto type3) {
132  caller(type1, type2, type3);
133  });
134  });
135  });
136  }
137 
138 
140  NONE = 0,
143  LZMA
144  };
145  std::string to_string (CompressorTypes);
146 
147 
151  LAGRANGE
152  };
153 
154 
155  enum CellTypes : std::uint8_t {
156  // Linear VTK cell types
157  VERTEX = 1,
158  /* POLY_VERTEX = 2, // not supported */
159  LINE = 3,
160  /* POLY_LINE = 4, // not supported */
161  TRIANGLE = 5,
162  /* TRIANGLE_STRIP = 6, // not supported */
163  POLYGON = 7,
164  /* PIXEL = 8, // not supported */
165  QUAD = 9,
166  TETRA = 10,
167  /* VOXEL = 11, // not supported */
169  WEDGE = 13,
170  PYRAMID = 14,
171  // Quadratic VTK cell types
177  // Arbitrary order Lagrange elements
185  };
186  GeometryType to_geometry (std::uint8_t);
187 
188 
190  class CellType
191  {
192  public:
193  CellType (GeometryType const& t, CellParametrization = LINEAR);
194 
196  std::uint8_t type () const
197  {
198  return type_;
199  }
200 
202  int permutation (int idx) const
203  {
204  return permutation_[idx];
205  }
206 
207  bool noPermutation () const
208  {
209  return noPermutation_;
210  }
211 
212  private:
213  std::uint8_t type_;
214  std::vector<int> permutation_;
215  bool noPermutation_ = true;
216  };
217 
218 
219  class FieldInfo
220  {
221  public:
222  template <class... Args>
223  explicit FieldInfo (std::string name, Args... args)
224  : name_(std::move(name))
225  , ncomps_(getArg<int,unsigned int,long,unsigned long>(args..., 1))
226  , rangeType_(getArg<Vtk::RangeTypes>(args..., Vtk::RangeTypes::AUTO))
227  , dataType_(getArg<Vtk::DataTypes>(args..., Vtk::DataTypes::FLOAT32))
228  {
229  if (rangeType_ == Vtk::RangeTypes::AUTO)
230  rangeType_ = rangeTypeOf(ncomps_);
231  }
232 
233  // Construct from dune-grid FieldInfo
234  FieldInfo (Dune::VTK::FieldInfo info)
235  : FieldInfo(info.name(), info.size(), rangeTypeOf(info.type()), dataTypeOf(info.precision()))
236  {}
237 
239  std::string const& name () const
240  {
241  return name_;
242  }
243 
245  int size () const
246  {
247  return ncomps_;
248  }
249 
252  {
253  return rangeType_;
254  }
255 
258  {
259  return dataType_;
260  }
261 
262  private:
263  std::string name_;
264  int ncomps_;
265  Vtk::RangeTypes rangeType_;
266  Vtk::DataTypes dataType_;
267  };
268 
269  } // end namespace Vtk
270 } // end namespace Dune
Macro for wrapping error checks and throwing exceptions.
#define VTK_ASSERT_MSG(cond, text)
check if condition cond holds; otherwise, throw a VtkError with a message.
Definition: errors.hh:19
Definition: datacollectorinterface.hh:9
Vtk::DataTypes dataTypeOf(Dune::VTK::Precision p)
Definition: types.cc:98
std::string to_string(Vtk::FormatTypes type)
Definition: types.cc:12
FormatTypes
Type used for representing the output format.
Definition: types.hh:21
@ BINARY
Definition: types.hh:23
@ APPENDED
Definition: types.hh:25
@ COMPRESSED
Definition: types.hh:24
@ ASCII
Definition: types.hh:22
void mapDataTypes(Vtk::DataTypes t, Caller caller)
Map a given enum DataType to a type passed to Caller as MetaType.
Definition: types.hh:90
CompressorTypes
Definition: types.hh:139
@ NONE
Definition: types.hh:140
@ ZLIB
Definition: types.hh:141
@ LZ4
Definition: types.hh:142
@ LZMA
Definition: types.hh:143
Vtk::FormatTypes formatTypeOf(Dune::VTK::OutputType o)
Map the dune-grid OutputType to FormatTypes.
Definition: types.cc:25
RangeTypes
Type used to determine whether to limit output components to e.g. 3 (vector), or 9 (tensor)
Definition: types.hh:34
CellTypes
Definition: types.hh:155
@ VERTEX
Definition: types.hh:157
@ QUADRATIC_QUAD
Definition: types.hh:174
@ LAGRANGE_TETRAHEDRON
Definition: types.hh:181
@ WEDGE
Definition: types.hh:169
@ QUADRATIC_HEXAHEDRON
Definition: types.hh:176
@ TETRA
Definition: types.hh:166
@ LAGRANGE_QUADRILATERAL
Definition: types.hh:180
@ QUADRATIC_TRIANGLE
Definition: types.hh:173
@ LAGRANGE_TRIANGLE
Definition: types.hh:179
@ LINE
Definition: types.hh:159
@ LAGRANGE_WEDGE
Definition: types.hh:183
@ HEXAHEDRON
Definition: types.hh:168
@ LAGRANGE_PYRAMID
Definition: types.hh:184
@ LAGRANGE_CURVE
Definition: types.hh:178
@ LAGRANGE_HEXAHEDRON
Definition: types.hh:182
@ TRIANGLE
Definition: types.hh:161
@ QUADRATIC_TETRA
Definition: types.hh:175
@ POLYGON
Definition: types.hh:163
@ QUAD
Definition: types.hh:165
@ PYRAMID
Definition: types.hh:170
@ QUADRATIC_EDGE
Definition: types.hh:172
DataTypes
Definition: types.hh:50
@ UINT64
Definition: types.hh:55
@ UINT8
Definition: types.hh:52
@ UINT16
Definition: types.hh:53
@ INT32
Definition: types.hh:54
@ INT64
Definition: types.hh:55
@ INT16
Definition: types.hh:53
@ UINT32
Definition: types.hh:54
@ FLOAT64
Definition: types.hh:57
@ UNKNOWN
Definition: types.hh:51
@ INT8
Definition: types.hh:52
@ FLOAT32
Definition: types.hh:56
decltype(auto) getArg(Arg0 &&arg0, Args &&... args)
Definition: arguments.hh:29
GeometryType to_geometry(std::uint8_t cell)
Definition: types.cc:145
Vtk::RangeTypes rangeTypeOf(Dune::VTK::FieldInfo::Type t)
Definition: types.cc:56
CellParametrization
Definition: types.hh:148
@ LINEAR
Definition: types.hh:149
@ QUADRATIC
Definition: types.hh:150
@ LAGRANGE
Definition: types.hh:151
Definition: types.hh:86
Mapping of Dune geometry types to VTK cell types.
Definition: types.hh:191
CellType(GeometryType const &t, CellParametrization=LINEAR)
Definition: types.cc:178
bool noPermutation() const
Definition: types.hh:207
std::uint8_t type() const
Return VTK Cell type.
Definition: types.hh:196
int permutation(int idx) const
Return a permutation of Dune elemenr vertices to conform to VTK element numbering.
Definition: types.hh:202
Definition: types.hh:220
int size() const
The number of components in the data field.
Definition: types.hh:245
std::string const & name() const
The name of the data field.
Definition: types.hh:239
FieldInfo(Dune::VTK::FieldInfo info)
Definition: types.hh:234
Vtk::RangeTypes rangeType() const
Return the category of the stored range.
Definition: types.hh:251
Vtk::DataTypes dataType() const
Return the data tpe of the data field.
Definition: types.hh:257
FieldInfo(std::string name, Args... args)
Definition: types.hh:223