49 #ifndef __INTREPID2_HGRAD_LINE_CN_FEM_HPP__
50 #define __INTREPID2_HGRAD_LINE_CN_FEM_HPP__
56 #include "Teuchos_LAPACK.hpp"
83 typedef struct Line<2> cell_topology_type;
87 template<EOperator opType>
89 template<
typename outputValueViewType,
90 typename inputPointViewType,
91 typename workViewType,
92 typename vinvViewType>
93 KOKKOS_INLINE_FUNCTION
95 getValues( outputValueViewType outputValues,
96 const inputPointViewType inputPoints,
98 const vinvViewType vinv,
99 const ordinal_type operatorDn = 0 );
102 template<
typename DeviceType, ordinal_type numPtsPerEval,
103 typename outputValueValueType,
class ...outputValueProperties,
104 typename inputPointValueType,
class ...inputPointProperties,
105 typename vinvValueType,
class ...vinvProperties>
107 getValues(
const typename DeviceType::execution_space& space,
108 Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
109 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
110 const Kokkos::DynRankView<vinvValueType, vinvProperties...> vinv,
111 const EOperator operatorType );
116 template<
typename outputValueViewType,
117 typename inputPointViewType,
118 typename vinvViewType,
119 typename workViewType,
121 ordinal_type numPtsEval>
123 outputValueViewType _outputValues;
124 const inputPointViewType _inputPoints;
125 const vinvViewType _vinv;
127 const ordinal_type _opDn;
129 KOKKOS_INLINE_FUNCTION
130 Functor( outputValueViewType outputValues_,
131 inputPointViewType inputPoints_,
134 const ordinal_type opDn_ = 0 )
135 : _outputValues(outputValues_), _inputPoints(inputPoints_),
136 _vinv(vinv_), _work(work_), _opDn(opDn_) {}
138 KOKKOS_INLINE_FUNCTION
139 void operator()(
const size_type iter)
const {
143 const auto ptRange = Kokkos::pair<ordinal_type,ordinal_type>(ptBegin, ptEnd);
144 const auto input = Kokkos::subview( _inputPoints, ptRange, Kokkos::ALL() );
146 typename workViewType::pointer_type ptr = _work.data() + _work.extent(0)*ptBegin*get_dimension_scalar(_work);
148 auto vcprop = Kokkos::common_view_alloc_prop(_work);
149 workViewType work(Kokkos::view_wrap(ptr,vcprop), (ptEnd-ptBegin)*_work.extent(0));
152 case OPERATOR_VALUE : {
153 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange );
158 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange, Kokkos::ALL() );
163 INTREPID2_TEST_FOR_ABORT(
true,
164 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_Cn_FEM::Functor) operator is not supported");
173 template<
typename DeviceType = void,
174 typename outputValueType = double,
175 typename pointValueType =
double>
177 :
public Basis<DeviceType,outputValueType,pointValueType> {
196 Kokkos::DynRankView<typename ScalarViewType::value_type,DeviceType>
vinv_;
197 EPointType pointType_;
202 const EPointType pointType = POINTTYPE_EQUISPACED);
211 const EOperator operatorType = OPERATOR_VALUE )
const override {
212 #ifdef HAVE_INTREPID2_DEBUG
213 Intrepid2::getValues_HGRAD_Args(outputValues,
219 constexpr ordinal_type numPtsPerEval = 1;
220 Impl::Basis_HGRAD_LINE_Cn_FEM::
221 getValues<DeviceType,numPtsPerEval>(space,
231 #ifdef HAVE_INTREPID2_DEBUG
233 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoords) != 2, std::invalid_argument,
234 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_Cn_FEM::getDofCoords) rank = 2 required for dofCoords array");
236 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoords.extent(0)) != this->
getCardinality(), std::invalid_argument,
237 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_Cn_FEM::getDofCoords) mismatch in number of dof and 0th dimension of dofCoords array");
239 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.extent(1) != this->
getBaseCellTopology().getDimension(), std::invalid_argument,
240 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_Cn_FEM::getDofCoords) incorrect reference cell (1st) dimension in dofCoords array");
242 Kokkos::deep_copy(dofCoords, this->
dofCoords_);
248 #ifdef HAVE_INTREPID2_DEBUG
250 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoeffs) != 1, std::invalid_argument,
251 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_Cn_FEM::getdofCoeffs) rank = 1 required for dofCoeffs array");
253 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoeffs.extent(0)) != this->
getCardinality(), std::invalid_argument,
254 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_Cn_FEM::getdofCoeffs) mismatch in number of dof and 0th dimension of dofCoeffs array");
256 Kokkos::deep_copy(dofCoeffs, 1.0);
262 return "Intrepid2_HGRAD_LINE_Cn_FEM";
268 Kokkos::deep_copy(vinv, this->
vinv_);
271 Kokkos::DynRankView<typename ScalarViewType::const_value_type,DeviceType>
272 getVandermondeInverse()
const {
277 getWorkSizePerPoint(
const EOperator operatorType)
const {
285 virtual HostBasisPtr<outputValueType,pointValueType>
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
See Intrepid2::Basis_HGRAD_LINE_Cn_FEM.
virtual const char * getName() const override
Returns basis name.
ordinal_type getCardinality() const
Returns cardinality of the basis.
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.
virtual void getValues(const ExecutionSpace &, OutputViewType, const PointViewType, const EOperator=OPERATOR_VALUE) const
Evaluation of a FEM basis on a reference cell.
Kokkos::DynRankView< typename ScalarViewType::value_type, DeviceType > vinv_
inverse of Generalized Vandermonde matrix, whose columns store the expansion coefficients of the noda...
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation https://trili...
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, DeviceType > ScalarViewType
View type for scalars.
ordinal_type basisDegree_
Degree of the largest complete polynomial space that can be represented by the basis.
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, DeviceType > ScalarViewType
View type for scalars.
virtual void getDofCoords(ScalarViewType dofCoords) const override
Returns spatial locations (coordinates) of degrees of freedom on the reference cell.
Kokkos::View< ordinal_type **, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray2DHost
View type for 2d host array.
Header file for the Intrepid2::Basis_HGRAD_LINE_Cn_FEM_JACOBI class.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
See Intrepid2::Basis_HGRAD_LINE_Cn_FEM.
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Coefficients for computing degrees of freedom for Lagrangian basis If P is an element of the space sp...
typename DeviceType::execution_space ExecutionSpace
(Kokkos) Execution space for basis.
Implementation of the locally H(grad)-compatible FEM basis of variable order on the [-1...
See Intrepid2::Basis_HGRAD_LINE_Cn_FEM.
Definition file for FEM basis functions of degree n for H(grad) functions on LINE.
Kokkos::View< ordinal_type ***, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray3DHost
View type for 3d host array.
Basis_HGRAD_LINE_Cn_FEM(const ordinal_type order, const EPointType pointType=POINTTYPE_EQUISPACED)
Constructor.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
virtual HostBasisPtr< outputValueType, pointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
Kokkos::DynRankView< scalarType, DeviceType > dofCoords_
Coordinates of degrees-of-freedom for basis functions defined in physical space.
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.