dune-vtk  0.2
lagrangedatacollector.hh
Go to the documentation of this file.
1 #pragma once
2 
3 #include <cassert>
4 #include <map>
5 #include <vector>
6 
7 #include <dune/geometry/referenceelements.hh>
8 #include <dune/grid/common/partitionset.hh>
9 #include <dune/grid/common/partitionset.hh>
10 #include <dune/vtk/types.hh>
12 
14 
15 namespace Dune
16 {
17  namespace Vtk
18  {
20  template <class GridView, int ORDER = -1>
22  : public UnstructuredDataCollectorInterface<GridView, LagrangeDataCollector<GridView,ORDER>, Partitions::All>
23  {
26 
27  public:
28  static_assert(ORDER != 0, "Order 0 not supported");
29  using Super::dim;
30  using Super::partition; // NOTE: lagrange data-collector currently implemented for the All partition only
31 
32  public:
33  LagrangeDataCollector (GridView const& gridView, int order = ORDER)
34  : Super(gridView)
35  , order_(order)
36  {
37  assert(order > 0 && "Order 0 not supported");
38  assert(ORDER < 0 || order == ORDER);
39  }
40 
42  void updateImpl ()
43  {
44  auto const& indexSet = gridView_.indexSet();
45 
46  pointSets_.clear();
47  for (auto gt : indexSet.types(0))
48  pointSets_.emplace(gt, order_);
49 
50  for (auto& pointSet : pointSets_)
51  pointSet.second.build(pointSet.first);
52 
53  numPoints_ = indexSet.size(dim);
54  for (auto const& pointSet : pointSets_) {
55  auto gt = pointSet.first;
56  auto refElem = referenceElement<double,dim>(gt);
57  numPoints_ += (pointSet.second.size() - refElem.size(dim)) * indexSet.size(gt);
58  }
59  }
60 
62  std::uint64_t numPointsImpl () const
63  {
64  return numPoints_;
65  }
66 
68 
72  template <class T>
73  std::vector<T> pointsImpl () const
74  {
75  std::vector<T> data(this->numPoints() * 3);
76  auto const& indexSet = gridView_.indexSet();
77 
78  std::size_t shift = indexSet.size(dim);
79 
80  for (auto const& element : elements(gridView_, partition)) {
81  auto geometry = element.geometry();
82  auto refElem = referenceElement<T,dim>(element.type());
83 
84  auto const& pointSet = pointSets_.at(element.type());
85  unsigned int vertexDOFs = refElem.size(dim);
86  unsigned int innerDOFs = pointSet.size() - vertexDOFs;
87 
88  for (std::size_t i = 0; i < pointSet.size(); ++i) {
89  auto const& p = pointSet[i];
90  if (i < vertexDOFs)
91  assert(p.localKey().codim() == dim);
92 
93  auto const& localKey = p.localKey();
94  std::size_t idx = 3 * (localKey.codim() == dim
95  ? indexSet.subIndex(element, localKey.subEntity(), dim)
96  : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
97 
98  auto v = geometry.global(p.point());
99  for (std::size_t j = 0; j < v.size(); ++j)
100  data[idx + j] = T(v[j]);
101  for (std::size_t j = v.size(); j < 3u; ++j)
102  data[idx + j] = T(0);
103  }
104  }
105  return data;
106  }
107 
109  std::uint64_t numCellsImpl () const
110  {
111  return gridView_.size(0);
112  }
113 
115 
119  Cells cellsImpl () const
120  {
121  Cells cells;
122  cells.connectivity.reserve(this->numPoints());
123  cells.offsets.reserve(this->numCells());
124  cells.types.reserve(this->numCells());
125 
126  auto const& indexSet = gridView_.indexSet();
127  std::size_t shift = indexSet.size(dim);
128 
129  std::int64_t old_o = 0;
130  for (auto const& element : elements(gridView_, partition)) {
131  auto refElem = referenceElement<double,dim>(element.type());
132  Vtk::CellType cellType(element.type(), Vtk::LAGRANGE);
133 
134  auto const& pointSet = pointSets_.at(element.type());
135  unsigned int vertexDOFs = refElem.size(dim);
136  unsigned int innerDOFs = pointSet.size() - vertexDOFs;
137 
138  for (std::size_t i = 0; i < pointSet.size(); ++i) {
139  auto const& p = pointSet[i];
140  auto const& localKey = p.localKey();
141  std::size_t idx = (localKey.codim() == dim
142  ? indexSet.subIndex(element, localKey.subEntity(), dim)
143  : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
144  cells.connectivity.push_back(idx);
145  }
146 
147  cells.offsets.push_back(old_o += pointSet.size());
148  cells.types.push_back(cellType.type());
149  }
150  return cells;
151  }
152 
154  template <class T, class GlobalFunction>
155  std::vector<T> pointDataImpl (GlobalFunction const& fct) const
156  {
157  int nComps = fct.numComponents();
158  std::vector<T> data(this->numPoints() * nComps);
159  auto const& indexSet = gridView_.indexSet();
160 
161  std::size_t shift = indexSet.size(dim);
162 
163  auto localFct = localFunction(fct);
164  for (auto const& element : elements(gridView_, partition)) {
165  localFct.bind(element);
166  auto refElem = referenceElement<T,dim>(element.type());
167 
168  auto const& pointSet = pointSets_.at(element.type());
169  unsigned int vertexDOFs = refElem.size(dim);
170  unsigned int innerDOFs = pointSet.size() - vertexDOFs;
171 
172  for (std::size_t i = 0; i < pointSet.size(); ++i) {
173  auto const& p = pointSet[i];
174  auto const& localKey = p.localKey();
175  std::size_t idx = nComps * (localKey.codim() == dim
176  ? indexSet.subIndex(element, localKey.subEntity(), dim)
177  : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
178 
179  for (int comp = 0; comp < nComps; ++comp)
180  data[idx + comp] = T(localFct.evaluate(comp, p.point()));
181  }
182  localFct.unbind();
183  }
184  return data;
185  }
186 
187  protected:
188  using Super::gridView_;
189 
190  unsigned int order_;
191  std::uint64_t numPoints_ = 0;
192 
194  std::map<GeometryType, PointSet> pointSets_;
195  };
196 
197  } // end namespace Vtk
198 } // end namespace Dune
Definition: datacollectorinterface.hh:9
@ LAGRANGE
Definition: types.hh:151
std::uint64_t numPoints() const
Return the number of points in (this partition of the) grid.
Definition: datacollectorinterface.hh:58
std::uint64_t numCells() const
Return the number of cells in (this partition of the) grid.
Definition: datacollectorinterface.hh:52
Implementation of DataCollector for lagrange cells.
Definition: lagrangedatacollector.hh:23
std::uint64_t numPointsImpl() const
Return number of lagrange nodes.
Definition: lagrangedatacollector.hh:62
void updateImpl()
Construct the point sets.
Definition: lagrangedatacollector.hh:42
Cells cellsImpl() const
Return cell types, offsets, and connectivity.
Definition: lagrangedatacollector.hh:119
std::vector< T > pointsImpl() const
Return a vector of point coordinates.
Definition: lagrangedatacollector.hh:73
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: lagrangedatacollector.hh:155
LagrangeDataCollector(GridView const &gridView, int order=ORDER)
Definition: lagrangedatacollector.hh:33
std::map< GeometryType, PointSet > pointSets_
Definition: lagrangedatacollector.hh:194
unsigned int order_
Definition: lagrangedatacollector.hh:190
std::uint64_t numPoints_
Definition: lagrangedatacollector.hh:191
std::uint64_t numCellsImpl() const
Return number of grid cells.
Definition: lagrangedatacollector.hh:109
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:191
std::uint8_t type() const
Return VTK Cell type.
Definition: types.hh:196
A set of lagrange points compatible with the numbering of VTK and Gmsh.
Definition: lagrangepoints.hh:30