Intrepid2
Intrepid2_HGRAD_TET_C1_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_TET_C1_FEM_DEF_HPP__
17#define __INTREPID2_HGRAD_TET_C1_FEM_DEF_HPP__
18
19namespace Intrepid2 {
20
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_HGRAD_TET_C1_FEM::Serial<opType>::
31 getValues( OutputViewType output,
32 const inputViewType input ) {
33 switch (opType) {
34 case OPERATOR_VALUE: {
35 const auto x = input(0);
36 const auto y = input(1);
37 const auto z = input(2);
38
39 // output is a rank-1 array with dimensions (basisCardinality_)
40 output.access(0) = 1.0 - x - y -z;
41 output.access(1) = x;
42 output.access(2) = y;
43 output.access(3) = z;
44 break;
45 }
46 case OPERATOR_GRAD: {
47 // output.access is a rank-2 array with dimensions (basisCardinality_,spaceDim)
48 output.access(0, 0) = -1.0;
49 output.access(0, 1) = -1.0;
50 output.access(0, 2) = -1.0;
51
52 output.access(1, 0) = 1.0;
53 output.access(1, 1) = 0.0;
54 output.access(1, 2) = 0.0;
55
56 output.access(2, 0) = 0.0;
57 output.access(2, 1) = 1.0;
58 output.access(2, 2) = 0.0;
59
60 output.access(3, 0) = 0.0;
61 output.access(3, 1) = 0.0;
62 output.access(3, 2) = 1.0;
63 break;
64 }
65 case OPERATOR_MAX: {
66 const ordinal_type jend = output.extent(1);
67 const ordinal_type iend = output.extent(0);
68
69 for (ordinal_type j=0;j<jend;++j)
70 for (ordinal_type i=0;i<iend;++i)
71 output.access(i, j) = 0.0;
72 break;
73 }
74 default: {
75 INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
76 opType != OPERATOR_GRAD &&
77 opType != OPERATOR_MAX,
78 ">>> ERROR: (Intrepid2::Basis_HGRAD_TET_C1_FEM::Serial::getValues) operator is not supported");
79 }
80 }
81 }
82
83 template<typename DT,
84 typename outputValueValueType, class ...outputValueProperties,
85 typename inputPointValueType, class ...inputPointProperties>
86 void
87 Basis_HGRAD_TET_C1_FEM::
88 getValues( const typename DT::execution_space& space,
89 Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
90 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
91 const EOperator operatorType ) {
92 typedef Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValueViewType;
93 typedef Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPointViewType;
94 typedef typename ExecSpace<typename inputPointViewType::execution_space,typename DT::execution_space>::ExecSpaceType ExecSpaceType;
95
96 // Number of evaluation points = dim 0 of inputPoints
97 const auto loopSize = inputPoints.extent(0);
98 Kokkos::RangePolicy<ExecSpaceType,Kokkos::Schedule<Kokkos::Static> > policy(space, 0, loopSize);
99
100 switch (operatorType) {
101
102 case OPERATOR_VALUE: {
103 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_VALUE> FunctorType;
104 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
105 break;
106 }
107 case OPERATOR_GRAD:
108 case OPERATOR_D1: {
109 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_GRAD> FunctorType;
110 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
111 break;
112 }
113 case OPERATOR_CURL: {
114 INTREPID2_TEST_FOR_EXCEPTION( operatorType == OPERATOR_CURL, std::invalid_argument,
115 ">>> ERROR (Basis_HGRAD_TET_C1_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D");
116 break;
117 }
118 case OPERATOR_DIV: {
119 INTREPID2_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
120 ">>> ERROR (Basis_HGRAD_TET_C1_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D");
121 break;
122 }
123 case OPERATOR_D2:
124 case OPERATOR_D3:
125 case OPERATOR_D4:
126 case OPERATOR_D5:
127 case OPERATOR_D6:
128 case OPERATOR_D7:
129 case OPERATOR_D8:
130 case OPERATOR_D9:
131 case OPERATOR_D10: {
132 typedef Functor<outputValueViewType,inputPointViewType,OPERATOR_MAX> FunctorType;
133 Kokkos::parallel_for( policy, FunctorType(outputValues, inputPoints) );
134 break;
135 }
136 default: {
137 INTREPID2_TEST_FOR_EXCEPTION( !( Intrepid2::isValidOperator(operatorType) ), std::invalid_argument,
138 ">>> ERROR (Basis_HGRAD_TET_C1_FEM): Invalid operator type");
139 }
140 }
141 }
142 }
143 // -------------------------------------------------------------------------------------
144
145 template<typename DT, typename OT, typename PT>
148 const ordinal_type spaceDim = 3;
149 this->basisCardinality_ = 4;
150 this->basisDegree_ = 1;
151 this->basisCellTopologyKey_ = shards::Tetrahedron<4>::key;
152 this->basisType_ = BASIS_FEM_DEFAULT;
153 this->basisCoordinates_ = COORDINATES_CARTESIAN;
154 this->functionSpace_ = FUNCTION_SPACE_HGRAD;
155
156 // initialize tags
157 {
158 // Basis-dependent intializations
159 const ordinal_type tagSize = 4; // size of DoF tag
160 const ordinal_type posScDim = 0; // position in the tag, counting from 0, of the subcell dim
161 const ordinal_type posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
162 const ordinal_type posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
163
164 // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
165 ordinal_type tags[16] = { 0, 0, 0, 1,
166 0, 1, 0, 1,
167 0, 2, 0, 1,
168 0, 3, 0, 1 };
169
170 // host tags
171 OrdinalTypeArray1DHost tagView(&tags[0], 16);
172
173 // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
174 //OrdinalTypeArray2DHost ordinalToTag;
175 //OrdinalTypeArray3DHost tagToOrdinal;
176 this->setOrdinalTagData(this->tagToOrdinal_,
177 this->ordinalToTag_,
178 tagView,
179 this->basisCardinality_,
180 tagSize,
181 posScDim,
182 posScOrd,
183 posDfOrd);
184
185 //this->tagToOrdinal_ = Kokkos::create_mirror_view(typename DT::memory_space(), tagToOrdinal);
186 //Kokkos::deep_copy(this->tagToOrdinal_, tagToOrdinal);
187
188 //this->ordinalToTag_ = Kokkos::create_mirror_view(typename DT::memory_space(), ordinalToTag);
189 //Kokkos::deep_copy(this->ordinalToTag_, ordinalToTag);
190 }
191
192 // dofCoords on host and create its mirror view to device
193 Kokkos::DynRankView<typename ScalarViewType::value_type,typename DT::execution_space::array_layout,Kokkos::HostSpace>
194 dofCoords("dofCoordsHost", this->basisCardinality_,spaceDim);
195
196 dofCoords(0,0) = 0.0; dofCoords(0,1) = 0.0; dofCoords(0,2) = 0.0;
197 dofCoords(1,0) = 1.0; dofCoords(1,1) = 0.0; dofCoords(1,2) = 0.0;
198 dofCoords(2,0) = 0.0; dofCoords(2,1) = 1.0; dofCoords(2,2) = 0.0;
199 dofCoords(3,0) = 0.0; dofCoords(3,1) = 0.0; dofCoords(3,2) = 1.0;
200
201 this->dofCoords_ = Kokkos::create_mirror_view(typename DT::memory_space(), dofCoords);
202 Kokkos::deep_copy(this->dofCoords_, dofCoords);
203 }
204
205 template<typename DT, typename OT, typename PT>
206 void
208 ordinal_type& perTeamSpaceSize,
209 ordinal_type& perThreadSpaceSize,
210 const PointViewType inputPoints,
211 const EOperator operatorType) const {
212 perTeamSpaceSize = 0;
213 perThreadSpaceSize = 0;
214 }
215
216 template<typename DT, typename OT, typename PT>
217 KOKKOS_INLINE_FUNCTION
218 void
220 OutputViewType outputValues,
221 const PointViewType inputPoints,
222 const EOperator operatorType,
223 const typename Kokkos::TeamPolicy<typename DT::execution_space>::member_type& team_member,
224 const typename DT::execution_space::scratch_memory_space & scratchStorage,
225 const ordinal_type subcellDim,
226 const ordinal_type subcellOrdinal) const {
227
228 INTREPID2_TEST_FOR_ABORT( !((subcellDim <= 0) && (subcellOrdinal == -1)),
229 ">>> ERROR: (Intrepid2::Basis_HGRAD_TET_C1_FEM::getValues), The capability of selecting subsets of basis functions has not been implemented yet.");
230
231 (void) scratchStorage; //avoid unused variable warning
232
233 const int numPoints = inputPoints.extent(0);
234
235 switch(operatorType) {
236 case OPERATOR_VALUE:
237 Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
238 auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
239 const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
241 });
242 break;
243 case OPERATOR_GRAD:
244 Kokkos::parallel_for (Kokkos::TeamThreadRange (team_member, numPoints), [=] (ordinal_type& pt) {
245 auto output = Kokkos::subview( outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
246 const auto input = Kokkos::subview( inputPoints, pt, Kokkos::ALL() );
247 Impl::Basis_HGRAD_TET_C1_FEM::Serial<OPERATOR_GRAD>::getValues( output, input);
248 });
249 break;
250 default: {}
251 }
252 }
253
254}// namespace Intrepid2
255#endif
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...
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.
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.