Intrepid2
Intrepid2_HGRAD_WEDGE_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_WEDGE_C2_FEM_DEF_HPP__
50#define __INTREPID2_HGRAD_WEDGE_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) = ((-1. + x + y)*(-1. + 2.*x + 2.*y)*(-1. + z)*z)/2.;
74 output.access(1) = (x*(-1. + 2.*x)*(-1. + z)*z)/2.;
75 output.access(2) = (y*(-1. + 2.*y)*(-1. + z)*z)/2.;
76 output.access(3) = ((-1. + x + y)*(-1. + 2.*x + 2.*y)*z*(1. + z))/2.;
77 output.access(4) = (x*(-1. + 2.*x)*z*(1. + z))/2.;
78 output.access(5) = (y*(-1. + 2.*y)*z*(1. + z))/2.;
79
80 output.access(6) = -2.*x*(-1. + x + y)*(-1. + z)*z;
81 output.access(7) = 2.*x*y*(-1. + z)*z;
82 output.access(8) = -2.*y*(-1. + x + y)*(-1. + z)*z;
83 output.access(9) = -((-1. + x + y)*(-1. + 2.*x + 2.*y)*(-1. + z)*(1. + z));
84 output.access(10) = -(x*(-1. + 2.*x)*(-1. + z)*(1. + z));
85 output.access(11) = -(y*(-1. + 2.*y)*(-1. + z)*(1. + z));
86 output.access(12) = -2.*x*(-1. + x + y)*z*(1. + z);
87 output.access(13) = 2.*x*y*z*(1. + z);
88 output.access(14) = -2.*y*(-1. + x + y)*z*(1. + z);
89 if constexpr (!serendipity) {
90 output.access(15) = 4.*x*(-1. + x + y)*(-1. + z)*(1. + z);
91 output.access(16) = -4.*x*y*(-1. + z)*(1. + z);
92 output.access(17) = 4.*y*(-1. + x + y)*(-1. + z)*(1. + z);
93 }
94 break;
95 }
96 case OPERATOR_GRAD: {
97 const auto x = input(0);
98 const auto y = input(1);
99 const auto z = input(2);
100
101 // output is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
102 output.access(0, 0) = ((-3 + 4*x + 4*y)*(-1 + z)*z)/2.;
103 output.access(0, 1) = ((-3 + 4*x + 4*y)*(-1 + z)*z)/2.;
104 output.access(0, 2) = ((-1 + x + y)*(-1 + 2*x + 2*y)*(-1 + 2*z))/2.;
105
106 output.access(1, 0) = ((-1 + 4*x)*(-1 + z)*z)/2.;
107 output.access(1, 1) = 0.;
108 output.access(1, 2) = (x*(-1 + 2*x)*(-1 + 2*z))/2.;
109
110 output.access(2, 0) = 0.;
111 output.access(2, 1) = ((-1 + 4*y)*(-1 + z)*z)/2.;
112 output.access(2, 2) = (y*(-1 + 2*y)*(-1 + 2*z))/2.;
113
114 output.access(3, 0) = ((-3 + 4*x + 4*y)*z*(1 + z))/2.;
115 output.access(3, 1) = ((-3 + 4*x + 4*y)*z*(1 + z))/2.;
116 output.access(3, 2) = ((-1 + x + y)*(-1 + 2*x + 2*y)*(1 + 2*z))/2.;
117
118 output.access(4, 0) = ((-1 + 4*x)*z*(1 + z))/2.;
119 output.access(4, 1) = 0.;
120 output.access(4, 2) = (x*(-1 + 2*x)*(1 + 2*z))/2.;
121
122 output.access(5, 0) = 0.;
123 output.access(5, 1) = ((-1 + 4*y)*z*(1 + z))/2.;
124 output.access(5, 2) = (y*(-1 + 2*y)*(1 + 2*z))/2.;
125
126 output.access(6, 0) = -2*(-1 + 2*x + y)*(-1 + z)*z;
127 output.access(6, 1) = -2*x*(-1 + z)*z;
128 output.access(6, 2) = 2*x*(-1 + x + y)*(1 - 2*z);
129
130 output.access(7, 0) = 2*y*(-1 + z)*z;
131 output.access(7, 1) = 2*x*(-1 + z)*z;
132 output.access(7, 2) = 2*x*y*(-1 + 2*z);
133
134 output.access(8, 0) = -2*y*(-1 + z)*z;
135 output.access(8, 1) = -2*(-1 + x + 2*y)*(-1 + z)*z;
136 output.access(8, 2) = 2*y*(-1 + x + y)*(1 - 2*z);
137
138 output.access(9, 0) = -(-3 + 4*x + 4*y)*(-1 + z*z);
139 output.access(9, 1) = -(-3 + 4*x + 4*y)*(-1 + z*z);
140 output.access(9, 2) = -2*(1 + 2*x*x - 3*y + 2*y*y + x*(-3 + 4*y))*z;
141
142 output.access(10, 0) = -(-1 + 4*x)*(-1 + z*z);
143 output.access(10, 1) = 0;
144 output.access(10, 2) = 2*(1 - 2*x)*x*z;
145
146 output.access(11, 0) = 0;
147 output.access(11, 1) = -(-1 + 4*y)*(-1 + z*z);
148 output.access(11, 2) = 2*(1 - 2*y)*y*z;
149
150 output.access(12, 0) = -2*(-1 + 2*x + y)*z*(1 + z);
151 output.access(12, 1) = -2*x*z*(1 + z);
152 output.access(12, 2) = -2*x*(-1 + x + y)*(1 + 2*z);
153
154 output.access(13, 0) = 2*y*z*(1 + z);
155 output.access(13, 1) = 2*x*z*(1 + z);
156 output.access(13, 2) = 2*x*y*(1 + 2*z);
157
158 output.access(14, 0) = -2*y*z*(1 + z);
159 output.access(14, 1) = -2*(-1 + x + 2*y)*z*(1 + z);
160 output.access(14, 2) = -2*y*(-1 + x + y)*(1 + 2*z);
161
162 if constexpr (!serendipity) {
163 output.access(15, 0) = 4*(-1 + 2*x + y)*(-1 + z*z);
164 output.access(15, 1) = 4*x*(-1 + z)*(1 + z);
165 output.access(15, 2) = 8*x*(-1 + x + y)*z;
166
167 output.access(16, 0) = -4*y*(-1 + z)*(1 + z);
168 output.access(16, 1) = -4*x*(-1 + z)*(1 + z);
169 output.access(16, 2) = -8*x*y*z;
170
171 output.access(17, 0) = 4*y*(-1 + z)*(1 + z);
172 output.access(17, 1) = 4*(-1 + x + 2*y)*(-1 + z*z);
173 output.access(17, 2) = 8*y*(-1 + x + y)*z;
174 }
175 break;
176 }
177 case OPERATOR_D2: {
178 const auto x = input(0);
179 const auto y = input(1);
180 const auto z = input(2);
181
182 output.access(0, 0) = 2.*(-1. + z)*z;
183 output.access(0, 1) = 2.*(-1. + z)*z;
184 output.access(0, 2) = ((-3. + 4.*x + 4.*y)*(-1. + 2.*z))/2.;
185 output.access(0, 3) = 2.*(-1. + z)*z;
186 output.access(0, 4) = ((-3. + 4.*x + 4.*y)*(-1. + 2.*z))/2.;
187 output.access(0, 5) = (-1. + x + y)*(-1. + 2.*x + 2.*y);
188
189 output.access(1, 0) = 2.*(-1. + z)*z;
190 output.access(1, 1) = 0.;
191 output.access(1, 2) = ((-1. + 4.*x)*(-1. + 2.*z))/2.;
192 output.access(1, 3) = 0.;
193 output.access(1, 4) = 0.;
194 output.access(1, 5) = x*(-1. + 2.*x);
195
196 output.access(2, 0) = 0.;
197 output.access(2, 1) = 0.;
198 output.access(2, 2) = 0.;
199 output.access(2, 3) = 2.*(-1. + z)*z;
200 output.access(2, 4) = ((-1. + 4.*y)*(-1. + 2.*z))/2.;
201 output.access(2, 5) = y*(-1. + 2.*y);
202
203 output.access(3, 0) = 2.*z*(1. + z);
204 output.access(3, 1) = 2.*z*(1. + z);
205 output.access(3, 2) = ((-3. + 4.*x + 4.*y)*(1. + 2.*z))/2.;
206 output.access(3, 3) = 2.*z*(1. + z);
207 output.access(3, 4) = ((-3. + 4.*x + 4.*y)*(1. + 2.*z))/2.;
208 output.access(3, 5) = (-1. + x + y)*(-1. + 2.*x + 2.*y);
209
210 output.access(4, 0) = 2.*z*(1. + z);
211 output.access(4, 1) = 0.;
212 output.access(4, 2) = ((-1. + 4.*x)*(1. + 2.*z))/2.;
213 output.access(4, 3) = 0.;
214 output.access(4, 4) = 0.;
215 output.access(4, 5) = x*(-1. + 2.*x);
216
217 output.access(5, 0) = 0.;
218 output.access(5, 1) = 0.;
219 output.access(5, 2) = 0.;
220 output.access(5, 3) = 2.*z*(1. + z);
221 output.access(5, 4) = ((-1. + 4.*y)*(1. + 2.*z))/2.;
222 output.access(5, 5) = y*(-1. + 2.*y);
223
224 output.access(6, 0) = -4.*(-1. + z)*z;
225 output.access(6, 1) = -2.*(-1. + z)*z;
226 output.access(6, 2) = -2.*(-1. + 2.*x + y)*(-1. + 2.*z);
227 output.access(6, 3) = 0.;
228 output.access(6, 4) = x*(2. - 4.*z);
229 output.access(6, 5) = -4.*x*(-1. + x + y);
230
231 output.access(7, 0) = 0.;
232 output.access(7, 1) = 2.*(-1. + z)*z;
233 output.access(7, 2) = 2.*y*(-1. + 2.*z);
234 output.access(7, 3) = 0.;
235 output.access(7, 4) = 2.*x*(-1. + 2.*z);
236 output.access(7, 5) = 4.*x*y;
237
238 output.access(8, 0) = 0.;
239 output.access(8, 1) = -2.*(-1. + z)*z;
240 output.access(8, 2) = y*(2. - 4.*z);
241 output.access(8, 3) = -4.*(-1. + z)*z;
242 output.access(8, 4) = -2.*(-1. + x + 2.*y)*(-1. + 2.*z);
243 output.access(8, 5) = -4.*y*(-1. + x + y);
244
245 output.access(9, 0) = 4. - 4.*z*z;
246 output.access(9, 1) = 4. - 4.*z*z;
247 output.access(9, 2) = -2.*(-3. + 4.*x + 4.*y)*z;
248 output.access(9, 3) = 4. - 4.*z*z;
249 output.access(9, 4) = -2.*(-3. + 4.*x + 4.*y)*z;
250 output.access(9, 5) = -2.*(-1. + x + y)*(-1. + 2.*x + 2.*y);
251
252 output.access(10, 0) = 4. - 4.*z*z;
253 output.access(10, 1) = 0.;
254 output.access(10, 2) = (2. - 8.*x)*z;
255 output.access(10, 3) = 0.;
256 output.access(10, 4) = 0.;
257 output.access(10, 5) = -2.*x*(-1. + 2.*x);
258
259 output.access(11, 0) = 0.;
260 output.access(11, 1) = 0.;
261 output.access(11, 2) = 0.;
262 output.access(11, 3) = 4. - 4.*z*z;
263 output.access(11, 4) = (2. - 8.*y)*z;
264 output.access(11, 5) = -2.*y*(-1. + 2.*y);
265
266 output.access(12, 0) = -4.*z*(1. + z);
267 output.access(12, 1) = -2.*z*(1. + z);
268 output.access(12, 2) = -2.*(-1. + 2.*x + y)*(1. + 2.*z);
269 output.access(12, 3) = 0.;
270 output.access(12, 4) = -2.*(x + 2.*x*z);
271 output.access(12, 5) = -4.*x*(-1. + x + y);
272
273 output.access(13, 0) = 0.;
274 output.access(13, 1) = 2.*z*(1. + z);
275 output.access(13, 2) = 2.*(y + 2.*y*z);
276 output.access(13, 3) = 0.;
277 output.access(13, 4) = 2.*(x + 2.*x*z);
278 output.access(13, 5) = 4.*x*y;
279
280 output.access(14, 0) = 0.;
281 output.access(14, 1) = -2.*z*(1. + z);
282 output.access(14, 2) = -2.*(y + 2.*y*z);
283 output.access(14, 3) = -4.*z*(1. + z);
284 output.access(14, 4) = -2.*(-1. + x + 2.*y)*(1. + 2.*z);
285 output.access(14, 5) = -4.*y*(-1. + x + y);
286
287 if constexpr (!serendipity) {
288 output.access(15, 0) = 8.*(-1. + z*z);
289 output.access(15, 1) = 4.*(-1. + z*z);
290 output.access(15, 2) = 8.*(-1. + 2.*x + y)*z;
291 output.access(15, 3) = 0.;
292 output.access(15, 4) = 8.*x*z;
293 output.access(15, 5) = 8.*x*(-1. + x + y);
294
295 output.access(16, 0) = 0.;
296 output.access(16, 1) = 4. - 4.*z*z;
297 output.access(16, 2) = -8.*y*z;
298 output.access(16, 3) = 0.;
299 output.access(16, 4) = -8.*x*z;
300 output.access(16, 5) = -8.*x*y;
301
302
303 output.access(17, 0) = 0.;
304 output.access(17, 1) = 4.*(-1. + z*z);
305 output.access(17, 2) = 8.*y*z;
306 output.access(17, 3) = 8.*(-1. + z*z);
307 output.access(17, 4) = 8.*(-1. + x + 2.*y)*z;
308 output.access(17, 5) = 8.*y*(-1. + x + y);
309 }
310 break;
311 }
312 case OPERATOR_D3: {
313 const auto x = input(0);
314 const auto y = input(1);
315 const auto z = input(2);
316
317 output.access(0, 0) = 0.;
318 output.access(0, 1) = 0.;
319 output.access(0, 2) = -2. + 4.*z;
320 output.access(0, 3) = 0.;
321 output.access(0, 4) = -2. + 4.*z;
322 output.access(0, 5) = -3. + 4.*x + 4.*y;
323 output.access(0, 6) = 0.;
324 output.access(0, 7) = -2. + 4.*z;
325 output.access(0, 8) = -3. + 4.*x + 4.*y;
326 output.access(0, 9) = 0.;
327
328 output.access(1, 0) = 0.;
329 output.access(1, 1) = 0.;
330 output.access(1, 2) = -2. + 4.*z;
331 output.access(1, 3) = 0.;
332 output.access(1, 4) = 0.;
333 output.access(1, 5) = -1 + 4.*x;
334 output.access(1, 6) = 0.;
335 output.access(1, 7) = 0.;
336 output.access(1, 8) = 0.;
337 output.access(1, 9) = 0.;
338
339 output.access(2, 0) = 0.;
340 output.access(2, 1) = 0.;
341 output.access(2, 2) = 0.;
342 output.access(2, 3) = 0.;
343 output.access(2, 4) = 0.;
344 output.access(2, 5) = 0.;
345 output.access(2, 6) = 0.;
346 output.access(2, 7) = -2. + 4.*z;
347 output.access(2, 8) = -1 + 4.*y;
348 output.access(2, 9) = 0.;
349
350 output.access(3, 0) = 0.;
351 output.access(3, 1) = 0.;
352 output.access(3, 2) = 2. + 4.*z;
353 output.access(3, 3) = 0.;
354 output.access(3, 4) = 2. + 4.*z;
355 output.access(3, 5) = -3. + 4.*x + 4.*y;
356 output.access(3, 6) = 0.;
357 output.access(3, 7) = 2. + 4.*z;
358 output.access(3, 8) = -3. + 4.*x + 4.*y;
359 output.access(3, 9) = 0.;
360
361 output.access(4, 0) = 0.;
362 output.access(4, 1) = 0.;
363 output.access(4, 2) = 2. + 4.*z;
364 output.access(4, 3) = 0.;
365 output.access(4, 4) = 0.;
366 output.access(4, 5) = -1 + 4.*x;
367 output.access(4, 6) = 0.;
368 output.access(4, 7) = 0.;
369 output.access(4, 8) = 0.;
370 output.access(4, 9) = 0.;
371
372 output.access(5, 0) = 0.;
373 output.access(5, 1) = 0.;
374 output.access(5, 2) = 0.;
375 output.access(5, 3) = 0.;
376 output.access(5, 4) = 0.;
377 output.access(5, 5) = 0.;
378 output.access(5, 6) = 0.;
379 output.access(5, 7) = 2. + 4.*z;
380 output.access(5, 8) = -1 + 4.*y;
381 output.access(5, 9) = 0.;
382
383 output.access(6, 0) = 0.;
384 output.access(6, 1) = 0.;
385 output.access(6, 2) = 4. - 8.*z;
386 output.access(6, 3) = 0.;
387 output.access(6, 4) = 2. - 4.*z;
388 output.access(6, 5) = -4.*(-1 + 2*x + y);
389 output.access(6, 6) = 0.;
390 output.access(6, 7) = 0.;
391 output.access(6, 8) = -4.*x;
392 output.access(6, 9) = 0.;
393
394 output.access(7, 0) = 0.;
395 output.access(7, 1) = 0.;
396 output.access(7, 2) = 0.;
397 output.access(7, 3) = 0.;
398 output.access(7, 4) = -2. + 4.*z;
399 output.access(7, 5) = 4.*y;
400 output.access(7, 6) = 0.;
401 output.access(7, 7) = 0.;
402 output.access(7, 8) = 4.*x;
403 output.access(7, 9) = 0.;
404
405 output.access(8, 0) = 0.;
406 output.access(8, 1) = 0.;
407 output.access(8, 2) = 0.;
408 output.access(8, 3) = 0.;
409 output.access(8, 4) = 2. - 4.*z;
410 output.access(8, 5) = -4.*y;
411 output.access(8, 6) = 0.;
412 output.access(8, 7) = 4. - 8.*z;
413 output.access(8, 8) = -4.*(-1 + x + 2*y);
414 output.access(8, 9) = 0.;
415
416 output.access(9, 0) = 0.;
417 output.access(9, 1) = 0.;
418 output.access(9, 2) = -8.*z;
419 output.access(9, 3) = 0.;
420 output.access(9, 4) = -8.*z;
421 output.access(9, 5) = 6. - 8.*x - 8.*y;
422 output.access(9, 6) = 0.;
423 output.access(9, 7) = -8.*z;
424 output.access(9, 8) = 6. - 8.*x - 8.*y;
425 output.access(9, 9) = 0.;
426
427 output.access(10, 0) = 0.;
428 output.access(10, 1) = 0.;
429 output.access(10, 2) = -8.*z;
430 output.access(10, 3) = 0.;
431 output.access(10, 4) = 0.;
432 output.access(10, 5) = 2. - 8.*x;
433 output.access(10, 6) = 0.;
434 output.access(10, 7) = 0.;
435 output.access(10, 8) = 0.;
436 output.access(10, 9) = 0.;
437
438 output.access(11, 0) = 0.;
439 output.access(11, 1) = 0.;
440 output.access(11, 2) = 0.;
441 output.access(11, 3) = 0.;
442 output.access(11, 4) = 0.;
443 output.access(11, 5) = 0.;
444 output.access(11, 6) = 0.;
445 output.access(11, 7) = -8.*z;
446 output.access(11, 8) = 2. - 8.*y;
447 output.access(11, 9) = 0.;
448
449 output.access(12, 0) = 0.;
450 output.access(12, 1) = 0.;
451 output.access(12, 2) = -4. - 8.*z;
452 output.access(12, 3) = 0.;
453 output.access(12, 4) = -2. - 4.*z;
454 output.access(12, 5) = -4.*(-1 + 2*x + y);
455 output.access(12, 6) = 0.;
456 output.access(12, 7) = 0.;
457 output.access(12, 8) = -4.*x;
458 output.access(12, 9) = 0.;
459
460 output.access(13, 0) = 0.;
461 output.access(13, 1) = 0.;
462 output.access(13, 2) = 0.;
463 output.access(13, 3) = 0.;
464 output.access(13, 4) = 2. + 4.*z;
465 output.access(13, 5) = 4.*y;
466 output.access(13, 6) = 0.;
467 output.access(13, 7) = 0.;
468 output.access(13, 8) = 4.*x;
469 output.access(13, 9) = 0.;
470
471 output.access(14, 0) = 0.;
472 output.access(14, 1) = 0.;
473 output.access(14, 2) = 0.;
474 output.access(14, 3) = 0.;
475 output.access(14, 4) = -2. - 4.*z;
476 output.access(14, 5) = -4.*y;
477 output.access(14, 6) = 0.;
478 output.access(14, 7) = -4. - 8.*z;
479 output.access(14, 8) = -4.*(-1 + x + 2*y);
480 output.access(14, 9) = 0.;
481
482 if constexpr (!serendipity) {
483 output.access(15, 0) = 0.;
484 output.access(15, 1) = 0.;
485 output.access(15, 2) = 16.*z;
486 output.access(15, 3) = 0.;
487 output.access(15, 4) = 8.*z;
488 output.access(15, 5) = 8.*(-1 + 2*x + y);
489 output.access(15, 6) = 0.;
490 output.access(15, 7) = 0.;
491 output.access(15, 8) = 8.*x;
492 output.access(15, 9) = 0.;
493
494 output.access(16, 0) = 0.;
495 output.access(16, 1) = 0.;
496 output.access(16, 2) = 0.;
497 output.access(16, 3) = 0.;
498 output.access(16, 4) = -8.*z;
499 output.access(16, 5) = -8.*y;
500 output.access(16, 6) = 0.;
501 output.access(16, 7) = 0.;
502 output.access(16, 8) = -8.*x;
503 output.access(16, 9) = 0.;
504
505 output.access(17, 0) = 0.;
506 output.access(17, 1) = 0.;
507 output.access(17, 2) = 0.;
508 output.access(17, 3) = 0.;
509 output.access(17, 4) = 8.*z;
510 output.access(17, 5) = 8.*y;
511 output.access(17, 6) = 0.;
512 output.access(17, 7) = 16.*z;
513 output.access(17, 8) = 8.*(-1 + x + 2*y);
514 output.access(17, 9) = 0.;
515 }
516 break;
517 }
518 case OPERATOR_D4: {
519 const ordinal_type jend = output.extent(1);
520 const ordinal_type iend = output.extent(0);
521
522 for (ordinal_type j=0;j<jend;++j)
523 for (ordinal_type i=0;i<iend;++i)
524 output.access(i, j) = 0.0;
525
526 output.access(0, 5) = 4.;
527 output.access(0, 8) = 4.;
528 output.access(0,12) = 4.;
529
530 output.access(1, 5) = 4.;
531
532 output.access(2,12) = 4.;
533
534 output.access(3, 5) = 4.;
535 output.access(3, 8) = 4.;
536 output.access(3,12) = 4.;
537
538 output.access(4, 5) = 4.0;
539
540 output.access(5,12) = 4.0;
541
542 output.access(6, 5) =-8.;
543 output.access(6, 8) =-4.;
544
545 output.access(7, 8) = 4.;
546
547 output.access(8, 8) =-4.;
548 output.access(8,12) =-8.;
549
550 output.access(9, 5) =-8.;
551 output.access(9, 8) =-8.;
552 output.access(9,12) =-8.;
553
554 output.access(10, 5) =-8.;
555
556 output.access(11,12) =-8.;
557
558 output.access(12, 5) =-8.;
559 output.access(12, 8) =-4.;
560
561 output.access(13, 8) = 4.;
562
563 output.access(14, 8) =-4;
564 output.access(14,12) =-8.;
565
566 if constexpr (!serendipity) {
567 output.access(15, 5) =16.;
568 output.access(15, 8) = 8.;
569
570 output.access(16, 8) =-8.;
571
572
573 output.access(17, 8) = 8.;
574 output.access(17,12) =16.;
575 }
576 break;
577 }
578 case OPERATOR_MAX : {
579 const ordinal_type jend = output.extent(1);
580 const ordinal_type iend = output.extent(0);
581
582 for (ordinal_type j=0;j<jend;++j)
583 for (ordinal_type i=0;i<iend;++i)
584 output.access(i, j) = 0.0;
585 break;
586 }
587
588 default: {
589 INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
590 opType != OPERATOR_GRAD &&
591 opType != OPERATOR_D2 &&
592 opType != OPERATOR_D3 &&
593 opType != OPERATOR_D4 &&
594 opType != OPERATOR_MAX,
595 ">>> ERROR: (Intrepid2::Basis_HGRAD_WEDGE_DEG2_FEM::Serial::getValues) operator is not supported");
596 }
597 }
598 }
599
600 template<bool serendipity>
601 template<typename DT,
602 typename outputValueValueType, class ...outputValueProperties,
603 typename inputPointValueType, class ...inputPointProperties>
604 void
606 getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
607 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
608 const EOperator operatorType ) {
609 typedef Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValueViewType;
610 typedef Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPointViewType;
611 typedef typename ExecSpace<typename inputPointViewType::execution_space,typename DT::execution_space>::ExecSpaceType ExecSpaceType;
612
613 // Number of evaluation points = dim 0 of inputPoints
614 const auto loopSize = inputPoints.extent(0);
615 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
616
617 switch (operatorType) {
618
619 case OPERATOR_VALUE: {
620 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_VALUE> FunctorType;
621 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
622 break;
623 }
624 case OPERATOR_GRAD:
625 case OPERATOR_D1: {
626 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_GRAD> FunctorType;
627 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
628 break;
629 }
630 case OPERATOR_CURL: {
631 INTREPID2_TEST_FOR_EXCEPTION( operatorType == OPERATOR_CURL, std::invalid_argument,
632 ">>> ERROR (Basis_HGRAD_WEDGE_DEG2_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D");
633 break;
634 }
635 case OPERATOR_DIV: {
636 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
637 ">>> ERROR (Basis_HGRAD_WEDGE_DEG2_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D");
638 break;
639 }
640 case OPERATOR_D2: {
641 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D2> FunctorType;
642 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
643 break;
644 }
645 case OPERATOR_D3: {
646 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D3> FunctorType;
647 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
648 break;
649 }
650 case OPERATOR_D4: {
651 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D4> FunctorType;
652 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
653 break;
654 }
655 case OPERATOR_D5:
656 case OPERATOR_D6:
657 case OPERATOR_D7:
658 case OPERATOR_D8:
659 case OPERATOR_D9:
660 case OPERATOR_D10: {
661 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_MAX> FunctorType;
662 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
663 break;
664 }
665 default: {
666 INTREPID2_TEST_FOR_EXCEPTION( !( Intrepid2::isValidOperator(operatorType) ), std::invalid_argument,
667 ">>> ERROR (Basis_HGRAD_WEDGE_DEG2_FEM): Invalid operator type");
668 }
669 }
670 }
671
672 }
673
674 // -------------------------------------------------------------------------------------
675
676 template<bool serendipity, typename DT, typename OT, typename PT>
679 this->basisCardinality_ = serendipity ? 15 : 18;
680 this->basisDegree_ = 2;
681 this->basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Wedge<6> >() );
682 this->basisType_ = BASIS_FEM_DEFAULT;
683 this->basisCoordinates_ = COORDINATES_CARTESIAN;
684 this->functionSpace_ = FUNCTION_SPACE_HGRAD;
685
686 // initialize tags
687 {
688 // Basis-dependent intializations
689 const ordinal_type tagSize = 4; // size of DoF tag
690 const ordinal_type posScDim = 0; // position in the tag, counting from 0, of the subcell dim
691 const ordinal_type posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
692 const ordinal_type posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
693
694 // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
695 ordinal_type tags[72] = { 0, 0, 0, 1,
696 0, 1, 0, 1,
697 0, 2, 0, 1,
698 0, 3, 0, 1,
699 0, 4, 0, 1,
700 0, 5, 0, 1,
701 1, 0, 0, 1,
702 1, 1, 0, 1,
703 1, 2, 0, 1,
704 1, 6, 0, 1,
705 1, 7, 0, 1,
706 1, 8, 0, 1,
707 1, 3, 0, 1,
708 1, 4, 0, 1,
709 1, 5, 0, 1,
710 // following entries not used for serendipity elements
711 2, 0, 0, 1,
712 2, 1, 0, 1,
713 2, 2, 0, 1
714 };
715
716 // host tags
717 OrdinalTypeArray1DHost tagView(&tags[0], serendipity ? 60 : 72);
718
719 // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
721 this->ordinalToTag_,
722 tagView,
723 this->basisCardinality_,
724 tagSize,
725 posScDim,
726 posScOrd,
727 posDfOrd);
728 }
729
730 // dofCoords on host and create its mirror view to device
731 Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
732 dofCoords("dofCoordsHost", this->basisCardinality_,this->basisCellTopology_.getDimension());
733
734 dofCoords(0,0) = 0.0; dofCoords(0,1) = 0.0; dofCoords(0,2) = -1.0;
735 dofCoords(1,0) = 1.0; dofCoords(1,1) = 0.0; dofCoords(1,2) = -1.0;
736 dofCoords(2,0) = 0.0; dofCoords(2,1) = 1.0; dofCoords(2,2) = -1.0;
737 dofCoords(3,0) = 0.0; dofCoords(3,1) = 0.0; dofCoords(3,2) = 1.0;
738 dofCoords(4,0) = 1.0; dofCoords(4,1) = 0.0; dofCoords(4,2) = 1.0;
739 dofCoords(5,0) = 0.0; dofCoords(5,1) = 1.0; dofCoords(5,2) = 1.0;
740
741 dofCoords(6,0) = 0.5; dofCoords(6,1) = 0.0; dofCoords(6,2) = -1.0;
742 dofCoords(7,0) = 0.5; dofCoords(7,1) = 0.5; dofCoords(7,2) = -1.0;
743 dofCoords(8,0) = 0.0; dofCoords(8,1) = 0.5; dofCoords(8,2) = -1.0;
744 dofCoords(9,0) = 0.0; dofCoords(9,1) = 0.0; dofCoords(9,2) = 0.0;
745 dofCoords(10,0)= 1.0; dofCoords(10,1)= 0.0; dofCoords(10,2)= 0.0;
746 dofCoords(11,0)= 0.0; dofCoords(11,1)= 1.0; dofCoords(11,2)= 0.0;
747
748 dofCoords(12,0)= 0.5; dofCoords(12,1)= 0.0; dofCoords(12,2)= 1.0;
749 dofCoords(13,0)= 0.5; dofCoords(13,1)= 0.5; dofCoords(13,2)= 1.0;
750 dofCoords(14,0)= 0.0; dofCoords(14,1)= 0.5; dofCoords(14,2)= 1.0;
751
752 if constexpr (!serendipity) {
753 dofCoords(15,0)= 0.5; dofCoords(15,1)= 0.0; dofCoords(15,2)= 0.0;
754 dofCoords(16,0)= 0.5; dofCoords(16,1)= 0.5; dofCoords(16,2)= 0.0;
755 dofCoords(17,0)= 0.0; dofCoords(17,1)= 0.5; dofCoords(17,2)= 0.0;
756 }
757
758 this->dofCoords_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoords);
759 Kokkos::deep_copy(this->dofCoords_, dofCoords);
760 }
761
762}// namespace Intrepid2
763#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_WEDGE_DEG2_FEM.