dune-localfunctions  2.4
pyramidp1localbasis.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_PYRAMID_P1_LOCALBASIS_HH
4 #define DUNE_PYRAMID_P1_LOCALBASIS_HH
5 
6 #include <dune/common/fmatrix.hh>
7 
9 
10 
11 namespace Dune
12 {
23  template<class D, class R>
25  {
26  public:
28  typedef LocalBasisTraits<D,3,Dune::FieldVector<D,3>,R,1,Dune::FieldVector<R,1>,
29  Dune::FieldMatrix<R,1,3> > Traits;
30 
32  unsigned int size () const
33  {
34  return 5;
35  }
36 
38  inline void evaluateFunction (const typename Traits::DomainType& in, // position
39  std::vector<typename Traits::RangeType>& out) const // return value
40  {
41  out.resize(5);
42 
43  if(in[0] > in[1])
44  {
45  out[0] = (1-in[0])*(1-in[1])-in[2]*(1-in[1]);
46  out[1] = in[0]*(1-in[1])-in[2]*in[1];
47  out[2] = (1-in[0])*in[1]-in[2]*in[1];
48  out[3] = in[0]*in[1]+in[2]*in[1];
49  }
50  else
51  {
52  out[0] = (1-in[0])*(1-in[1])-in[2]*(1-in[0]);
53  out[1] = in[0]*(1-in[1])-in[2]*in[0];
54  out[2] = (1-in[0])*in[1]-in[2]*in[0];
55  out[3] = in[0]*in[1]+in[2]*in[0];
56  }
57 
58 
59  out[4] = in[2];
60 
61 
62  }
63 
65  inline void
66  evaluateJacobian (const typename Traits::DomainType& in, // position
67  std::vector<typename Traits::JacobianType>& out) const // return value
68  {
69  out.resize(5);
70 
71  if(in[0] > in[1])
72  {
73  out[0][0][0] = -1 + in[1]; out[0][0][1] = -1 + in[0] + in[2]; out[0][0][2] = -1 + in[1];
74  out[1][0][0] = 1 - in[1]; out[1][0][1] = -in[0] - in[2]; out[1][0][2] = -in[1];
75  out[2][0][0] = -in[1]; out[2][0][1] = 1 - in[0] - in[2]; out[2][0][2] = -in[1];
76  out[3][0][0] = in[1]; out[3][0][1] = in[0]+in[2]; out[3][0][2] = in[1];
77  }
78  else
79  {
80  out[0][0][0] = -1 + in[1] + in[2]; out[0][0][1] = -1 + in[0]; out[0][0][2] = -1 + in[0];
81  out[1][0][0] = 1 - in[1] - in[2]; out[1][0][1] = -in[0]; out[1][0][2] = -in[0];
82  out[2][0][0] = -in[1] - in[2]; out[2][0][1] = 1 - in[0]; out[2][0][2] = -in[0];
83  out[3][0][0] = in[1] + in[2]; out[3][0][1] = in[0]; out[3][0][2] = in[0];
84 
85  }
86 
87  out[4][0][0] = 0; out[4][0][1] = 0; out[4][0][2] = 1;
88  }
89 
91  unsigned int order () const
92  {
93  return 1;
94  }
95  };
96 }
97 #endif
unsigned int size() const
number of shape functions
Definition: pyramidp1localbasis.hh:32
Linear Lagrange shape functions on the pyramid.
Definition: pyramidp1localbasis.hh:24
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: pyramidp1localbasis.hh:38
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
unsigned int order() const
Polynomial order of the shape functions.
Definition: pyramidp1localbasis.hh:91
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, 3 > > Traits
export type traits for function signature
Definition: pyramidp1localbasis.hh:29
D DomainType
domain type
Definition: localbasis.hh:49
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: pyramidp1localbasis.hh:66