Intrepid2
Static Public Member Functions | List of all members
Intrepid2::FunctionSpaceTools< ExecSpaceType > Class Template Reference

Defines expert-level interfaces for the evaluation of functions and operators in physical space (supported for FE, FV, and FD methods) and FE reference space; in addition, provides several function transformation utilities. More...

#include <Intrepid2_FunctionSpaceTools.hpp>

Static Public Member Functions

template<typename outputValueType , class... outputProperties, typename inputValueType , class... inputProperties>
static void HGRADtransformVALUE (Kokkos::DynRankView< outputValueType, outputProperties...> output, const Kokkos::DynRankView< inputValueType, inputProperties...> input)
 Transformation of a (scalar) value field in the H-grad space, defined at points on a reference cell, stored in the user-provided container input and indexed by (F,P), into the output container output, defined on cells in physical space and indexed by (C,F,P). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianInverseValueType , class... jacobianInverseProperties, typename inputValValueType , class... inputValProperties>
static void HGRADtransformGRAD (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianInverseValueType, jacobianInverseProperties...> jacobianInverse, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a gradient field in the H-grad space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianInverseValueType , class... jacobianInverseProperties, typename inputValValueType , class... inputValProperties>
static void HCURLtransformVALUE (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianInverseValueType, jacobianInverseProperties...> jacobianInverse, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a (vector) value field in the H-curl space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianValueType , class... jacobianProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
static void HCURLtransformCURL (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianValueType, jacobianProperties...> jacobian, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...> jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a 3D curl field in the H-curl space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
static void HCURLtransformCURL (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...> jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a 2D curl field in the H-curl space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianValueType , class... jacobianProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
static void HGRADtransformCURL (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianValueType, jacobianProperties...> jacobian, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...> jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a 2D curl field in the H-grad space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianValueType , class... jacobianProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
static void HDIVtransformVALUE (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianValueType, jacobianProperties...> jacobian, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...> jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a (vector) value field in the H-div space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
static void HDIVtransformDIV (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...> jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a divergence field in the H-div space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P). More...
 
template<typename outputValValueType , class... outputValProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
static void HVOLtransformVALUE (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...> jacobianDet, const Kokkos::DynRankView< inputValValueType, inputValProperties...> inputVals)
 Transformation of a (scalar) value field in the H-vol space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P). More...
 
template<typename outputValueValueType , class... outputValueProperties, typename leftValueValueType , class... leftValueProperties, typename rightValueValueType , class... rightValueProperties>
static void integrate (Kokkos::DynRankView< outputValueValueType, outputValueProperties...> outputValues, const Kokkos::DynRankView< leftValueValueType, leftValueProperties...> leftValues, const Kokkos::DynRankView< rightValueValueType, rightValueProperties...> rightValues, const bool sumInto=false)
 Contracts leftValues and rightValues arrays on the point and possibly space dimensions and stores the result in outputValues; this is a generic, high-level integration routine that calls either FunctionSpaceTools::operatorIntegral, or FunctionSpaceTools::functionalIntegral, or FunctionSpaceTools::dataIntegral methods, depending on the rank of the outputValues array. More...
 
template<typename outputValValueType , class... outputValProperties, typename inputDetValueType , class... inputDetPropertes, typename inputWeightValueType , class... inputWeightPropertes>
static bool computeCellMeasure (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< inputDetValueType, inputDetPropertes...> inputDet, const Kokkos::DynRankView< inputWeightValueType, inputWeightPropertes...> inputWeights)
 Returns the weighted integration measures outputVals with dimensions (C,P) used for the computation of cell integrals, by multiplying absolute values of the user-provided cell Jacobian determinants inputDet with dimensions (C,P) with the user-provided integration weights inputWeights with dimensions (P). More...
 
template<typename outputValValueType , class... outputValProperties, typename inputJacValueType , class... inputJacProperties, typename inputWeightValueType , class... inputWeightPropertes, typename scratchValueType , class... scratchProperties>
static void computeFaceMeasure (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< inputJacValueType, inputJacProperties...> inputJac, const Kokkos::DynRankView< inputWeightValueType, inputWeightPropertes...> inputWeights, const int whichFace, const shards::CellTopology parentCell, const Kokkos::DynRankView< scratchValueType, scratchProperties...> scratch)
 Returns the weighted integration measures outputVals with dimensions (C,P) used for the computation of face integrals, based on the provided cell Jacobian array inputJac with dimensions (C,P,D,D) and the provided integration weights inputWeights with dimensions (P). More...
 
template<typename outputValValueType , class... outputValProperties, typename inputJacValueType , class... inputJacProperties, typename inputWeightValueType , class... inputWeightProperties, typename scratchValueType , class... scratchProperties>
static void computeEdgeMeasure (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< inputJacValueType, inputJacProperties...> inputJac, const Kokkos::DynRankView< inputWeightValueType, inputWeightProperties...> inputWeights, const int whichEdge, const shards::CellTopology parentCell, const Kokkos::DynRankView< scratchValueType, scratchProperties...> scratch)
 Returns the weighted integration measures outVals with dimensions (C,P) used for the computation of edge integrals, based on the provided cell Jacobian array inputJac with dimensions (C,P,D,D) and the provided integration weights inWeights with dimensions (P). More...
 
template<typename outputValValueType , class... outputValProperties, typename inputMeasureValueType , class... inputMeasureProperties, typename inputValValueType , class... inputValProperteis>
static void multiplyMeasure (Kokkos::DynRankView< outputValValueType, outputValProperties...> outputVals, const Kokkos::DynRankView< inputMeasureValueType, inputMeasureProperties...> inputMeasure, const Kokkos::DynRankView< inputValValueType, inputValProperteis...> inputVals)
 Multiplies fields inputVals by weighted measures inputMeasure and returns the field array outputVals; this is a simple redirection to the call FunctionSpaceTools::scalarMultiplyDataField. More...
 
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataPropertes, typename inputFieldValueType , class... inputFieldProperties>
static void scalarMultiplyDataField (Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...> outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataPropertes...> inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...> inputFields, const bool reciprocal=false)
 Scalar multiplication of data and fields; please read the description below. More...
 
template<typename outputDataValuetype , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
static void scalarMultiplyDataData (Kokkos::DynRankView< outputDataValuetype, outputDataProperties...> outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...> inputDataRight, const bool reciprocal=false)
 Scalar multiplication of data and data; please read the description below. More...
 
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataProperties, typename inputFieldValueType , class... inputFieldProperties>
static void dotMultiplyDataField (Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...> outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties...> inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...> inputFields)
 Dot product of data and fields; please read the description below. More...
 
template<typename outputDataValueType , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
static void dotMultiplyDataData (Kokkos::DynRankView< outputDataValueType, outputDataProperties...> outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...> inputDataRight)
 Dot product of data and data; please read the description below. More...
 
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataProperties, typename inputFieldValueType , class... inputFieldProperties>
static void vectorMultiplyDataField (Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...> outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties...> inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...> inputFields)
 Cross or outer product of data and fields; please read the description below. More...
 
template<typename outputDataValueType , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
static void vectorMultiplyDataData (Kokkos::DynRankView< outputDataValueType, outputDataProperties...> outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...> inputDataRight)
 Cross or outer product of data and data; please read the description below. More...
 
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataProperties, typename inputFieldValueType , class... inputFieldProperties>
static void tensorMultiplyDataField (Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...> outputFields, const Kokkos::DynRankView< inputDataValueType, inputDataProperties...> inputData, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...> inputFields, const char transpose= 'N')
 Matrix-vector or matrix-matrix product of data and fields; please read the description below. More...
 
template<typename outputDataValueType , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
static void tensorMultiplyDataData (Kokkos::DynRankView< outputDataValueType, outputDataProperties...> outputData, const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...> inputDataLeft, const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...> inputDataRight, const char transpose= 'N')
 Matrix-vector or matrix-matrix product of data and data; please read the description below. More...
 
template<typename inoutOperatorValueType , class... inoutOperatorProperties, typename fieldSignValueType , class... fieldSignProperties>
static void applyLeftFieldSigns (Kokkos::DynRankView< inoutOperatorValueType, inoutOperatorProperties...> inoutOperator, const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties...> fieldSigns)
 Applies left (row) signs, stored in the user-provided container fieldSigns and indexed by (C,L), to the operator inoutOperator indexed by (C,L,R). More...
 
template<typename inoutOperatorValueType , class... inoutOperatorProperties, typename fieldSignValueType , class... fieldSignProperties>
static void applyRightFieldSigns (Kokkos::DynRankView< inoutOperatorValueType, inoutOperatorProperties...> inoutOperator, const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties...> fieldSigns)
 Applies right (column) signs, stored in the user-provided container fieldSigns and indexed by (C,R), to the operator inoutOperator indexed by (C,L,R). More...
 
template<typename inoutFunctionValueType , class... inoutFunctionProperties, typename fieldSignValueType , class... fieldSignProperties>
static void applyFieldSigns (Kokkos::DynRankView< inoutFunctionValueType, inoutFunctionProperties...> inoutFunction, const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties...> fieldSigns)
 Applies field signs, stored in the user-provided container fieldSigns and indexed by (C,F), to the function inoutFunction indexed by (C,F), (C,F,P), (C,F,P,D1) or (C,F,P,D1,D2). More...
 
template<typename outputPointValueType , class... outputPointProperties, typename inputCoeffValueType , class... inputCoeffProperties, typename inputFieldValueType , class... inputFieldProperties>
static void evaluate (Kokkos::DynRankView< outputPointValueType, outputPointProperties...> outputPointVals, const Kokkos::DynRankView< inputCoeffValueType, inputCoeffProperties...> inputCoeffs, const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...> inputFields)
 Computes point values outPointVals of a discrete function specified by the basis inFields and coefficients inCoeffs. More...
 

Detailed Description

template<typename ExecSpaceType = void>
class Intrepid2::FunctionSpaceTools< ExecSpaceType >

Defines expert-level interfaces for the evaluation of functions and operators in physical space (supported for FE, FV, and FD methods) and FE reference space; in addition, provides several function transformation utilities.

Definition at line 78 of file Intrepid2_FunctionSpaceTools.hpp.

Member Function Documentation

template<typename SpT >
template<typename inoutFunctionValueType , class... inoutFunctionProperties, typename fieldSignValueType , class... fieldSignProperties>
void Intrepid2::FunctionSpaceTools< SpT >::applyFieldSigns ( Kokkos::DynRankView< inoutFunctionValueType, inoutFunctionProperties...>  inoutFunction,
const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties...>  fieldSigns 
)
static

Applies field signs, stored in the user-provided container fieldSigns and indexed by (C,F), to the function inoutFunction indexed by (C,F), (C,F,P), (C,F,P,D1) or (C,F,P,D1,D2).

Returns

\[ \mbox{inoutFunction}(c,f,*) = \mbox{fieldSigns}(c,f)*\mbox{inoutFunction}(c,f,*) \]

See Section Pullbacks for discussion of field signs.

C - num. integration domains
F - num. fields
P - num. integration points
D1 - spatial dimension
D2 - spatial dimension
Parameters
inoutFunction[in/out] - Input / output function array.
fieldSigns[in] - Right field signs.

Definition at line 933 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename inoutOperatorValueType , class... inoutOperatorProperties, typename fieldSignValueType , class... fieldSignProperties>
void Intrepid2::FunctionSpaceTools< SpT >::applyLeftFieldSigns ( Kokkos::DynRankView< inoutOperatorValueType, inoutOperatorProperties...>  inoutOperator,
const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties...>  fieldSigns 
)
static

Applies left (row) signs, stored in the user-provided container fieldSigns and indexed by (C,L), to the operator inoutOperator indexed by (C,L,R).

Mathematically, this method computes the matrix-matrix product

\[ \mathbf{K}^{c} = \mbox{diag}(\sigma^c_0,\ldots,\sigma^c_{L-1}) \mathbf{K}^c \]

where $\mathbf{K}^{c} \in \mathbf{R}^{L\times R}$ is array of matrices indexed by cell number c and stored in the rank-3 array inoutOperator, and $\{\sigma^c_l\}_{l=0}^{L-1}$ is array of left field signs indexed by cell number c and stored in the rank-2 container fieldSigns; see Section Pullbacks for discussion of field signs. This operation is required for operators generated by HCURL and HDIV-conforming vector-valued finite element basis functions; see Sections Pullbacks and Section Evaluation of finite element operators and functionals for applications of this method.

C - num. integration domains
L - num. left fields
R - num. right fields
Parameters
inoutOperator[in/out] - Input / output operator array.
fieldSigns[in] - Left field signs.

Definition at line 811 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename inoutOperatorValueType , class... inoutOperatorProperties, typename fieldSignValueType , class... fieldSignProperties>
void Intrepid2::FunctionSpaceTools< SpT >::applyRightFieldSigns ( Kokkos::DynRankView< inoutOperatorValueType, inoutOperatorProperties...>  inoutOperator,
const Kokkos::DynRankView< fieldSignValueType, fieldSignProperties...>  fieldSigns 
)
static

Applies right (column) signs, stored in the user-provided container fieldSigns and indexed by (C,R), to the operator inoutOperator indexed by (C,L,R).

Mathematically, this method computes the matrix-matrix product

\[ \mathbf{K}^{c} = \mathbf{K}^c \mbox{diag}(\sigma^c_0,\ldots,\sigma^c_{R-1}) \]

where $\mathbf{K}^{c} \in \mathbf{R}^{L\times R}$ is array of matrices indexed by cell number c and stored in the rank-3 container inoutOperator, and $\{\sigma^c_r\}_{r=0}^{R-1}$ is array of right field signs indexed by cell number c and stored in the rank-2 container fieldSigns; see Section Pullbacks for discussion of field signs. This operation is required for operators generated by HCURL and HDIV-conforming vector-valued finite element basis functions; see Sections Pullbacks and Section Evaluation of finite element operators and functionals for applications of this method.

C - num. integration domains
L - num. left fields
R - num. right fields
Parameters
inoutOperator[in/out] - Input / output operator array.
fieldSigns[in] - Right field signs.

Definition at line 870 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename inputDetValueType , class... inputDetProperties, typename inputWeightValueType , class... inputWeightProperties>
bool Intrepid2::FunctionSpaceTools< SpT >::computeCellMeasure ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< inputDetValueType, inputDetProperties...>  inputDet,
const Kokkos::DynRankView< inputWeightValueType, inputWeightProperties...>  inputWeights 
)
static

Returns the weighted integration measures outputVals with dimensions (C,P) used for the computation of cell integrals, by multiplying absolute values of the user-provided cell Jacobian determinants inputDet with dimensions (C,P) with the user-provided integration weights inputWeights with dimensions (P).

Returns a rank-2 array (C, P) array such that

\[ \mbox{outputVals}(c,p) = |\mbox{det}(DF_{c}(\widehat{x}_p))|\omega_{p} \,, \]

where $\{(\widehat{x}_p,\omega_p)\}$ is a cubature rule defined on a reference cell (a set of integration points and their associated weights; see Intrepid2::Cubature::getCubature for getting cubature rules on reference cells).

Warning
The user is responsible for providing input arrays with consistent data: the determinants in inputDet should be evaluated at integration points on the reference cell corresponding to the weights in inputWeights.
Remarks
See Intrepid2::CellTools::setJacobian for computation of DF and Intrepid2::CellTools::setJacobianDet for computation of its determinant.
C - num. integration domains dim0 in all containers
P - num. integration points dim1 in all containers
Parameters
outputVals[out] - Output array with weighted cell measures.
inputDet[in] - Input array containing determinants of cell Jacobians.
inputWeights[in] - Input integration weights.

Definition at line 438 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename inputJacValueType , class... inputJacProperties, typename inputWeightValueType , class... inputWeightProperties, typename scratchValueType , class... scratchProperties>
void Intrepid2::FunctionSpaceTools< SpT >::computeEdgeMeasure ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< inputJacValueType, inputJacProperties...>  inputJac,
const Kokkos::DynRankView< inputWeightValueType, inputWeightProperties...>  inputWeights,
const int  whichEdge,
const shards::CellTopology  parentCell,
const Kokkos::DynRankView< scratchValueType, scratchProperties...>  scratch 
)
static

Returns the weighted integration measures outVals with dimensions (C,P) used for the computation of edge integrals, based on the provided cell Jacobian array inputJac with dimensions (C,P,D,D) and the provided integration weights inWeights with dimensions (P).

Returns a rank-2 array (C, P) array such that

\[ \mbox{outVals}(c,p) = \left\|\frac{d \Phi_c(\widehat{x}_p)}{d s}\right\|\omega_{p} \,, \]

where:

  • $\{(\widehat{x}_p,\omega_p)\}$ is a cubature rule defined on reference edge $\widehat{\mathcal{E}}$, with ordinal whichEdge relative to the specified parent reference cell;
  • $ \Phi_c : R \mapsto \mathcal{E} $ is parameterization of the physical edge corresponding to $\widehat{\mathcal{E}}$; see Section Parametrization of physical 1- and 2-subcells.
Warning
The user is responsible for providing input arrays with consistent data: the Jacobians in inputJac should be evaluated at integration points on the reference edge corresponding to the weights in inputWeights. Additionally, the user is responsible for providing a scratch array of the same dimension as the inputJac array.
Remarks
Cubature rules on reference edges are defined by a two-step process:
See Intrepid2::CellTools::setJacobian for computation of DF and Intrepid2::CellTools::setJacobianDet for computation of its determinant.
C - num. integration domains dim0 in all input containers
P - num. integration points dim1 in all input containers
D - spatial dimension dim2 and dim3 in Jacobian container
Parameters
outputVals[out] - Output array with weighted edge measures.
inputJac[in] - Input array containing cell Jacobians.
inputWeights[in] - Input integration weights.
whichEdge[in] - Index of the edge subcell relative to the parent cell; defines the domain of integration.
parentCell[in] - Parent cell topology.
scratch[in] - Scratch space, sized like inputJac

Definition at line 525 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::CellTools< ExecSpaceType >::getPhysicalEdgeTangents(), Intrepid2::ArrayTools< ExecSpaceType >::scalarMultiplyDataData(), and Intrepid2::RealSpaceTools< ExecSpaceType >::vectorNorm().

Referenced by Intrepid2::CubatureControlVolumeBoundary< ExecSpaceType, pointValueType, weightValueType >::getCubature().

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename inputJacValueType , class... inputJacProperties, typename inputWeightValueType , class... inputWeightProperties, typename scratchValueType , class... scratchProperties>
void Intrepid2::FunctionSpaceTools< SpT >::computeFaceMeasure ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< inputJacValueType, inputJacProperties...>  inputJac,
const Kokkos::DynRankView< inputWeightValueType, inputWeightProperties...>  inputWeights,
const int  whichFace,
const shards::CellTopology  parentCell,
const Kokkos::DynRankView< scratchValueType, scratchProperties...>  scratch 
)
static

Returns the weighted integration measures outputVals with dimensions (C,P) used for the computation of face integrals, based on the provided cell Jacobian array inputJac with dimensions (C,P,D,D) and the provided integration weights inputWeights with dimensions (P).

Returns a rank-2 array (C, P) array such that

\[ \mbox{outputVals}(c,p) = \left\|\frac{\partial\Phi_c(\widehat{x}_p)}{\partial u}\times \frac{\partial\Phi_c(\widehat{x}_p)}{\partial v}\right\|\omega_{p} \,, \]

where:

  • $\{(\widehat{x}_p,\omega_p)\}$ is a cubature rule defined on reference face $\widehat{\mathcal{F}}$, with ordinal whichFace relative to the specified parent reference cell;
  • $ \Phi_c : R \mapsto \mathcal{F} $ is parameterization of the physical face corresponding to $\widehat{\mathcal{F}}$; see Section Parametrization of physical 1- and 2-subcells.
Warning
The user is responsible for providing input arrays with consistent data: the Jacobians in inputJac should be evaluated at integration points on the reference face corresponding to the weights in inputWeights. Additionally, the user is responsible for providing a scratch array of the same dimension as the inputJac array.
Remarks
Cubature rules on reference faces are defined by a two-step process:
  • A cubature rule is defined on the parametrization domain R of the face (R is the standard 2-simplex {(0,0),(1,0),(0,1)} or the standard 2-cube [-1,1] X [-1,1]).
  • The points are mapped to a reference face using Intrepid2::CellTools::mapToReferenceSubcell
See Intrepid2::CellTools::setJacobian for computation of DF and Intrepid2::CellTools::setJacobianDet for computation of its determinant.
C - num. integration domains dim0 in all input containers
P - num. integration points dim1 in all input containers
D - spatial dimension dim2 and dim3 in Jacobian and scratch containers
Parameters
outputVals[out] - Output array with weighted face measures.
inputJac[in] - Input array containing cell Jacobians.
inputWeights[in] - Input integration weights.
whichFace[in] - Index of the face subcell relative to the parent cell; defines the domain of integration.
parentCell[in] - Parent cell topology.
scratch[in] - Scratch space, sized like inputJac

Definition at line 478 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::CellTools< ExecSpaceType >::getPhysicalFaceNormals(), Intrepid2::ArrayTools< ExecSpaceType >::scalarMultiplyDataData(), and Intrepid2::RealSpaceTools< ExecSpaceType >::vectorNorm().

Referenced by Intrepid2::CubatureControlVolumeBoundary< ExecSpaceType, pointValueType, weightValueType >::getCubature().

template<typename SpT >
template<typename outputDataValueType , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
void Intrepid2::FunctionSpaceTools< SpT >::dotMultiplyDataData ( Kokkos::DynRankView< outputDataValueType, outputDataProperties...>  outputData,
const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...>  inputDataLeft,
const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...>  inputDataRight 
)
static

Dot product of data and data; please read the description below.

There are two use cases:

  • dot product of a rank-2, 3 or 4 container inputDataRight with dimensions (C,P) (C,P,D1) or (C,P,D1,D2), representing the values of a scalar, vector or a tensor set of data, by the values in a rank-2, 3 or 4 container inputDataLeft indexed by (C,P), (C,P,D1), or (C,P,D1,D2) representing the values of scalar, vector or tensor data, OR
  • dot product of a rank-2, 3 or 4 container inputDataRight with dimensions (P), (P,D1) or (P,D1,D2), representing the values of scalar, vector or tensor data, by the values in a rank-2 container inputDataLeft indexed by (C,P), (C,P,D1) or (C,P,D1,D2), representing the values of scalar, vector, or tensor data; the output value container outputData is indexed by (C,P), regardless of which of the two use cases is considered.

For input fields containers without a dimension index, this operation reduces to scalar multiplication.

C - num. integration domains
P - num. integration points
D1 - first spatial (tensor) dimension index
D2 - second spatial (tensor) dimension index
Parameters
outputData[out] - Output (dot product) data array.
inputDataLeft[in] - Left input data array.
inputDataRight[in] - Right input data array.

Definition at line 639 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::dotMultiplyDataData().

template<typename SpT >
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataProperties, typename inputFieldValueType , class... inputFieldProperties>
void Intrepid2::FunctionSpaceTools< SpT >::dotMultiplyDataField ( Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...>  outputFields,
const Kokkos::DynRankView< inputDataValueType, inputDataProperties...>  inputData,
const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...>  inputFields 
)
static

Dot product of data and fields; please read the description below.

There are two use cases:

  • dot product of a rank-3, 4 or 5 container inputFields with dimensions (C,F,P) (C,F,P,D1) or (C,F,P,D1,D2), representing the values of a set of scalar, vector or tensor fields, by the values in a rank-2, 3 or 4 container inputData indexed by (C,P), (C,P,D1), or (C,P,D1,D2) representing the values of scalar, vector or tensor data, OR
  • dot product of a rank-2, 3 or 4 container inputFields with dimensions (F,P), (F,P,D1) or (F,P,D1,D2), representing the values of a scalar, vector or tensor field, by the values in a rank-2 container inputData indexed by (C,P), (C,P,D1) or (C,P,D1,D2), representing the values of scalar, vector or tensor data; the output value container outputFields is indexed by (C,F,P), regardless of which of the two use cases is considered.

For input fields containers without a dimension index, this operation reduces to scalar multiplication.

C - num. integration domains
F - num. fields
P - num. integration points
D1 - first spatial (tensor) dimension index
D2 - second spatial (tensor) dimension index
Parameters
outputFields[out] - Output (dot product) fields array.
inputData[in] - Data array.
inputFields[in] - Input fields array.

Definition at line 623 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::dotMultiplyDataField().

template<typename SpT >
template<typename outputPointValueType , class... outputPointProperties, typename inputCoeffValueType , class... inputCoeffProperties, typename inputFieldValueType , class... inputFieldProperties>
void Intrepid2::FunctionSpaceTools< SpT >::evaluate ( Kokkos::DynRankView< outputPointValueType, outputPointProperties...>  outputPointVals,
const Kokkos::DynRankView< inputCoeffValueType, inputCoeffProperties...>  inputCoeffs,
const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...>  inputFields 
)
static

Computes point values outPointVals of a discrete function specified by the basis inFields and coefficients inCoeffs.

The array inFields with dimensions (C,F,P), (C,F,P,D1), or (C,F,P,D1,D2) represents the signed, transformed field (basis) values at points in REFERENCE frame; the outPointVals array with dimensions (C,P), (C,P,D1), or (C,P,D1,D2), respectively, represents values of a discrete function at points in PHYSICAL frame. The array inCoeffs dimensioned (C,F) supplies the coefficients for the field (basis) array.

Returns rank-2,3 or 4 array such that

\[ outPointValues(c,p,*) = \sum_{f=0}^{F-1} \sigma_{c,f} u_{c,f}(x_p) \]

where $\{u_{c,f}\}_{f=0}^{F-1} $ is scalar, vector or tensor valued finite element basis defined on physical cell $\mathcal{C}$ and $\{\sigma_{c,f}\}_{f=0}^{F-1} $ are the field signs of the basis functions; see Section Pullbacks. This method implements the last step in a four step process; please see Section Evaluation of finite element fields for details about the first three steps that prepare the necessary data for this method.

C - num. integration domains
F - num. fields
P - num. integration points
D1 - spatial dimension
D2 - spatial dimension
Parameters
outputPointVals[out] - Output point values of a discrete function.
inputCoeffs[in] - Coefficients associated with the fields (basis) array.
inputFields[in] - Field (basis) values.

Definition at line 1003 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianValueType , class... jacobianProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HCURLtransformCURL ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianValueType, jacobianProperties...>  jacobian,
const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...>  jacobianDet,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a 3D curl field in the H-curl space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D).

Computes pullback of curls of HCURL functions $\Phi^*(\widehat{\bf u}_f) = \left(J^{-1}_{c} DF_{c}\cdot\nabla\times\widehat{\bf u}_f\right)\circ F^{-1}_{c}$ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the curls of the vector function set $\{\widehat{\bf u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p,*) = \nabla\times\widehat{\bf u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left(J^{-1}_{c} DF_{c}\cdot\nabla\times\widehat{\bf u}_f\right)\circ F^{-1}_{c}(x_{c,p}) = J^{-1}_{c}(\widehat{x}_p) DF_{c}(\widehat{x}_p)\cdot\nabla\times\widehat{\bf u}_f(\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
| D | space dim | 0 <= D < spatial dimension |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobian[in] - Input array containing cell Jacobians.
jacobianDet[in] - Input array containing cell Jacobian determinants.
inputVals[in] - Input array of reference HCURL curls.

Definition at line 189 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HCURLtransformCURL ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...>  jacobianDet,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a 2D curl field in the H-curl space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P).

Computes pullback of curls of HCURL functions $\Phi^*(\widehat{\bf u}_f) = \left(J^{-1}_{c}\nabla\times\widehat{\bf u}_{f}\right) \circ F^{-1}_{c} $ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the 2d curls of the vector function set $\{\widehat{\bf u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p) = \nabla\times\widehat{\bf u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left(J^{-1}_{c}\nabla\times\widehat{\bf u}_{f}\right) \circ F^{-1}_{c} (x_{c,p}) = J^{-1}_{c}(\widehat{x}_p) \nabla\times\widehat{\bf u}_{f} (\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
| D | space dim | 0 <= D < spatial dimension |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobianDet[in] - Input array containing cell Jacobian determinants.
inputVals[in] - Input array of reference HCURL curls.

Definition at line 207 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianInverseValueType , class... jacobianInverseProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HCURLtransformVALUE ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianInverseValueType, jacobianInverseProperties...>  jacobianInverse,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a (vector) value field in the H-curl space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D).

Computes pullback of HCURL functions $\Phi^*(\widehat{\bf u}_f) = \left((DF_c)^{-{\sf T}}\cdot\widehat{\bf u}_f\right)\circ F^{-1}_{c}$ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the values of the vector function set $\{\widehat{\bf u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p,*) = \widehat{\bf u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left((DF_c)^{-{\sf T}}\cdot\widehat{\bf u}_f\right)\circ F^{-1}_{c}(x_{c,p}) = (DF_c)^{-{\sf T}}(\widehat{x}_p)\cdot\widehat{\bf u}_f(\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of native basis |
| P | point | 0 <= P < num. integration points |
| D | space dim | 0 <= D < spatial dimension |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobianInverse[in] - Input array containing cell Jacobian inverses.
inputVals[in] - Input array of reference HCURL values.

Definition at line 174 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::matvecProductDataField().

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HDIVtransformDIV ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...>  jacobianDet,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a divergence field in the H-div space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P).

Computes pullback of the divergence of HDIV functions $\Phi^*(\widehat{\bf u}_f) = \left(J^{-1}_{c}\nabla\cdot\widehat{\bf u}_{f}\right) \circ F^{-1}_{c} $ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the divergencies of the vector function set $\{\widehat{\bf u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p) = \nabla\cdot\widehat{\bf u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left(J^{-1}_{c}\nabla\cdot\widehat{\bf u}_{f}\right) \circ F^{-1}_{c} (x_{c,p}) = J^{-1}_{c}(\widehat{x}_p) \nabla\cdot\widehat{\bf u}_{f} (\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobianDet[in] - Input array containing cell Jacobian determinants.
inputVals[in] - Input array of reference HDIV divergences.

Definition at line 266 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianValueType , class... jacobianProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HDIVtransformVALUE ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianValueType, jacobianProperties...>  jacobian,
const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...>  jacobianDet,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a (vector) value field in the H-div space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D).

Computes pullback of HDIV functions $\Phi^*(\widehat{\bf u}_f) = \left(J^{-1}_{c} DF_{c}\cdot\widehat{\bf u}_f\right)\circ F^{-1}_{c} $ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the values of the vector function set $\{\widehat{\bf u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p,*) = \widehat{\bf u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left(J^{-1}_{c} DF_{c}\cdot \widehat{\bf u}_f\right)\circ F^{-1}_{c}(x_{c,p}) = J^{-1}_{c}(\widehat{x}_p) DF_{c}(\widehat{x}_p)\cdot\widehat{\bf u}_f(\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
| D | space dim | 0 <= D < spatial dimension |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobian[in] - Input array containing cell Jacobians.
jacobianDet[in] - Input array containing cell Jacobian determinants.
inputVals[in] - Input array of reference HDIV values.

Definition at line 250 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::matvecProductDataField(), and Intrepid2::ArrayTools< ExecSpaceType >::scalarMultiplyDataField().

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianValueType , class... jacobianProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HGRADtransformCURL ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianValueType, jacobianProperties...>  jacobian,
const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...>  jacobianDet,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a 2D curl field in the H-grad space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D).

Computes pullback of curls of 2D HGRAD functions $\Phi^*(\widehat{\bf u}_f) = \left(J^{-1}_{c} DF_{c}\cdot\nabla\times\widehat{\bf u}_f\right)\circ F^{-1}_{c}$ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the curls of the vector function set $\{\widehat{\bf u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p,*) = \nabla\times\widehat{\bf u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left(J^{-1}_{c} DF_{c}\cdot\nabla\times\widehat{\bf u}_f\right)\circ F^{-1}_{c}(x_{c,p}) = J^{-1}_{c}(\widehat{x}_p) DF_{c}(\widehat{x}_p)\cdot\nabla\times\widehat{\bf u}_f(\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
| D | space dim | 0 <= D < spatial dimension |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobian[in] - Input array containing cell Jacobians.
jacobianDet[in] - Input array containing cell Jacobian determinants.
inputVals[in] - Input array of reference HDIV values.

Definition at line 228 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianInverseValueType , class... jacobianInverseProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HGRADtransformGRAD ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianInverseValueType, jacobianInverseProperties...>  jacobianInverse,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a gradient field in the H-grad space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P,D), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P,D).

Computes pullback of gradients of HGRAD functions $\Phi^*(\nabla\widehat{u}_f) = \left((DF_c)^{-{\sf T}}\cdot\nabla\widehat{u}_f\right)\circ F^{-1}_{c}$ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the gradients of the function set $\{\widehat{u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p,*) = \nabla\widehat{u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left((DF_c)^{-{\sf T}}\cdot\nabla\widehat{u}_f\right)\circ F^{-1}_{c}(x_{c,p}) = (DF_c)^{-{\sf T}}(\widehat{x}_p)\cdot\nabla\widehat{u}_f(\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
| D | space dim | 0 <= D < spatial dimension |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobianInverse[in] - Input array containing cell Jacobian inverses.
inputVals[in] - Input array of reference HGRAD gradients.

Definition at line 124 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputValueType , class... outputProperties, typename inputValueType , class... inputProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HGRADtransformVALUE ( Kokkos::DynRankView< outputValueType, outputProperties...>  output,
const Kokkos::DynRankView< inputValueType, inputProperties...>  input 
)
static

Transformation of a (scalar) value field in the H-grad space, defined at points on a reference cell, stored in the user-provided container input and indexed by (F,P), into the output container output, defined on cells in physical space and indexed by (C,F,P).

Computes pullback of HGRAD functions $\Phi^*(\widehat{u}_f) = \widehat{u}_f\circ F^{-1}_{c} $ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the values of the function set $\{\widehat{u}_f\}_{f=0}^{F}$ at the reference points:

\[ input(f,p) = \widehat{u}_f(\widehat{x}_p) \,. \]

The method returns

\[ output(c,f,p) = \widehat{u}_f\circ F^{-1}_{c}(x_{c,p}) = \widehat{u}_f(\widehat{x}_p) = input(f,p) \qquad 0\le c < C \,, \]

i.e., it simply replicates the values in the user-provided container to every cell. See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
|------|----------------------|--------------------------------------------------|
Parameters
output[out] - Output array with transformed values
input[in] - Input array of reference HGRAD values.

Definition at line 60 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::RealSpaceTools< ExecSpaceType >::clone(), and Intrepid2::ArrayTools< ExecSpaceType >::cloneFields().

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename jacobianDetValueType , class... jacobianDetProperties, typename inputValValueType , class... inputValProperties>
void Intrepid2::FunctionSpaceTools< SpT >::HVOLtransformVALUE ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< jacobianDetValueType, jacobianDetProperties...>  jacobianDet,
const Kokkos::DynRankView< inputValValueType, inputValProperties...>  inputVals 
)
static

Transformation of a (scalar) value field in the H-vol space, defined at points on a reference cell, stored in the user-provided container inputVals and indexed by (F,P), into the output container outputVals, defined on cells in physical space and indexed by (C,F,P).

Computes pullback of HVOL functions $\Phi^*(\widehat{u}_f) = \left(J^{-1}_{c}\widehat{u}_{f}\right) \circ F^{-1}_{c} $ for points in one or more physical cells that are images of a given set of points in the reference cell:

\[ \{ x_{c,p} \}_{p=0}^P = \{ F_{c} (\widehat{x}_p) \}_{p=0}^{P}\qquad 0\le c < C \,. \]

In this case $ F^{-1}_{c}(x_{c,p}) = \widehat{x}_p $ and the user-provided container should contain the values of the functions in the set $\{\widehat{\bf u}_f\}_{f=0}^{F}$ at the reference points:

\[ inputVals(f,p) = \widehat{u}_f(\widehat{x}_p) \,. \]

The method returns

\[ outputVals(c,f,p,*) = \left(J^{-1}_{c}\widehat{u}_{f}\right) \circ F^{-1}_{c} (x_{c,p}) = J^{-1}_{c}(\widehat{x}_p) \widehat{u}_{f} (\widehat{x}_p) \qquad 0\le c < C \,. \]

See Section Pullbacks for more details about pullbacks.

|------|----------------------|--------------------------------------------------|
| | Index | Dimension |
|------|----------------------|--------------------------------------------------|
| C | cell | 0 <= C < num. integration domains |
| F | field | 0 <= F < dim. of the basis |
| P | point | 0 <= P < num. integration points |
|------|----------------------|--------------------------------------------------|
Parameters
outputVals[out] - Output array with transformed values
jacobianDet[in] - Input array containing cell Jacobian determinants.
inputVals[in] - Input array of reference HVOL values.

Definition at line 280 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::scalarMultiplyDataField().

template<typename SpT >
template<typename outputValueValueType , class... outputValueProperties, typename leftValueValueType , class... leftValueProperties, typename rightValueValueType , class... rightValueProperties>
void Intrepid2::FunctionSpaceTools< SpT >::integrate ( Kokkos::DynRankView< outputValueValueType, outputValueProperties...>  outputValues,
const Kokkos::DynRankView< leftValueValueType, leftValueProperties...>  leftValues,
const Kokkos::DynRankView< rightValueValueType, rightValueProperties...>  rightValues,
const bool  sumInto = false 
)
static

Contracts leftValues and rightValues arrays on the point and possibly space dimensions and stores the result in outputValues; this is a generic, high-level integration routine that calls either FunctionSpaceTools::operatorIntegral, or FunctionSpaceTools::functionalIntegral, or FunctionSpaceTools::dataIntegral methods, depending on the rank of the outputValues array.

Parameters
outputValues[out] - Output array.
leftValues[in] - Left input array.
rightValues[in] - Right input array.
sumInto[in] - If TRUE, sum into given output array, otherwise overwrite it. Default: FALSE.

Definition at line 294 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::contractDataDataScalar(), Intrepid2::ArrayTools< ExecSpaceType >::contractDataDataTensor(), Intrepid2::ArrayTools< ExecSpaceType >::contractDataDataVector(), Intrepid2::ArrayTools< ExecSpaceType >::contractDataFieldScalar(), Intrepid2::ArrayTools< ExecSpaceType >::contractDataFieldTensor(), Intrepid2::ArrayTools< ExecSpaceType >::contractDataFieldVector(), Intrepid2::ArrayTools< ExecSpaceType >::contractFieldFieldScalar(), Intrepid2::ArrayTools< ExecSpaceType >::contractFieldFieldTensor(), and Intrepid2::ArrayTools< ExecSpaceType >::contractFieldFieldVector().

template<typename SpT >
template<typename outputValValueType , class... outputValProperties, typename inputMeasureValueType , class... inputMeasureProperties, typename inputValValueType , class... inputValProperteis>
void Intrepid2::FunctionSpaceTools< SpT >::multiplyMeasure ( Kokkos::DynRankView< outputValValueType, outputValProperties...>  outputVals,
const Kokkos::DynRankView< inputMeasureValueType, inputMeasureProperties...>  inputMeasure,
const Kokkos::DynRankView< inputValValueType, inputValProperteis...>  inputVals 
)
static

Multiplies fields inputVals by weighted measures inputMeasure and returns the field array outputVals; this is a simple redirection to the call FunctionSpaceTools::scalarMultiplyDataField.

Parameters
outputVals[out] - Output array with scaled field values.
inputMeasure[in] - Input array containing weighted measures.
inputVals[in] - Input fields.

Definition at line 571 of file Intrepid2_FunctionSpaceToolsDef.hpp.

template<typename SpT >
template<typename outputDataValuetype , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
void Intrepid2::FunctionSpaceTools< SpT >::scalarMultiplyDataData ( Kokkos::DynRankView< outputDataValuetype, outputDataProperties...>  outputData,
const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...>  inputDataLeft,
const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...>  inputDataRight,
const bool  reciprocal = false 
)
static

Scalar multiplication of data and data; please read the description below.

There are two use cases:

  • multiplies a rank-2, 3, or 4 container inputDataRight with dimensions (C,P), (C,P,D1) or (C,P,D1,D2), representing the values of a set of scalar, vector or tensor data, by the values in a rank-2 container inputDataLeft indexed by (C,P), representing the values of scalar data, OR
  • multiplies a rank-1, 2, or 3 container inputDataRight with dimensions (P), (P,D1) or (P,D1,D2), representing the values of scalar, vector or tensor data, by the values in a rank-2 container inputDataLeft indexed by (C,P), representing the values of scalar data; the output value container outputData is indexed by (C,P), (C,P,D1) or (C,P,D1,D2), regardless of which of the two use cases is considered.
C - num. integration domains
P - num. integration points
D1 - first spatial (tensor) dimension index
D2 - second spatial (tensor) dimension index
Note
The arguments inputDataLeft, inputDataRight can be changed! This enables in-place multiplication.
Parameters
outputData[out] - Output data array.
inputDataLeft[in] - Left (multiplying) data array.
inputDataRight[in] - Right (being multiplied) data array.
reciprocal[in] - If TRUE, divides input fields by the data (instead of multiplying). Default: FALSE.

Definition at line 605 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::scalarMultiplyDataData().

template<typename SpT >
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataProperties, typename inputFieldValueType , class... inputFieldProperties>
void Intrepid2::FunctionSpaceTools< SpT >::scalarMultiplyDataField ( Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...>  outputFields,
const Kokkos::DynRankView< inputDataValueType, inputDataProperties...>  inputData,
const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...>  inputFields,
const bool  reciprocal = false 
)
static

Scalar multiplication of data and fields; please read the description below.

There are two use cases:

  • multiplies a rank-3, 4, or 5 container inputFields with dimensions (C,F,P), (C,F,P,D1) or (C,F,P,D1,D2), representing the values of a set of scalar, vector or tensor fields, by the values in a rank-2 container inputData indexed by (C,P), representing the values of scalar data, OR
  • multiplies a rank-2, 3, or 4 container inputFields with dimensions (F,P), (F,P,D1) or (F,P,D1,D2), representing the values of a scalar, vector or a tensor field, by the values in a rank-2 container inputData indexed by (C,P), representing the values of scalar data; the output value container outputFields is indexed by (C,F,P), (C,F,P,D1) or (C,F,P,D1,D2), regardless of which of the two use cases is considered.
C - num. integration domains
F - num. fields
P - num. integration points
D1 - first spatial (tensor) dimension index
D2 - second spatial (tensor) dimension index
Note
The argument inputFields can be changed! This enables in-place multiplication.
Parameters
outputFields[out] - Output (product) fields array.
inputData[in] - Data (multiplying) array.
inputFields[in] - Input (being multiplied) fields array.
reciprocal[in] - If TRUE, divides input fields by the data (instead of multiplying). Default: FALSE.

Definition at line 587 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::scalarMultiplyDataField().

template<typename SpT >
template<typename outputDataValueType , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
void Intrepid2::FunctionSpaceTools< SpT >::tensorMultiplyDataData ( Kokkos::DynRankView< outputDataValueType, outputDataProperties...>  outputData,
const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...>  inputDataLeft,
const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...>  inputDataRight,
const char  transpose = 'N' 
)
static

Matrix-vector or matrix-matrix product of data and data; please read the description below.

There are four use cases:

  • matrix-vector product of a rank-3 container inputDataRight with dimensions (C,P,D), representing the values of a set of vector data, on the left by the values in a rank-2, 3, or 4 container inputDataLeft indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data, OR
  • matrix-vector product of a rank-2 container inputDataRight with dimensions (P,D), representing the values of vector data, on the left by the values in a rank-2, 3, or 4 container inputDataLeft indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data, OR
  • matrix-matrix product of a rank-4 container inputDataRight with dimensions (C,P,D,D), representing the values of a set of tensor data, on the left by the values in a rank-2, 3, or 4 container inputDataLeft indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data, OR
  • matrix-matrix product of a rank-3 container inputDataRight with dimensions (P,D,D), representing the values of tensor data, on the left by the values in a rank-2, 3, or 4 container inputDataLeft indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data; for matrix-vector products, the output value container outputData is indexed by (C,P,D); for matrix-matrix products, the output value container outputData is indexed by (C,P,D1,D2).
Remarks
The rank of inputDataLeft implicitly defines the type of tensor data:
  • rank = 2 corresponds to a constant diagonal tensor $ diag(a,\ldots,a) $
  • rank = 3 corresponds to a nonconstant diagonal tensor $ diag(a_1,\ldots,a_d) $
  • rank = 4 corresponds to a full tensor $ \{a_{ij}\}$
Note
It is assumed that all tensors are square!
The method is defined for spatial dimensions D = 1, 2, 3
C - num. integration domains
P - num. integration points
D - spatial dimension
Parameters
outputData[out] - Output (matrix-vector product) data array.
inputDataLeft[in] - Left input data array.
inputDataRight[in] - Right input data array.
transpose[in] - If 'T', use transposed tensor; if 'N', no transpose. Default: 'N'.

Definition at line 751 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::matmatProductDataData(), and Intrepid2::ArrayTools< ExecSpaceType >::matvecProductDataData().

template<typename SpT >
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataProperties, typename inputFieldValueType , class... inputFieldProperties>
void Intrepid2::FunctionSpaceTools< SpT >::tensorMultiplyDataField ( Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...>  outputFields,
const Kokkos::DynRankView< inputDataValueType, inputDataProperties...>  inputData,
const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...>  inputFields,
const char  transpose = 'N' 
)
static

Matrix-vector or matrix-matrix product of data and fields; please read the description below.

There are four use cases:

  • matrix-vector product of a rank-4 container inputFields with dimensions (C,F,P,D), representing the values of a set of vector fields, on the left by the values in a rank-2, 3, or 4 container inputData indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data, OR
  • matrix-vector product of a rank-3 container inputFields with dimensions (F,P,D), representing the values of a vector field, on the left by the values in a rank-2, 3, or 4 container inputData indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data, OR
  • matrix-matrix product of a rank-5 container inputFields with dimensions (C,F,P,D,D), representing the values of a set of tensor fields, on the left by the values in a rank-2, 3, or 4 container inputData indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data, OR
  • matrix-matrix product of a rank-4 container inputFields with dimensions (F,P,D,D), representing the values of a tensor field, on the left by the values in a rank-2, 3, or 4 container inputData indexed by (C,P), (C,P,D) or (C,P,D,D), respectively, representing the values of tensor data; for matrix-vector products, the output value container outputFields is indexed by (C,F,P,D); for matrix-matrix products the output value container outputFields is indexed by (C,F,P,D,D).
Remarks
The rank of inputData implicitly defines the type of tensor data:
  • rank = 2 corresponds to a constant diagonal tensor $ diag(a,\ldots,a) $
  • rank = 3 corresponds to a nonconstant diagonal tensor $ diag(a_1,\ldots,a_d) $
  • rank = 4 corresponds to a full tensor $ \{a_{ij}\}$
Note
It is assumed that all tensors are square!
The method is defined for spatial dimensions D = 1, 2, 3
C - num. integration domains
F - num. fields
P - num. integration points
D - spatial dimension
Parameters
outputFields[out] - Output (matrix-vector or matrix-matrix product) fields array.
inputData[in] - Data array.
inputFields[in] - Input fields array.
transpose[in] - If 'T', use transposed left data tensor; if 'N', no transpose. Default: 'N'.

Definition at line 717 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::matmatProductDataField(), and Intrepid2::ArrayTools< ExecSpaceType >::matvecProductDataField().

template<typename SpT >
template<typename outputDataValueType , class... outputDataProperties, typename inputDataLeftValueType , class... inputDataLeftProperties, typename inputDataRightValueType , class... inputDataRightProperties>
void Intrepid2::FunctionSpaceTools< SpT >::vectorMultiplyDataData ( Kokkos::DynRankView< outputDataValueType, outputDataProperties...>  outputData,
const Kokkos::DynRankView< inputDataLeftValueType, inputDataLeftProperties...>  inputDataLeft,
const Kokkos::DynRankView< inputDataRightValueType, inputDataRightProperties...>  inputDataRight 
)
static

Cross or outer product of data and data; please read the description below.

There are four use cases:

  • cross product of a rank-3 container inputDataRight with dimensions (C,P,D), representing the values of a set of vector data, on the left by the values in a rank-3 container inputDataLeft indexed by (C,P,D) representing the values of vector data, OR
  • cross product of a rank-2 container inputDataRight with dimensions (P,D), representing the values of vector data, on the left by the values in a rank-3 container inputDataLeft indexed by (C,P,D), representing the values of vector data, OR
  • outer product of a rank-3 container inputDataRight with dimensions (C,P,D), representing the values of a set of vector data, on the left by the values in a rank-3 container inputDataLeft indexed by (C,P,D) representing the values of vector data, OR
  • outer product of a rank-2 container inputDataRight with dimensions (P,D), representing the values of vector data, on the left by the values in a rank-3 container inputDataLeft indexed by (C,P,D), representing the values of vector data; for cross products, the output value container outputData is indexed by (C,P,D) in 3D (vector output) and by (C,P) in 2D (scalar output); for outer products, the output value container outputData is indexed by (C,P,D,D).
C - num. integration domains
P - num. integration points
D - spatial dimension, must be 2 or 3
Parameters
outputData[out] - Output (cross or outer product) data array.
inputDataLeft[in] - Left input data array.
inputDataRight[in] - Right input data array.

Definition at line 686 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::crossProductDataData(), and Intrepid2::ArrayTools< ExecSpaceType >::outerProductDataData().

template<typename SpT >
template<typename outputFieldValueType , class... outputFieldProperties, typename inputDataValueType , class... inputDataProperties, typename inputFieldValueType , class... inputFieldProperties>
void Intrepid2::FunctionSpaceTools< SpT >::vectorMultiplyDataField ( Kokkos::DynRankView< outputFieldValueType, outputFieldProperties...>  outputFields,
const Kokkos::DynRankView< inputDataValueType, inputDataProperties...>  inputData,
const Kokkos::DynRankView< inputFieldValueType, inputFieldProperties...>  inputFields 
)
static

Cross or outer product of data and fields; please read the description below.

There are four use cases:

  • cross product of a rank-4 container inputFields with dimensions (C,F,P,D), representing the values of a set of vector fields, on the left by the values in a rank-3 container inputData indexed by (C,P,D), representing the values of vector data, OR
  • cross product of a rank-3 container inputFields with dimensions (F,P,D), representing the values of a vector field, on the left by the values in a rank-3 container inputData indexed by (C,P,D), representing the values of vector data, OR
  • outer product of a rank-4 container inputFields with dimensions (C,F,P,D), representing the values of a set of vector fields, on the left by the values in a rank-3 container inputData indexed by (C,P,D), representing the values of vector data, OR
  • outer product of a rank-3 container inputFields with dimensions (F,P,D), representing the values of a vector field, on the left by the values in a rank-3 container inputData indexed by (C,P,D), representing the values of vector data; for cross products, the output value container outputFields is indexed by (C,F,P,D) in 3D (vector output) and by (C,F,P) in 2D (scalar output); for outer products, the output value container outputFields is indexed by (C,F,P,D,D).
C - num. integration domains
F - num. fields
P - num. integration points
D - spatial dimension, must be 2 or 3
Parameters
outputFields[out] - Output (cross or outer product) fields array.
inputData[in] - Data array.
inputFields[in] - Input fields array.

Definition at line 655 of file Intrepid2_FunctionSpaceToolsDef.hpp.

References Intrepid2::ArrayTools< ExecSpaceType >::crossProductDataField(), and Intrepid2::ArrayTools< ExecSpaceType >::outerProductDataField().


The documentation for this class was generated from the following files: