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