Intrepid2
Intrepid2_HGRAD_QUAD_C2_FEMDef.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Intrepid2 Package
4 //
5 // Copyright 2007 NTESS and the Intrepid2 contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
16 #ifndef __INTREPID2_HGRAD_QUAD_C2_FEM_DEF_HPP__
17 #define __INTREPID2_HGRAD_QUAD_C2_FEM_DEF_HPP__
18 
19 namespace Intrepid2 {
20 
21  // -------------------------------------------------------------------------------------
22 
23  namespace Impl {
24 
25  template<bool serendipity>
26  template<EOperator opType>
27  template<typename OutputViewType,
28  typename inputViewType>
29  KOKKOS_INLINE_FUNCTION
30  void
31  Basis_HGRAD_QUAD_DEG2_FEM<serendipity>::Serial<opType>::
32  getValues( OutputViewType output,
33  const inputViewType input ) {
34  switch (opType) {
35  case OPERATOR_VALUE : {
36  const auto x = input(0);
37  const auto y = input(1);
38 
39  // output is a rank-1 array with dimensions (basisCardinality_)
40  if constexpr (!serendipity) {
41  output.access(0) = x*(x - 1.0)*y*(y - 1.0)/4.0;
42  output.access(1) = x*(x + 1.0)*y*(y - 1.0)/4.0;
43  output.access(2) = x*(x + 1.0)*y*(y + 1.0)/4.0;
44  output.access(3) = x*(x - 1.0)*y*(y + 1.0)/4.0;
45  // edge midpoints basis functions
46  output.access(4) = (1.0 - x)*(1.0 + x)*y*(y - 1.0)/2.0;
47  output.access(5) = x*(x + 1.0)*(1.0 - y)*(1.0 + y)/2.0;
48  output.access(6) = (1.0 - x)*(1.0 + x)*y*(y + 1.0)/2.0;
49  output.access(7) = x*(x - 1.0)*(1.0 - y)*(1.0 + y)/2.0;
50 
51  // quad bubble basis function
52  output.access(8) = (1.0 - x)*(1.0 + x)*(1.0 - y)*(1.0 + y);
53 
54  } else { //serendipity
55 
56  output.access(0) = 0.25*(1.0 - x)*(1.0 - y)*(-x - y - 1.0);
57  output.access(1) = 0.25*(1.0 + x)*(1.0 - y)*( x - y - 1.0);
58  output.access(2) = 0.25*(1.0 + x)*(1.0 + y)*( x + y - 1.0);
59  output.access(3) = 0.25*(1.0 - x)*(1.0 + y)*(-x + y - 1.0);
60 
61  output.access(4) = 0.5*(1.0 - x*x)*(1.0 - y);
62  output.access(5) = 0.5*(1.0 + x)*(1.0 - y*y);
63  output.access(6) = 0.5*(1.0 - x*x)*(1.0 + y);
64  output.access(7) = 0.5*(1.0 - x)*(1.0 - y*y);
65  }
66 
67  break;
68  }
69  case OPERATOR_D1 :
70  case OPERATOR_GRAD : {
71  const auto x = input(0);
72  const auto y = input(1);
73 
74  if constexpr (!serendipity) {
75 
76  output.access(0, 0) = (-0.25 + 0.5*x)*(-1. + y)*y;
77  output.access(0, 1) = (-1.0 + x)*x*(-0.25 + 0.5*y);
78 
79  output.access(1, 0) = (0.25 + 0.5*x)*(-1. + y)*y;
80  output.access(1, 1) = x*(1. + x)*(-0.25 + 0.5*y);
81 
82  output.access(2, 0) = (0.25 + 0.5*x)*y*(1. + y);
83  output.access(2, 1) = x*(1. + x)*(0.25 + 0.5*y);
84 
85  output.access(3, 0) = (-0.25 + 0.5*x)*y*(1. + y);
86  output.access(3, 1) = (-1. + x)*x*(0.25 + 0.5*y);
87 
88  output.access(4, 0) = x*(1.0 - y)*y;
89  output.access(4, 1) = 0.5*(1.0 - x)*(1.0 + x)*(-1.0 + 2.0*y);
90 
91  output.access(5, 0) = 0.5*(1.0 - y)*(1.0 + y)*(1.0 + 2.0*x);
92  output.access(5, 1) =-x*(1.0 + x)*y;
93 
94  output.access(6, 0) =-y*(1.0 + y)*x;
95  output.access(6, 1) = 0.5*(1.0 - x)*(1.0 + x)*(1.0 + 2.0*y);
96 
97  output.access(7, 0) = 0.5*(1.0 - y)*(1.0+ y)*(-1.0 + 2.0*x);
98  output.access(7, 1) = (1.0 - x)*x*y;
99 
100  output.access(8, 0) =-2.0*(1.0 - y)*(1.0 + y)*x;
101  output.access(8, 1) =-2.0*(1.0 - x)*(1.0 + x)*y;
102 
103  } else { //serendipity
104 
105  output.access(0, 0) = -0.25*(1.0-y)*(-x-y-1.0) - 0.25*(1.0-x)*(1.0-y);
106  output.access(0, 1) = -0.25*(1.0-x)*(-x-y-1.0) - 0.25*(1.0-x)*(1.0-y);
107 
108  output.access(1, 0) = 0.25*(1.0-y)*( x-y-1.0) + 0.25*(1.0+x)*(1.0-y);
109  output.access(1, 1) = -0.25*(1.0+x)*( x-y-1.0) - 0.25*(1.0+x)*(1.0-y);
110 
111  output.access(2, 0) = 0.25*(1.0+y)*( x+y-1.0) + 0.25*(1.0+x)*(1.0+y);
112  output.access(2, 1) = 0.25*(1.0+x)*( x+y-1.0) + 0.25*(1.0+x)*(1.0+y);
113 
114  output.access(3, 0) = -0.25*(1.0+y)*(-x+y-1.0) - 0.25*(1.0-x)*(1.0+y);
115  output.access(3, 1) = 0.25*(1.0-x)*(-x+y-1.0) + 0.25*(1.0-x)*(1.0+y);
116 
117  output.access(4, 0) = -x*(1.0-y);
118  output.access(4, 1) = -0.5*(1.0-x*x);
119 
120  output.access(5, 0) = 0.5*(1.0-y*y);
121  output.access(5, 1) = -y*(1.0+x);
122 
123  output.access(6, 0) = -x*(1.0+y);
124  output.access(6, 1) = 0.5*(1.0-x*x);
125 
126  output.access(7, 0) = -0.5*(1.0-y*y);
127  output.access(7, 1) = -y*(1.0-x);
128  }
129  break;
130  }
131  case OPERATOR_CURL : {
132  const auto x = input(0);
133  const auto y = input(1);
134 
135  // output.access is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
136  // CURL(u) = (u_y, -u_x), is rotated GRAD
137 
138  if constexpr (!serendipity) {
139  output.access(0, 1) =-(-0.25 + 0.5*x)*(-1. + y)*y;
140  output.access(0, 0) = (-1.0 + x)*x*(-0.25 + 0.5*y);
141 
142  output.access(1, 1) =-(0.25 + 0.5*x)*(-1. + y)*y;
143  output.access(1, 0) = x*(1. + x)*(-0.25 + 0.5*y);
144 
145  output.access(2, 1) =-(0.25 + 0.5*x)*y*(1. + y);
146  output.access(2, 0) = x*(1. + x)*(0.25 + 0.5*y);
147 
148  output.access(3, 1) =-(-0.25 + 0.5*x)*y*(1. + y);
149  output.access(3, 0) = (-1. + x)*x*(0.25 + 0.5*y);
150 
151  output.access(4, 1) =-x*(1.0 - y)*y;
152  output.access(4, 0) = 0.5*(1.0 - x)*(1.0 + x)*(-1.0 + 2.0*y);
153 
154  output.access(5, 1) =-0.5*(1.0 - y)*(1.0 + y)*(1.0 + 2.0*x);
155  output.access(5, 0) =-x*(1.0 + x)*y;
156 
157  output.access(6, 1) = y*(1.0 + y)*x;
158  output.access(6, 0) = 0.5*(1.0 - x)*(1.0 + x)*(1.0 + 2.0*y);
159 
160  output.access(7, 1) =-0.5*(1.0 - y)*(1.0 + y)*(-1.0 + 2.0*x);
161  output.access(7, 0) = (1.0 - x)*x*y;
162 
163  output.access(8, 1) = 2.0*(1.0 - y)*(1.0 + y)*x;
164  output.access(8, 0) =-2.0*(1.0 - x)*(1.0 + x)*y;
165 
166  } else { //serendipity
167  output.access(0, 1) = 0.25*(1.0-y)*(-x-y-1.0) + 0.25*(1.0-x)*(1.0-y);
168  output.access(0, 0) = -0.25*(1.0-x)*(-x-y-1.0) - 0.25*(1.0-x)*(1.0-y);
169 
170  output.access(1, 1) = -0.25*(1.0-y)*( x-y-1.0) - 0.25*(1.0+x)*(1.0-y);
171  output.access(1, 0) = -0.25*(1.0+x)*( x-y-1.0) - 0.25*(1.0+x)*(1.0-y);
172 
173  output.access(2, 1) = -0.25*(1.0+y)*( x+y-1.0) - 0.25*(1.0+x)*(1.0+y);
174  output.access(2, 0) = 0.25*(1.0+x)*( x+y-1.0) + 0.25*(1.0+x)*(1.0+y);
175 
176  output.access(3, 1) = 0.25*(1.0+y)*(-x+y-1.0) + 0.25*(1.0-x)*(1.0+y);
177  output.access(3, 0) = 0.25*(1.0-x)*(-x+y-1.0) + 0.25*(1.0-x)*(1.0+y);
178 
179  output.access(4, 1) = x*(1.0-y);
180  output.access(4, 0) = -0.5*(1.0-x*x);
181 
182  output.access(5, 1) = -0.5*(1.0-y*y);
183  output.access(5, 0) = -y*(1.0+x);
184 
185  output.access(6, 1) = x*(1.0+y);
186  output.access(6, 0) = 0.5*(1.0-x*x);
187 
188  output.access(7, 1) = 0.5*(1.0-y*y);
189  output.access(7, 0) = -y*(1.0-x);
190  }
191  break;
192  }
193  case OPERATOR_D2 : {
194  const auto x = input(0);
195  const auto y = input(1);
196  // output.access is a rank-3 array with dimensions (basisCardinality_, D2Cardinality=3)
197 
198  if constexpr (!serendipity) {
199 
200  output.access(0, 0) = 0.5*(-1.0 + y)*y;
201  output.access(0, 1) = 0.25 - 0.5*y + x*(-0.5 + 1.*y);
202  output.access(0, 2) = 0.5*(-1.0 + x)*x;
203 
204  output.access(1, 0) = 0.5*(-1.0 + y)*y;
205  output.access(1, 1) =-0.25 + 0.5*y + x*(-0.5 + 1.*y);
206  output.access(1, 2) = 0.5*x*(1.0 + x);
207 
208  output.access(2, 0) = 0.5*y*(1.0 + y);
209  output.access(2, 1) = 0.25 + 0.5*y + x*(0.5 + 1.*y);
210  output.access(2, 2) = 0.5*x*(1.0 + x);
211 
212  output.access(3, 0) = 0.5*y*(1.0 + y);
213  output.access(3, 1) =-0.25 - 0.5*y + x*(0.5 + 1.*y);
214  output.access(3, 2) = 0.5*(-1.0 + x)*x;
215 
216  output.access(4, 0) = (1.0 - y)*y;
217  output.access(4, 1) = x*(1. - 2.*y);
218  output.access(4, 2) = (1.0 - x)*(1.0 + x);
219 
220  output.access(5, 0) = (1.0 - y)*(1.0 + y);
221  output.access(5, 1) = x*(0. - 2.*y) - 1.*y;
222  output.access(5, 2) =-x*(1.0 + x);
223 
224  output.access(6, 0) =-y*(1.0 + y);
225  output.access(6, 1) = x*(-1. - 2.*y);
226  output.access(6, 2) = (1.0 - x)*(1.0 + x);
227 
228  output.access(7, 0) = (1.0 - y)*(1.0 + y);
229  output.access(7, 1) = x*(0. - 2.*y) + 1.*y;
230  output.access(7, 2) = (1.0 - x)*x;
231 
232  output.access(8, 0) =-2.0 + 2.0*y*y;
233  output.access(8, 1) = 4*x*y;
234  output.access(8, 2) =-2.0 + 2.0*x*x;
235 
236  } else { //serendipity
237 
238  output.access(0, 0) = 0.5*(1.0 - y);
239  output.access(0, 1) = 0.25*(1.0 - 2.0*x - 2.0*y);
240  output.access(0, 2) = 0.5*(1.0 - x);
241 
242  output.access(1, 0) = 0.5*(1.0 - y);
243  output.access(1, 1) = -0.25*(1.0 + 2.0*x - 2.0*y);
244  output.access(1, 2) = 0.5*(1.0 + x);
245 
246  output.access(2, 0) = 0.5*(1.0 + y);
247  output.access(2, 1) = 0.25*(1.0 + 2.0*x + 2.0*y);
248  output.access(2, 2) = 0.5*(1.0 + x);
249 
250  output.access(3, 0) = 0.5*(1.0 + y);
251  output.access(3, 1) = -0.25*(1.0 - 2.0*x + 2.0*y);
252  output.access(3, 2) = 0.5*(1.0 - x);
253 
254  output.access(4, 0) = -(1.0 - y);
255  output.access(4, 1) = x;
256  output.access(4, 2) = 0.0;
257 
258  output.access(5, 0) = 0.0;
259  output.access(5, 1) = -y;
260  output.access(5, 2) = -(1.0 + x);
261 
262  output.access(6, 0) = -(1.0 + y);
263  output.access(6, 1) = -x;
264  output.access(6, 2) = 0.0;
265 
266  output.access(7, 0) = 0.0;
267  output.access(7, 1) = y;
268  output.access(7, 2) = -(1.0 - x);
269 
270  }
271  break;
272  }
273  case OPERATOR_D3 : {
274  if constexpr (!serendipity) {
275  const auto x = input(0);
276  const auto y = input(1);
277  output.access(0, 0) = 0.0;
278  output.access(0, 1) =-0.5 + y;
279  output.access(0, 2) =-0.5 + x;
280  output.access(0, 3) = 0.0;
281 
282  output.access(1, 0) = 0.0;
283  output.access(1, 1) =-0.5 + y;
284  output.access(1, 2) = 0.5 + x;
285  output.access(1, 3) = 0.0;
286 
287  output.access(2, 0) = 0.0;
288  output.access(2, 1) = 0.5 + y;
289  output.access(2, 2) = 0.5 + x;
290  output.access(2, 3) = 0.0;
291 
292  output.access(3, 0) = 0.0;
293  output.access(3, 1) = 0.5 + y;
294  output.access(3, 2) =-0.5 + x;
295  output.access(3, 3) = 0.0;
296 
297  output.access(4, 0) = 0.0;
298  output.access(4, 1) = 1.0 - 2.0*y;
299  output.access(4, 2) =-2.0*x;
300  output.access(4, 3) = 0.0;
301 
302  output.access(5, 0) = 0.0;
303  output.access(5, 1) =-2.0*y;
304  output.access(5, 2) =-1.0 - 2.0*x;
305  output.access(5, 3) = 0.0;
306 
307  output.access(6, 0) = 0.0;
308  output.access(6, 1) =-1.0 - 2.0*y;
309  output.access(6, 2) =-2.0*x;
310  output.access(6, 3) = 0.0;
311 
312  output.access(7, 0) = 0.0;
313  output.access(7, 1) =-2.0*y;
314  output.access(7, 2) = 1.0 - 2.0*x;
315  output.access(7, 3) = 0.0;
316 
317  output.access(8, 0) = 0.0;
318  output.access(8, 1) = 4.0*y;
319  output.access(8, 2) = 4.0*x;
320  output.access(8, 3) = 0.0;
321 
322  } else { //serendipity
323 
324  output.access(0, 0) = 0.0;
325  output.access(0, 1) =-0.5;
326  output.access(0, 2) =-0.5;
327  output.access(0, 3) = 0.0;
328 
329  output.access(1, 0) = 0.0;
330  output.access(1, 1) =-0.5;
331  output.access(1, 2) = 0.5;
332  output.access(1, 3) = 0.0;
333 
334  output.access(2, 0) = 0.0;
335  output.access(2, 1) = 0.5;
336  output.access(2, 2) = 0.5;
337  output.access(2, 3) = 0.0;
338 
339  output.access(3, 0) = 0.0;
340  output.access(3, 1) = 0.5;
341  output.access(3, 2) =-0.5;
342  output.access(3, 3) = 0.0;
343 
344  output.access(4, 0) = 0.0;
345  output.access(4, 1) = 1.0;
346  output.access(4, 2) = 0.0;
347  output.access(4, 3) = 0.0;
348 
349  output.access(5, 0) = 0.0;
350  output.access(5, 1) = 0.0;
351  output.access(5, 2) =-1.0;
352  output.access(5, 3) = 0.0;
353 
354  output.access(6, 0) = 0.0;
355  output.access(6, 1) =-1.0;
356  output.access(6, 2) = 0.0;
357  output.access(6, 3) = 0.0;
358 
359  output.access(7, 0) = 0.0;
360  output.access(7, 1) = 0.0;
361  output.access(7, 2) = 1.0;
362  output.access(7, 3) = 0.0;
363  }
364  break;
365  }
366  case OPERATOR_D4 : {
367 
368  const ordinal_type jend = output.extent(1);
369  const ordinal_type iend = output.extent(0);
370 
371  for (ordinal_type j=0;j<jend;++j)
372  for (ordinal_type i=0;i<iend;++i)
373  output.access(i, j) = 0.0;
374 
375  if constexpr (!serendipity) {
376  output.access(0, 2) = 1.0;
377  output.access(1, 2) = 1.0;
378  output.access(2, 2) = 1.0;
379  output.access(3, 2) = 1.0;
380 
381  output.access(4, 2) =-2.0;
382  output.access(5, 2) =-2.0;
383  output.access(6, 2) =-2.0;
384  output.access(7, 2) =-2.0;
385 
386  output.access(8, 2) = 4.0;
387  }
388  break;
389  }
390  case OPERATOR_MAX : {
391  const ordinal_type jend = output.extent(1);
392  const ordinal_type iend = output.extent(0);
393 
394  for (ordinal_type j=0;j<jend;++j)
395  for (ordinal_type i=0;i<iend;++i)
396  output.access(i, j) = 0.0;
397  break;
398  }
399  default: {
400  INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
401  opType != OPERATOR_GRAD &&
402  opType != OPERATOR_CURL &&
403  opType != OPERATOR_D1 &&
404  opType != OPERATOR_D2 &&
405  opType != OPERATOR_D3 &&
406  opType != OPERATOR_D4 &&
407  opType != OPERATOR_MAX,
408  ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_C2_FEM::Serial::getValues) operator is not supported");
409 
410  }
411  }
412  }
413 
414  template<bool serendipity>
415  template<typename DT,
416  typename outputValueValueType, class ...outputValueProperties,
417  typename inputPointValueType, class ...inputPointProperties>
418  void
419  Basis_HGRAD_QUAD_DEG2_FEM<serendipity>::
420  getValues( const typename DT::execution_space& space,
421  Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
422  const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
423  const EOperator operatorType ) {
424  typedef Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValueViewType;
425  typedef Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPointViewType;
426  typedef typename ExecSpace<typename inputPointViewType::execution_space,typename DT::execution_space>::ExecSpaceType ExecSpaceType;
427 
428  // Number of evaluation points = dim 0 of inputPoints
429  const auto loopSize = inputPoints.extent(0);
430  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(space, 0, loopSize);
431 
432  switch (operatorType) {
433 
434  case OPERATOR_VALUE: {
435  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_VALUE> FunctorType;
436  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
437  break;
438  }
439  case OPERATOR_GRAD:
440  case OPERATOR_D1: {
441  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_GRAD> FunctorType;
442  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
443  break;
444  }
445  case OPERATOR_CURL: {
446  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_CURL> FunctorType;
447  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
448  break;
449  }
450  case OPERATOR_DIV: {
451  INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
452  ">>> ERROR (Basis_HGRAD_QUAD_C2_FEM): DIV is invalid operator for rank-0 (scalar) functions in 2D");
453  break;
454  }
455  case OPERATOR_D2: {
456  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D2> FunctorType;
457  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
458  break;
459  }
460  case OPERATOR_D3: {
461  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, D3Cardinality=4)
462  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D3> FunctorType;
463  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
464  break;
465  }
466  case OPERATOR_D4: {
467  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, D4Cardinality=5)
468  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_D4> FunctorType;
469  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
470  break;
471  }
472  case OPERATOR_D5:
473  case OPERATOR_D6:
474  case OPERATOR_D7:
475  case OPERATOR_D8:
476  case OPERATOR_D9:
477  case OPERATOR_D10: {
478  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_MAX> FunctorType;
479  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
480  break;
481  }
482  default: {
483  INTREPID2_TEST_FOR_EXCEPTION( !( Intrepid2::isValidOperator(operatorType) ), std::invalid_argument,
484  ">>> ERROR (Basis_HGRAD_QUAD_C2_FEM): Invalid operator type");
485  }
486  }
487  }
488 
489  }
490  // -------------------------------------------------------------------------------------
491 
492 
493  template<bool serendipity, typename DT, typename OT, typename PT>
494  Basis_HGRAD_QUAD_DEG2_FEM<serendipity, DT,OT,PT>::
495  Basis_HGRAD_QUAD_DEG2_FEM() {
496  const ordinal_type spaceDim = 2;
497  this->basisCardinality_ = serendipity ? 8 : 9;
498  this->basisDegree_ = 2;
499  this->basisCellTopologyKey_ = shards::Quadrilateral<4>::key;
500  this->basisType_ = BASIS_FEM_DEFAULT;
501  this->basisCoordinates_ = COORDINATES_CARTESIAN;
502  this->functionSpace_ = FUNCTION_SPACE_HGRAD;
503 
504  {
505  // Basis-dependent intializations
506  const ordinal_type tagSize = 4; // size of DoF tag, i.e., number of fields in the tag
507  const ordinal_type posScDim = 0; // position in the tag, counting from 0, of the subcell dim
508  const ordinal_type posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
509  const ordinal_type posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
510 
511  // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
512  ordinal_type tags[36] = { 0, 0, 0, 1,
513  0, 1, 0, 1,
514  0, 2, 0, 1,
515  0, 3, 0, 1,
516  // edge midpoints
517  1, 0, 0, 1,
518  1, 1, 0, 1,
519  1, 2, 0, 1,
520  1, 3, 0, 1,
521  // quad center, not used for serendipity elements
522  2, 0, 0, 1};
523 
524  //host view
525  OrdinalTypeArray1DHost tagView(&tags[0], serendipity ? 32 : 36);
526 
527  // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
528  this->setOrdinalTagData(this->tagToOrdinal_,
529  this->ordinalToTag_,
530  tagView,
531  this->basisCardinality_,
532  tagSize,
533  posScDim,
534  posScOrd,
535  posDfOrd);
536  }
537 
538  // dofCoords on host and create its mirror view to device
539  Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
540  dofCoords("dofCoordsHost", this->basisCardinality_,spaceDim);
541 
542  dofCoords(0,0) = -1.0; dofCoords(0,1) = -1.0;
543  dofCoords(1,0) = 1.0; dofCoords(1,1) = -1.0;
544  dofCoords(2,0) = 1.0; dofCoords(2,1) = 1.0;
545  dofCoords(3,0) = -1.0; dofCoords(3,1) = 1.0;
546 
547  dofCoords(4,0) = 0.0; dofCoords(4,1) = -1.0;
548  dofCoords(5,0) = 1.0; dofCoords(5,1) = 0.0;
549  dofCoords(6,0) = 0.0; dofCoords(6,1) = 1.0;
550  dofCoords(7,0) = -1.0; dofCoords(7,1) = 0.0;
551 
552  if constexpr (!serendipity) {
553  dofCoords(8,0) = 0.0; dofCoords(8,1) = 0.0;
554  }
555 
556  this->dofCoords_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoords);
557  Kokkos::deep_copy(this->dofCoords_, dofCoords);
558  }
559 
560  template<bool serendipity, typename DT, typename OT, typename PT>
561  void
562  Basis_HGRAD_QUAD_DEG2_FEM<serendipity,DT,OT,PT>::getScratchSpaceSize(
563  ordinal_type& perTeamSpaceSize,
564  ordinal_type& perThreadSpaceSize,
565  const PointViewType inputPoints,
566  const EOperator operatorType) const {
567  perTeamSpaceSize = 0;
568  perThreadSpaceSize = 0;
569  }
570 
571  template<bool serendipity, typename DT, typename OT, typename PT>
572  KOKKOS_INLINE_FUNCTION
573  void
574  Basis_HGRAD_QUAD_DEG2_FEM<serendipity,DT,OT,PT>::getValues(
575  OutputViewType outputValues,
576  const PointViewType inputPoints,
577  const EOperator operatorType,
578  const typename Kokkos::TeamPolicy<typename DT::execution_space>::member_type& team_member,
579  const typename DT::execution_space::scratch_memory_space & scratchStorage,
580  const ordinal_type subcellDim,
581  const ordinal_type subcellOrdinal) const {
582 
583  INTREPID2_TEST_FOR_ABORT( !((subcellDim <= 0) && (subcellOrdinal == -1)),
584  ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM::getValues), The capability of selecting subsets of basis functions has not been implemented yet.");
585 
586  (void) scratchStorage; //avoid unused variable warning
587 
588  const int numPoints = inputPoints.extent(0);
589 
590  switch(operatorType) {
591  case OPERATOR_VALUE:
592  Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
593  auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
594  const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
595  using SerialValue = typename Impl::Basis_HGRAD_QUAD_DEG2_FEM<serendipity>::template Serial<OPERATOR_VALUE>;
596  SerialValue::getValues( output, input);
597  });
598  break;
599  case OPERATOR_GRAD:
600  Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
601  auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
602  const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
603  using SerialGrad = typename Impl::Basis_HGRAD_QUAD_DEG2_FEM<serendipity>::template Serial<OPERATOR_GRAD>;
604  SerialGrad::getValues( output, input);
605  });
606  break;
607  case OPERATOR_CURL:
608  Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
609  auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
610  const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
611  using SerialCurl = typename Impl::Basis_HGRAD_QUAD_DEG2_FEM<serendipity>::template Serial<OPERATOR_CURL>;
612  SerialCurl::getValues( output, input);
613  });
614  break;
615  default: {
616  INTREPID2_TEST_FOR_ABORT( true, ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM::getValues), Operator Type not supported.");
617  }
618  }
619  }
620 
621 }// namespace Intrepid2
622 #endif
Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM::OrdinalTypeArray1DHost
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
Definition: Intrepid2_Basis.hpp:125
Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM::getScratchSpaceSize
virtual void getScratchSpaceSize(ordinal_type &perTeamSpaceSize, ordinal_type &perThreadSpaceSize, const PointViewType inputPointsconst, const EOperator operatorType=OPERATOR_VALUE) const override
Return the size of the scratch space, in bytes, needed for using the team-level implementation of get...
Definition: Intrepid2_HGRAD_QUAD_C2_FEMDef.hpp:562
Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM::PointViewType
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
Definition: Intrepid2_Basis.hpp:338
Intrepid2::Impl::Basis_HGRAD_QUAD_DEG2_FEM
See Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM.
Definition: Intrepid2_HGRAD_QUAD_C2_FEM.hpp:75
Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM::Basis_HGRAD_QUAD_DEG2_FEM
Basis_HGRAD_QUAD_DEG2_FEM()
Constructor.
Definition: Intrepid2_HGRAD_QUAD_C2_FEMDef.hpp:495
Intrepid2::Basis_HGRAD_QUAD_DEG2_FEM::getValues
virtual void getValues(const ExecutionSpace &space, OutputViewType outputValues, const PointViewType inputPoints, const EOperator operatorType=OPERATOR_VALUE) const override
Evaluation of a FEM basis on a reference cell.
Definition: Intrepid2_HGRAD_QUAD_C2_FEM.hpp:178
Generated by   doxygen 1.8.5