Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Public Member Functions | Protected Member Functions | Private Member Functions | List of all members
Nektar::StdRegions::StdPyrExp Class Reference

#include <StdPyrExp.h>

Inheritance diagram for Nektar::StdRegions::StdPyrExp:
[legend]

Public Member Functions

 StdPyrExp ()
 
 StdPyrExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc)
 
 StdPyrExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc, NekDouble *coeffs, NekDouble *phys)
 
 StdPyrExp (const StdPyrExp &T)
 
 ~StdPyrExp ()
 
- Public Member Functions inherited from Nektar::StdRegions::StdExpansion3D
 StdExpansion3D ()
 
 StdExpansion3D (int numcoeffs, const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc)
 
 StdExpansion3D (const StdExpansion3D &T)
 
virtual ~StdExpansion3D ()
 
void PhysTensorDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray_d1, Array< OneD, NekDouble > &outarray_d2, Array< OneD, NekDouble > &outarray_d3)
 Calculate the 3D derivative in the local tensor/collapsed coordinate at the physical points. More...
 
void BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
 
void IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
 
int GetNedges () const
 return the number of edges in 3D expansion More...
 
int GetEdgeNcoeffs (const int i) const
 This function returns the number of expansion coefficients belonging to the i-th edge. More...
 
void GetEdgeInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards)
 
- Public Member Functions inherited from Nektar::StdRegions::StdExpansion
 StdExpansion ()
 Default Constructor. More...
 
 StdExpansion (const int numcoeffs, const int numbases, const LibUtilities::BasisKey &Ba=LibUtilities::NullBasisKey, const LibUtilities::BasisKey &Bb=LibUtilities::NullBasisKey, const LibUtilities::BasisKey &Bc=LibUtilities::NullBasisKey)
 Constructor. More...
 
 StdExpansion (const StdExpansion &T)
 Copy Constructor. More...
 
virtual ~StdExpansion ()
 Destructor. More...
 
int GetNumBases () const
 This function returns the number of 1D bases used in the expansion. More...
 
const Array< OneD, const LibUtilities::BasisSharedPtr > & GetBase () const
 This function gets the shared point to basis. More...
 
const LibUtilities::BasisSharedPtrGetBasis (int dir) const
 This function gets the shared point to basis in the dir direction. More...
 
int GetNcoeffs (void) const
 This function returns the total number of coefficients used in the expansion. More...
 
int GetTotPoints () const
 This function returns the total number of quadrature points used in the element. More...
 
LibUtilities::BasisType GetBasisType (const int dir) const
 This function returns the type of basis used in the dir direction. More...
 
int GetBasisNumModes (const int dir) const
 This function returns the number of expansion modes in the dir direction. More...
 
int EvalBasisNumModesMax (void) const
 This function returns the maximum number of expansion modes over all local directions. More...
 
LibUtilities::PointsType GetPointsType (const int dir) const
 This function returns the type of quadrature points used in the dir direction. More...
 
int GetNumPoints (const int dir) const
 This function returns the number of quadrature points in the dir direction. More...
 
const Array< OneD, const NekDouble > & GetPoints (const int dir) const
 This function returns a pointer to the array containing the quadrature points in dir direction. More...
 
int GetNverts () const
 This function returns the number of vertices of the expansion domain. More...
 
int GetTraceNcoeffs (const int i) const
 This function returns the number of expansion coefficients belonging to the i-th trace. More...
 
int GetTraceIntNcoeffs (const int i) const
 
int GetTraceNumPoints (const int i) const
 This function returns the number of quadrature points belonging to the i-th trace. More...
 
const LibUtilities::BasisKey GetTraceBasisKey (const int i, int k=-1) const
 This function returns the basis key belonging to the i-th trace. More...
 
LibUtilities::PointsKey GetTracePointsKey (const int i, int k=-1) const
 This function returns the basis key belonging to the i-th trace. More...
 
int NumBndryCoeffs (void) const
 
int NumDGBndryCoeffs (void) const
 
const LibUtilities::PointsKey GetNodalPointsKey () const
 This function returns the type of expansion Nodal point type if defined. More...
 
int GetNtraces () const
 Returns the number of trace elements connected to this element. More...
 
LibUtilities::ShapeType DetShapeType () const
 This function returns the shape of the expansion domain. More...
 
std::shared_ptr< StdExpansionGetStdExp (void) const
 
std::shared_ptr< StdExpansionGetLinStdExp (void) const
 
int GetShapeDimension () const
 
bool IsBoundaryInteriorExpansion ()
 
bool IsNodalNonTensorialExp ()
 
void BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs the Backward transformation from coefficient space to physical space. More...
 
void FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs the Forward transformation from physical space to coefficient space. More...
 
void FwdTrans_BndConstrained (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
NekDouble Integral (const Array< OneD, const NekDouble > &inarray)
 This function integrates the specified function over the domain. More...
 
void FillMode (const int mode, Array< OneD, NekDouble > &outarray)
 This function fills the array outarray with the mode-th mode of the expansion. More...
 
void IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 this function calculates the inner product of a given function f with the different modes of the expansion More...
 
void IProductWRTBase (const Array< OneD, const NekDouble > &base, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, int coll_check)
 
void IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDirectionalDerivBase (const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
int GetElmtId ()
 Get the element id of this expansion when used in a list by returning value of m_elmt_id. More...
 
void SetElmtId (const int id)
 Set the element id of this expansion when used in a list by returning value of m_elmt_id. More...
 
void GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2=NullNekDouble1DArray, Array< OneD, NekDouble > &coords_3=NullNekDouble1DArray)
 this function returns the physical coordinates of the quadrature points of the expansion More...
 
void GetCoord (const Array< OneD, const NekDouble > &Lcoord, Array< OneD, NekDouble > &coord)
 given the coordinates of a point of the element in the local collapsed coordinate system, this function calculates the physical coordinates of the point More...
 
DNekMatSharedPtr GetStdMatrix (const StdMatrixKey &mkey)
 
DNekBlkMatSharedPtr GetStdStaticCondMatrix (const StdMatrixKey &mkey)
 
void NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
 
void NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
void NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, const Array< OneD, const NekDouble > &Fz, Array< OneD, NekDouble > &outarray)
 
void NormVectorIProductWRTBase (const Array< OneD, const Array< OneD, NekDouble > > &Fvec, Array< OneD, NekDouble > &outarray)
 
DNekScalBlkMatSharedPtr GetLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
void DropLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
int CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset)
 
NekDouble StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
 
int GetCoordim ()
 
void GetBoundaryMap (Array< OneD, unsigned int > &outarray)
 
void GetInteriorMap (Array< OneD, unsigned int > &outarray)
 
int GetVertexMap (const int localVertexId, bool useCoeffPacking=false)
 
void GetTraceToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
 
void GetTraceInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eForwards)
 
void GetTraceNumModes (const int tid, int &numModes0, int &numModes1, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2)
 
void MultiplyByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
DNekMatSharedPtr CreateGeneralMatrix (const StdMatrixKey &mkey)
 this function generates the mass matrix \(\mathbf{M}[i][j] = \int \phi_i(\mathbf{x}) \phi_j(\mathbf{x}) d\mathbf{x}\) More...
 
void GeneralMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void ReduceOrderCoeffs (int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdMatrixKey &mkey)
 
void ExponentialFilter (Array< OneD, NekDouble > &array, const NekDouble alpha, const NekDouble exponent, const NekDouble cutoff)
 
void LaplacianMatrixOp (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDerivMatrixOp (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDirectionalDerivMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassLevelCurvatureMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionDiffusionReactionMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
 
void HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
DNekMatSharedPtr GenMatrix (const StdMatrixKey &mkey)
 
void PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
 
void PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void PhysDeriv_s (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_ds)
 
void PhysDeriv_n (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_dn)
 
void PhysDirectionalDeriv (const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &direction, Array< OneD, NekDouble > &outarray)
 
void StdPhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
 
void StdPhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
NekDouble PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble PhysEvaluate (const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals)
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble PhysEvaluateBasis (const Array< OneD, const NekDouble > &coords, int mode)
 This function evaluates the basis function mode mode at a point coords of the domain. More...
 
void LocCoordToLocCollapsed (const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta)
 Convert local cartesian coordinate xi into local collapsed coordinates eta. More...
 
void LocCollapsedToLocCoord (const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi)
 Convert local collapsed coordinates eta into local cartesian coordinate xi. More...
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, const Array< OneD, const NekDouble > &Fz, Array< OneD, NekDouble > &outarray)
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const Array< OneD, NekDouble > > &Fvec, Array< OneD, NekDouble > &outarray)
 
virtual DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual void v_DropLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
NekDouble Linf (const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
 Function to evaluate the discrete \( L_\infty\) error \( |\epsilon|_\infty = \max |u - u_{exact}|\) where \( u_{exact}\) is given by the array sol. More...
 
NekDouble L2 (const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
 Function to evaluate the discrete \( L_2\) error, \( | \epsilon |_{2} = \left [ \int^1_{-1} [u - u_{exact}]^2 dx \right]^{1/2} d\xi_1 \) where \( u_{exact}\) is given by the array sol. More...
 
NekDouble H1 (const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
 Function to evaluate the discrete \( H^1\) error, \( | \epsilon |^1_{2} = \left [ \int^1_{-1} [u - u_{exact}]^2 + \nabla(u - u_{exact})\cdot\nabla(u - u_{exact})\cdot dx \right]^{1/2} d\xi_1 \) where \( u_{exact}\) is given by the array sol. More...
 
const LibUtilities::PointsKeyVector GetPointsKeys () const
 
DNekMatSharedPtr BuildInverseTransformationMatrix (const DNekScalMatSharedPtr &m_transformationmatrix)
 
void PhysInterpToSimplexEquiSpaced (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, int npset=-1)
 This function performs an interpolation from the physical space points provided at input into an array of equispaced points which are not the collapsed coordinate. So for a tetrahedron you will only get a tetrahedral number of values. More...
 
void GetSimplexEquiSpacedConnectivity (Array< OneD, int > &conn, bool standard=true)
 This function provides the connectivity of local simplices (triangles or tets) to connect the equispaced data points provided by PhysInterpToSimplexEquiSpaced. More...
 
void EquiSpacedToCoeffs (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs a projection/interpolation from the equispaced points sometimes used in post-processing onto the coefficient space. More...
 
template<class T >
std::shared_ptr< T > as ()
 
void IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
 
void GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat)
 

Protected Member Functions

virtual void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
 Calculate the derivative of the physical points. More...
 
virtual void v_PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 Calculate the derivative of the physical points in a given direction. More...
 
virtual void v_StdPhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
 
virtual void v_StdPhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 Backward transformation is evaluated at the quadrature points. More...
 
virtual void v_BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
 
virtual void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 Forward transform from physical quadrature space stored in inarray and evaluate the expansion coefficients and store in outarray. More...
 
virtual void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 Inner product of inarray over region with respect to the expansion basis m_base[0]->GetBdata(),m_base[1]->GetBdata(), m_base[2]->GetBdata() and return in outarray. More...
 
virtual void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
 
virtual void v_IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
 
virtual void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_LocCoordToLocCollapsed (const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta)
 
virtual void v_LocCollapsedToLocCoord (const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi)
 
virtual void v_GetCoords (Array< OneD, NekDouble > &xi_x, Array< OneD, NekDouble > &xi_y, Array< OneD, NekDouble > &xi_z)
 
virtual void v_FillMode (const int mode, Array< OneD, NekDouble > &outarray)
 
virtual void v_GetTraceNumModes (const int fid, int &numModes0, int &numModes1, Orientation faceOrient=eDir1FwdDir1_Dir2FwdDir2)
 
NekDouble v_PhysEvaluateBasis (const Array< OneD, const NekDouble > &coords, int mode) final
 
virtual int v_GetNverts () const
 
virtual int v_GetNedges () const
 
virtual int v_GetNtraces () const
 
virtual LibUtilities::ShapeType v_DetShapeType () const
 
virtual int v_NumBndryCoeffs () const
 
virtual int v_NumDGBndryCoeffs () const
 
virtual int v_GetTraceNcoeffs (const int i) const
 
virtual int v_GetTraceIntNcoeffs (const int i) const
 
virtual int v_GetTraceNumPoints (const int i) const
 
virtual int v_GetEdgeNcoeffs (const int i) const
 
virtual int v_CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset)
 
virtual const LibUtilities::BasisKey v_GetTraceBasisKey (const int i, const int k) const
 
virtual int v_GetVertexMap (int localVertexId, bool useCoeffPacking=false)
 
virtual void v_GetInteriorMap (Array< OneD, unsigned int > &outarray)
 
virtual void v_GetBoundaryMap (Array< OneD, unsigned int > &outarray)
 
virtual void v_GetTraceToElementMap (const int fid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation faceOrient, int P, int Q)
 
virtual void v_GetEdgeInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2)
 
virtual void v_GetTraceInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2)
 
virtual DNekMatSharedPtr v_GenMatrix (const StdMatrixKey &mkey)
 
virtual DNekMatSharedPtr v_CreateStdMatrix (const StdMatrixKey &mkey)
 
virtual void v_SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdMatrixKey &mkey)
 
