49 #ifndef __INTREPID2_HCURL_TRI_IN_FEM_HPP__
50 #define __INTREPID2_HCURL_TRI_IN_FEM_HPP__
57 #include "Teuchos_LAPACK.hpp"
81 #define CardinalityHCurlTri(order) (order*(order+2))
90 typedef struct Triangle<3> cell_topology_type;
95 template<EOperator opType>
97 template<
typename outputValueViewType,
98 typename inputPointViewType,
99 typename workViewType,
100 typename vinvViewType>
101 KOKKOS_INLINE_FUNCTION
103 getValues( outputValueViewType outputValues,
104 const inputPointViewType inputPoints,
106 const vinvViewType vinv );
108 KOKKOS_INLINE_FUNCTION
110 getWorkSizePerPoint(ordinal_type order) {
111 auto cardinality = CardinalityHCurlTri(order);
116 return 5*cardinality;
118 return getDkCardinality<opType,2>()*cardinality;
123 template<
typename DeviceType, ordinal_type numPtsPerEval,
124 typename outputValueValueType,
class ...outputValueProperties,
125 typename inputPointValueType,
class ...inputPointProperties,
126 typename vinvValueType,
class ...vinvProperties>
129 const typename DeviceType::execution_space& space,
130 Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
131 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
132 const Kokkos::DynRankView<vinvValueType, vinvProperties...> vinv,
133 const EOperator operatorType);
138 template<
typename outputValueViewType,
139 typename inputPointViewType,
140 typename vinvViewType,
141 typename workViewType,
143 ordinal_type numPtsEval>
145 outputValueViewType _outputValues;
146 const inputPointViewType _inputPoints;
147 const vinvViewType _coeffs;
150 KOKKOS_INLINE_FUNCTION
151 Functor( outputValueViewType outputValues_,
152 inputPointViewType inputPoints_,
153 vinvViewType coeffs_,
155 : _outputValues(outputValues_), _inputPoints(inputPoints_),
156 _coeffs(coeffs_), _work(work_) {}
158 KOKKOS_INLINE_FUNCTION
159 void operator()(
const size_type iter)
const {
163 const auto ptRange = Kokkos::pair<ordinal_type,ordinal_type>(ptBegin, ptEnd);
164 const auto input = Kokkos::subview( _inputPoints, ptRange, Kokkos::ALL() );
167 typename workViewType::pointer_type ptr = _work.data() + _work.extent(0)*ptBegin*get_dimension_scalar(_work);
169 auto vcprop = Kokkos::common_view_alloc_prop(_work);
170 workViewType work(Kokkos::view_wrap(ptr,vcprop), (ptEnd-ptBegin)*_work.extent(0));
173 case OPERATOR_VALUE : {
174 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange, Kokkos::ALL() );
178 case OPERATOR_CURL: {
179 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), ptRange);
184 INTREPID2_TEST_FOR_ABORT(
true,
185 ">>> ERROR: (Intrepid2::Basis_HCURL_TRI_In_FEM::Functor) operator is not supported");
194 template<
typename DeviceType = void,
195 typename outputValueType = double,
196 typename pointValueType =
double>
198 :
public Basis<DeviceType,outputValueType,pointValueType> {
209 const EPointType pointType = POINTTYPE_EQUISPACED);
227 const EOperator operatorType = OPERATOR_VALUE)
const override {
228 #ifdef HAVE_INTREPID2_DEBUG
229 Intrepid2::getValues_HCURL_Args(outputValues,
236 Impl::Basis_HCURL_TRI_In_FEM::
237 getValues<DeviceType,numPtsPerEval>(
248 #ifdef HAVE_INTREPID2_DEBUG
250 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoords) != 2, std::invalid_argument,
251 ">>> ERROR: (Intrepid2::Basis_HCURL_TRI_In_FEM::getDofCoords) rank = 2 required for dofCoords array");
253 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoords.extent(0)) != this->
getCardinality(), std::invalid_argument,
254 ">>> ERROR: (Intrepid2::Basis_HCURL_TRI_In_FEM::getDofCoords) mismatch in number of dof and 0th dimension of dofCoords array");
256 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.extent(1) != this->
getBaseCellTopology().getDimension(), std::invalid_argument,
257 ">>> ERROR: (Intrepid2::Basis_HCURL_TRI_In_FEM::getDofCoords) incorrect reference cell (1st) dimension in dofCoords array");
259 Kokkos::deep_copy(dofCoords, this->
dofCoords_);
265 #ifdef HAVE_INTREPID2_DEBUG
267 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoeffs) != 2, std::invalid_argument,
268 ">>> ERROR: (Intrepid2::Basis_HCURL_TRI_In_FEM::getDofCoeffs) rank = 2 required for dofCoeffs array");
270 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoeffs.extent(0)) != this->
getCardinality(), std::invalid_argument,
271 ">>> ERROR: (Intrepid2::Basis_HCURL_TRI_In_FEM::getDofCoeffs) mismatch in number of dof and 0th dimension of dofCoeffs array");
273 INTREPID2_TEST_FOR_EXCEPTION( dofCoeffs.extent(1) != this->
getBaseCellTopology().getDimension(), std::invalid_argument,
274 ">>> ERROR: (Intrepid2::Basis_HCURL_TRI_In_FEM::getDofCoeffs) incorrect reference cell (1st) dimension in dofCoeffs array");
276 Kokkos::deep_copy(dofCoeffs, this->
dofCoeffs_);
282 Kokkos::deep_copy(coeffs, this->
coeffs_);
288 return "Intrepid2_HCURL_TRI_In_FEM";
306 BasisPtr<DeviceType,outputValueType,pointValueType>
308 if(subCellDim == 1) {
309 return Teuchos::rcp(
new
313 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
316 BasisPtr<typename Kokkos::HostSpace::device_type,outputValueType,pointValueType>
324 Kokkos::DynRankView<scalarType,DeviceType>
coeffs_;
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.
ScalarTraits< pointValueType >::scalar_type scalarType
Scalar type for point values.
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.
Header file for the Intrepid2::Basis_HGRAD_TRI_Cn_FEM_ORTH class.
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...
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.
Basis_HCURL_TRI_In_FEM(const ordinal_type order, const EPointType pointType=POINTTYPE_EQUISPACED)
Constructor.
typename DeviceType::execution_space ExecutionSpace
(Kokkos) Execution space for basis.
virtual bool requireOrientation() const override
True if orientation is required.
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...
ordinal_type basisDegree_
Degree of the largest complete polynomial space that can be represented by the basis.
virtual void getDofCoords(ScalarViewType dofCoords) const override
Returns spatial locations (coordinates) of degrees of freedom on the reference cell.
EPointType pointType_
type of lattice used for creating the DoF coordinates
Header file for the Intrepid2::Basis_HVOL_LINE_Cn_FEM class.
See Intrepid2::Basis_HCURL_TRI_In_FEM.
Kokkos::View< ordinal_type **, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray2DHost
View type for 2d host array.
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, DeviceType > ScalarViewType
View type for scalars.
static constexpr ordinal_type MaxNumPtsPerBasisEval
The maximum number of points to eval in serial mode.
Kokkos::DynRankView< scalarType, DeviceType > coeffs_
expansion coefficients of the nodal basis in terms of the orthgonal one
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
BasisPtr< DeviceType, outputValueType, pointValueType > getSubCellRefBasis(const ordinal_type subCellDim, const ordinal_type subCellOrd) const override
returns the basis associated to a subCell.
See Intrepid2::Basis_HCURL_TRI_In_FEM.
Implementation of the locally HVOL-compatible FEM basis of variable order on the [-1,1] reference line cell, using Lagrange polynomials.
Kokkos::View< ordinal_type ***, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray3DHost
View type for 3d host array.
See Intrepid2::Basis_HCURL_TRI_In_FEM.
Definition file for FEM basis functions of degree n for H(curl) functions on TRI. ...
Implementation of the default H(curl)-compatible Nedelec (first kind) basis of arbitrary degree on Tr...
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
Kokkos::DynRankView< scalarType, DeviceType > dofCoeffs_
Coefficients for computing degrees of freedom for Lagrangian basis If P is an element of the space sp...
virtual const char * getName() const override
Returns basis name.
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.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.