Caribou
VTKReader.h
1 #pragma once
2 
3 #include <string>
4 #include <iostream>
5 #include <utility>
6 #include <functional>
7 #include <unordered_map>
8 
9 #include <Caribou/config.h>
10 #include <Caribou/Topology/Mesh.h>
11 #include <Caribou/Topology/Domain.h>
12 
13 #include <vtkSmartPointer.h>
14 #include <vtkUnstructuredGridReader.h>
15 #include <vtkCellType.h>
16 
17 namespace caribou::topology::io {
18 
19 template<UNSIGNED_INTEGER_TYPE Dimension>
20 class VTKReader {
21  static_assert(Dimension == 1 or Dimension == 2 or Dimension == 3, "The VTKReader can only read 1D, 2D or 3D fields.");
22 public:
23 
24  using ElementsIndices = Eigen::Matrix<UNSIGNED_INTEGER_TYPE, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>;
25  using DomainBuilder = std::function<BaseDomain* (Mesh<Dimension> &, const ElementsIndices &)>;
26  using MeshType = Mesh<Dimension>;
27 
29  static auto Read(const std::string & filepath) -> VTKReader;
30 
32  void print (std::ostream &out) const;
33 
35  [[nodiscard]]
36  auto mesh() const -> MeshType;
37 
39  template<typename Element>
40  auto register_element_type(const VTKCellType & vtkCellType) -> VTKReader & {
41  p_domain_builders[vtkCellType] = [](MeshType & m, const ElementsIndices & indices) {
42  return m.template add_domain<Element>("domain_"+std::to_string(m.number_of_domains()+1), indices);
43  };
44  return *this;
45  }
46 
84  auto register_element_type(const VTKCellType & vtkCellType, DomainBuilder builder) -> VTKReader & {
85  p_domain_builders[vtkCellType] = builder;
86  return *this;
87  }
88 
89 private:
90  VTKReader(std::string filepath, vtkSmartPointer<vtkUnstructuredGridReader> reader, std::array<UNSIGNED_INTEGER_TYPE, Dimension> axes);
91 
92  const std::string p_filepath;
93  vtkSmartPointer<vtkUnstructuredGridReader> p_reader;
94  std::array<UNSIGNED_INTEGER_TYPE, Dimension> p_axes;
95  std::unordered_map<VTKCellType, DomainBuilder> p_domain_builders;
96 };
97 
98 extern template class VTKReader<1>;
99 extern template class VTKReader<2>;
100 extern template class VTKReader<3>;
101 
102 template<UNSIGNED_INTEGER_TYPE Dimension>
103 auto operator<<(std::ostream& os, const caribou::topology::io::VTKReader<Dimension> & t) -> std::ostream&
104 {
105  t.print(os);
106  return os;
107 }
108 
109 }
caribou::topology::io::VTKReader
Definition: VTKReader.h:20
caribou::topology::io::VTKReader::print
void print(std::ostream &out) const
Print information about the current vtk file.
Definition: VTKReader.cpp:159
caribou::topology::io::VTKReader::mesh
auto mesh() const -> MeshType
Build the actual unstructured mesh from the vtk file.
Definition: VTKReader.cpp:95
caribou::topology::Mesh
Definition: Mesh.h:16
caribou::topology::BaseDomain
A domain is subspace of a mesh containing a set of points and the topological relation between them.
Definition: BaseDomain.h:13
caribou::topology::io::VTKReader::Read
static auto Read(const std::string &filepath) -> VTKReader
Build a new VTKReader instance by reading a vtk file.
Definition: VTKReader.cpp:75
caribou::topology::Mesh::number_of_domains
auto number_of_domains() const -> UNSIGNED_INTEGER_TYPE final
Definition: Mesh.h:96
caribou::topology::io::VTKReader::register_element_type
auto register_element_type(const VTKCellType &vtkCellType) -> VTKReader &
Register an element type to the given VTK cell type.
Definition: VTKReader.h:40
caribou::topology::io::VTKReader::register_element_type
auto register_element_type(const VTKCellType &vtkCellType, DomainBuilder builder) -> VTKReader &
Register an element type to the given VTK cell type and give it a domain builder callback function.
Definition: VTKReader.h:84