OperatorTermBase.hpp 3.75 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
    static constexpr int dim = Element::dimension;
27
    
28
29
    using QuadratureRule = Dune::QuadratureRule<double, dim>;
    using PointList = std::vector<Dune::QuadraturePoint<double, dim>>;
30
      
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
  public:
    virtual void init(Element const& element, 
                      PointList const& points) = 0;
                      
    virtual double evalZot(size_t iq, 
        Dune::FieldVector<double,1> const& test, 
        Dune::FieldVector<double,1> const trial = 1.0) const = 0;
        
    virtual double evalFot1(size_t iq, 
        Dune::FieldVector<double,1> const& test, 
        Dune::FieldVector<double,dim> const& grad_trial) const = 0;
        
    virtual double evalFot2(size_t iq, 
        Dune::FieldVector<double,dim> const& grad_test, 
        Dune::FieldVector<double,1> const trial = 1.0) const = 0;
        
    virtual double evalSot(size_t iq, 
        Dune::FieldVector<double,dim> const& grad_test, 
        Dune::FieldVector<double,dim> const& grad_trial) const = 0;
                          
    virtual int getDegree() const = 0;
  };
53
54
55
  
  
  
56
57
58
59
60
  /// Base class for all operators based on expressions
  template <class MeshView, class Term, class Traits = tag::none>
  class GenericOperatorTerm 
      : public OperatorTerm<MeshView>
      , public OperatorEvaluation
61
  {
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
    using Super   = OperatorTerm<MeshView>;
    using Element = typename Super::Element;
    using PointList = typename Super::PointList;
    
    static constexpr int dim = Element::dimension;

  public:
    GenericOperatorTerm(Term const& term)
      : term(term)
    {}

    virtual void init(Element const& element, 
                      PointList const& points) override
    {
      term.init(element, points);
77
      
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
      // cache term evaluation
      values.resize(points.size());
      for (size_t iq = 0; iq < points.size(); ++iq)
        values[iq] = term[iq];
    }
                      
    virtual double evalZot(size_t iq, 
        Dune::FieldVector<double,1> const& test, 
        Dune::FieldVector<double,1> const trial = 1.0) const override
    {
      return this->evalZotImpl(_cat{}, _traits{}, values[iq], test, trial);
    }
        
    virtual double evalFot1(size_t iq, 
        Dune::FieldVector<double,1> const& test, 
        Dune::FieldVector<double,dim> const& grad_trial) const override
    {
      return this->evalFotImpl(_cat{}, _traits{}, values[iq], grad_trial, test);
    }
        
    virtual double evalFot2(size_t iq, 
        Dune::FieldVector<double,dim> const& grad_test, 
        Dune::FieldVector<double,1> const trial = 1.0) const override
    {
      return this->evalFotImpl(_cat{}, _traits{}, values[iq], grad_test, trial);
    }
        
    virtual double evalSot(size_t iq, 
        Dune::FieldVector<double,dim> const& grad_test, 
        Dune::FieldVector<double,dim> const& grad_trial) const override
    {
      return this->evalSotImpl(_cat{}, _traits{}, values[iq], grad_test, grad_trial);
    }
                            
    virtual int getDegree() const override
    {
      return term.getDegree();
    }
116

117
118
119
120
121
122
123
124
125
126
127
  private:
    Term term;
    
    using value_type = std::decay_t< decltype( std::declval<Term>()[std::declval<size_t>()] ) >;
    using _cat    = ValueCategory_t<value_type>;
    using _traits = Traits;
    
    std::vector<value_type> values;
  };
  
} // end namespace AMDiS