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