Intrepid2
Intrepid2_HGRAD_HEX_C2_FEMDef.hpp
Go to the documentation of this file.
1// @HEADER
2// ************************************************************************
3//
4// Intrepid2 Package
5// Copyright (2007) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38// Mauro Perego (mperego@sandia.gov)
39//
40// ************************************************************************
41// @HEADER
42
48
49#ifndef __INTREPID2_HGRAD_HEX_C2_FEM_DEF_HPP__
50#define __INTREPID2_HGRAD_HEX_C2_FEM_DEF_HPP__
51
52namespace Intrepid2 {
53
54 // -------------------------------------------------------------------------------------
55 namespace Impl {
56
57 template<bool serendipity>
58 template<EOperator opType>
59 template<typename OutputViewType,
60 typename inputViewType>
61 KOKKOS_INLINE_FUNCTION
62 void
64 getValues( OutputViewType output,
65 const inputViewType input ) {
66 switch (opType) {
67 case OPERATOR_VALUE : {
68 const auto x = input(0);
69 const auto y = input(1);
70 const auto z = input(2);
71
72 // output is a rank-2 array with dimensions (basisCardinality_, dim0)
73 output.access( 0) = 0.125*(-1. + x)*x*(-1. + y)*y*(-1. + z)*z;
74 output.access( 1) = 0.125*x*(1.+ x)*(-1. + y)*y*(-1. + z)*z;
75 output.access( 2) = 0.125*x*(1.+ x)*y*(1.+ y)*(-1. + z)*z;
76 output.access( 3) = 0.125*(-1. + x)*x*y*(1.+ y)*(-1. + z)*z;
77 output.access( 4) = 0.125*(-1. + x)*x*(-1. + y)*y*z*(1.+ z);
78 output.access( 5) = 0.125*x*(1.+ x)*(-1. + y)*y*z*(1.+ z);
79 output.access( 6) = 0.125*x*(1.+ x)*y*(1.+ y)*z*(1.+ z);
80 output.access( 7) = 0.125*(-1. + x)*x*y*(1.+ y)*z*(1.+ z);
81 output.access( 8) = 0.25*(1. - x)*(1. + x)*(-1. + y)*y*(-1. + z)*z;
82 output.access( 9) = 0.25*x*(1.+ x)*(1. - y)*(1. + y)*(-1. + z)*z;
83 output.access(10) = 0.25*(1. - x)*(1. + x)*y*(1.+ y)*(-1. + z)*z;
84 output.access(11) = 0.25*(-1. + x)*x*(1. - y)*(1. + y)*(-1. + z)*z;
85 output.access(12) = 0.25*(-1. + x)*x*(-1. + y)*y*(1. - z)*(1. + z);
86 output.access(13) = 0.25*x*(1.+ x)*(-1. + y)*y*(1. - z)*(1. + z);
87 output.access(14) = 0.25*x*(1.+ x)*y*(1.+ y)*(1. - z)*(1. + z);
88 output.access(15) = 0.25*(-1. + x)*x*y*(1.+ y)*(1. - z)*(1. + z);
89 output.access(16) = 0.25*(1. - x)*(1. + x)*(-1. + y)*y*z*(1.+ z);
90 output.access(17) = 0.25*x*(1.+ x)*(1. - y)*(1. + y)*z*(1.+ z);
91 output.access(18) = 0.25*(1. - x)*(1. + x)*y*(1.+ y)*z*(1.+ z);
92 output.access(19) = 0.25*(-1. + x)*x*(1. - y)*(1. + y)*z*(1.+ z);
93 if constexpr (!serendipity) {
94 output.access(20) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
95 output.access(21) = 0.5*(1. - x)*(1. + x)*(1. - y)*(1. + y)*(-1. + z)*z;
96 output.access(22) = 0.5*(1. - x)*(1. + x)*(1. - y)*(1. + y)*z*(1.+ z);
97 output.access(23) = 0.5*(-1. + x)*x*(1. - y)*(1. + y)*(1. - z)*(1. + z);
98 output.access(24) = 0.5*x*(1.+ x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
99 output.access(25) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y*(1. - z)*(1. + z);
100 output.access(26) = 0.5*(1. - x)*(1. + x)*y*(1.+ y)*(1. - z)*(1. + z);
101 }
102 break;
103 }
104 case OPERATOR_GRAD :
105 case OPERATOR_D1 : {
106 const auto x = input(0);
107 const auto y = input(1);
108 const auto z = input(2);
109
110 // output.access is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
111 output.access(0, 0) = (-0.125 + 0.25*x)*(-1. + y)*y*(-1. + z)*z;
112 output.access(0, 1) = (-1. + x)*x*(-0.125 + 0.25*y)*(-1. + z)*z;
113 output.access(0, 2) = (-1. + x)*x*(-1. + y)*y*(-0.125 + 0.25*z);
114
115 output.access(1, 0) = (0.125 + 0.25*x)*(-1. + y)*y*(-1. + z)*z;
116 output.access(1, 1) = x*(1. + x)*(-0.125 + 0.25*y)*(-1. + z)*z;
117 output.access(1, 2) = x*(1. + x)*(-1. + y)*y*(-0.125 + 0.25*z);
118
119 output.access(2, 0) = (0.125 + 0.25*x)*y*(1. + y)*(-1. + z)*z;
120 output.access(2, 1) = x*(1. + x)*(0.125 + 0.25*y)*(-1. + z)*z;
121 output.access(2, 2) = x*(1. + x)*y*(1. + y)*(-0.125 + 0.25*z);
122
123 output.access(3, 0) = (-0.125 + 0.25*x)*y*(1. + y)*(-1. + z)*z;
124 output.access(3, 1) = (-1. + x)*x*(0.125 + 0.25*y)*(-1. + z)*z;
125 output.access(3, 2) = (-1. + x)*x*y*(1. + y)*(-0.125 + 0.25*z);
126
127 output.access(4, 0) = (-0.125 + 0.25*x)*(-1. + y)*y*z*(1. + z);
128 output.access(4, 1) = (-1. + x)*x*(-0.125 + 0.25*y)*z*(1. + z);
129 output.access(4, 2) = (-1. + x)*x*(-1. + y)*y*(0.125 + 0.25*z);
130
131 output.access(5, 0) = (0.125 + 0.25*x)*(-1. + y)*y*z*(1. + z);
132 output.access(5, 1) = x*(1. + x)*(-0.125 + 0.25*y)*z*(1. + z);
133 output.access(5, 2) = x*(1. + x)*(-1. + y)*y*(0.125 + 0.25*z);
134
135 output.access(6, 0) = (0.125 + 0.25*x)*y*(1. + y)*z*(1. + z);
136 output.access(6, 1) = x*(1. + x)*(0.125 + 0.25*y)*z*(1. + z);
137 output.access(6, 2) = x*(1. + x)*y*(1. + y)*(0.125 + 0.25*z);
138
139 output.access(7, 0) = (-0.125 + 0.25*x)*y*(1. + y)*z*(1. + z);
140 output.access(7, 1) = (-1. + x)*x*(0.125 + 0.25*y)*z*(1. + z);
141 output.access(7, 2) = (-1. + x)*x*y*(1. + y)*(0.125 + 0.25*z);
142
143 output.access(8, 0) = -0.5*x*(-1. + y)*y*(-1. + z)*z;
144 output.access(8, 1) = (1. - x)*(1. + x)*(-0.25 + 0.5*y)*(-1. + z)*z;
145 output.access(8, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(-0.25 + 0.5*z);
146
147 output.access(9, 0) = (0.25 + 0.5*x)*(1. - y)*(1. + y)*(-1. + z)*z;
148 output.access(9, 1) = x*(1. + x)*(-0.5*y)*(-1. + z)*z;
149 output.access(9, 2) = x*(1. + x)*(1. - y)*(1. + y)*(-0.25 + 0.5*z);
150
151 output.access(10, 0) = -0.5*x*y*(1. + y)*(-1. + z)*z;
152 output.access(10, 1) = (1. - x)*(1. + x)*(0.25 + 0.5*y)*(-1. + z)*z;
153 output.access(10, 2) = (1. - x)*(1. + x)*y*(1. + y)*(-0.25 + 0.5*z);
154
155 output.access(11, 0) = (-0.25 + 0.5*x)*(1. - y)*(1. + y)*(-1. + z)*z;
156 output.access(11, 1) = (-1. + x)*x*(-0.5*y)*(-1. + z)*z;
157 output.access(11, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(-0.25 + 0.5*z);
158
159 output.access(12, 0) = (-0.25 + 0.5*x)*(-1. + y)*y*(1. - z)*(1. + z);
160 output.access(12, 1) = (-1. + x)*x*(-0.25 + 0.5*y)*(1. - z)*(1. + z);
161 output.access(12, 2) = (-1. + x)*x*(-1. + y)*y*(-0.5*z);
162
163 output.access(13, 0) = (0.25 + 0.5*x)*(-1. + y)*y*(1. - z)*(1. + z);
164 output.access(13, 1) = x*(1. + x)*(-0.25 + 0.5*y)*(1. - z)*(1. + z);
165 output.access(13, 2) = x*(1. + x)*(-1. + y)*y*(-0.5*z);
166
167 output.access(14, 0) = (0.25 + 0.5*x)*y*(1. + y)*(1. - z)*(1. + z);
168 output.access(14, 1) = x*(1. + x)*(0.25 + 0.5*y)*(1. - z)*(1. + z);
169 output.access(14, 2) = x*(1. + x)*y*(1. + y)*(-0.5*z);
170
171 output.access(15, 0) = (-0.25 + 0.5*x)*y*(1. + y)*(1. - z)*(1. + z);
172 output.access(15, 1) = (-1. + x)*x*(0.25 + 0.5*y)*(1. - z)*(1. + z);
173 output.access(15, 2) = (-1. + x)*x*y*(1. + y)*(-0.5*z);
174
175 output.access(16, 0) = -0.5*x*(-1. + y)*y*z*(1. + z);
176 output.access(16, 1) = (1. - x)*(1. + x)*(-0.25 + 0.5*y)*z*(1. + z);
177 output.access(16, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(0.25 + 0.5*z);
178
179 output.access(17, 0) = (0.25 + 0.5*x)*(1. - y)*(1. + y)*z*(1. + z);
180 output.access(17, 1) = x*(1. + x)*(-0.5*y)*z*(1. + z);
181 output.access(17, 2) = x*(1. + x)*(1. - y)*(1. + y)*(0.25 + 0.5*z);
182
183 output.access(18, 0) = -0.5*x*y*(1. + y)*z*(1. + z);
184 output.access(18, 1) = (1. - x)*(1. + x)*(0.25 + 0.5*y)*z*(1. + z);
185 output.access(18, 2) = (1. - x)*(1. + x)*y*(1. + y)*(0.25 + 0.5*z);
186
187 output.access(19, 0) = (-0.25 + 0.5*x)*(1. - y)*(1. + y)*z*(1. + z);
188 output.access(19, 1) = (-1. + x)*x*(-0.5*y)*z*(1. + z);
189 output.access(19, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(0.25 + 0.5*z);
190
191 if constexpr (!serendipity) {
192 output.access(20, 0) = -2.*x*(1. - y)*(1. + y)*(1. - z)*(1. + z);
193 output.access(20, 1) = (1. - x)*(1. + x)*(-2.*y)*(1. - z)*(1. + z);
194 output.access(20, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(-2.*z);
195
196 output.access(21, 0) = -x*(1. - y)*(1. + y)*(-1. + z)*z;
197 output.access(21, 1) = (1. - x)*(1. + x)*(-y)*(-1. + z)*z;
198 output.access(21, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(-0.5 + z);
199
200 output.access(22, 0) = -x*(1. - y)*(1. + y)*z*(1. + z);
201 output.access(22, 1) = (1. - x)*(1. + x)*(-y)*z*(1. + z);
202 output.access(22, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(0.5 + z);
203
204 output.access(23, 0) = (-0.5 + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
205 output.access(23, 1) = (-1. + x)*x*(-y)*(1. - z)*(1. + z);
206 output.access(23, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(-z);
207
208 output.access(24, 0) = (0.5 + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
209 output.access(24, 1) = x*(1. + x)*(-y)*(1. - z)*(1. + z);
210 output.access(24, 2) = x*(1. + x)*(1. - y)*(1. + y)*(-z);
211
212 output.access(25, 0) = -x*(-1. + y)*y*(1. - z)*(1. + z);
213 output.access(25, 1) = (1. - x)*(1. + x)*(-0.5 + y)*(1. - z)*(1. + z);
214 output.access(25, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(-z);
215
216 output.access(26, 0) = -x*y*(1. + y)*(1. - z)*(1. + z);
217 output.access(26, 1) = (1. - x)*(1. + x)*(0.5 + y)*(1. - z)*(1. + z);
218 output.access(26, 2) = (1. - x)*(1. + x)*y*(1. + y)*(-z);
219 }
220 break;
221 }
222 case OPERATOR_D2 : {
223 const auto x = input(0);
224 const auto y = input(1);
225 const auto z = input(2);
226
227 // output.access is a rank-3 array with dimensions (basisCardinality_, dim0, D2Cardinality = 6)
228 output.access(0, 0) = 0.25*(-1. + y)*y*(-1. + z)*z;
229 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;
230 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);
231 output.access(0, 3) = 0.25*(-1. + x)*x*(-1. + z)*z;
232 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);
233 output.access(0, 5) = 0.25*(-1. + x)*x*(-1. + y)*y;
234
235 output.access(1, 0) = 0.25*(-1. + y)*y*(-1. + z)*z;
236 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;
237 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);
238 output.access(1, 3) = 0.25*x*(1 + x)*(-1. + z)*z;
239 output.access(1, 4) = x*(1. + x)*(0.125 + y*(-0.25 + 0.5*z) - 0.25*z);
240 output.access(1, 5) = 0.25*x*(1 + x)*(-1. + y)*y;
241
242 output.access(2, 0) = 0.25*y*(1 + y)*(-1. + z)*z;
243 output.access(2, 1) = (0.125 + x*(0.25 + 0.5*y) + 0.25*y)*(-1. + z)*z;
244 output.access(2, 2) = y*(1. + y)*(-0.125 + x*(-0.25 + 0.5*z) + 0.25*z);
245 output.access(2, 3) = 0.25*x*(1 + x)*(-1. + z)*z;
246 output.access(2, 4) = x*(1. + x)*(-0.125 + y*(-0.25 + 0.5*z) + 0.25*z);
247 output.access(2, 5) = 0.25*x*(1 + x)*y*(1 + y);
248
249 output.access(3, 0) = 0.25*y*(1 + y)*(-1. + z)*z;
250 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;
251 output.access(3, 2) = y*(1. + y)*(0.125 + x*(-0.25 + 0.5*z) - 0.25*z);
252 output.access(3, 3) = 0.25*(-1. + x)*x*(-1. + z)*z;
253 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);
254 output.access(3, 5) = 0.25*(-1. + x)*x*y*(1 + y);
255
256 output.access(4, 0) = 0.25*(-1. + y)*y*z*(1 + z);
257 output.access(4, 1) = (0.125 + x*(-0.25 + 0.5*y) - 0.25*y)*z*(1. + z);
258 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);
259 output.access(4, 3) = 0.25*(-1. + x)*x*z*(1 + z);
260 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);
261 output.access(4, 5) = 0.25*(-1. + x)*x*(-1. + y)*y;
262
263 output.access(5, 0) = 0.25*(-1. + y)*y*z*(1 + z);
264 output.access(5, 1) = (-0.125 + x*(-0.25 + 0.5*y) + 0.25*y)*z*(1. + z);
265 output.access(5, 2) = (-1. + y)*y*(0.125 + x*(0.25 + 0.5*z) + 0.25*z);
266 output.access(5, 3) = 0.25*x*(1 + x)*z*(1 + z);
267 output.access(5, 4) = x*(1. + x)*(-0.125 + y*(0.25 + 0.5*z) - 0.25*z);
268 output.access(5, 5) = 0.25*x*(1 + x)*(-1. + y)*y;
269
270 output.access(6, 0) = 0.25*y*(1 + y)*z*(1 + z);
271 output.access(6, 1) = (0.125 + x*(0.25 + 0.5*y) + 0.25*y)*z*(1. + z);
272 output.access(6, 2) = y*(1. + y)*(0.125 + x*(0.25 + 0.5*z) + 0.25*z);
273 output.access(6, 3) = 0.25*x*(1 + x)*z*(1 + z);
274 output.access(6, 4) = x*(1. + x)*(0.125 + y*(0.25 + 0.5*z) + 0.25*z);
275 output.access(6, 5) = 0.25*x*(1 + x)*y*(1 + y);
276
277 output.access(7, 0) = 0.25*y*(1 + y)*z*(1 + z);
278 output.access(7, 1) = (-0.125 + x*(0.25 + 0.5*y) - 0.25*y)*z*(1. + z);
279 output.access(7, 2) = y*(1. + y)*(-0.125 + x*(0.25 + 0.5*z) - 0.25*z);
280 output.access(7, 3) = 0.25*(-1. + x)*x*z*(1 + z);
281 output.access(7, 4) = (-1. + x)*x*(0.125 + y*(0.25 + 0.5*z) + 0.25*z);
282 output.access(7, 5) = 0.25*(-1. + x)*x*y*(1 + y);
283
284 output.access(8, 0) = -0.5*(-1. + y)*y*(-1. + z)*z;
285 output.access(8, 1) = (0. + x*(-0.5 + y))*z + (x*(0.5 - y) )*(z*z);
286 output.access(8, 2) = (y*y)*(x*(0.5 - z) ) + y*(x*(-0.5 + z));
287 output.access(8, 3) = 0.5*(1. - x)*(1. + x)*(-1. + z)*z;
288 output.access(8, 4) = 0.25 + (x*x)*(-0.25 + y*(0.5 - z) + 0.5*z) - 0.5*z + y*(-0.5 + z);
289 output.access(8, 5) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y;
290
291 output.access(9, 0) = 0.5*(1. - y)*(1. + y)*(-1. + z)*z;
292 output.access(9, 1) = (0.5*y + x*(y))*z + (x*(-y) - 0.5*y)*(z*z);
293 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);
294 output.access(9, 3) = -0.5*x*(1 + x)*(-1. + z)*z;
295 output.access(9, 4) = x*(y*(0.5 - z) ) + (x*x)*(y*(0.5 - z) );
296 output.access(9, 5) = 0.5*x*(1 + x)*(1. - y)*(1. + y);
297
298 output.access(10, 0) = -0.5*y*(1 + y)*(-1. + z)*z;
299 output.access(10, 1) = (0. + x*(0.5 + y))*z + (x*(-0.5 - y) )*(z*z);
300 output.access(10, 2) = y*(x*(0.5 - z) ) + (y*y)*(x*(0.5 - z) );
301 output.access(10, 3) = 0.5*(1. - x)*(1. + x)*(-1. + z)*z;
302 output.access(10, 4) = -0.25 + (x*x)*(0.25 + y*(0.5 - z) - 0.5*z) + 0.5*z + y*(-0.5 + z);
303 output.access(10, 5) = 0.5*(1. - x)*(1. + x)*y*(1 + y);
304
305 output.access(11, 0) = 0.5*(1. - y)*(1. + y)*(-1. + z)*z;
306 output.access(11, 1) = (-0.5*y + x*(y))*z + (x*(-y) + 0.5*y)*(z*z);
307 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);
308 output.access(11, 3) = -0.5*(-1. + x)*x*(-1. + z)*z;
309 output.access(11, 4) = (x*x)*(y*(0.5 - z) ) + x*(y*(-0.5 + z));
310 output.access(11, 5) = 0.5*(-1. + x)*x*(1. - y)*(1. + y);
311
312 output.access(12, 0) = 0.5*(-1. + y)*y*(1. - z)*(1. + z);
313 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)));
314 output.access(12, 2) = (y*y)*(x*(-z) + 0.5*z) + y*(-0.5*z + x*(z));
315 output.access(12, 3) = 0.5*(-1. + x)*x*(1. - z)*(1. + z);
316 output.access(12, 4) = (x*x)*(y*(-z) + 0.5*z) + x*(-0.5*z + y*(z));
317 output.access(12, 5) = -0.5*(-1. + x)*x*(-1. + y)*y;
318
319 output.access(13, 0) = 0.5*(-1. + y)*y*(1. - z)*(1. + z);
320 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)));
321 output.access(13, 2) = (y*y)*(x*(-z) - 0.5*z) + y*(0.5*z + x*(z));
322 output.access(13, 3) = 0.5*x*(1 + x)*(1. - z)*(1. + z);
323 output.access(13, 4) = x*(y*(-z) + 0.5*z) + (x*x)*(y*(-z) + 0.5*z);
324 output.access(13, 5) = -0.5*x*(1 + x)*(-1. + y)*y;
325
326 output.access(14, 0) = 0.5*y*(1 + y)*(1. - z)*(1. + z);
327 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)));
328 output.access(14, 2) = y*(x*(-z) - 0.5*z) + (y*y)*(x*(-z) - 0.5*z);
329 output.access(14, 3) = 0.5*x*(1 + x)*(1. - z)*(1. + z);
330 output.access(14, 4) = x*(y*(-z) - 0.5*z) + (x*x)*(y*(-z) - 0.5*z);
331 output.access(14, 5) = -0.5*x*(1 + x)*y*(1 + y);
332
333 output.access(15, 0) = 0.5*y*(1 + y)*(1. - z)*(1. + z);
334 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)));
335 output.access(15, 2) = y*(x*(-z) + 0.5*z) + (y*y)*(x*(-z) + 0.5*z);
336 output.access(15, 3) = 0.5*(-1. + x)*x*(1. - z)*(1. + z);
337 output.access(15, 4) = (x*x)*(y*(-z) - 0.5*z) + x*(0.5*z + y*(z));
338 output.access(15, 5) = -0.5*(-1. + x)*x*y*(1 + y);
339
340 output.access(16, 0) = -0.5*(-1. + y)*y*z*(1 + z);
341 output.access(16, 1) = (x*(0.5 - y) )*z + (x*(0.5 - y) )*(z*z);
342 output.access(16, 2) = (y*y)*(x*(-0.5 - z) ) + y*(x*(0.5 + z));
343 output.access(16, 3) = 0.5*(1. - x)*(1. + x)*z*(1 + z);
344 output.access(16, 4) = -0.25 + (x*x)*(0.25 + y*(-0.5 - z) + 0.5*z) - 0.5*z + y*(0.5 + z);
345 output.access(16, 5) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y;
346
347 output.access(17, 0) = 0.5*(1. - y)*(1. + y)*z*(1 + z);
348 output.access(17, 1) = (x*(-y) - 0.5*y)*z + (x*(-y) - 0.5*y)*(z*z);
349 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);
350 output.access(17, 3) = -0.5*x*(1 + x)*z*(1 + z);
351 output.access(17, 4) = x*(y*(-0.5 - z) ) + (x*x)*(y*(-0.5 - z) );
352 output.access(17, 5) = 0.5*x*(1 + x)*(1. - y)*(1. + y);
353
354 output.access(18, 0) = -0.5*y*(1 + y)*z*(1 + z);
355 output.access(18, 1) = (x*(-0.5 - y) )*z + (x*(-0.5 - y) )*(z*z);
356 output.access(18, 2) = y*(x*(-0.5 - z) ) + (y*y)*(x*(-0.5 - z) );
357 output.access(18, 3) = 0.5*(1. - x)*(1. + x)*z*(1 + z);
358 output.access(18, 4) = 0.25 + (x*x)*(-0.25 + y*(-0.5 - z) - 0.5*z) + 0.5*z + y*(0.5 + z);
359 output.access(18, 5) = 0.5*(1. - x)*(1. + x)*y*(1 + y);
360
361 output.access(19, 0) = 0.5*(1. - y)*(1. + y)*z*(1 + z);
362 output.access(19, 1) = (x*(-y) + 0.5*y)*z + (x*(-y) + 0.5*y)*(z*z);
363 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);
364 output.access(19, 3) = -0.5*(-1. + x)*x*z*(1 + z);
365 output.access(19, 4) = (x*x)*(y*(-0.5 - z) ) + x*(y*(0.5 + z));
366 output.access(19, 5) = 0.5*(-1. + x)*x*(1. - y)*(1. + y);
367
368 if constexpr (!serendipity) {
369 output.access(20, 0) = -2.*(1. - y)*(1. + y)*(1. - z)*(1. + z);
370 output.access(20, 1) = -4.*x*y*(-1. + z*z);
371 output.access(20, 2) = x*((y*y)*(-4.*z) + 4.*z);
372 output.access(20, 3) = -2.*(1. - x)*(1. + x)*(1. - z)*(1. + z);
373 output.access(20, 4) = (x*x)*(y*(-4.*z) ) + y*(4.*z);
374 output.access(20, 5) = -2.*(1. - x)*(1. + x)*(1. - y)*(1. + y);
375
376 output.access(21, 0) = -(1. - y)*(1. + y)*(-1. + z)*z;
377 output.access(21, 1) = (x*(-2.*y) )*z + (0. + x*(2.*y))*(z*z);
378 output.access(21, 2) = x*(1. - 2.*z + (y*y)*(-1. + 2.*z));
379 output.access(21, 3) = -(1. - x)*(1. + x)*(-1. + z)*z;
380 output.access(21, 4) = y*(1. - 2.*z) + (x*x)*(y*(-1. + 2.*z));
381 output.access(21, 5) = (1. - x)*(1. + x)*(1. - y)*(1. + y);
382
383 output.access(22, 0) = -(1. - y)*(1. + y)*z*(1 + z);
384 output.access(22, 1) = (0. + x*(2.*y))*z + (0. + x*(2.*y))*(z*z);
385 output.access(22, 2) = x*(-1. - 2.*z + (y*y)*(1. + 2.*z));
386 output.access(22, 3) = -(1. - x)*(1. + x)*z*(1 + z);
387 output.access(22, 4) = y*(-1. - 2.*z) + (x*x)*(y*(1. + 2.*z));
388 output.access(22, 5) = (1. - x)*(1. + x)*(1. - y)*(1. + y);
389
390 output.access(23, 0) = (1. - y)*(1. + y)*(1. - z)*(1. + z);
391 output.access(23, 1) = (-1. + 2.*x)*y*(-1. + z*z);
392 output.access(23, 2) = (-1. + 2.*x)*(-1. + y*y)*z;
393 output.access(23, 3) =-(-1. + x)*x*(1. - z)*(1. + z);
394 output.access(23, 4) = 2.*(-1. + x)*x*y*z;
395 output.access(23, 5) =-(-1. + x)*x*(1. - y)*(1. + y);
396
397 output.access(24, 0) = (1. - y)*(1. + y)*(1. - z)*(1. + z);
398 output.access(24, 1) = (1. + 2.*x)*y*(-1. + z*z);
399 output.access(24, 2) = (1. + 2.*x)*(-1. + y*y)*z;
400 output.access(24, 3) = x*(1. + x)*(-1. + z)*(1. + z);
401 output.access(24, 4) = 2.*x*(1. + x)*y*z;
402 output.access(24, 5) = x*(1. + x)*(-1. + y)*(1. + y);
403
404 output.access(25, 0) = -(-1. + y)*y*(1. - z)*(1. + z);
405 output.access(25, 1) = x*(-1. + 2.*y)*(-1. + z*z);
406 output.access(25, 2) = 2.*x*(-1. + y)*y*z;
407 output.access(25, 3) = (1. - x)*(1. + x)*(1. - z)*(1. + z);
408 output.access(25, 4) = (-1. + x*x)*(-1. + 2.*y)*z;
409 output.access(25, 5) =-(1. - x)*(1. + x)*(-1. + y)*y;
410
411 output.access(26, 0) = y*(1. + y)*(-1. + z)*(1. + z);
412 output.access(26, 1) = x*(1. + 2.*y)*(-1. + z*z);
413 output.access(26, 2) = 2.*x*y*(1. + y)*z;
414 output.access(26, 3) = (-1. + x)*(1. + x)*(-1. + z)*(1. + z);
415 output.access(26, 4) = (-1. + x*x)*(1. + 2.*y)*z;
416 output.access(26, 5) = (-1. + x)*(1. + x)*y*(1. + y);
417 }
418 break;
419 }
420 case OPERATOR_D3 : {
421 const auto x = input(0);
422 const auto y = input(1);
423 const auto z = input(2);
424
425 output.access(0, 0) = 0.;
426 output.access(0, 1) = ((-1.+ 2.*y)*(-1.+ z)*z)/4.;
427 output.access(0, 2) = ((-1.+ y)*y*(-1.+ 2.*z))/4.;
428 output.access(0, 3) = ((-1.+ 2.*x)*(-1.+ z)*z)/4.;
429 output.access(0, 4) = ((-1.+ 2.*x)*(-1.+ 2.*y)*(-1.+ 2.*z))/8.;
430 output.access(0, 5) = ((-1.+ 2.*x)*(-1.+ y)*y)/4.;
431 output.access(0, 6) = 0.;
432 output.access(0, 7) = ((-1.+ x)*x*(-1.+ 2.*z))/4.;
433 output.access(0, 8) = ((-1.+ x)*x*(-1.+ 2.*y))/4.;
434 output.access(0, 9) = 0.;
435
436 output.access(1, 0) = 0.;
437 output.access(1, 1) = ((-1.+ 2.*y)*(-1.+ z)*z)/4.;
438 output.access(1, 2) = ((-1.+ y)*y*(-1.+ 2.*z))/4.;
439 output.access(1, 3) = ((1.+ 2.*x)*(-1.+ z)*z)/4.;
440 output.access(1, 4) = ((1.+ 2.*x)*(-1.+ 2.*y)*(-1.+ 2.*z))/8.;
441 output.access(1, 5) = ((1.+ 2.*x)*(-1.+ y)*y)/4.;
442 output.access(1, 6) = 0.;
443 output.access(1, 7) = (x*(1.+ x)*(-1.+ 2.*z))/4.;
444 output.access(1, 8) = (x*(1.+ x)*(-1.+ 2.*y))/4.;
445 output.access(1, 9) = 0.;
446
447 output.access(2, 0) = 0.;
448 output.access(2, 1) = ((1.+ 2.*y)*(-1.+ z)*z)/4.;
449 output.access(2, 2) = (y*(1.+ y)*(-1.+ 2.*z))/4.;
450 output.access(2, 3) = ((1.+ 2.*x)*(-1.+ z)*z)/4.;
451 output.access(2, 4) = ((1.+ 2.*x)*(1.+ 2.*y)*(-1.+ 2.*z))/8.;
452 output.access(2, 5) = ((1.+ 2.*x)*y*(1.+ y))/4.;
453 output.access(2, 6) = 0.;
454 output.access(2, 7) = (x*(1.+ x)*(-1.+ 2.*z))/4.;
455 output.access(2, 8) = (x*(1.+ x)*(1.+ 2.*y))/4.;
456 output.access(2, 9) = 0.;
457
458 output.access(3, 0) = 0.;
459 output.access(3, 1) = ((1.+ 2.*y)*(-1.+ z)*z)/4.;
460 output.access(3, 2) = (y*(1.+ y)*(-1.+ 2.*z))/4.;
461 output.access(3, 3) = ((-1.+ 2.*x)*(-1.+ z)*z)/4.;
462 output.access(3, 4) = ((-1.+ 2.*x)*(1.+ 2.*y)*(-1.+ 2.*z))/8.;
463 output.access(3, 5) = ((-1.+ 2.*x)*y*(1.+ y))/4.;
464 output.access(3, 6) = 0.;
465 output.access(3, 7) = ((-1.+ x)*x*(-1.+ 2.*z))/4.;
466 output.access(3, 8) = ((-1.+ x)*x*(1.+ 2.*y))/4.;
467 output.access(3, 9) = 0.;
468
469 output.access(4, 0) = 0.;
470 output.access(4, 1) = ((-1.+ 2.*y)*z*(1.+ z))/4.;
471 output.access(4, 2) = ((-1.+ y)*y*(1.+ 2.*z))/4.;
472 output.access(4, 3) = ((-1.+ 2.*x)*z*(1.+ z))/4.;
473 output.access(4, 4) = ((-1.+ 2.*x)*(-1.+ 2.*y)*(1.+ 2.*z))/8.;
474 output.access(4, 5) = ((-1.+ 2.*x)*(-1.+ y)*y)/4.;
475 output.access(4, 6) = 0.;
476 output.access(4, 7) = ((-1.+ x)*x*(1.+ 2.*z))/4.;
477 output.access(4, 8) = ((-1.+ x)*x*(-1.+ 2.*y))/4.;
478 output.access(4, 9) = 0.;
479
480 output.access(5, 0) = 0.;
481 output.access(5, 1) = ((-1.+ 2.*y)*z*(1.+ z))/4.;
482 output.access(5, 2) = ((-1.+ y)*y*(1.+ 2.*z))/4.;
483 output.access(5, 3) = ((1.+ 2.*x)*z*(1.+ z))/4.;
484 output.access(5, 4) = ((1.+ 2.*x)*(-1.+ 2.*y)*(1.+ 2.*z))/8.;
485 output.access(5, 5) = ((1.+ 2.*x)*(-1.+ y)*y)/4.;
486 output.access(5, 6) = 0.;
487 output.access(5, 7) = (x*(1.+ x)*(1.+ 2.*z))/4.;
488 output.access(5, 8) = (x*(1.+ x)*(-1.+ 2.*y))/4.;
489 output.access(5, 9) = 0.;
490
491 output.access(6, 0) = 0.;
492 output.access(6, 1) = ((1.+ 2.*y)*z*(1.+ z))/4.;
493 output.access(6, 2) = (y*(1.+ y)*(1.+ 2.*z))/4.;
494 output.access(6, 3) = ((1.+ 2.*x)*z*(1.+ z))/4.;
495 output.access(6, 4) = ((1.+ 2.*x)*(1.+ 2.*y)*(1.+ 2.*z))/8.;
496 output.access(6, 5) = ((1.+ 2.*x)*y*(1.+ y))/4.;
497 output.access(6, 6) = 0.;
498 output.access(6, 7) = (x*(1.+ x)*(1.+ 2.*z))/4.;
499 output.access(6, 8) = (x*(1.+ x)*(1.+ 2.*y))/4.;
500 output.access(6, 9) = 0.;
501
502 output.access(7, 0) = 0.;
503 output.access(7, 1) = ((1.+ 2.*y)*z*(1.+ z))/4.;
504 output.access(7, 2) = (y*(1.+ y)*(1.+ 2.*z))/4.;
505 output.access(7, 3) = ((-1.+ 2.*x)*z*(1.+ z))/4.;
506 output.access(7, 4) = ((-1.+ 2.*x)*(1.+ 2.*y)*(1.+ 2.*z))/8.;
507 output.access(7, 5) = ((-1.+ 2.*x)*y*(1.+ y))/4.;
508 output.access(7, 6) = 0.;
509 output.access(7, 7) = ((-1.+ x)*x*(1.+ 2.*z))/4.;
510 output.access(7, 8) = ((-1.+ x)*x*(1.+ 2.*y))/4.;
511 output.access(7, 9) = 0.;
512
513 output.access(8, 0) = 0.;
514 output.access(8, 1) = -((-1.+ 2.*y)*(-1.+ z)*z)/2.;
515 output.access(8, 2) = -((-1.+ y)*y*(-1.+ 2.*z))/2.;
516 output.access(8, 3) = -(x*(-1.+ z)*z);
517 output.access(8, 4) = -(x*(-1.+ 2.*y)*(-1.+ 2.*z))/2.;
518 output.access(8, 5) = -(x*(-1.+ y)*y);
519 output.access(8, 6) = 0.;
520 output.access(8, 7) = -((-1.+ (x*x))*(-1.+ 2.*z))/2.;
521 output.access(8, 8) = -((-1.+ (x*x))*(-1.+ 2.*y))/2.;
522 output.access(8, 9) = 0.;
523
524 output.access(9, 0) = 0.;
525 output.access(9, 1) = -(y*(-1.+ z)*z);
526 output.access(9, 2) = -((-1.+ (y*y))*(-1.+ 2.*z))/2.;
527 output.access(9, 3) = -((1.+ 2.*x)*(-1.+ z)*z)/2.;
528 output.access(9, 4) = -((1.+ 2.*x)*y*(-1.+ 2.*z))/2.;
529 output.access(9, 5) = -((1.+ 2.*x)*(-1.+ (y*y)))/2.;
530 output.access(9, 6) = 0.;
531 output.access(9, 7) = -(x*(1.+ x)*(-1.+ 2.*z))/2.;
532 output.access(9, 8) = -(x*(1.+ x)*y);
533 output.access(9, 9) = 0.;
534
535 output.access(10, 0) = 0.;
536 output.access(10, 1) = -((1.+ 2.*y)*(-1.+ z)*z)/2.;
537 output.access(10, 2) = -(y*(1.+ y)*(-1.+ 2.*z))/2.;
538 output.access(10, 3) = -(x*(-1.+ z)*z);
539 output.access(10, 4) = -(x*(1.+ 2.*y)*(-1.+ 2.*z))/2.;
540 output.access(10, 5) = -(x*y*(1.+ y));
541 output.access(10, 6) = 0.;
542 output.access(10, 7) = -((-1.+ (x*x))*(-1.+ 2.*z))/2.;
543 output.access(10, 8) = -((-1.+ (x*x))*(1.+ 2.*y))/2.;
544 output.access(10, 9) = 0.;
545
546 output.access(11, 0) = 0.;
547 output.access(11, 1) = -(y*(-1.+ z)*z);
548 output.access(11, 2) = -((-1.+ (y*y))*(-1.+ 2.*z))/2.;
549 output.access(11, 3) = -((-1.+ 2.*x)*(-1.+ z)*z)/2.;
550 output.access(11, 4) = -((-1.+ 2.*x)*y*(-1.+ 2.*z))/2.;
551 output.access(11, 5) = -((-1.+ 2.*x)*(-1.+ (y*y)))/2.;
552 output.access(11, 6) = 0.;
553 output.access(11, 7) = -((-1.+ x)*x*(-1.+ 2.*z))/2.;
554 output.access(11, 8) = -((-1.+ x)*x*y);
555 output.access(11, 9) = 0.;
556
557 output.access(12, 0) = 0.;
558 output.access(12, 1) = -((-1.+ 2.*y)*(-1.+ (z*z)))/2.;
559 output.access(12, 2) = -((-1.+ y)*y*z);
560 output.access(12, 3) = -((-1.+ 2.*x)*(-1.+ (z*z)))/2.;
561 output.access(12, 4) = -((-1.+ 2.*x)*(-1.+ 2.*y)*z)/2.;
562 output.access(12, 5) = -((-1.+ 2.*x)*(-1.+ y)*y)/2.;
563 output.access(12, 6) = 0.;
564 output.access(12, 7) = -((-1.+ x)*x*z);
565 output.access(12, 8) = -((-1.+ x)*x*(-1.+ 2.*y))/2.;
566 output.access(12, 9) = 0.;
567
568 output.access(13, 0) = 0.;
569 output.access(13, 1) = -((-1.+ 2.*y)*(-1.+ (z*z)))/2.;
570 output.access(13, 2) = -((-1.+ y)*y*z);
571 output.access(13, 3) = -((1.+ 2.*x)*(-1.+ (z*z)))/2.;
572 output.access(13, 4) = -((1.+ 2.*x)*(-1.+ 2.*y)*z)/2.;
573 output.access(13, 5) = -((1.+ 2.*x)*(-1.+ y)*y)/2.;
574 output.access(13, 6) = 0.;
575 output.access(13, 7) = -(x*(1.+ x)*z);
576 output.access(13, 8) = -(x*(1.+ x)*(-1.+ 2.*y))/2.;
577 output.access(13, 9) = 0.;
578
579 output.access(14, 0) = 0.;
580 output.access(14, 1) = -((1.+ 2.*y)*(-1.+ (z*z)))/2.;
581 output.access(14, 2) = -(y*(1.+ y)*z);
582 output.access(14, 3) = -((1.+ 2.*x)*(-1.+ (z*z)))/2.;
583 output.access(14, 4) = -((1.+ 2.*x)*(1.+ 2.*y)*z)/2.;
584 output.access(14, 5) = -((1.+ 2.*x)*y*(1.+ y))/2.;
585 output.access(14, 6) = 0.;
586 output.access(14, 7) = -(x*(1.+ x)*z);
587 output.access(14, 8) = -(x*(1.+ x)*(1.+ 2.*y))/2.;
588 output.access(14, 9) = 0.;
589
590 output.access(15, 0) = 0.;
591 output.access(15, 1) = -((1.+ 2.*y)*(-1.+ (z*z)))/2.;
592 output.access(15, 2) = -(y*(1.+ y)*z);
593 output.access(15, 3) = -((-1.+ 2.*x)*(-1.+ (z*z)))/2.;
594 output.access(15, 4) = -((-1.+ 2.*x)*(1.+ 2.*y)*z)/2.;
595 output.access(15, 5) = -((-1.+ 2.*x)*y*(1.+ y))/2.;
596 output.access(15, 6) = 0.;
597 output.access(15, 7) = -((-1.+ x)*x*z);
598 output.access(15, 8) = -((-1.+ x)*x*(1.+ 2.*y))/2.;
599 output.access(15, 9) = 0.;
600
601 output.access(16, 0) = 0.;
602 output.access(16, 1) = -((-1.+ 2.*y)*z*(1.+ z))/2.;
603 output.access(16, 2) = -((-1.+ y)*y*(1.+ 2.*z))/2.;
604 output.access(16, 3) = -(x*z*(1.+ z));
605 output.access(16, 4) = -(x*(-1.+ 2.*y)*(1.+ 2.*z))/2.;
606 output.access(16, 5) = -(x*(-1.+ y)*y);
607 output.access(16, 6) = 0.;
608 output.access(16, 7) = -((-1.+ (x*x))*(1.+ 2.*z))/2.;
609 output.access(16, 8) = -((-1.+ (x*x))*(-1.+ 2.*y))/2.;
610 output.access(16, 9) = 0.;
611
612 output.access(17, 0) = 0.;
613 output.access(17, 1) = -(y*z*(1.+ z));
614 output.access(17, 2) = -((-1.+ (y*y))*(1.+ 2.*z))/2.;
615 output.access(17, 3) = -((1.+ 2.*x)*z*(1.+ z))/2.;
616 output.access(17, 4) = -((1.+ 2.*x)*y*(1.+ 2.*z))/2.;
617 output.access(17, 5) = -((1.+ 2.*x)*(-1.+ (y*y)))/2.;
618 output.access(17, 6) = 0.;
619 output.access(17, 7) = -(x*(1.+ x)*(1.+ 2.*z))/2.;
620 output.access(17, 8) = -(x*(1.+ x)*y);
621 output.access(17, 9) = 0.;
622
623 output.access(18, 0) = 0.;
624 output.access(18, 1) = -((1.+ 2.*y)*z*(1.+ z))/2.;
625 output.access(18, 2) = -(y*(1.+ y)*(1.+ 2.*z))/2.;
626 output.access(18, 3) = -(x*z*(1.+ z));
627 output.access(18, 4) = -(x*(1.+ 2.*y)*(1.+ 2.*z))/2.;
628 output.access(18, 5) = -(x*y*(1.+ y));
629 output.access(18, 6) = 0.;
630 output.access(18, 7) = -((-1.+ (x*x))*(1.+ 2.*z))/2.;
631 output.access(18, 8) = -((-1.+ (x*x))*(1.+ 2.*y))/2.;
632 output.access(18, 9) = 0.;
633
634 output.access(19, 0) = 0.;
635 output.access(19, 1) = -(y*z*(1.+ z));
636 output.access(19, 2) = -((-1.+ (y*y))*(1.+ 2.*z))/2.;
637 output.access(19, 3) = -((-1.+ 2.*x)*z*(1.+ z))/2.;
638 output.access(19, 4) = -((-1.+ 2.*x)*y*(1.+ 2.*z))/2.;
639 output.access(19, 5) = -((-1.+ 2.*x)*(-1.+ (y*y)))/2.;
640 output.access(19, 6) = 0.;
641 output.access(19, 7) = -((-1.+ x)*x*(1.+ 2.*z))/2.;
642 output.access(19, 8) = -((-1.+ x)*x*y);
643 output.access(19, 9) = 0.;
644
645 if constexpr(!serendipity) {
646 output.access(20, 0) = 0.;
647 output.access(20, 1) = -4*y*(-1.+ (z*z));
648 output.access(20, 2) = -4*(-1.+ (y*y))*z;
649 output.access(20, 3) = -4*x*(-1.+ (z*z));
650 output.access(20, 4) = -8*x*y*z;
651 output.access(20, 5) = -4*x*(-1.+ (y*y));
652 output.access(20, 6) = 0.;
653 output.access(20, 7) = -4*(-1.+ (x*x))*z;
654 output.access(20, 8) = -4*(-1.+ (x*x))*y;
655 output.access(20, 9) = 0.;
656
657 output.access(21, 0) = 0.;
658 output.access(21, 1) = 2.*y*(-1.+ z)*z;
659 output.access(21, 2) = (-1.+ (y*y))*(-1.+ 2.*z);
660 output.access(21, 3) = 2.*x*(-1.+ z)*z;
661 output.access(21, 4) = 2.*x*y*(-1.+ 2.*z);
662 output.access(21, 5) = 2.*x*(-1.+ (y*y));
663 output.access(21, 6) = 0.;
664 output.access(21, 7) = (-1.+ (x*x))*(-1.+ 2.*z);
665 output.access(21, 8) = 2.*(-1.+ (x*x))*y;
666 output.access(21, 9) = 0.;
667
668 output.access(22, 0) = 0.;
669 output.access(22, 1) = 2.*y*z*(1.+ z);
670 output.access(22, 2) = (-1.+ (y*y))*(1.+ 2.*z);
671 output.access(22, 3) = 2.*x*z*(1.+ z);
672 output.access(22, 4) = 2.*x*y*(1.+ 2.*z);
673 output.access(22, 5) = 2.*x*(-1.+ (y*y));
674 output.access(22, 6) = 0.;
675 output.access(22, 7) = (-1.+ (x*x))*(1.+ 2.*z);
676 output.access(22, 8) = 2.*(-1.+ (x*x))*y;
677 output.access(22, 9) = 0.;
678
679 output.access(23, 0) = 0.;
680 output.access(23, 1) = 2.*y*(-1.+ (z*z));
681 output.access(23, 2) = 2.*(-1.+ (y*y))*z;
682 output.access(23, 3) = (-1.+ 2.*x)*(-1.+ (z*z));
683 output.access(23, 4) = 2.*(-1.+ 2.*x)*y*z;
684 output.access(23, 5) = (-1.+ 2.*x)*(-1.+ (y*y));
685 output.access(23, 6) = 0.;
686 output.access(23, 7) = 2.*(-1.+ x)*x*z;
687 output.access(23, 8) = 2.*(-1.+ x)*x*y;
688 output.access(23, 9) = 0.;
689
690 output.access(24, 0) = 0.;
691 output.access(24, 1) = 2.*y*(-1.+ (z*z));
692 output.access(24, 2) = 2.*(-1.+ (y*y))*z;
693 output.access(24, 3) = (1.+ 2.*x)*(-1.+ (z*z));
694 output.access(24, 4) = 2.*(1.+ 2.*x)*y*z;
695 output.access(24, 5) = (1.+ 2.*x)*(-1.+ (y*y));
696 output.access(24, 6) = 0.;
697 output.access(24, 7) = 2.*x*(1.+ x)*z;
698 output.access(24, 8) = 2.*x*(1.+ x)*y;
699 output.access(24, 9) = 0.;
700
701 output.access(25, 0) = 0.;
702 output.access(25, 1) = (-1.+ 2.*y)*(-1.+ (z*z));
703 output.access(25, 2) = 2.*(-1.+ y)*y*z;
704 output.access(25, 3) = 2.*x*(-1.+ (z*z));
705 output.access(25, 4) = 2.*x*(-1.+ 2.*y)*z;
706 output.access(25, 5) = 2.*x*(-1.+ y)*y;
707 output.access(25, 6) = 0.;
708 output.access(25, 7) = 2.*(-1.+ (x*x))*z;
709 output.access(25, 8) = (-1.+ (x*x))*(-1.+ 2.*y);
710 output.access(25, 9) = 0.;
711
712 output.access(26, 0) = 0.;
713 output.access(26, 1) = (1.+ 2.*y)*(-1.+ (z*z));
714 output.access(26, 2) = 2.*y*(1.+ y)*z;
715 output.access(26, 3) = 2.*x*(-1.+ (z*z));
716 output.access(26, 4) = 2.*x*(1.+ 2.*y)*z;
717 output.access(26, 5) = 2.*x*y*(1.+ y);
718 output.access(26, 6) = 0.;
719 output.access(26, 7) = 2.*(-1.+ (x*x))*z;
720 output.access(26, 8) = (-1.+ (x*x))*(1.+ 2.*y);
721 output.access(26, 9) = 0.;
722 }
723 break;
724 }
725 case OPERATOR_D4 : {
726 // Non-zero entries have Dk (derivative cardinality) indices {3,4,5,7,8,12}, all other entries are 0.
727 // Intitialize array by zero and then fill only non-zero entries.
728 const ordinal_type jend = output.extent(1);
729 const ordinal_type iend = output.extent(0);
730
731 for (ordinal_type j=0;j<jend;++j)
732 for (ordinal_type i=0;i<iend;++i)
733 output.access(i, j) = 0.0;
734
735 const auto x = input(0);
736 const auto y = input(1);
737 const auto z = input(2);
738 output.access(0, 3) = ((-1.+ z)*z)/2.;
739 output.access(0, 4) = ((-1.+ 2.*y)*(-1.+ 2.*z))/4.;
740 output.access(0, 5) = ((-1.+ y)*y)/2.;
741 output.access(0, 7) = ((-1.+ 2.*x)*(-1.+ 2.*z))/4.;
742 output.access(0, 8) = ((-1.+ 2.*x)*(-1.+ 2.*y))/4.;
743 output.access(0, 12)= ((-1.+ x)*x)/2.;
744
745 output.access(1, 3) = ((-1.+ z)*z)/2.;
746 output.access(1, 4) = ((-1.+ 2.*y)*(-1.+ 2.*z))/4.;
747 output.access(1, 5) = ((-1.+ y)*y)/2.;
748 output.access(1, 7) = ((1. + 2.*x)*(-1.+ 2.*z))/4.;
749 output.access(1, 8) = ((1. + 2.*x)*(-1.+ 2.*y))/4.;
750 output.access(1, 12)= (x*(1. + x))/2.;
751
752 output.access(2, 3) = ((-1.+ z)*z)/2.;
753 output.access(2, 4) = ((1. + 2.*y)*(-1.+ 2.*z))/4.;
754 output.access(2, 5) = (y*(1. + y))/2.;
755 output.access(2, 7) = ((1. + 2.*x)*(-1.+ 2.*z))/4.;
756 output.access(2, 8) = ((1. + 2.*x)*(1. + 2.*y))/4.;
757 output.access(2, 12)= (x*(1. + x))/2.;
758
759 output.access(3, 3) = ((-1.+ z)*z)/2.;
760 output.access(3, 4) = ((1. + 2.*y)*(-1.+ 2.*z))/4.;
761 output.access(3, 5) = (y*(1. + y))/2.;
762 output.access(3, 7) = ((-1.+ 2.*x)*(-1.+ 2.*z))/4.;
763 output.access(3, 8) = ((-1.+ 2.*x)*(1. + 2.*y))/4.;
764 output.access(3, 12)= ((-1.+ x)*x)/2.;
765
766 output.access(4, 3) = (z*(1. + z))/2.;
767 output.access(4, 4) = ((-1.+ 2.*y)*(1. + 2.*z))/4.;
768 output.access(4, 5) = ((-1.+ y)*y)/2.;
769 output.access(4, 7) = ((-1.+ 2.*x)*(1. + 2.*z))/4.;
770 output.access(4, 8) = ((-1.+ 2.*x)*(-1.+ 2.*y))/4.;
771 output.access(4, 12)= ((-1.+ x)*x)/2.;
772
773 output.access(5, 3) = (z*(1. + z))/2.;
774 output.access(5, 4) = ((-1.+ 2.*y)*(1. + 2.*z))/4.;
775 output.access(5, 5) = ((-1.+ y)*y)/2.;
776 output.access(5, 7) = ((1. + 2.*x)*(1. + 2.*z))/4.;
777 output.access(5, 8) = ((1. + 2.*x)*(-1.+ 2.*y))/4.;
778 output.access(5, 12)= (x*(1. + x))/2.;
779
780 output.access(6, 3) = (z*(1. + z))/2.;
781 output.access(6, 4) = ((1. + 2.*y)*(1. + 2.*z))/4.;
782 output.access(6, 5) = (y*(1. + y))/2.;
783 output.access(6, 7) = ((1. + 2.*x)*(1. + 2.*z))/4.;
784 output.access(6, 8) = ((1. + 2.*x)*(1. + 2.*y))/4.;
785 output.access(6, 12)= (x*(1. + x))/2.;
786
787 output.access(7, 3) = (z*(1. + z))/2.;
788 output.access(7, 4) = ((1. + 2.*y)*(1. + 2.*z))/4.;
789 output.access(7, 5) = (y*(1. + y))/2.;
790 output.access(7, 7) = ((-1.+ 2.*x)*(1. + 2.*z))/4.;
791 output.access(7, 8) = ((-1.+ 2.*x)*(1. + 2.*y))/4.;
792 output.access(7, 12)= ((-1.+ x)*x)/2.;
793
794 output.access(8, 3) = -((-1.+ z)*z);
795 output.access(8, 4) = -0.5 + y + z - 2.*y*z;
796 output.access(8, 5) = -((-1.+ y)*y);
797 output.access(8, 7) = x - 2.*x*z;
798 output.access(8, 8) = x - 2.*x*y;
799 output.access(8, 12)= 1. - x*x;
800
801 output.access(9, 3) = -((-1.+ z)*z);
802 output.access(9, 4) = y - 2.*y*z;
803 output.access(9, 5) = 1 - y*y;
804 output.access(9, 7) = 0.5 + x - z - 2.*x*z;
805 output.access(9, 8) = -((1. + 2.*x)*y);
806 output.access(9, 12)= -(x*(1. + x));
807
808 output.access(10, 3) = -((-1.+ z)*z);
809 output.access(10, 4) = 0.5 + y - z - 2.*y*z;
810 output.access(10, 5) = -(y*(1. + y));
811 output.access(10, 7) = x - 2.*x*z;
812 output.access(10, 8) = -(x*(1. + 2.*y));
813 output.access(10, 12)= 1. - x*x;
814
815 output.access(11, 3) = -((-1.+ z)*z);
816 output.access(11, 4) = y - 2.*y*z;
817 output.access(11, 5) = 1. - y*y;
818 output.access(11, 7) = -0.5 + x + z - 2.*x*z;
819 output.access(11, 8) = y - 2.*x*y;
820 output.access(11, 12)= -((-1.+ x)*x);
821
822 output.access(12, 3) = 1. - z*z;
823 output.access(12, 4) = z - 2.*y*z;
824 output.access(12, 5) = -((-1.+ y)*y);
825 output.access(12, 7) = z - 2.*x*z;
826 output.access(12, 8) = -0.5 + x + y - 2.*x*y;
827 output.access(12, 12)= -((-1.+ x)*x);
828
829 output.access(13, 3) = 1. - z*z;
830 output.access(13, 4) = z - 2.*y*z;
831 output.access(13, 5) = -((-1.+ y)*y);
832 output.access(13, 7) = -((1. + 2.*x)*z);
833 output.access(13, 8) = 0.5 + x - y - 2.*x*y;
834 output.access(13, 12)= -(x*(1. + x));
835
836 output.access(14, 3) = 1. - z*z;
837 output.access(14, 4) = -((1. + 2.*y)*z);
838 output.access(14, 5) = -(y*(1. + y));
839 output.access(14, 7) = -((1. + 2.*x)*z);
840 output.access(14, 8) = -((1. + 2.*x)*(1. + 2.*y))/2.;
841 output.access(14, 12)= -(x*(1. + x));
842
843 output.access(15, 3) = 1. - z*z;
844 output.access(15, 4) = -((1. + 2.*y)*z);
845 output.access(15, 5) = -(y*(1. + y));
846 output.access(15, 7) = z - 2.*x*z;
847 output.access(15, 8) = 0.5 + y - x*(1. + 2.*y);
848 output.access(15, 12)= -((-1.+ x)*x);
849
850 output.access(16, 3) = -(z*(1. + z));
851 output.access(16, 4) = 0.5 + z - y*(1. + 2.*z);
852 output.access(16, 5) = -((-1.+ y)*y);
853 output.access(16, 7) = -(x*(1. + 2.*z));
854 output.access(16, 8) = x - 2.*x*y;
855 output.access(16, 12)= 1. - x*x;
856
857 output.access(17, 3) = -(z*(1. + z));
858 output.access(17, 4) = -(y*(1. + 2.*z));
859 output.access(17, 5) = 1. - y*y;
860 output.access(17, 7) = -((1. + 2.*x)*(1. + 2.*z))/2.;
861 output.access(17, 8) = -((1. + 2.*x)*y);
862 output.access(17, 12)= -(x*(1. + x));
863
864 output.access(18, 3) = -(z*(1. + z));
865 output.access(18, 4) = -((1. + 2.*y)*(1. + 2.*z))/2.;
866 output.access(18, 5) = -(y*(1. + y));
867 output.access(18, 7) = -(x*(1. + 2.*z));
868 output.access(18, 8) = -(x*(1. + 2.*y));
869 output.access(18, 12)= 1. - x*x;
870
871 output.access(19, 3) = -(z*(1. + z));
872 output.access(19, 4) = -(y*(1. + 2.*z));
873 output.access(19, 5) = 1. - y*y;
874 output.access(19, 7) = 0.5 + z - x*(1. + 2.*z);
875 output.access(19, 8) = y - 2.*x*y;
876 output.access(19, 12)= -((-1.+ x)*x);
877
878 if constexpr (!serendipity) {
879 output.access(20, 3) = 4. - 4.*z*z;
880 output.access(20, 4) = -8.*y*z;
881 output.access(20, 5) = 4. - 4.*y*y;
882 output.access(20, 7) = -8.*x*z;
883 output.access(20, 8) = -8.*x*y;
884 output.access(20, 12)= 4. - 4.*x*x;
885
886 output.access(21, 3) = 2.*(-1.+ z)*z;
887 output.access(21, 4) = 2.*y*(-1.+ 2.*z);
888 output.access(21, 5) = 2.*(-1.+ y*y);
889 output.access(21, 7) = 2.*x*(-1.+ 2.*z);
890 output.access(21, 8) = 4.*x*y;
891 output.access(21, 12)= 2.*(-1.+ x*x);
892
893 output.access(22, 3) = 2.*z*(1. + z);
894 output.access(22, 4) = 2.*(y + 2.*y*z);
895 output.access(22, 5) = 2.*(-1.+ y*y);
896 output.access(22, 7) = 2.*(x + 2.*x*z);
897 output.access(22, 8) = 4.*x*y;
898 output.access(22, 12)= 2.*(-1.+ x*x);
899
900 output.access(23, 3) = 2.*(-1.+ z*z);
901 output.access(23, 4) = 4.*y*z;
902 output.access(23, 5) = 2.*(-1.+ y*y);
903 output.access(23, 7) = 2.*(-1.+ 2.*x)*z;
904 output.access(23, 8) = 2.*(-1.+ 2.*x)*y;
905 output.access(23, 12)= 2.*(-1.+ x)*x;
906
907 output.access(24, 3) = 2.*(-1.+ z*z);
908 output.access(24, 4) = 4.*y*z;
909 output.access(24, 5) = 2.*(-1.+ y*y);
910 output.access(24, 7) = 2.*(z + 2.*x*z);
911 output.access(24, 8) = 2.*(y + 2.*x*y);
912 output.access(24, 12)= 2.*x*(1. + x);
913
914 output.access(25, 3) = 2.*(-1.+ z*z);
915 output.access(25, 4) = 2.*(-1.+ 2.*y)*z;
916 output.access(25, 5) = 2.*(-1.+ y)*y;
917 output.access(25, 7) = 4.*x*z;
918 output.access(25, 8) = 2.*x*(-1.+ 2.*y);
919 output.access(25, 12)= 2.*(-1.+ x*x);
920
921 output.access(26, 3) = 2.*(-1.+ z*z);
922 output.access(26, 4) = 2.*(z + 2.*y*z);
923 output.access(26, 5) = 2.*y*(1. + y);
924 output.access(26, 7) = 4.*x*z;
925 output.access(26, 8) = 2.*(x + 2.*x*y);
926 output.access(26, 12)= 2.*(-1.+ x*x);
927 }
928 break;
929 }
930 case OPERATOR_MAX : {
931 const ordinal_type jend = output.extent(1);
932 const ordinal_type iend = output.extent(0);
933
934 for (ordinal_type j=0;j<jend;++j)
935 for (ordinal_type i=0;i<iend;++i)
936 output.access(i, j) = 0.0;
937 break;
938 }
939 default: {
940 INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
941 opType != OPERATOR_GRAD &&
942 opType != OPERATOR_CURL &&
943 opType != OPERATOR_D1 &&
944 opType != OPERATOR_D2 &&
945 opType != OPERATOR_D3 &&
946 opType != OPERATOR_D4 &&
947 opType != OPERATOR_MAX,
948 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_DEG2_FEM::Serial::getValues) operator is not supported");
949
950 }
951 }
952 }
953
954 template<bool serendipity>
955 template<typename DT,
956 typename outputValueValueType, class ...outputValueProperties,
957 typename inputPointValueType, class ...inputPointProperties>
958 void
960 getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
961 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
962 const EOperator operatorType ) {
963 typedef Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValueViewType;
964 typedef Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPointViewType;
965 typedef typename ExecSpace<typename inputPointViewType::execution_space,typename DT::execution_space>::ExecSpaceType ExecSpaceType;
966
967 // Number of evaluation points = dim 0 of inputPoints
968 const auto loopSize = inputPoints.extent(0);
969 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
970
971 switch (operatorType) {
972
973 case OPERATOR_VALUE: {
974 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_VALUE> FunctorType;
975 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
976 break;
977 }
978 case OPERATOR_GRAD:
979 case OPERATOR_D1: {
980 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_GRAD> FunctorType;
981 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
982 break;
983 }
984 case OPERATOR_CURL:
985 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument,
986 ">>> ERROR (Basis_HGRAD_HEX_DEG2_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D");
987 break;
988
989 case OPERATOR_DIV:
990 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
991 ">>> ERROR (Basis_HGRAD_HEX_DEG2_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D");
992 break;
993
994 case OPERATOR_D2: {
995 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D2> FunctorType;
996 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
997 break;
998 }
999 case OPERATOR_D3: {
1000 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D3> FunctorType;
1001 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
1002 break;
1003 }
1004 case OPERATOR_D4: {
1005 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D4> FunctorType;
1006 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
1007 break;
1008 }
1009 case OPERATOR_D5:
1010 case OPERATOR_D6: {
1011 INTREPID2_TEST_FOR_EXCEPTION( true, std::invalid_argument,
1012 ">>> ERROR (Basis_HGRAD_HEX_DEG2_FEM): operator not supported");
1013 // break; commented out becuase this always throws
1014 }
1015 case OPERATOR_D7:
1016 case OPERATOR_D8:
1017 case OPERATOR_D9:
1018 case OPERATOR_D10: {
1019 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_MAX> FunctorType;
1020 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
1021 break;
1022 }
1023 default: {
1024 INTREPID2_TEST_FOR_EXCEPTION( !( Intrepid2::isValidOperator(operatorType) ), std::invalid_argument,
1025 ">>> ERROR (Basis_HGRAD_HEX_DEG2_FEM): Invalid operator type");
1026 }
1027 }
1028 }
1029
1030 }
1031
1032 // -------------------------------------------------------------------------------------
1033
1034 template<bool serendipity, typename DT, typename OT, typename PT>
1037 this -> basisCardinality_ = serendipity ? 20 : 27;
1038 this -> basisDegree_ = 2;
1039 this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Hexahedron<8> >() );
1040 this -> basisType_ = BASIS_FEM_DEFAULT;
1041 this -> basisCoordinates_ = COORDINATES_CARTESIAN;
1042 this->functionSpace_ = FUNCTION_SPACE_HGRAD;
1043
1044 // initialize tags
1045 {
1046 // Basis-dependent intializations
1047 const ordinal_type tagSize = 4; // size of DoF tag, i.e., number of fields in the tag
1048 const ordinal_type posScDim = 0; // position in the tag, counting from 0, of the subcell dim
1049 const ordinal_type posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
1050 const ordinal_type posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
1051
1052 // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
1053 ordinal_type tags[108] = { 0, 0, 0, 1, // Nodes 0 to 7 follow vertex order of the topology
1054 0, 1, 0, 1,
1055 0, 2, 0, 1,
1056 0, 3, 0, 1,
1057 0, 4, 0, 1,
1058 0, 5, 0, 1,
1059 0, 6, 0, 1,
1060 0, 7, 0, 1,
1061 1, 0, 0, 1, // Node 8 -> edge 0
1062 1, 1, 0, 1, // Node 9 -> edge 1
1063 1, 2, 0, 1, // Node 10 -> edge 2
1064 1, 3, 0, 1, // Node 11 -> edge 3
1065 1, 8, 0, 1, // Node 12 -> edge 8
1066 1, 9, 0, 1, // Node 13 -> edge 9
1067 1,10, 0, 1, // Node 14 -> edge 10
1068 1,11, 0, 1, // Node 15 -> edge 11
1069 1, 4, 0, 1, // Node 16 -> edge 4
1070 1, 5, 0, 1, // Node 17 -> edge 5
1071 1, 6, 0, 1, // Node 18 -> edge 6
1072 1, 7, 0, 1, // Node 19 -> edge 7
1073
1074 // following entries not used for serendipity elements
1075 3, 0, 0, 1, // Node 20 -> Hexahedron
1076 2, 4, 0, 1, // Node 21 -> face 4
1077 2, 5, 0, 1, // Node 22 -> face 5
1078 2, 3, 0, 1, // Node 23 -> face 3
1079 2, 1, 0, 1, // Node 24 -> face 1
1080 2, 0, 0, 1, // Node 25 -> face 0
1081 2, 2, 0, 1, // Node 26 -> face 2
1082 };
1083
1084 // host tags
1085 OrdinalTypeArray1DHost tagView(&tags[0], serendipity ? 80 : 108);
1086
1087 // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
1088 this->setOrdinalTagData(this -> tagToOrdinal_,
1089 this -> ordinalToTag_,
1090 tagView,
1091 this -> basisCardinality_,
1092 tagSize,
1093 posScDim,
1094 posScOrd,
1095 posDfOrd);
1096 }
1097 // dofCoords on host and create its mirror view to device
1098 Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
1099 dofCoords("dofCoordsHost", this->basisCardinality_,this->basisCellTopology_.getDimension());
1100
1101 dofCoords(0,0) = -1.0; dofCoords(0,1) = -1.0; dofCoords(0,2) = -1.0;
1102 dofCoords(1,0) = 1.0; dofCoords(1,1) = -1.0; dofCoords(1,2) = -1.0;
1103 dofCoords(2,0) = 1.0; dofCoords(2,1) = 1.0; dofCoords(2,2) = -1.0;
1104 dofCoords(3,0) = -1.0; dofCoords(3,1) = 1.0; dofCoords(3,2) = -1.0;
1105 dofCoords(4,0) = -1.0; dofCoords(4,1) = -1.0; dofCoords(4,2) = 1.0;
1106 dofCoords(5,0) = 1.0; dofCoords(5,1) = -1.0; dofCoords(5,2) = 1.0;
1107 dofCoords(6,0) = 1.0; dofCoords(6,1) = 1.0; dofCoords(6,2) = 1.0;
1108 dofCoords(7,0) = -1.0; dofCoords(7,1) = 1.0; dofCoords(7,2) = 1.0;
1109
1110 dofCoords(8,0) = 0.0; dofCoords(8,1) = -1.0; dofCoords(8,2) = -1.0;
1111 dofCoords(9,0) = 1.0; dofCoords(9,1) = 0.0; dofCoords(9,2) = -1.0;
1112 dofCoords(10,0) = 0.0; dofCoords(10,1) = 1.0; dofCoords(10,2) = -1.0;
1113 dofCoords(11,0) = -1.0; dofCoords(11,1) = 0.0; dofCoords(11,2) = -1.0;
1114 dofCoords(12,0) = -1.0; dofCoords(12,1) = -1.0; dofCoords(12,2) = 0.0;
1115 dofCoords(13,0) = 1.0; dofCoords(13,1) = -1.0; dofCoords(13,2) = 0.0;
1116 dofCoords(14,0) = 1.0; dofCoords(14,1) = 1.0; dofCoords(14,2) = 0.0;
1117 dofCoords(15,0) = -1.0; dofCoords(15,1) = 1.0; dofCoords(15,2) = 0.0;
1118 dofCoords(16,0) = 0.0; dofCoords(16,1) = -1.0; dofCoords(16,2) = 1.0;
1119 dofCoords(17,0) = 1.0; dofCoords(17,1) = 0.0; dofCoords(17,2) = 1.0;
1120 dofCoords(18,0) = 0.0; dofCoords(18,1) = 1.0; dofCoords(18,2) = 1.0;
1121 dofCoords(19,0) = -1.0; dofCoords(19,1) = 0.0; dofCoords(19,2) = 1.0;
1122
1123 if constexpr (!serendipity) {
1124 dofCoords(20,0) = 0.0; dofCoords(20,1) = 0.0; dofCoords(20,2) = 0.0;
1125
1126 dofCoords(21,0) = 0.0; dofCoords(21,1) = 0.0; dofCoords(21,2) = -1.0;
1127 dofCoords(22,0) = 0.0; dofCoords(22,1) = 0.0; dofCoords(22,2) = 1.0;
1128 dofCoords(23,0) = -1.0; dofCoords(23,1) = 0.0; dofCoords(23,2) = 0.0;
1129 dofCoords(24,0) = 1.0; dofCoords(24,1) = 0.0; dofCoords(24,2) = 0.0;
1130 dofCoords(25,0) = 0.0; dofCoords(25,1) = -1.0; dofCoords(25,2) = 0.0;
1131 dofCoords(26,0) = 0.0; dofCoords(26,1) = 1.0; dofCoords(26,2) = 0.0;
1132 }
1133
1134 this->dofCoords_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoords);
1135 Kokkos::deep_copy(this->dofCoords_, dofCoords);
1136
1137 }
1138
1139}// namespace Intrepid2
1140#endif
void setOrdinalTagData(OrdinalTypeView3D &tagToOrdinal, OrdinalTypeView2D &ordinalToTag, const OrdinalTypeView1D tags, const ordinal_type basisCard, const ordinal_type tagSize, const ordinal_type posScDim, const ordinal_type posScOrd, const ordinal_type posDfOrd)
Kokkos::DynRankView< scalarType, DeviceType > dofCoords_
See Intrepid2::Basis_HGRAD_HEX_DEG2_FEM.