datacollectorinterface.hh 2.97 KB
Newer Older
1
2
#pragma once

3
#include <dune/vtk/vtktypes.hh>
4

5
namespace Dune {
6

Praetorius, Simon's avatar
Praetorius, Simon committed
7
template <class GridView, class Derived>
8
9
10
class DataCollectorInterface
{
public:
Praetorius, Simon's avatar
Praetorius, Simon committed
11
12
13
14
  DataCollectorInterface (GridView const& gridView)
    : gridView_(gridView)
  {}

15
  /// Update the DataCollector on the current GridView
Praetorius, Simon's avatar
Praetorius, Simon committed
16
17
18
19
  void update ()
  {
    asDerived().updateImpl();
  }
20

21
22
23
24
25
26
  /// Return the number of overlapping elements
  int ghostLevel () const
  {
    return asDerived().ghostLevelImpl();
  }

27
28
29
30
31
32
  /// \brief Return the number of points in the grid
  std::uint64_t numPoints () const
  {
    return asDerived().numPointsImpl();
  }

33
34
35
36
37
38
39
40
  /// \brief Return a flat vector of point coordinates
  /**
   * All coordinates are extended to 3 components and concatenated.
   * [p0_x, p0_y, p0_z, p1_x, p1_y, p1_z, ...]
   * If the GridView::dimensionworld < 3, the remaining components are
   * set to 0
   **/
  template <class T>
Praetorius, Simon's avatar
Praetorius, Simon committed
41
42
43
44
  std::vector<T> points () const
  {
    return asDerived().template pointsImpl<T>();
  }
45
46
47
48
49
50
51
52
53
54

  /// \brief Return a flat vector of function values evaluated at the grid points.
  /**
   * In case of a vector valued function, flat the vector entries:
   * [fct(p0)_0, fct(p0)_1, fct(p0)_2, fct(p1)_0, ...]
   * where the vector dimension must be 3 (possible extended by 0s)
   * In case of tensor valued function, flat the tensor row-wise:
   * [fct(p0)_00, fct(p0)_01, fct(p0)_02, fct(p0)_10, fct(p0)_11, fct(p0)_12, fct(p0)_20...]
   * where the tensor dimension must be 3x3 (possible extended by 0s)
   **/
55
56
  template <class T, class VtkFunction>
  std::vector<T> pointData (VtkFunction const& fct) const
Praetorius, Simon's avatar
Praetorius, Simon committed
57
58
59
  {
    return asDerived().template pointDataImpl<T>(fct);
  }
60
61

  /// \brief Return a flat vector of function values evaluated at the grid cells. \see pointData.
62
63
  template <class T, class VtkFunction>
  std::vector<T> cellData (VtkFunction const& fct) const
Praetorius, Simon's avatar
Praetorius, Simon committed
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
  {
    return asDerived().template cellDataImpl<T>(fct);
  }

protected: // cast to derived type

  Derived& asDerived ()
  {
    return static_cast<Derived&>(*this);
  }

  const Derived& asDerived () const
  {
    return static_cast<const Derived&>(*this);
  }

80
public: // default implementations
Praetorius, Simon's avatar
Praetorius, Simon committed
81

82
83
84
85
86
87
88
89
90
  void updateImpl ()
  {
    /* do nothing */
  }

  int ghostLevelImpl () const
  {
    return gridView_.overlapSize(0);
  }
Praetorius, Simon's avatar
Praetorius, Simon committed
91
92

  // Evaluate `fct` in center of cell
93
94
  template <class T, class VtkFunction>
  std::vector<T> cellDataImpl (VtkFunction const& fct) const
Praetorius, Simon's avatar
Praetorius, Simon committed
95
  {
96
    std::vector<T> data(gridView_.size(0) * fct.ncomps());
Praetorius, Simon's avatar
Praetorius, Simon committed
97
98
    auto const& indexSet = gridView_.indexSet();
    auto localFct = localFunction(fct);
99
    for (auto const& e : elements(gridView_, Partitions::all)) {
Praetorius, Simon's avatar
Praetorius, Simon committed
100
      localFct.bind(e);
101
      auto refElem = referenceElement<T,GridView::dimension>(e.type());
Praetorius, Simon's avatar
Praetorius, Simon committed
102
103
      std::size_t idx = fct.ncomps() * indexSet.index(e);
      for (int comp = 0; comp < fct.ncomps(); ++comp)
104
        data[idx + comp] = T(localFct.evaluate(comp, refElem.position(0,0)));
Praetorius, Simon's avatar
Praetorius, Simon committed
105
106
107
108
109
110
111
      localFct.unbind();
    }
    return data;
  }

protected:
  GridView gridView_;
112
113
};

114
} // end namespace Dune