48 #ifndef __INTREPID2_HVOL_TET_CN_FEM_HPP__
49 #define __INTREPID2_HVOL_TET_CN_FEM_HPP__
54 #include "Teuchos_LAPACK.hpp"
85 template<EOperator opType>
87 template<
typename outputValueViewType,
88 typename inputPointViewType,
89 typename workViewType,
90 typename vinvViewType>
91 KOKKOS_INLINE_FUNCTION
93 getValues( outputValueViewType outputValues,
94 const inputPointViewType inputPoints,
96 const vinvViewType vinv );
99 KOKKOS_INLINE_FUNCTION
101 getWorkSizePerPoint(ordinal_type order) {
102 auto cardinality = getPnCardinality<3>(order);
107 return 7*cardinality;
109 return getDkCardinality<opType,3>()*cardinality;
114 template<
typename ExecSpaceType, ordinal_type numPtsPerEval,
115 typename outputValueValueType,
class ...outputValueProperties,
116 typename inputPointValueType,
class ...inputPointProperties,
117 typename vinvValueType,
class ...vinvProperties>
119 getValues( Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
120 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
121 const Kokkos::DynRankView<vinvValueType, vinvProperties...> vinv,
122 const EOperator operatorType);
127 template<
typename outputValueViewType,
128 typename inputPointViewType,
129 typename vinvViewType,
130 typename workViewType,
132 ordinal_type numPtsEval>
134 outputValueViewType _outputValues;
135 const inputPointViewType _inputPoints;
136 const vinvViewType _vinv;
139 KOKKOS_INLINE_FUNCTION
140 Functor( outputValueViewType outputValues_,
141 inputPointViewType inputPoints_,
144 : _outputValues(outputValues_), _inputPoints(inputPoints_),
145 _vinv(vinv_), _work(work_) {}
147 KOKKOS_INLINE_FUNCTION
148 void operator()(
const size_type iter)
const {
152 const auto ptRange = Kokkos::pair<ordinal_type,ordinal_type>(ptBegin, ptEnd);
153 const auto input = Kokkos::subview( _inputPoints, ptRange, Kokkos::ALL() );
155 typename workViewType::pointer_type ptr = _work.data() + _work.extent(0)*ptBegin*get_dimension_scalar(_work);
157 auto vcprop = Kokkos::common_view_alloc_prop(_work);
158 workViewType work(Kokkos::view_wrap(ptr,vcprop), (ptEnd-ptBegin)*_work.extent(0));
161 case OPERATOR_VALUE : {
162 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange );
171 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange, Kokkos::ALL() );
176 INTREPID2_TEST_FOR_ABORT(
true,
177 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::Functor) operator is not supported");
186 template<
typename ExecSpaceType = void,
187 typename outputValueType = double,
188 typename pointValueType =
double>
190 :
public Basis<ExecSpaceType,outputValueType,pointValueType> {
199 const EPointType pointType = POINTTYPE_EQUISPACED);
213 const EOperator operatorType = OPERATOR_VALUE)
const {
214 #ifdef HAVE_INTREPID2_DEBUG
215 Intrepid2::getValues_HVOL_Args(outputValues,
222 Impl::Basis_HVOL_TET_Cn_FEM::
223 getValues<ExecSpaceType,numPtsPerEval>( outputValues,
232 #ifdef HAVE_INTREPID2_DEBUG
234 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.rank() != 2, std::invalid_argument,
235 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getDofCoords) rank = 2 required for dofCoords array");
237 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoords.extent(0)) != this->
getCardinality(), std::invalid_argument,
238 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getDofCoords) mismatch in number of dof and 0th dimension of dofCoords array");
240 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.extent(1) != this->
getBaseCellTopology().getDimension(), std::invalid_argument,
241 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getDofCoords) incorrect reference cell (1st) dimension in dofCoords array");
243 Kokkos::deep_copy(dofCoords, this->
dofCoords_);
249 #ifdef HAVE_INTREPID2_DEBUG
251 INTREPID2_TEST_FOR_EXCEPTION( dofCoeffs.rank() != 1, std::invalid_argument,
252 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getdofCoeffs) rank = 1 required for dofCoeffs array");
254 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoeffs.extent(0)) != this->
getCardinality(), std::invalid_argument,
255 ">>> ERROR: (Intrepid2::Basis_HVOL_TET_Cn_FEM::getdofCoeffs) mismatch in number of dof and 0th dimension of dofCoeffs array");
257 Kokkos::deep_copy(dofCoeffs, 1.0);
261 getVandermondeInverse( scalarViewType vinv )
const {
263 Kokkos::deep_copy(vinv, this->
vinv_);
269 return "Intrepid2_HVOL_TET_Cn_FEM";
282 Kokkos::DynRankView<scalarType,ExecSpaceType>
vinv_;
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.
Kokkos::DynRankView< scalarType, ExecSpaceType > vinv_
inverse of Generalized Vandermonde matrix, whose columns store the expansion coefficients of the noda...
virtual void getDofCoords(scalarViewType dofCoords) const
Returns spatial locations (coordinates) of degrees of freedom on the reference cell.
ordinal_type getCardinality() const
Returns cardinality of the basis.
Kokkos::View< ordinal_type ***, typename ExecSpaceType::array_layout, Kokkos::HostSpace > ordinal_type_array_3d_host
View type for 3d host array.
virtual void getDofCoeffs(scalarViewType dofCoeffs) const
Coefficients for computing degrees of freedom for Lagrangian basis If P is an element of the space sp...
virtual void getValues(outputViewType outputValues, const pointViewType inputPoints, const EOperator operatorType=OPERATOR_VALUE) const
Evaluation of a FEM basis on a reference cell.
Basis_HVOL_TET_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.
See Intrepid2::Basis_HVOL_TET_Cn_FEM.
Kokkos::DynRankView< pointValueType, Kokkos::LayoutStride, ExecSpaceType > pointViewType
View type for input points.
See Intrepid2::Basis_HVOL_TET_Cn_FEM.
virtual const char * getName() const
Returns basis name.
virtual bool requireOrientation() const
True if orientation is required.
Definition file for FEM basis functions of degree n for H(vol) functions on TET.
See Intrepid2::Basis_HVOL_TET_Cn_FEM.
ScalarTraits< pointValueType >::scalar_type scalarType
Scalar type for point values.
Implementation of the default HVOL-compatible Lagrange basis of arbitrary degree on Tetrahedron cell...
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.