dune-localfunctions  2.4.1-rc2
rannacherturek2dlocalbasis.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_RANNACHER_TUREK_2D_LOCALBASIS_HH
4 #define DUNE_RANNACHER_TUREK_2D_LOCALBASIS_HH
5 
6 #include <vector>
7 
8 #include <dune/common/fvector.hh>
9 #include <dune/common/fmatrix.hh>
10 
12 
13 namespace Dune
14 {
15 
16  template< class D, class R >
18  {
20  R, 1, FieldVector< R, 1 >,
21  FieldMatrix< R, 1, 2 > > Traits;
22 
24  unsigned int size () const
25  {
26  return 4;
27  }
28 
30  inline void evaluateFunction ( const typename Traits::DomainType &in,
31  std::vector< typename Traits::RangeType > &out ) const
32  {
33  out.resize(4);
34  typename Traits::DomainFieldType qbase = in[0]*in[0]-in[1]*in[1];
35  out[0] = .75 - 2*in[0] + in[1] + qbase;
36  out[1] = -.25 + in[1] + qbase;
37  out[2] = .75 + in[0] - 2*in[1] - qbase;
38  out[3] = -.25 + in[0] - qbase;
39  }
40 
42  inline void evaluateJacobian ( const typename Traits::DomainType &in,
43  std::vector< typename Traits::JacobianType > &out ) const
44  {
45  out.resize(4);
46 
47  // see http://www.dune-project.org/doc/doxygen/html/classDune_1_1C1LocalBasisInterface.html#d6f8368f8aa43439cc7ef10419f6e2ea
48  // out[i][j][k] = d_k \phi^i_j , where \phi^i_j is the j'th component of the i'th shape function.
49 
50  out[0][0][0] = -2 + 2*in[0]; out[0][0][1] = 1 - 2*in[1];
51  out[1][0][0] = 2*in[0]; out[1][0][1] = 1 - 2*in[1];
52  out[2][0][0] = 1 - 2*in[0]; out[2][0][1] = -2 + 2*in[1];
53  out[3][0][0] = 1 - 2*in[0]; out[3][0][1] = 2*in[1];
54  }
55 
57  unsigned int order () const
58  {
59  // must be 2 here since it contains x^2 and x^2
60  return 2;
61  }
62  };
63 
64 } //namespace Dune
65 
66 #endif // #ifndef DUNE_RANNACHER_TUREK_2D_LOCALBASIS_HH
unsigned int order() const
polynomial order of the shape functions
Definition: rannacherturek2dlocalbasis.hh:57
DF DomainFieldType
Export type for domain field.
Definition: localbasis.hh:40
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
evaluate jacobian of all shape functions
Definition: rannacherturek2dlocalbasis.hh:42
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
evaluate all shape functions
Definition: rannacherturek2dlocalbasis.hh:30
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
LocalBasisTraits< D, 2, FieldVector< D, 2 >, R, 1, FieldVector< R, 1 >, FieldMatrix< R, 1, 2 > > Traits
Definition: rannacherturek2dlocalbasis.hh:21
unsigned int size() const
number of shape functions
Definition: rannacherturek2dlocalbasis.hh:24
Definition: rannacherturek2dlocalbasis.hh:17
D DomainType
domain type
Definition: localbasis.hh:49