53 #ifndef Intrepid2_DerivedBasis_HIV_QUAD_h
54 #define Intrepid2_DerivedBasis_HIV_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 LineGradBasis = HGRAD_LINE;
80 using LineHVolBasis = HVOL_LINE;
82 using BasisBase =
typename HGRAD_LINE::BasisBase;
93 TensorBasis(Teuchos::rcp(new LineHVolBasis(polyOrder_x-1,pointType)),
94 Teuchos::rcp(new LineGradBasis(polyOrder_y,pointType)))
96 this->functionSpace_ = FUNCTION_SPACE_HDIV;
97 this->setShardsTopologyAndTags();
104 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
105 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
106 const EOperator DIV = Intrepid2::OPERATOR_DIV;
110 const double weight = -1.0;
111 if (operatorType == VALUE)
113 std::vector< std::vector<EOperator> > ops(2);
114 ops[0] = std::vector<EOperator>{};
115 ops[1] = std::vector<EOperator>{VALUE,VALUE};
116 std::vector<double> weights {0.0,weight};
119 else if (operatorType == DIV)
122 std::vector< std::vector<EOperator> > ops(1);
123 ops[0] = std::vector<EOperator>{VALUE,GRAD};
124 std::vector<double> weights {weight};
129 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
142 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
143 const PointViewType inputPoints1,
const PointViewType inputPoints2,
144 bool tensorPoints)
const override
148 const double weight = -1.0;
150 Intrepid2::EOperator op1, op2;
151 if (operatorType == Intrepid2::OPERATOR_VALUE)
153 op1 = Intrepid2::OPERATOR_VALUE;
154 op2 = Intrepid2::OPERATOR_VALUE;
157 auto outputValuesComponent1 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
158 auto outputValuesComponent2 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
162 inputPoints2, op2, tensorPoints, weight);
164 Kokkos::deep_copy(outputValuesComponent1,0.0);
166 else if (operatorType == Intrepid2::OPERATOR_DIV)
170 op1 = Intrepid2::OPERATOR_VALUE;
171 op2 = Intrepid2::OPERATOR_GRAD;
175 inputPoints2, op2, tensorPoints, weight);
179 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
195 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
196 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),1);
197 Kokkos::deep_copy(dofCoeffs1,0.0);
200 Kokkos::parallel_for( Kokkos::RangePolicy<ExecutionSpace>(0, dofCoeffs2.extent(0)),
201 KOKKOS_LAMBDA (
const int i){ dofCoeffs2(i) *= -1.0; });
206 template<
class HGRAD_LINE,
class HVOL_LINE>
210 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
211 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
212 using PointValueType =
typename HGRAD_LINE::PointValueType;
214 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
215 using PointViewType =
typename HGRAD_LINE::PointViewType ;
216 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
218 using LineGradBasis = HGRAD_LINE;
219 using LineHVolBasis = HVOL_LINE;
221 using BasisBase =
typename HGRAD_LINE::BasisBase;
232 TensorBasis(Teuchos::rcp(new LineGradBasis(polyOrder_x,pointType)),
233 Teuchos::rcp(new LineHVolBasis(polyOrder_y-1,pointType)))
235 this->functionSpace_ = FUNCTION_SPACE_HDIV;
236 this->setShardsTopologyAndTags();
243 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
244 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
245 const EOperator DIV = Intrepid2::OPERATOR_DIV;
247 const double weight = 1.0;
248 if (operatorType == VALUE)
250 std::vector< std::vector<EOperator> > ops(2);
251 ops[0] = std::vector<EOperator>{VALUE,VALUE};
252 ops[1] = std::vector<EOperator>{};
253 std::vector<double> weights {weight, 0.0};
256 else if (operatorType == DIV)
259 std::vector< std::vector<EOperator> > ops(1);
260 ops[0] = std::vector<EOperator>{GRAD,VALUE};
261 std::vector<double> weights {weight};
266 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
279 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
280 const PointViewType inputPoints1,
const PointViewType inputPoints2,
281 bool tensorPoints)
const override
283 Intrepid2::EOperator op1, op2;
284 if (operatorType == Intrepid2::OPERATOR_VALUE)
286 op1 = Intrepid2::OPERATOR_VALUE;
287 op2 = Intrepid2::OPERATOR_VALUE;
290 auto outputValuesComponent1 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
291 auto outputValuesComponent2 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
295 inputPoints2, op2, tensorPoints);
297 Kokkos::deep_copy(outputValuesComponent2, 0.0);
299 else if (operatorType == Intrepid2::OPERATOR_DIV)
303 op1 = Intrepid2::OPERATOR_GRAD;
304 op2 = Intrepid2::OPERATOR_VALUE;
308 inputPoints2, op2, tensorPoints);
312 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
328 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
329 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),1);
331 Kokkos::deep_copy(dofCoeffs2,0.0);
336 template<
class HGRAD_LINE,
class HVOL_LINE>
345 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
346 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
347 using PointValueType =
typename HGRAD_LINE::PointValueType;
349 using BasisBase =
typename HGRAD_LINE::BasisBase;
353 ordinal_type order_x_;
354 ordinal_type order_y_;
355 EPointType pointType_;
365 Teuchos::rcp(new
Family2(polyOrder_x, polyOrder_y, pointType)))
367 this->functionSpace_ = FUNCTION_SPACE_HDIV;
369 std::ostringstream basisName;
371 name_ = basisName.str();
373 order_x_ = polyOrder_x;
374 order_y_ = polyOrder_y;
375 pointType_ = pointType;
397 return name_.c_str();
410 Teuchos::RCP<BasisBase>
412 if(subCellDim == 1) {
416 return Teuchos::rcp(
new HVOL_LINE(order_x_-1, pointType_) );
419 return Teuchos::rcp(
new HVOL_LINE(order_y_-1, pointType_) );
423 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
430 virtual HostBasisPtr<OutputValueType, PointValueType>
434 auto hostBasis = Teuchos::rcp(
new HostBasis(order_x_, order_y_, pointType_));
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...
Implementation of bases that are tensor products of two or three component bases. ...
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_Derived_HDIV_Family1_QUAD(int polyOrder_x, int polyOrder_y, const EPointType pointType)
Constructor.
Teuchos::RCP< BasisBase > getSubCellRefBasis(const ordinal_type subCellDim, const ordinal_type subCellOrd) const override
returns the basis associated to a subCell.
Free functions, callable from device code, that implement various polynomials useful in basis definit...
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...
Implementation of a basis that is the direct sum of two other bases.
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...
Basis_Derived_HDIV_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.
Header file to include all Sacado headers that are required if using Intrepid2 with Sacado types...
virtual bool requireOrientation() const override
True if orientation is required.
virtual const char * getName() const override
Returns basis name.
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 void getDofCoeffs(typename BasisBase::ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom on the reference cell.
Basis_Derived_HDIV_Family2_QUAD(int polyOrder_x, int polyOrder_y, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
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 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_Derived_HDIV_QUAD(int polyOrder, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Basis defined as the tensor product of two component bases.
virtual HostBasisPtr< OutputValueType, PointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...