57 #ifndef Intrepid2_DerivedBasis_HDIV_HEX_h
58 #define Intrepid2_DerivedBasis_HDIV_HEX_h
60 #include <Kokkos_DynRankView.hpp>
70 template<
class HGRAD_LINE,
class HVOL_LINE>
76 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
77 using PointViewType =
typename HGRAD_LINE::PointViewType ;
78 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
80 using LineGradBasis = HGRAD_LINE;
81 using LineHVolBasis = HVOL_LINE;
93 TensorBasis3(Teuchos::rcp( new LineGradBasis(polyOrder_x,pointType)),
94 Teuchos::rcp( new LineHVolBasis(polyOrder_y-1,pointType)),
95 Teuchos::rcp( new LineHVolBasis(polyOrder_z-1,pointType)),
98 this->functionSpace_ = FUNCTION_SPACE_HDIV;
105 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
106 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
107 const EOperator DIV = Intrepid2::OPERATOR_DIV;
109 const double weight = 1.0;
110 if (operatorType == VALUE)
112 std::vector< std::vector<EOperator> > ops(3);
113 ops[0] = std::vector<EOperator>{VALUE,VALUE,VALUE};
114 ops[1] = std::vector<EOperator>{};
115 ops[2] = std::vector<EOperator>{};
116 std::vector<double> weights {weight,0.0,0.0};
119 else if (operatorType == DIV)
122 std::vector< std::vector<EOperator> > ops(1);
123 ops[0] = std::vector<EOperator>{GRAD,VALUE,VALUE};
124 std::vector<double> weights {weight};
129 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
143 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
144 const PointViewType inputPoints1,
const PointViewType inputPoints2,
const PointViewType inputPoints3,
145 bool tensorPoints)
const override
147 Intrepid2::EOperator op1, op2, op3;
148 if (operatorType == Intrepid2::OPERATOR_VALUE)
150 op1 = Intrepid2::OPERATOR_VALUE;
151 op2 = Intrepid2::OPERATOR_VALUE;
152 op3 = Intrepid2::OPERATOR_VALUE;
155 auto outputValuesComponent1 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
156 auto outputValuesComponent23 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),std::make_pair(1,3));
161 inputPoints3, op3, tensorPoints);
163 Kokkos::deep_copy(outputValuesComponent23,0.0);
165 else if (operatorType == Intrepid2::OPERATOR_DIV)
170 op1 = Intrepid2::OPERATOR_GRAD;
171 op2 = Intrepid2::OPERATOR_VALUE;
172 op3 = Intrepid2::OPERATOR_VALUE;
178 inputPoints3, op3, tensorPoints, weight);
182 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
198 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
199 auto dofCoeffs23 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),std::make_pair(1,3));
201 Kokkos::deep_copy(dofCoeffs23,0.0);
205 template<
class HGRAD_LINE,
class HVOL_LINE>
211 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
212 using PointViewType =
typename HGRAD_LINE::PointViewType ;
213 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
215 using LineGradBasis = HGRAD_LINE;
216 using LineHVolBasis = HVOL_LINE;
228 TensorBasis3(Teuchos::rcp( new LineHVolBasis(polyOrder_x-1,pointType) ),
229 Teuchos::rcp( new LineGradBasis(polyOrder_y,pointType) ),
230 Teuchos::rcp( new LineHVolBasis(polyOrder_z-1,pointType) ),
233 this->functionSpace_ = FUNCTION_SPACE_HDIV;
240 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
241 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
242 const EOperator DIV = Intrepid2::OPERATOR_DIV;
244 const double weight = 1.0;
245 if (operatorType == VALUE)
247 std::vector< std::vector<EOperator> > ops(3);
248 ops[0] = std::vector<EOperator>{};
249 ops[1] = std::vector<EOperator>{VALUE,VALUE,VALUE};
250 ops[2] = std::vector<EOperator>{};
251 std::vector<double> weights {0.0,weight,0.0};
254 else if (operatorType == DIV)
257 std::vector< std::vector<EOperator> > ops(1);
258 ops[0] = std::vector<EOperator>{VALUE,GRAD,VALUE};
259 std::vector<double> weights {weight};
264 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
278 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
279 const PointViewType inputPoints1,
const PointViewType inputPoints2,
const PointViewType inputPoints3,
280 bool tensorPoints)
const override
282 Intrepid2::EOperator op1, op2, op3;
283 if (operatorType == Intrepid2::OPERATOR_VALUE)
285 op1 = Intrepid2::OPERATOR_VALUE;
286 op2 = Intrepid2::OPERATOR_VALUE;
287 op3 = Intrepid2::OPERATOR_VALUE;
290 auto outputValuesComponent_x = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
291 auto outputValuesComponent_y = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
292 auto outputValuesComponent_z = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
295 Kokkos::deep_copy(outputValuesComponent_x,0.0);
301 inputPoints3, op3, tensorPoints, weight);
304 Kokkos::deep_copy(outputValuesComponent_z,0.0);
306 else if (operatorType == Intrepid2::OPERATOR_DIV)
309 op1 = Intrepid2::OPERATOR_VALUE;
310 op2 = Intrepid2::OPERATOR_GRAD;
311 op3 = Intrepid2::OPERATOR_VALUE;
317 inputPoints3, op3, tensorPoints, weight);
321 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
337 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
338 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),1);
339 auto dofCoeffs3 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),2);
340 Kokkos::deep_copy(dofCoeffs1,0.0);
342 Kokkos::deep_copy(dofCoeffs3,0.0);
346 template<
class HGRAD_LINE,
class HVOL_LINE>
351 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
352 using PointViewType =
typename HGRAD_LINE::PointViewType ;
353 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
355 using LineGradBasis = HGRAD_LINE;
356 using LineHVolBasis = HVOL_LINE;
368 TensorBasis3(Teuchos::rcp( new LineHVolBasis(polyOrder_x-1,pointType) ),
369 Teuchos::rcp( new LineHVolBasis(polyOrder_y-1,pointType) ),
370 Teuchos::rcp( new LineGradBasis(polyOrder_z,pointType) ),
373 this->functionSpace_ = FUNCTION_SPACE_HDIV;
380 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
381 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
382 const EOperator DIV = Intrepid2::OPERATOR_DIV;
384 const double weight = 1.0;
385 if (operatorType == VALUE)
387 std::vector< std::vector<EOperator> > ops(3);
388 ops[0] = std::vector<EOperator>{};
389 ops[1] = std::vector<EOperator>{};
390 ops[2] = std::vector<EOperator>{VALUE,VALUE,VALUE};
391 std::vector<double> weights {0.0,0.0,weight};
394 else if (operatorType == DIV)
397 std::vector< std::vector<EOperator> > ops(1);
398 ops[0] = std::vector<EOperator>{VALUE,VALUE,GRAD};
399 std::vector<double> weights {weight};
404 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
418 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
419 const PointViewType inputPoints1,
const PointViewType inputPoints2,
const PointViewType inputPoints3,
420 bool tensorPoints)
const override
422 Intrepid2::EOperator op1, op2, op3;
423 if (operatorType == Intrepid2::OPERATOR_VALUE)
425 op1 = Intrepid2::OPERATOR_VALUE;
426 op2 = Intrepid2::OPERATOR_VALUE;
427 op3 = Intrepid2::OPERATOR_VALUE;
430 auto outputValuesComponent_xy = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),std::make_pair(0,2));
431 auto outputValuesComponent_z = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
434 Kokkos::deep_copy(outputValuesComponent_xy,0.0);
440 inputPoints3, op3, tensorPoints);
442 else if (operatorType == Intrepid2::OPERATOR_DIV)
447 op1 = Intrepid2::OPERATOR_VALUE;
448 op2 = Intrepid2::OPERATOR_VALUE;
449 op3 = Intrepid2::OPERATOR_GRAD;
455 inputPoints3, op3, tensorPoints, weight);
459 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
475 auto dofCoeffs12 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),std::make_pair(0,2));
476 auto dofCoeffs3 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),2);
477 Kokkos::deep_copy(dofCoeffs12,0.0);
486 template<
class HGRAD_LINE,
class HVOL_LINE>
503 Teuchos::rcp( new
Family1(polyOrder_x, polyOrder_y, polyOrder_z, pointType))) {
504 this->functionSpace_ = FUNCTION_SPACE_HDIV;
508 template<
class HGRAD_LINE,
class HVOL_LINE>
517 ordinal_type order_x_;
518 ordinal_type order_y_;
519 ordinal_type order_z_;
520 EPointType pointType_;
523 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
524 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
525 using PointValueType =
typename HGRAD_LINE::PointValueType;
527 using BasisBase =
typename HGRAD_LINE::BasisBase;
538 Teuchos::rcp(new
Family2 (polyOrder_x, polyOrder_y, polyOrder_z, pointType))) {
539 this->functionSpace_ = FUNCTION_SPACE_HDIV;
541 std::ostringstream basisName;
543 name_ = basisName.str();
545 order_x_ = polyOrder_x;
546 order_y_ = polyOrder_y;
547 order_z_ = polyOrder_z;
548 pointType_ = pointType;
569 return name_.c_str();
582 Teuchos::RCP<BasisBase>
587 if(subCellDim == 2) {
590 return Teuchos::rcp(
new QuadBasis(order_x_-1, order_z_-1, pointType_) );
592 return Teuchos::rcp(
new QuadBasis(order_y_-1,order_z_-1, pointType_) );
594 return Teuchos::rcp(
new QuadBasis(order_x_-1, order_z_-1, pointType_) );
596 return Teuchos::rcp(
new QuadBasis(order_z_-1, order_y_-1, pointType_) );
598 return Teuchos::rcp(
new QuadBasis(order_y_-1, order_x_-1, pointType_) );
600 return Teuchos::rcp(
new QuadBasis(order_x_-1, order_y_-1, pointType_) );
604 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
611 virtual HostBasisPtr<OutputValueType, PointValueType>
615 auto hostBasis = Teuchos::rcp(
new HostBasis(order_x_, order_y_, order_z_, pointType_));
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, const PointViewType inputPoints3, bool tensorPoints) const override
multi-component getValues() method (required/called by TensorBasis3)
Implementation of H(vol) basis on the quadrilateral that is templated on H(vol) on the line...
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 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_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
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 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 inputPoints12, const PointViewType inputPoints3, bool tensorPoints) const override
Evaluation of a tensor FEM basis on a reference cell.
virtual const char * getName() const override
Returns basis name.
Implementation of a basis that is the direct sum of two other bases.
A basis that is the direct sum of two other bases.
virtual const char * getName() const override
Returns basis name.
virtual bool requireOrientation() const override
True if orientation is required.
Basis_Derived_HDIV_HEX(int polyOrder, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Basis_Derived_HDIV_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Header file to include all Sacado headers that are required if using Intrepid2 with Sacado types...
Basis_Derived_HDIV_Family3_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
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...
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, const PointViewType inputPoints3, bool tensorPoints) const override
multi-component getValues() method (required/called by TensorBasis3)
virtual void getValues(OutputViewType outputValues, const EOperator operatorType, const PointViewType inputPoints1, const PointViewType inputPoints2, const PointViewType inputPoints3, bool tensorPoints) const override
multi-component getValues() method (required/called by TensorBasis3)
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_Family1_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Basis_Derived_HDIV_Family3_Family1_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType)
Constructor.
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...
Teuchos::RCP< BasisBase > getSubCellRefBasis(const ordinal_type subCellDim, const ordinal_type subCellOrd) const override
returns the basis associated to a subCell.