54 #ifndef Intrepid2_DerivedBasis_HCURL_HEX_h
55 #define Intrepid2_DerivedBasis_HCURL_HEX_h
57 #include <Kokkos_DynRankView.hpp>
67 template<
class HGRAD_LINE,
class HVOL_LINE>
72 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
73 using PointViewType =
typename HGRAD_LINE::PointViewType ;
74 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
76 using LineGradBasis = HGRAD_LINE;
77 using LineVolBasis = HVOL_LINE;
89 TensorBasis3(Teuchos::rcp(new LineVolBasis (polyOrder_x-1,pointType)),
90 Teuchos::rcp(new LineGradBasis(polyOrder_y,pointType)),
91 Teuchos::rcp(new LineGradBasis(polyOrder_z,pointType)),
94 this->functionSpace_ = FUNCTION_SPACE_HCURL;
101 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
102 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
104 if (operatorType == Intrepid2::OPERATOR_VALUE)
107 std::vector< std::vector<EOperator> > ops(3);
108 ops[0] = std::vector<EOperator>{VALUE,VALUE,VALUE};
109 ops[1] = std::vector<EOperator>{};
110 ops[2] = std::vector<EOperator>{};
111 std::vector<double> weights {1.0, 0.0, 0.0};
114 else if (operatorType == Intrepid2::OPERATOR_CURL)
121 std::vector< std::vector<EOperator> > ops(3);
122 ops[0] = std::vector<EOperator>{};
123 ops[1] = std::vector<EOperator>{VALUE,VALUE,GRAD};
124 ops[2] = std::vector<EOperator>{VALUE,GRAD,VALUE};
126 std::vector<double> weights {0.0, 1.0, -1.0};
131 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
145 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
146 const PointViewType inputPoints1,
const PointViewType inputPoints2,
const PointViewType inputPoints3,
147 bool tensorPoints)
const override
149 Intrepid2::EOperator op1, op2, op3;
150 if (operatorType == Intrepid2::OPERATOR_VALUE)
152 op1 = Intrepid2::OPERATOR_VALUE;
153 op2 = Intrepid2::OPERATOR_VALUE;
154 op3 = Intrepid2::OPERATOR_VALUE;
157 auto outputValuesComponent1 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
158 auto outputValuesComponent23 = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),std::make_pair(1,3));
163 inputPoints3, op3, tensorPoints);
165 Kokkos::deep_copy(outputValuesComponent23,0.0);
167 else if (operatorType == Intrepid2::OPERATOR_CURL)
170 auto outputValuesComponent_x = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
171 auto outputValuesComponent_y = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
172 auto outputValuesComponent_z = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
175 Kokkos::deep_copy(outputValuesComponent_x, 0.0);
178 op1 = Intrepid2::OPERATOR_VALUE;
179 op2 = Intrepid2::OPERATOR_VALUE;
180 op3 = Intrepid2::OPERATOR_GRAD;
186 inputPoints3, op3, tensorPoints, weight);
189 op1 = Intrepid2::OPERATOR_VALUE;
190 op2 = Intrepid2::OPERATOR_GRAD;
191 op3 = Intrepid2::OPERATOR_VALUE;
196 inputPoints3, op3, tensorPoints, weight);
200 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
216 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
217 auto dofCoeffs23 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),std::make_pair(1,3));
219 Kokkos::deep_copy(dofCoeffs23,0.0);
223 template<
class HGRAD_LINE,
class HVOL_LINE>
228 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
229 using PointViewType =
typename HGRAD_LINE::PointViewType ;
230 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
232 using LineGradBasis = HGRAD_LINE;
233 using LineVolBasis = HVOL_LINE;
245 TensorBasis3(Teuchos::rcp( new LineGradBasis(polyOrder_x,pointType)),
246 Teuchos::rcp( new LineVolBasis (polyOrder_y-1,pointType)),
247 Teuchos::rcp( new LineGradBasis(polyOrder_z,pointType)),
250 this->functionSpace_ = FUNCTION_SPACE_HCURL;
257 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
258 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
259 const EOperator CURL = Intrepid2::OPERATOR_CURL;
260 if (operatorType == VALUE)
263 std::vector< std::vector<EOperator> > ops(3);
264 ops[0] = std::vector<EOperator>{};
265 ops[1] = std::vector<EOperator>{VALUE,VALUE,VALUE};
266 ops[2] = std::vector<EOperator>{};
267 std::vector<double> weights {0.0, 1.0, 0.0};
270 else if (operatorType == CURL)
276 std::vector< std::vector<EOperator> > ops(3);
277 ops[0] = std::vector<EOperator>{VALUE,VALUE,GRAD};
278 ops[1] = std::vector<EOperator>{};
279 ops[2] = std::vector<EOperator>{GRAD,VALUE,VALUE};
281 std::vector<double> weights {-1.0, 0.0, 1.0};
286 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
300 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
301 const PointViewType inputPoints1,
const PointViewType inputPoints2,
const PointViewType inputPoints3,
302 bool tensorPoints)
const override
304 Intrepid2::EOperator op1, op2, op3;
305 if (operatorType == Intrepid2::OPERATOR_VALUE)
307 op1 = Intrepid2::OPERATOR_VALUE;
308 op2 = Intrepid2::OPERATOR_VALUE;
309 op3 = Intrepid2::OPERATOR_VALUE;
312 auto outputValuesComponent_x = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
313 auto outputValuesComponent_y = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
314 auto outputValuesComponent_z = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
317 Kokkos::deep_copy(outputValuesComponent_x,0.0);
322 inputPoints3, op3, tensorPoints);
324 Kokkos::deep_copy(outputValuesComponent_z,0.0);
326 else if (operatorType == Intrepid2::OPERATOR_CURL)
329 auto outputValuesComponent_x = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
330 auto outputValuesComponent_y = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
331 auto outputValuesComponent_z = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
334 op1 = Intrepid2::OPERATOR_VALUE;
335 op2 = Intrepid2::OPERATOR_VALUE;
336 op3 = Intrepid2::OPERATOR_GRAD;
338 double weight = -1.0;
342 inputPoints3, op3, tensorPoints, weight);
345 Kokkos::deep_copy(outputValuesComponent_y, 0.0);
348 op1 = Intrepid2::OPERATOR_GRAD;
349 op2 = Intrepid2::OPERATOR_VALUE;
350 op3 = Intrepid2::OPERATOR_VALUE;
355 inputPoints3, op3, tensorPoints, weight);
359 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
375 auto dofCoeffs1 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),0);
376 auto dofCoeffs2 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),1);
377 auto dofCoeffs3 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),2);
378 Kokkos::deep_copy(dofCoeffs1,0.0);
380 Kokkos::deep_copy(dofCoeffs3,0.0);
386 template<
class HGRAD_LINE,
class HVOL_LINE>
390 using OutputViewType =
typename HGRAD_LINE::OutputViewType;
391 using PointViewType =
typename HGRAD_LINE::PointViewType ;
392 using ScalarViewType =
typename HGRAD_LINE::ScalarViewType;
394 using LineGradBasis = HGRAD_LINE;
395 using LineVolBasis = HVOL_LINE;
407 TensorBasis3(Teuchos::rcp(new LineGradBasis(polyOrder_x,pointType)),
408 Teuchos::rcp(new LineGradBasis(polyOrder_y,pointType)),
409 Teuchos::rcp(new LineVolBasis (polyOrder_z-1,pointType)),
412 this->functionSpace_ = FUNCTION_SPACE_HCURL;
419 const EOperator VALUE = Intrepid2::OPERATOR_VALUE;
420 const EOperator GRAD = Intrepid2::OPERATOR_GRAD;
421 if (operatorType == Intrepid2::OPERATOR_VALUE)
424 std::vector< std::vector<EOperator> > ops(3);
425 ops[0] = std::vector<EOperator>{};
426 ops[1] = std::vector<EOperator>{};
427 ops[2] = std::vector<EOperator>{VALUE,VALUE,VALUE};
428 std::vector<double> weights {0.0, 0.0, 1.0};
431 else if (operatorType == Intrepid2::OPERATOR_CURL)
437 std::vector< std::vector<EOperator> > ops(3);
438 ops[0] = std::vector<EOperator>{VALUE,GRAD,VALUE};
439 ops[1] = std::vector<EOperator>{GRAD,VALUE,VALUE};
440 ops[2] = std::vector<EOperator>{};
442 std::vector<double> weights {1.0, -1.0, 0.0};
447 INTREPID2_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
"Unsupported operator type");
461 virtual void getValues(OutputViewType outputValues,
const EOperator operatorType,
462 const PointViewType inputPoints1,
const PointViewType inputPoints2,
const PointViewType inputPoints3,
463 bool tensorPoints)
const override
465 Intrepid2::EOperator op1, op2, op3;
466 if (operatorType == Intrepid2::OPERATOR_VALUE)
468 op1 = Intrepid2::OPERATOR_VALUE;
469 op2 = Intrepid2::OPERATOR_VALUE;
470 op3 = Intrepid2::OPERATOR_VALUE;
473 auto outputValuesComponent_xy = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),std::make_pair(0,2));
474 auto outputValuesComponent_z = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
477 Kokkos::deep_copy(outputValuesComponent_xy,0.0);
482 inputPoints3, op3, tensorPoints);
484 else if (operatorType == Intrepid2::OPERATOR_CURL)
487 auto outputValuesComponent_x = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),0);
488 auto outputValuesComponent_y = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),1);
489 auto outputValuesComponent_z = Kokkos::subview(outputValues,Kokkos::ALL(),Kokkos::ALL(),2);
492 op1 = Intrepid2::OPERATOR_VALUE;
493 op2 = Intrepid2::OPERATOR_GRAD;
494 op3 = Intrepid2::OPERATOR_VALUE;
500 inputPoints3, op3, tensorPoints, weight);
502 op1 = Intrepid2::OPERATOR_GRAD;
503 op2 = Intrepid2::OPERATOR_VALUE;
504 op3 = Intrepid2::OPERATOR_VALUE;
509 inputPoints3, op3, tensorPoints, weight);
512 Kokkos::deep_copy(outputValuesComponent_z, 0.0);
516 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"operator not yet supported");
532 auto dofCoeffs12 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),std::make_pair(0,2));
533 auto dofCoeffs3 = Kokkos::subview(dofCoeffs,Kokkos::ALL(),2);
534 Kokkos::deep_copy(dofCoeffs12,0.0);
539 template<
class HGRAD_LINE,
class HVOL_LINE>
556 Teuchos::rcp(new
Family2(polyOrder_x, polyOrder_y, polyOrder_z, pointType))) {}
559 template<
class HGRAD_LINE,
class HVOL_LINE>
568 ordinal_type order_x_;
569 ordinal_type order_y_;
570 ordinal_type order_z_;
571 EPointType pointType_;
575 using ExecutionSpace =
typename HGRAD_LINE::ExecutionSpace;
576 using OutputValueType =
typename HGRAD_LINE::OutputValueType;
577 using PointValueType =
typename HGRAD_LINE::PointValueType;
579 using BasisBase =
typename HGRAD_LINE::BasisBase;
590 Teuchos::rcp(new
Family3 (polyOrder_x, polyOrder_y, polyOrder_z, pointType))) {
591 this->functionSpace_ = FUNCTION_SPACE_HCURL;
593 std::ostringstream basisName;
595 name_ = basisName.str();
597 order_x_ = polyOrder_x;
598 order_y_ = polyOrder_y;
599 order_z_ = polyOrder_z;
600 pointType_ = pointType;
621 return name_.c_str();
634 Teuchos::RCP<BasisBase>
637 using LineBasis = HVOL_LINE;
640 if(subCellDim == 1) {
646 return Teuchos::rcp(
new LineBasis(order_x_-1, pointType_) );
651 return Teuchos::rcp(
new LineBasis(order_y_-1, pointType_) );
656 return Teuchos::rcp(
new LineBasis(order_z_-1, pointType_) );
658 }
else if(subCellDim == 2) {
661 return Teuchos::rcp(
new QuadBasis(order_x_, order_z_, pointType_) );
663 return Teuchos::rcp(
new QuadBasis(order_y_,order_z_, pointType_) );
665 return Teuchos::rcp(
new QuadBasis(order_x_, order_z_, pointType_) );
667 return Teuchos::rcp(
new QuadBasis(order_z_, order_y_, pointType_) );
669 return Teuchos::rcp(
new QuadBasis(order_y_, order_x_, pointType_) );
671 return Teuchos::rcp(
new QuadBasis(order_x_, order_y_, pointType_) );
675 INTREPID2_TEST_FOR_EXCEPTION(
true,std::invalid_argument,
"Input parameters out of bounds");
681 virtual HostBasisPtr<OutputValueType, PointValueType>
685 auto hostBasis = Teuchos::rcp(
new HostBasis(order_x_, order_y_, order_z_, pointType_));
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.
virtual bool requireOrientation() const override
True if orientation is required.
Basis_Derived_HCURL_Family1_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Basis_Derived_HCURL_Family1_Family2_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType)
Constructor.
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...
Implementation of bases that are tensor products of two or three component bases. ...
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(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 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_Family2_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
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.
Basis_Derived_HCURL_Family3_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
Basis_Derived_HCURL_HEX(int polyOrder_x, int polyOrder_y, int polyOrder_z, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
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 getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis 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)
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...
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)
Basis_Derived_HCURL_HEX(int polyOrder, const EPointType pointType=POINTTYPE_DEFAULT)
Constructor.
virtual HostBasisPtr< OutputValueType, PointValueType > getHostBasis() const override
Creates and returns a Basis object whose DeviceType template argument is Kokkos::HostSpace::device_ty...
virtual void getDofCoeffs(ScalarViewType dofCoeffs) const override
Fills in coefficients of degrees of freedom for Lagrangian basis on the reference cell...