dune-localfunctions  2.4
dualq1.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 #ifndef DUNE_DUAL_Q1_LOCALFINITEELEMENT_HH
4 #define DUNE_DUAL_Q1_LOCALFINITEELEMENT_HH
5 
6 #include <array>
7 
8 #include <dune/common/fvector.hh>
9 #include <dune/common/fmatrix.hh>
10 
11 #include <dune/geometry/type.hh>
12 #include <dune/geometry/quadraturerules.hh>
13 
18 
19 namespace Dune
20 {
21 
27  template<class D, class R, int dim>
29  {
30  public:
35 
39  {
40  gt.makeCube(dim);
41 
42  // dual basis functions are linear combinations of Lagrange elements
43  // compute these coefficients here because the basis and the local interpolation needs them
44 
45  const Dune::QuadratureRule<D,dim>& quad = Dune::QuadratureRules<D,dim>::rule(gt, 2*dim);
46 
47  const int size = 1 <<dim;
48 
49  // assemble mass matrix on the reference element
50  Dune::FieldMatrix<R, size, size> massMat;
51  massMat = 0;
52 
53  // and the integrals of the lagrange shape functions
54  std::vector<Dune::FieldVector<R,1> > integral(size);
55  for (int i=0; i<size; i++)
56  integral[i] = 0;
57 
58  for(size_t pt=0; pt<quad.size(); pt++) {
59 
60  const Dune::FieldVector<D ,dim>& pos = quad[pt].position();
61  std::vector<typename Traits::LocalBasisType::Traits::RangeType> q1Values(size);
62 
63  // evaluate q1 basis functions
64  for (int i=0; i<size; i++) {
65 
66  q1Values[i] = 1;
67 
68  for (int j=0; j<dim; j++)
69  // if j-th bit of i is set multiply with in[j], else with 1-in[j]
70  q1Values[i] *= (i & (1<<j)) ? pos[j] : 1-pos[j];
71  }
72 
73  double weight = quad[pt].weight();
74 
75  for (int k=0; k<size; k++) {
76  integral[k] += q1Values[k]*weight;
77 
78  for (int l=0; l<=k; l++)
79  massMat[k][l] += weight*(q1Values[k]*q1Values[l]);
80  }
81  }
82 
83  // make matrix symmetric
84  for (int i=0; i<size-1; i++)
85  for (int j=i+1; j<size; j++)
86  massMat[i][j] = massMat[j][i];
87 
88  //solve for the coefficients
89  std::array<Dune::FieldVector<R, size>, size> coefficients;
90 
91  for (int i=0; i<size; i++) {
92 
93  Dune::FieldVector<R, size> rhs(0);
94  rhs[i] = integral[i];
95 
96  coefficients[i] = 0;
97  massMat.solve(coefficients[i] ,rhs);
98 
99  }
100 
101  basis.setCoefficients(coefficients);
102  interpolation.setCoefficients(coefficients);
103  }
104 
107  const typename Traits::LocalBasisType& localBasis () const
108  {
109  return basis;
110  }
111 
115  {
116  return coefficients;
117  }
118 
122  {
123  return interpolation;
124  }
125 
127  unsigned int size () const
128  {
129  return basis.size();
130  }
131 
134  GeometryType type () const
135  {
136  return gt;
137  }
138 
140  {
141  return new DualQ1LocalFiniteElement(*this);
142  }
143 
144  private:
146  DualQ1LocalCoefficients<dim> coefficients;
148  GeometryType gt;
149  };
150 }
151 
152 #endif
LocalFiniteElementTraits< DualQ1LocalBasis< D, R, dim >, DualQ1LocalCoefficients< dim >, DualQ1LocalInterpolation< dim, DualQ1LocalBasis< D, R, dim > > > Traits
Definition: dualq1.hh:34
unsigned int size() const
Number of shape functions in this finite element.
Definition: dualq1.hh:127
DualQ1LocalFiniteElement * clone() const
Definition: dualq1.hh:139
Definition: dualq1localinterpolation.hh:17
GeometryType type() const
Definition: dualq1.hh:134
DualQ1LocalFiniteElement()
Definition: dualq1.hh:38
The local dual Q1 finite element on cubes.
Definition: dualq1.hh:28
const Traits::LocalBasisType & localBasis() const
Definition: dualq1.hh:107
traits helper struct
Definition: localfiniteelementtraits.hh:10
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:22
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
Layout map for dual Q1 elements.
Definition: dualq1localcoefficients.hh:22
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:18
const Traits::LocalInterpolationType & localInterpolation() const
Definition: dualq1.hh:121
LB LocalBasisType
Definition: localfiniteelementtraits.hh:14
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: dualq1.hh:114
Dual Lagrange shape functions of order 1 on the reference cube.
Definition: dualq1localbasis.hh:25