virtual void v_ReduceOrderCoeffs (int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion3D
virtual NekDouble v_PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
virtual NekDouble v_PhysEvaluate (const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals)
 
virtual void v_LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
 
virtual void v_HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
 
virtual NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray)
 Integrates the specified function over the domain. More...
 
virtual void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat)
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion
DNekMatSharedPtr CreateStdMatrix (const StdMatrixKey &mkey)
 
DNekBlkMatSharedPtr CreateStdStaticCondMatrix (const StdMatrixKey &mkey)
 Create the static condensation of a matrix when using a boundary interior decomposition. More...
 
void BwdTrans_MatOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDirectionalDerivBase_SumFac (const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void GeneralMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree_Kernel (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp)
 
void LaplacianMatrixOp_MatFree_GenericImpl (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDerivMatrixOp_MatFree (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDirectionalDerivMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassLevelCurvatureMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionDiffusionReactionMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
 
void HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void HelmholtzMatrixOp_MatFree_GenericImpl (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
virtual void v_SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_SetCoeffsToOrientation (Array< OneD, NekDouble > &coeffs, StdRegions::Orientation dir)
 
virtual NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
 
template<int DIR>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals)
 This function performs the barycentric interpolation of the polynomial stored in coord at a point physvals using barycentric interpolation weights in direction. More...
 
template<int DIR>
NekDouble BaryEvaluateBasis (const NekDouble &coord, const int &mode)
 

Private Member Functions

int GetMode (int I, int J, int K)
 Compute the mode number in the expansion for a particular tensorial combination. More...
 

Additional Inherited Members

- Protected Attributes inherited from Nektar::StdRegions::StdExpansion
Array< OneD, LibUtilities::BasisSharedPtrm_base
 
int m_elmt_id
 
int m_ncoeffs
 
LibUtilities::NekManager< StdMatrixKey, DNekMat, StdMatrixKey::opLessm_stdMatrixManager
 
LibUtilities::NekManager< StdMatrixKey, DNekBlkMat, StdMatrixKey::opLessm_stdStaticCondMatrixManager
 

Detailed Description

Definition at line 81 of file StdPyrExp.h.

Constructor & Destructor Documentation

◆ StdPyrExp() [1/4]

Nektar::StdRegions::StdPyrExp::StdPyrExp ( )

Definition at line 47 of file StdPyrExp.cpp.

48  {
49  }

◆ StdPyrExp() [2/4]

Nektar::StdRegions::StdPyrExp::StdPyrExp ( const LibUtilities::BasisKey Ba,
const LibUtilities::BasisKey Bb,
const LibUtilities::BasisKey Bc 
)

Definition at line 51 of file StdPyrExp.cpp.

55  Ba.GetNumModes(),
56  Bb.GetNumModes(),
57  Bc.GetNumModes()),
58  3, Ba, Bb, Bc),
60  Ba.GetNumModes(),
61  Bb.GetNumModes(),
62  Bc.GetNumModes()),
63  Ba, Bb, Bc)
64  {
65 
66  ASSERTL0(Ba.GetNumModes() <= Bc.GetNumModes(), "order in 'a' direction is higher "
67  "than order in 'c' direction");
68  ASSERTL0(Bb.GetNumModes() <= Bc.GetNumModes(), "order in 'b' direction is higher "
69  "than order in 'c' direction");
70  ASSERTL1(Bc.GetBasisType() == LibUtilities::eModifiedPyr_C ||
71  Bc.GetBasisType() == LibUtilities::eOrthoPyr_C,
72  "Expected basis type in 'c' direction to be ModifiedPyr_C or OrthoPyr_C");
73 
74  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:216
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:250
StdExpansion()
Default Constructor.
int getNumberOfCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:240
@ eModifiedPyr_C
Principle Modified Functions .
Definition: BasisType.h:52
@ eOrthoPyr_C
Principle Orthogonal Functions .
Definition: BasisType.h:51

References ASSERTL0, ASSERTL1, Nektar::LibUtilities::eModifiedPyr_C, Nektar::LibUtilities::eOrthoPyr_C, Nektar::LibUtilities::BasisKey::GetBasisType(), and Nektar::LibUtilities::BasisKey::GetNumModes().

◆ StdPyrExp() [3/4]

Nektar::StdRegions::StdPyrExp::StdPyrExp ( const LibUtilities::BasisKey Ba,
const LibUtilities::BasisKey Bb,
const LibUtilities::BasisKey Bc,
NekDouble coeffs,
NekDouble phys 
)

◆ StdPyrExp() [4/4]

Nektar::StdRegions::StdPyrExp::StdPyrExp ( const StdPyrExp T)

Definition at line 76 of file StdPyrExp.cpp.

77  : StdExpansion (T),
79  {
80  }

◆ ~StdPyrExp()

Nektar::StdRegions::StdPyrExp::~StdPyrExp ( )

Definition at line 84 of file StdPyrExp.cpp.

85  {
86  }

Member Function Documentation

◆ GetMode()

int Nektar::StdRegions::StdPyrExp::GetMode ( int  I,
int  J,
int  K 
)
private

Compute the mode number in the expansion for a particular tensorial combination.

Modes are numbered with the r index travelling fastest, followed by q and then p, and each q-r plane is of size

(R+1-p)*(Q+1) - l(l+1)/2 where l = max(0,Q-p)

For example, when P=2, Q=3 and R=4 the indexing inside each q-r plane (with r increasing upwards and q to the right) is:

p = 0: p = 1: p = 2:

4 3 8 17 21 2 7 11 16 20 24 29 32 35 1 6 10 13 15 19 23 26 28 31 34 37 0 5 9 12 14 18 22 25 27 30 33 36

Note that in this element, we must have that \( P,Q \leq R\).

Definition at line 1895 of file StdPyrExp.cpp.

1896  {
1897  const int Q = m_base[1]->GetNumModes()-1;
1898  const int R = m_base[2]->GetNumModes()-1;
1899 
1900  int i,l;
1901  int cnt = 0;
1902 
1903  // Traverse to q-r plane number I
1904  for (i = 0; i < I; ++i)
1905  {
1906  // Size of triangle part
1907  l = max(0,Q-i);
1908 
1909  // Size of rectangle part
1910  cnt += (R+1-i)*(Q+1) - l*(l+1)/2;
1911  }
1912 
1913  // Traverse to q column J (Pretend this is a face of width J)
1914  l = max(0,J-1-I);
1915  cnt += (R+1-I)*J - l*(l+1)/2;
1916 
1917  // Traverse up stacks to K
1918  cnt += K;
1919 
1920  return cnt;
1921  }
Array< OneD, LibUtilities::BasisSharedPtr > m_base

References Nektar::StdRegions::StdExpansion::m_base.

Referenced by v_GetBoundaryMap(), v_GetEdgeInteriorToElementMap(), v_GetInteriorMap(), v_GetTraceInteriorToElementMap(), v_GetTraceToElementMap(), and v_GetVertexMap().

◆ v_BwdTrans()

void Nektar::StdRegions::StdPyrExp::v_BwdTrans ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Backward transformation is evaluated at the quadrature points.

\( u^{\delta} (\xi_{1i}, \xi_{2j}, \xi_{3k}) = \sum_{m(pqr)} \hat u_{pqr} \phi_{pqr} (\xi_{1i}, \xi_{2j}, \xi_{3k})\)

Backward transformation is three dimensional tensorial expansion

\( u (\xi_{1i}, \xi_{2j}, \xi_{3k}) = \sum_{p=0}^{Q_x} \psi_p^a (\xi_{1i}) \lbrace { \sum_{q=0}^{Q_y} \psi_{q}^a (\xi_{2j}) \lbrace { \sum_{r=0}^{Q_z} \hat u_{pqr} \psi_{pqr}^c (\xi_{3k}) \rbrace} \rbrace}. \)

And sumfactorizing step of the form is as:\ \ \( f_{pqr} (\xi_{3k}) = \sum_{r=0}^{Q_z} \hat u_{pqr} \psi_{pqr}^c (\xi_{3k}),\\ g_{p} (\xi_{2j}, \xi_{3k}) = \sum_{r=0}^{Q_y} \psi_{p}^a (\xi_{2j}) f_{pqr} (\xi_{3k}),\\ u(\xi_{1i}, \xi_{2j}, \xi_{3k}) = \sum_{p=0}^{Q_x} \psi_{p}^a (\xi_{1i}) g_{p} (\xi_{2j}, \xi_{3k}). \)

Implements Nektar::StdRegions::StdExpansion.

Definition at line 252 of file StdPyrExp.cpp.

254  {
255  if (m_base[0]->Collocation() &&
256  m_base[1]->Collocation() &&
257  m_base[2]->Collocation())
258  {
260  m_base[1]->GetNumPoints()*
261  m_base[2]->GetNumPoints(),
262  inarray, 1, outarray, 1);
263  }
264  else
265  {
266  StdPyrExp::v_BwdTrans_SumFac(inarray,outarray);
267  }
268  }
int GetNumPoints(const int dir) const
This function returns the number of quadrature points in the dir direction.
Definition: StdExpansion.h:221
virtual void v_BwdTrans_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdPyrExp.cpp:273
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1199

References Nektar::StdRegions::StdExpansion::GetNumPoints(), Nektar::StdRegions::StdExpansion::m_base, v_BwdTrans_SumFac(), and Vmath::Vcopy().

Referenced by v_FillMode().

◆ v_BwdTrans_SumFac()

void Nektar::StdRegions::StdPyrExp::v_BwdTrans_SumFac ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Sum-factorisation implementation of the BwdTrans operation.

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 273 of file StdPyrExp.cpp.

276  {
277  int nquad0 = m_base[0]->GetNumPoints();
278  int nquad1 = m_base[1]->GetNumPoints();
279  int nquad2 = m_base[2]->GetNumPoints();
280  int order0 = m_base[0]->GetNumModes();
281  int order1 = m_base[1]->GetNumModes();
282 
283  Array<OneD, NekDouble> wsp(nquad2*order0*order1+
284  nquad2*nquad1*nquad0);
285 
286  v_BwdTrans_SumFacKernel(m_base[0]->GetBdata(),
287  m_base[1]->GetBdata(),
288  m_base[2]->GetBdata(),
289  inarray,outarray,wsp,
290  true,true,true);
291  }
virtual void v_BwdTrans_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
Definition: StdPyrExp.cpp:293

References Nektar::StdRegions::StdExpansion::m_base, and v_BwdTrans_SumFacKernel().

Referenced by v_BwdTrans().

◆ v_BwdTrans_SumFacKernel()

void Nektar::StdRegions::StdPyrExp::v_BwdTrans_SumFacKernel ( const Array< OneD, const NekDouble > &  base0,
const Array< OneD, const NekDouble > &  base1,
const Array< OneD, const NekDouble > &  base2,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
Array< OneD, NekDouble > &  wsp,
bool  doCheckCollDir0,
bool  doCheckCollDir1,
bool  doCheckCollDir2 
)
protectedvirtual

Implements Nektar::StdRegions::StdExpansion3D.

Definition at line 293 of file StdPyrExp.cpp.

