dune-vtk  0.2
quadraticdatacollector.hh
Go to the documentation of this file.
1 #pragma once
2 
3 #include <vector>
4 
5 #include <dune/geometry/referenceelements.hh>
6 #include <dune/grid/common/partitionset.hh>
7 #include <dune/vtk/types.hh>
8 
10 
11 namespace Dune
12 {
13  namespace Vtk
14  {
16  template <class GridView>
18  : public UnstructuredDataCollectorInterface<GridView, QuadraticDataCollector<GridView>, Partitions::All>
19  {
22 
23  public:
24  using Super::dim;
25  using Super::partition; // NOTE: quadratic data-collector currently implemented for the All partition only
26 
27  public:
28  QuadraticDataCollector (GridView const& gridView)
29  : Super(gridView)
30  {}
31 
33  std::uint64_t numPointsImpl () const
34  {
35  return gridView_.size(dim) + gridView_.size(dim-1);
36  }
37 
39 
43  template <class T>
44  std::vector<T> pointsImpl () const
45  {
46  std::vector<T> data(this->numPoints() * 3);
47  auto const& indexSet = gridView_.indexSet();
48  for (auto const& element : elements(gridView_, partition)) {
49  auto geometry = element.geometry();
50  auto refElem = referenceElement<T,dim>(element.type());
51 
52  // vertices
53  for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
54  std::size_t idx = 3 * indexSet.subIndex(element, i, dim);
55  auto v = geometry.global(refElem.position(i,dim));
56  for (std::size_t j = 0; j < v.size(); ++j)
57  data[idx + j] = T(v[j]);
58  for (std::size_t j = v.size(); j < 3u; ++j)
59  data[idx + j] = T(0);
60  }
61  // edge centers
62  for (unsigned int i = 0; i < element.subEntities(dim-1); ++i) {
63  std::size_t idx = 3 * (indexSet.subIndex(element, i, dim-1) + gridView_.size(dim));
64  auto v = geometry.global(refElem.position(i,dim-1));
65  for (std::size_t j = 0; j < v.size(); ++j)
66  data[idx + j] = T(v[j]);
67  for (std::size_t j = v.size(); j < 3u; ++j)
68  data[idx + j] = T(0);
69  }
70  }
71  return data;
72  }
73 
75  std::uint64_t numCellsImpl () const
76  {
77  return gridView_.size(0);
78  }
79 
81 
85  Cells cellsImpl () const
86  {
87  Cells cells;
88  cells.connectivity.reserve(this->numPoints());
89  cells.offsets.reserve(this->numCells());
90  cells.types.reserve(this->numCells());
91 
92  std::int64_t old_o = 0;
93  auto const& indexSet = gridView_.indexSet();
94  for (auto const& c : elements(gridView_, partition)) {
95  Vtk::CellType cellType(c.type(), Vtk::CellType::QUADRATIC);
96  for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
97  int k = cellType.permutation(j);
98  std::int64_t point_idx = indexSet.subIndex(c,k,dim);
99  cells.connectivity.push_back(point_idx);
100  }
101  for (unsigned int j = 0; j < c.subEntities(dim-1); ++j) {
102  int k = cellType.permutation(c.subEntities(dim) + j);
103  std::int64_t point_idx = (indexSet.subIndex(c,k,dim-1) + gridView_.size(dim));
104  cells.connectivity.push_back(point_idx);
105  }
106  cells.offsets.push_back(old_o += c.subEntities(dim)+c.subEntities(dim-1));
107  cells.types.push_back(cellType.type());
108  }
109  return cells;
110  }
111 
113  template <class T, class GlobalFunction>
114  std::vector<T> pointDataImpl (GlobalFunction const& fct) const
115  {
116  std::vector<T> data(this->numPoints() * fct.numComponents());
117  auto const& indexSet = gridView_.indexSet();
118  auto localFct = localFunction(fct);
119  for (auto const& e : elements(gridView_, partition)) {
120  localFct.bind(e);
122  auto refElem = referenceElement(e.geometry());
123  for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
124  int k = cellType.permutation(j);
125  std::size_t idx = fct.numComponents() * indexSet.subIndex(e, k, dim);
126  for (int comp = 0; comp < fct.numComponents(); ++comp)
127  data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim)));
128  }
129  for (unsigned int j = 0; j < e.subEntities(dim-1); ++j) {
130  int k = cellType.permutation(e.subEntities(dim) + j);
131  std::size_t idx = fct.numComponents() * (indexSet.subIndex(e, k, dim-1) + gridView_.size(dim));
132  for (int comp = 0; comp < fct.numComponents(); ++comp)
133  data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim-1)));
134  }
135  localFct.unbind();
136  }
137  return data;
138  }
139 
140  protected:
141  using Super::gridView_;
142  };
143 
144  } // end namespace Vtk
145 } // end namespace Dune
Definition: writer.hh:13
std::uint64_t numPoints() const
Return the number of points in (this partition of the) grid.
Definition: datacollectorinterface.hh:58
GridViewType GridView
Definition: datacollectorinterface.hh:25
std::uint64_t numCells() const
Return the number of cells in (this partition of the) grid.
Definition: datacollectorinterface.hh:52
Implementation of DataCollector for quadratic cells, with continuous data.
Definition: quadraticdatacollector.hh:19
QuadraticDataCollector(GridView const &gridView)
Definition: quadraticdatacollector.hh:28
Cells cellsImpl() const
Return cell types, offsets, and connectivity.
Definition: quadraticdatacollector.hh:85
std::vector< T > pointDataImpl(GlobalFunction const &fct) const
Evaluate the fct at element vertices and edge centers in the same order as the point coords.
Definition: quadraticdatacollector.hh:114
std::vector< T > pointsImpl() const
Return a vector of point coordinates.
Definition: quadraticdatacollector.hh:44
std::uint64_t numPointsImpl() const
Return number of vertices + number of edge.
Definition: quadraticdatacollector.hh:33
std::uint64_t numCellsImpl() const
Return number of grid cells.
Definition: quadraticdatacollector.hh:75
Definition: unstructureddatacollector.hh:14
std::vector< std::int64_t > offsets
Definition: unstructureddatacollector.hh:16
std::vector< std::int64_t > connectivity
Definition: unstructureddatacollector.hh:17
std::vector< std::uint8_t > types
Definition: unstructureddatacollector.hh:15
Definition: unstructureddatacollector.hh:23
@ dim
Definition: datacollectorinterface.hh:28
static constexpr auto partition
The partitionset to collect data from.
Definition: datacollectorinterface.hh:23
Cells cells() const
Return cell types, offsets, and connectivity.
Definition: unstructureddatacollector.hh:36
Mapping of Dune geometry types to VTK cell types.
Definition: types.hh:160
@ QUADRATIC
Definition: types.hh:163
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