3 #ifndef DUNE_DUAL_Q1_LOCALFINITEELEMENT_HH
4 #define DUNE_DUAL_Q1_LOCALFINITEELEMENT_HH
8 #include <dune/common/fvector.hh>
9 #include <dune/common/fmatrix.hh>
11 #include <dune/geometry/type.hh>
12 #include <dune/geometry/quadraturerules.hh>
27 template<
class D,
class R,
int dim>
45 const Dune::QuadratureRule<D,dim>& quad = Dune::QuadratureRules<D,dim>::rule(gt, 2*dim);
47 const int size = 1 <<dim;
50 Dune::FieldMatrix<R, size, size> massMat;
54 std::vector<Dune::FieldVector<R,1> > integral(size);
55 for (
int i=0; i<
size; i++)
58 for(
size_t pt=0; pt<quad.size(); pt++) {
60 const Dune::FieldVector<D ,dim>& pos = quad[pt].position();
61 std::vector<typename Traits::LocalBasisType::Traits::RangeType> q1Values(size);
64 for (
int i=0; i<
size; i++) {
68 for (
int j=0; j<dim; j++)
70 q1Values[i] *= (i & (1<<j)) ? pos[j] : 1-pos[j];
73 double weight = quad[pt].weight();
75 for (
int k=0; k<
size; k++) {
76 integral[k] += q1Values[k]*weight;
78 for (
int l=0; l<=k; l++)
79 massMat[k][l] += weight*(q1Values[k]*q1Values[l]);
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];
89 std::array<Dune::FieldVector<R, size>, size> coefficients;
91 for (
int i=0; i<
size; i++) {
93 Dune::FieldVector<R, size> rhs(0);
97 massMat.solve(coefficients[i] ,rhs);
101 basis.setCoefficients(coefficients);
102 interpolation.setCoefficients(coefficients);
123 return interpolation;
146 DualQ1LocalCoefficients<dim> coefficients;
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