303  {
304  boost::ignore_unused(doCheckCollDir0, doCheckCollDir1,
305  doCheckCollDir2);
306 
307  int nquad0 = m_base[0]->GetNumPoints();
308  int nquad1 = m_base[1]->GetNumPoints();
309  int nquad2 = m_base[2]->GetNumPoints();
310 
311  int order0 = m_base[0]->GetNumModes();
312  int order1 = m_base[1]->GetNumModes();
313  int order2 = m_base[2]->GetNumModes();
314 
315  Array<OneD, NekDouble > tmp = wsp;
316  Array<OneD, NekDouble > tmp1 = tmp + nquad2*order0*order1;
317 
318  int i, j, mode, mode1, cnt;
319 
320  // Perform summation over '2' direction
321  mode = mode1 = cnt = 0;
322  for(i = 0; i < order0; ++i)
323  {
324  for(j = 0; j < order1; ++j, ++cnt)
325  {
326  int ijmax = max(i,j);
327  Blas::Dgemv('N', nquad2, order2-ijmax,
328  1.0, base2.get()+mode*nquad2, nquad2,
329  inarray.get()+mode1, 1,
330  0.0, tmp.get()+cnt*nquad2, 1);
331  mode += order2-ijmax;
332  mode1 += order2-ijmax;
333  }
334  //increment mode in case order1!=order2
335  for(j = order1; j < order2-i; ++j)
336  {
337  int ijmax = max(i,j);
338  mode += order2-ijmax;
339  }
340  }
341 
342  // fix for modified basis by adding split of top singular
343  // vertex mode - currently (1+c)/2 x (1-b)/2 x (1-a)/2
344  // component is evaluated
346  {
347 
348  // Not sure why we could not use basis as 1.0
349  // top singular vertex - (1+c)/2 x (1+b)/2 x (1-a)/2 component
350  Blas::Daxpy(nquad2,inarray[1],base2.get()+nquad2,1,
351  &tmp[0]+nquad2,1);
352 
353  // top singular vertex - (1+c)/2 x (1-b)/2 x (1+a)/2 component
354  Blas::Daxpy(nquad2,inarray[1],base2.get()+nquad2,1,
355  &tmp[0]+order1*nquad2,1);
356 
357  // top singular vertex - (1+c)/2 x (1+b)/2 x (1+a)/2 component
358  Blas::Daxpy(nquad2,inarray[1],base2.get()+nquad2,1,
359  &tmp[0]+order1*nquad2+nquad2,1);
360  }
361 
362  // Perform summation over '1' direction
363  mode = 0;
364  for(i = 0; i < order0; ++i)
365  {
366  Blas::Dgemm('N', 'T', nquad1, nquad2, order1,
367  1.0, base1.get(), nquad1,
368  tmp.get()+mode*nquad2, nquad2,
369  0.0, tmp1.get()+i*nquad1*nquad2, nquad1);
370  mode += order1;
371  }
372 
373  // Perform summation over '0' direction
374  Blas::Dgemm('N', 'T', nquad0, nquad1*nquad2, order0,
375  1.0, base0.get(), nquad0,
376  tmp1.get(), nquad1*nquad2,
377  0.0, outarray.get(), nquad0);
378  }
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
Definition: StdExpansion.h:158
static void Dgemv(const char &trans, const int &m, const int &n, const double &alpha, const double *a, const int &lda, const double *x, const int &incx, const double &beta, double *y, const int &incy)
BLAS level 2: Matrix vector multiply y = A x where A[m x n].
Definition: Blas.hpp:265
static void Dgemm(const char &transa, const char &transb, const int &m, const int &n, const int &k, const double &alpha, const double *a, const int &lda, const double *b, const int &ldb, const double &beta, double *c, const int &ldc)
BLAS level 3: Matrix-matrix multiply C = A x B where op(A)[m x k], op(B)[k x n], C[m x n] DGEMM perfo...
Definition: Blas.hpp:394
static void Daxpy(const int &n, const double &alpha, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: y = alpha x plus y.
Definition: Blas.hpp:167
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:48

References Blas::Daxpy(), Blas::Dgemm(), Blas::Dgemv(), Nektar::LibUtilities::eModified_A, Nektar::StdRegions::StdExpansion::GetBasisType(), and Nektar::StdRegions::StdExpansion::m_base.

Referenced by v_BwdTrans_SumFac().

◆ v_CalcNumberOfCoefficients()

int Nektar::StdRegions::StdPyrExp::v_CalcNumberOfCoefficients ( const std::vector< unsigned int > &  nummodes,
int &  modes_offset 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1072 of file StdPyrExp.cpp.

1075  {
1077  nummodes[modes_offset],
1078  nummodes[modes_offset+1],
1079  nummodes[modes_offset+2]);
1080 
1081  modes_offset += 3;
1082  return nmodes;
1083  }

References Nektar::LibUtilities::StdPyrData::getNumberOfCoefficients().

◆ v_CreateStdMatrix()

DNekMatSharedPtr Nektar::StdRegions::StdPyrExp::v_CreateStdMatrix ( const StdMatrixKey mkey)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 1865 of file StdPyrExp.cpp.

1866  {
1867  return v_GenMatrix(mkey);
1868  }
virtual DNekMatSharedPtr v_GenMatrix(const StdMatrixKey &mkey)
Definition: StdPyrExp.cpp:1860

References v_GenMatrix().

◆ v_DetShapeType()

LibUtilities::ShapeType Nektar::StdRegions::StdPyrExp::v_DetShapeType ( ) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 906 of file StdPyrExp.cpp.

907  {
908  return LibUtilities::ePyramid;
909  }

References Nektar::LibUtilities::ePyramid.

◆ v_FillMode()

void Nektar::StdRegions::StdPyrExp::v_FillMode ( const int  mode,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 780 of file StdPyrExp.cpp.

781  {
782  Array<OneD, NekDouble> tmp(m_ncoeffs, 0.0);
783  tmp[mode] = 1.0;
784  v_BwdTrans(tmp, outarray);
785  }
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Backward transformation is evaluated at the quadrature points.
Definition: StdPyrExp.cpp:252

References Nektar::StdRegions::StdExpansion::m_ncoeffs, and v_BwdTrans().

◆ v_FwdTrans()

void Nektar::StdRegions::StdPyrExp::v_FwdTrans ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Forward transform from physical quadrature space stored in inarray and evaluate the expansion coefficients and store in outarray.

Inputs:

  • inarray: array of physical quadrature points to be transformed

Outputs:

  • outarray: updated array of expansion coefficients.

Implements Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 393 of file StdPyrExp.cpp.

395  {
396  v_IProductWRTBase(inarray,outarray);
397 
398  // get Mass matrix inverse
399  StdMatrixKey imasskey(eInvMass,DetShapeType(),*this);
400  DNekMatSharedPtr imatsys = GetStdMatrix(imasskey);
401 
402 
403  // copy inarray in case inarray == outarray
404  DNekVec in (m_ncoeffs, outarray);
405  DNekVec out(m_ncoeffs, outarray, eWrapper);
406 
407  out = (*imatsys)*in;
408  }
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
Definition: StdExpansion.h:617
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
Definition: StdExpansion.h:376
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Inner product of inarray over region with respect to the expansion basis m_base[0]->GetBdata(),...
Definition: StdPyrExp.cpp:429
NekVector< NekDouble > DNekVec
Definition: NekTypeDefs.hpp:48
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:69

References Nektar::StdRegions::StdExpansion::DetShapeType(), Nektar::StdRegions::eInvMass, Nektar::eWrapper, Nektar::StdRegions::StdExpansion::GetStdMatrix(), Nektar::StdRegions::StdExpansion::m_ncoeffs, and v_IProductWRTBase().

◆ v_GenMatrix()

DNekMatSharedPtr Nektar::StdRegions::StdPyrExp::v_GenMatrix ( const StdMatrixKey mkey)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 1860 of file StdPyrExp.cpp.

1861  {
1862  return CreateGeneralMatrix(mkey);
1863  }
DNekMatSharedPtr CreateGeneralMatrix(const StdMatrixKey &mkey)
this function generates the mass matrix

References Nektar::StdRegions::StdExpansion::CreateGeneralMatrix().

Referenced by v_CreateStdMatrix().

◆ v_GetBoundaryMap()

void Nektar::StdRegions::StdPyrExp::v_GetBoundaryMap ( Array< OneD, unsigned int > &  outarray)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1184 of file StdPyrExp.cpp.

1185  {
1188  "BasisType is not a boundary interior form");
1191  "BasisType is not a boundary interior form");
1194  "BasisType is not a boundary interior form");
1195 
1196  int P = m_base[0]->GetNumModes() - 1, p;
1197  int Q = m_base[1]->GetNumModes() - 1, q;
1198  int R = m_base[2]->GetNumModes() - 1, r;
1199  int idx = 0;
1200 
1201  int nBnd = NumBndryCoeffs();
1202 
1203  if (maparray.size() != nBnd)
1204  {
1205  maparray = Array<OneD, unsigned int>(nBnd);
1206  }
1207 
1208  // Loop over all boundary modes (in ascending order).
1209  for (p = 0; p <= P; ++p)
1210  {
1211  // First two q-r planes are entirely boundary modes.
1212  if (p <= 1)
1213  {
1214  for (q = 0; q <= Q; ++q)
1215  {
1216  int maxpq = max(p,q);
1217  for (r = 0; r <= R-maxpq; ++r)
1218  {
1219  maparray[idx++] = GetMode(p,q,r);
1220  }
1221  }
1222  }
1223  else
1224  {
1225  // Remaining q-r planes contain boundary modes on the two
1226  // front and back sides and edges 0 2.
1227  for (q = 0; q <= Q; ++q)
1228  {
1229  if (q <= 1)
1230  {
1231  for (r = 0; r <= R-p; ++r)
1232  {
1233  maparray[idx++] = GetMode(p,q,r);
1234  }
1235  }
1236  else
1237  {
1238  maparray[idx++] = GetMode(p,q,0);
1239  }
1240  }
1241  }
1242  }
1243  }
int GetMode(int I, int J, int K)
Compute the mode number in the expansion for a particular tensorial combination.
Definition: StdPyrExp.cpp:1895
@ eGLL_Lagrange
Lagrange for SEM basis .
Definition: BasisType.h:54
P
Definition: main.py:133

References ASSERTL1, Nektar::LibUtilities::eGLL_Lagrange, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModifiedPyr_C, Nektar::StdRegions::StdExpansion::GetBasisType(), GetMode(), Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::NumBndryCoeffs(), CellMLToNektar.cellml_metadata::p, and main::P.

◆ v_GetCoords()

void Nektar::StdRegions::StdPyrExp::v_GetCoords ( Array< OneD, NekDouble > &  xi_x,
Array< OneD, NekDouble > &  xi_y,
Array< OneD, NekDouble > &  xi_z 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 752 of file StdPyrExp.cpp.

755  {
756  Array<OneD, const NekDouble> etaBar_x = m_base[0]->GetZ();
757  Array<OneD, const NekDouble> eta_y = m_base[1]->GetZ();
758  Array<OneD, const NekDouble> eta_z = m_base[2]->GetZ();
759  int Qx = GetNumPoints(0);
760  int Qy = GetNumPoints(1);
761  int Qz = GetNumPoints(2);
762 
763  // Convert collapsed coordinates into cartesian coordinates: eta --> xi
764  for (int k = 0; k < Qz; ++k )
765  {
766  for (int j = 0; j < Qy; ++j)
767  {
768  for (int i = 0; i < Qx; ++i)
769  {
770  int s = i + Qx*(j + Qy*k);
771 
772  xi_z[s] = eta_z[k];
773  xi_y[s] = (1.0 + eta_y[j]) * (1.0 - eta_z[k]) / 2.0 - 1.0;
774  xi_x[s] = (1.0 + etaBar_x[i]) * (1.0 - eta_z[k]) / 2.0 - 1.0;
775  }
776  }
777  }
778  }

References Nektar::StdRegions::StdExpansion::GetNumPoints(), and Nektar::StdRegions::StdExpansion::m_base.

◆ v_GetEdgeInteriorToElementMap()

void Nektar::StdRegions::StdPyrExp::v_GetEdgeInteriorToElementMap ( const int  tid,
Array< OneD, unsigned int > &  maparray,
Array< OneD, int > &  signarray,
const Orientation  traceOrient = eDir1FwdDir1_Dir2FwdDir2 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion3D.

Definition at line 1577 of file StdPyrExp.cpp.

1582  {
1583  int i;
1584  bool signChange;
1585  const int P = m_base[0]->GetNumModes() - 1;
1586  const int Q = m_base[1]->GetNumModes() - 1;
1587  const int R = m_base[2]->GetNumModes() - 1;
1588  const int nEdgeIntCoeffs = v_GetEdgeNcoeffs(eid) - 2;
1589 
1590  if (maparray.size() != nEdgeIntCoeffs)
1591  {
1592  maparray = Array<OneD, unsigned int>(nEdgeIntCoeffs);
1593  }
1594 
1595  if(signarray.size() != nEdgeIntCoeffs)
1596  {
1597  signarray = Array<OneD, int>(nEdgeIntCoeffs,1);
1598  }
1599  else
1600  {
1601  fill(signarray.get(), signarray.get()+nEdgeIntCoeffs, 1);
1602  }
1603 
1604  // If edge is oriented backwards, change sign of modes which have
1605  // degree 2n+1, n >= 1.
1606  signChange = edgeOrient == eBackwards;
1607 
1608  switch (eid)
1609  {
1610  case 0:
1611  for (i = 2; i <= P; ++i)
1612  {
1613  maparray[i-2] = GetMode(i,0,0);
1614  }
1615  break;
1616 
1617  case 1:
1618  for (i = 2; i <= Q; ++i)
1619  {
1620  maparray[i-2] = GetMode(1,i,0);
1621  }
1622  break;
1623  case 2:
1624  for (i = 2; i <= P; ++i)
1625  {
1626  maparray[i-2] = GetMode(i,1,0);
1627  }
1628  break;
1629 
1630  case 3:
1631  for (i = 2; i <= Q; ++i)
1632  {
1633  maparray[i-2] = GetMode(0,i,0);
1634  }
1635  break;
1636  case 4:
1637  for (i = 2; i <= R; ++i)
1638  {
1639  maparray[i-2] = GetMode(0,0,i);
1640  }
1641  break;
1642 
1643  case 5:
1644  for (i = 1; i <= R-1; ++i)
1645  {
1646  maparray[i-1] = GetMode(1,0,i);
1647  }
1648  break;
1649  case 6:
1650  for (i = 1; i <= R-1; ++i)
1651  {
1652  maparray[i-1] = GetMode(1,1,i);
1653  }
1654  break;
1655 
1656  case 7:
1657  for (i = 1; i <= R-1; ++i)
1658  {
1659  maparray[i-1] = GetMode(0,1,i);
1660  }
1661  break;
1662  default:
1663  ASSERTL0(false, "Edge not defined.");
1664  break;
1665  }
1666 
1667  if (signChange)
1668  {
1669  for (i = 1; i < nEdgeIntCoeffs; i += 2)
1670  {
1671  signarray[i] = -1;
1672  }
1673  }
1674  }
virtual int v_GetEdgeNcoeffs(const int i) const
Definition: StdPyrExp.cpp:1016

References ASSERTL0, Nektar::StdRegions::eBackwards, GetMode(), Nektar::StdRegions::StdExpansion::m_base, main::P, and v_GetEdgeNcoeffs().

◆ v_GetEdgeNcoeffs()

int Nektar::StdRegions::StdPyrExp::v_GetEdgeNcoeffs ( const int  i) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion3D.

Definition at line 1016 of file StdPyrExp.cpp.

1017  {
1018  ASSERTL2(i >= 0 && i <= 7, "edge id is out of range");
1019 
1020  if (i == 0 || i == 2)
1021  {
1022  return GetBasisNumModes(0);
1023  }
1024  else if (i == 1 || i == 3)
1025  {
1026  return GetBasisNumModes(1);
1027  }
1028  else
1029  {
1030  return GetBasisNumModes(2);
1031  }
1032  }
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed to...
Definition: ErrorUtil.hpp:274
int GetBasisNumModes(const int dir) const
This function returns the number of expansion modes in the dir direction.
Definition: StdExpansion.h:171

References ASSERTL2, and Nektar::StdRegions::StdExpansion::GetBasisNumModes().

Referenced by v_GetEdgeInteriorToElementMap().

◆ v_GetInteriorMap()

void Nektar::StdRegions::StdPyrExp::v_GetInteriorMap ( Array< OneD, unsigned int > &  outarray)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1144 of file StdPyrExp.cpp.

1145  {
1148  "BasisType is not a boundary interior form");
1151  "BasisType is not a boundary interior form");
1154  "BasisType is not a boundary interior form");
1155 
1156 
1157  int P = m_base[0]->GetNumModes() - 1, p;
1158  int Q = m_base[1]->GetNumModes() - 1, q;
1159  int R = m_base[2]->GetNumModes() - 1, r;
1160 
1161  int nIntCoeffs = m_ncoeffs - NumBndryCoeffs();
1162 
1163  if(outarray.size()!=nIntCoeffs)
1164  {
1165  outarray = Array<OneD, unsigned int>(nIntCoeffs);
1166  }
1167 
1168  int idx = 0;
1169 
1170  // Loop over all interior modes.
1171  for (p = 2; p <= P; ++p)
1172  {
1173  for (q = 2; q <= Q; ++q)
1174  {
1175  int maxpq = max(p,q);
1176  for (r = 1; r <= R-maxpq; ++r)
1177  {
1178  outarray[idx++] = GetMode(p,q,r);
1179  }
1180  }
1181  }
1182  }

