LocalView.hpp 4.39 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
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
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
#pragma once

#include <tuple>
#include <optional>
#include <vector>

#include <dune/common/concept.hh>
#include <dune/functions/functionspacebases/concepts.hh>

#include <amdis/functions/NodeCache.hpp>
#include <amdis/functions/Nodes.hpp>

// NOTE: This is a variant of dune-functions DefaultLocalView

namespace AMDiS
{
  /// \brief The restriction of a finite element basis to a single element
  template <class GB>
  class LocalView
  {
    using PrefixPath = Dune::TypeTree::HybridTreePath<>;

    // Node index set provided by PreBasis
    using NodeIndexSet = NodeIndexSet_t<typename GB::PreBasis, PrefixPath>;

  public:
    /// The global FE basis that this is a view on
    using GlobalBasis = GB;

    /// The grid view the global FE basis lives on
    using GridView = typename GlobalBasis::GridView;

    /// Type of the grid element we are bound to
    using Element = typename GridView::template Codim<0>::Entity;

    /// The type used for sizes
    using size_type = std::size_t;

    /// Tree of local finite elements / local shape function sets
    using Tree = Node_t<typename GlobalBasis::PreBasis, PrefixPath>;

    /// Cached basis-tree
    using TreeCache = NodeCache_t<Tree>;

    /// Type used for global numbering of the basis vectors
    using MultiIndex = typename NodeIndexSet::MultiIndex;

  private:
    template <class NIS>
    using hasIndices = decltype(std::declval<NIS>().indices(std::declval<std::vector<typename NIS::MultiIndex>>().begin()));

  public:
    /// \brief Construct local view for a given global finite element basis
    LocalView (GlobalBasis const& globalBasis)
      : globalBasis_(&globalBasis)
      , tree_(makeNode(globalBasis_->preBasis(), PrefixPath{}))
      , treeCache_(makeNodeCache(tree_))
      , nodeIndexSet_(makeNodeIndexSet(globalBasis_->preBasis(), PrefixPath{}))
    {
      static_assert(Dune::models<Dune::Functions::Concept::BasisTree<GridView>, Tree>());
      initializeTree(tree_);
    }

    /// \brief Bind the view to a grid element
    /**
     * Having to bind the view to an element before being able to actually access any of
     * its data members offers to centralize some expensive setup code in the `bind`
     * method, which can save a lot of run-time.
     */
    void bind (Element const& element)
    {
      element_ = element;
      bindTree(tree_, *element_);
      nodeIndexSet_.bind(tree_);
      indices_.resize(size());

      if constexpr (Dune::Std::is_detected<hasIndices,NodeIndexSet>{})
        nodeIndexSet_.indices(indices_.begin());
      else
        for (size_type i = 0; i < size(); ++i)
          indices_[i] = nodeIndexSet_.index(i);
    }

    /// \brief Return if the view is bound to a grid element
    bool isBound () const
    {
      return bool(element_);
    }

    /// \brief Return the grid element that the view is bound to
    Element const& element () const
    {
      return *element_;
    }

    /// \brief Unbind from the current element
    /**
     * Calling this method should only be a hint that the view can be unbound.
     */
    void unbind ()
    {
      nodeIndexSet_.unbind();
      element_.reset();
    }

    /// \brief Return the local ansatz tree associated to the bound entity
    Tree const& tree () const
    {
      return tree_;
    }

    /// \brief Cached version of the local ansatz tree
    TreeCache const& treeCache () const
    {
      return treeCache_;
    }

    /// \brief Total number of degrees of freedom on this element
    size_type size () const
    {
      return tree_.size();
    }

    /// \brief Maximum local size for any element on the GridView
    /**
     * This is the maximal size needed for local matrices and local vectors.
     */
    size_type maxSize () const
    {
      return globalBasis_->preBasis().maxNodeSize();
    }

    /// \brief Maps from subtree index set [0..size-1] to a globally unique multi index
    /// in global basis
    MultiIndex index (size_type i) const
    {
      return indices_[i];
    }

    /// \brief Return the global basis that we are a view on
    GlobalBasis const& globalBasis () const
    {
      return *globalBasis_;
    }

    /// \brief Return this local-view
    LocalView const& rootLocalView () const
    {
      return *this;
    }

  protected:
    GlobalBasis const* globalBasis_;
    std::optional<Element> element_;
    Tree tree_;
    TreeCache treeCache_;
    NodeIndexSet nodeIndexSet_;
    std::vector<MultiIndex> indices_;
  };

} // end namespace AMDiS