dune-localfunctions  2.2.1
hierarchicalprismp2localbasis.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_HIERARCHICAL_PRISM_P2_LOCALBASIS_HH
3 #define DUNE_HIERARCHICAL_PRISM_P2_LOCALBASIS_HH
4 
9 #include <dune/common/fvector.hh>
10 #include <dune/common/fmatrix.hh>
11 
13 
14 namespace Dune
15 {
16  template<class D, class R>
18  {
19  public:
21  typedef LocalBasisTraits<D,3,Dune::FieldVector<D,3>,R,1,Dune::FieldVector<R,1>, Dune::FieldMatrix<R,1,3> > Traits;
22 
24 unsigned int size () const
25 {
26 return 18;
27 }
28 
30 void evaluateFunction (const typename Traits::DomainType& in,
31  std::vector<typename Traits::RangeType> & out) const
32 {
33 out.resize(18);
34 
35  out[0]=(1.0-in[0]-in[1])*(1.0-in[2]);
36  out[1]= in[0]*(1-in[2]);
37  out[2]=in[1]*(1-in[2]);
38  out[3]=in[2]*(1.0-in[0]-in[1]);
39  out[4]=in[0]*in[2];
40  out[5]=in[1]*in[2];
41 
42  //edges
43  out[6]=2*(1.0-in[0]-in[1])*(0.5-in[0]-in[1])*(4*in[2]-4*in[2]*in[2]);
44  out[7]=2*in[0]*(-0.5+in[0])*(4*in[2]-4*in[2]*in[2]);
45  out[8]=2*in[1]*(-0.5+in[1])*(4*in[2]-4*in[2]*in[2]);
46  out[9]=4*in[0]*(1-in[0]-in[1])*(1-3*in[2]+2*in[2]*in[2]);
47  out[10]=4*in[1]*(1-in[0]-in[1])*(1-3*in[2]+2*in[2]*in[2]);
48  out[11]=4*in[0]*in[1]*(1-3*in[2]+2*in[2]*in[2]);
49  out[12]=4*in[0]*(1-in[0]-in[1])*(-in[2]+2*in[2]*in[2]);
50  out[13]=4*in[1]*(1-in[0]-in[1])*(-in[2]+2*in[2]*in[2]);
51  out[14]=4*in[0]*in[1]*(-in[2]+2*in[2]*in[2]);
52 
53  //faces
54  out[15]=4*in[0]*(1-in[0]-in[1])*(4*in[2]-4*in[2]*in[2]);
55  out[16]=4*in[1]*(1-in[0]-in[1])*(4*in[2]-4*in[2]*in[2]);
56  out[17]=4*in[0]*in[1]*(4*in[2]-4*in[2]*in[2]);
57 }
58 
59 
60 
62  void evaluateJacobian (const typename Traits::DomainType& in, //position
63  std::vector<typename Traits::JacobianType>& out) const //return value
64 {
65 out.resize(18);
66 
67  //vertices
68  out[0][0][0] = in[2]-1;
69  out[0][0][1] = in[2]-1;
70  out[0][0][2] = in[0]+in[1]-1;
71 
72  out[1][0][0] = 1-in[2];
73  out[1][0][1] = 0;
74  out[1][0][2] =-in[0];
75 
76  out[2][0][0] = 0;
77  out[2][0][1] = 1-in[2];
78  out[2][0][2] = -in[1];
79 
80  out[3][0][0] = -in[2];
81  out[3][0][1] = -in[2];
82  out[3][0][2] = 1-in[0]-in[1];
83 
84  out[4][0][0] = in[2];
85  out[4][0][1] = 0;
86  out[4][0][2] = in[0];
87 
88  out[5][0][0] = 0;
89  out[5][0][1] = in[2];
90  out[5][0][2] = in[1];
91 
92  //edges
93  out[6][0][0] = (-3+4*in[0]+4*in[1])*(4*in[2]-4*in[2]*in[2]);
94  out[6][0][1] = (-3+4*in[0]+4*in[1])*(4*in[2]-4*in[2]*in[2]);
95  out[6][0][2] = 2*(1-in[0]-in[1])*(0.5-in[0]-in[1])*(4-8*in[2]);
96 
97  out[7][0][0] = (-1+4*in[0])*(4*in[2]-4*in[2]*in[2]);
98  out[7][0][1] = 0;
99  out[7][0][2] = 2*in[0]*(-0.5+in[0])*(4-8*in[2]);
100 
101  out[8][0][0] = 0;
102  out[8][0][1] = (-1+4*in[1])*(4*in[2]-4*in[2]*in[2]);
103  out[8][0][2] = 2*in[1]*(-0.5+in[1])*(4-8*in[2]);
104 
105  out[9][0][0] = (4-8*in[0]-4*in[1])*(1-3*in[2]+2*in[2]*in[2]);
106  out[9][0][1] = -4*in[0]*(1-3*in[2]+2*in[2]*in[2]);
107  out[9][0][2] = 4*in[0]*(1-in[0]-in[1])*(-3+4*in[2]);
108 
109  out[10][0][0] = (-4*in[1])*(1-3*in[2]+2*in[2]*in[2]);
110  out[10][0][1] = (4-4*in[0]-8*in[1])*(1-3*in[2]+2*in[2]*in[2]);
111  out[10][0][2] = 4*in[1]*(1-in[0]-in[1])*(-3+4*in[2]);
112 
113  out[11][0][0] = 4*in[1]*(1-3*in[2]+2*in[2]*in[2]);
114  out[11][0][1] = 4*in[0]*(1-3*in[2]+2*in[2]*in[2]);
115  out[11][0][2] = 4*in[0]*in[1]*(-3+4*in[2]);
116 
117  out[12][0][0] = (4-8*in[0]-4*in[1])*(-in[2]+2*in[2]*in[2]);
118  out[12][0][1] = (-4*in[0])*(-in[2]+2*in[2]*in[2]);
119  out[12][0][2] = 4*in[0]*(1-in[0]-in[1])*(-1+4*in[2]);
120 
121  out[13][0][0] = -4*in[1]*(-in[2]+2*in[2]*in[2]);
122  out[13][0][1] = (4-4*in[0]-8*in[1])*(-in[2]+2*in[2]*in[2]);
123  out[13][0][2] = 4*in[1]*(1-in[0]-in[1])*(-1+4*in[2]);
124 
125  out[14][0][0] = 4*in[1]*(-in[2]+2*in[2]*in[2]);
126  out[14][0][1] = 4*in[0]*(-in[2]+2*in[2]*in[2]);
127  out[14][0][2] = 4*in[0]*in[1]*(-1+4*in[2]);
128 
129  //faces
130  out[15][0][0] = (4-8*in[0]-4*in[1])*(4*in[2]-4*in[2]*in[2]);
131  out[15][0][1] = -4*in[0]*(4*in[2]-4*in[2]*in[2]);
132  out[15][0][2] = 4*in[0]*(1-in[0]-in[1])*(4-8*in[2]);
133 
134  out[16][0][0] = -4*in[1]*(4*in[2]-4*in[2]*in[2]);
135  out[16][0][1] = (4-4*in[0]-8*in[1])*(4*in[2]-4*in[2]*in[2]);
136  out[16][0][2] = 4*in[1]*(1-in[0]-in[1])*(4-8*in[2]);
137 
138  out[17][0][0] = 4*in[1]*(4*in[2]-4*in[2]*in[2]);
139  out[17][0][1] = 4*in[0]*(4*in[2]-4*in[2]*in[2]);
140  out[17][0][2] = 4*in[0]*in[1]*(4-8*in[2]);
141 
142  }
145 unsigned int order() const
146 {
147 return 2;
148 }
149 
150  };
151 }
152 #endif