References ASSERTL1, Nektar::LibUtilities::eGLL_Lagrange, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModifiedPyr_C, Nektar::StdRegions::StdExpansion::GetBasisType(), GetMode(), Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::m_ncoeffs, Nektar::StdRegions::StdExpansion::NumBndryCoeffs(), CellMLToNektar.cellml_metadata::p, and main::P.

◆ v_GetNedges()

int Nektar::StdRegions::StdPyrExp::v_GetNedges ( void  ) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion3D.

Definition at line 896 of file StdPyrExp.cpp.

897  {
898  return 8;
899  }

◆ v_GetNtraces()

int Nektar::StdRegions::StdPyrExp::v_GetNtraces ( ) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 901 of file StdPyrExp.cpp.

902  {
903  return 5;
904  }

◆ v_GetNverts()

int Nektar::StdRegions::StdPyrExp::v_GetNverts ( ) const
protectedvirtual

Implements Nektar::StdRegions::StdExpansion.

Definition at line 891 of file StdPyrExp.cpp.

892  {
893  return 5;
894  }

◆ v_GetTraceBasisKey()

const LibUtilities::BasisKey Nektar::StdRegions::StdPyrExp::v_GetTraceBasisKey ( const int  i,
const int  k 
) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1034 of file StdPyrExp.cpp.

1036  {
1037  ASSERTL2(i >= 0 && i <= 4, "face id is out of range");
1038  ASSERTL2(k >= 0 && k <= 1, "basis key id is out of range");
1039 
1040  switch(i)
1041  {
1042  case 0:
1043  {
1044  return EvaluateQuadFaceBasisKey(k,
1045  m_base[k]->GetBasisType(),
1046  m_base[k]->GetNumPoints(),
1047  m_base[k]->GetNumModes());
1048 
1049  }
1050  case 1:
1051  case 3:
1052  {
1053  return EvaluateTriFaceBasisKey(k,
1054  m_base[2*k]->GetBasisType(),
1055  m_base[2*k]->GetNumPoints(),
1056  m_base[2*k]->GetNumModes());
1057  }
1058  case 2:
1059  case 4:
1060  {
1061  return EvaluateTriFaceBasisKey(k,
1062  m_base[k+1]->GetBasisType(),
1063  m_base[k+1]->GetNumPoints(),
1064  m_base[k+1]->GetNumModes());
1065  }
1066  }
1067 
1068  // Should never get here.
1070  }
static const BasisKey NullBasisKey(eNoBasisType, 0, NullPointsKey)
Defines a null basis with no type or points.
LibUtilities::BasisKey EvaluateTriFaceBasisKey(const int facedir, const LibUtilities::BasisType faceDirBasisType, const int numpoints, const int nummodes)
LibUtilities::BasisKey EvaluateQuadFaceBasisKey(const int facedir, const LibUtilities::BasisType faceDirBasisType, const int numpoints, const int nummodes)

References ASSERTL2, Nektar::StdRegions::EvaluateQuadFaceBasisKey(), Nektar::StdRegions::EvaluateTriFaceBasisKey(), Nektar::StdRegions::StdExpansion::GetBasisType(), Nektar::StdRegions::StdExpansion::GetNumPoints(), Nektar::StdRegions::StdExpansion::m_base, and Nektar::LibUtilities::NullBasisKey().

◆ v_GetTraceInteriorToElementMap()

