DOFVectorView.hpp 5.05 KB
Newer Older
1
2
#pragma once

3
#include <amdis/GridFunctions.hpp>
4
#include <amdis/functions/Interpolate.hpp>
5
#include <amdis/gridfunctions/DiscreteFunction.hpp>
6

7
8
namespace AMDiS
{
9
10
11
12
13
14
15
  namespace tag
  {
    struct average {};
    struct assign {};

  } // end namespace tag

16
  /// A mutable view on the subspace of a DOFVector, \relates DiscreteFunction
17
  template <class GB, class VT, class TP>
18
  class DOFVectorView
19
      : public DiscreteFunction<GB, VT, TP>
20
  {
21
    using Self = DOFVectorView;
22
    using Super = DiscreteFunction<GB, VT, TP>;
23

24
25
    using GlobalBasis = GB;
    using TreePath = TP;
26
27

  public:
28
    /// Constructor. Stores a pointer to the mutable `dofvector`.
29
30
31
    template <class PreTreePath>
    DOFVectorView(DOFVector<GB,VT>& dofVector, PreTreePath const& preTreePath)
      : Super(dofVector, makeTreePath(preTreePath))
32
      , mutableDofVector_(&dofVector)
33
34
    {}

35
36
37
38
39
    /// Constructor forwards to the treePath constructor, with empty TreePath
    DOFVectorView(DOFVector<GB,VT>& dofVector)
      : DOFVectorView(dofVector, Dune::TypeTree::hybridTreePath())
    {}

40
41
42
  public:
    /// \brief Interpolation of GridFunction to DOFVector, assuming that there is no
    /// reference to this DOFVector in the expression.
43
44
45
    /**
     * **Example:**
     * ```
Praetorius, Simon's avatar
Praetorius, Simon committed
46
     * auto v = makeDOFVectorView(prob.solutionVector(),0);
47
     * v.interpolate_noalias([](auto const& x) { return x[0]; });
48
49
     * ```
     **/
50
51
    template <class Expr, class Tag = tag::average>
    void interpolate_noalias(Expr&& expr, Tag strategy = {})
52
    {
53
      auto const& basis = *this->basis();
54
      auto const& treePath = this->treePath();
55

56
      auto&& gridFct = makeGridFunction(FWD(expr), basis.gridView());
57

58
59
60
61
62
63
64
65
66
67
68
69
70
      if (std::is_same<Tag, tag::average>::value) {
        thread_local std::vector<std::uint8_t> counter;
        counter.clear();
        AMDiS::interpolate(basis, treePath, coefficients(), counter, FWD(gridFct), std::true_type{});

        auto& coeff = coefficients().vector();
        for (std::size_t i = 0; i < counter.size(); ++i) {
          if (counter[i] > 0)
            coeff[i] /= double(counter[i]);
        }
      } else {
        AMDiS::interpolate(basis, treePath, coefficients(), FWD(gridFct));
      }
71
72
    }

73
74
75
76
    /// \brief Interpolation of GridFunction to DOFVector
    /**
     * **Example:**
     * ```
Praetorius, Simon's avatar
Praetorius, Simon committed
77
     * auto v = makeDOFVectorView(prob.solutionVector(),0);
78
79
80
81
82
     * v.interpolate(v + [](auto const& x) { return x[0]; });
     * ```
     * Allows to have a reference to the DOFVector in the expression, e.g. as
     * \ref DiscreteFunction or \ref gradientAtQP() of a DiscreteFunction.
     **/
83
84
    template <class Expr, class Tag = tag::average>
    void interpolate(Expr&& expr, Tag strategy = {})
85
    {
86
      // create temporary copy of data
87
      DOFVector<GB,VT> tmp(coefficients());
88
      Self tmpView{tmp, this->treePath()};
89
      tmpView.interpolate_noalias(FWD(expr), strategy);
90

91
92
      // move data from temporary vector into stored DOFVector
      coefficients().vector() = std::move(tmp.vector());
93
94
    }

95
96
97
    /// \brief Interpolation of GridFunction to DOFVector, alias to \ref interpolate()
    template <class Expr>
    DOFVectorView& operator<<(Expr&& expr)
98
    {
99
100
      interpolate(expr);
      return *this;
101
102
    }

103
    /// \brief interpolate `(*this) + expr` to DOFVector
Praetorius, Simon's avatar
Praetorius, Simon committed
104
    template <class Expr>
105
    DOFVectorView& operator+=(Expr&& expr)
106
    {
107
      interpolate((*this) + expr);
108
109
      return *this;
    }
110

111
    /// \brief interpolate `(*this) - expr` to DOFVector
112
    template <class Expr>
113
    DOFVectorView& operator-=(Expr&& expr)
114
    {
115
116
      interpolate((*this) - expr);
      return *this;
117
118
    }

119
    /// Return the mutable DOFVector
120
    DOFVector<GB,VT>& coefficients() { return *mutableDofVector_; }
121

122
    /// Return the const DOFVector
Praetorius, Simon's avatar
Praetorius, Simon committed
123
    using Super::coefficients;
124

125
  protected:
126
    DOFVector<GB,VT>* mutableDofVector_;
127
128
  };

129

130
#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
131
132
133
134
135
  // Deduction guide for DOFVectorView class
  template <class GlobalBasis, class ValueType, class PreTreePath>
  DOFVectorView(DOFVector<GlobalBasis, ValueType>& dofVector, PreTreePath const& preTreePath)
    -> DOFVectorView<GlobalBasis, ValueType, TYPEOF(makeTreePath(preTreePath))>;

136
137
138
139
140
141
  // Deduction guide for DOFVectorView class
  template <class GlobalBasis, class ValueType>
  DOFVectorView(DOFVector<GlobalBasis, ValueType>& dofVector)
    -> DOFVectorView<GlobalBasis, ValueType, Dune::TypeTree::HybridTreePath<>>;
#endif

142
  /// A Generator for a mutable \ref DOFVectorView
143
144
  template <class GlobalBasis, class ValueType, class PreTreePath>
  auto makeDOFVectorView(DOFVector<GlobalBasis, ValueType>& dofVector, PreTreePath const& preTreePath)
145
  {
146
147
    auto treePath = makeTreePath(preTreePath);
    return DOFVectorView<GlobalBasis, ValueType, decltype(treePath)>{dofVector, treePath};
148
149
  }

150
  /// A Generator for a mutable \ref DOFVectorView
151
152
  template <class GlobalBasis, class ValueType>
  auto makeDOFVectorView(DOFVector<GlobalBasis, ValueType>& dofVector)
153
154
  {
    auto treePath = Dune::TypeTree::hybridTreePath();
155
    return DOFVectorView<GlobalBasis, ValueType, Dune::TypeTree::HybridTreePath<>>{dofVector, treePath};
156
157
  }

158
} // end namespace AMDiS