64 #ifndef Intrepid2_DerivedBasis_HDIV_WEDGE_h
65 #define Intrepid2_DerivedBasis_HDIV_WEDGE_h
67 #include <Kokkos_DynRankView.hpp>
77 template<
class HDIV_TRI,
class HVOL_LINE>
83 using OutputViewType =
typename HVOL_LINE::OutputViewType;
84 using PointViewType =
typename HVOL_LINE::PointViewType ;
85 using ScalarViewType =
typename HVOL_LINE::ScalarViewType;
87 using TriDivBasis = HDIV_TRI;
88 using LineHVolBasis = HVOL_LINE;
90 using BasisBase =
typename HVOL_LINE::BasisBase;
93 using DeviceType =
typename BasisBase::DeviceType;
94 using ExecutionSpace =
typename BasisBase::ExecutionSpace;
95 using OutputValueType =
typename BasisBase::OutputValueType;
96 using PointValueType =
typename BasisBase::PointValueType;
105 TensorBasis(Teuchos::rcp( new TriDivBasis(polyOrder_xy, pointType) ),
106 Teuchos::rcp( new LineHVolBasis(polyOrder_z-1, pointType) ))
108 this->functionSpace_ = FUNCTION_SPACE_HDIV;
109 this->setShardsTopologyAndTags();
116 if (operatorType == OPERATOR_DIV)
120 std::vector< std::vector<EOperator> > ops(1);
121 ops[0] = std::vector<EOperator>{OPERATOR_DIV,OPERATOR_VALUE};
122 std::vector<double> weights {1.0};
124 return opDecomposition;
126 else if (OPERATOR_VALUE == operatorType)
129 std::vector< std::vector<EOperator> > ops(2);
130 ops[0] = std::vector<EOperator>{OPERATOR_VALUE,OPERATOR_VALUE};
131 ops[1] = std::vector<EOperator>{};
132 std::vector<double> weights {1.0,0.0};
137 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
150 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
151 const PointViewType inputPoints1,
const PointViewType inputPoints2,
152 bool tensorPoints)
const override
155 if (operatorType == OPERATOR_VALUE)
157 op1 = OPERATOR_VALUE;
158 op2 = OPERATOR_VALUE;
161 auto outputValuesComponent12 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),std::pair<int,int>{0,2});
162 auto outputValuesComponent3 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
165 Kokkos::deep_copy(outputValuesComponent3, 0.0);
168 inputPoints2, op2, tensorPoints);
171 else if (operatorType == OPERATOR_DIV)
176 op2 = OPERATOR_VALUE;
180 inputPoints2, op2, tensorPoints);
184 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
200 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(), std::make_pair(0,2));
201 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(), 2);
203 Kokkos::deep_copy(dofCoeffs2,0.0);
207 template<
class HVOL_TRI,
class HGRAD_LINE>
212 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
213 using PointViewType =
typename HGRAD_LINE::PointViewType ;
214 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
216 using BasisBase =
typename HGRAD_LINE::BasisBase;
218 using DeviceType =
typename BasisBase::DeviceType;
219 using ExecutionSpace =
typename BasisBase::ExecutionSpace;
220 using OutputValueType =
typename BasisBase::OutputValueType;
221 using PointValueType =
typename BasisBase::PointValueType;
223 using TriVolBasis = HVOL_TRI;
224 using LineGradBasis = HGRAD_LINE;
235 TensorBasis(Teuchos::rcp( new TriVolBasis(polyOrder_xy-1,pointType)),
236 Teuchos::rcp( new LineGradBasis(polyOrder_z,pointType)))
238 this->functionSpace_ = FUNCTION_SPACE_HDIV;
239 this->setShardsTopologyAndTags();
248 const EOperator & VALUE = OPERATOR_VALUE;
249 const EOperator & GRAD = OPERATOR_GRAD;
250 const EOperator & DIV = OPERATOR_DIV;
251 if (operatorType == VALUE)
254 std::vector< std::vector<EOperator> > ops(2);
255 ops[0] = std::vector<EOperator>{};
256 ops[1] = std::vector<EOperator>{VALUE,VALUE};
257 std::vector<double> weights {0.,1.0};
259 return opDecomposition;
261 else if (operatorType == DIV)
264 std::vector< std::vector<EOperator> > ops(1);
265 ops[0] = std::vector<EOperator>{VALUE,GRAD};
266 std::vector<double> weights {1.0};
268 return opDecomposition;
272 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
283 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
284 const PointViewType inputPoints1,
const PointViewType inputPoints2,
285 bool tensorPoints)
const override
288 if (operatorType == OPERATOR_VALUE)
290 op1 = OPERATOR_VALUE;
291 op2 = OPERATOR_VALUE;
294 auto outputValuesComponent12 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),std::pair<int,int>{0,2});
295 auto outputValuesComponent3 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
298 Kokkos::deep_copy(outputValuesComponent12, 0.0);
301 inputPoints2, op2, tensorPoints);
304 else if (operatorType == OPERATOR_DIV)
307 op1 = OPERATOR_VALUE;
312 inputPoints2, op2, tensorPoints);
316 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
332 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(), std::make_pair(0,2));
333 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(), 2);
334 Kokkos::deep_copy(dofCoeffs1,0.0);
340 template<
class HDIV_TRI,
class HVOL_TRI,
class HGRAD_LINE,
class HVOL_LINE>
348 using BasisBase =
typename HGRAD_LINE::BasisBase;
352 ordinal_type order_xy_;
353 ordinal_type order_z_;
354 EPointType pointType_;
357 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
358 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
359 using PointValueType =
typename HGRAD_LINE::PointValueType;
369 Teuchos::rcp( new
Family2(polyOrder_xy, polyOrder_z, pointType) ))
371 this->functionSpace_ = FUNCTION_SPACE_HDIV;
373 std::ostringstream basisName;
375 name_ = basisName.str();
377 order_xy_ = polyOrder_xy;
378 order_z_ = polyOrder_z;
379 pointType_ = pointType;
402 return name_.c_str();
415 Teuchos::RCP<BasisBase>
419 using TriBasis = HVOL_TRI;
421 if(subCellDim == 2) {
424 return Teuchos::rcp(
new QuadBasis(order_xy_-1, order_z_-1, pointType_) );
426 return Teuchos::rcp(
new QuadBasis(order_xy_-1, order_z_-1, pointType_) );
428 return Teuchos::rcp(
new QuadBasis(order_z_-1, order_xy_-1, pointType_) );
430 return Teuchos::rcp(
new TriBasis(order_xy_-1, pointType_) );
432 return Teuchos::rcp(
new TriBasis(order_xy_-1, pointType_) );
434 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"subCellOrd is out of bounds");
437 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"subCellDim is out of bounds");
444 virtual HostBasisPtr<OutputValueType, PointValueType>
448 auto hostBasis = Teuchos::rcp(
new HostBasis(order_xy_, order_z_, pointType_));
virtual HostBasisPtr< OutputValueType, PointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
virtual bool requireOrientation() const override
True if orientation is required.
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)
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...
Implementation of H(vol) basis on the quadrilateral that is templated on H(vol) on the line...
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...
Basis_Derived_HDIV_Family1_WEDGE(int polyOrder_xy, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
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.
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_HDIV_WEDGE(int polyOrder_xy, int polyOrder_z, 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 void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...
Header file to include all Sacado headers that are required if using Intrepid2 with Sacado types...
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_HDIV_Family2_WEDGE(int polyOrder_xy, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
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.
Basis_Derived_HDIV_WEDGE(int polyOrder, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Basis defined as the tensor product of two 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)