49 #ifndef __INTREPID2_HGRAD_LINE_C1_FEM_HPP__
50 #define __INTREPID2_HGRAD_LINE_C1_FEM_HPP__
85 typedef struct Line<2> cell_topology_type;
89 template<EOperator opType>
91 template<
typename OutputViewType,
92 typename inputViewType>
93 KOKKOS_INLINE_FUNCTION
95 getValues( OutputViewType output,
96 const inputViewType input );
100 template<
typename DeviceType,
101 typename outputValueValueType,
class ...outputValueProperties,
102 typename inputPointValueType,
class ...inputPointProperties>
104 getValues(
const typename DeviceType::execution_space& space,
105 Kokkos::DynRankView<outputValueValueType,outputValueProperties...> outputValues,
106 const Kokkos::DynRankView<inputPointValueType, inputPointProperties...> inputPoints,
107 const EOperator operatorType);
112 template<
typename outputValueViewType,
113 typename inputPointViewType,
116 outputValueViewType _outputValues;
117 const inputPointViewType _inputPoints;
119 KOKKOS_INLINE_FUNCTION
120 Functor( outputValueViewType outputValues_,
121 inputPointViewType inputPoints_ )
122 : _outputValues(outputValues_), _inputPoints(inputPoints_) {}
124 KOKKOS_INLINE_FUNCTION
125 void operator()(
const ordinal_type pt)
const {
127 case OPERATOR_VALUE : {
128 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), pt );
129 const auto input = Kokkos::subview( _inputPoints, pt, Kokkos::ALL() );
134 case OPERATOR_MAX : {
135 auto output = Kokkos::subview( _outputValues, Kokkos::ALL(), pt, Kokkos::ALL() );
136 const auto input = Kokkos::subview( _inputPoints, pt, Kokkos::ALL() );
141 INTREPID2_TEST_FOR_ABORT( opType != OPERATOR_VALUE &&
142 opType != OPERATOR_GRAD &&
143 opType != OPERATOR_MAX,
144 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_C1_FEM::Serial::getValues) operator is not supported");
153 template<
typename DeviceType = void,
154 typename outputValueType = double,
155 typename pointValueType =
double>
157 :
public Basis<DeviceType,outputValueType,pointValueType> {
181 const EOperator operatorType = OPERATOR_VALUE )
const override {
182 #ifdef HAVE_INTREPID2_DEBUG
184 Intrepid2::getValues_HGRAD_Args(outputValues,
190 Impl::Basis_HGRAD_LINE_C1_FEM::
191 getValues<DeviceType>(space,
200 #ifdef HAVE_INTREPID2_DEBUG
202 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoords) != 2, std::invalid_argument,
203 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_C1_FEM::getDofCoords) rank = 2 required for dofCoords array");
205 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoords.extent(0)) != this->
basisCardinality_, std::invalid_argument,
206 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_C1_FEM::getDofCoords) mismatch in number of dof and 0th dimension of dofCoords array");
208 INTREPID2_TEST_FOR_EXCEPTION( dofCoords.extent(1) != this->
basisCellTopology_.getDimension(), std::invalid_argument,
209 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_C1_FEM::getDofCoords) incorrect reference cell (1st) dimension in dofCoords array");
211 Kokkos::deep_copy(dofCoords, this->
dofCoords_);
217 #ifdef HAVE_INTREPID2_DEBUG
219 INTREPID2_TEST_FOR_EXCEPTION( rank(dofCoeffs) != 1, std::invalid_argument,
220 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_C1_FEM::getdofCoeffs) rank = 1 required for dofCoeffs array");
222 INTREPID2_TEST_FOR_EXCEPTION( static_cast<ordinal_type>(dofCoeffs.extent(0)) != this->
getCardinality(), std::invalid_argument,
223 ">>> ERROR: (Intrepid2::Basis_HGRAD_LINE_C1_FEM::getdofCoeffs) mismatch in number of dof and 0th dimension of dofCoeffs array");
225 Kokkos::deep_copy(dofCoeffs, 1.0);
231 return "Intrepid2_HGRAD_LINE_C1_FEM";
234 BasisPtr<typename Kokkos::HostSpace::device_type,outputValueType,pointValueType>
See Intrepid2::Basis_HGRAD_LINE_C1_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...
Kokkos::View< ordinal_type *, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray1DHost
View type for 1d host array.
Definition file for FEM basis functions of degree 1 for H(grad) functions on a Line.
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.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
An abstract base class that defines interface for concrete basis implementations for Finite Element (...
virtual void getDofCoords(ScalarViewType dofCoords) const override
Returns spatial locations (coordinates) of degrees of freedom on the reference cell.
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.
Kokkos::DynRankView< PointValueType, Kokkos::LayoutStride, DeviceType > PointViewType
View type for input points.
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.
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...
Kokkos::View< ordinal_type **, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray2DHost
View type for 2d host array.
Kokkos::DynRankView< OutputValueType, Kokkos::LayoutStride, DeviceType > OutputViewType
View type for basis value output.
typename DeviceType::execution_space ExecutionSpace
(Kokkos) Execution space for basis.
Kokkos::DynRankView< scalarType, Kokkos::LayoutStride, DeviceType > ScalarViewType
View type for scalars.
ordinal_type basisCardinality_
Cardinality of the basis, i.e., the number of basis functions/degrees-of-freedom. ...
See Intrepid2::Basis_HGRAD_LINE_C1_FEM.
Kokkos::View< ordinal_type ***, typename ExecutionSpace::array_layout, Kokkos::HostSpace > OrdinalTypeArray3DHost
View type for 3d host array.
virtual const char * getName() const override
Returns basis name.
Implementation of the default H(grad)-compatible FEM basis of degree 1 on Line cell.
shards::CellTopology basisCellTopology_
Base topology of the cells for which the basis is defined. See the Shards package for definition of b...
Basis_HGRAD_LINE_C1_FEM()
Constructor.
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.
See Intrepid2::Basis_HGRAD_LINE_C1_FEM.