void Nektar::StdRegions::StdPyrExp::v_GetTraceInteriorToElementMap ( const int  tid,
Array< OneD, unsigned int > &  maparray,
Array< OneD, int > &  signarray,
const Orientation  traceOrient = eDir1FwdDir1_Dir2FwdDir2 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1676 of file StdPyrExp.cpp.

1681  {
1682  const int P = m_base[0]->GetNumModes() - 1;
1683  const int Q = m_base[1]->GetNumModes() - 1;
1684  const int R = m_base[2]->GetNumModes() - 1;
1685  const int nFaceIntCoeffs = v_GetTraceIntNcoeffs(fid);
1686  int p, q, r, idx = 0;
1687  int nummodesA = 0;
1688  int nummodesB = 0;
1689  int i, j;
1690 
1691  if (maparray.size() != nFaceIntCoeffs)
1692  {
1693  maparray = Array<OneD, unsigned int>(nFaceIntCoeffs);
1694  }
1695 
1696  if (signarray.size() != nFaceIntCoeffs)
1697  {
1698  signarray = Array<OneD, int>(nFaceIntCoeffs, 1);
1699  }
1700  else
1701  {
1702  fill(signarray.get(), signarray.get()+nFaceIntCoeffs, 1);
1703  }
1704 
1705  // Set up an array indexing for quad faces, since the ordering may
1706  // need to be transposed depending on orientation.
1707  Array<OneD, int> arrayindx(nFaceIntCoeffs);
1708  if (fid == 0)
1709  {
1710  nummodesA = P-1;
1711  nummodesB = Q-1;
1712 
1713  for (i = 0; i < nummodesB; i++)
1714  {
1715  for (j = 0; j < nummodesA; j++)
1716  {
1717  if (faceOrient < 9)
1718  {
1719  arrayindx[i*nummodesA+j] = i*nummodesA+j;
1720  }
1721  else
1722  {
1723  arrayindx[i*nummodesA+j] = j*nummodesB+i;
1724  }
1725  }
1726  }
1727  }
1728 
1729  switch (fid)
1730  {
1731  case 0: // Bottom quad
1732  for (q = 2; q <= Q; ++q)
1733  {
1734  for (p = 2; p <= P; ++p)
1735  {
1736  maparray[arrayindx[(q-2)*nummodesA+(p-2)]] = GetMode(p,q,0);
1737  }
1738  }
1739  break;
1740  case 1: // Front triangle
1741  for (p = 2; p <= P; ++p)
1742  {
1743  for (r = 1; r <= R-p; ++r)
1744  {
1745  if ((int)faceOrient == 7)
1746  {
1747  signarray[idx] = p % 2 ? -1 : 1;
1748  }
1749  maparray[idx++] = GetMode(p,0,r);
1750  }
1751  }
1752  break;
1753  case 2: // Right triangle
1754  for (q = 2; q <= Q; ++q)
1755  {
1756  for (r = 1; r <= R-q; ++r)
1757  {
1758  if ((int)faceOrient == 7)
1759  {
1760  signarray[idx] = q % 2 ? -1 : 1;
1761  }
1762  maparray[idx++] = GetMode(1, q, r);
1763  }
1764  }
1765  break;
1766 
1767  case 3: // Rear triangle
1768  for (p = 2; p <= P; ++p)
1769  {
1770  for (r = 1; r <= R-p; ++r)
1771  {
1772  if ((int)faceOrient == 7)
1773  {
1774  signarray[idx] = p % 2 ? -1 : 1;
1775  }
1776  maparray[idx++] = GetMode(p, 1, r);
1777  }
1778  }
1779  break;
1780 
1781  case 4: // Left triangle
1782  for (q = 2; q <= Q; ++q)
1783  {
1784  for (r = 1; r <= R-q; ++r)
1785  {
1786  if ((int)faceOrient == 7)
1787  {
1788  signarray[idx] = q % 2 ? -1 : 1;
1789  }
1790  maparray[idx++] = GetMode(0, q, r);
1791  }
1792  }
1793  break;
1794  default:
1795  ASSERTL0(false, "Face interior map not available.");
1796  }
1797 
1798  // Triangular faces are processed in the above switch loop; for
1799  // remaining quad faces, set up orientation if necessary.
1800  if (fid > 0)
1801  {
1802  return;
1803  }
1804 
1805  if (faceOrient == 6 || faceOrient == 8 ||
1806  faceOrient == 11 || faceOrient == 12)
1807  {
1808  if (faceOrient < 9)
1809  {
1810  for (i = 1; i < nummodesB; i += 2)
1811  {
1812  for (j = 0; j < nummodesA; j++)
1813  {
1814  signarray[arrayindx[i*nummodesA+j]] *= -1;
1815  }
1816  }
1817  }
1818  else
1819  {
1820  for (i = 0; i < nummodesB; i++)
1821  {
1822  for (j = 1; j < nummodesA; j += 2)
1823  {
1824  signarray[arrayindx[i*nummodesA+j]] *= -1;
1825  }
1826  }
1827  }
1828  }
1829 
1830  if (faceOrient == 7 || faceOrient == 8 ||
1831  faceOrient == 10 || faceOrient == 12)
1832  {
1833  if (faceOrient < 9)
1834  {
1835  for (i = 0; i < nummodesB; i++)
1836  {
1837  for (j = 1; j < nummodesA; j += 2)
1838  {
1839  signarray[arrayindx[i*nummodesA+j]] *= -1;
1840  }
1841  }
1842  }
1843  else
1844  {
1845  for (i = 1; i < nummodesB; i += 2)
1846  {
1847  for (j = 0; j < nummodesA; j++)
1848  {
1849  signarray[arrayindx[i*nummodesA+j]] *= -1;
1850  }
1851  }
1852  }
1853  }
1854  }
virtual int v_GetTraceIntNcoeffs(const int i) const
Definition: StdPyrExp.cpp:973

References ASSERTL0, GetMode(), Nektar::StdRegions::StdExpansion::m_base, CellMLToNektar.cellml_metadata::p, main::P, and v_GetTraceIntNcoeffs().

◆ v_GetTraceIntNcoeffs()

int Nektar::StdRegions::StdPyrExp::v_GetTraceIntNcoeffs ( const int  i) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 973 of file StdPyrExp.cpp.

974  {
975  ASSERTL2(i >= 0 && i <= 4, "face id is out of range");
976 
977  int P = m_base[0]->GetNumModes()-1;
978  int Q = m_base[1]->GetNumModes()-1;
979  int R = m_base[2]->GetNumModes()-1;
980 
981  if (i == 0)
982  {
983  return (P-1)*(Q-1);
984  }
985  else if (i == 1 || i == 3)
986  {
987  return (P-1) * (2*(R-1) - (P-1) - 1) / 2;
988  }
989  else
990  {
991  return (Q-1) * (2*(R-1) - (Q-1) - 1) / 2;
992  }
993  }

References ASSERTL2, Nektar::StdRegions::StdExpansion::m_base, and main::P.

Referenced by v_GetTraceInteriorToElementMap().

◆ v_GetTraceNcoeffs()

int Nektar::StdRegions::StdPyrExp::v_GetTraceNcoeffs ( const int  i) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 953 of file StdPyrExp.cpp.

954  {
955  ASSERTL2(i >= 0 && i <= 4, "face id is out of range");
956 
957  if (i == 0)
958  {
959  return GetBasisNumModes(0)*GetBasisNumModes(1);
960  }
961  else if (i == 1 || i == 3)
962  {
963  int P = GetBasisNumModes(0)-1, Q = GetBasisNumModes(2)-1;
964  return Q+1 + (P*(1 + 2*Q - P))/2;
965  }
966  else
967  {
968  int P = GetBasisNumModes(1)-1, Q = GetBasisNumModes(2)-1;
969  return Q+1 + (P*(1 + 2*Q - P))/2;
970  }
971  }

References ASSERTL2, Nektar::StdRegions::StdExpansion::GetBasisNumModes(), and main::P.

◆ v_GetTraceNumModes()

void Nektar::StdRegions::StdPyrExp::v_GetTraceNumModes ( const int  fid,
int &  numModes0,
int &  numModes1,
Orientation  faceOrient = eDir1FwdDir1_Dir2FwdDir2 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 787 of file StdPyrExp.cpp.

792  {
793  int nummodes [3] = {m_base[0]->GetNumModes(),
794  m_base[1]->GetNumModes(),
795  m_base[2]->GetNumModes()};
796  switch(fid)
797  {
798  // quad
799  case 0:
800  {
801  numModes0 = nummodes[0];
802  numModes1 = nummodes[1];
803  }
804  break;
805  case 1:
806  case 3:
807  {
808  numModes0 = nummodes[0];
809  numModes1 = nummodes[2];
810  }
811  break;
812  case 2:
813  case 4:
814  {
815  numModes0 = nummodes[1];
816  numModes1 = nummodes[2];
817  }
818  break;
819  }
820 
821  if ( faceOrient >= 9 )
822  {
823  std::swap(numModes0, numModes1);
824  }
825  }

References Nektar::StdRegions::StdExpansion::m_base.

◆ v_GetTraceNumPoints()

int Nektar::StdRegions::StdPyrExp::v_GetTraceNumPoints ( const int  i) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 995 of file StdPyrExp.cpp.

996  {
997  ASSERTL2(i >= 0 && i <= 4, "face id is out of range");
998 
999  if (i == 0)
1000  {
1001  return m_base[0]->GetNumPoints()*
1002  m_base[1]->GetNumPoints();
1003  }
1004  else if (i == 1 || i == 3)
1005  {
1006  return m_base[0]->GetNumPoints()*
1007  m_base[2]->GetNumPoints();
1008  }
1009  else
1010  {
1011  return m_base[1]->GetNumPoints()*
1012  m_base[2]->GetNumPoints();
1013  }
1014  }

References ASSERTL2, and Nektar::StdRegions::StdExpansion::m_base.

◆ v_GetTraceToElementMap()

void Nektar::StdRegions::StdPyrExp::v_GetTraceToElementMap ( const int  fid,
Array< OneD, unsigned int > &  maparray,
Array< OneD, int > &  signarray,
Orientation  faceOrient,
int  P,
int  Q 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1245 of file StdPyrExp.cpp.

1252  {
1254  "Method only implemented if BasisType is identical"
1255  "in x and y directions");
1258  "Method only implemented for Modified_A BasisType"
1259  "(x and y direction) and ModifiedPyr_C BasisType (z "
1260  "direction)");
1261 
1262  int i, j, k, p, q, r, nFaceCoeffs, idx = 0;
1263  int nummodesA=0, nummodesB=0;
1264 
1265  int order0 = m_base[0]->GetNumModes();
1266  int order1 = m_base[1]->GetNumModes();
1267  int order2 = m_base[2]->GetNumModes();
1268 
1269  switch (fid)
1270  {
1271  case 0:
1272  nummodesA = order0;
1273  nummodesB = order1;
1274  break;
1275  case 1:
1276  case 3:
1277  nummodesA = order0;
1278  nummodesB = order2;
1279  break;
1280  case 2:
1281  case 4:
1282  nummodesA = order1;
1283  nummodesB = order2;
1284  break;
1285  default:
1286  ASSERTL0(false,"fid must be between 0 and 4");
1287  }
1288 
1289  bool CheckForZeroedModes = false;
1290 
1291  if (P == -1)
1292  {
1293  P = nummodesA;
1294  Q = nummodesB;
1295  nFaceCoeffs = GetTraceNcoeffs(fid);
1296  }
1297  else if (fid > 0)
1298  {
1299  nFaceCoeffs = P*(2*Q-P+1)/2;
1300  CheckForZeroedModes = true;
1301  }
1302  else
1303  {
1304  nFaceCoeffs = P*Q;
1305  CheckForZeroedModes = true;
1306  }
1307 
1308  // Allocate the map array and sign array; set sign array to ones (+)
1309  if (maparray.size() != nFaceCoeffs)
1310  {
1311  maparray = Array<OneD, unsigned int>(nFaceCoeffs);
1312  }
1313 
1314  if (signarray.size() != nFaceCoeffs)
1315  {
1316  signarray = Array<OneD, int>(nFaceCoeffs,1);
1317  }
1318  else
1319  {
1320  fill(signarray.get(), signarray.get() + nFaceCoeffs, 1);
1321  }
1322 
1323  // Set up an array indexing for quads, since the ordering may need
1324  // to be transposed.
1325  Array<OneD, int> arrayindx(nFaceCoeffs,-1);
1326 
1327  if (fid == 0)
1328  {
1329  for (i = 0; i < Q; i++)
1330  {
1331  for (j = 0; j < P; j++)
1332  {
1333  if (faceOrient < 9)
1334  {
1335  arrayindx[i*P+j] = i*P+j;
1336  }
1337  else
1338  {
1339  arrayindx[i*P+j] = j*Q+i;
1340  }
1341  }
1342  }
1343  }
1344 
1345  // Set up ordering inside each 2D face. Also for triangular faces,
1346  // populate signarray.
1347  switch (fid)
1348  {
1349  case 0: // Bottom quad
1350 
1351  for (q = 0; q < Q; ++q)
1352  {
1353  for (p = 0; p < P; ++p)
1354  {
1355  maparray[arrayindx[q*P+p]] = GetMode(p,q,0);
1356  }
1357  }
1358  break;
1359 
1360  case 1: // Front triangle
1361  for (p = 0; p < P; ++p)
1362  {
1363  for (r = 0; r < Q-p; ++r)
1364  {
1365  if ((int)faceOrient == 7 && p > 1)
1366  {
1367  signarray[idx] = p % 2 ? -1 : 1;
1368  }
1369  maparray[idx++] = GetMode(p,0,r);
1370  }
1371  }
1372  break;
1373 
1374  case 2: // Right triangle
1375  maparray[idx++] = GetMode(1,0,0);
1376  maparray[idx++] = GetMode(0,0,1);
1377  for (r = 1; r < Q-1; ++r)
1378  {
1379  maparray[idx++] = GetMode(1,0,r);
1380  }
1381 
1382  for (q = 1; q < P; ++q)
1383  {
1384  for (r = 0; r < Q-q; ++r)
1385  {
1386  if ((int)faceOrient == 7 && q > 1)
1387  {
1388  signarray[idx] = q % 2 ? -1 : 1;
1389  }
1390  maparray[idx++] = GetMode(1,q,r);
1391  }
1392  }
1393  break;
1394 
1395  case 3: // Rear triangle
1396  maparray[idx++] = GetMode(0,1,0);
1397  maparray[idx++] = GetMode(0,0,1);
1398  for (r = 1; r < Q-1; ++r)
1399  {
1400  maparray[idx++] = GetMode(0,1,r);
1401  }
1402 
1403  for (p = 1; p < P; ++p)
1404  {
1405  for (r = 0; r < Q-p; ++r)
1406  {
1407  if ((int)faceOrient == 7 && p > 1)
1408  {
1409  signarray[idx] = p % 2 ? -1 : 1;
1410  }
1411  maparray[idx++] = GetMode(p, 1, r);
1412  }
1413  }
1414  break;
1415 
1416  case 4: // Left triangle
1417  for (q = 0; q < P; ++q)
1418  {
1419  for (r = 0; r < Q-q; ++r)
1420  {
1421  if ((int)faceOrient == 7 && q > 1)
1422  {
1423  signarray[idx] = q % 2 ? -1 : 1;
1424  }
1425  maparray[idx++] = GetMode(0,q,r);
1426  }
1427  }
1428  break;
1429 
1430  default:
1431  ASSERTL0(false, "Face to element map unavailable.");
1432  }
1433 
1434  if (fid > 0)
1435  {
1436  if(CheckForZeroedModes)
1437  {
1438  // zero signmap and set maparray to zero if elemental
1439  // modes are not as large as face modesl
1440  int idx = 0;
1441  for (j = 0; j < P; ++j)
1442  {
1443  idx += Q-j;
1444  for (k = Q-j; k < Q-j; ++k)
1445  {
1446  signarray[idx] = 0.0;
1447  maparray[idx++] = maparray[0];
1448  }
1449  }
1450 
1451  for (j = P; j < P; ++j)
1452  {
1453  for (k = 0; k < Q-j; ++k)
1454  {
1455  signarray[idx] = 0.0;
1456  maparray[idx++] = maparray[0];
1457  }
1458  }
1459  }
1460 
1461  // Triangles only have one possible orientation (base
1462  // direction reversed); swap edge modes.
1463  if ((int)faceOrient == 7)
1464  {
1465  swap(maparray[0], maparray[Q]);
1466  for (i = 1; i < Q-1; ++i)
1467  {
1468  swap(maparray[i+1], maparray[Q+i]);
1469  }
1470  }
1471  }
1472  else
1473  {
1474  if(CheckForZeroedModes)
1475  {
1476  // zero signmap and set maparray to zero if elemental
1477  // modes are not as large as face modesl
1478  for (j = 0; j < P; ++j)
1479  {
1480  for (k = Q; k < Q; ++k)
1481  {
1482  signarray[arrayindx[j+k*P]] = 0.0;
1483  maparray[arrayindx[j+k*P]] = maparray[0];
1484  }
1485  }
1486 
1487  for (j = P; j < P; ++j)
1488  {
1489  for (k = 0; k < Q; ++k)
1490  {
1491  signarray[arrayindx[j+k*P]] = 0.0;
1492  maparray[arrayindx[j+k*P]] = maparray[0];
1493  }
1494  }
1495  }
1496 
1497  // The code below is exactly the same as that taken from
1498  // StdHexExp and reverses the 'b' and 'a' directions as
1499  // appropriate (1st and 2nd if statements respectively) in
1500  // quadrilateral faces.
1501  if (faceOrient == 6 || faceOrient == 8 ||
1502  faceOrient == 11 || faceOrient == 12)
1503  {
1504  if (faceOrient < 9)
1505  {
1506  for (i = 3; i < Q; i += 2)
1507  {
1508  for (j = 0; j < P; j++)
1509  {
1510  signarray[arrayindx[i*P+j]] *= -1;
1511  }
1512  }
1513 
1514  for (i = 0; i < P; i++)
1515  {
1516  swap(maparray [i], maparray [i+P]);
1517  swap(signarray[i], signarray[i+P]);
1518  }
1519  }
1520  else
1521  {
1522  for (i = 0; i < Q; i++)
1523  {
1524  for (j = 3; j < P; j += 2)
1525  {
1526  signarray[arrayindx[i*P+j]] *= -1;
1527  }
1528  }
1529 
1530  for (i = 0; i < Q; i++)
1531  {
1532  swap (maparray [i], maparray [i+Q]);
1533  swap (signarray[i], signarray[i+Q]);
1534  }
1535  }
1536  }
1537 
1538  if (faceOrient == 7 || faceOrient == 8 ||
1539  faceOrient == 10 || faceOrient == 12)
1540  {
1541  if (faceOrient < 9)
1542  {
1543  for (i = 0; i < Q; i++)
1544  {
1545  for (j = 3; j < P; j += 2)
1546  {
1547  signarray[arrayindx[i*P+j]] *= -1;
1548  }
1549  }
1550 
1551  for(i = 0; i < Q; i++)
1552  {
1553  swap(maparray [i*P], maparray [i*P+1]);
1554  swap(signarray[i*P], signarray[i*P+1]);
1555  }
1556  }
1557  else
1558  {
1559  for (i = 3; i < Q; i += 2)
1560  {
1561  for (j = 0; j < P; j++)
1562  {
1563  signarray[arrayindx[i*P+j]] *= -1;
1564  }
1565  }
1566 
1567  for (i = 0; i < P; i++)
1568  {
1569  swap(maparray [i*Q], maparray [i*Q+1]);
1570  swap(signarray[i*Q], signarray[i*Q+1]);
1571  }
1572  }
1573  }
1574  }
1575  }
int GetTraceNcoeffs(const int i) const
This function returns the number of expansion coefficients belonging to the i-th trace.
Definition: StdExpansion.h:265

References ASSERTL0, ASSERTL1, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModifiedPyr_C, Nektar::StdRegions::StdExpansion::GetBasisType(), GetMode(), Nektar::StdRegions::StdExpansion::GetTraceNcoeffs(), Nektar::StdRegions::StdExpansion::m_base, CellMLToNektar.cellml_metadata::p, and main::P.

◆ v_GetVertexMap()

int Nektar::StdRegions::StdPyrExp::v_GetVertexMap ( int  localVertexId,
bool  useCoeffPacking = false 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1085 of file StdPyrExp.cpp.

1086  {
1090  "Mapping not defined for this type of basis");
1091 
1092  int l = 0;
1093 
1094  if(useCoeffPacking == true) // follow packing of coefficients i.e q,r,p
1095  {
1096  switch (vId)
1097  {
1098  case 0:
1099  l = GetMode(0,0,0);
1100  break;
1101  case 1:
1102  l = GetMode(0,0,1);
1103  break;
1104  case 2:
1105  l = GetMode(0,1,0);
1106  break;
1107  case 3:
1108  l = GetMode(1,0,0);
1109  break;
1110  case 4:
1111  l = GetMode(1,1,0);
1112  break;
1113  default:
1114  ASSERTL0(false, "local vertex id must be between 0 and 4");
1115  }
1116  }
1117  else
1118  {
1119  switch (vId)
1120  {
1121  case 0:
1122  l = GetMode(0,0,0);
1123  break;
1124  case 1:
1125  l = GetMode(1,0,0);
1126  break;
1127  case 2:
1128  l = GetMode(1,1,0);
1129  break;
1130  case 3:
1131  l = GetMode(0,1,0);
1132  break;
1133  case 4:
1134  l = GetMode(0,0,1);
1135  break;
1136  default:
1137  ASSERTL0(false, "local vertex id must be between 0 and 4");
1138  }
1139  }
1140 
1141  return l;
1142  }

References ASSERTL0, ASSERTL1, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModifiedPyr_C, Nektar::StdRegions::StdExpansion::GetBasisType(), and GetMode().

◆ v_IProductWRTBase()

void Nektar::StdRegions::StdPyrExp::v_IProductWRTBase ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Inner product of inarray over region with respect to the expansion basis m_base[0]->GetBdata(),m_base[1]->GetBdata(), m_base[2]->GetBdata() and return in outarray.

Wrapper call to StdPyrExp::IProductWRTBase

Input:

  • inarray: array of function evaluated at the physical collocation points

Output:

  • outarray: array of inner product with respect to each basis over region

Implements Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 429 of file StdPyrExp.cpp.

432  {
433  if (m_base[0]->Collocation() &&
434  m_base[1]->Collocation() &&
435  m_base[2]->Collocation())
436  {
437  v_MultiplyByStdQuadratureMetric(inarray, outarray);
438  }
439  else
440  {
441  StdPyrExp::v_IProductWRTBase_SumFac(inarray,outarray);
442  }
443  }
void v_MultiplyByStdQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdPyrExp.cpp:1924
virtual void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
Definition: StdPyrExp.cpp:445

References Nektar::StdRegions::StdExpansion::m_base, v_IProductWRTBase_SumFac(), and v_MultiplyByStdQuadratureMetric().

Referenced by v_FwdTrans().

◆ v_IProductWRTBase_SumFac()

void Nektar::StdRegions::StdPyrExp::v_IProductWRTBase_SumFac ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
bool  multiplybyweights = true 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 445 of file StdPyrExp.cpp.

449  {
450 
451  int nquad1 = m_base[1]->GetNumPoints();
452  int nquad2 = m_base[2]->GetNumPoints();
453  int order0 = m_base[0]->GetNumModes();
454  int order1 = m_base[1]->GetNumModes();
455 
456  Array<OneD, NekDouble> wsp(order0*nquad2*(nquad1 + order1));
457 
458  if(multiplybyweights)
459  {
460  Array<OneD, NekDouble> tmp(inarray.size());
461 
462  v_MultiplyByStdQuadratureMetric(inarray, tmp);
463 
464  v_IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(),
465  m_base[1]->GetBdata(),
466  m_base[2]->GetBdata(),
467  tmp,outarray,wsp,
468  true,true,true);
469  }
470  else
471  {
472  v_IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(),
473  m_base[1]->GetBdata(),
474  m_base[2]->GetBdata(),
475  inarray,outarray,wsp,
476  true,true,true);
477  }
478  }
virtual void v_IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
Definition: StdPyrExp.cpp:480

