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