51 #ifndef Intrepid2_DerivedBasis_HGRAD_QUAD_h 52 #define Intrepid2_DerivedBasis_HGRAD_QUAD_h 58 template<
class HGRAD_LINE>
64 ordinal_type order_x_, order_y_;
67 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
68 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
69 using PointValueType =
typename HGRAD_LINE::PointValueType;
71 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
72 using PointViewType =
typename HGRAD_LINE::PointViewType ;
73 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
75 using BasisBase =
typename HGRAD_LINE::BasisBase;
76 using LineBasis = HGRAD_LINE;
86 TensorBasis(Teuchos::rcp( new LineBasis(polyOrder_x, pointType)),
87 Teuchos::rcp( new LineBasis(polyOrder_y, pointType)))
89 this->functionSpace_ = FUNCTION_SPACE_HGRAD;
91 std::ostringstream basisName;
93 name_ = basisName.str();
95 order_x_= polyOrder_x;
96 order_y_ = polyOrder_x;
97 pointType_ = pointType;
110 return (this->getDofCount(1,0) > 1);
117 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
118 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
120 if (operatorType == VALUE)
124 else if (operatorType == GRAD)
129 std::vector< std::vector<EOperator> > ops;
130 ops.push_back(std::vector<EOperator>{GRAD, VALUE});
131 ops.push_back(std::vector<EOperator>{VALUE, GRAD});
133 std::vector<double> weights(ops.size(), 1.0);
139 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
143 using BasisBase::getValues;
153 const PointViewType inputPoints1,
const PointViewType inputPoints2,
154 bool tensorPoints)
const override 157 if (operatorType == Intrepid2::OPERATOR_VALUE)
159 op1 = Intrepid2::OPERATOR_VALUE;
160 op2 = Intrepid2::OPERATOR_VALUE;
164 inputPoints2, op2, tensorPoints);
166 else if (operatorType == Intrepid2::OPERATOR_GRAD)
175 auto outputValuesComponent1 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
176 auto outputValuesComponent2 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
179 op1 = Intrepid2::OPERATOR_GRAD;
180 op2 = Intrepid2::OPERATOR_VALUE;
184 inputPoints2, op2, tensorPoints);
187 op1 = Intrepid2::OPERATOR_VALUE;
188 op2 = Intrepid2::OPERATOR_GRAD;
192 inputPoints2, op2, tensorPoints);
196 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
208 return name_.c_str();
220 Teuchos::RCP<BasisBase>
222 if(subCellDim == 1) {
226 return Teuchos::rcp(
new LineBasis(order_x_, pointType_) );
229 return Teuchos::rcp(
new LineBasis(order_y_, pointType_) );
233 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
244 auto hostBasis = Teuchos::rcp(
new HostBasis(order_x_, order_y_, pointType_));
Implementation of bases that are tensor products of two or three component bases. ...
virtual HostBasisPtr< OutputValueType, PointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
Basis_Derived_HGRAD_QUAD(int polyOrder_x, int polyOrder_y, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
virtual void getValues(OutputViewType outputValues, const EOperator operatorType, const PointViewType inputPoints1, const PointViewType inputPoints2, bool tensorPoints) const override
multi-component getValues() method (required/called by TensorBasis)
For a multi-component tensor basis, specifies the operators to be applied to the components to produc...
BasisPtr< typename Kokkos::HostSpace::device_type, OutputType, PointType > HostBasisPtr
Pointer to a Basis whose device type is on the host (Kokkos::HostSpace::device_type), allowing host access to input and output views, and ensuring host execution of basis evaluation.
virtual const char * getName() const override
Returns basis name.
EOperator
Enumeration of primitive operators available in Intrepid. Primitive operators act on reconstructed fu...
virtual OperatorTensorDecomposition getSimpleOperatorDecomposition(const EOperator operatorType) const override
Returns a simple decomposition of the specified operator: what operator(s) should be applied to basis...
virtual const char * getName() const override
Returns basis name.
Basis_Derived_HGRAD_QUAD(int polyOrder, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
EPointType
Enumeration of types of point distributions in Intrepid.
virtual bool requireOrientation() const override
True if orientation is required.
virtual void getValues(BasisValues< OutputValueType, DeviceType > outputValues, const TensorPoints< PointValueType, DeviceType > inputPoints, const EOperator operatorType=OPERATOR_VALUE) const override
Evaluation of a FEM basis on a reference cell, using point and output value containers that allow pre...
Basis defined as the tensor product of two component bases.
Teuchos::RCP< BasisBase > getSubCellRefBasis(const ordinal_type subCellDim, const ordinal_type subCellOrd) const override
returns the basis associated to a subCell.