Intrepid2
Intrepid2_HCURL_WEDGE_I1_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_HCURL_WEDGE_I1_FEM_DEF_HPP__
17 #define __INTREPID2_HCURL_WEDGE_I1_FEM_DEF_HPP__
18 
19 namespace Intrepid2 {
20 
21  // -------------------------------------------------------------------------------------
22  namespace Impl {
23 
24  template<EOperator opType>
25  template<typename OutputViewType,
26  typename inputViewType>
27  KOKKOS_INLINE_FUNCTION
28  void
29  Basis_HCURL_WEDGE_I1_FEM::Serial<opType>::
30  getValues( OutputViewType output,
31  const inputViewType input ) {
32  switch (opType) {
33  case OPERATOR_VALUE: {
34  const auto x = input(0);
35  const auto y = input(1);
36  const auto z = input(2);
37 
38  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
39  output.access(0, 0) = (1.0 - z)*(1.0 - y);
40  output.access(0, 1) = x*(1.0 - z);
41  output.access(0, 2) = 0.0;
42 
43  output.access(1, 0) = y*(z - 1.0);
44  output.access(1, 1) = x*(1.0 - z);
45  output.access(1, 2) = 0.0;
46 
47  output.access(2, 0) = y*(z - 1.0);
48  output.access(2, 1) = (1.0 - x)*(z - 1.0);
49  output.access(2, 2) = 0.0;
50 
51  output.access(3, 0) = (1.0 - y)*(1.0 + z);
52  output.access(3, 1) = x*(1.0 + z);
53  output.access(3, 2) = 0.0;
54 
55  output.access(4, 0) =-y*(1.0 + z);
56  output.access(4, 1) = x*(1.0 + z);
57  output.access(4, 2) = 0.0;
58 
59  output.access(5, 0) = -y*(1.0 + z);
60  output.access(5, 1) = (x - 1.0)*(1.0 + z);
61  output.access(5, 2) = 0.0;
62 
63  output.access(6, 0) = 0.0;
64  output.access(6, 1) = 0.0;
65  output.access(6, 2) = (1.0 - x - y);
66 
67  output.access(7, 0) = 0.0;
68  output.access(7, 1) = 0.0;
69  output.access(7, 2) = x;
70 
71  output.access(8, 0) = 0.0;
72  output.access(8, 1) = 0.0;
73  output.access(8, 2) = y;
74  break;
75  }
76  case OPERATOR_CURL: {
77  const auto x = input(0);
78  const auto y = input(1);
79  const auto z = input(2);
80  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
81  output.access(0, 0) = x;
82  output.access(0, 1) = (y - 1.0);
83  output.access(0, 2) = 2*(1.0 - z);
84 
85  output.access(1, 0) = x;
86  output.access(1, 1) = y;
87  output.access(1, 2) = 2.0*(1.0 - z);
88 
89  output.access(2, 0) = (x - 1.0);
90  output.access(2, 1) = y;
91  output.access(2, 2) = 2.0*(1.0 - z);
92 
93  output.access(3, 0) = -x;
94  output.access(3, 1) = (1.0 - y);
95  output.access(3, 2) = 2.0*(1.0 + z);
96 
97  output.access(4, 0) = -x;
98  output.access(4, 1) = -y;
99  output.access(4, 2) = 2*(1.0 + z);
100 
101  output.access(5, 0) = (1.0 - x);
102  output.access(5, 1) = -y;
103  output.access(5, 2) = 2.0*(1.0 + z);
104 
105  output.access(6, 0) =-1.0;
106  output.access(6, 1) = 1.0;
107  output.access(6, 2) = 0.0;
108 
109  output.access(7, 0) = 0.0;
110  output.access(7, 1) =-1.0;
111  output.access(7, 2) = 0.0;
112 
113  output.access(8, 0) = 1.0;
114  output.access(8, 1) = 0.0;
115  output.access(8, 2) = 0.0;
116  break;
117  }
118  default: {
119  INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
120  opType != OPERATOR_CURL,
121  ">>> ERROR: (Intrepid2::Basis_HCURL_WEDGE_I1_FEM::Serial::getValues) operator is not supported");
122 
123  }
124  }
125  }
126 
127  template<typename DT,
128  typename outputValueValueType, class ...outputValueProperties,
129  typename inputPointValueType, class ...inputPointProperties>
130  void Basis_HCURL_WEDGE_I1_FEM::
131  getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
132  const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
133  const EOperator operatorType ) {
134  typedef Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValueViewType;
135  typedef Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPointViewType;
136  typedef typename ExecSpace<typename inputPointViewType::execution_space,typename DT::execution_space>::ExecSpaceType ExecSpaceType;
137 
138  // Number of evaluation points = dim 0 of inputPoints
139  const auto loopSize = inputPoints.extent(0);
140  Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
141 
142  switch (operatorType) {
143  case OPERATOR_VALUE: {
144  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_VALUE> FunctorType;
145  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
146  break;
147  }
148  case OPERATOR_CURL: {
149  typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_CURL> FunctorType;
150  Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
151  break;
152  }
153  case OPERATOR_DIV: {
154  INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
155  ">>> ERROR (Basis_HCURL_WEDGE_I1_FEM): DIV is invalid operator for HCURL Basis Functions");
156  break;
157  }
158  case OPERATOR_GRAD: {
159  INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_GRAD), std::invalid_argument,
160  ">>> ERROR (Basis_HCURL_WEDGE_I1_FEM): GRAD is invalid operator for HCURL Basis Functions");
161  break;
162  }
163  case OPERATOR_D1:
164  case OPERATOR_D2:
165  case OPERATOR_D3:
166  case OPERATOR_D4:
167  case OPERATOR_D5:
168  case OPERATOR_D6:
169  case OPERATOR_D7:
170  case OPERATOR_D8:
171  case OPERATOR_D9:
172  case OPERATOR_D10: {
173  INTREPID2_TEST_FOR_EXCEPTION( ( (operatorType == OPERATOR_D1) ||
174  (operatorType == OPERATOR_D2) ||
175  (operatorType == OPERATOR_D3) ||
176  (operatorType == OPERATOR_D4) ||
177  (operatorType == OPERATOR_D5) ||
178  (operatorType == OPERATOR_D6) ||
179  (operatorType == OPERATOR_D7) ||
180  (operatorType == OPERATOR_D8) ||
181  (operatorType == OPERATOR_D9) ||
182  (operatorType == OPERATOR_D10) ),
183  std::invalid_argument,
184  ">>> ERROR (Basis_HCURL_WEDGE_I1_FEM): Invalid operator type");
185  break;
186  }
187  default: {
188  INTREPID2_TEST_FOR_EXCEPTION( !(Intrepid2::isValidOperator(operatorType)),
189  std::invalid_argument,
190  ">>> ERROR (Basis_HCURL_WEDGE_I1_FEM): Invalid operator type");
191  }
192  }
193  }
194 
195  }
196 
197  template<typename DT, typename OT, typename PT>
198  Basis_HCURL_WEDGE_I1_FEM<DT,OT,PT>::
199  Basis_HCURL_WEDGE_I1_FEM() {
200  const ordinal_type spaceDim = 3;
201  this->basisCardinality_ = 9;
202  this->basisDegree_ = 1;
203  this->basisCellTopologyKey_ = shards::Wedge<6>::key;
204  this->basisType_ = BASIS_FEM_DEFAULT;
205  this->basisCoordinates_ = COORDINATES_CARTESIAN;
206  this->functionSpace_ = FUNCTION_SPACE_HCURL;
207 
208 
209  // initialize tags
210  {
211  // Basis-dependent intializations
212  const ordinal_type tagSize = 4; // size of DoF tag
213  const ordinal_type posScDim = 0; // position in the tag, counting from 0, of the subcell dim
214  const ordinal_type posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
215  const ordinal_type posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
216 
217  // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
218  ordinal_type tags[36] = {
219  1, 0, 0, 1,
220  1, 1, 0, 1,
221  1, 2, 0, 1,
222  1, 3, 0, 1,
223  1, 4, 0, 1,
224  1, 5, 0, 1,
225  1, 6, 0, 1,
226  1, 7, 0, 1,
227  1, 8, 0, 1 };
228 
229  //host tags
230  OrdinalTypeArray1DHost tagView(&tags[0],36);
231 
232  // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
233  this->setOrdinalTagData(this->tagToOrdinal_,
234  this->ordinalToTag_,
235  tagView,
236  this->basisCardinality_,
237  tagSize,
238  posScDim,
239  posScOrd,
240  posDfOrd);
241  }
242  // dofCoords on host and create its mirror view to device
243  Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
244  dofCoords("dofCoordsHost", this->basisCardinality_,spaceDim);
245 
246  dofCoords(0,0) = 0.5; dofCoords(0,1) = 0.0; dofCoords(0,2) = -1.0;
247  dofCoords(1,0) = 0.5; dofCoords(1,1) = 0.5; dofCoords(1,2) = -1.0;
248  dofCoords(2,0) = 0.0; dofCoords(2,1) = 0.5; dofCoords(2,2) = -1.0;
249  dofCoords(3,0) = 0.5; dofCoords(3,1) = 0.0; dofCoords(3,2) = 1.0;
250  dofCoords(4,0) = 0.5; dofCoords(4,1) = 0.5; dofCoords(4,2) = 1.0;
251  dofCoords(5,0) = 0.0; dofCoords(5,1) = 0.5; dofCoords(5,2) = 1.0;
252  dofCoords(6,0) = 0.0; dofCoords(6,1) = 0.0; dofCoords(6,2) = 0.0;
253  dofCoords(7,0) = 1.0; dofCoords(7,1) = 0.0; dofCoords(7,2) = 0.0;
254  dofCoords(8,0) = 0.0; dofCoords(8,1) = 1.0; dofCoords(8,2) = 0.0;
255 
256  this->dofCoords_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoords);
257  Kokkos::deep_copy(this->dofCoords_, dofCoords);
258 
259  // dofCoeffs on host and create its mirror view to device
260  Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
261  dofCoeffs("dofCoeffsHost", this->basisCardinality_,spaceDim);
262 
263  // for HCURL_WEDGE_I1 dofCoeffs are the tangents on the hexahedron edges
264  dofCoeffs(0,0) = 0.5; dofCoeffs(0,1) = 0.0; dofCoeffs(0,2) = 0.0;
265  dofCoeffs(1,0) = -0.5; dofCoeffs(1,1) = 0.5; dofCoeffs(1,2) = 0.0;
266  dofCoeffs(2,0) = 0.0; dofCoeffs(2,1) = -0.5; dofCoeffs(2,2) = 0.0;
267  dofCoeffs(3,0) = 0.5; dofCoeffs(3,1) = 0.0; dofCoeffs(3,2) = 0.0;
268  dofCoeffs(4,0) = -0.5; dofCoeffs(4,1) = 0.5; dofCoeffs(4,2) = 0.0;
269  dofCoeffs(5,0) = 0.0; dofCoeffs(5,1) = -0.5; dofCoeffs(5,2) = 0.0;
270  dofCoeffs(6,0) = 0.0; dofCoeffs(6,1) = 0.0; dofCoeffs(6,2) = 1.0;
271  dofCoeffs(7,0) = 0.0; dofCoeffs(7,1) = 0.0; dofCoeffs(7,2) = 1.0;
272  dofCoeffs(8,0) = 0.0; dofCoeffs(8,1) = 0.0; dofCoeffs(8,2) = 1.0;
273 
274  this->dofCoeffs_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoeffs);
275  Kokkos::deep_copy(this->dofCoeffs_, dofCoeffs);
276 
277  }
278 
279  template<typename DT, typename OT, typename PT>
280  void
281  Basis_HCURL_WEDGE_I1_FEM<DT,OT,PT>::getScratchSpaceSize(
282  ordinal_type& perTeamSpaceSize,
283  ordinal_type& perThreadSpaceSize,
284  const PointViewType inputPoints,
285  const EOperator operatorType) const {
286  perTeamSpaceSize = 0;
287  perThreadSpaceSize = 0;
288  }
289 
290  template<typename DT, typename OT, typename PT>
291  KOKKOS_INLINE_FUNCTION
292  void
293  Basis_HCURL_WEDGE_I1_FEM<DT,OT,PT>::getValues(
294  OutputViewType outputValues,
295  const PointViewType inputPoints,
296  const EOperator operatorType,
297  const typename Kokkos::TeamPolicy<typename DT::execution_space>::member_type& team_member,
298  const typename DT::execution_space::scratch_memory_space & scratchStorage,
299  const ordinal_type subcellDim,
300  const ordinal_type subcellOrdinal) const {
301 
302  INTREPID2_TEST_FOR_ABORT( !((subcellDim == -1) && (subcellOrdinal == -1)),
303  ">>> ERROR: (Intrepid2::Basis_HCURL_WEDGE_I1_FEM::getValues), The capability of selecting subsets of basis functions has not been implemented yet.");
304 
305  (void) scratchStorage; //avoid unused variable warning
306 
307  const int numPoints = inputPoints.extent(0);
308 
309  switch(operatorType) {
310  case OPERATOR_VALUE:
311  Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
312  auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
313  const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
314  Impl::Basis_HCURL_WEDGE_I1_FEM::Serial<OPERATOR_VALUE>::getValues( output, input);
315  });
316  break;
317  case OPERATOR_CURL:
318  Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
319  auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
320  const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
321  Impl::Basis_HCURL_WEDGE_I1_FEM::Serial<OPERATOR_CURL>::getValues( output, input);
322  });
323  break;
324  default: {
325  INTREPID2_TEST_FOR_ABORT( true, ">>> ERROR: (Intrepid2::Basis_HCURL_WEDGE_I1_FEM::getValues), Operator Type not supported.");
326  }
327  }
328  }
329 
330 }// namespace Intrepid2
331 
332 #endif
Intrepid2::Basis_HCURL_WEDGE_I1_FEM
Implementation of the default H(curl)-compatible FEM basis of degree 1 on Wedge cell.
Definition: Intrepid2_HCURL_WEDGE_I1_FEM.hpp:153
Intrepid2::Basis_HCURL_WEDGE_I1_FEM::Basis_HCURL_WEDGE_I1_FEM
Basis_HCURL_WEDGE_I1_FEM()
Constructor.
Definition: Intrepid2_HCURL_WEDGE_I1_FEMDef.hpp:199
Generated by   doxygen 1.8.5