3 #ifndef DUNE_PYRAMID_P2_LOCALBASIS_HH
4 #define DUNE_PYRAMID_P2_LOCALBASIS_HH
6 #include <dune/common/fmatrix.hh>
26 template<
class D,
class R>
42 std::vector<typename Traits::RangeType>& out)
const
47 const R x = 2.0*in[0] + in[2] - 1.0;
48 const R y = 2.0*in[1] + in[2] - 1.0;
54 out[0] = 0.25*(x + z)*(x + z - 1)*(y - z - 1)*(y - z);
55 out[1] = -0.25*(x + z)*(y - z)*((x + z + 1)*(-y + z + 1) - 4*z) - z*(x - y);
56 out[2] = 0.25*(x + z)*(y - z)*(y - z + 1)*(x + z - 1);
57 out[3] = 0.25*(y - z)*(x + z)*(y - z + 1)*(x + z + 1);
61 out[5] = -0.5*(y - z + 1)*(x + z - 1)*((y - 1)*(x + 1) + z*(x - y + z + 1));
62 out[6] = -0.5*(y - z + 1)*(((x + z + 1)*(y - 1)*x - z) + z*(2*y + 1));
63 out[7] = -0.5*(x + z - 1)*(((y - z - 1)*(x + 1)*y - z) + z*(2*x + 1));
64 out[8] = -0.5*(y - z + 1)*(x + z - 1)*(x + 1)*y;
67 out[9] = z*(x + z - 1)*(y - z - 1);
68 out[10] = -z*((x + z + 1)*(y - z - 1) + 4*z);
69 out[11] = -z*(y - z + 1)*(x + z - 1);
70 out[12] = z*(y - z + 1)*(x + z + 1);
73 out[13] = (y - z + 1)*(x + z - 1)*((y - 1)*(x + 1) + z*(x - y + z + 1));
78 out[0] = 0.25*(y + z)*(y + z - 1)*(x - z - 1)*(x - z);
79 out[1] = -0.25*(x - z)*(y + z)*(x - z + 1)*(-y - z + 1);
80 out[2] = 0.25*(x - z)*(y + z)*((x - z - 1)*(y + z + 1) + 4*z) + z*(x - y);
81 out[3] = 0.25*(y + z)*(x - z)*(x - z + 1)*(y + z + 1);
85 out[5] = -0.5*(y + z - 1)*(((x - z - 1)*(y + 1)*x - z) + z*(2*y + 1));
86 out[6] = -0.5*(x - z + 1)*(y + z - 1)*(y + 1)*x;
87 out[7] = -0.5*(x - z + 1)*(y + z - 1)*(x - 1)*y;
88 out[8] = -0.5*(x - z + 1)*(((y + z + 1)*(x - 1)*y - z) + z*(2*x + 1));
91 out[9] = z*(y + z - 1)*(x - z - 1);
92 out[10] = -z*(x - z + 1)*(y + z - 1);
93 out[11] = -z*((y + z + 1)*(x - z - 1) + 4*z);
94 out[12] = z*(x - z + 1)*(y + z + 1);
97 out[13] = (x - z + 1)*(y + z - 1)*((y + 1)*(x - 1) - z*(x - y - z - 1));
104 std::vector<typename Traits::JacobianType>& out)
const
109 const R x = 2.0*in[0] + in[2] - 1.0;
110 const R y = 2.0*in[1] + in[2] - 1.0;
118 out[0][0][0] = 0.5*(y - z - 1)*(y - z)*(2*x + 2*z - 1);
119 out[0][0][1] = 0.5*(x + z)*(x + z - 1)*(2*y - 2*z - 1);
120 out[0][0][2] = 0.5*(out[0][0][0] + out[0][0][1])
121 + 0.25*((2*x + 2*z - 1)*(y - z - 1)*(y - z)
122 + (x + z)*(x + z - 1)*(-2*y + 2*z + 1));
124 out[1][0][0] = 2*(-0.25*((y - z)*((x + z + 1)*(-y + z + 1) - 4*z)
125 + (x + z)*(y - z)*(-y + z + 1)) - z);
126 out[1][0][1] = 2*(-0.25*((x + z)*((x + z + 1)*(-y + z + 1) - 4*z)
127 + (x + z)*(y - z)*(-(x + z + 1))) + z);
128 out[1][0][2] = 0.5*(out[1][0][0] + out[1][0][1])
129 - 0.25*((y - z)*((x + z + 1)*(-y + z + 1) - 4*z)
130 - (x + z)*((x + z + 1)*(-y + z + 1) - 4*z)
131 + (x + z)*(y - z)*(x - y + 2*z - 2))
134 out[2][0][0] = 0.5*(y - z)*(y - z + 1)*(2*x + 2*z - 1);
135 out[2][0][1] = 0.5*(x + z)*(2*y - 2*z + 1)*(x + z - 1);
136 out[2][0][2] = 0.5*(out[2][0][0] + out[2][0][1])
137 + 0.25*((y - x - 2*z)*(y - z + 1)*(x + z - 1)
138 + (x + z)*(y - z)*(y - x - 2*z + 2));
140 out[3][0][0] = 0.5*(y - z)*(2*x + 2*z + 1)*(y - z + 1);
141 out[3][0][1] = 0.5*(2*y - 2*z + 1)*(x + z)*(x + z + 1);
142 out[3][0][2] = 0.5*(out[3][0][0] + out[3][0][1])
143 + 0.25*((y - x - 2*z)*(y - z + 1)*(x + z + 1)
144 + (y - z)*(x + z)*(y - x - 2*z));
148 out[4][0][2] = 4*z - 1;
151 out[5][0][0] = -((y - z + 1)*((y - 1)*(x + 1) + z*(x - y + z + 1))
152 + (y - z + 1)*(x + z - 1)*((y - 1) + z));
153 out[5][0][1] = -((x + z - 1)*((y - 1)*(x + 1) + z*(x - y + z + 1))
154 + (y - z + 1)*(x + z - 1)*((x + 1) - z));
155 out[5][0][2] = 0.5*(out[5][0][0] + out[5][0][1])
156 - 0.5*((-x + y - 2*z + 2)*((y - 1)*(x + 1) + z*(x - y + z + 1))
157 + (y - z + 1)*(x + z - 1)*(x - y + 2*z + 1));
159 out[6][0][0] = -(y - z + 1)*(2*x + z + 1)*(y - 1);
160 out[6][0][1] = -(((x + z + 1)*(y - 1)*x - z) + z*(2*y + 1)
161 + (y - z + 1)*((x + z + 1)*x + 2*z));
162 out[6][0][2] = 0.5*(out[6][0][0] + out[6][0][1])
163 - 0.5*(-(((x + z + 1)*(y - 1)*x - z) + z*(2*y + 1))
164 + (y - z + 1)*(((y - 1)*x - 1) + 2*y + 1));
166 out[7][0][0] = -(((y - z - 1)*(x + 1)*y - z) + z*(2*x + 1)
167 + (x + z - 1)*((y - z - 1)*y + 2*z));
168 out[7][0][1] = -(x + z - 1)*(2*y - z - 1)*(x + 1);
169 out[7][0][2] = 0.5*(out[7][0][0] + out[7][0][1])
170 - 0.5*(((y - z - 1)*(x + 1)*y - z) + z*(2*x + 1)
171 + (x + z - 1)*((-(x + 1)*y - 1) + 2*x + 1));
173 out[8][0][0] = -(y - z + 1)*(2*x + z)*y;
174 out[8][0][1] = -(2*y - z + 1)*(x + z - 1)*(x + 1);
175 out[8][0][2] = 0.5*(out[8][0][0] + out[8][0][1])
176 - 0.5*(-x + y - 2*z + 2)*(x + 1)*y;
179 out[9][0][0] = 2*z*(y - z - 1);
180 out[9][0][1] = 2*z*(x + z - 1);
181 out[9][0][2] = 0.5*(out[9][0][0] + out[9][0][1])
182 + (x + z - 1)*(y - z - 1) + z*(-x + y - 2*z);
184 out[10][0][0] = -2*z*(y - z - 1);
185 out[10][0][1] = -2*z*(x + z + 1);
186 out[10][0][2] = 0.5*(out[10][0][0] + out[10][0][1])
187 - ((x + z + 1)*(y - z - 1) + 4*z)
188 - z*(-x + y - 2*z + 2);
190 out[11][0][0] = -2*z*(y - z + 1);
191 out[11][0][1] = -2*z*(x + z - 1);
192 out[11][0][2] = 0.5*(out[11][0][0] + out[11][0][1])
193 - (y - z + 1)*(x + z - 1) - z*(-x + y - 2*z + 2);
195 out[12][0][0] = 2*z*(y - z + 1);
196 out[12][0][1] = 2*z*(x + z + 1);
197 out[12][0][2] = 0.5*(out[12][0][0] + out[12][0][1])
198 + (y - z + 1)*(x + z + 1) + z*(-x + y - 2*z);
201 out[13][0][0] = 2*((y - z + 1)*((y - 1)*(x + 1) + z*(x - y + z + 1))
202 + (y - z + 1)*(x + z - 1)*(y - 1 + z));
203 out[13][0][1] = 2*((x + z - 1)*((y - 1)*(x + 1) + z*(x - y + z + 1))
204 + (y - z + 1)*(x + z - 1)*(x + 1 - z));
205 out[13][0][2] = 0.5*(out[13][0][0] + out[13][0][1])
206 + ((-x + y - 2*z + 2)*((y - 1)*(x + 1) + z*(x - y + z + 1))
207 + (y - z + 1)*(x + z - 1)*(x - y + 2*z + 1));
212 out[0][0][0] = 0.5*(y + z)*(y + z - 1)*(2*x - 2*z - 1);
213 out[0][0][1] = 0.5*(2*y + 2*z - 1)*(x - z - 1)*(x - z);
214 out[0][0][2] = 0.5*(out[0][0][0] + out[0][0][1])
215 + 0.25*((2*y + 2*z - 1)*(x - z - 1)*(x - z)
216 + (y + z)*(y + z - 1)*(-2*x + 2*z + 1));
218 out[1][0][0] = -0.5*(y + z)*(2*x - 2*z + 1)*(-y - z + 1);
219 out[1][0][1] = -0.5*(x - z)*(x - z + 1)*(-2*y - 2*z + 1);
220 out[1][0][2] = 0.5*(out[1][0][0] + out[1][0][1])
221 - 0.25*((x - y - 2*z)*(x - z + 1)*(-y - z + 1)
222 + (x - z)*(y + z)*(-x + y + 2*z - 2));
224 out[2][0][0] = 0.5*((y + z)*((x - z - 1)*(y + z + 1) + 4*z)
225 + (x - z)*(y + z)*(y + z + 1) + 4*z);
226 out[2][0][1] = 0.5*((x - z)*((x - z - 1)*(y + z + 1) + 4*z)
227 + (x - z)*(y + z)*(x - z - 1) - 4*z);
228 out[2][0][2] = 0.5*(out[2][0][0] + out[2][0][1])
229 + 0.25*((x - y - 2*z)*((x - z - 1)*(y + z + 1) + 4*z)
230 + (x - z)*(y + z)*(x - y - 2*z + 2) + 4*(x - y));
232 out[3][0][0] = 0.5*(y + z)*(2*x - 2*z + 1)*(y + z + 1);
233 out[3][0][1] = 0.5*(x - z)*(x - z + 1)*(2*y + 2*z + 1);
234 out[3][0][2] = 0.5*(out[3][0][0] + out[3][0][1])
235 + 0.25*((x - y - 2*z)*(x - z + 1)*(y + z + 1)
236 + (y + z)*(x - z)*(x - y - 2*z));
240 out[4][0][2] = 4*z - 1;
243 out[5][0][0] = -(y + z - 1)*(2*x - z - 1)*(y + 1);
244 out[5][0][1] = -(((x - z - 1)*(y + 1)*x - z) + z*(2*y + 1)
245 + (y + z - 1)*((x - z - 1)*x + 2*z));
246 out[5][0][2] = 0.5*(out[5][0][0] + out[5][0][1])
247 - 0.5*((((x - z - 1)*(y + 1)*x - z) + z*(2*y + 1))
248 + (y + z - 1)*((-(y + 1)*x - 1) + 2*y + 1));
250 out[6][0][0] = -(2*x - z + 1)*(y + z - 1)*(y + 1);
251 out[6][0][1] = -(x - z + 1)*(2*y + z)*x;
252 out[6][0][2] = 0.5*(out[6][0][0] + out[6][0][1])
253 - 0.5*(x - y - 2*z + 2)*(y + 1)*x;
255 out[7][0][0] = -(2*x - z)*(y + z - 1)*y;
256 out[7][0][1] = -(x - z + 1)*(2*y + z - 1)*(x - 1);
257 out[7][0][2] = 0.5*(out[7][0][0] + out[7][0][1])
258 - 0.5*(x - y - 2*z + 2)*(x - 1)*y;
260 out[8][0][0] = -(((y + z + 1)*(x - 1)*y - z) + z*(2*x + 1)
261 + (x - z + 1)*((y + z + 1)*y + 2*z));
262 out[8][0][1] = -(x - z + 1)*(2*y + z + 1)*(x - 1);
263 out[8][0][2] = 0.5*(out[8][0][0] + out[8][0][1])
264 - 0.5*(-(((y + z + 1)*(x - 1)*y - z) + z*(2*x + 1))
265 + (x - z + 1)*(((x - 1)*y - 1) + 2*x + 1));
268 out[9][0][0] = 2*z*(y + z - 1);
269 out[9][0][1] = 2*z*(x - z - 1);
270 out[9][0][2] = 0.5*(out[9][0][0] + out[9][0][1])
271 + (y + z - 1)*(x - z - 1) + z*(x - y - 2*z);
273 out[10][0][0] = -2*z*(y + z - 1);
274 out[10][0][1] = -2*z*(x - z + 1);
275 out[10][0][2] = 0.5*(out[10][0][0] + out[10][0][1])
276 - (x - z + 1)*(y + z - 1) - z*(x - y - 2*z + 2);
278 out[11][0][0] = -2*z*(y + z + 1);
279 out[11][0][1] = -2*z*(x - z - 1);
280 out[11][0][2] = 0.5*(out[11][0][0] + out[11][0][1])
281 - ((y + z + 1)*(x - z - 1) + 4*z) - z*(x - y - 2*z + 2);
283 out[12][0][0] = 2*z*(y + z + 1);
284 out[12][0][1] = 2*z*(x - z + 1);
285 out[12][0][2] = 0.5*(out[12][0][0] + out[12][0][1])
286 + (x - z + 1)*(y + z + 1) + z*(x - y - 2*z);
289 out[13][0][0] = 2*((y + z - 1)*((y + 1)*(x - 1) - z*(x - y - z - 1))
290 + (x - z + 1)*(y + z - 1)*(y + 1 - z));
291 out[13][0][1] = 2*((x - z + 1)*((y + 1)*(x - 1) - z*(x - y - z - 1))
292 + (x - z + 1)*(y + z - 1)*(x - 1 + z));
293 out[13][0][2] = 0.5*(out[13][0][0] + out[13][0][1])
294 + (x - y - 2*z + 2)*((y + 1)*(x - 1) - z*(x - y - z - 1))
295 + (x - z + 1)*(y + z - 1)*(-(x - y - 2*z - 1));