53 #ifndef Intrepid2_DerivedBasis_HCURL_QUAD_h
54 #define Intrepid2_DerivedBasis_HCURL_QUAD_h
56 #include <Kokkos_DynRankView.hpp>
66 template<
class HGRAD_LINE,
class HVOL_LINE>
71 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
72 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
73 using PointValueType =
typename HGRAD_LINE::PointValueType;
75 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
76 using PointViewType =
typename HGRAD_LINE::PointViewType ;
77 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
79 using BasisBase =
typename HGRAD_LINE::BasisBase;
81 using LineGradBasis = HGRAD_LINE;
82 using LineHVolBasis = HVOL_LINE;
93 TensorBasis(Teuchos::rcp( new LineHVolBasis(polyOrder_x-1,pointType)),
94 Teuchos::rcp( new LineGradBasis(polyOrder_y,pointType)))
96 this->functionSpace_ = FUNCTION_SPACE_HCURL;
97 this->setShardsTopologyAndTags();
104 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
105 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
106 const EOperator CURL = Intrepid2::OPERATOR_CURL;
107 if (operatorType == VALUE)
110 std::vector< std::vector<EOperator> > ops(2);
111 ops[0] = std::vector<EOperator>{VALUE,VALUE};
112 ops[1] = std::vector<EOperator>{};
113 std::vector<double> weights {1.0, 0.0};
116 else if (operatorType == CURL)
124 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
137 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
138 const PointViewType inputPoints1,
const PointViewType inputPoints2,
139 bool tensorPoints)
const override
141 Intrepid2::EOperator op1, op2;
142 if (operatorType == Intrepid2::OPERATOR_VALUE)
144 op1 = Intrepid2::OPERATOR_VALUE;
145 op2 = Intrepid2::OPERATOR_VALUE;
148 OutputViewType outputValuesComponent1 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
149 OutputViewType outputValuesComponent2 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
153 inputPoints2, op2, tensorPoints);
155 Kokkos::deep_copy(outputValuesComponent2,0);
157 else if (operatorType == Intrepid2::OPERATOR_CURL)
161 op1 = Intrepid2::OPERATOR_VALUE;
162 op2 = Intrepid2::OPERATOR_GRAD;
164 double weight = -1.0;
167 inputPoints2, op2, tensorPoints, weight);
171 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
187 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
188 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),1);
190 Kokkos::deep_copy(dofCoeffs2,0.0);
194 template<
class HGRAD_LINE,
class HVOL_LINE>
200 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
201 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
202 using PointValueType =
typename HGRAD_LINE::PointValueType;
204 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
205 using PointViewType =
typename HGRAD_LINE::PointViewType ;
206 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
208 using LineGradBasis = HGRAD_LINE;
209 using LineHVolBasis = HVOL_LINE;
211 using BasisBase =
typename HGRAD_LINE::BasisBase;
222 TensorBasis(Teuchos::rcp( new LineGradBasis(polyOrder_x,pointType) ),
223 Teuchos::rcp( new LineHVolBasis(polyOrder_y-1,pointType) ))
225 this->functionSpace_ = FUNCTION_SPACE_HCURL;
226 this->setShardsTopologyAndTags();
233 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
234 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
235 const EOperator CURL = Intrepid2::OPERATOR_CURL;
236 if (operatorType == VALUE)
239 std::vector< std::vector<EOperator> > ops(2);
240 ops[0] = std::vector<EOperator>{};
241 ops[1] = std::vector<EOperator>{VALUE,VALUE};
242 std::vector<double> weights {0.0, 1.0};
245 else if (operatorType == CURL)
253 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
266 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
267 const PointViewType inputPoints1,
const PointViewType inputPoints2,
268 bool tensorPoints)
const override
270 Intrepid2::EOperator op1, op2;
271 if (operatorType == Intrepid2::OPERATOR_VALUE)
273 op1 = Intrepid2::OPERATOR_VALUE;
274 op2 = Intrepid2::OPERATOR_VALUE;
277 auto outputValuesComponent1 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
278 auto outputValuesComponent2 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
281 Kokkos::deep_copy(outputValuesComponent1, 0.0);
284 inputPoints2, op2, tensorPoints);
287 else if (operatorType == Intrepid2::OPERATOR_CURL)
291 op1 = Intrepid2::OPERATOR_GRAD;
292 op2 = Intrepid2::OPERATOR_VALUE;
296 inputPoints2, op2, tensorPoints);
300 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
316 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
317 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),1);
318 Kokkos::deep_copy(dofCoeffs1,0.0);
323 template<
class HGRAD_LINE,
class HVOL_LINE>
331 using BasisBase =
typename HGRAD_LINE::BasisBase;
335 ordinal_type order_x_;
336 ordinal_type order_y_;
337 EPointType pointType_;
340 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
341 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
342 using PointValueType =
typename HGRAD_LINE::PointValueType;
352 Teuchos::rcp( new
Family2(polyOrder_x, polyOrder_y, pointType) ))
354 this->functionSpace_ = FUNCTION_SPACE_HCURL;
356 std::ostringstream basisName;
358 name_ = basisName.str();
360 order_x_ = polyOrder_x;
361 order_y_ = polyOrder_y;
362 pointType_ = pointType;
385 return name_.c_str();
398 Teuchos::RCP<BasisBase>
400 if(subCellDim == 1) {
404 return Teuchos::rcp(
new HVOL_LINE(order_x_-1, pointType_) );
407 return Teuchos::rcp(
new HVOL_LINE(order_y_-1, pointType_) );
411 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
418 virtual HostBasisPtr<OutputValueType, PointValueType>
422 auto hostBasis = Teuchos::rcp(
new HostBasis(order_x_, order_y_, pointType_));
virtual OperatorTensorDecomposition getSimpleOperatorDecomposition(const EOperator &operatorType) const override
Returns a simple decomposition of the specified operator: what operator(s) should be applied to basis...
Basis_Derived_HCURL_QUAD(int polyOrder_x, int polyOrder_y, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Implementation of bases that are tensor products of two or three component bases. ...
Teuchos::RCP< BasisBase > getSubCellRefBasis(const ordinal_type subCellDim, const ordinal_type subCellOrd) const override
returns the basis associated to a subCell.
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...
Free functions, callable from device code, that implement various polynomials useful in basis definit...
virtual bool requireOrientation() const override
True if orientation is required.
For a multi-component tensor basis, specifies the operators to be applied to the components to produc...
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...
virtual OperatorTensorDecomposition getSimpleOperatorDecomposition(const EOperator &operatorType) const override
Returns a simple decomposition of the specified operator: what operator(s) should be applied to basis...
Implementation of a basis that is the direct sum of two other bases.
Basis_Derived_HCURL_Family1_QUAD(int polyOrder_x, int polyOrder_y, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
A basis that is the direct sum of two other bases.
virtual const char * getName() const override
Returns basis name.
virtual HostBasisPtr< OutputValueType, PointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
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)
Header file to include all Sacado headers that are required if using Intrepid2 with Sacado types...
virtual const char * getName() const override
Returns basis name.
virtual void getDofCoeffs(typename BasisBase::ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom on the reference cell.
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)
Basis defined as the tensor product of two component bases.
Basis_Derived_HCURL_Family2_QUAD(int polyOrder_x, int polyOrder_y, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Basis_Derived_HCURL_QUAD(int polyOrder, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.