50 #ifndef __INTREPID2_BASIS_HPP__
51 #define __INTREPID2_BASIS_HPP__
53 #include "Intrepid2_ConfigDefs.hpp"
58 #include "Kokkos_Vector.hpp"
59 #include "Shards_CellTopology.hpp"
92 template<
typename ExecSpaceType = void,
93 typename outputValueType = double,
94 typename pointValueType =
double>
112 using ordinal_view_type INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalViewType instead",
"OrdinalViewType") =
OrdinalViewType INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalViewType instead");
117 using ebasis_view_type INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use EBasisViewType instead",
"EBasisViewType") =
EBasisViewType INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use EBasisViewType instead");
122 using ecoordinates_view_type INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use ECoordinatesViewType instead",
"ECoordinatesViewType") =
ECoordinatesViewType INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use ECoordinatesViewType instead");
130 using ordinal_type_array_1d_host INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArray1DHost instead",
"OrdinalTypeArray1DHost") =
OrdinalTypeArray1DHost INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArray1DHost instead");
135 using ordinal_type_array_2d_host INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArray2DHost instead",
"OrdinalTypeArray2DHost") =
OrdinalTypeArray2DHost INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArray2DHost instead");
140 using ordinal_type_array_3d_host INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArray3DHost instead",
"OrdinalTypeArray3DHost") =
OrdinalTypeArray3DHost INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArray3DHost instead");
145 using ordinal_type_array_stride_1d_host INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArrayStride1DHost instead",
"OrdinalTypeArrayStride1DHost") =
OrdinalTypeArrayStride1DHost INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArrayStride1DHost instead");
150 using ordinal_type_array_1d INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArray1D instead",
"OrdinalTypeArray1D") =
OrdinalTypeArray1D INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArray1D instead");
155 using ordinal_type_array_2d INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArray2D instead",
"OrdinalTypeArray2D") =
OrdinalTypeArray2D INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArray2D instead");
160 using ordinal_type_array_3d INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArray3D instead",
"OrdinalTypeArray3D") =
OrdinalTypeArray3D INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArray3D instead");
165 using ordinal_type_array_stride_1d INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OrdinalTypeArrayStride1D instead",
"OrdinalTypeArrayStride1D") =
OrdinalTypeArrayStride1D INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OrdinalTypeArrayStride1D instead");
169 typedef typename ScalarTraits<pointValueType>::scalar_type
scalarType;
239 template<
typename OrdinalTypeView3D,
240 typename OrdinalTypeView2D,
241 typename OrdinalTypeView1D>
243 OrdinalTypeView2D &ordinalToTag,
244 const OrdinalTypeView1D tags,
245 const ordinal_type basisCard,
246 const ordinal_type tagSize,
247 const ordinal_type posScDim,
248 const ordinal_type posScOrd,
249 const ordinal_type posDfOrd ) {
251 ordinalToTag = OrdinalTypeView2D(
"ordinalToTag", basisCard, tagSize);
254 Kokkos::deep_copy( ordinalToTag, -1 );
257 for (ordinal_type i=0;i<basisCard;++i)
258 for (ordinal_type j=0;j<tagSize;++j)
259 ordinalToTag(i, j) = tags(i*tagSize + j);
263 for (ordinal_type i=0;i<basisCard;++i)
264 if (maxScDim < tags(i*tagSize + posScDim))
265 maxScDim = tags(i*tagSize + posScDim);
269 for (ordinal_type i=0;i<basisCard;++i)
270 if (maxScOrd < tags(i*tagSize + posScOrd))
271 maxScOrd = tags(i*tagSize + posScOrd);
275 for (ordinal_type i=0;i<basisCard;++i)
276 if (maxDfOrd < tags(i*tagSize + posDfOrd))
277 maxDfOrd = tags(i*tagSize + posDfOrd);
281 tagToOrdinal = OrdinalTypeView3D(
"tagToOrdinal", maxScDim, maxScOrd, maxDfOrd);
284 Kokkos::deep_copy( tagToOrdinal, -1 );
287 for (ordinal_type i=0;i<basisCard;++i)
288 tagToOrdinal(tags(i*tagSize), tags(i*tagSize+1), tags(i*tagSize+2)) = i;
317 virtual~
Basis() =
default;
330 using OutputViewType = Kokkos::DynRankView<OutputValueType,Kokkos::LayoutStride,ExecSpaceType>;
331 using outputViewType INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use OutputViewType instead",
"OutputViewType") =
OutputViewType INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use OutputViewType instead");
335 using PointViewType = Kokkos::DynRankView<PointValueType,Kokkos::LayoutStride,ExecSpaceType>;
336 using pointViewType INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use PointViewType instead",
"PointViewType") =
PointViewType INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use PointViewType instead");
340 using ScalarViewType = Kokkos::DynRankView<scalarType,Kokkos::LayoutStride,ExecSpaceType>;
341 using scalarViewType INTREPID2_DEPRECATED_TYPENAME_REPLACEMENT(
"use ScalarViewType instead",
"ScalarViewType") =
ScalarViewType INTREPID2_DEPRECATED_TYPENAME_TRAILING_ATTRIBUTE(
"use ScalarViewType instead");
365 const EOperator = OPERATOR_VALUE )
const {
366 INTREPID2_TEST_FOR_EXCEPTION(
true, std::logic_error,
367 ">>> ERROR (Basis::getValues): this method (FEM) is not supported or should be over-riden accordingly by derived classes.");
394 const EOperator = OPERATOR_VALUE )
const {
395 INTREPID2_TEST_FOR_EXCEPTION(
true, std::logic_error,
396 ">>> ERROR (Basis::getValues): this method (FVM) is not supported or should be over-riden accordingly by derived classes.");
406 INTREPID2_TEST_FOR_EXCEPTION(
true, std::logic_error,
407 ">>> ERROR (Basis::getDofCoords): this method is not supported or should be over-riden accordingly by derived classes.");
421 INTREPID2_TEST_FOR_EXCEPTION(
true, std::logic_error,
422 ">>> ERROR (Basis::getDofCoeffs): this method is not supported or should be over-riden accordingly by derived classes.");
434 INTREPID2_TEST_FOR_EXCEPTION(
basisType_ != BASIS_FEM_HIERARCHICAL, std::logic_error,
435 ">>> ERROR (Basis::getFieldOrdinalsForDegree): this method is not supported for non-hierarchical bases.");
437 int requestedDegreeLength = degrees.extent_int(0);
438 INTREPID2_TEST_FOR_EXCEPTION(degreeEntryLength != requestedDegreeLength, std::invalid_argument,
"length of degrees does not match the entries in fieldOrdinalPolynomialDegree_");
439 Kokkos::vector<int> fieldOrdinalsVector;
443 for (
int d=0; d<degreeEntryLength; d++)
447 if (matches) fieldOrdinalsVector.push_back(basisOrdinal);
450 for (
int i=0; i<fieldOrdinalsVector.size(); i++)
452 fieldOrdinals(i) = fieldOrdinalsVector[i];
454 return fieldOrdinals;
469 INTREPID2_TEST_FOR_EXCEPTION(
basisType_ != BASIS_FEM_HIERARCHICAL, std::logic_error,
470 ">>> ERROR (Basis::getFieldOrdinalsForDegree): this method is not supported for non-hierarchical bases.");
472 for (
int d=0; d<degrees.size(); d++)
474 degreesView(d) = degrees[d];
477 Kokkos::vector<int> fieldOrdinalsVector(fieldOrdinalsView.extent_int(0));
478 for (
int i=0; i<fieldOrdinalsView.extent_int(0); i++)
480 fieldOrdinalsVector[i] = fieldOrdinalsView(i);
482 return fieldOrdinalsVector;
493 INTREPID2_TEST_FOR_EXCEPTION(
basisType_ != BASIS_FEM_HIERARCHICAL, std::logic_error,
494 ">>> ERROR (Basis::getPolynomialDegreeOfField): this method is not supported for non-hierarchical bases.");
495 INTREPID2_TEST_FOR_EXCEPTION(fieldOrdinal < 0, std::invalid_argument,
"field ordinal must be non-negative");
496 INTREPID2_TEST_FOR_EXCEPTION(fieldOrdinal >=
fieldOrdinalPolynomialDegree_.extent_int(0), std::invalid_argument,
"field ordinal out of bounds");
510 INTREPID2_TEST_FOR_EXCEPTION(
basisType_ != BASIS_FEM_HIERARCHICAL, std::logic_error,
511 ">>> ERROR (Basis::getPolynomialDegreeOfFieldAsVector): this method is not supported for non-hierarchical bases.");
513 Kokkos::vector<int> polynomialDegree(polynomialDegreeView.extent_int(0));
515 for (
int d=0; d<polynomialDegree.size(); d++)
517 polynomialDegree[d] = polynomialDegreeView(d);
519 return polynomialDegree;
526 INTREPID2_TEST_FOR_EXCEPTION(
basisType_ != BASIS_FEM_HIERARCHICAL, std::logic_error,
527 ">>> ERROR (Basis::getPolynomialDegreeLength): this method is not supported for non-hierarchical bases.");
538 return "Intrepid2_Basis";
616 const ordinal_type subcOrd )
const {
617 if ( subcDim >= 0 && subcDim <
static_cast<ordinal_type
>(
tagToOrdinal_.extent(0)) &&
618 subcOrd >= 0 && subcOrd <
static_cast<ordinal_type
>(
tagToOrdinal_.extent(1)) )
621 if (firstDofOrdinal == -1)
return static_cast<ordinal_type
>(0);
623 return static_cast<ordinal_type
>(this->
getDofTag(firstDofOrdinal)[3]);
628 return static_cast<ordinal_type
>(0);
642 const ordinal_type subcOrd,
643 const ordinal_type subcDofOrd )
const {
645 #ifdef HAVE_INTREPID2_DEBUG
646 INTREPID2_TEST_FOR_EXCEPTION( subcDim < 0 || subcDim >=
static_cast<ordinal_type
>(
tagToOrdinal_.extent(0)), std::out_of_range,
647 ">>> ERROR (Basis::getDofOrdinal): subcDim is out of range");
648 INTREPID2_TEST_FOR_EXCEPTION( subcOrd < 0 || subcOrd >=
static_cast<ordinal_type
>(
tagToOrdinal_.extent(1)), std::out_of_range,
649 ">>> ERROR (Basis::getDofOrdinal): subcOrd is out of range");
650 INTREPID2_TEST_FOR_EXCEPTION( subcDofOrd < 0 || subcDofOrd >=
static_cast<ordinal_type
>(
tagToOrdinal_.extent(2)), std::out_of_range,
651 ">>> ERROR (Basis::getDofOrdinal): subcDofOrd is out of range");
653 ordinal_type r_val = -1;
654 if ( subcDim <
static_cast<ordinal_type
>(
tagToOrdinal_.extent(0)) &&
655 subcOrd <
static_cast<ordinal_type
>(
tagToOrdinal_.extent(1)) &&
656 subcDofOrd <
static_cast<ordinal_type
>(
tagToOrdinal_.extent(2)) )
658 #ifdef HAVE_INTREPID2_DEBUG
659 INTREPID2_TEST_FOR_EXCEPTION( r_val == -1, std::runtime_error,
660 ">>> ERROR (Basis::getDofOrdinal): Invalid DoF tag is found.");
684 #ifdef HAVE_INTREPID2_DEBUG
685 INTREPID2_TEST_FOR_EXCEPTION( dofOrd < 0 || dofOrd >=
static_cast<ordinal_type
>(
ordinalToTag_.extent(0)), std::out_of_range,
686 ">>> ERROR (Basis::getDofTag): dofOrd is out of range");
688 return Kokkos::subview(
ordinalToTag_, dofOrd, Kokkos::ALL());
730 KOKKOS_INLINE_FUNCTION
731 ordinal_type getFieldRank(
const EFunctionSpace spaceType);
768 KOKKOS_INLINE_FUNCTION
769 ordinal_type getOperatorRank(
const EFunctionSpace spaceType,
770 const EOperator operatorType,
771 const ordinal_type spaceDim);
778 KOKKOS_INLINE_FUNCTION
779 ordinal_type getOperatorOrder(
const EOperator operatorType);
781 template<EOperator operatorType>
782 KOKKOS_INLINE_FUNCTION
783 constexpr ordinal_type getOperatorOrder();
808 template<ordinal_type spaceDim>
809 KOKKOS_INLINE_FUNCTION
810 ordinal_type getDkEnumeration(
const ordinal_type xMult,
811 const ordinal_type yMult = -1,
812 const ordinal_type zMult = -1);
825 template<ordinal_type spaceDim>
826 KOKKOS_INLINE_FUNCTION
827 ordinal_type getPnEnumeration(
const ordinal_type p,
828 const ordinal_type q = 0,
829 const ordinal_type r = 0);
851 template<
typename value_type>
852 KOKKOS_INLINE_FUNCTION
853 void getJacobyRecurrenceCoeffs (
857 const ordinal_type alpha,
858 const ordinal_type beta ,
859 const ordinal_type n);
898 KOKKOS_INLINE_FUNCTION
899 ordinal_type getDkCardinality(
const EOperator operatorType,
900 const ordinal_type spaceDim);
902 template<EOperator operatorType, ordinal_type spaceDim>
903 KOKKOS_INLINE_FUNCTION
904 constexpr ordinal_type getDkCardinality();
917 template<ordinal_type spaceDim>
918 KOKKOS_INLINE_FUNCTION
919 ordinal_type getPnCardinality (ordinal_type n);
921 template<ordinal_type spaceDim, ordinal_type n>
922 KOKKOS_INLINE_FUNCTION
923 constexpr ordinal_type getPnCardinality ();
942 template<
typename outputValueViewType,
943 typename inputPointViewType>
944 void getValues_HGRAD_Args(
const outputValueViewType outputValues,
945 const inputPointViewType inputPoints,
946 const EOperator operatorType,
947 const shards::CellTopology cellTopo,
948 const ordinal_type basisCard );
960 template<
typename outputValueViewType,
961 typename inputPointViewType>
962 void getValues_HCURL_Args(
const outputValueViewType outputValues,
963 const inputPointViewType inputPoints,
964 const EOperator operatorType,
965 const shards::CellTopology cellTopo,
966 const ordinal_type basisCard );
978 template<
typename outputValueViewType,
979 typename inputPointViewType>
980 void getValues_HDIV_Args(
const outputValueViewType outputValues,
981 const inputPointViewType inputPoints,
982 const EOperator operatorType,
983 const shards::CellTopology cellTopo,
984 const ordinal_type basisCard );
996 template<
typename outputValueViewType,
997 typename inputPointViewType>
998 void getValues_HVOL_Args(
const outputValueViewType outputValues,
999 const inputPointViewType inputPoints,
1000 const EOperator operatorType,
1001 const shards::CellTopology cellTopo,
1002 const ordinal_type basisCard );