dune-localfunctions  2.2.1
brezzidouglasmarini12dlocalinterpolation.hh
Go to the documentation of this file.
1 #ifndef DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI12DLOCALINTERPOLATION_HH
2 #define DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI12DLOCALINTERPOLATION_HH
3 
4 #include <vector>
5 
6 #include <dune/geometry/quadraturerules.hh>
7 
8 namespace Dune
9 {
18  template<class LB>
20  {
21 
22 public:
25  {
26  sign0 = sign1 = sign2 = 1.0;
27  }
28 
34  BDM12DLocalInterpolation (unsigned int s)
35  {
36  sign0 = sign1 = sign2 = 1.0;
37  if (s & 1)
38  {
39  sign0 = -1.0;
40  }
41  if (s & 2)
42  {
43  sign1 = -1.0;
44  }
45  if (s & 4)
46  {
47  sign2 = -1.0;
48  }
49 
50  n0[0] = 0.0;
51  n0[1] = -1.0;
52  n1[0] = -1.0;
53  n1[1] = 0.0;
54  n2[0] = 1.0/sqrt(2.0);
55  n2[1] = 1.0/sqrt(2.0);
56  c0 = 0.5*n0[0] - 1.0*n0[1];
57  c1 = -1.0*n1[0] + 0.5*n1[1];
58  c2 = 0.5*n2[0] + 0.5*n2[1];
59  }
60 
69  template<typename F, typename C>
70  void interpolate (const F& f, std::vector<C>& out) const
71  {
72  // f gives v*outer normal at a point on the edge!
73  typedef typename LB::Traits::RangeFieldType Scalar;
74  typename F::Traits::RangeType y;
75 
76  out.resize(6);
77  fill(out.begin(), out.end(), 0.0);
78 
79  const int qOrder = 4;
80  const Dune::QuadratureRule<Scalar,1>& rule = Dune::QuadratureRules<Scalar,1>::rule(Dune::GeometryType(Dune::GeometryType::simplex,1), qOrder);
81 
82  for (typename Dune::QuadratureRule<Scalar,1>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
83  {
84  Scalar qPos = it->position();
85  typename LB::Traits::DomainType localPos;
86 
87  localPos[0] = qPos;
88  localPos[1] = 0.0;
89  f.evaluate(localPos, y);
90  out[0] += (y[0]*n0[0] + y[1]*n0[1])*it->weight()*sign0/c0;
91  out[3] += (y[0]*n0[0] + y[1]*n0[1])*(2.0*qPos - 1.0)*it->weight()/c0;
92 
93  localPos[0] = 0.0;
94  localPos[1] = qPos;
95  f.evaluate(localPos, y);
96  out[1] += (y[0]*n1[0] + y[1]*n1[1])*it->weight()*sign1/c1;
97  out[4] += (y[0]*n1[0] + y[1]*n1[1])*(1.0 - 2.0*qPos)*it->weight()/c1;
98 
99  localPos[0] = 1.0 - qPos;
100  localPos[1] = qPos;
101  f.evaluate(localPos, y);
102  out[2] += (y[0]*n2[0] + y[1]*n2[1])*it->weight()*sign2/c2;
103  out[5] += (y[0]*n2[0] + y[1]*n2[1])*(2.0*qPos - 1.0)*it->weight()/c2;
104  }
105  }
106 
107 private:
108  typename LB::Traits::RangeFieldType sign0,sign1,sign2;
109  typename LB::Traits::DomainType n0,n1,n2;
110  typename LB::Traits::RangeFieldType c0,c1,c2;
111  };
112 }
113 
114 #endif // DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI12DLOCALINTERPOLATION_HH