References Nektar::StdRegions::StdExpansion::m_base, v_IProductWRTBase_SumFacKernel(), and v_MultiplyByStdQuadratureMetric().

Referenced by v_IProductWRTBase().

◆ v_IProductWRTBase_SumFacKernel()

void Nektar::StdRegions::StdPyrExp::v_IProductWRTBase_SumFacKernel ( const Array< OneD, const NekDouble > &  base0,
const Array< OneD, const NekDouble > &  base1,
const Array< OneD, const NekDouble > &  base2,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
Array< OneD, NekDouble > &  wsp,
bool  doCheckCollDir0,
bool  doCheckCollDir1,
bool  doCheckCollDir2 
)
protectedvirtual

Implements Nektar::StdRegions::StdExpansion3D.

Definition at line 480 of file StdPyrExp.cpp.

490  {
491  boost::ignore_unused(doCheckCollDir0, doCheckCollDir1,
492  doCheckCollDir2);
493 
494  int nquad0 = m_base[0]->GetNumPoints();
495  int nquad1 = m_base[1]->GetNumPoints();
496  int nquad2 = m_base[2]->GetNumPoints();
497 
498  int order0 = m_base[0]->GetNumModes();
499  int order1 = m_base[1]->GetNumModes();
500  int order2 = m_base[2]->GetNumModes();
501 
502  ASSERTL1(wsp.size() >= nquad1*nquad2*order0 +
503  nquad2*order0*order1,
504  "Insufficient workspace size");
505 
506  Array<OneD, NekDouble > tmp1 = wsp;
507  Array<OneD, NekDouble > tmp2 = wsp + nquad1*nquad2*order0;
508 
509  int i,j, mode,mode1, cnt;
510 
511  // Inner product with respect to the '0' direction
512  Blas::Dgemm('T', 'N', nquad1*nquad2, order0, nquad0,
513  1.0, inarray.get(), nquad0,
514  base0.get(), nquad0,
515  0.0, tmp1.get(), nquad1*nquad2);
516 
517  // Inner product with respect to the '1' direction
518  for(mode=i=0; i < order0; ++i)
519  {
520  Blas::Dgemm('T', 'N', nquad2, order1, nquad1,
521  1.0, tmp1.get()+i*nquad1*nquad2, nquad1,
522  base1.get(), nquad1,
523  0.0, tmp2.get()+mode*nquad2, nquad2);
524  mode += order1;
525  }
526 
527 
528  // Inner product with respect to the '2' direction
529  mode = mode1 = cnt = 0;
530  for(i = 0; i < order0; ++i)
531  {
532  for(j = 0; j < order1; ++j, ++cnt)
533  {
534  int ijmax = max(i,j);
535 
536  Blas::Dgemv('T', nquad2, order2-ijmax,
537  1.0, base2.get()+mode*nquad2, nquad2,
538  tmp2.get()+cnt*nquad2, 1,
539  0.0, outarray.get()+mode1, 1);
540  mode += order2-ijmax;
541  mode1 += order2-ijmax;
542  }
543 
544  //increment mode in case order1!=order2
545  for(j = order1; j < order2; ++j)
546  {
547  int ijmax = max(i,j);
548  mode += order2-ijmax;
549  }
550  }
551 
552  // fix for modified basis for top singular vertex component
553  // Already have evaluated (1+c)/2 (1-b)/2 (1-a)/2
555  {
556  // add in (1+c)/2 (1+b)/2 (1-a)/2 component
557  outarray[1] += Blas::Ddot(nquad2,base2.get()+nquad2,1,
558  &tmp2[nquad2],1);
559 
560  // add in (1+c)/2 (1-b)/2 (1+a)/2 component
561  outarray[1] += Blas::Ddot(nquad2,base2.get()+nquad2,1,
562  &tmp2[nquad2*order1],1);
563 
564  // add in (1+c)/2 (1+b)/2 (1+a)/2 component
565  outarray[1] += Blas::Ddot(nquad2,base2.get()+nquad2,1,
566  &tmp2[nquad2*order1+nquad2],1);
567  }
568  }
static double Ddot(const int &n, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: output = .
Definition: Blas.hpp:197

References ASSERTL1, Blas::Ddot(), Blas::Dgemm(), Blas::Dgemv(), Nektar::LibUtilities::eModified_A, Nektar::StdRegions::StdExpansion::GetBasisType(), and Nektar::StdRegions::StdExpansion::m_base.

Referenced by v_IProductWRTBase_SumFac().

◆ v_IProductWRTDerivBase()

void Nektar::StdRegions::StdPyrExp::v_IProductWRTDerivBase ( const int  dir,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 570 of file StdPyrExp.cpp.

574  {
575  StdPyrExp::v_IProductWRTDerivBase_SumFac(dir,inarray,outarray);
576  }
virtual void v_IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdPyrExp.cpp:584

References v_IProductWRTDerivBase_SumFac().

◆ v_IProductWRTDerivBase_SumFac()

void Nektar::StdRegions::StdPyrExp::v_IProductWRTDerivBase_SumFac ( const int  dir,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual
Parameters
inarrayFunction evaluated at physical collocation points.
outarrayInner product with respect to each basis function over the element.

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 584 of file StdPyrExp.cpp.

588  {
589  int i;
590  int nquad0 = m_base[0]->GetNumPoints();
591  int nquad1 = m_base[1]->GetNumPoints();
592  int nquad2 = m_base[2]->GetNumPoints();
593  int nqtot = nquad0*nquad1*nquad2;
594 
595  Array<OneD, NekDouble> gfac0(nquad0);
596  Array<OneD, NekDouble> gfac1(nquad1);
597  Array<OneD, NekDouble> gfac2(nquad2);
598  Array<OneD, NekDouble> tmp0 (nqtot);
599 
600 
601  int order0 = m_base[0]->GetNumModes();
602  int order1 = m_base[1]->GetNumModes();
603 
604  Array<OneD, NekDouble> wsp(nquad1*nquad2*order0 +
605  nquad2*order0*order1);
606 
607  const Array<OneD, const NekDouble>& z0 = m_base[0]->GetZ();
608  const Array<OneD, const NekDouble>& z1 = m_base[1]->GetZ();
609  const Array<OneD, const NekDouble>& z2 = m_base[2]->GetZ();
610 
611  // set up geometric factor: (1+z0)/2
612  for(i = 0; i < nquad0; ++i)
613  {
614  gfac0[i] = 0.5*(1+z0[i]);
615  }
616 
617  // set up geometric factor: (1+z1)/2
618  for(i = 0; i < nquad1; ++i)
619  {
620  gfac1[i] = 0.5*(1+z1[i]);
621  }
622 
623  // Set up geometric factor: 2/(1-z2)
624  for(i = 0; i < nquad2; ++i)
625  {
626  gfac2[i] = 2.0/(1-z2[i]);
627  }
628 
629  // Derivative in first/second direction is always scaled as follows
630  const int nq01 = nquad0*nquad1;
631  for(i = 0; i < nquad2; ++i)
632  {
633  Vmath::Smul(nq01, gfac2[i],
634  &inarray[0] + i*nq01, 1,
635  &tmp0 [0] + i*nq01, 1);
636  }
637 
639 
640  switch(dir)
641  {
642  case 0:
643  {
644  IProductWRTBase_SumFacKernel(m_base[0]->GetDbdata(),
645  m_base[1]->GetBdata (),
646  m_base[2]->GetBdata (),
647  tmp0, outarray, wsp,
648  false, true, true);
649  break;
650  }
651  case 1:
652  {
653  IProductWRTBase_SumFacKernel(m_base[0]->GetBdata (),
654  m_base[1]->GetDbdata(),
655  m_base[2]->GetBdata (),
656  tmp0, outarray, wsp,
657  true, false, true);
658  break;
659  }
660  case 2:
661  {
662  Array<OneD, NekDouble> tmp3(m_ncoeffs);
663  Array<OneD, NekDouble> tmp4(m_ncoeffs);
664 
665  // Scale eta_1 derivative by gfac0
666  for(i = 0; i < nquad1*nquad2; ++i)
667  {
668  Vmath::Vmul(nquad0, tmp0 .get() + i*nquad0, 1,
669  gfac0.get(), 1,
670  tmp0 .get() + i*nquad0, 1);
671  }
672  IProductWRTBase_SumFacKernel(m_base[0]->GetDbdata(),
673  m_base[1]->GetBdata(),
674  m_base[2]->GetBdata(),
675  tmp0, tmp3, wsp,
676  false, true, true);
677 
678  // Scale eta_2 derivative by gfac1*gfac2
679  for(i = 0; i < nquad2; ++i)
680  {
681  Vmath::Smul(nq01, gfac2[i],
682  &inarray[0] + i*nq01, 1,
683  &tmp0 [0] + i*nq01, 1);
684  }
685  for(i = 0; i < nquad1*nquad2; ++i)
686  {
687  Vmath::Smul(nquad0, gfac1[i%nquad1],
688  &tmp0[0] + i*nquad0, 1,
689  &tmp0[0] + i*nquad0, 1);
690  }
691 
693  IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(),
694  m_base[1]->GetDbdata(),
695  m_base[2]->GetBdata(),
696  tmp0, tmp4, wsp,
697  true, false, true);
698 
699  v_MultiplyByStdQuadratureMetric(inarray,tmp0);
700  IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(),
701  m_base[1]->GetBdata(),
702  m_base[2]->GetDbdata(),
703  tmp0,outarray,wsp,
704  true, true, false);
705 
706  Vmath::Vadd(m_ncoeffs,&tmp3[0],1,&outarray[0],1,&outarray[0],1);
707  Vmath::Vadd(m_ncoeffs,&tmp4[0],1,&outarray[0],1,&outarray[0],1);
708  break;
709  }
710  default:
711  {
712  ASSERTL1(false, "input dir is out of range");
713  break;
714  }
715  }
716  }
void IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:192
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:322
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.cpp:225

References ASSERTL1, Nektar::StdRegions::StdExpansion3D::IProductWRTBase_SumFacKernel(), Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::m_ncoeffs, Vmath::Smul(), v_MultiplyByStdQuadratureMetric(), Vmath::Vadd(), and Vmath::Vmul().

Referenced by v_IProductWRTDerivBase().

◆ v_LocCollapsedToLocCoord()

void Nektar::StdRegions::StdPyrExp::v_LocCollapsedToLocCoord ( const Array< OneD, const NekDouble > &  eta,
Array< OneD, NekDouble > &  xi 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 743 of file StdPyrExp.cpp.

746  {
747  xi[0] = (1.0 + eta[0]) * (1.0 - eta[2]) * 0.5 - 1.0;
748  xi[1] = (1.0 + eta[1]) * (1.0 - eta[2]) * 0.5 - 1.0;
749  xi[2] = eta[2];
750  }

◆ v_LocCoordToLocCollapsed()

void Nektar::StdRegions::StdPyrExp::v_LocCoordToLocCollapsed ( const Array< OneD, const NekDouble > &  xi,
Array< OneD, NekDouble > &  eta 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 722 of file StdPyrExp.cpp.

725  {
726  NekDouble d2 = 1.0 - xi[2];
727  if(fabs(d2) < NekConstants::kNekZeroTol)
728  {
729  if(d2>=0.)
730  {
732  }
733  else
734  {
736  }
737  }
738  eta[2] = xi[2]; // eta_z = xi_z
739  eta[1] = 2.0*(1.0 + xi[1])/d2 - 1.0;
740  eta[0] = 2.0*(1.0 + xi[0])/d2 - 1.0;
741  }
static const NekDouble kNekZeroTol
double NekDouble

References Nektar::NekConstants::kNekZeroTol.

◆ v_MultiplyByStdQuadratureMetric()

void Nektar::StdRegions::StdPyrExp::v_MultiplyByStdQuadratureMetric ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1924 of file StdPyrExp.cpp.

1927  {
1928  int i, j;
1929 
1930  int nquad0 = m_base[0]->GetNumPoints();
1931  int nquad1 = m_base[1]->GetNumPoints();
1932  int nquad2 = m_base[2]->GetNumPoints();
1933 
1934  const Array<OneD, const NekDouble>& w0 = m_base[0]->GetW();
1935  const Array<OneD, const NekDouble>& w1 = m_base[1]->GetW();
1936  const Array<OneD, const NekDouble>& w2 = m_base[2]->GetW();
1937 
1938  const Array<OneD, const NekDouble>& z2 = m_base[2]->GetZ();
1939 
1940  // Multiply by integration constants in x-direction
1941  for(i = 0; i < nquad1*nquad2; ++i)
1942  {
1943  Vmath::Vmul(nquad0, inarray.get()+i*nquad0, 1,
1944  w0.get(), 1, outarray.get()+i*nquad0,1);
1945  }
1946 
1947  // Multiply by integration constants in y-direction
1948  for(j = 0; j < nquad2; ++j)
1949  {
1950  for(i = 0; i < nquad1; ++i)
1951  {
1952  Blas::Dscal(nquad0,w1[i], &outarray[0]+i*nquad0 +
1953  j*nquad0*nquad1,1);
1954  }
1955  }
1956 
1957  // Multiply by integration constants in z-direction; need to
1958  // incorporate factor [(1-eta_3)/2]^2 into weights, but only if
1959  // using GLL quadrature points.
1960  switch(m_base[2]->GetPointsType())
1961  {
1962  // (2,0) Jacobi inner product.
1964  for(i = 0; i < nquad2; ++i)
1965  {
1966  Blas::Dscal(nquad0*nquad1, 0.25*w2[i],
1967  &outarray[0]+i*nquad0*nquad1, 1);
1968  }
1969  break;
1970 
1971  default:
1972  for(i = 0; i < nquad2; ++i)
1973  {
1974  Blas::Dscal(nquad0*nquad1,0.25*(1-z2[i])*(1-z2[i])*w2[i],
1975  &outarray[0]+i*nquad0*nquad1,1);
1976  }
1977  break;
1978  }
1979  }
LibUtilities::PointsType GetPointsType(const int dir) const
This function returns the type of quadrature points used in the dir direction.
Definition: StdExpansion.h:208
static void Dscal(const int &n, const double &alpha, double *x, const int &incx)
BLAS level 1: x = alpha x.
Definition: Blas.hpp:182
@ eGaussRadauMAlpha2Beta0
Gauss Radau pinned at x=-1, .
Definition: PointsType.h:59

References Blas::Dscal(), Nektar::LibUtilities::eGaussRadauMAlpha2Beta0, Nektar::StdRegions::StdExpansion::GetPointsType(), Nektar::StdRegions::StdExpansion::m_base, and Vmath::Vmul().

Referenced by v_IProductWRTBase(), v_IProductWRTBase_SumFac(), and v_IProductWRTDerivBase_SumFac().

◆ v_NumBndryCoeffs()

int Nektar::StdRegions::StdPyrExp::v_NumBndryCoeffs ( ) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 911 of file StdPyrExp.cpp.

912  {
915  "BasisType is not a boundary interior form");
918  "BasisType is not a boundary interior form");
921  "BasisType is not a boundary interior form");
922 
923  int P = m_base[0]->GetNumModes();
924  int Q = m_base[1]->GetNumModes();
925  int R = m_base[2]->GetNumModes();
926 
929  }
int getNumberOfBndCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:267

