49 #ifndef __INTREPID2_HGRAD_HEX_C2_FEM_DEF_HPP__
50 #define __INTREPID2_HGRAD_HEX_C2_FEM_DEF_HPP__
57 template<EOperator opType>
58 template<
typename outputViewType,
59 typename inputViewType>
60 KOKKOS_INLINE_FUNCTION
62 Basis_HGRAD_HEX_C2_FEM::Serial<opType>::
63 getValues( outputViewType output,
64 const inputViewType input ) {
66 case OPERATOR_VALUE : {
67 const auto x = input(0);
68 const auto y = input(1);
69 const auto z = input(2);
72 output.access( 0) = 0.125*(-1. + x)*x*(-1. + y)*y*(-1. + z)*z;
73 output.access( 1) = 0.125*x*(1.+ x)*(-1. + y)*y*(-1. + z)*z;
74 output.access( 2) = 0.125*x*(1.+ x)*y*(1.+ y)*(-1. + z)*z;
75 output.access( 3) = 0.125*(-1. + x)*x*y*(1.+ y)*(-1. + z)*z;
76 output.access( 4) = 0.125*(-1. + x)*x*(-1. + y)*y*z*(1.+ z);
77 output.access( 5) = 0.125*x*(1.+ x)*(-1. + y)*y*z*(1.+ z);
78 output.access( 6) = 0.125*x*(1.+ x)*y*(1.+ y)*z*(1.+ z);
79 output.access( 7) = 0.125*(-1. + x)*x*y*(1.+ y)*z*(1.+ z);
80 output.access( 8) = 0.25*(1. - x)*(1. + x)*(-1. + y)*y*(-1. + z)*z;
81 output.access( 9) = 0.25*x*(1.+ x)*(1. - y)*(1. + y)*(-1. + z)*z;
82 output.access(10) = 0.25*(1. - x)*(1. + x)*y*(1.+ y)*(-1. + z)*z;
83 output.access(11) = 0.25*(-1. + x)*x*(1. - y)*(1. + y)*(-1. + z)*z;
84 output.access(12) = 0.25*(-1. + x)*x*(-1. + y)*y*(1. - z)*(1. + z);
85 output.access(13) = 0.25*x*(1.+ x)*(-1. + y)*y*(1. - z)*(1. + z);
86 output.access(14) = 0.25*x*(1.+ x)*y*(1.+ y)*(1. - z)*(1. + z);
87 output.access(15) = 0.25*(-1. + x)*x*y*(1.+ y)*(1. - z)*(1. + z);
88 output.access(16) = 0.25*(1. - x)*(1. + x)*(-1. + y)*y*z*(1.+ z);
89 output.access(17) = 0.25*x*(1.+ x)*(1. - y)*(1. + y)*z*(1.+ z);
90 output.access(18) = 0.25*(1. - x)*(1. + x)*y*(1.+ y)*z*(1.+ z);
91 output.access(19) = 0.25*(-1. + x)*x*(1. - y)*(1. + y)*z*(1.+ z);
92 output.access(20) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
93 output.access(21) = 0.5*(1. - x)*(1. + x)*(1. - y)*(1. + y)*(-1. + z)*z;
94 output.access(22) = 0.5*(1. - x)*(1. + x)*(1. - y)*(1. + y)*z*(1.+ z);
95 output.access(23) = 0.5*(-1. + x)*x*(1. - y)*(1. + y)*(1. - z)*(1. + z);
96 output.access(24) = 0.5*x*(1.+ x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
97 output.access(25) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y*(1. - z)*(1. + z);
98 output.access(26) = 0.5*(1. - x)*(1. + x)*y*(1.+ y)*(1. - z)*(1. + z);
103 const auto x = input(0);
104 const auto y = input(1);
105 const auto z = input(2);
108 output.access(0, 0) = (-0.125 + 0.25*x)*(-1. + y)*y*(-1. + z)*z;
109 output.access(0, 1) = (-1. + x)*x*(-0.125 + 0.25*y)*(-1. + z)*z;
110 output.access(0, 2) = (-1. + x)*x*(-1. + y)*y*(-0.125 + 0.25*z);
112 output.access(1, 0) = (0.125 + 0.25*x)*(-1. + y)*y*(-1. + z)*z;
113 output.access(1, 1) = x*(1. + x)*(-0.125 + 0.25*y)*(-1. + z)*z;
114 output.access(1, 2) = x*(1. + x)*(-1. + y)*y*(-0.125 + 0.25*z);
116 output.access(2, 0) = (0.125 + 0.25*x)*y*(1. + y)*(-1. + z)*z;
117 output.access(2, 1) = x*(1. + x)*(0.125 + 0.25*y)*(-1. + z)*z;
118 output.access(2, 2) = x*(1. + x)*y*(1. + y)*(-0.125 + 0.25*z);
120 output.access(3, 0) = (-0.125 + 0.25*x)*y*(1. + y)*(-1. + z)*z;
121 output.access(3, 1) = (-1. + x)*x*(0.125 + 0.25*y)*(-1. + z)*z;
122 output.access(3, 2) = (-1. + x)*x*y*(1. + y)*(-0.125 + 0.25*z);
124 output.access(4, 0) = (-0.125 + 0.25*x)*(-1. + y)*y*z*(1. + z);
125 output.access(4, 1) = (-1. + x)*x*(-0.125 + 0.25*y)*z*(1. + z);
126 output.access(4, 2) = (-1. + x)*x*(-1. + y)*y*(0.125 + 0.25*z);
128 output.access(5, 0) = (0.125 + 0.25*x)*(-1. + y)*y*z*(1. + z);
129 output.access(5, 1) = x*(1. + x)*(-0.125 + 0.25*y)*z*(1. + z);
130 output.access(5, 2) = x*(1. + x)*(-1. + y)*y*(0.125 + 0.25*z);
132 output.access(6, 0) = (0.125 + 0.25*x)*y*(1. + y)*z*(1. + z);
133 output.access(6, 1) = x*(1. + x)*(0.125 + 0.25*y)*z*(1. + z);
134 output.access(6, 2) = x*(1. + x)*y*(1. + y)*(0.125 + 0.25*z);
136 output.access(7, 0) = (-0.125 + 0.25*x)*y*(1. + y)*z*(1. + z);
137 output.access(7, 1) = (-1. + x)*x*(0.125 + 0.25*y)*z*(1. + z);
138 output.access(7, 2) = (-1. + x)*x*y*(1. + y)*(0.125 + 0.25*z);
140 output.access(8, 0) = -0.5*x*(-1. + y)*y*(-1. + z)*z;
141 output.access(8, 1) = (1. - x)*(1. + x)*(-0.25 + 0.5*y)*(-1. + z)*z;
142 output.access(8, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(-0.25 + 0.5*z);
144 output.access(9, 0) = (0.25 + 0.5*x)*(1. - y)*(1. + y)*(-1. + z)*z;
145 output.access(9, 1) = x*(1. + x)*(-0.5*y)*(-1. + z)*z;
146 output.access(9, 2) = x*(1. + x)*(1. - y)*(1. + y)*(-0.25 + 0.5*z);
148 output.access(10, 0) = -0.5*x*y*(1. + y)*(-1. + z)*z;
149 output.access(10, 1) = (1. - x)*(1. + x)*(0.25 + 0.5*y)*(-1. + z)*z;
150 output.access(10, 2) = (1. - x)*(1. + x)*y*(1. + y)*(-0.25 + 0.5*z);
152 output.access(11, 0) = (-0.25 + 0.5*x)*(1. - y)*(1. + y)*(-1. + z)*z;
153 output.access(11, 1) = (-1. + x)*x*(-0.5*y)*(-1. + z)*z;
154 output.access(11, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(-0.25 + 0.5*z);
156 output.access(12, 0) = (-0.25 + 0.5*x)*(-1. + y)*y*(1. - z)*(1. + z);
157 output.access(12, 1) = (-1. + x)*x*(-0.25 + 0.5*y)*(1. - z)*(1. + z);
158 output.access(12, 2) = (-1. + x)*x*(-1. + y)*y*(-0.5*z);
160 output.access(13, 0) = (0.25 + 0.5*x)*(-1. + y)*y*(1. - z)*(1. + z);
161 output.access(13, 1) = x*(1. + x)*(-0.25 + 0.5*y)*(1. - z)*(1. + z);
162 output.access(13, 2) = x*(1. + x)*(-1. + y)*y*(-0.5*z);
164 output.access(14, 0) = (0.25 + 0.5*x)*y*(1. + y)*(1. - z)*(1. + z);
165 output.access(14, 1) = x*(1. + x)*(0.25 + 0.5*y)*(1. - z)*(1. + z);
166 output.access(14, 2) = x*(1. + x)*y*(1. + y)*(-0.5*z);
168 output.access(15, 0) = (-0.25 + 0.5*x)*y*(1. + y)*(1. - z)*(1. + z);
169 output.access(15, 1) = (-1. + x)*x*(0.25 + 0.5*y)*(1. - z)*(1. + z);
170 output.access(15, 2) = (-1. + x)*x*y*(1. + y)*(-0.5*z);
172 output.access(16, 0) = -0.5*x*(-1. + y)*y*z*(1. + z);
173 output.access(16, 1) = (1. - x)*(1. + x)*(-0.25 + 0.5*y)*z*(1. + z);
174 output.access(16, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(0.25 + 0.5*z);
176 output.access(17, 0) = (0.25 + 0.5*x)*(1. - y)*(1. + y)*z*(1. + z);
177 output.access(17, 1) = x*(1. + x)*(-0.5*y)*z*(1. + z);
178 output.access(17, 2) = x*(1. + x)*(1. - y)*(1. + y)*(0.25 + 0.5*z);
180 output.access(18, 0) = -0.5*x*y*(1. + y)*z*(1. + z);
181 output.access(18, 1) = (1. - x)*(1. + x)*(0.25 + 0.5*y)*z*(1. + z);
182 output.access(18, 2) = (1. - x)*(1. + x)*y*(1. + y)*(0.25 + 0.5*z);
184 output.access(19, 0) = (-0.25 + 0.5*x)*(1. - y)*(1. + y)*z*(1. + z);
185 output.access(19, 1) = (-1. + x)*x*(-0.5*y)*z*(1. + z);
186 output.access(19, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(0.25 + 0.5*z);
188 output.access(20, 0) = -2.*x*(1. - y)*(1. + y)*(1. - z)*(1. + z);
189 output.access(20, 1) = (1. - x)*(1. + x)*(-2.*y)*(1. - z)*(1. + z);
190 output.access(20, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(-2.*z);
192 output.access(21, 0) = -x*(1. - y)*(1. + y)*(-1. + z)*z;
193 output.access(21, 1) = (1. - x)*(1. + x)*(-y)*(-1. + z)*z;
194 output.access(21, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(-0.5 + z);
196 output.access(22, 0) = -x*(1. - y)*(1. + y)*z*(1. + z);
197 output.access(22, 1) = (1. - x)*(1. + x)*(-y)*z*(1. + z);
198 output.access(22, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(0.5 + z);
200 output.access(23, 0) = (-0.5 + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
201 output.access(23, 1) = (-1. + x)*x*(-y)*(1. - z)*(1. + z);
202 output.access(23, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(-z);
204 output.access(24, 0) = (0.5 + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
205 output.access(24, 1) = x*(1. + x)*(-y)*(1. - z)*(1. + z);
206 output.access(24, 2) = x*(1. + x)*(1. - y)*(1. + y)*(-z);
208 output.access(25, 0) = -x*(-1. + y)*y*(1. - z)*(1. + z);
209 output.access(25, 1) = (1. - x)*(1. + x)*(-0.5 + y)*(1. - z)*(1. + z);
210 output.access(25, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(-z);
212 output.access(26, 0) = -x*y*(1. + y)*(1. - z)*(1. + z);
213 output.access(26, 1) = (1. - x)*(1. + x)*(0.5 + y)*(1. - z)*(1. + z);
214 output.access(26, 2) = (1. - x)*(1. + x)*y*(1. + y)*(-z);
218 const auto x = input(0);
219 const auto y = input(1);
220 const auto z = input(2);
223 output.access(0, 0) = 0.25*(-1. + y)*y*(-1. + z)*z;
224 output.access(0, 1) = (-0.125 + y*(0.25 - 0.25*z) + x*(0.25 + y*(-0.5 + 0.5*z) - 0.25*z) + 0.125*z)*z;
225 output.access(0, 2) = y*(-0.125 + x*(0.25 + y*(-0.25 + 0.5*z) - 0.5*z) + y*(0.125 - 0.25*z) + 0.25*z);
226 output.access(0, 3) = 0.25*(-1. + x)*x*(-1. + z)*z;
227 output.access(0, 4) = x*(-0.125 + y*(0.25 - 0.5*z) + x*(0.125 + y*(-0.25 + 0.5*z) - 0.25*z) + 0.25*z);
228 output.access(0, 5) = 0.25*(-1. + x)*x*(-1. + y)*y;
230 output.access(1, 0) = 0.25*(-1. + y)*y*(-1. + z)*z;
231 output.access(1, 1) = (0.125 + x*(0.25 + y*(-0.5 + 0.5*z) - 0.25*z) + y*(-0.25 + 0.25*z) - 0.125*z)*z;
232 output.access(1, 2) = y*(0.125 + x*(0.25 + y*(-0.25 + 0.5*z) - 0.5*z) + y*(-0.125 + 0.25*z) - 0.25*z);
233 output.access(1, 3) = 0.25*x*(1 + x)*(-1. + z)*z;
234 output.access(1, 4) = x*(1. + x)*(0.125 + y*(-0.25 + 0.5*z) - 0.25*z);
235 output.access(1, 5) = 0.25*x*(1 + x)*(-1. + y)*y;
237 output.access(2, 0) = 0.25*y*(1 + y)*(-1. + z)*z;
238 output.access(2, 1) = (0.125 + x*(0.25 + 0.5*y) + 0.25*y)*(-1. + z)*z;
239 output.access(2, 2) = y*(1. + y)*(-0.125 + x*(-0.25 + 0.5*z) + 0.25*z);
240 output.access(2, 3) = 0.25*x*(1 + x)*(-1. + z)*z;
241 output.access(2, 4) = x*(1. + x)*(-0.125 + y*(-0.25 + 0.5*z) + 0.25*z);
242 output.access(2, 5) = 0.25*x*(1 + x)*y*(1 + y);
244 output.access(3, 0) = 0.25*y*(1 + y)*(-1. + z)*z;
245 output.access(3, 1) = (0.125 + y*(0.25 - 0.25*z) + x*(-0.25 + y*(-0.5 + 0.5*z) + 0.25*z) - 0.125*z)*z;
246 output.access(3, 2) = y*(1. + y)*(0.125 + x*(-0.25 + 0.5*z) - 0.25*z);
247 output.access(3, 3) = 0.25*(-1. + x)*x*(-1. + z)*z;
248 output.access(3, 4) = x*(0.125 + y*(0.25 - 0.5*z) + x*(-0.125 + y*(-0.25 + 0.5*z) + 0.25*z) - 0.25*z);
249 output.access(3, 5) = 0.25*(-1. + x)*x*y*(1 + y);
251 output.access(4, 0) = 0.25*(-1. + y)*y*z*(1 + z);
252 output.access(4, 1) = (0.125 + x*(-0.25 + 0.5*y) - 0.25*y)*z*(1. + z);
253 output.access(4, 2) = y*(0.125 + x*(-0.25 + y*(0.25 + 0.5*z) - 0.5*z) + y*(-0.125 - 0.25*z) + 0.25*z);
254 output.access(4, 3) = 0.25*(-1. + x)*x*z*(1 + z);
255 output.access(4, 4) = x*(0.125 + y*(-0.25 - 0.5*z) + x*(-0.125 + y*(0.25 + 0.5*z) - 0.25*z) + 0.25*z);
256 output.access(4, 5) = 0.25*(-1. + x)*x*(-1. + y)*y;
258 output.access(5, 0) = 0.25*(-1. + y)*y*z*(1 + z);
259 output.access(5, 1) = (-0.125 + x*(-0.25 + 0.5*y) + 0.25*y)*z*(1. + z);
260 output.access(5, 2) = (-1. + y)*y*(0.125 + x*(0.25 + 0.5*z) + 0.25*z);
261 output.access(5, 3) = 0.25*x*(1 + x)*z*(1 + z);
262 output.access(5, 4) = x*(1. + x)*(-0.125 + y*(0.25 + 0.5*z) - 0.25*z);
263 output.access(5, 5) = 0.25*x*(1 + x)*(-1. + y)*y;
265 output.access(6, 0) = 0.25*y*(1 + y)*z*(1 + z);
266 output.access(6, 1) = (0.125 + x*(0.25 + 0.5*y) + 0.25*y)*z*(1. + z);
267 output.access(6, 2) = y*(1. + y)*(0.125 + x*(0.25 + 0.5*z) + 0.25*z);
268 output.access(6, 3) = 0.25*x*(1 + x)*z*(1 + z);
269 output.access(6, 4) = x*(1. + x)*(0.125 + y*(0.25 + 0.5*z) + 0.25*z);
270 output.access(6, 5) = 0.25*x*(1 + x)*y*(1 + y);
272 output.access(7, 0) = 0.25*y*(1 + y)*z*(1 + z);
273 output.access(7, 1) = (-0.125 + x*(0.25 + 0.5*y) - 0.25*y)*z*(1. + z);
274 output.access(7, 2) = y*(1. + y)*(-0.125 + x*(0.25 + 0.5*z) - 0.25*z);
275 output.access(7, 3) = 0.25*(-1. + x)*x*z*(1 + z);
276 output.access(7, 4) = (-1. + x)*x*(0.125 + y*(0.25 + 0.5*z) + 0.25*z);
277 output.access(7, 5) = 0.25*(-1. + x)*x*y*(1 + y);
279 output.access(8, 0) = -0.5*(-1. + y)*y*(-1. + z)*z;
280 output.access(8, 1) = (0. + x*(-0.5 + y))*z + (x*(0.5 - y) )*(z*z);
281 output.access(8, 2) = (y*y)*(x*(0.5 - z) ) + y*(x*(-0.5 + z));
282 output.access(8, 3) = 0.5*(1. - x)*(1. + x)*(-1. + z)*z;
283 output.access(8, 4) = 0.25 + (x*x)*(-0.25 + y*(0.5 - z) + 0.5*z) - 0.5*z + y*(-0.5 + z);
284 output.access(8, 5) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y;
286 output.access(9, 0) = 0.5*(1. - y)*(1. + y)*(-1. + z)*z;
287 output.access(9, 1) = (0.5*y + x*(y))*z + (x*(-y) - 0.5*y)*(z*z);
288 output.access(9, 2) = -0.25 + (y*y)*(0.25 - 0.5*z) + 0.5*z + x*(-0.5 + (y*y)*(0.5 - z) + z);
289 output.access(9, 3) = -0.5*x*(1 + x)*(-1. + z)*z;
290 output.access(9, 4) = x*(y*(0.5 - z) ) + (x*x)*(y*(0.5 - z) );
291 output.access(9, 5) = 0.5*x*(1 + x)*(1. - y)*(1. + y);
293 output.access(10, 0) = -0.5*y*(1 + y)*(-1. + z)*z;
294 output.access(10, 1) = (0. + x*(0.5 + y))*z + (x*(-0.5 - y) )*(z*z);
295 output.access(10, 2) = y*(x*(0.5 - z) ) + (y*y)*(x*(0.5 - z) );
296 output.access(10, 3) = 0.5*(1. - x)*(1. + x)*(-1. + z)*z;
297 output.access(10, 4) = -0.25 + (x*x)*(0.25 + y*(0.5 - z) - 0.5*z) + 0.5*z + y*(-0.5 + z);
298 output.access(10, 5) = 0.5*(1. - x)*(1. + x)*y*(1 + y);
300 output.access(11, 0) = 0.5*(1. - y)*(1. + y)*(-1. + z)*z;
301 output.access(11, 1) = (-0.5*y + x*(y))*z + (x*(-y) + 0.5*y)*(z*z);
302 output.access(11, 2) = 0.25 + (y*y)*(-0.25 + 0.5*z) - 0.5*z + x*(-0.5 + (y*y)*(0.5 - z) + z);
303 output.access(11, 3) = -0.5*(-1. + x)*x*(-1. + z)*z;
304 output.access(11, 4) = (x*x)*(y*(0.5 - z) ) + x*(y*(-0.5 + z));
305 output.access(11, 5) = 0.5*(-1. + x)*x*(1. - y)*(1. + y);
307 output.access(12, 0) = 0.5*(-1. + y)*y*(1. - z)*(1. + z);
308 output.access(12, 1) = 0.25 - 0.25*(z*z) + y*(-0.5 + 0.5*(z*z)) + x*(-0.5 + 0.5*(z*z) + y*(1. - (z*z)));
309 output.access(12, 2) = (y*y)*(x*(-z) + 0.5*z) + y*(-0.5*z + x*(z));
310 output.access(12, 3) = 0.5*(-1. + x)*x*(1. - z)*(1. + z);
311 output.access(12, 4) = (x*x)*(y*(-z) + 0.5*z) + x*(-0.5*z + y*(z));
312 output.access(12, 5) = -0.5*(-1. + x)*x*(-1. + y)*y;
314 output.access(13, 0) = 0.5*(-1. + y)*y*(1. - z)*(1. + z);
315 output.access(13, 1) = -0.25 + 0.25*(z*z) + y*(0.5 - 0.5*(z*z)) + x*(-0.5 + 0.5*(z*z) + y*(1. - (z*z)));
316 output.access(13, 2) = (y*y)*(x*(-z) - 0.5*z) + y*(0.5*z + x*(z));
317 output.access(13, 3) = 0.5*x*(1 + x)*(1. - z)*(1. + z);
318 output.access(13, 4) = x*(y*(-z) + 0.5*z) + (x*x)*(y*(-z) + 0.5*z);
319 output.access(13, 5) = -0.5*x*(1 + x)*(-1. + y)*y;
321 output.access(14, 0) = 0.5*y*(1 + y)*(1. - z)*(1. + z);
322 output.access(14, 1) = 0.25 - 0.25*(z*z) + y*(0.5 - 0.5*(z*z)) + x*(0.5 - 0.5*(z*z) + y*(1. - (z*z)));
323 output.access(14, 2) = y*(x*(-z) - 0.5*z) + (y*y)*(x*(-z) - 0.5*z);
324 output.access(14, 3) = 0.5*x*(1 + x)*(1. - z)*(1. + z);
325 output.access(14, 4) = x*(y*(-z) - 0.5*z) + (x*x)*(y*(-z) - 0.5*z);
326 output.access(14, 5) = -0.5*x*(1 + x)*y*(1 + y);
328 output.access(15, 0) = 0.5*y*(1 + y)*(1. - z)*(1. + z);
329 output.access(15, 1) = -0.25 + 0.25*(z*z) + y*(-0.5 + 0.5*(z*z)) + x*(0.5 - 0.5*(z*z) + y*(1. - (z*z)));
330 output.access(15, 2) = y*(x*(-z) + 0.5*z) + (y*y)*(x*(-z) + 0.5*z);
331 output.access(15, 3) = 0.5*(-1. + x)*x*(1. - z)*(1. + z);
332 output.access(15, 4) = (x*x)*(y*(-z) - 0.5*z) + x*(0.5*z + y*(z));
333 output.access(15, 5) = -0.5*(-1. + x)*x*y*(1 + y);
335 output.access(16, 0) = -0.5*(-1. + y)*y*z*(1 + z);
336 output.access(16, 1) = (x*(0.5 - y) )*z + (x*(0.5 - y) )*(z*z);
337 output.access(16, 2) = (y*y)*(x*(-0.5 - z) ) + y*(x*(0.5 + z));
338 output.access(16, 3) = 0.5*(1. - x)*(1. + x)*z*(1 + z);
339 output.access(16, 4) = -0.25 + (x*x)*(0.25 + y*(-0.5 - z) + 0.5*z) - 0.5*z + y*(0.5 + z);
340 output.access(16, 5) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y;
342 output.access(17, 0) = 0.5*(1. - y)*(1. + y)*z*(1 + z);
343 output.access(17, 1) = (x*(-y) - 0.5*y)*z + (x*(-y) - 0.5*y)*(z*z);
344 output.access(17, 2) = 0.25 + (y*y)*(-0.25 - 0.5*z) + 0.5*z + x*(0.5 + (y*y)*(-0.5 - z) + z);
345 output.access(17, 3) = -0.5*x*(1 + x)*z*(1 + z);
346 output.access(17, 4) = x*(y*(-0.5 - z) ) + (x*x)*(y*(-0.5 - z) );
347 output.access(17, 5) = 0.5*x*(1 + x)*(1. - y)*(1. + y);
349 output.access(18, 0) = -0.5*y*(1 + y)*z*(1 + z);
350 output.access(18, 1) = (x*(-0.5 - y) )*z + (x*(-0.5 - y) )*(z*z);
351 output.access(18, 2) = y*(x*(-0.5 - z) ) + (y*y)*(x*(-0.5 - z) );
352 output.access(18, 3) = 0.5*(1. - x)*(1. + x)*z*(1 + z);
353 output.access(18, 4) = 0.25 + (x*x)*(-0.25 + y*(-0.5 - z) - 0.5*z) + 0.5*z + y*(0.5 + z);
354 output.access(18, 5) = 0.5*(1. - x)*(1. + x)*y*(1 + y);
356 output.access(19, 0) = 0.5*(1. - y)*(1. + y)*z*(1 + z);
357 output.access(19, 1) = (x*(-y) + 0.5*y)*z + (x*(-y) + 0.5*y)*(z*z);
358 output.access(19, 2) = -0.25 + (y*y)*(0.25 + 0.5*z) - 0.5*z + x*(0.5 + (y*y)*(-0.5 - z) + z);
359 output.access(19, 3) = -0.5*(-1. + x)*x*z*(1 + z);
360 output.access(19, 4) = (x*x)*(y*(-0.5 - z) ) + x*(y*(0.5 + z));
361 output.access(19, 5) = 0.5*(-1. + x)*x*(1. - y)*(1. + y);
363 output.access(20, 0) = -2.*(1. - y)*(1. + y)*(1. - z)*(1. + z);
364 output.access(20, 1) = -4.*x*y*(-1. + z*z);
365 output.access(20, 2) = x*((y*y)*(-4.*z) + 4.*z);
366 output.access(20, 3) = -2.*(1. - x)*(1. + x)*(1. - z)*(1. + z);
367 output.access(20, 4) = (x*x)*(y*(-4.*z) ) + y*(4.*z);
368 output.access(20, 5) = -2.*(1. - x)*(1. + x)*(1. - y)*(1. + y);
370 output.access(21, 0) = -(1. - y)*(1. + y)*(-1. + z)*z;
371 output.access(21, 1) = (x*(-2.*y) )*z + (0. + x*(2.*y))*(z*z);
372 output.access(21, 2) = x*(1. - 2.*z + (y*y)*(-1. + 2.*z));
373 output.access(21, 3) = -(1. - x)*(1. + x)*(-1. + z)*z;
374 output.access(21, 4) = y*(1. - 2.*z) + (x*x)*(y*(-1. + 2.*z));
375 output.access(21, 5) = (1. - x)*(1. + x)*(1. - y)*(1. + y);
377 output.access(22, 0) = -(1. - y)*(1. + y)*z*(1 + z);
378 output.access(22, 1) = (0. + x*(2.*y))*z + (0. + x*(2.*y))*(z*z);
379 output.access(22, 2) = x*(-1. - 2.*z + (y*y)*(1. + 2.*z));
380 output.access(22, 3) = -(1. - x)*(1. + x)*z*(1 + z);
381 output.access(22, 4) = y*(-1. - 2.*z) + (x*x)*(y*(1. + 2.*z));
382 output.access(22, 5) = (1. - x)*(1. + x)*(1. - y)*(1. + y);
384 output.access(23, 0) = (1. - y)*(1. + y)*(1. - z)*(1. + z);
385 output.access(23, 1) = (-1. + 2.*x)*y*(-1. + z*z);
386 output.access(23, 2) = (-1. + 2.*x)*(-1. + y*y)*z;
387 output.access(23, 3) =-(-1. + x)*x*(1. - z)*(1. + z);
388 output.access(23, 4) = 2.*(-1. + x)*x*y*z;
389 output.access(23, 5) =-(-1. + x)*x*(1. - y)*(1. + y);
391 output.access(24, 0) = (1. - y)*(1. + y)*(1. - z)*(1. + z);
392 output.access(24, 1) = (1. + 2.*x)*y*(-1. + z*z);
393 output.access(24, 2) = (1. + 2.*x)*(-1. + y*y)*z;
394 output.access(24, 3) = x*(1. + x)*(-1. + z)*(1. + z);
395 output.access(24, 4) = 2.*x*(1. + x)*y*z;
396 output.access(24, 5) = x*(1. + x)*(-1. + y)*(1. + y);
398 output.access(25, 0) = -(-1. + y)*y*(1. - z)*(1. + z);
399 output.access(25, 1) = x*(-1. + 2.*y)*(-1. + z*z);
400 output.access(25, 2) = 2.*x*(-1. + y)*y*z;
401 output.access(25, 3) = (1. - x)*(1. + x)*(1. - z)*(1. + z);
402 output.access(25, 4) = (-1. + x*x)*(-1. + 2.*y)*z;
403 output.access(25, 5) =-(1. - x)*(1. + x)*(-1. + y)*y;
405 output.access(26, 0) = y*(1. + y)*(-1. + z)*(1. + z);
406 output.access(26, 1) = x*(1. + 2.*y)*(-1. + z*z);
407 output.access(26, 2) = 2.*x*y*(1. + y)*z;
408 output.access(26, 3) = (-1. + x)*(1. + x)*(-1. + z)*(1. + z);
409 output.access(26, 4) = (-1. + x*x)*(1. + 2.*y)*z;
410 output.access(26, 5) = (-1. + x)*(1. + x)*y*(1. + y);
414 const auto x = input(0);
415 const auto y = input(1);
416 const auto z = input(2);
418 output.access(0, 0) = 0.;
419 output.access(0, 1) = ((-1.+ 2.*y)*(-1.+ z)*z)/4.;
420 output.access(0, 2) = ((-1.+ y)*y*(-1.+ 2.*z))/4.;
421 output.access(0, 3) = ((-1.+ 2.*x)*(-1.+ z)*z)/4.;
422 output.access(0, 4) = ((-1.+ 2.*x)*(-1.+ 2.*y)*(-1.+ 2.*z))/8.;
423 output.access(0, 5) = ((-1.+ 2.*x)*(-1.+ y)*y)/4.;
424 output.access(0, 6) = 0.;
425 output.access(0, 7) = ((-1.+ x)*x*(-1.+ 2.*z))/4.;
426 output.access(0, 8) = ((-1.+ x)*x*(-1.+ 2.*y))/4.;
427 output.access(0, 9) = 0.;
429 output.access(1, 0) = 0.;
430 output.access(1, 1) = ((-1.+ 2.*y)*(-1.+ z)*z)/4.;
431 output.access(1, 2) = ((-1.+ y)*y*(-1.+ 2.*z))/4.;
432 output.access(1, 3) = ((1.+ 2.*x)*(-1.+ z)*z)/4.;
433 output.access(1, 4) = ((1.+ 2.*x)*(-1.+ 2.*y)*(-1.+ 2.*z))/8.;
434 output.access(1, 5) = ((1.+ 2.*x)*(-1.+ y)*y)/4.;
435 output.access(1, 6) = 0.;
436 output.access(1, 7) = (x*(1.+ x)*(-1.+ 2.*z))/4.;
437 output.access(1, 8) = (x*(1.+ x)*(-1.+ 2.*y))/4.;
438 output.access(1, 9) = 0.;
440 output.access(2, 0) = 0.;
441 output.access(2, 1) = ((1.+ 2.*y)*(-1.+ z)*z)/4.;
442 output.access(2, 2) = (y*(1.+ y)*(-1.+ 2.*z))/4.;
443 output.access(2, 3) = ((1.+ 2.*x)*(-1.+ z)*z)/4.;
444 output.access(2, 4) = ((1.+ 2.*x)*(1.+ 2.*y)*(-1.+ 2.*z))/8.;
445 output.access(2, 5) = ((1.+ 2.*x)*y*(1.+ y))/4.;
446 output.access(2, 6) = 0.;
447 output.access(2, 7) = (x*(1.+ x)*(-1.+ 2.*z))/4.;
448 output.access(2, 8) = (x*(1.+ x)*(1.+ 2.*y))/4.;
449 output.access(2, 9) = 0.;
451 output.access(3, 0) = 0.;
452 output.access(3, 1) = ((1.+ 2.*y)*(-1.+ z)*z)/4.;
453 output.access(3, 2) = (y*(1.+ y)*(-1.+ 2.*z))/4.;
454 output.access(3, 3) = ((-1.+ 2.*x)*(-1.+ z)*z)/4.;
455 output.access(3, 4) = ((-1.+ 2.*x)*(1.+ 2.*y)*(-1.+ 2.*z))/8.;
456 output.access(3, 5) = ((-1.+ 2.*x)*y*(1.+ y))/4.;
457 output.access(3, 6) = 0.;
458 output.access(3, 7) = ((-1.+ x)*x*(-1.+ 2.*z))/4.;
459 output.access(3, 8) = ((-1.+ x)*x*(1.+ 2.*y))/4.;
460 output.access(3, 9) = 0.;
462 output.access(4, 0) = 0.;
463 output.access(4, 1) = ((-1.+ 2.*y)*z*(1.+ z))/4.;
464 output.access(4, 2) = ((-1.+ y)*y*(1.+ 2.*z))/4.;
465 output.access(4, 3) = ((-1.+ 2.*x)*z*(1.+ z))/4.;
466 output.access(4, 4) = ((-1.+ 2.*x)*(-1.+ 2.*y)*(1.+ 2.*z))/8.;
467 output.access(4, 5) = ((-1.+ 2.*x)*(-1.+ y)*y)/4.;
468 output.access(4, 6) = 0.;
469 output.access(4, 7) = ((-1.+ x)*x*(1.+ 2.*z))/4.;
470 output.access(4, 8) = ((-1.+ x)*x*(-1.+ 2.*y))/4.;
471 output.access(4, 9) = 0.;
473 output.access(5, 0) = 0.;
474 output.access(5, 1) = ((-1.+ 2.*y)*z*(1.+ z))/4.;
475 output.access(5, 2) = ((-1.+ y)*y*(1.+ 2.*z))/4.;
476 output.access(5, 3) = ((1.+ 2.*x)*z*(1.+ z))/4.;
477 output.access(5, 4) = ((1.+ 2.*x)*(-1.+ 2.*y)*(1.+ 2.*z))/8.;
478 output.access(5, 5) = ((1.+ 2.*x)*(-1.+ y)*y)/4.;
479 output.access(5, 6) = 0.;
480 output.access(5, 7) = (x*(1.+ x)*(1.+ 2.*z))/4.;
481 output.access(5, 8) = (x*(1.+ x)*(-1.+ 2.*y))/4.;
482 output.access(5, 9) = 0.;
484 output.access(6, 0) = 0.;
485 output.access(6, 1) = ((1.+ 2.*y)*z*(1.+ z))/4.;
486 output.access(6, 2) = (y*(1.+ y)*(1.+ 2.*z))/4.;
487 output.access(6, 3) = ((1.+ 2.*x)*z*(1.+ z))/4.;
488 output.access(6, 4) = ((1.+ 2.*x)*(1.+ 2.*y)*(1.+ 2.*z))/8.;
489 output.access(6, 5) = ((1.+ 2.*x)*y*(1.+ y))/4.;
490 output.access(6, 6) = 0.;
491 output.access(6, 7) = (x*(1.+ x)*(1.+ 2.*z))/4.;
492 output.access(6, 8) = (x*(1.+ x)*(1.+ 2.*y))/4.;
493 output.access(6, 9) = 0.;
495 output.access(7, 0) = 0.;
496 output.access(7, 1) = ((1.+ 2.*y)*z*(1.+ z))/4.;
497 output.access(7, 2) = (y*(1.+ y)*(1.+ 2.*z))/4.;
498 output.access(7, 3) = ((-1.+ 2.*x)*z*(1.+ z))/4.;
499 output.access(7, 4) = ((-1.+ 2.*x)*(1.+ 2.*y)*(1.+ 2.*z))/8.;
500 output.access(7, 5) = ((-1.+ 2.*x)*y*(1.+ y))/4.;
501 output.access(7, 6) = 0.;
502 output.access(7, 7) = ((-1.+ x)*x*(1.+ 2.*z))/4.;
503 output.access(7, 8) = ((-1.+ x)*x*(1.+ 2.*y))/4.;
504 output.access(7, 9) = 0.;
506 output.access(8, 0) = 0.;
507 output.access(8, 1) = -((-1.+ 2.*y)*(-1.+ z)*z)/2.;
508 output.access(8, 2) = -((-1.+ y)*y*(-1.+ 2.*z))/2.;
509 output.access(8, 3) = -(x*(-1.+ z)*z);
510 output.access(8, 4) = -(x*(-1.+ 2.*y)*(-1.+ 2.*z))/2.;
511 output.access(8, 5) = -(x*(-1.+ y)*y);
512 output.access(8, 6) = 0.;
513 output.access(8, 7) = -((-1.+ (x*x))*(-1.+ 2.*z))/2.;
514 output.access(8, 8) = -((-1.+ (x*x))*(-1.+ 2.*y))/2.;
515 output.access(8, 9) = 0.;
517 output.access(9, 0) = 0.;
518 output.access(9, 1) = -(y*(-1.+ z)*z);
519 output.access(9, 2) = -((-1.+ (y*y))*(-1.+ 2.*z))/2.;
520 output.access(9, 3) = -((1.+ 2.*x)*(-1.+ z)*z)/2.;
521 output.access(9, 4) = -((1.+ 2.*x)*y*(-1.+ 2.*z))/2.;
522 output.access(9, 5) = -((1.+ 2.*x)*(-1.+ (y*y)))/2.;
523 output.access(9, 6) = 0.;
524 output.access(9, 7) = -(x*(1.+ x)*(-1.+ 2.*z))/2.;
525 output.access(9, 8) = -(x*(1.+ x)*y);
526 output.access(9, 9) = 0.;
528 output.access(10, 0) = 0.;
529 output.access(10, 1) = -((1.+ 2.*y)*(-1.+ z)*z)/2.;
530 output.access(10, 2) = -(y*(1.+ y)*(-1.+ 2.*z))/2.;
531 output.access(10, 3) = -(x*(-1.+ z)*z);
532 output.access(10, 4) = -(x*(1.+ 2.*y)*(-1.+ 2.*z))/2.;
533 output.access(10, 5) = -(x*y*(1.+ y));
534 output.access(10, 6) = 0.;
535 output.access(10, 7) = -((-1.+ (x*x))*(-1.+ 2.*z))/2.;
536 output.access(10, 8) = -((-1.+ (x*x))*(1.+ 2.*y))/2.;
537 output.access(10, 9) = 0.;
539 output.access(11, 0) = 0.;
540 output.access(11, 1) = -(y*(-1.+ z)*z);
541 output.access(11, 2) = -((-1.+ (y*y))*(-1.+ 2.*z))/2.;
542 output.access(11, 3) = -((-1.+ 2.*x)*(-1.+ z)*z)/2.;
543 output.access(11, 4) = -((-1.+ 2.*x)*y*(-1.+ 2.*z))/2.;
544 output.access(11, 5) = -((-1.+ 2.*x)*(-1.+ (y*y)))/2.;
545 output.access(11, 6) = 0.;
546 output.access(11, 7) = -((-1.+ x)*x*(-1.+ 2.*z))/2.;
547 output.access(11, 8) = -((-1.+ x)*x*y);
548 output.access(11, 9) = 0.;
550 output.access(12, 0) = 0.;
551 output.access(12, 1) = -((-1.+ 2.*y)*(-1.+ (z*z)))/2.;
552 output.access(12, 2) = -((-1.+ y)*y*z);
553 output.access(12, 3) = -((-1.+ 2.*x)*(-1.+ (z*z)))/2.;
554 output.access(12, 4) = -((-1.+ 2.*x)*(-1.+ 2.*y)*z)/2.;
555 output.access(12, 5) = -((-1.+ 2.*x)*(-1.+ y)*y)/2.;
556 output.access(12, 6) = 0.;
557 output.access(12, 7) = -((-1.+ x)*x*z);
558 output.access(12, 8) = -((-1.+ x)*x*(-1.+ 2.*y))/2.;
559 output.access(12, 9) = 0.;
561 output.access(13, 0) = 0.;
562 output.access(13, 1) = -((-1.+ 2.*y)*(-1.+ (z*z)))/2.;
563 output.access(13, 2) = -((-1.+ y)*y*z);
564 output.access(13, 3) = -((1.+ 2.*x)*(-1.+ (z*z)))/2.;
565 output.access(13, 4) = -((1.+ 2.*x)*(-1.+ 2.*y)*z)/2.;
566 output.access(13, 5) = -((1.+ 2.*x)*(-1.+ y)*y)/2.;
567 output.access(13, 6) = 0.;
568 output.access(13, 7) = -(x*(1.+ x)*z);
569 output.access(13, 8) = -(x*(1.+ x)*(-1.+ 2.*y))/2.;
570 output.access(13, 9) = 0.;
572 output.access(14, 0) = 0.;
573 output.access(14, 1) = -((1.+ 2.*y)*(-1.+ (z*z)))/2.;
574 output.access(14, 2) = -(y*(1.+ y)*z);
575 output.access(14, 3) = -((1.+ 2.*x)*(-1.+ (z*z)))/2.;
576 output.access(14, 4) = -((1.+ 2.*x)*(1.+ 2.*y)*z)/2.;
577 output.access(14, 5) = -((1.+ 2.*x)*y*(1.+ y))/2.;
578 output.access(14, 6) = 0.;
579 output.access(14, 7) = -(x*(1.+ x)*z);
580 output.access(14, 8) = -(x*(1.+ x)*(1.+ 2.*y))/2.;
581 output.access(14, 9) = 0.;
583 output.access(15, 0) = 0.;
584 output.access(15, 1) = -((1.+ 2.*y)*(-1.+ (z*z)))/2.;
585 output.access(15, 2) = -(y*(1.+ y)*z);
586 output.access(15, 3) = -((-1.+ 2.*x)*(-1.+ (z*z)))/2.;
587 output.access(15, 4) = -((-1.+ 2.*x)*(1.+ 2.*y)*z)/2.;
588 output.access(15, 5) = -((-1.+ 2.*x)*y*(1.+ y))/2.;
589 output.access(15, 6) = 0.;
590 output.access(15, 7) = -((-1.+ x)*x*z);
591 output.access(15, 8) = -((-1.+ x)*x*(1.+ 2.*y))/2.;
592 output.access(15, 9) = 0.;
594 output.access(16, 0) = 0.;
595 output.access(16, 1) = -((-1.+ 2.*y)*z*(1.+ z))/2.;
596 output.access(16, 2) = -((-1.+ y)*y*(1.+ 2.*z))/2.;
597 output.access(16, 3) = -(x*z*(1.+ z));
598 output.access(16, 4) = -(x*(-1.+ 2.*y)*(1.+ 2.*z))/2.;
599 output.access(16, 5) = -(x*(-1.+ y)*y);
600 output.access(16, 6) = 0.;
601 output.access(16, 7) = -((-1.+ (x*x))*(1.+ 2.*z))/2.;
602 output.access(16, 8) = -((-1.+ (x*x))*(-1.+ 2.*y))/2.;
603 output.access(16, 9) = 0.;
605 output.access(17, 0) = 0.;
606 output.access(17, 1) = -(y*z*(1.+ z));
607 output.access(17, 2) = -((-1.+ (y*y))*(1.+ 2.*z))/2.;
608 output.access(17, 3) = -((1.+ 2.*x)*z*(1.+ z))/2.;
609 output.access(17, 4) = -((1.+ 2.*x)*y*(1.+ 2.*z))/2.;
610 output.access(17, 5) = -((1.+ 2.*x)*(-1.+ (y*y)))/2.;
611 output.access(17, 6) = 0.;
612 output.access(17, 7) = -(x*(1.+ x)*(1.+ 2.*z))/2.;
613 output.access(17, 8) = -(x*(1.+ x)*y);
614 output.access(17, 9) = 0.;
616 output.access(18, 0) = 0.;
617 output.access(18, 1) = -((1.+ 2.*y)*z*(1.+ z))/2.;
618 output.access(18, 2) = -(y*(1.+ y)*(1.+ 2.*z))/2.;
619 output.access(18, 3) = -(x*z*(1.+ z));
620 output.access(18, 4) = -(x*(1.+ 2.*y)*(1.+ 2.*z))/2.;
621 output.access(18, 5) = -(x*y*(1.+ y));
622 output.access(18, 6) = 0.;
623 output.access(18, 7) = -((-1.+ (x*x))*(1.+ 2.*z))/2.;
624 output.access(18, 8) = -((-1.+ (x*x))*(1.+ 2.*y))/2.;
625 output.access(18, 9) = 0.;
627 output.access(19, 0) = 0.;
628 output.access(19, 1) = -(y*z*(1.+ z));
629 output.access(19, 2) = -((-1.+ (y*y))*(1.+ 2.*z))/2.;
630 output.access(19, 3) = -((-1.+ 2.*x)*z*(1.+ z))/2.;
631 output.access(19, 4) = -((-1.+ 2.*x)*y*(1.+ 2.*z))/2.;
632 output.access(19, 5) = -((-1.+ 2.*x)*(-1.+ (y*y)))/2.;
633 output.access(19, 6) = 0.;
634 output.access(19, 7) = -((-1.+ x)*x*(1.+ 2.*z))/2.;
635 output.access(19, 8) = -((-1.+ x)*x*y);
636 output.access(19, 9) = 0.;
638 output.access(20, 0) = 0.;
639 output.access(20, 1) = -4*y*(-1.+ (z*z));
640 output.access(20, 2) = -4*(-1.+ (y*y))*z;
641 output.access(20, 3) = -4*x*(-1.+ (z*z));
642 output.access(20, 4) = -8*x*y*z;
643 output.access(20, 5) = -4*x*(-1.+ (y*y));
644 output.access(20, 6) = 0.;
645 output.access(20, 7) = -4*(-1.+ (x*x))*z;
646 output.access(20, 8) = -4*(-1.+ (x*x))*y;
647 output.access(20, 9) = 0.;
649 output.access(21, 0) = 0.;
650 output.access(21, 1) = 2.*y*(-1.+ z)*z;
651 output.access(21, 2) = (-1.+ (y*y))*(-1.+ 2.*z);
652 output.access(21, 3) = 2.*x*(-1.+ z)*z;
653 output.access(21, 4) = 2.*x*y*(-1.+ 2.*z);
654 output.access(21, 5) = 2.*x*(-1.+ (y*y));
655 output.access(21, 6) = 0.;
656 output.access(21, 7) = (-1.+ (x*x))*(-1.+ 2.*z);
657 output.access(21, 8) = 2.*(-1.+ (x*x))*y;
658 output.access(21, 9) = 0.;
660 output.access(22, 0) = 0.;
661 output.access(22, 1) = 2.*y*z*(1.+ z);
662 output.access(22, 2) = (-1.+ (y*y))*(1.+ 2.*z);
663 output.access(22, 3) = 2.*x*z*(1.+ z);
664 output.access(22, 4) = 2.*x*y*(1.+ 2.*z);
665 output.access(22, 5) = 2.*x*(-1.+ (y*y));
666 output.access(22, 6) = 0.;
667 output.access(22, 7) = (-1.+ (x*x))*(1.+ 2.*z);
668 output.access(22, 8) = 2.*(-1.+ (x*x))*y;
669 output.access(22, 9) = 0.;
671 output.access(23, 0) = 0.;
672 output.access(23, 1) = 2.*y*(-1.+ (z*z));
673 output.access(23, 2) = 2.*(-1.+ (y*y))*z;
674 output.access(23, 3) = (-1.+ 2.*x)*(-1.+ (z*z));
675 output.access(23, 4) = 2.*(-1.+ 2.*x)*y*z;
676 output.access(23, 5) = (-1.+ 2.*x)*(-1.+ (y*y));
677 output.access(23, 6) = 0.;
678 output.access(23, 7) = 2.*(-1.+ x)*x*z;
679 output.access(23, 8) = 2.*(-1.+ x)*x*y;
680 output.access(23, 9) = 0.;
682 output.access(24, 0) = 0.;
683 output.access(24, 1) = 2.*y*(-1.+ (z*z));
684 output.access(24, 2) = 2.*(-1.+ (y*y))*z;
685 output.access(24, 3) = (1.+ 2.*x)*(-1.+ (z*z));
686 output.access(24, 4) = 2.*(1.+ 2.*x)*y*z;
687 output.access(24, 5) = (1.+ 2.*x)*(-1.+ (y*y));
688 output.access(24, 6) = 0.;
689 output.access(24, 7) = 2.*x*(1.+ x)*z;
690 output.access(24, 8) = 2.*x*(1.+ x)*y;
691 output.access(24, 9) = 0.;
693 output.access(25, 0) = 0.;
694 output.access(25, 1) = (-1.+ 2.*y)*(-1.+ (z*z));
695 output.access(25, 2) = 2.*(-1.+ y)*y*z;
696 output.access(25, 3) = 2.*x*(-1.+ (z*z));
697 output.access(25, 4) = 2.*x*(-1.+ 2.*y)*z;
698 output.access(25, 5) = 2.*x*(-1.+ y)*y;
699 output.access(25, 6) = 0.;
700 output.access(25, 7) = 2.*(-1.+ (x*x))*z;
701 output.access(25, 8) = (-1.+ (x*x))*(-1.+ 2.*y);
702 output.access(25, 9) = 0.;
704 output.access(26, 0) = 0.;
705 output.access(26, 1) = (1.+ 2.*y)*(-1.+ (z*z));
706 output.access(26, 2) = 2.*y*(1.+ y)*z;
707 output.access(26, 3) = 2.*x*(-1.+ (z*z));
708 output.access(26, 4) = 2.*x*(1.+ 2.*y)*z;
709 output.access(26, 5) = 2.*x*y*(1.+ y);
710 output.access(26, 6) = 0.;
711 output.access(26, 7) = 2.*(-1.+ (x*x))*z;
712 output.access(26, 8) = (-1.+ (x*x))*(1.+ 2.*y);
713 output.access(26, 9) = 0.;
719 const ordinal_type jend = output.extent(1);
720 const ordinal_type iend = output.extent(0);
722 for (ordinal_type j=0;j<jend;++j)
723 for (ordinal_type i=0;i<iend;++i)
724 output.access(i, j) = 0.0;
726 const auto x = input(0);
727 const auto y = input(1);
728 const auto z = input(2);
729 output.access(0, 3) = ((-1.+ z)*z)/2.;
730 output.access(0, 4) = ((-1.+ 2.*y)*(-1.+ 2.*z))/4.;
731 output.access(0, 5) = ((-1.+ y)*y)/2.;
732 output.access(0, 7) = ((-1.+ 2.*x)*(-1.+ 2.*z))/4.;
733 output.access(0, 8) = ((-1.+ 2.*x)*(-1.+ 2.*y))/4.;
734 output.access(0, 12)= ((-1.+ x)*x)/2.;
736 output.access(1, 3) = ((-1.+ z)*z)/2.;
737 output.access(1, 4) = ((-1.+ 2.*y)*(-1.+ 2.*z))/4.;
738 output.access(1, 5) = ((-1.+ y)*y)/2.;
739 output.access(1, 7) = ((1. + 2.*x)*(-1.+ 2.*z))/4.;
740 output.access(1, 8) = ((1. + 2.*x)*(-1.+ 2.*y))/4.;
741 output.access(1, 12)= (x*(1. + x))/2.;
743 output.access(2, 3) = ((-1.+ z)*z)/2.;
744 output.access(2, 4) = ((1. + 2.*y)*(-1.+ 2.*z))/4.;
745 output.access(2, 5) = (y*(1. + y))/2.;
746 output.access(2, 7) = ((1. + 2.*x)*(-1.+ 2.*z))/4.;
747 output.access(2, 8) = ((1. + 2.*x)*(1. + 2.*y))/4.;
748 output.access(2, 12)= (x*(1. + x))/2.;
750 output.access(3, 3) = ((-1.+ z)*z)/2.;
751 output.access(3, 4) = ((1. + 2.*y)*(-1.+ 2.*z))/4.;
752 output.access(3, 5) = (y*(1. + y))/2.;
753 output.access(3, 7) = ((-1.+ 2.*x)*(-1.+ 2.*z))/4.;
754 output.access(3, 8) = ((-1.+ 2.*x)*(1. + 2.*y))/4.;
755 output.access(3, 12)= ((-1.+ x)*x)/2.;
757 output.access(4, 3) = (z*(1. + z))/2.;
758 output.access(4, 4) = ((-1.+ 2.*y)*(1. + 2.*z))/4.;
759 output.access(4, 5) = ((-1.+ y)*y)/2.;
760 output.access(4, 7) = ((-1.+ 2.*x)*(1. + 2.*z))/4.;
761 output.access(4, 8) = ((-1.+ 2.*x)*(-1.+ 2.*y))/4.;
762 output.access(4, 12)= ((-1.+ x)*x)/2.;
764 output.access(5, 3) = (z*(1. + z))/2.;
765 output.access(5, 4) = ((-1.+ 2.*y)*(1. + 2.*z))/4.;
766 output.access(5, 5) = ((-1.+ y)*y)/2.;
767 output.access(5, 7) = ((1. + 2.*x)*(1. + 2.*z))/4.;
768 output.access(5, 8) = ((1. + 2.*x)*(-1.+ 2.*y))/4.;
769 output.access(5, 12)= (x*(1. + x))/2.;
771 output.access(6, 3) = (z*(1. + z))/2.;
772 output.access(6, 4) = ((1. + 2.*y)*(1. + 2.*z))/4.;
773 output.access(6, 5) = (y*(1. + y))/2.;
774 output.access(6, 7) = ((1. + 2.*x)*(1. + 2.*z))/4.;
775 output.access(6, 8) = ((1. + 2.*x)*(1. + 2.*y))/4.;
776 output.access(6, 12)= (x*(1. + x))/2.;
778 output.access(7, 3) = (z*(1. + z))/2.;
779 output.access(7, 4) = ((1. + 2.*y)*(1. + 2.*z))/4.;
780 output.access(7, 5) = (y*(1. + y))/2.;
781 output.access(7, 7) = ((-1.+ 2.*x)*(1. + 2.*z))/4.;
782 output.access(7, 8) = ((-1.+ 2.*x)*(1. + 2.*y))/4.;
783 output.access(7, 12)= ((-1.+ x)*x)/2.;
785 output.access(8, 3) = -((-1.+ z)*z);
786 output.access(8, 4) = -0.5 + y + z - 2.*y*z;
787 output.access(8, 5) = -((-1.+ y)*y);
788 output.access(8, 7) = x - 2.*x*z;
789 output.access(8, 8) = x - 2.*x*y;
790 output.access(8, 12)= 1. - x*x;
792 output.access(9, 3) = -((-1.+ z)*z);
793 output.access(9, 4) = y - 2.*y*z;
794 output.access(9, 5) = 1 - y*y;
795 output.access(9, 7) = 0.5 + x - z - 2.*x*z;
796 output.access(9, 8) = -((1. + 2.*x)*y);
797 output.access(9, 12)= -(x*(1. + x));
799 output.access(10, 3) = -((-1.+ z)*z);
800 output.access(10, 4) = 0.5 + y - z - 2.*y*z;
801 output.access(10, 5) = -(y*(1. + y));
802 output.access(10, 7) = x - 2.*x*z;
803 output.access(10, 8) = -(x*(1. + 2.*y));
804 output.access(10, 12)= 1. - x*x;
806 output.access(11, 3) = -((-1.+ z)*z);
807 output.access(11, 4) = y - 2.*y*z;
808 output.access(11, 5) = 1. - y*y;
809 output.access(11, 7) = -0.5 + x + z - 2.*x*z;
810 output.access(11, 8) = y - 2.*x*y;
811 output.access(11, 12)= -((-1.+ x)*x);
813 output.access(12, 3) = 1. - z*z;
814 output.access(12, 4) = z - 2.*y*z;
815 output.access(12, 5) = -((-1.+ y)*y);
816 output.access(12, 7) = z - 2.*x*z;
817 output.access(12, 8) = -0.5 + x + y - 2.*x*y;
818 output.access(12, 12)= -((-1.+ x)*x);
820 output.access(13, 3) = 1. - z*z;
821 output.access(13, 4) = z - 2.*y*z;
822 output.access(13, 5) = -((-1.+ y)*y);
823 output.access(13, 7) = -((1. + 2.*x)*z);
824 output.access(13, 8) = 0.5 + x - y - 2.*x*y;
825 output.access(13, 12)= -(x*(1. + x));
827 output.access(14, 3) = 1. - z*z;
828 output.access(14, 4) = -((1. + 2.*y)*z);
829 output.access(14, 5) = -(y*(1. + y));
830 output.access(14, 7) = -((1. + 2.*x)*z);
831 output.access(14, 8) = -((1. + 2.*x)*(1. + 2.*y))/2.;
832 output.access(14, 12)= -(x*(1. + x));
834 output.access(15, 3) = 1. - z*z;
835 output.access(15, 4) = -((1. + 2.*y)*z);
836 output.access(15, 5) = -(y*(1. + y));
837 output.access(15, 7) = z - 2.*x*z;
838 output.access(15, 8) = 0.5 + y - x*(1. + 2.*y);
839 output.access(15, 12)= -((-1.+ x)*x);
841 output.access(16, 3) = -(z*(1. + z));
842 output.access(16, 4) = 0.5 + z - y*(1. + 2.*z);
843 output.access(16, 5) = -((-1.+ y)*y);
844 output.access(16, 7) = -(x*(1. + 2.*z));
845 output.access(16, 8) = x - 2.*x*y;
846 output.access(16, 12)= 1. - x*x;
848 output.access(17, 3) = -(z*(1. + z));
849 output.access(17, 4) = -(y*(1. + 2.*z));
850 output.access(17, 5) = 1. - y*y;
851 output.access(17, 7) = -((1. + 2.*x)*(1. + 2.*z))/2.;
852 output.access(17, 8) = -((1. + 2.*x)*y);
853 output.access(17, 12)= -(x*(1. + x));
855 output.access(18, 3) = -(z*(1. + z));
856 output.access(18, 4) = -((1. + 2.*y)*(1. + 2.*z))/2.;
857 output.access(18, 5) = -(y*(1. + y));
858 output.access(18, 7) = -(x*(1. + 2.*z));
859 output.access(18, 8) = -(x*(1. + 2.*y));
860 output.access(18, 12)= 1. - x*x;
862 output.access(19, 3) = -(z*(1. + z));
863 output.access(19, 4) = -(y*(1. + 2.*z));
864 output.access(19, 5) = 1. - y*y;
865 output.access(19, 7) = 0.5 + z - x*(1. + 2.*z);
866 output.access(19, 8) = y - 2.*x*y;
867 output.access(19, 12)= -((-1.+ x)*x);
869 output.access(20, 3) = 4. - 4.*z*z;
870 output.access(20, 4) = -8.*y*z;
871 output.access(20, 5) = 4. - 4.*y*y;
872 output.access(20, 7) = -8.*x*z;
873 output.access(20, 8) = -8.*x*y;
874 output.access(20, 12)= 4. - 4.*x*x;
876 output.access(21, 3) = 2.*(-1.+ z)*z;
877 output.access(21, 4) = 2.*y*(-1.+ 2.*z);
878 output.access(21, 5) = 2.*(-1.+ y*y);
879 output.access(21, 7) = 2.*x*(-1.+ 2.*z);
880 output.access(21, 8) = 4.*x*y;
881 output.access(21, 12)= 2.*(-1.+ x*x);
883 output.access(22, 3) = 2.*z*(1. + z);
884 output.access(22, 4) = 2.*(y + 2.*y*z);
885 output.access(22, 5) = 2.*(-1.+ y*y);
886 output.access(22, 7) = 2.*(x + 2.*x*z);
887 output.access(22, 8) = 4.*x*y;
888 output.access(22, 12)= 2.*(-1.+ x*x);
890 output.access(23, 3) = 2.*(-1.+ z*z);
891 output.access(23, 4) = 4.*y*z;
892 output.access(23, 5) = 2.*(-1.+ y*y);
893 output.access(23, 7) = 2.*(-1.+ 2.*x)*z;
894 output.access(23, 8) = 2.*(-1.+ 2.*x)*y;
895 output.access(23, 12)= 2.*(-1.+ x)*x;
897 output.access(24, 3) = 2.*(-1.+ z*z);
898 output.access(24, 4) = 4.*y*z;
899 output.access(24, 5) = 2.*(-1.+ y*y);
900 output.access(24, 7) = 2.*(z + 2.*x*z);
901 output.access(24, 8) = 2.*(y + 2.*x*y);
902 output.access(24, 12)= 2.*x*(1. + x);
904 output.access(25, 3) = 2.*(-1.+ z*z);
905 output.access(25, 4) = 2.*(-1.+ 2.*y)*z;
906 output.access(25, 5) = 2.*(-1.+ y)*y;
907 output.access(25, 7) = 4.*x*z;
908 output.access(25, 8) = 2.*x*(-1.+ 2.*y);
909 output.access(25, 12)= 2.*(-1.+ x*x);
911 output.access(26, 3) = 2.*(-1.+ z*z);
912 output.access(26, 4) = 2.*(z + 2.*y*z);
913 output.access(26, 5) = 2.*y*(1. + y);
914 output.access(26, 7) = 4.*x*z;
915 output.access(26, 8) = 2.*(x + 2.*x*y);
916 output.access(26, 12)= 2.*(-1.+ x*x);
920 case OPERATOR_MAX : {
921 const ordinal_type jend = output.extent(1);
922 const ordinal_type iend = output.extent(0);
924 for (ordinal_type j=0;j<jend;++j)
925 for (ordinal_type i=0;i<iend;++i)
926 output.access(i, j) = 0.0;
930 INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
931 opType != OPERATOR_GRAD &&
932 opType != OPERATOR_CURL &&
933 opType != OPERATOR_D1 &&
934 opType != OPERATOR_D2 &&
935 opType != OPERATOR_D3 &&
936 opType != OPERATOR_D4 &&
937 opType != OPERATOR_MAX,
938 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_C2_FEM::Serial::getValues) operator is not supported");
944 template<
typename SpT,
945 typename outputValueValueType,
class ...outputValueProperties,
946 typename inputPointValueType,
class ...inputPointProperties>
948 Basis_HGRAD_HEX_C2_FEM::
949 getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
950 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
951 const EOperator operatorType ) {
952 typedef Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValueViewType;
953 typedef Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPointViewType;
954 typedef typename ExecSpace<typename inputPointViewType::execution_space,SpT>::ExecSpaceType ExecSpaceType;
957 const auto loopSize = inputPoints.extent(0);
958 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
960 switch (operatorType) {
962 case OPERATOR_VALUE: {
963 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_VALUE> FunctorType;
964 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
969 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_GRAD> FunctorType;
970 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
974 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument,
975 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D");
979 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
980 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D");
984 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D2> FunctorType;
985 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
989 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D3> FunctorType;
990 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
994 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D4> FunctorType;
995 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
1000 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
1001 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): operator not supported");
1007 case OPERATOR_D10: {
1008 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_MAX> FunctorType;
1009 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
1013 INTREPID2_TEST_FOR_EXCEPTION( !( Intrepid2::isValidOperator(operatorType) ), std::invalid_argument,
1014 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): Invalid operator type");
1023 template<
typename SpT,
typename OT,
typename PT>
1026 this -> basisCardinality_ = 27;
1027 this -> basisDegree_ = 2;
1028 this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Hexahedron<8> >() );
1029 this -> basisType_ = BASIS_FEM_DEFAULT;
1030 this -> basisCoordinates_ = COORDINATES_CARTESIAN;
1035 const ordinal_type tagSize = 4;
1036 const ordinal_type posScDim = 0;
1037 const ordinal_type posScOrd = 1;
1038 const ordinal_type posDfOrd = 2;
1041 ordinal_type tags[108] = { 0, 0, 0, 1,
1074 this->setOrdinalTagData(
this -> tagToOrdinal_,
1075 this -> ordinalToTag_,
1077 this -> basisCardinality_,
1084 Kokkos::DynRankView<typename scalarViewType::value_type,typename SpT::array_layout,Kokkos::HostSpace>
1085 dofCoords(
"dofCoordsHost", this->basisCardinality_,this->basisCellTopology_.getDimension());
1087 dofCoords(0,0) = -1.0; dofCoords(0,1) = -1.0; dofCoords(0,2) = -1.0;
1088 dofCoords(1,0) = 1.0; dofCoords(1,1) = -1.0; dofCoords(1,2) = -1.0;
1089 dofCoords(2,0) = 1.0; dofCoords(2,1) = 1.0; dofCoords(2,2) = -1.0;
1090 dofCoords(3,0) = -1.0; dofCoords(3,1) = 1.0; dofCoords(3,2) = -1.0;
1091 dofCoords(4,0) = -1.0; dofCoords(4,1) = -1.0; dofCoords(4,2) = 1.0;
1092 dofCoords(5,0) = 1.0; dofCoords(5,1) = -1.0; dofCoords(5,2) = 1.0;
1093 dofCoords(6,0) = 1.0; dofCoords(6,1) = 1.0; dofCoords(6,2) = 1.0;
1094 dofCoords(7,0) = -1.0; dofCoords(7,1) = 1.0; dofCoords(7,2) = 1.0;
1096 dofCoords(8,0) = 0.0; dofCoords(8,1) = -1.0; dofCoords(8,2) = -1.0;
1097 dofCoords(9,0) = 1.0; dofCoords(9,1) = 0.0; dofCoords(9,2) = -1.0;
1098 dofCoords(10,0) = 0.0; dofCoords(10,1) = 1.0; dofCoords(10,2) = -1.0;
1099 dofCoords(11,0) = -1.0; dofCoords(11,1) = 0.0; dofCoords(11,2) = -1.0;
1100 dofCoords(12,0) = -1.0; dofCoords(12,1) = -1.0; dofCoords(12,2) = 0.0;
1101 dofCoords(13,0) = 1.0; dofCoords(13,1) = -1.0; dofCoords(13,2) = 0.0;
1102 dofCoords(14,0) = 1.0; dofCoords(14,1) = 1.0; dofCoords(14,2) = 0.0;
1103 dofCoords(15,0) = -1.0; dofCoords(15,1) = 1.0; dofCoords(15,2) = 0.0;
1104 dofCoords(16,0) = 0.0; dofCoords(16,1) = -1.0; dofCoords(16,2) = 1.0;
1105 dofCoords(17,0) = 1.0; dofCoords(17,1) = 0.0; dofCoords(17,2) = 1.0;
1106 dofCoords(18,0) = 0.0; dofCoords(18,1) = 1.0; dofCoords(18,2) = 1.0;
1107 dofCoords(19,0) = -1.0; dofCoords(19,1) = 0.0; dofCoords(19,2) = 1.0;
1109 dofCoords(20,0) = 0.0; dofCoords(20,1) = 0.0; dofCoords(20,2) = 0.0;
1111 dofCoords(21,0) = 0.0; dofCoords(21,1) = 0.0; dofCoords(21,2) = -1.0;
1112 dofCoords(22,0) = 0.0; dofCoords(22,1) = 0.0; dofCoords(22,2) = 1.0;
1113 dofCoords(23,0) = -1.0; dofCoords(23,1) = 0.0; dofCoords(23,2) = 0.0;
1114 dofCoords(24,0) = 1.0; dofCoords(24,1) = 0.0; dofCoords(24,2) = 0.0;
1115 dofCoords(25,0) = 0.0; dofCoords(25,1) = -1.0; dofCoords(25,2) = 0.0;
1116 dofCoords(26,0) = 0.0; dofCoords(26,1) = 1.0; dofCoords(26,2) = 0.0;
1118 this->dofCoords_ = Kokkos::create_mirror_view(
typename SpT::memory_space(), dofCoords);
1119 Kokkos::deep_copy(this->dofCoords_, dofCoords);
Basis_HGRAD_HEX_C2_FEM()
Constructor.
Kokkos::View< ordinal_type *,typename ExecSpaceType::array_layout, Kokkos::HostSpace > ordinal_type_array_1d_host
View type for 1d host array.