3#ifndef DUNE_PDELAB_GRIDFUNCTIONSPACE_DUNEFUNCTIONSLOCALFUNCTIONSPACE_HH
4#define DUNE_PDELAB_GRIDFUNCTIONSPACE_DUNEFUNCTIONSLOCALFUNCTIONSPACE_HH
11#include <dune/geometry/referenceelements.hh>
13#include <dune/localfunctions/common/interfaceswitch.hh>
14#include <dune/localfunctions/common/localkey.hh>
16#include <dune/typetree/typetree.hh>
18#include <dune/pdelab/gridfunctionspace/tags.hh>
28 namespace Experimental {
30 template<
typename LFS>
32 :
public TypeTree::LeafNode
35 const auto& finiteElement()
const
37 return static_cast<const LFS*
>(
this)->tree().finiteElement();
40 template<
typename Tree>
43 using FiniteElement =
typename Tree::FiniteElement;
44 using FiniteElementType = FiniteElement;
48 template<
typename GFS,
typename TreePath = TypeTree::Hybr
idTreePath<>>
49 class LocalFunctionSpace
50 :
public LeafLFSMixin<LocalFunctionSpace<GFS,TreePath>>
55 using Basis =
typename GFS::Basis;
56 using LocalView =
typename Basis::LocalView;
57 using Tree = TypeTree::ChildForTreePath<typename LocalView::Tree,TreePath>;
58 using DOFIndex =
typename GFS::Ordering::Traits::DOFIndex;
60 template<
typename LFS,
typename C,
typename Tag,
bool fast>
61 friend class Dune::PDELab::LFSIndexCacheBase;
64 :
public LeafLFSMixin<LocalFunctionSpace<GFS,TreePath>>::template Traits<Tree>
67 using GridFunctionSpace = GFS;
68 using GridView =
typename GFS::Traits::GridView;
69 using SizeType = std::size_t;
70 using DOFIndex =
typename GFS::Ordering::Traits::DOFIndex;
71 using ConstraintsType =
typename GFS::Traits::ConstraintsType;
75 using size_type = std::size_t;
77 LocalFunctionSpace(std::shared_ptr<const GFS> gfs, TreePath tree_path = TreePath(), size_type offset = 0)
79 , _local_view(gfs->basis().localView())
80 , _tree_path(tree_path)
81 , _tree(TypeTree::
child(_local_view.tree(),tree_path))
84 size_type subSpaceDepth()
const
90 size_type
size ()
const
92 return _local_view.size();
95 size_type maxSize ()
const
98 return _local_view.maxSize();
102 size_type localIndex (size_type index)
const
104 return _tree.localIndex(index);
108 DOFIndex dofIndex(size_type index)
const
112 auto localKey = _local_view.tree().finiteElement().localCoefficients().localKey(index);
114 const auto& indexSet = _gfs->basis().gridView().indexSet();
117 auto gt = refElement.type(localKey.subEntity(), localKey.codim());
120 auto indexOnEntity = indexSet.subIndex(_local_view.element(),
121 localKey.subEntity(),
126 GFS::Ordering::Traits::DOFIndexAccessor::store(result,
gt,indexOnEntity,localKey.index());
131 auto containerIndex(size_type index)
const
133 MultiIndex<std::size_t,1> result;
134 result.set({_local_view.index(_tree.localIndex(index))});
139 const GFS& gridFunctionSpace()
const
144 void bind(
const typename GFS::Traits::EntitySet::template Codim<0>::Entity& e)
149 const typename Traits::ConstraintsType&
constraints()
const
151 return _gfs->constraints();
154 const Tree& tree()
const
161 typename GFS::Ordering::Traits::ContainerIndex containerIndex(
const DOFIndex& i)
const
163 return _gfs->ordering().containerIndex(i);
166 std::shared_ptr<const GFS> _gfs;
167 LocalView _local_view;
174 template<
typename DFBasis,
typename V,
typename CE=NoConstra
ints>
175 class GridFunctionSpace;
181 template<
typename DFBasis,
typename V,
typename CE,
typename TAG>
182 class LocalFunctionSpace<Experimental::GridFunctionSpace<DFBasis,V,CE>,TAG>
183 :
public Experimental::LocalFunctionSpace<Experimental::GridFunctionSpace<DFBasis,V,CE>>
186 using GFS = Experimental::GridFunctionSpace<DFBasis,V,CE>;
190 LocalFunctionSpace(std::shared_ptr<const GFS> gfs)
191 : Experimental::LocalFunctionSpace<GFS>(gfs)
194 LocalFunctionSpace(
const GFS& gfs)
200 template<
typename DFBasis,
typename V,
typename CE>
201 class LocalFunctionSpace<Experimental::GridFunctionSpace<DFBasis,V,CE>,AnySpaceTag>
202 :
public Experimental::LocalFunctionSpace<Experimental::GridFunctionSpace<DFBasis,V,CE>>
205 using GFS = Experimental::GridFunctionSpace<DFBasis,V,CE>;
209 LocalFunctionSpace(std::shared_ptr<const GFS> gfs)
210 : Experimental::LocalFunctionSpace<GFS>(gfs)
213 LocalFunctionSpace(
const GFS& gfs)
bool gt(const T &first, const T &second, typename EpsilonType< T >::Type epsilon)
test if first greater than second
Definition: float_cmp.cc:158
void constraints(const GFS &gfs, CG &cg, const bool verbose=false)
construct constraints
Definition: constraints.hh:749
Dune namespace.
Definition: alignedallocator.hh:13
std::shared_ptr< T > stackobject_to_shared_ptr(T &t)
Create a shared_ptr for a stack-allocated object.
Definition: shared_ptr.hh:72
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75
This file implements several utilities related to std::shared_ptr.
Standard Dune debug streams.
static const ReferenceElement & general(const GeometryType &type)
get general reference elements
Definition: referenceelements.hh:156