49 #ifndef __INTREPID2_HGRAD_HEX_CN_FEM_HPP__
50 #define __INTREPID2_HGRAD_HEX_CN_FEM_HPP__
63 typedef struct Hexahedron<8> cell_topology_type;
67 template<EOperator opType>
69 template<
typename outputValueViewType,
70 typename inputPointViewType,
71 typename workViewType,
72 typename vinvViewType>
73 KOKKOS_INLINE_FUNCTION
75 getValues( outputValueViewType outputValues,
76 const inputPointViewType inputPoints,
78 const vinvViewType vinv,
79 const ordinal_type operatorDn = 0 );
82 template<
typename DeviceType, ordinal_type numPtsPerEval,
83 typename outputValueValueType,
class ...outputValueProperties,
84 typename inputPointValueType,
class ...inputPointProperties,
85 typename vinvValueType,
class ...vinvProperties>
87 getValues(
const typename DeviceType::execution_space& space,
88 Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
89 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
90 const Kokkos::DynRankView<vinvValueType, vinvProperties...> vinv,
91 const EOperator operatorType );
96 template<
typename outputValueViewType,
97 typename inputPointViewType,
98 typename vinvViewType,
99 typename workViewType,
101 ordinal_type numPtsEval>
103 outputValueViewType _outputValues;
104 const inputPointViewType _inputPoints;
105 const vinvViewType _vinv;
107 const ordinal_type _opDn;
109 KOKKOS_INLINE_FUNCTION
110 Functor( outputValueViewType outputValues_,
111 inputPointViewType inputPoints_,
114 const ordinal_type opDn_ = 0 )
115 : _outputValues(outputValues_), _inputPoints(inputPoints_),
116 _vinv(vinv_), _work(work_), _opDn(opDn_) {}
118 KOKKOS_INLINE_FUNCTION
119 void operator()(
const size_type iter)
const {
123 const auto ptRange = Kokkos::pair<ordinal_type,ordinal_type>(ptBegin, ptEnd);
124 const auto input = Kokkos::subview( _inputPoints, ptRange, Kokkos::ALL() );
126 typename workViewType::pointer_type ptr = _work.data() + _work.extent(0)*ptBegin*get_dimension_scalar(_work);
128 auto vcprop = Kokkos::common_view_alloc_prop(_work);
129 workViewType work(Kokkos::view_wrap(ptr,vcprop), (ptEnd-ptBegin)*_work.extent(0));
132 case OPERATOR_VALUE : {
133 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange );
139 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange, Kokkos::ALL() );
144 INTREPID2_TEST_FOR_ABORT(
true,
145 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_Cn_FEM::Functor) operator is not supported");
165 template<
typename DeviceType = void,
166 typename outputValueType = double,
167 typename pointValueType =
double>
169 :
public Basis<DeviceType,outputValueType,pointValueType> {
184 Kokkos::DynRankView<typename ScalarViewType::value_type,DeviceType>
vinv_;
194 const EPointType pointType = POINTTYPE_EQUISPACED);
203 const EOperator operatorType = OPERATOR_VALUE )
const override {
204 #ifdef HAVE_INTREPID2_DEBUG
205 Intrepid2::getValues_HGRAD_Args(outputValues,
212 Impl::Basis_HGRAD_HEX_Cn_FEM::
213 getValues<DeviceType,numPtsPerEval>(space,
224 #ifdef HAVE_INTREPID2_DEBUG
226 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoords) != 2, std::invalid_argument,
227 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_Cn_FEM::getDofCoords) rank = 2 required for dofCoords array");
229 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoords.extent(0)) != this->
getCardinality(), std::invalid_argument,
230 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_Cn_FEM::getDofCoords) mismatch in number of dof and 0th dimension of dofCoords array");
232 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.extent(1) != this->
getBaseCellTopology().getDimension(), std::invalid_argument,
233 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_Cn_FEM::getDofCoords) incorrect reference cell (1st) dimension in dofCoords array");
235 Kokkos::deep_copy(dofCoords, this->
dofCoords_);
241 #ifdef HAVE_INTREPID2_DEBUG
243 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoeffs) != 1, std::invalid_argument,
244 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_Cn_FEM::getdofCoeffs) rank = 1 required for dofCoeffs array");
246 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoeffs.extent(0)) != this->
getCardinality(), std::invalid_argument,
247 ">>> ERROR: (Intrepid2::Basis_HGRAD_HEX_Cn_FEM::getdofCoeffs) mismatch in number of dof and 0th dimension of dofCoeffs array");
249 Kokkos::deep_copy(dofCoeffs, 1.0);
255 return "Intrepid2_HGRAD_HEX_Cn_FEM";
264 Kokkos::DynRankView<typename ScalarViewType::const_value_type,DeviceType>
265 getVandermondeInverse() {
270 getWorkSizePerPoint(
const EOperator operatorType) {
282 BasisPtr<DeviceType,outputValueType,pointValueType>
284 if(subCellDim == 1) {
285 return Teuchos::rcp(
new
288 }
else if(subCellDim == 2) {
289 return Teuchos::rcp(
new
293 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
296 BasisPtr<typename Kokkos::HostSpace::device_type,outputValueType,pointValueType>
typename DeviceType::execution_space ExecutionSpace
(Kokkos) Execution space for basis.
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
ordinal_type getCardinality() const
Returns cardinality of the basis.
virtual void getValues(const ExecutionSpace &, OutputViewType, const PointViewType, const EOperator=OPERATOR_VALUE) const
Evaluation of a FEM basis on a reference cell.
See Intrepid2::Basis_HGRAD_HEX_Cn_FEM.
See Intrepid2::Basis_HGRAD_HEX_Cn_FEM.
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.
See Intrepid2::Basis_HGRAD_HEX_Cn_FEM.
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.
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...
Basis_HGRAD_HEX_Cn_FEM(const ordinal_type order, const EPointType pointType=POINTTYPE_EQUISPACED)
Constructor.
virtual const char * getName() const override
Returns basis name.
ordinal_type basisDegree_
Degree of the largest complete polynomial space that can be represented by the basis.
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Hexahedron cell...
Definition file for basis function of degree n for H(grad) functions on HEX cells.
virtual bool requireOrientation() const override
True if orientation is required.
virtual void getDofCoords(ScalarViewType dofCoords) const override
Returns spatial locations (coordinates) of degrees of freedom on the reference cell.
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
Kokkos::View< ordinal_type **, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray2DHost
View type for 2d host array.
static constexpr ordinal_type MaxNumPtsPerBasisEval
The maximum number of points to eval in serial mode.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
Header file for the Intrepid2::Basis_HGRAD_QUAD_Cn_FEM class.
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, DeviceType > ScalarViewType
View type for scalars.
BasisPtr< DeviceType, outputValueType, pointValueType > getSubCellRefBasis(const ordinal_type subCellDim, const ordinal_type subCellOrd) const override
returns the basis associated to a subCell.
Implementation of the locally H(grad)-compatible FEM basis of variable order on the [-1...
Kokkos::View< ordinal_type ***, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray3DHost
View type for 3d host array.
BasisPtr< typename Kokkos::HostSpace::device_type, outputValueType, pointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
EPointType pointType_
type of lattice used for creating the DoF coordinates
Implementation of the default H(grad)-compatible FEM basis of degree n on Quadrilateral cell Implemen...
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.
Kokkos::DynRankView< typename ScalarViewType::value_type, DeviceType > vinv_
inverse of Generalized Vandermonde matrix (isotropic order)
Header file for the abstract base class Intrepid2::Basis.
typename DeviceType::execution_space ExecutionSpace
(Kokkos) Execution space for basis.