dune-pdelab  2.0.0
bdm1cube2dfem.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_FINITEELEMENTMAP_BDM1CUBE2DFEM_HH
3 #define DUNE_PDELAB_FINITEELEMENTMAP_BDM1CUBE2DFEM_HH
4 
5 #include <vector>
6 #include <dune/localfunctions/brezzidouglasmarini/brezzidouglasmarini1cube2d.hh>
7 #include "finiteelementmap.hh"
8 
9 namespace Dune {
10  namespace PDELab {
11 
14  template<typename GV, typename D, typename R>
17  LocalFiniteElementMapTraits<Dune::BDM1Cube2DLocalFiniteElement<D,R> >,
18  BDM1Cube2DLocalFiniteElementMap<GV,D,R> >
19  {
20  typedef Dune::BDM1Cube2DLocalFiniteElement<D,R> FE;
21  typedef typename GV::IndexSet IndexSet;
22 
23  public:
26 
29  : gv(gv_), is(gv_.indexSet()), orient(gv_.size(0))
30  {
31  // create all variants
32  for (int i = 0; i < 16; i++)
33  {
34  variant[i] = FE(i);
35  }
36 
37  // compute orientation for all elements
38  typedef typename GV::Traits::template Codim<0>::Iterator ElementIterator;
39  typedef typename GV::IntersectionIterator IntersectionIterator;
40 
41  // loop once over the grid
42  for (ElementIterator it = gv.template begin<0>(); it != gv.template end<0>(); ++it)
43  {
44  unsigned int myId = is.template index<0>(*it);
45  orient[myId] = 0;
46 
47  IntersectionIterator endit = gv.iend(*it);
48  for (IntersectionIterator iit = gv.ibegin(*it); iit != endit; ++iit)
49  {
50  if (iit->neighbor()
51  && is.template index<0>(*(iit->outside())) > myId)
52  {
53  orient[myId] |= 1 << iit->indexInInside();
54  }
55  }
56  }
57  }
58 
60  template<class EntityType>
61  const typename Traits::FiniteElementType& find(const EntityType& e) const
62  {
63  return variant[orient[is.template index<0>(e)]];
64  }
65 
66  bool fixedSize() const
67  {
68  return true;
69  }
70 
71  std::size_t size(GeometryType gt) const
72  {
73  switch (gt.dim())
74  {
75  case 1:
76  return 2;
77  default:
78  return 0;
79  }
80  }
81 
82  std::size_t maxLocalSize() const
83  {
84  return 8;
85  }
86 
87  private:
88  GV gv;
89  FE variant[16];
90  const IndexSet& is;
91  std::vector<unsigned char> orient;
92  };
93  } // end namespace PDELab
94 } // end namespace Dune
95 
96 #endif // DUNE_PDELAB_FINITEELEMENTMAP_BDM1CUBE2DFEM_HH
const Traits::FiniteElementType & find(const EntityType &e) const
get local basis functions for entity
Definition: bdm1cube2dfem.hh:61
interface for a finite element map
Definition: finiteelementmap.hh:42
std::size_t maxLocalSize() const
Definition: bdm1cube2dfem.hh:82
T FiniteElementType
Type of finite element from local functions.
Definition: finiteelementmap.hh:30
bool fixedSize() const
Definition: bdm1cube2dfem.hh:66
BDM1Cube2DLocalFiniteElementMap(const GV &gv_)
Use when Imp has a standard constructor.
Definition: bdm1cube2dfem.hh:28
Definition: bdm1cube2dfem.hh:15
std::size_t size(GeometryType gt) const
Definition: bdm1cube2dfem.hh:71
LocalFiniteElementMapTraits< FE > Traits
export type of the signature
Definition: bdm1cube2dfem.hh:25
collect types exported by a finite element map
Definition: finiteelementmap.hh:38
const E & e
Definition: interpolate.hh:172