References ASSERTL1, Nektar::LibUtilities::eGLL_Lagrange, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModifiedPyr_C, Nektar::StdRegions::StdExpansion::GetBasisType(), Nektar::LibUtilities::StdPyrData::getNumberOfBndCoefficients(), Nektar::StdRegions::StdExpansion::m_base, and main::P.

◆ v_NumDGBndryCoeffs()

int Nektar::StdRegions::StdPyrExp::v_NumDGBndryCoeffs ( ) const
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 931 of file StdPyrExp.cpp.

932  {
935  "BasisType is not a boundary interior form");
938  "BasisType is not a boundary interior form");
941  "BasisType is not a boundary interior form");
942 
943  int P = m_base[0]->GetNumModes()-1;
944  int Q = m_base[1]->GetNumModes()-1;
945  int R = m_base[2]->GetNumModes()-1;
946 
947  return (P+1)*(Q+1) // 1 rect. face on base
948  + 2*(R+1) + P*(1 + 2*R - P) // 2 tri. (P,R) faces
949  + 2*(R+1) + Q*(1 + 2*R - Q); // 2 tri. (Q,R) faces
950  }

References ASSERTL1, Nektar::LibUtilities::eGLL_Lagrange, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModifiedPyr_C, Nektar::StdRegions::StdExpansion::GetBasisType(), Nektar::StdRegions::StdExpansion::m_base, and main::P.

◆ v_PhysDeriv() [1/2]

void Nektar::StdRegions::StdPyrExp::v_PhysDeriv ( const Array< OneD, const NekDouble > &  u_physical,
Array< OneD, NekDouble > &  out_dxi1,
Array< OneD, NekDouble > &  out_dxi2,
Array< OneD, NekDouble > &  out_dxi3 
)
protectedvirtual

Calculate the derivative of the physical points.

The derivative is evaluated at the nodal physical points. Derivatives with respect to the local Cartesian coordinates.

\(\begin{Bmatrix} \frac {\partial} {\partial \xi_1} \\ \frac {\partial} {\partial \xi_2} \\ \frac {\partial} {\partial \xi_3} \end{Bmatrix} = \begin{Bmatrix} \frac 2 {(1-\eta_3)} \frac \partial {\partial \bar \eta_1} \\ \frac {\partial} {\partial \xi_2} \ \ \frac {(1 + \bar \eta_1)} {(1 - \eta_3)} \frac \partial {\partial \bar \eta_1} + \frac {\partial} {\partial \eta_3} \end{Bmatrix}\)

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 106 of file StdPyrExp.cpp.

111  {
112  // PhysDerivative implementation based on Spen's book page 152.
113  int Qx = m_base[0]->GetNumPoints();
114  int Qy = m_base[1]->GetNumPoints();
115  int Qz = m_base[2]->GetNumPoints();
116 
117  Array<OneD, NekDouble> dEta_bar1(Qx*Qy*Qz,0.0);
118  Array<OneD, NekDouble> dXi2 (Qx*Qy*Qz,0.0);
119  Array<OneD, NekDouble> dEta3 (Qx*Qy*Qz,0.0);
120  PhysTensorDeriv(u_physical, dEta_bar1, dXi2, dEta3);
121 
122  Array<OneD, const NekDouble> eta_x, eta_y, eta_z;
123  eta_x = m_base[0]->GetZ();
124  eta_y = m_base[1]->GetZ();
125  eta_z = m_base[2]->GetZ();
126 
127  int i, j, k, n;
128 
129  if (out_dxi1.size() > 0)
130  {
131  for (k = 0, n = 0; k < Qz; ++k)
132  {
133  NekDouble fac = 2.0/(1.0 - eta_z[k]);
134  for (j = 0; j < Qy; ++j)
135  {
136  for (i = 0; i < Qx; ++i, ++n)
137  {
138  out_dxi1[n] = fac * dEta_bar1[n];
139  }
140  }
141  }
142  }
143 
144  if (out_dxi2.size() > 0)
145  {
146  for (k = 0, n = 0; k < Qz; ++k)
147  {
148  NekDouble fac = 2.0/(1.0 - eta_z[k]);
149  for (j = 0; j < Qy; ++j)
150  {
151  for (i = 0; i < Qx; ++i, ++n)
152  {
153  out_dxi2[n] = fac * dXi2[n];
154  }
155  }
156  }
157  }
158 
159  if (out_dxi3.size() > 0)
160  {
161  for (k = 0, n = 0; k < Qz; ++k)
162  {
163  NekDouble fac = 1.0/(1.0 - eta_z[k]);
164  for (j = 0; j < Qy; ++j)
165  {
166  NekDouble fac1 = (1.0+eta_y[j]);
167  for (i = 0; i < Qx; ++i, ++n)
168  {
169  out_dxi3[n] = (1.0+eta_x[i])*fac*dEta_bar1[n] +
170  fac1*fac*dXi2[n] + dEta3[n];
171  }
172  }
173  }
174  }
175  }
void PhysTensorDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray_d1, Array< OneD, NekDouble > &outarray_d2, Array< OneD, NekDouble > &outarray_d3)
Calculate the 3D derivative in the local tensor/collapsed coordinate at the physical points.

References Nektar::StdRegions::StdExpansion::m_base, and Nektar::StdRegions::StdExpansion3D::PhysTensorDeriv().

Referenced by v_PhysDeriv(), and v_StdPhysDeriv().

◆ v_PhysDeriv() [2/2]

void Nektar::StdRegions::StdPyrExp::v_PhysDeriv ( const int  dir,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  out_d0 
)
protectedvirtual

Calculate the derivative of the physical points in a given direction.

See also
StdRegions::StdExpansion::PhysDeriv

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 177 of file StdPyrExp.cpp.

