Intrepid2
Intrepid2_HDIV_TET_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_HDIV_TET_I1_FEM_DEF_HPP__
17#define __INTREPID2_HDIV_TET_I1_FEM_DEF_HPP__
18
19namespace 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_HDIV_TET_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 // output is a rank-2 array with dimensions (basisCardinality_)
39 output.access(0, 0) = x;
40 output.access(0, 1) = (y - 1.0);
41 output.access(0, 2) = z;
42
43 output.access(1, 0) = x;
44 output.access(1, 1) = y;
45 output.access(1, 2) = z;
46
47 output.access(2, 0) = (x - 1.0);
48 output.access(2, 1) = y;
49 output.access(2, 2) = z;
50
51 output.access(3, 0) = x;
52 output.access(3, 1) = y;
53 output.access(3, 2) = (z - 1.0);
54 break;
55 }
56 case OPERATOR_DIV: {
57 // output is a rank-3 array with dimensions (basisCardinality_, spaceDim)
58 output.access(0) = 3;
59 output.access(1) = 3;
60 output.access(2) = 3;
61 output.access(3) = 3;
62 break;
63 }
64 default: {
65 INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
66 opType != OPERATOR_DIV,
67 ">>> ERROR: (Intrepid2::Basis_HDIV_TET_I1_FEM::Serial::getValues) operator is not supported");
68 }
69 }
70 }
71
72
73 template<typename DT,
74 typename outputValueValueType, class ...outputValueProperties,
75 typename inputPointValueType, class ...inputPointProperties>
76 void
77 Basis_HDIV_TET_I1_FEM::
78 getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
79 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
80 const EOperator operatorType ) {
81 typedef Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValueViewType;
82 typedef Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPointViewType;
83 typedef typename ExecSpace<typename inputPointViewType::execution_space,typename DT::execution_space>::ExecSpaceType ExecSpaceType;
84
85 // Number of evaluation points = dim 0 of inputPoints
86 const auto loopSize = inputPoints.extent(0);
87 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(0, loopSize);
88
89 switch (operatorType) {
90 case OPERATOR_VALUE: {
91 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_VALUE> FunctorType;
92 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
93 break;
94 }
95 case OPERATOR_DIV: {
96 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_DIV> FunctorType;
97 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
98 break;
99 }
100 case OPERATOR_CURL: {
101 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument,
102 ">>> ERROR (Basis_HDIV_TET_I1_FEM): CURL is invalid operator for HDIV Basis Functions");
103 break;
104 }
105
106 case OPERATOR_GRAD: {
107 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_GRAD), std::invalid_argument,
108 ">>> ERROR (Basis_HDIV_TET_I1_FEM): GRAD is invalid operator for HDIV Basis Functions");
109 break;
110 }
111 case OPERATOR_D1:
112 case OPERATOR_D2:
113 case OPERATOR_D3:
114 case OPERATOR_D4:
115 case OPERATOR_D5:
116 case OPERATOR_D6:
117 case OPERATOR_D7:
118 case OPERATOR_D8:
119 case OPERATOR_D9:
120 case OPERATOR_D10: {
121 INTREPID2_TEST_FOR_EXCEPTION( ( (operatorType == OPERATOR_D1) &&
122 (operatorType == OPERATOR_D2) &&
123 (operatorType == OPERATOR_D3) &&
124 (operatorType == OPERATOR_D4) &&
125 (operatorType == OPERATOR_D5) &&
126 (operatorType == OPERATOR_D6) &&
127 (operatorType == OPERATOR_D7) &&
128 (operatorType == OPERATOR_D8) &&
129 (operatorType == OPERATOR_D9) &&
130 (operatorType == OPERATOR_D10) ),
131 std::invalid_argument,
132 ">>> ERROR (Basis_HDIV_TET_I1_FEM): Invalid operator type");
133 break;
134 }
135 default: {
136 INTREPID2_TEST_FOR_EXCEPTION( ( (operatorType != OPERATOR_VALUE) &&
137 (operatorType != OPERATOR_GRAD) &&
138 (operatorType != OPERATOR_CURL) &&
139 (operatorType != OPERATOR_DIV) &&
140 (operatorType != OPERATOR_D1) &&
141 (operatorType != OPERATOR_D2) &&
142 (operatorType != OPERATOR_D3) &&
143 (operatorType != OPERATOR_D4) &&
144 (operatorType != OPERATOR_D5) &&
145 (operatorType != OPERATOR_D6) &&
146 (operatorType != OPERATOR_D7) &&
147 (operatorType != OPERATOR_D8) &&
148 (operatorType != OPERATOR_D9) &&
149 (operatorType != OPERATOR_D10) ),
150 std::invalid_argument,
151 ">>> ERROR (Basis_HDIV_TET_I1_FEM): Invalid operator type");
152 }
153 }
154 }
155
156 }
157
158 // -------------------------------------------------------------------------------------
159
160 template<typename DT, typename OT, typename PT>
163 this->basisCardinality_ = 4;
164 this->basisDegree_ = 1;
165 constexpr ordinal_type spaceDim = 3;
166 this->basisCellTopologyKey_ = shards::Tetrahedron<4>::key;
167 this->basisType_ = BASIS_FEM_DEFAULT;
168 this->basisCoordinates_ = COORDINATES_CARTESIAN;
169 this->functionSpace_ = FUNCTION_SPACE_HDIV;
170
171 // initialize tags
172 {
173 // Basis-dependent initializations
174 const ordinal_type tagSize = 4; // size of DoF tag, i.e., number of fields in the tag
175 const ordinal_type posScDim = 0; // position in the tag, counting from 0, of the subcell dim
176 const ordinal_type posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
177 const ordinal_type posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
178
179 // An array with local DoF tags assigned to the basis functions, in the order of their local enumeration
180 ordinal_type tags[16] = { 2, 0, 0, 1,
181 2, 1, 0, 1,
182 2, 2, 0, 1,
183 2, 3, 0, 1};
184
185 // host tags
186 OrdinalTypeArray1DHost tagView(&tags[0], 16);
187
188 // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
189 //OrdinalTypeArray2DHost ordinalToTag;
190 //OrdinalTypeArray3DHost tagToOrdinal;
191 this->setOrdinalTagData(this->tagToOrdinal_,
192 this->ordinalToTag_,
193 tagView,
194 this->basisCardinality_,
195 tagSize,
196 posScDim,
197 posScOrd,
198 posDfOrd);
199
200 //this->tagToOrdinal_ = Kokkos::create_mirror_view(typename DT::memory_space(), tagToOrdinal);
201 //Kokkos::deep_copy(this->tagToOrdinal_, tagToOrdinal);
202
203 //this->ordinalToTag_ = Kokkos::create_mirror_view(typename DT::memory_space(), ordinalToTag);
204 //Kokkos::deep_copy(this->ordinalToTag_, ordinalToTag);
205 }
206
207 // dofCoords on host and create its mirror view to device
208 Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
209 dofCoords("dofCoordsHost", this->basisCardinality_,spaceDim);
210
211 dofCoords(0,0) = 1.0/3.0; dofCoords(0,1) = 0.0; dofCoords(0,2) = 1.0/3.0;
212 dofCoords(1,0) = 1.0/3.0; dofCoords(1,1) = 1.0/3.0; dofCoords(1,2) = 1.0/3.0;
213 dofCoords(2,0) = 0.0; dofCoords(2,1) = 1.0/3.0; dofCoords(2,2) = 1.0/3.0;
214 dofCoords(3,0) = 1.0/3.0; dofCoords(3,1) = 1.0/3.0; dofCoords(3,2) = 0.0;
215
216 this->dofCoords_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoords);
217 Kokkos::deep_copy(this->dofCoords_, dofCoords);
218
219 // dofCoeffs on host and create its mirror view to device
220 Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
221 dofCoeffs("dofCoeffsHost", this->basisCardinality_,spaceDim);
222
223 // dofCoeffs are normals to faces, having magnitude equal to faces' measures
224 dofCoeffs(0,0) = 0.0; dofCoeffs(0,1) = -1.0; dofCoeffs(0,2) = 0.0;
225 dofCoeffs(1,0) = 1.0; dofCoeffs(1,1) = 1.0; dofCoeffs(1,2) = 1.0;
226 dofCoeffs(2,0) = -1.0; dofCoeffs(2,1) = 0.0; dofCoeffs(2,2) = 0.0;
227 dofCoeffs(3,0) = 0.0; dofCoeffs(3,1) = 0.0; dofCoeffs(3,2) = -1.0;
228
229 this->dofCoeffs_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoeffs);
230 Kokkos::deep_copy(this->dofCoeffs_, dofCoeffs);
231 }
232
233 template<typename DT, typename OT, typename PT>
234 void
236 ordinal_type& perTeamSpaceSize,
237 ordinal_type& perThreadSpaceSize,
238 const PointViewType inputPoints,
239 const EOperator operatorType) const {
240 perTeamSpaceSize = 0;
241 perThreadSpaceSize = 0;
242 }
243
244 template<typename DT, typename OT, typename PT>
245 KOKKOS_INLINE_FUNCTION
246 void
248 OutputViewType outputValues,
249 const PointViewType inputPoints,
250 const EOperator operatorType,
251 const typename Kokkos::TeamPolicy<typename DT::execution_space>::member_type& team_member,
252 const typename DT::execution_space::scratch_memory_space & scratchStorage,
253 const ordinal_type subcellDim,
254 const ordinal_type subcellOrdinal) const {
255
256 INTREPID2_TEST_FOR_ABORT( !((subcellDim == -1) && (subcellOrdinal == -1)),
257 ">>> ERROR: (Intrepid2::Basis_HDIV_TET_I1_FEM::getValues), The capability of selecting subsets of basis functions has not been implemented yet.");
258
259 (void) scratchStorage; //avoid unused variable warning
260
261 const int numPoints = inputPoints.extent(0);
262
263 switch(operatorType) {
264 case OPERATOR_VALUE:
265 Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
266 auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
267 const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
269 });
270 break;
271 case OPERATOR_DIV:
272 Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
273 auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
274 const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
275 Impl::Basis_HDIV_TET_I1_FEM::Serial<OPERATOR_DIV>::getValues( output, input);
276 });
277 break;
278 default: {
279 INTREPID2_TEST_FOR_ABORT( true, ">>> ERROR: (Intrepid2::Basis_HDIV_TET_I1_FEM::getValues), Operator Type not supported.");
280 }
281 }
282 }
283
284}// namespace Intrepid2
285
286#endif
287
Implementation of the default H(div)-compatible FEM basis of degree 1 on a Tetrahedron cell.
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...
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.