dune-localfunctions  2.7.0
q1.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 
4 #ifndef DUNE_Q1_LOCALFINITEELEMENT_HH
5 #define DUNE_Q1_LOCALFINITEELEMENT_HH
6 
7 #include <dune/geometry/type.hh>
8 
12 
13 namespace Dune
14 {
15 
23  // The test test-q1.cc triggers compiler bugs in gcc-6 and earlier when Q1LocalFiniteElement
24  // is simply redefined as LagrangeCubeLocalFiniteElement: it grabs more and more main memory,
25  // and eventually stalls the machine. Since I didn't find the real cause for this let's just
26  // keep the relevant parts of Q1LocalFiniteElement until we can retire the problematic gcc versions.
27 #if !defined __GNUC__ || __GNUC__ > 6
28  template<class D, class R, int dim>
30 #else
31  template <int dim>
32  class Q1LocalCoefficients
33  {
34  public:
36  Q1LocalCoefficients () : li(1<<dim)
37  {
38  for (std::size_t i=0; i<(1<<dim); i++)
39  li[i] = LocalKey(i,dim,0);
40  }
41 
43  std::size_t size () const
44  {
45  return 1<<dim;
46  }
47 
49  const LocalKey& localKey (std::size_t i) const
50  {
51  return li[i];
52  }
53 
54  private:
55  std::vector<LocalKey> li;
56  };
57 
58  template<class D, class R, int dim>
60  {
61  public:
62  // user-defined default constructor is required for clang 3.8,
63  // see https://gitlab.dune-project.org/core/dune-localfunctions/merge_requests/60
66 
69  typedef LocalFiniteElementTraits<Impl::LagrangeCubeLocalBasis<D,R,dim,1>,Q1LocalCoefficients<dim>,
70  Impl::LagrangeCubeLocalInterpolation<Impl::LagrangeCubeLocalBasis<D,R,dim,1> > > Traits;
71 
74  const typename Traits::LocalBasisType& localBasis () const
75  {
76  return basis;
77  }
78 
81  const typename Traits::LocalCoefficientsType& localCoefficients () const
82  {
83  return coefficients;
84  }
85 
88  const typename Traits::LocalInterpolationType& localInterpolation () const
89  {
90  return interpolation;
91  }
92 
94  unsigned int size () const
95  {
96  return basis.size();
97  }
98 
101  static constexpr GeometryType type ()
102  {
103  return GeometryTypes::cube(dim);
104  }
105 
106  private:
107  Impl::LagrangeCubeLocalBasis<D,R,dim,1> basis;
108  Q1LocalCoefficients<dim> coefficients;
109  Impl::LagrangeCubeLocalInterpolation<Impl::LagrangeCubeLocalBasis<D,R,dim,1> > interpolation;
110  };
111 #endif
112 
114 
119  template<class Geometry, class RF>
122 #if !defined __GNUC__ || __GNUC__ > 6
124  typename Geometry::ctype, RF, Geometry::mydimension, 1
125  >,
126 #else
128  typename Geometry::ctype, RF, Geometry::mydimension
129  >,
130 #endif
131  Geometry
132  >
133  {
134 #if !defined __GNUC__ || __GNUC__ > 6
136  typename Geometry::ctype, RF, Geometry::mydimension, 1
137  > LFE;
138 #else
139  typedef Q1LocalFiniteElement<
140  typename Geometry::ctype, RF, Geometry::mydimension
141  > LFE;
142 #endif
144 
145  static const LFE lfe;
146 
147  public:
150  };
151 
152  template<class Geometry, class RF>
153  const typename Q1FiniteElementFactory<Geometry, RF>::LFE
154  Q1FiniteElementFactory<Geometry, RF>::lfe;
155 }
156 
157 #endif
localfiniteelementtraits.hh
Dune
Definition: bdfmcube.hh:15
Dune::Q1LocalFiniteElement
LagrangeCubeLocalFiniteElement< D, R, dim, 1 > Q1LocalFiniteElement
The local Q1 finite element on cubes.
Definition: q1.hh:29
Dune::LagrangeCubeLocalFiniteElement::localCoefficients
const Traits::LocalCoefficientsType & localCoefficients() const
Returns the assignment of the degrees of freedom to the element subentities.
Definition: lagrangecube.hh:733
Dune::LagrangeCubeLocalFiniteElement
Lagrange finite element for cubes with arbitrary compile-time dimension and polynomial order.
Definition: lagrangecube.hh:708
Dune::LocalFiniteElementTraits::LocalInterpolationType
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:22
Dune::LagrangeCubeLocalFiniteElement::Traits
LocalFiniteElementTraits< Impl::LagrangeCubeLocalBasis< D, R, dim, k >, Impl::LagrangeCubeLocalCoefficients< dim, k >, Impl::LagrangeCubeLocalInterpolation< Impl::LagrangeCubeLocalBasis< D, R, dim, k > > > Traits
Export number types, dimensions, etc.
Definition: lagrangecube.hh:715
Dune::LocalFiniteElementTraits::LocalCoefficientsType
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:18
Dune::LagrangeCubeLocalFiniteElement::localBasis
const Traits::LocalBasisType & localBasis() const
Returns the local basis, i.e., the set of shape functions.
Definition: lagrangecube.hh:726
Dune::ScalarLocalToGlobalFiniteElementAdaptorFactory
Factory for ScalarLocalToGlobalFiniteElementAdaptor objects.
Definition: localtoglobaladaptors.hh:242
Dune::LagrangeCubeLocalFiniteElement::type
static constexpr GeometryType type()
The reference element that the local finite element is defined on.
Definition: lagrangecube.hh:753
lagrangecube.hh
Dune::LagrangeCubeLocalFiniteElement::localInterpolation
const Traits::LocalInterpolationType & localInterpolation() const
Returns object that evaluates degrees of freedom.
Definition: lagrangecube.hh:740
Dune::Q1FiniteElementFactory::Q1FiniteElementFactory
Q1FiniteElementFactory()
default constructor
Definition: q1.hh:149
Dune::LagrangeCubeLocalFiniteElement::size
static constexpr std::size_t size()
The number of shape functions.
Definition: lagrangecube.hh:746
localtoglobaladaptors.hh
Dune::LocalFiniteElementTraits::LocalBasisType
LB LocalBasisType
Definition: localfiniteelementtraits.hh:14
Dune::Q1FiniteElementFactory
Factory for global-valued Q1 elements.
Definition: q1.hh:120
Dune::LocalKey
Describe position of one degree of freedom.
Definition: localkey.hh:20