4 #ifndef DUNE_PDELAB_FINITEELEMENTMAP_FINITEELEMENTMAP_HH
5 #define DUNE_PDELAB_FINITEELEMENTMAP_FINITEELEMENTMAP_HH
7 #include <dune/common/deprecated.hh>
10 #include <dune/geometry/referenceelements.hh>
41 template<
class T,
class Imp>
54 template<
class EntityType>
55 const typename Traits::FiniteElementType&
56 find (
const EntityType&
e)
const =
delete;
82 static constexpr std::size_t
size(GeometryType gt) =
delete;
87 static constexpr
bool hasDOFs(
int codim) =
delete;
97 template<
class Imp,
int dim_>
100 SimpleLocalFiniteElementMap<Imp,dim_> >
115 template<
class EntityType>
147 template<
typename GV,
typename FE,
typename Imp>
150 LocalFiniteElementMapTraits<FE>,
154 typedef typename GV::IndexSet IndexSet;
155 static const int dim = GV::dimension;
166 : gv(gv_), orient(gv_.
size(0))
168 using ct =
typename GV::Grid::ctype;
170 ReferenceElements<ct, dim>::general(FE().type());
172 auto &idSet = gv.grid().globalIdSet();
175 variant.resize(1 << refElem.size(dim-1));
176 for (std::size_t i=0; i<variant.size(); i++)
180 auto& indexSet = gv.indexSet();
183 for (
const auto& element : elements(gv))
185 auto elemid = indexSet.index(element);
188 std::vector<typename GV::Grid::GlobalIdSet::IdType> vid(refElem.size(dim));
189 for(std::size_t i = 0; i < vid.size(); ++i)
190 vid[i] = idSet.subId(element, i, dim);
193 for(
int i = 0; i < refElem.size(dim-1); ++i) {
194 auto v0 = refElem.subEntity(i, dim-1, 0, dim);
195 auto v1 = refElem.subEntity(i, dim-1, 1, dim);
197 if((v0 > v1) != (vid[v0] > vid[v1]))
198 orient[elemid] |= 1 << i;
204 template<
class EntityType>
207 return variant[orient[gv.indexSet().index(
e)]];
212 std::vector<FE> variant;
213 std::vector<unsigned char> orient;
218 template<
typename GV,
typename FE,
typename Imp, std::
size_t Variants>
221 LocalFiniteElementMapTraits<FE>,
224 typedef FE FiniteElement;
225 typedef typename GV::IndexSet IndexSet;
236 : gv(gv_), is(gv_.indexSet()), orient(gv_.
size(0))
239 for (std::size_t i = 0; i < Variants; i++)
241 variant[i] = FiniteElement(i);
246 for(
const auto& cell : elements(gv))
248 auto myId = is.index(cell);
251 for (
const auto& intersection : intersections(gv,cell))
253 if (intersection.neighbor()
254 && is.index(intersection.outside()) > myId)
256 orient[myId] |= 1 << intersection.indexInInside();
263 template<
class EntityType>
266 return variant[orient[is.index(
e)]];
271 std::array<FiniteElement,Variants> variant;
273 std::vector<unsigned char> orient;
PDELab-specific exceptions.
const Entity & e
Definition: localfunctionspace.hh:123
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Base class for all PDELab exceptions.
Definition: exceptions.hh:19
general FiniteElementMap exception
Definition: finiteelementmap.hh:19
FiniteElementMap exception concerning the computation of the FiniteElementMap size.
Definition: finiteelementmap.hh:21
FiniteElementMap exception raised when trying to obtain a finite element for an unsupported GeometryT...
Definition: finiteelementmap.hh:23
collect types exported by a finite element map
Definition: finiteelementmap.hh:28
T FiniteElement
Type of finite element from local functions.
Definition: finiteelementmap.hh:33
T FiniteElementType
Type of finite element from local functions.
Definition: finiteelementmap.hh:30
collect types exported by a finite element map
Definition: finiteelementmap.hh:38
interface for a finite element map
Definition: finiteelementmap.hh:43
static constexpr bool hasDOFs(int codim)=delete
return if FiniteElementMap has degrees of freedom for given codimension
static constexpr bool fixedSize()=delete
a FiniteElementMap is fixedSize iif the size of the local functions space for each GeometryType is fi...
static constexpr std::size_t maxLocalSize()=delete
compute an upper bound for the local number of DOFs.
T Traits
Export traits.
Definition: finiteelementmap.hh:46
static constexpr int dimension
dimension of the domain this FEM is defined on.
Definition: finiteelementmap.hh:74
static constexpr std::size_t size(GeometryType gt)=delete
if the FiniteElementMap is fixedSize, the size methods computes the number of DOFs for given Geometry...
const Traits::FiniteElementType & find(const EntityType &e) const =delete
Return local basis for the given entity.
simple implementation where all entities have the same finite element
Definition: finiteelementmap.hh:101
SimpleLocalFiniteElementMap(const Imp &imp_)
Constructor where an instance of Imp can be provided.
Definition: finiteelementmap.hh:111
LocalFiniteElementMapTraits< Imp > Traits
export type of the signature
Definition: finiteelementmap.hh:104
const Traits::FiniteElementType & find(const EntityType &e) const
get local basis functions for entity
Definition: finiteelementmap.hh:116
static constexpr int dimension
Definition: finiteelementmap.hh:121
SimpleLocalFiniteElementMap()
Use when Imp has a standard constructor.
Definition: finiteelementmap.hh:107
implementation for finite elements requiring oriented edges
Definition: finiteelementmap.hh:153
const Traits::FiniteElementType & find(const EntityType &e) const
get local basis functions for entity
Definition: finiteelementmap.hh:205
EdgeS0LocalFiniteElementMap(const GV &gv_)
construct EdgeSLocalFiniteElementMap
Definition: finiteelementmap.hh:165
static constexpr int dimension
The dimension of the finite elements returned by this map.
Definition: finiteelementmap.hh:162
LocalFiniteElementMapTraits< FE > Traits
export type of the signature
Definition: finiteelementmap.hh:159
Definition: finiteelementmap.hh:223
static constexpr int dimension
The dimension of the finite elements returned by this map.
Definition: finiteelementmap.hh:232
const Traits::FiniteElementType & find(const EntityType &e) const
get local basis functions for entity
Definition: finiteelementmap.hh:264
RTLocalFiniteElementMap(const GV &gv_)
Use when Imp has a standard constructor.
Definition: finiteelementmap.hh:235
LocalFiniteElementMapTraits< FE > Traits
export type of the signature
Definition: finiteelementmap.hh:229