49 #ifndef Intrepid2_LegendreBasis_HVOL_TRI_h
50 #define Intrepid2_LegendreBasis_HVOL_TRI_h
52 #include <Kokkos_DynRankView.hpp>
54 #include <Intrepid2_config.h>
68 template<
class DeviceType,
class OutputScalar,
class PointScalar,
69 class OutputFieldType,
class InputPointsType>
72 using ExecutionSpace =
typename DeviceType::execution_space;
73 using ScratchSpace =
typename ExecutionSpace::scratch_memory_space;
74 using OutputScratchView = Kokkos::View<OutputScalar*,ScratchSpace,Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
75 using PointScratchView = Kokkos::View<PointScalar*, ScratchSpace,Kokkos::MemoryTraits<Kokkos::Unmanaged>>;
77 using TeamPolicy = Kokkos::TeamPolicy<ExecutionSpace>;
78 using TeamMember =
typename TeamPolicy::member_type;
82 OutputFieldType output_;
83 InputPointsType inputPoints_;
86 int numFields_, numPoints_;
88 size_t fad_size_output_;
91 : opType_(opType), output_(output), inputPoints_(inputPoints),
92 polyOrder_(polyOrder),
93 fad_size_output_(getScalarDimensionForView(output))
95 numFields_ = output.extent_int(0);
96 numPoints_ = output.extent_int(1);
97 INTREPID2_TEST_FOR_EXCEPTION(numPoints_ != inputPoints.extent_int(0), std::invalid_argument,
"point counts need to match!");
98 INTREPID2_TEST_FOR_EXCEPTION(numFields_ != (polyOrder_+1)*(polyOrder_+2)/2, std::invalid_argument,
"output field size does not match basis cardinality");
101 KOKKOS_INLINE_FUNCTION
102 void operator()(
const TeamMember & teamMember )
const
104 auto pointOrdinal = teamMember.league_rank();
105 OutputScratchView legendre_field_values_at_point, jacobi_values_at_point;
106 if (fad_size_output_ > 0) {
107 legendre_field_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_ + 1, fad_size_output_);
108 jacobi_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_ + 1, fad_size_output_);
111 legendre_field_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_ + 1);
112 jacobi_values_at_point = OutputScratchView(teamMember.team_shmem(), polyOrder_ + 1);
115 const auto & x = inputPoints_(pointOrdinal,0);
116 const auto & y = inputPoints_(pointOrdinal,1);
119 const PointScalar lambda[3] = {1. - x - y, x, y};
127 const PointScalar tLegendre = lambda[0] + lambda[1];
128 Polynomials::shiftedScaledLegendreValues(legendre_field_values_at_point, polyOrder_, lambda[1], tLegendre);
130 int fieldOrdinalOffset = 0;
131 const int max_ij_sum = polyOrder_;
132 for (
int ij_sum=0; ij_sum<=max_ij_sum; ij_sum++)
134 for (
int i=0; i<=ij_sum; i++)
136 const int j = ij_sum - i;
137 const auto & legendreValue = legendre_field_values_at_point(i);
138 const double alpha = i*2.0+1;
140 const PointScalar tJacobi = 1.0;
141 Polynomials::shiftedScaledJacobiValues(jacobi_values_at_point, alpha, polyOrder_, lambda[2], tJacobi);
143 const auto & jacobiValue = jacobi_values_at_point(j);
144 output_(fieldOrdinalOffset,pointOrdinal) = legendreValue * jacobiValue;
145 fieldOrdinalOffset++;
153 device_assert(
false);
160 size_t team_shmem_size (
int team_size)
const
164 size_t shmem_size = 0;
165 if (fad_size_output_ > 0)
166 shmem_size += 2 * OutputScratchView::shmem_size(polyOrder_ + 1, fad_size_output_);
168 shmem_size += 2 * OutputScratchView::shmem_size(polyOrder_ + 1);
184 template<
typename DeviceType,
185 typename OutputScalar = double,
186 typename PointScalar =
double>
188 :
public Basis<DeviceType,OutputScalar,PointScalar>
205 EPointType pointType_;
215 polyOrder_(polyOrder),
216 pointType_(pointType)
218 INTREPID2_TEST_FOR_EXCEPTION(pointType!=POINTTYPE_DEFAULT,std::invalid_argument,
"PointType not supported");
222 this->
basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Triangle<> >() );
227 const int degreeLength = 1;
231 int fieldOrdinalOffset = 0;
233 const int max_ij_sum = polyOrder;
234 for (
int ij_sum=0; ij_sum<=max_ij_sum; ij_sum++)
236 for (
int i=0; i<=ij_sum; i++)
238 const int j = ij_sum - i;
241 fieldOrdinalOffset++;
244 INTREPID2_TEST_FOR_EXCEPTION(fieldOrdinalOffset != this->
basisCardinality_, std::invalid_argument,
"Internal error: basis enumeration is incorrect");
251 const ordinal_type tagSize = 4;
252 const ordinal_type posScDim = 0;
253 const ordinal_type posScOrd = 1;
254 const ordinal_type posDfOrd = 2;
257 const int faceDim = 2;
259 for (ordinal_type i=0;i<cardinality;++i) {
260 tagView(i*tagSize+0) = faceDim;
261 tagView(i*tagSize+1) = 0;
262 tagView(i*tagSize+2) = i;
263 tagView(i*tagSize+3) = cardinality;
284 return "Intrepid2_LegendreBasis_HVOL_TRI";
317 const EOperator operatorType = OPERATOR_VALUE )
const override
319 auto numPoints = inputPoints.extent_int(0);
323 FunctorType functor(operatorType, outputValues, inputPoints, polyOrder_);
325 const int outputVectorSize = getVectorSizeForHierarchicalParallelism<OutputScalar>();
326 const int pointVectorSize = getVectorSizeForHierarchicalParallelism<PointScalar>();
327 const int vectorSize = std::max(outputVectorSize,pointVectorSize);
328 const int teamSize = 1;
330 auto policy = Kokkos::TeamPolicy<ExecutionSpace>(numPoints,teamSize,vectorSize);
331 Kokkos::parallel_for(
"Hierarchical_HVOL_TRI_Functor", policy, functor);
342 BasisPtr<DeviceType,OutputScalar,PointScalar>
345 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
352 virtual BasisPtr<typename Kokkos::HostSpace::device_type, OutputScalar, PointScalar>
354 using HostDeviceType =
typename Kokkos::HostSpace::device_type;
356 return Teuchos::rcp(
new HostBasisType(polyOrder_, pointType_) );
H(grad) basis on the line based on integrated Legendre polynomials.
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
const char * getName() const override
Returns basis name.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
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.
Functor for computing values for the LegendreBasis_HVOL_TRI class.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
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...
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...
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::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
BasisPtr< DeviceType, OutputScalar, PointScalar > getSubCellRefBasis(const ordinal_type subCellDim, const ordinal_type subCellOrd) const override
returns the basis associated to a subCell.
virtual void getValues(OutputViewType outputValues, const PointViewType inputPoints, const EOperator operatorType=OPERATOR_VALUE) const override
Evaluation of a FEM basis on a reference cell.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
ordinal_type basisCardinality_
Cardinality of the basis, i.e., the number of basis functions/degrees-of-freedom. ...
OrdinalTypeArray3DHost tagToOrdinal_
DoF tag to ordinal lookup table.
LegendreBasis_HVOL_TRI(int polyOrder, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
OrdinalTypeArray2DHost ordinalToTag_
"true" if tagToOrdinal_ and ordinalToTag_ have been initialized
virtual bool requireOrientation() const override
True if orientation is required.
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...
OrdinalTypeArray2DHost fieldOrdinalPolynomialDegree_
Polynomial degree for each degree of freedom. Only defined for hierarchical bases right now...
Basis defining Legendre basis on the line, a polynomial subspace of H(vol) on the line: extension to ...
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.
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.