49 #ifndef __INTREPID2_HGRAD_QUAD_CN_FEM_HPP__
50 #define __INTREPID2_HGRAD_QUAD_CN_FEM_HPP__
71 template<EOperator opType>
73 template<
typename outputValueViewType,
74 typename inputPointViewType,
75 typename workViewType,
76 typename vinvViewType>
77 KOKKOS_INLINE_FUNCTION
79 getValues( outputValueViewType outputValues,
80 const inputPointViewType inputPoints,
82 const vinvViewType vinv,
83 const ordinal_type operatorDn = 0 );
86 template<
typename ExecSpaceType, ordinal_type numPtsPerEval,
87 typename outputValueValueType,
class ...outputValueProperties,
88 typename inputPointValueType,
class ...inputPointProperties,
89 typename vinvValueType,
class ...vinvProperties>
91 getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
92 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
93 const Kokkos::DynRankView<vinvValueType, vinvProperties...> vinv,
94 const EOperator operatorType);
100 template<
typename outputValueViewType,
101 typename inputPointViewType,
102 typename vinvViewType,
103 typename workViewType,
105 ordinal_type numPtsEval>
107 outputValueViewType _outputValues;
108 const inputPointViewType _inputPoints;
109 const vinvViewType _vinv;
111 const ordinal_type _opDn;
113 KOKKOS_INLINE_FUNCTION
114 Functor( outputValueViewType outputValues_,
115 inputPointViewType inputPoints_,
118 const ordinal_type opDn_ = 0 )
119 : _outputValues(outputValues_), _inputPoints(inputPoints_),
120 _vinv(vinv_), _work(work_), _opDn(opDn_) {}
122 KOKKOS_INLINE_FUNCTION
123 void operator()(
const size_type iter)
const {
127 const auto ptRange = Kokkos::pair<ordinal_type,ordinal_type>(ptBegin, ptEnd);
128 const auto input = Kokkos::subview( _inputPoints, ptRange, Kokkos::ALL() );
130 typename workViewType::pointer_type ptr = _work.data() + _work.extent(0)*ptBegin*get_dimension_scalar(_work);
132 auto vcprop = Kokkos::common_view_alloc_prop(_work);
133 workViewType work(Kokkos::view_wrap(ptr,vcprop), (ptEnd-ptBegin)*_work.extent(0));
136 case OPERATOR_VALUE : {
137 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange );
143 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange, Kokkos::ALL() );
148 INTREPID2_TEST_FOR_ABORT(
true,
149 ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_Cn_FEM::Functor) operator is not supported");
163 template<
typename ExecSpaceType = void,
164 typename outputValueType = double,
165 typename pointValueType =
double>
167 :
public Basis<ExecSpaceType,outputValueType,pointValueType> {
179 Kokkos::DynRankView<typename scalarViewType::value_type,ExecSpaceType>
vinv_;
185 const EPointType pointType = POINTTYPE_EQUISPACED);
193 const EOperator operatorType = OPERATOR_VALUE )
const {
194 #ifdef HAVE_INTREPID2_DEBUG
195 Intrepid2::getValues_HGRAD_Args(outputValues,
202 Impl::Basis_HGRAD_QUAD_Cn_FEM::
203 getValues<ExecSpaceType,numPtsPerEval>( outputValues,
212 #ifdef HAVE_INTREPID2_DEBUG
214 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.rank() != 2, std::invalid_argument,
215 ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_Cn_FEM::getDofCoords) rank = 2 required for dofCoords array");
217 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoords.extent(0)) != this->
getCardinality(), std::invalid_argument,
218 ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_Cn_FEM::getDofCoords) mismatch in number of dof and 0th dimension of dofCoords array");
220 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.extent(1) != this->
getBaseCellTopology().getDimension(), std::invalid_argument,
221 ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_Cn_FEM::getDofCoords) incorrect reference cell (1st) dimension in dofCoords array");
223 Kokkos::deep_copy(dofCoords, this->
dofCoords_);
229 #ifdef HAVE_INTREPID2_DEBUG
231 INTREPID2_TEST_FOR_EXCEPTION( dofCoeffs.rank() != 1, std::invalid_argument,
232 ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_Cn_FEM::getdofCoeffs) rank = 1 required for dofCoeffs array");
234 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoeffs.extent(0)) != this->
getCardinality(), std::invalid_argument,
235 ">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_Cn_FEM::getdofCoeffs) mismatch in number of dof and 0th dimension of dofCoeffs array");
237 Kokkos::deep_copy(dofCoeffs, 1.0);
243 return "Intrepid2_HGRAD_QUAD_Cn_FEM";
252 Kokkos::DynRankView<typename scalarViewType::const_value_type,ExecSpaceType>
253 getVandermondeInverse()
const {
258 getWorkSizePerPoint(
const EOperator operatorType)
const {
Header file for the Intrepid2::Basis_HGRAD_LINE_Cn_FEM class.
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, ExecSpaceType > scalarViewType
View type for scalars.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
Kokkos::DynRankView< outputValueType, Kokkos::LayoutStride, ExecSpaceType > outputViewType
View type for basis value output.
shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation https://trili...
Kokkos::DynRankView< scalarType, ExecSpaceType > dofCoords_
Coordinates of degrees-of-freedom for basis functions defined in physical space.
ordinal_type getCardinality() const
Returns cardinality of the basis.
virtual bool requireOrientation() const
True if orientation is required.
Kokkos::View< ordinal_type ***, typename ExecSpaceType::array_layout, Kokkos::HostSpace > ordinal_type_array_3d_host
View type for 3d host array.
ordinal_type basisDegree_
Degree of the largest complete polynomial space that can be represented by the basis.
See Intrepid2::Basis_HGRAD_QUAD_Cn_FEM.
Basis_HGRAD_QUAD_Cn_FEM(const ordinal_type order, const EPointType pointType=POINTTYPE_EQUISPACED)
Constructor.
static constexpr ordinal_type MaxNumPtsPerBasisEval
The maximum number of points to eval in serial mode.
Kokkos::View< ordinal_type *,typename ExecSpaceType::array_layout, Kokkos::HostSpace > ordinal_type_array_1d_host
View type for 1d host array.
virtual const char * getName() const
Returns basis name.
Definition file for FEM basis functions of degree n for H(grad) functions on QUAD cells...
Kokkos::DynRankView< pointValueType, Kokkos::LayoutStride, ExecSpaceType > pointViewType
View type for input points.
See Intrepid2::Basis_HGRAD_QUAD_Cn_FEM work is a rank 1 view having the same value_type of inputPoint...
See Intrepid2::Basis_HGRAD_QUAD_Cn_FEM.
virtual void getValues(outputViewType outputValues, const pointViewType inputPoints, const EOperator operatorType=OPERATOR_VALUE) const
Evaluation of a FEM basis on a reference cell.
virtual void getDofCoords(scalarViewType dofCoords) const
Returns spatial locations (coordinates) of degrees of freedom on the reference cell.
virtual void getDofCoeffs(scalarViewType dofCoeffs) const
Coefficients for computing degrees of freedom for Lagrangian basis If P is an element of the space sp...
Implementation of the default H(grad)-compatible FEM basis of degree n on Quadrilateral cell Implemen...
Kokkos::DynRankView< typename scalarViewType::value_type, ExecSpaceType > vinv_
inverse of Generalized Vandermonde matrix (isotropic order)
Header file for the abstract base class Intrepid2::Basis.
Kokkos::View< ordinal_type **,typename ExecSpaceType::array_layout, Kokkos::HostSpace > ordinal_type_array_2d_host
View type for 2d host array.