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

A class providing static members to perform Lagrangian interpolation on a finite element. More...

Static Public Member Functions

template<typename basisCoeffsViewType , typename funcViewType , typename BasisType , typename ortViewType >
static void getBasisCoeffs (basisCoeffsViewType basisCoeffs, const funcViewType functionAtDofCoords, const BasisType *cellBasis, const ortViewType orts)
 Computes the basis weights of the function interpolation. More...
 

Detailed Description

template<typename DeviceType>
class Intrepid2::LagrangianInterpolation< DeviceType >

A class providing static members to perform Lagrangian interpolation on a finite element.

classes The Lagrangian interpolant $\mathcal I_h$ is defined as

\[ \mathcal I_h(f) = \sum_i \alpha^f_i \phi_i, \quad \alpha_i^f := L_i(f), \]

where $\{\phi_i\}$ is the basis of the finite element, $\alpha_i^f$ are the basisCoeffs and $L_i$ are the DOFs of the basis. This class works only for Lagrangian elements that provide a set of orthonormal DOFs defined as

\[ L_i(f) := f(\mathbf x_i) \cdot \beta_i, \quad L_i(\phi_j) = \delta_{ij}, \]

where $\beta_i$ are referred to as dofCoeffs, and $\mathbf x_i$ are the coordinates of the basis nodes.

In order to perform the interpolation, evaluate the function $f$ at the set of points $\{\mathbf x_i\}$ computed using the basis method getDoCoords and then obtain the basis coefficients $\alpha_i^f$ by calling the function getBasisCoeffs.

Remarks
The interpolation is performed at the reference element. Therefore, the function $f$, which is contravariant, needs to mapped to the reference space, using the inverse operation of a pullback, before calling getBasisCoeffs.

Definition at line 130 of file Intrepid2_LagrangianInterpolation.hpp.

Member Function Documentation

template<typename DeviceType >
template<typename basisCoeffsViewType , typename funcViewType , typename BasisType , typename ortViewType >
void Intrepid2::LagrangianInterpolation< DeviceType >::getBasisCoeffs ( basisCoeffsViewType  basisCoeffs,
const funcViewType  functionAtDofCoords,
const BasisType *  cellBasis,
const ortViewType  orts 
)
static

Computes the basis weights of the function interpolation.

C - num. cells
F - num. fields
D - spatial dimension
Parameters
basisCoeffs[out] - rank-2 view (C,F) that will contain the basis coefficients of the interpolation.
functionAtDofCoords[in] - variable rank view that contains the function evaluated at reference DOF coordinates.
cellBasis[in] - pointer to the basis for the interpolation
cellOrientations[in] - rank-1 view (C) containing the Orientation objects at each cell
Remarks
The output views need to be pre-allocated. dofCoeffs and functionAtDofCoords have rank 2, (C,F) for scalar basis and 3, (C,F,D) for vector basis. functValsAtDofCoords contains the function evaluated at the reference dofCoords and contravariantly transformed to the reference element. The reference dofCoords are obtained with the method cellBasis->getDofCoords(...), NOT with the getOrientedDofCoords() method

Definition at line 342 of file Intrepid2_LagrangianInterpolationDef.hpp.

References Intrepid2::RealSpaceTools< DeviceType >::clone(), Intrepid2::ArrayTools< DeviceType >::dotMultiplyDataData(), and Intrepid2::OrientationTools< DeviceType >::modifyBasisByOrientationInverse().


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