3#ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_DEFAULTLOCALVIEW_HH
4#define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_DEFAULTLOCALVIEW_HH
10#include <dune/common/concept.hh>
11#include <dune/common/hybridutilities.hh>
25template<
typename PreBasis>
26using PreBasisHasIndicesConcept =
decltype(std::declval<PreBasis>().indices(std::declval<typename PreBasis::Node>(), std::declval<std::vector<typename PreBasis::MultiIndex>>().begin()),
true);
31template<
typename PreBasis>
32using preBasisHasIndices = Std::is_detected<PreBasisHasIndicesConcept, PreBasis>;
43template<
typename PreBasis,
typename Node,
typename It>
44It preBasisIndices(
const PreBasis& preBasis,
const Node& node, It multiIndices)
46 if constexpr (preBasisHasIndices<PreBasis>{})
47 return preBasis.indices(node, multiIndices);
50 auto indexSet = preBasis().makeIndexSet();
52 return indexSet.indices(multiIndices);
64class DefaultNodeIndexSet
68 using size_type = std::size_t;
70 using MultiIndex =
typename PreBasis::MultiIndex;
71 using Node =
typename PreBasis::Node;
73 DefaultNodeIndexSet(
const PreBasis& preBasis) :
78 void bind(
const Node& node) {
86 size_type size()
const {
87 assert(node_ !=
nullptr);
92 It indices(It it)
const {
93 assert(node_ !=
nullptr);
94 return preBasis_->indices(*node_, it);
98 const PreBasis* preBasis_;
119 using Element =
typename GridView::template Codim<0>::Entity;
125 using Tree =
typename GlobalBasis::PreBasis::Node;
128 using MultiIndex =
typename GlobalBasis::PreBasis::MultiIndex;
133 template<
typename NodeIndexSet_>
134 using hasIndices =
decltype(std::declval<NodeIndexSet_>().indices(std::declval<std::vector<typename NodeIndexSet_::MultiIndex>>().begin()));
138 struct DummyNodeIndexSet {};
140 static auto makeIndexSet(
const typename GlobalBasis::PreBasis& preBasis)
142 if constexpr (Impl::preBasisHasIndices<typename GlobalBasis::PreBasis>{})
143 return DummyNodeIndexSet{};
145 return preBasis.makeIndexSet();
148 using NodeIndexSet = std::decay_t<decltype(makeIndexSet(std::declval<typename GlobalBasis::PreBasis>()))>;
158 static_assert(models<Concept::BasisTree<GridView>,
Tree>(),
"Tree type passed to DefaultLocalView does not model the BasisNode concept.");
172 if constexpr (Impl::preBasisHasIndices<typename GlobalBasis::PreBasis>{})
177 if constexpr (Std::is_detected<hasIndices,NodeIndexSet>{})
206 if constexpr (not Impl::preBasisHasIndices<typename GlobalBasis::PreBasis>{})
Definition: polynomial.hh:10
void bindTree(Tree &tree, const Entity &entity, std::size_t offset=0)
Definition: nodes.hh:251
void initializeTree(Tree &tree, std::size_t treeIndexOffset=0)
Definition: nodes.hh:258
The restriction of a finite element basis to a single element.
Definition: defaultlocalview.hh:109
void unbind()
Unbind from the current element.
Definition: defaultlocalview.hh:204
typename GlobalBasis::GridView GridView
The grid view the global FE basis lives on.
Definition: defaultlocalview.hh:116
std::optional< Element > element_
Definition: defaultlocalview.hh:258
typename GridView::template Codim< 0 >::Entity Element
Type of the grid element we are bound to.
Definition: defaultlocalview.hh:119
const Tree & tree() const
Return the local ansatz tree associated to the bound entity.
Definition: defaultlocalview.hh:215
void bind(const Element &e)
Bind the view to a grid element.
Definition: defaultlocalview.hh:167
const Element & element() const
Return the grid element that the view is bound to.
Definition: defaultlocalview.hh:195
MultiIndex index(size_type i) const
Maps from subtree index set [0..size-1] to a globally unique multi index in global basis.
Definition: defaultlocalview.hh:239
size_type size() const
Total number of degrees of freedom on this element.
Definition: defaultlocalview.hh:222
GB GlobalBasis
The global FE basis that this is a view on.
Definition: defaultlocalview.hh:113
Tree tree_
Definition: defaultlocalview.hh:259
size_type maxSize() const
Maximum local size for any element on the GridView.
Definition: defaultlocalview.hh:233
std::size_t size_type
The type used for sizes.
Definition: defaultlocalview.hh:122
bool isBound() const
Return if the view is bound to a grid element.
Definition: defaultlocalview.hh:187
const DefaultLocalView & rootLocalView() const
Definition: defaultlocalview.hh:251
NodeIndexSet nodeIndexSet_
Definition: defaultlocalview.hh:260
typename GlobalBasis::PreBasis::MultiIndex MultiIndex
Type used for global numbering of the basis vectors.
Definition: defaultlocalview.hh:128
typename GlobalBasis::PreBasis::Node Tree
Tree of local finite elements / local shape function sets.
Definition: defaultlocalview.hh:125
std::vector< MultiIndex > indices_
Definition: defaultlocalview.hh:261
DefaultLocalView(const GlobalBasis &globalBasis)
Construct local view for a given global finite element basis.
Definition: defaultlocalview.hh:153
const GlobalBasis * globalBasis_
Definition: defaultlocalview.hh:257
const GlobalBasis & globalBasis() const
Return the global basis that we are a view on.
Definition: defaultlocalview.hh:246