49 #ifndef Intrepid2_IntegratedLegendreBasis_HGRAD_LINE_h
50 #define Intrepid2_IntegratedLegendreBasis_HGRAD_LINE_h
52 #include <Kokkos_DynRankView.hpp>
54 #include <Intrepid2_config.h>
67 template<
class DeviceType,
class OutputScalar,
class PointScalar,
68 class OutputFieldType,
class InputPointsType>
71 using ExecutionSpace =
typename DeviceType::execution_space;
72 using ScratchSpace =
typename ExecutionSpace::scratch_memory_space;
73 using OutputScratchView = Kokkos::View<OutputScalar*,ScratchSpace,Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
74 using PointScratchView = Kokkos::View<PointScalar*, ScratchSpace,Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
76 using TeamPolicy = Kokkos::TeamPolicy<ExecutionSpace>;
77 using TeamMember =
typename TeamPolicy::member_type;
81 OutputFieldType output_;
82 InputPointsType inputPoints_;
85 bool defineVertexFunctions_;
86 int numFields_, numPoints_;
88 size_t fad_size_output_;
91 int polyOrder,
bool defineVertexFunctions)
92 : opType_(opType), output_(output), inputPoints_(inputPoints),
93 polyOrder_(polyOrder), defineVertexFunctions_(defineVertexFunctions),
94 fad_size_output_(getScalarDimensionForView(output))
96 numFields_ = output.extent_int(0);
97 numPoints_ = output.extent_int(1);
98 INTREPID2_TEST_FOR_EXCEPTION(numPoints_ != inputPoints.extent_int(0), std::invalid_argument,
"point counts need to match!");
99 INTREPID2_TEST_FOR_EXCEPTION(numFields_ != polyOrder_+1, std::invalid_argument,
"output field size does not match basis cardinality");
102 KOKKOS_INLINE_FUNCTION
103 void operator()(
const TeamMember & teamMember )
const
105 auto pointOrdinal = teamMember.league_rank();
106 OutputScratchView field_values_at_point;
107 if (fad_size_output_ > 0) {
108 field_values_at_point = OutputScratchView(teamMember.team_shmem(), numFields_, fad_size_output_);
111 field_values_at_point = OutputScratchView(teamMember.team_shmem(), numFields_);
114 const auto & input_x = inputPoints_(pointOrdinal,0);
115 const bool taking_derivative = (opType_ != OPERATOR_VALUE);
116 const bool callingShiftedScaledLegendre = (opType_ == OPERATOR_VALUE) || (opType_ == OPERATOR_GRAD) || (opType_ == OPERATOR_D1);
119 const PointScalar x = callingShiftedScaledLegendre ? PointScalar((input_x + 1.0)/2.0) : PointScalar(input_x);
120 const double legendreScaling = 1.0;
121 const double outputScaling = taking_derivative ? 0.5 : 1.0;
128 Polynomials::shiftedScaledIntegratedLegendreValues(field_values_at_point, polyOrder_, x, legendreScaling);
133 if (defineVertexFunctions_)
135 field_values_at_point(0) = 1. - x;
136 field_values_at_point(1) = x;
143 Polynomials::shiftedScaledIntegratedLegendreValues_dx(field_values_at_point, polyOrder_, x, legendreScaling);
148 if (defineVertexFunctions_)
150 field_values_at_point(0) = -1.0;
151 field_values_at_point(1) = 1.0;
164 auto derivativeOrder = getOperatorOrder(opType_) - 1;
165 Polynomials::legendreDerivativeValues(field_values_at_point, polyOrder_, x, derivativeOrder);
170 OutputScalar Pn_minus_one = field_values_at_point(1);
171 for (
int fieldOrdinal=2; fieldOrdinal<numFields_; fieldOrdinal++)
173 OutputScalar Pn = field_values_at_point(fieldOrdinal);
174 field_values_at_point(fieldOrdinal) = Pn_minus_one;
178 if (numFields_ >= 1) field_values_at_point(0) = 0.0;
179 if (numFields_ >= 2) field_values_at_point(1) = 0.0;
187 device_assert(
false);
191 for (
int fieldOrdinal=0; fieldOrdinal<numFields_; fieldOrdinal++)
194 output_.access(fieldOrdinal,pointOrdinal,0) = outputScaling * field_values_at_point(fieldOrdinal);
201 size_t team_shmem_size (
int team_size)
const
204 size_t shmem_size = 0;
205 if (fad_size_output_ > 0)
206 shmem_size += OutputScratchView::shmem_size(numFields_, fad_size_output_);
208 shmem_size += OutputScratchView::shmem_size(numFields_);
231 template<
typename DeviceType,
232 typename OutputScalar = double,
233 typename PointScalar = double,
234 bool defineVertexFunctions =
true,
235 bool useMinusOneToOneReferenceElement =
true>
237 :
public Basis<DeviceType,OutputScalar,PointScalar>
254 bool defineVertexFunctions_;
255 EPointType pointType_;
269 polyOrder_(polyOrder),
270 pointType_(pointType)
272 INTREPID2_TEST_FOR_EXCEPTION(pointType!=POINTTYPE_DEFAULT,std::invalid_argument,
"PointType not supported");
276 this->
basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Line<2> >() );
281 const int degreeLength = 1;
292 if (defineVertexFunctions)
303 const ordinal_type tagSize = 4;
304 const ordinal_type posScDim = 0;
305 const ordinal_type posScOrd = 1;
306 const ordinal_type posDfOrd = 2;
310 if (defineVertexFunctions) {
312 const ordinal_type v0 = 0;
313 tagView(v0*tagSize+0) = 0;
314 tagView(v0*tagSize+1) = 0;
315 tagView(v0*tagSize+2) = 0;
316 tagView(v0*tagSize+3) = 1;
318 const ordinal_type v1 = 1;
319 tagView(v1*tagSize+0) = 0;
320 tagView(v1*tagSize+1) = 1;
321 tagView(v1*tagSize+2) = 0;
322 tagView(v1*tagSize+3) = 1;
324 const ordinal_type iend = cardinality - 2;
325 for (ordinal_type i=0;i<iend;++i) {
326 const auto e = i + 2;
327 tagView(e*tagSize+0) = 1;
328 tagView(e*tagSize+1) = 0;
329 tagView(e*tagSize+2) = i;
330 tagView(e*tagSize+3) = iend;
334 for (ordinal_type i=0;i<cardinality;++i) {
335 tagView(i*tagSize+0) = 1;
336 tagView(i*tagSize+1) = 0;
337 tagView(i*tagSize+2) = i;
338 tagView(i*tagSize+3) = cardinality;
347 this->basisCardinality_,
360 return "Intrepid2_IntegratedLegendreBasis_HGRAD_LINE";
393 const EOperator operatorType = OPERATOR_VALUE )
const override
395 auto numPoints = inputPoints.extent_int(0);
399 FunctorType functor(operatorType, outputValues, inputPoints, polyOrder_, defineVertexFunctions);
401 const int outputVectorSize = getVectorSizeForHierarchicalParallelism<OutputScalar>();
402 const int pointVectorSize = getVectorSizeForHierarchicalParallelism<PointScalar>();
403 const int vectorSize = std::max(outputVectorSize,pointVectorSize);
404 const int teamSize = 1;
406 auto policy = Kokkos::TeamPolicy<ExecutionSpace>(numPoints,teamSize,vectorSize);
407 Kokkos::parallel_for(
"Hierarchical_HGRAD_LINE_Functor", policy, functor);
414 virtual BasisPtr<typename Kokkos::HostSpace::device_type, OutputScalar, PointScalar>
416 using HostDeviceType =
typename Kokkos::HostSpace::device_type;
418 return Teuchos::rcp(
new HostBasisType(polyOrder_, pointType_) );
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
virtual bool requireOrientation() const override
True if orientation is required.
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
virtual void getValues(const ExecutionSpace &, OutputViewType, const PointViewType, const EOperator=OPERATOR_VALUE) const
Evaluation of a FEM basis on a reference cell.
EBasis basisType_
Type of the basis.
virtual void getValues(OutputViewType outputValues, const PointViewType inputPoints, const EOperator operatorType=OPERATOR_VALUE) const override
Evaluation of a FEM basis on a reference cell.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, DeviceType > ScalarViewType
View type for scalars.
Free functions, callable from device code, that implement various polynomials useful in basis definit...
IntegratedLegendreBasis_HGRAD_LINE(int polyOrder, EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
EFunctionSpace functionSpace_
The function space in which the basis is defined.
ordinal_type basisDegree_
Degree of the largest complete polynomial space that can be represented by the basis.
Header function for Intrepid2::Util class and other utility functions.
Kokkos::View< ordinal_type **, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray2DHost
View type for 2d host array.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
Basis defining integrated Legendre basis on the line, a polynomial subspace of H(grad) on the line...
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
ordinal_type basisCardinality_
Cardinality of the basis, i.e., the number of basis functions/degrees-of-freedom. ...
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
OrdinalTypeArray3DHost tagToOrdinal_
DoF tag to ordinal lookup table.
const char * getName() const override
Returns basis name.
OrdinalTypeArray2DHost ordinalToTag_
"true" if tagToOrdinal_ and ordinalToTag_ have been initialized
ECoordinates basisCoordinates_
The coordinate system for which the basis is defined.
shards::CellTopology basisCellTopology_
Base topology of the cells for which the basis is defined. See the Shards package for definition of b...
OrdinalTypeArray2DHost fieldOrdinalH1PolynomialDegree_
H^1 polynomial degree for each degree of freedom. Only defined for hierarchical bases right now...
virtual BasisPtr< typename Kokkos::HostSpace::device_type, OutputScalar, PointScalar > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
OrdinalTypeArray2DHost fieldOrdinalPolynomialDegree_
Polynomial degree for each degree of freedom. Only defined for hierarchical bases right now...
void setOrdinalTagData(OrdinalTypeView3D &tagToOrdinal, OrdinalTypeView2D &ordinalToTag, const OrdinalTypeView1D tags, const ordinal_type basisCard, const ordinal_type tagSize, const ordinal_type posScDim, const ordinal_type posScOrd, const ordinal_type posDfOrd)
Fills ordinalToTag_ and tagToOrdinal_ by basis-specific tag data.
Functor for computing values for the IntegratedLegendreBasis_HGRAD_LINE class.
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
Header file for the abstract base class Intrepid2::Basis.
typename DeviceType::execution_space ExecutionSpace
(Kokkos) Execution space for basis.