DUNE PDELab (git)

defaultlocalview.hh
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3
4// SPDX-FileCopyrightText: Copyright © DUNE Project contributors, see file AUTHORS.md
5// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception OR LGPL-3.0-or-later
6
7#ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_DEFAULTLOCALVIEW_HH
8#define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_DEFAULTLOCALVIEW_HH
9
10
11#include <tuple>
12#include <optional>
13
15#include <dune/common/hybridutilities.hh>
17
18#include <dune/functions/common/overflowarray.hh>
19#include <dune/functions/common/multiindex.hh>
20#include <dune/functions/functionspacebases/concepts.hh>
21
22
23
24namespace Dune {
25namespace Functions {
26
27
28
30template<class GB>
32{
33public:
34
36 using GlobalBasis = GB;
37
39 using GridView = typename GlobalBasis::GridView;
40
42 using Element = typename GridView::template Codim<0>::Entity;
43
45 using size_type = std::size_t;
46
48 using Tree = typename GlobalBasis::PreBasis::Node;
49
50protected:
51
52 using PreBasis = typename GlobalBasis::PreBasis;
53
54 // Type used to store the multi indices of the basis vectors.
55 // In contrast to MultiIndex this always has dynamic size.
56 // It's guaranteed, that you can always cast it to MultiIndex
57 using MultiIndexStorage =
58 std::conditional_t<(PreBasis::minMultiIndexSize == PreBasis::maxMultiIndexSize),
61
62public:
63
65 using MultiIndex =
66 std::conditional_t<(PreBasis::minMultiIndexSize == PreBasis::maxMultiIndexSize),
69
70
73 globalBasis_(&globalBasis),
74 tree_(globalBasis_->preBasis().makeNode())
75 {
76 static_assert(models<Concept::BasisTree<GridView>, Tree>(), "Tree type passed to DefaultLocalView does not model the BasisNode concept.");
77 initializeTree(tree_);
78 }
79
85 void bind(const Element& e)
86 {
87 element_ = e;
88 bindTree(tree_, *element_);
89 indices_.resize(size());
90 globalBasis_->preBasis().indices(tree_, indices_.begin());
91 }
92
95 bool bound() const
96 {
97 return static_cast<bool>(element_);
98 }
99
104 const Element& element() const
105 {
106 return *element_;
107 }
108
113 void unbind()
114 {
115 element_.reset();
116 }
117
122 const Tree& tree() const
123 {
124 return tree_;
125 }
126
130 {
131 return tree_.size();
132 }
133
141 {
142 return globalBasis_->preBasis().maxNodeSize();
143 }
144
146 const MultiIndex& index(size_type i) const
147 {
148 return indices_[i];
149 }
150
154 {
155 return *globalBasis_;
156 }
157
158 const DefaultLocalView& rootLocalView() const
159 {
160 return *this;
161 }
162
163protected:
164 const GlobalBasis* globalBasis_;
165 std::optional<Element> element_;
166 Tree tree_;
167 std::vector<MultiIndexStorage> indices_;
168};
169
170
171
172} // end namespace Functions
173} // end namespace Dune
174
175
176
177#endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_DEFAULTLOCALVIEW_HH
The restriction of a finite element basis to a single element.
Definition: defaultlocalview.hh:32
void unbind()
Unbind from the current element.
Definition: defaultlocalview.hh:113
bool bound() const
Return if the view is bound to a grid element.
Definition: defaultlocalview.hh:95
typename GlobalBasis::GridView GridView
The grid view the global FE basis lives on.
Definition: defaultlocalview.hh:39
const 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:146
typename GridView::template Codim< 0 >::Entity Element
Type of the grid element we are bound to.
Definition: defaultlocalview.hh:42
const Tree & tree() const
Return the local ansatz tree associated to the bound entity.
Definition: defaultlocalview.hh:122
void bind(const Element &e)
Bind the view to a grid element.
Definition: defaultlocalview.hh:85
const Element & element() const
Return the grid element that the view is bound to.
Definition: defaultlocalview.hh:104
size_type size() const
Total number of degrees of freedom on this element.
Definition: defaultlocalview.hh:129
GB GlobalBasis
The global FE basis that this is a view on.
Definition: defaultlocalview.hh:36
size_type maxSize() const
Maximum local size for any element on the GridView.
Definition: defaultlocalview.hh:140
std::size_t size_type
The type used for sizes.
Definition: defaultlocalview.hh:45
std::conditional_t<(PreBasis::minMultiIndexSize==PreBasis::maxMultiIndexSize), StaticMultiIndex< size_type, PreBasis::maxMultiIndexSize >, Dune::ReservedVector< size_type, PreBasis::multiIndexBufferSize > > MultiIndex
Type used for global numbering of the basis vectors.
Definition: defaultlocalview.hh:68
typename GlobalBasis::PreBasis::Node Tree
Tree of local finite elements / local shape function sets.
Definition: defaultlocalview.hh:48
DefaultLocalView(const GlobalBasis &globalBasis)
Construct local view for a given global finite element basis.
Definition: defaultlocalview.hh:72
const GlobalBasis & globalBasis() const
Return the global basis that we are a view on.
Definition: defaultlocalview.hh:153
A dynamically sized array-like class with overflow.
Definition: overflowarray.hh:49
A statically sized MultiIndex type.
Definition: multiindex.hh:29
A Vector class with statically reserved memory.
Definition: reservedvector.hh:47
Infrastructure for concepts.
Dune namespace.
Definition: alignedallocator.hh:13
An stl-compliant random-access container which stores everything on the stack.
Static tag representing a codimension.
Definition: dimension.hh:24
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 13, 23:29, 2024)