OperatorTermBase.hpp 3.84 KB
Newer Older
1
2
#pragma once

3
4
5
6
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif

7
8
#include <vector>
#include <type_traits>
9

10
#include <dune/geometry/quadraturerules.hh>
11

12
13
#include <dune/amdis/OperatorEvaluation.hpp>
#include <dune/amdis/common/ValueCategory.hpp>
14

15
namespace AMDiS
16
{
17

18
  /// Abstract base class of all operator terms
19
20
21
22
23
24
  template <class MeshView>
  class OperatorTerm
  {
  protected:
    using Codim0  = typename MeshView::template Codim<0>;
    using Element = typename Codim0::Entity;
25

26
27
    static constexpr int dim = MeshView::dimension;
    static constexpr int dow = MeshView::dimensionworld;
28

29
30
    using QuadratureRule = Dune::QuadratureRule<double, dim>;
    using PointList = std::vector<Dune::QuadraturePoint<double, dim>>;
31

32
  public:
Praetorius, Simon's avatar
Praetorius, Simon committed
33
    virtual void init(Element const& element, PointList const& points) = 0;
34

35
    virtual double evalZot(std::size_t iq,
36
        Dune::FieldVector<double,1> const& test,
37
        Dune::FieldVector<double,1> const trial = 1.0) const = 0;
38

39
    virtual double evalFot1(std::size_t iq,
40
        Dune::FieldVector<double,1> const& test,
41
        Dune::FieldVector<double,dow> const& grad_trial) const = 0;
42

43
    virtual double evalFot2(std::size_t iq,
44
        Dune::FieldVector<double,dow> const& grad_test,
45
        Dune::FieldVector<double,1> const trial = 1.0) const = 0;
46

47
    virtual double evalSot(std::size_t iq,
48
49
        Dune::FieldVector<double,dow> const& grad_test,
        Dune::FieldVector<double,dow> const& grad_trial) const = 0;
50
51

    virtual int getDegree(Dune::GeometryType const& t) const = 0;
52
  };
53
54
55



56
57
  /// Base class for all operators based on expressions
  template <class MeshView, class Term, class Traits = tag::none>
58
  class GenericOperatorTerm
59
      : public OperatorTerm<MeshView>
Praetorius, Simon's avatar
Praetorius, Simon committed
60
      , public OperatorEvaluation // NOTE: better: use static functions
61
  {
62
63
64
    using Super   = OperatorTerm<MeshView>;
    using Element = typename Super::Element;
    using PointList = typename Super::PointList;
65

66
67
    static constexpr int dim = MeshView::dimension;
    static constexpr int dow = MeshView::dimensionworld;
68
69

  public:
70
    GenericOperatorTerm(Term const& term, Traits traits = {})
71
      : term(term)
72
      , traits(traits)
73
74
    {}

Praetorius, Simon's avatar
Praetorius, Simon committed
75
    virtual void init(Element const& element, PointList const& points) override
76
77
    {
      term.init(element, points);
78

79
80
      // cache term evaluation
      values.resize(points.size());
81
      for (std::size_t iq = 0; iq < points.size(); ++iq)
82
83
        values[iq] = term[iq];
    }
84

85
    virtual double evalZot(std::size_t iq,
86
        Dune::FieldVector<double,1> const& test,
87
88
        Dune::FieldVector<double,1> const trial = 1.0) const override
    {
89
      return this->evalZotImpl(_cat{}, traits, values[iq], test, trial);
90
    }
91

92
    virtual double evalFot1(std::size_t iq,
93
        Dune::FieldVector<double,1> const& test,
94
        Dune::FieldVector<double,dow> const& grad_trial) const override
95
    {
96
      return this->evalFotImpl(_cat{}, traits, values[iq], grad_trial, test);
97
    }
98

99
    virtual double evalFot2(std::size_t iq,
100
        Dune::FieldVector<double,dow> const& grad_test,
101
102
        Dune::FieldVector<double,1> const trial = 1.0) const override
    {
103
      return this->evalFotImpl(_cat{}, traits, values[iq], grad_test, trial);
104
    }
105

106
    virtual double evalSot(std::size_t iq,
107
108
        Dune::FieldVector<double,dow> const& grad_test,
        Dune::FieldVector<double,dow> const& grad_trial) const override
109
    {
110
      return this->evalSotImpl(_cat{}, traits, values[iq], grad_test, grad_trial);
111
    }
112
113

    virtual int getDegree(Dune::GeometryType const& t) const override
114
    {
115
      return term.getDegree(t);
116
    }
117

118
119
  private:
    Term term;
120
    Traits traits;
121

122
    using value_type = std::decay_t< decltype( std::declval<Term>()[std::size_t(0)] ) >;
Praetorius, Simon's avatar
Praetorius, Simon committed
123
    using _cat = ValueCategory_t<value_type>;
124

Praetorius, Simon's avatar
Praetorius, Simon committed
125
    std::vector<value_type> values; // NOTE: maybe caching is not necessary here, since cached already in the term
126
  };
127

128
} // end namespace AMDiS