180  {
181  switch(dir)
182  {
183  case 0:
184  {
185  v_PhysDeriv(inarray, outarray, NullNekDouble1DArray,
187  break;
188  }
189 
190  case 1:
191  {
192  v_PhysDeriv(inarray, NullNekDouble1DArray, outarray,
194  break;
195  }
196 
197  case 2:
198  {
200  NullNekDouble1DArray, outarray);
201  break;
202  }
203 
204  default:
205  {
206  ASSERTL1(false,"input dir is out of range");
207  }
208  break;
209  }
210  }
virtual void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
Calculate the derivative of the physical points.
Definition: StdPyrExp.cpp:106
static Array< OneD, NekDouble > NullNekDouble1DArray

References ASSERTL1, Nektar::NullNekDouble1DArray, and v_PhysDeriv().

◆ v_PhysEvaluateBasis()

NekDouble Nektar::StdRegions::StdPyrExp::v_PhysEvaluateBasis ( const Array< OneD, const NekDouble > &  coords,
int  mode 
)
finalprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 827 of file StdPyrExp.cpp.

830  {
831  Array<OneD, NekDouble> coll(3);
832  LocCoordToLocCollapsed(coords, coll);
833 
834  const int nm0 = m_base[0]->GetNumModes();
835  const int nm1 = m_base[1]->GetNumModes();
836  const int nm2 = m_base[2]->GetNumModes();
837 
838  int mode0 = 0, mode1 = 0, mode2 = 0, cnt = 0;
839 
840  bool found = false;
841  for (mode0 = 0; mode0 < nm0; ++mode0)
842  {
843  for (mode1 = 0; mode1 < nm1; ++mode1)
844  {
845  int maxpq = max(mode0, mode1);
846  for (mode2 = 0; mode2 < nm2 - maxpq; ++mode2, ++cnt)
847  {
848  if (cnt == mode)
849  {
850  found = true;
851  break;
852  }
853  }
854 
855  if (found)
856  {
857  break;
858  }
859  }
860 
861  if (found)
862  {
863  break;
864  }
865 
866  for (int j = nm1; j < nm2; ++j)
867  {
868  int ijmax = max(mode0, j);
869  mode2 += nm2 - ijmax;
870  }
871  }
872 
873  if (mode == 1 &&
875  {
876  return StdExpansion::BaryEvaluateBasis<2>(coll[2], 1);
877  }
878  else
879  {
880  return
881  StdExpansion::BaryEvaluateBasis<0>(coll[0], mode0) *
882  StdExpansion::BaryEvaluateBasis<1>(coll[1], mode1) *
883  StdExpansion::BaryEvaluateBasis<2>(coll[2], mode2);
884  }
885  }
void LocCoordToLocCollapsed(const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta)
Convert local cartesian coordinate xi into local collapsed coordinates eta.
Definition: StdExpansion.h:982

References Nektar::LibUtilities::eModified_A, Nektar::StdRegions::StdExpansion::GetBasisType(), Nektar::StdRegions::StdExpansion::LocCoordToLocCollapsed(), and Nektar::StdRegions::StdExpansion::m_base.

◆ v_ReduceOrderCoeffs()

void Nektar::StdRegions::StdPyrExp::v_ReduceOrderCoeffs ( int  numMin,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 2117 of file StdPyrExp.cpp.

2121  {
2122  int nquad0 = m_base[0]->GetNumPoints();
2123  int nquad1 = m_base[1]->GetNumPoints();
2124  int nquad2 = m_base[2]->GetNumPoints();
2125  int nqtot = nquad0*nquad1*nquad2;
2126  int nmodes0 = m_base[0]->GetNumModes();
2127  int nmodes1 = m_base[1]->GetNumModes();
2128  int nmodes2 = m_base[2]->GetNumModes();
2129  int numMax = nmodes0;
2130 
2131  Array<OneD, NekDouble> coeff (m_ncoeffs);
2132  Array<OneD, NekDouble> coeff_tmp1(m_ncoeffs, 0.0);
2133  Array<OneD, NekDouble> phys_tmp (nqtot, 0.0);
2134  Array<OneD, NekDouble> tmp, tmp2, tmp3, tmp4;
2135 
2136 
2137  const LibUtilities::PointsKey Pkey0 = m_base[0]->GetPointsKey();
2138  const LibUtilities::PointsKey Pkey1 = m_base[1]->GetPointsKey();
2139  const LibUtilities::PointsKey Pkey2 = m_base[2]->GetPointsKey();
2140 
2141  LibUtilities::BasisKey bortho0(
2142  LibUtilities::eOrtho_A, nmodes0, Pkey0);
2143  LibUtilities::BasisKey bortho1(
2144  LibUtilities::eOrtho_A, nmodes1, Pkey1);
2145  LibUtilities::BasisKey bortho2(
2146  LibUtilities::eOrthoPyr_C, nmodes2, Pkey2);
2147 
2148  int cnt = 0;
2149  int u = 0;
2150  int i = 0;
2151  StdRegions::StdPyrExpSharedPtr OrthoPyrExp;
2152 
2154  ::AllocateSharedPtr(bortho0, bortho1, bortho2);
2155 
2156  BwdTrans(inarray,phys_tmp);
2157  OrthoPyrExp->FwdTrans(phys_tmp, coeff);
2158 
2159  // filtering
2160  for (u = 0; u < numMin; ++u)
2161  {
2162  for (i = 0; i < numMin; ++i)
2163  {
2164 
2165  int maxui = max(u,i);
2166  Vmath::Vcopy(numMin - maxui, tmp = coeff + cnt, 1,
2167  tmp2 = coeff_tmp1 + cnt, 1);
2168  cnt += nmodes2 - maxui;
2169  }
2170 
2171  for (i = numMin; i < nmodes1; ++i)
2172  {
2173  int maxui = max(u,i);
2174  cnt += numMax - maxui;
2175  }
2176  }
2177 
2178  OrthoPyrExp->BwdTrans(coeff_tmp1, phys_tmp);
2179  StdPyrExp::FwdTrans(phys_tmp, outarray);
2180  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function performs the Backward transformation from coefficient space to physical space.
Definition: StdExpansion.h:433
void FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function performs the Forward transformation from physical space to coefficient space.
@ eOrtho_A
Principle Orthogonal Functions .
Definition: BasisType.h:45
std::shared_ptr< StdPyrExp > StdPyrExpSharedPtr
Definition: StdPyrExp.h:278

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::StdRegions::StdExpansion::BwdTrans(), Nektar::LibUtilities::eOrtho_A, Nektar::LibUtilities::eOrthoPyr_C, Nektar::StdRegions::StdExpansion::FwdTrans(), Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::m_ncoeffs, and Vmath::Vcopy().

◆ v_StdPhysDeriv() [1/2]

void Nektar::StdRegions::StdPyrExp::v_StdPhysDeriv ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  out_d0,
Array< OneD, NekDouble > &  out_d1,
Array< OneD, NekDouble > &  out_d2 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 212 of file StdPyrExp.cpp.

216  {
217  StdPyrExp::v_PhysDeriv(inarray, out_d0, out_d1, out_d2);
218  }

References v_PhysDeriv().

◆ v_StdPhysDeriv() [2/2]

void Nektar::StdRegions::StdPyrExp::v_StdPhysDeriv ( const int  dir,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 220 of file StdPyrExp.cpp.

223  {
224  StdPyrExp::v_PhysDeriv(dir, inarray, outarray);
225  }

References v_PhysDeriv().

◆ v_SVVLaplacianFilter()

void Nektar::StdRegions::StdPyrExp::v_SVVLaplacianFilter ( Array< OneD, NekDouble > &  array,
const StdMatrixKey mkey 
)
protectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Reimplemented in Nektar::LocalRegions::PyrExp.

Definition at line 1982 of file StdPyrExp.cpp.

1984  {
1985  // Generate an orthonogal expansion
1986  int qa = m_base[0]->GetNumPoints();
1987  int qb = m_base[1]->GetNumPoints();
1988  int qc = m_base[2]->GetNumPoints();
1989  int nmodes_a = m_base[0]->GetNumModes();
1990  int nmodes_b = m_base[1]->GetNumModes();
1991  int nmodes_c = m_base[2]->GetNumModes();
1992  // Declare orthogonal basis.
1993  LibUtilities::PointsKey pa(qa,m_base[0]->GetPointsType());
1994  LibUtilities::PointsKey pb(qb,m_base[1]->GetPointsType());
1995  LibUtilities::PointsKey pc(qc,m_base[2]->GetPointsType());
1996 
1997  LibUtilities::BasisKey Ba(LibUtilities::eOrtho_A,nmodes_a,pa);
1998  LibUtilities::BasisKey Bb(LibUtilities::eOrtho_A,nmodes_b,pb);
1999  LibUtilities::BasisKey Bc(LibUtilities::eOrthoPyr_C,nmodes_c,pc);
2000  StdPyrExp OrthoExp(Ba,Bb,Bc);
2001 
2002  Array<OneD, NekDouble> orthocoeffs(OrthoExp.GetNcoeffs());
2003  int i,j,k,cnt = 0;
2004 
2005  // project onto modal space.
2006  OrthoExp.FwdTrans(array,orthocoeffs);
2007 
2008  if(mkey.ConstFactorExists(eFactorSVVPowerKerDiffCoeff))
2009  {
2010  // Rodrigo's power kernel
2011  NekDouble cutoff = mkey.GetConstFactor(eFactorSVVCutoffRatio);
2012  NekDouble SvvDiffCoeff =
2013  mkey.GetConstFactor(eFactorSVVPowerKerDiffCoeff)*
2014  mkey.GetConstFactor(eFactorSVVDiffCoeff);
2015 
2016  for(int i = 0; i < nmodes_a; ++i)
2017  {
2018  for(int j = 0; j < nmodes_b; ++j)
2019  {
2020  int maxij = max(i,j);
2021  NekDouble fac1 = std::max(
2022  pow((1.0*i)/(nmodes_a-1),cutoff*nmodes_a),
2023  pow((1.0*j)/(nmodes_b-1),cutoff*nmodes_b));
2024 
2025  for(int k = 0; k < nmodes_c-maxij; ++k)
2026  {
2027  NekDouble fac = std::max(fac1,
2028  pow((1.0*k)/(nmodes_c-1),cutoff*nmodes_c));
2029 
2030  orthocoeffs[cnt] *= SvvDiffCoeff * fac;
2031  cnt++;
2032  }
2033  }
2034  }
2035  }
2036  else if(mkey.ConstFactorExists(eFactorSVVDGKerDiffCoeff)) // Rodrigo/Mansoor's DG Kernel
2037  {
2038  NekDouble SvvDiffCoeff =
2039  mkey.GetConstFactor(eFactorSVVDGKerDiffCoeff)*
2040  mkey.GetConstFactor(eFactorSVVDiffCoeff);
2041 
2042  int max_abc = max(nmodes_a-kSVVDGFiltermodesmin,
2043  nmodes_b-kSVVDGFiltermodesmin);
2044  max_abc = max(max_abc, nmodes_c-kSVVDGFiltermodesmin);
2045  // clamp max_abc
2046  max_abc = max(max_abc,0);
2047  max_abc = min(max_abc,kSVVDGFiltermodesmax-kSVVDGFiltermodesmin);
2048 
2049  for(int i = 0; i < nmodes_a; ++i)
2050  {
2051  for(int j = 0; j < nmodes_b; ++j)
2052  {
2053  int maxij = max(i,j);
2054 
2055  for(int k = 0; k < nmodes_c-maxij; ++k)
2056  {
2057  int maxijk = max(maxij,k);
2058  maxijk = min(maxijk,kSVVDGFiltermodesmax-1);
2059 
2060  orthocoeffs[cnt] *= SvvDiffCoeff *
2061  kSVVDGFilter[max_abc][maxijk];
2062  cnt++;
2063  }
2064  }
2065  }
2066  }
2067  else
2068  {
2069  //SVV filter paramaters (how much added diffusion relative
2070  // to physical one and fraction of modes from which you
2071  // start applying this added diffusion)
2072  //
2073  NekDouble SvvDiffCoeff = mkey.GetConstFactor(StdRegions::eFactorSVVDiffCoeff);
2074  NekDouble SVVCutOff = mkey.GetConstFactor(StdRegions::eFactorSVVCutoffRatio);
2075 
2076  //Defining the cut of mode
2077  int cutoff_a = (int) (SVVCutOff*nmodes_a);
2078  int cutoff_b = (int) (SVVCutOff*nmodes_b);
2079  int cutoff_c = (int) (SVVCutOff*nmodes_c);
2080  //To avoid the fac[j] from blowing up
2081  NekDouble epsilon = 1;
2082 
2083  int nmodes = min(min(nmodes_a,nmodes_b),nmodes_c);
2084  NekDouble cutoff = min(min(cutoff_a,cutoff_b),cutoff_c);
2085 
2086  for(i = 0; i < nmodes_a; ++i)//P
2087  {
2088  for(j = 0; j < nmodes_b; ++j) //Q
2089  {
2090  int maxij = max(i,j);
2091  for(k = 0; k < nmodes_c-maxij; ++k) //R
2092  {
2093  if(j + k >= cutoff || i + k >= cutoff)
2094  {
2095  orthocoeffs[cnt] *=
2096  (SvvDiffCoeff*exp(-(i+k-nmodes)*(i+k-nmodes)/
2097  ((NekDouble)((i+k-cutoff+epsilon)*
2098  (i+k-cutoff+epsilon))))*
2099  exp(-(j-nmodes)*(j-nmodes)/
2100  ((NekDouble)((j-cutoff+epsilon)*
2101  (j-cutoff+epsilon)))));
2102  }
2103  else
2104  {
2105  orthocoeffs[cnt] *= 0.0;
2106  }
2107  cnt++;
2108  }
2109  }
2110  }
2111  }
2112 
2113  // backward transform to physical space
2114  OrthoExp.BwdTrans(orthocoeffs,array);
2115  }
const int kSVVDGFiltermodesmin
Definition: StdRegions.hpp:388
const int kSVVDGFiltermodesmax
Definition: StdRegions.hpp:389
const NekDouble kSVVDGFilter[9][11]
Definition: StdRegions.hpp:391

References Nektar::StdRegions::StdExpansion::BwdTrans(), Nektar::StdRegions::StdMatrixKey::ConstFactorExists(), Nektar::StdRegions::eFactorSVVCutoffRatio, Nektar::StdRegions::eFactorSVVDGKerDiffCoeff, Nektar::StdRegions::eFactorSVVDiffCoeff, Nektar::StdRegions::eFactorSVVPowerKerDiffCoeff, Nektar::LibUtilities::eOrtho_A, Nektar::LibUtilities::eOrthoPyr_C, Nektar::StdRegions::StdExpansion::FwdTrans(), Nektar::StdRegions::StdMatrixKey::GetConstFactor(), Nektar::StdRegions::StdExpansion::GetNcoeffs(), Nektar::StdRegions::StdExpansion::GetPointsType(), Nektar::StdRegions::kSVVDGFilter, Nektar::StdRegions::kSVVDGFiltermodesmax, Nektar::StdRegions::kSVVDGFiltermodesmin, and Nektar::StdRegions::StdExpansion::m_base.