dune-localfunctions  2.2.1
raviartthomas0q2dall.hh
Go to the documentation of this file.
1 #ifndef DUNE_RT0Q2DALL_HH
2 #define DUNE_RT0Q2DALL_HH
3 
4 #include <cstddef>
5 #include <vector>
6 
7 #include <dune/common/fmatrix.hh>
8 
11 
12 namespace Dune
13 {
22  template<class D, class R>
24  {
25  public:
26  typedef LocalBasisTraits<D,2,Dune::FieldVector<D,2>,R,2,Dune::FieldVector<R,2>,
27  Dune::FieldMatrix<R,2,2> > Traits;
28 
31  {
32  sign0 = sign1 = sign2 = sign3 = 1.0;
33  }
34 
36  RT0Q2DLocalBasis (unsigned int s)
37  {
38  sign0 = sign1 = sign2 = sign3 = 1.0;
39  if (s&1) sign0 = -1.0;
40  if (s&2) sign1 = -1.0;
41  if (s&4) sign2 = -1.0;
42  if (s&8) sign3 = -1.0;
43  }
44 
46  unsigned int size () const
47  {
48  return 4;
49  }
50 
52  inline void evaluateFunction (const typename Traits::DomainType& in,
53  std::vector<typename Traits::RangeType>& out) const
54  {
55  out.resize(4);
56  out[0][0] = sign0*(in[0]-1.0); out[0][1]=0.0;
57  out[1][0] = sign1*(in[0]); out[1][1]=0.0;
58  out[2][0] = 0.0; out[2][1]=sign2*(in[1]-1.0);
59  out[3][0] = 0.0; out[3][1]=sign3*(in[1]);
60  }
61 
63  inline void
64  evaluateJacobian (const typename Traits::DomainType& in, // position
65  std::vector<typename Traits::JacobianType>& out) const // return value
66  {
67  out.resize(4);
68  out[0][0][0] = sign0; out[0][0][1] = 0;
69  out[0][1][0] = 0; out[0][1][1] = 0;
70 
71  out[1][0][0] = sign1; out[1][0][1] = 0;
72  out[1][1][0] = 0; out[1][1][1] = 0;
73 
74  out[2][0][0] = 0; out[2][0][1] = 0;
75  out[2][1][0] = 0; out[2][1][1] = sign2;
76 
77  out[3][0][0] = 0; out[3][0][1] = 0;
78  out[3][1][0] = 0; out[3][1][1] = sign3;
79  }
80 
82  unsigned int order () const
83  {
84  return 1;
85  }
86 
87  private:
88  R sign0, sign1, sign2, sign3;
89  };
90 
91 
99  template<class LB>
101  {
102  public:
103 
106  {
107  sign0 = sign1 = sign2 = sign3 = 1.0;
108  }
109 
111  RT0Q2DLocalInterpolation (unsigned int s)
112  {
113  sign0 = sign1 = sign2 = sign3 = 1.0;
114  if (s&1) sign0 *= -1.0;
115  if (s&2) sign1 *= -1.0;
116  if (s&4) sign2 *= -1.0;
117  if (s&8) sign3 *= -1.0;
118 
119  m0[0] = 0.0; m0[1] = 0.5;
120  m1[0] = 1.0; m1[1] = 0.5;
121  m2[0] = 0.5; m2[1] = 0.0;
122  m3[0] = 0.5; m3[1] = 1.0;
123 
124  n0[0] = -1.0; n0[1] = 0.0;
125  n1[0] = 1.0; n1[1] = 0.0;
126  n2[0] = 0.0; n2[1] = -1.0;
127  n3[0] = 0.0; n3[1] = 1.0;
128  }
129 
130  template<typename F, typename C>
131  void interpolate (const F& f, std::vector<C>& out) const
132  {
133  // f gives v*outer normal at a point on the edge!
134  typename F::Traits::RangeType y;
135 
136  out.resize(4);
137 
138  f.evaluate(m0,y); out[0] = (y[0]*n0[0]+y[1]*n0[1])*sign0;
139  f.evaluate(m1,y); out[1] = (y[0]*n1[0]+y[1]*n1[1])*sign1;
140  f.evaluate(m2,y); out[2] = (y[0]*n2[0]+y[1]*n2[1])*sign2;
141  f.evaluate(m3,y); out[3] = (y[0]*n3[0]+y[1]*n3[1])*sign3;
142  }
143 
144  private:
145  typename LB::Traits::RangeFieldType sign0,sign1,sign2,sign3;
146  typename LB::Traits::DomainType m0,m1,m2,m3;
147  typename LB::Traits::DomainType n0,n1,n2,n3;
148  };
149 
157  {
158  public:
161  {
162  for (std::size_t i=0; i<4; i++)
163  li[i] = LocalKey(i,1,0);
164  }
165 
167  std::size_t size () const
168  {
169  return 4;
170  }
171 
173  const LocalKey& localKey (std::size_t i) const
174  {
175  return li[i];
176  }
177 
178  private:
179  std::vector<LocalKey> li;
180  };
181 
182 }
183 #endif