Nektar++
Public Member Functions | Protected Member Functions | Private Member Functions | Private Attributes | List of all members
Nektar::LocalRegions::TetExp Class Reference

#include <TetExp.h>

Inheritance diagram for Nektar::LocalRegions::TetExp:
[legend]

Public Member Functions

 TetExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc, const SpatialDomains::TetGeomSharedPtr &geom)
 Constructor using BasisKey class for quadrature points and order definition. More...
 
 TetExp (const TetExp &T)
 Copy Constructor. More...
 
 ~TetExp () override=default
 
- Public Member Functions inherited from Nektar::StdRegions::StdTetExp
 StdTetExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc)
 
 StdTetExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc, NekDouble *coeffs, NekDouble *phys)
 
 StdTetExp (const StdTetExp &T)=default
 
 ~StdTetExp () override=default
 
LibUtilities::ShapeType DetShapeType () const
 
- Public Member Functions inherited from Nektar::StdRegions::StdExpansion3D
 StdExpansion3D (int numcoeffs, const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc)
 
 StdExpansion3D ()=default
 
 StdExpansion3D (const StdExpansion3D &T)=default
 
 ~StdExpansion3D () override=default
 
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, bool UseGLL=false) 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 () const
 
std::shared_ptr< StdExpansionGetLinStdExp (void) const
 
int GetShapeDimension () const
 
bool IsBoundaryInteriorExpansion () const
 
bool IsNodalNonTensorialExp ()
 
void NodalToModal (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
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 FwdTransBndConstrained (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 GetTraceCoeffMap (const unsigned int traceid, Array< OneD, unsigned int > &maparray)
 
void GetElmtTraceToTraceMap (const unsigned 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 LinearAdvectionMatrixOp (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, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs)
 This function evaluates the first derivative of the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble PhysEvaluate (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs, std::array< NekDouble, 6 > &secondOrderDerivs)
 
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...
 
void PhysInterp (std::shared_ptr< StdExpansion > fromExp, const Array< OneD, const NekDouble > &fromData, Array< OneD, NekDouble > &toData)
 interpolate from one set of quadrature points available from FromExp to the set of quadrature points in the current expansion. If the points are the same this routine will just copy the data More...
 
virtual int v_CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset)
 
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)
 
- Public Member Functions inherited from Nektar::LocalRegions::Expansion3D
 Expansion3D (SpatialDomains::Geometry3DSharedPtr pGeom)
 
 ~Expansion3D () override=default
 
void SetTraceToGeomOrientation (Array< OneD, NekDouble > &inout)
 Align trace orientation with the geometry orientation. More...
 
void SetFaceToGeomOrientation (const int face, Array< OneD, NekDouble > &inout)
 Align face orientation with the geometry orientation. More...
 
void AddHDGHelmholtzFaceTerms (const NekDouble tau, const int edge, Array< OneD, NekDouble > &facePhys, const StdRegions::VarCoeffMap &dirForcing, Array< OneD, NekDouble > &outarray)
 
void AddNormTraceInt (const int dir, Array< OneD, ExpansionSharedPtr > &FaceExp, Array< OneD, Array< OneD, NekDouble > > &faceCoeffs, Array< OneD, NekDouble > &outarray)
 
void AddNormTraceInt (const int dir, Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &FaceExp, Array< OneD, NekDouble > &outarray, const StdRegions::VarCoeffMap &varcoeffs)
 
void AddFaceBoundaryInt (const int face, ExpansionSharedPtr &FaceExp, Array< OneD, NekDouble > &facePhys, Array< OneD, NekDouble > &outarray, const StdRegions::VarCoeffMap &varcoeffs=StdRegions::NullVarCoeffMap)
 
SpatialDomains::Geometry3DSharedPtr GetGeom3D () const
 
void v_ReOrientTracePhysMap (const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1) override
 
void v_NormVectorIProductWRTBase (const Array< OneD, const Array< OneD, NekDouble > > &Fvec, Array< OneD, NekDouble > &outarray) override
 
Array< OneD, unsigned int > GetEdgeInverseBoundaryMap (int eid)
 
Array< OneD, unsigned int > GetTraceInverseBoundaryMap (int fid, StdRegions::Orientation faceOrient=StdRegions::eNoOrientation, int P1=-1, int P2=-1)
 
void GetInverseBoundaryMaps (Array< OneD, unsigned int > &vmap, Array< OneD, Array< OneD, unsigned int > > &emap, Array< OneD, Array< OneD, unsigned int > > &fmap)
 
DNekScalMatSharedPtr CreateMatrix (const MatrixKey &mkey)
 
- Public Member Functions inherited from Nektar::LocalRegions::Expansion
 Expansion (SpatialDomains::GeometrySharedPtr pGeom)
 
 Expansion (const Expansion &pSrc)
 
 ~Expansion () override
 
void SetTraceExp (const int traceid, ExpansionSharedPtr &f)
 
ExpansionSharedPtr GetTraceExp (const int traceid)
 
DNekScalMatSharedPtr GetLocMatrix (const LocalRegions::MatrixKey &mkey)
 
void DropLocMatrix (const LocalRegions::MatrixKey &mkey)
 
DNekScalMatSharedPtr GetLocMatrix (const StdRegions::MatrixType mtype, const StdRegions::ConstFactorMap &factors=StdRegions::NullConstFactorMap, const StdRegions::VarCoeffMap &varcoeffs=StdRegions::NullVarCoeffMap)
 
SpatialDomains::GeometrySharedPtr GetGeom () const
 
void Reset ()
 
IndexMapValuesSharedPtr CreateIndexMap (const IndexMapKey &ikey)
 
DNekScalBlkMatSharedPtr CreateStaticCondMatrix (const MatrixKey &mkey)
 
const SpatialDomains::GeomFactorsSharedPtrGetMetricInfo () const
 
DNekMatSharedPtr BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
 
DNekMatSharedPtr BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd)
 
void ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
 
void AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
void AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
void AddFaceNormBoundaryInt (const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
void DGDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &coeffs, Array< OneD, NekDouble > &outarray)
 
NekDouble VectorFlux (const Array< OneD, Array< OneD, NekDouble > > &vec)
 
void NormalTraceDerivFactors (Array< OneD, Array< OneD, NekDouble > > &factors, Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors)
 
IndexMapValuesSharedPtr GetIndexMap (const IndexMapKey &ikey)
 
void AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
ExpansionSharedPtr GetLeftAdjacentElementExp () const
 
ExpansionSharedPtr GetRightAdjacentElementExp () const
 
int GetLeftAdjacentElementTrace () const
 
int GetRightAdjacentElementTrace () const
 
void SetAdjacentElementExp (int traceid, ExpansionSharedPtr &e)
 
StdRegions::Orientation GetTraceOrient (int trace)
 
void SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void DivideByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 Divided by the metric jacobi and quadrature weights. More...
 
void GetTraceQFactors (const int trace, Array< OneD, NekDouble > &outarray)
 Extract the metric factors to compute the contravariant fluxes along edge edge and stores them into outarray following the local edge orientation (i.e. anticlockwise convention). More...
 
void GetTracePhysVals (const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient=StdRegions::eNoOrientation)
 
void GetTracePhysMap (const int edge, Array< OneD, int > &outarray)
 
void ReOrientTracePhysMap (const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1)
 
const NormalVectorGetTraceNormal (const int id)
 
void ComputeTraceNormal (const int id)
 
const Array< OneD, const NekDouble > & GetPhysNormals (void)
 
void SetPhysNormals (Array< OneD, const NekDouble > &normal)
 
void SetUpPhysNormals (const int trace)
 
void AddRobinMassMatrix (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
 
void TraceNormLen (const int traceid, NekDouble &h, NekDouble &p)
 
void AddRobinTraceContribution (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
 
const Array< OneD, const NekDouble > & GetElmtBndNormDirElmtLen (const int nbnd) const
 
void StdDerivBaseOnTraceMat (Array< OneD, DNekMatSharedPtr > &DerivMat)
 

Protected Member Functions

NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray) override
 Integrate the physical point list inarray over region. More...
 
void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2) override
 Differentiate inarray in the three coordinate directions. More...
 
void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Forward transform from physical quadrature space stored in inarray and evaluate the expansion coefficients and store in (this)->_coeffs. More...
 
void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculate the inner product of inarray with respect to the basis B=m_base0*m_base1*m_base2 and put into outarray: More...
 
void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculates the inner product \( I_{pqr} = (u, \partial_{x_i} \phi_{pqr}) \). More...
 
void v_AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals) override
 
NekDouble v_PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals) override
 
NekDouble v_PhysEvalFirstDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) override
 
void v_GetCoord (const Array< OneD, const NekDouble > &Lcoords, Array< OneD, NekDouble > &coords) override
 Get the coordinates "coords" at the local coordinates "Lcoords". More...
 
void v_GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
 
LibUtilities::ShapeType v_DetShapeType () const override
 Return Shape of region, using ShapeType enum list. More...
 
StdRegions::StdExpansionSharedPtr v_GetStdExp (void) const override
 
StdRegions::StdExpansionSharedPtr v_GetLinStdExp (void) const override
 
void v_ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int mode_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType) override
 
void v_GetTracePhysMap (const int face, Array< OneD, int > &outarray) override
 Returns the physical values at the quadrature points of a face. More...
 
void v_ComputeTraceNormal (const int face) override
 Compute the normal of a triangular face. More...
 
void v_HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_LaplacianMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_LaplacianMatrixOp (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdRegions::StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_CreateStdMatrix (const StdRegions::StdMatrixKey &mkey) override
 
DNekScalMatSharedPtr v_GetLocMatrix (const MatrixKey &mkey) override
 
void v_DropLocMatrix (const MatrixKey &mkey) override
 
DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix (const MatrixKey &mkey) override
 
void v_DropLocStaticCondMatrix (const MatrixKey &mkey) override
 
void SetUpInverseTransformationMatrix (const DNekMatSharedPtr &m_transformationmatrix, DNekMatSharedPtr m_inversetransformationmatrix, DNekMatSharedPtr m_inversetransposedtransformationmatrix)
 
void v_ComputeLaplacianMetric () override
 
void v_NormalTraceDerivFactors (Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors, Array< OneD, Array< OneD, NekDouble > > &d2factors) override
 : This method gets all of the factors which are required as part of the Gradient Jump Penalty stabilisation and involves the product of the normal and geometric factors along the element trace. More...
 
- Protected Member Functions inherited from Nektar::StdRegions::StdTetExp
void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_dx, Array< OneD, NekDouble > &out_dy, Array< OneD, NekDouble > &out_dz) override
 Calculate the derivative of the physical points. More...
 
void v_PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_StdPhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2) override
 
void v_StdPhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
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) override
 
void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
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) override
 
void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_LocCoordToLocCollapsed (const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta) override
 
void v_LocCollapsedToLocCoord (const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi) override
 
void v_GetCoords (Array< OneD, NekDouble > &coords_x, Array< OneD, NekDouble > &coords_y, Array< OneD, NekDouble > &coords_z) override
 
void v_FillMode (const int mode, Array< OneD, NekDouble > &outarray) override
 
NekDouble v_PhysEvaluateBasis (const Array< OneD, const NekDouble > &coords, int mode) final
 
NekDouble v_PhysEvalFirstDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) override
 
void v_GetTraceNumModes (const int fid, int &numModes0, int &numModes1, Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2) override
 
int v_GetNverts () const override
 
int v_GetNedges () const override
 
int v_GetNtraces () const override
 
LibUtilities::ShapeType v_DetShapeType () const override
 
int v_NumBndryCoeffs () const override
 
int v_NumDGBndryCoeffs () const override
 
int v_GetTraceNcoeffs (const int i) const override
 
int v_GetTraceIntNcoeffs (const int i) const override
 
int v_GetTraceNumPoints (const int i) const override
 
int v_GetEdgeNcoeffs (const int i) const override
 
LibUtilities::PointsKey v_GetTracePointsKey (const int i, const int j) const override
 
int v_CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset) override
 
const LibUtilities::BasisKey v_GetTraceBasisKey (const int i, const int k, bool UseGLL=false) const override
 
bool v_IsBoundaryInteriorExpansion () const override
 
int v_GetVertexMap (int localVertexId, bool useCoeffPacking=false) override
 
void v_GetInteriorMap (Array< OneD, unsigned int > &outarray) override
 
void v_GetBoundaryMap (Array< OneD, unsigned int > &outarray) override
 
void v_GetTraceCoeffMap (const unsigned int fid, Array< OneD, unsigned int > &maparray) override
 
void v_GetElmtTraceToTraceMap (const unsigned int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1) override
 
void v_GetEdgeInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2) override
 
void v_GetTraceInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2) override
 
DNekMatSharedPtr v_GenMatrix (const StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_CreateStdMatrix (const StdMatrixKey &mkey) override
 
void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdMatrixKey &mkey) override
 
void v_ReduceOrderCoeffs (int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_GetSimplexEquiSpacedConnectivity (Array< OneD, int > &conn, bool standard=true) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion3D
NekDouble v_PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals) override
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble v_PhysEvaluateInterp (const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals) override
 
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)=0
 
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)=0
 
void v_LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray) override
 Integrates the specified function over the domain. More...
 
virtual int v_GetNedges (void) const
 
virtual int v_GetEdgeNcoeffs (const int i) const
 
NekDouble BaryTensorDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs)
 
virtual void v_GetEdgeInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards)
 
void v_GetTraceToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient, int P, int Q) override
 
void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat) override
 
void v_PhysInterp (std::shared_ptr< StdExpansion > fromExp, const Array< OneD, const NekDouble > &fromData, Array< OneD, NekDouble > &toData) override
 
- 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_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 LinearAdvectionMatrixOp_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 NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
 
virtual void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat)
 
virtual void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals, NekDouble &deriv, NekDouble &deriv2)
 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)
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals)
 Helper function to pass an unused value by reference into BaryEvaluate. More...
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals, NekDouble &deriv)
 
- Protected Member Functions inherited from Nektar::LocalRegions::Expansion3D
void v_DGDeriv (const int dir, const Array< OneD, const NekDouble > &incoeffs, Array< OneD, ExpansionSharedPtr > &FaceExp, Array< OneD, Array< OneD, NekDouble > > &faceCoeffs, Array< OneD, NekDouble > &out_d) override
 Evaluate coefficients of weak deriviative in the direction dir given the input coefficicents incoeffs and the imposed boundary values in EdgeExp (which will have its phys space updated). More...
 
DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
void v_AddFaceNormBoundaryInt (const int face, const ExpansionSharedPtr &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray) override
 
void v_AddRobinMassMatrix (const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat) override
 
StdRegions::Orientation v_GetTraceOrient (int face) override
 
void v_GetTracePhysVals (const int face, const StdRegions::StdExpansionSharedPtr &FaceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient) override
 Extract the physical values along face face from inarray into outarray following the local face orientation and point distribution defined by defined in FaceExp. More...
 
void v_GenTraceExp (const int traceid, ExpansionSharedPtr &exp) override
 
void GetPhysFaceVarCoeffsFromElement (const int face, ExpansionSharedPtr &FaceExp, const Array< OneD, const NekDouble > &varcoeff, Array< OneD, NekDouble > &outarray)
 
DNekMatSharedPtr v_BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType) override
 
DNekMatSharedPtr v_BuildInverseTransformationMatrix (const DNekScalMatSharedPtr &transformationmatrix) override
 Build inverse and inverse transposed transformation matrix: \(\mathbf{R^{-1}}\) and \(\mathbf{R^{-T}}\). More...
 
DNekMatSharedPtr v_BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd) override
 
void v_TraceNormLen (const int traceid, NekDouble &h, NekDouble &p) override
 
- Protected Member Functions inherited from Nektar::LocalRegions::Expansion
void ComputeLaplacianMetric ()
 
void ComputeQuadratureMetric ()
 
void ComputeGmatcdotMF (const Array< TwoD, const NekDouble > &df, const Array< OneD, const NekDouble > &direction, Array< OneD, Array< OneD, NekDouble > > &dfdir)
 
Array< OneD, NekDoubleGetMF (const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
 
Array< OneD, NekDoubleGetMFDiv (const int dir, const StdRegions::VarCoeffMap &varcoeffs)
 
Array< OneD, NekDoubleGetMFMag (const int dir, const StdRegions::VarCoeffMap &varcoeffs)
 
void v_MultiplyByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_DivideByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_ComputeLaplacianMetric ()
 
int v_GetCoordim () const override
 
void v_GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
 
virtual DNekScalMatSharedPtr v_GetLocMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual void v_DropLocMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual DNekMatSharedPtr v_BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
 
virtual DNekMatSharedPtr v_BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd)
 
virtual void v_ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
 
virtual void v_AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddFaceNormBoundaryInt (const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
virtual void v_DGDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &coeffs, Array< OneD, NekDouble > &outarray)
 
virtual NekDouble v_VectorFlux (const Array< OneD, Array< OneD, NekDouble > > &vec)
 
virtual void v_NormalTraceDerivFactors (Array< OneD, Array< OneD, NekDouble > > &factors, Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors)
 
virtual void v_AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual StdRegions::Orientation v_GetTraceOrient (int trace)
 
void v_SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_GetTraceQFactors (const int trace, Array< OneD, NekDouble > &outarray)
 
virtual void v_GetTracePhysVals (const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient)
 
virtual void v_GetTracePhysMap (const int edge, Array< OneD, int > &outarray)
 
virtual void v_ReOrientTracePhysMap (const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1=-1)
 
virtual void v_ComputeTraceNormal (const int id)
 
virtual const Array< OneD, const NekDouble > & v_GetPhysNormals ()
 
virtual void v_SetPhysNormals (Array< OneD, const NekDouble > &normal)
 
virtual void v_SetUpPhysNormals (const int id)
 
virtual void v_AddRobinMassMatrix (const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
 
virtual void v_AddRobinTraceContribution (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
 
virtual void v_TraceNormLen (const int traceid, NekDouble &h, NekDouble &p)
 
virtual void v_GenTraceExp (const int traceid, ExpansionSharedPtr &exp)
 

Private Member Functions

void GeneralMatrixOp_MatOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
 
void v_LaplacianMatrixOp_MatFree_Kernel (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp) override
 

Private Attributes

LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLessm_matrixManager
 
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLessm_staticCondMatrixManager
 

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
 
- Protected Attributes inherited from Nektar::LocalRegions::Expansion3D
std::map< int, NormalVectorm_faceNormals
 
- Protected Attributes inherited from Nektar::LocalRegions::Expansion
LibUtilities::NekManager< IndexMapKey, IndexMapValues, IndexMapKey::opLessm_indexMapManager
 
std::map< int, ExpansionWeakPtrm_traceExp
 
SpatialDomains::GeometrySharedPtr m_geom
 
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
 
MetricMap m_metrics
 
std::map< int, NormalVectorm_traceNormals
 
ExpansionWeakPtr m_elementLeft
 
ExpansionWeakPtr m_elementRight
 
int m_elementTraceLeft = -1
 
int m_elementTraceRight = -1
 
std::map< int, Array< OneD, NekDouble > > m_elmtBndNormDirElmtLen
 the element length in each element boundary(Vertex, edge or face) normal direction calculated based on the local m_metricinfo times the standard element length (which is 2.0) More...
 

Detailed Description

Defines a Tetrahedral local expansion.

Definition at line 48 of file TetExp.h.

Constructor & Destructor Documentation

◆ TetExp() [1/2]

Nektar::LocalRegions::TetExp::TetExp ( const LibUtilities::BasisKey Ba,
const LibUtilities::BasisKey Bb,
const LibUtilities::BasisKey Bc,
const SpatialDomains::TetGeomSharedPtr geom 
)

Constructor using BasisKey class for quadrature points and order definition.

Parameters
BaBasis key for first coordinate.
BbBasis key for second coordinate.
BcBasis key for third coordinate.

Definition at line 57 of file TetExp.cpp.

62 Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes()),
63 3, Ba, Bb, Bc),
65 Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes()),
66 Ba, Bb, Bc),
67 StdTetExp(Ba, Bb, Bc), Expansion(geom), Expansion3D(geom),
69 std::bind(&Expansion3D::CreateMatrix, this, std::placeholders::_1),
70 std::string("TetExpMatrix")),
72 this, std::placeholders::_1),
73 std::string("TetExpStaticCondMatrix"))
74{
75}
Expansion3D(SpatialDomains::Geometry3DSharedPtr pGeom)
Definition: Expansion3D.h:59
DNekScalMatSharedPtr CreateMatrix(const MatrixKey &mkey)
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
Definition: Expansion.cpp:272
Expansion(SpatialDomains::GeometrySharedPtr pGeom)
Definition: Expansion.cpp:43
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLess > m_staticCondMatrixManager
Definition: TetExp.h:202
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager
Definition: TetExp.h:200
StdExpansion()
Default Constructor.
StdTetExp(const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc)
Definition: StdTetExp.cpp:42
int getNumberOfCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:187

◆ TetExp() [2/2]

Nektar::LocalRegions::TetExp::TetExp ( const TetExp T)

Copy Constructor.

Definition at line 80 of file TetExp.cpp.

81 : StdRegions::StdExpansion(T), StdRegions::StdExpansion3D(T),
82 StdRegions::StdTetExp(T), Expansion(T), Expansion3D(T),
83 m_matrixManager(T.m_matrixManager),
84 m_staticCondMatrixManager(T.m_staticCondMatrixManager)
85{
86}

◆ ~TetExp()

Nektar::LocalRegions::TetExp::~TetExp ( )
overridedefault

Member Function Documentation

◆ GeneralMatrixOp_MatOp()

void Nektar::LocalRegions::TetExp::GeneralMatrixOp_MatOp ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
const StdRegions::StdMatrixKey mkey 
)
private

Definition at line 1068 of file TetExp.cpp.

1071{
1073
1074 if (inarray.data() == outarray.data())
1075 {
1076 Array<OneD, NekDouble> tmp(m_ncoeffs);
1077 Vmath::Vcopy(m_ncoeffs, inarray.data(), 1, tmp.data(), 1);
1078
1079 Blas::Dgemv('N', m_ncoeffs, m_ncoeffs, mat->Scale(),
1080 (mat->GetOwnedMatrix())->GetPtr().data(), m_ncoeffs,
1081 tmp.data(), 1, 0.0, outarray.data(), 1);
1082 }
1083 else
1084 {
1085 Blas::Dgemv('N', m_ncoeffs, m_ncoeffs, mat->Scale(),
1086 (mat->GetOwnedMatrix())->GetPtr().data(), m_ncoeffs,
1087 inarray.data(), 1, 0.0, outarray.data(), 1);
1088 }
1089}
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:84
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 = alpha A x plus beta y where A[m x n].
Definition: Blas.hpp:211
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References Blas::Dgemv(), Nektar::LocalRegions::Expansion::GetLocMatrix(), Nektar::StdRegions::StdExpansion::m_ncoeffs, and Vmath::Vcopy().

◆ SetUpInverseTransformationMatrix()

void Nektar::LocalRegions::TetExp::SetUpInverseTransformationMatrix ( const DNekMatSharedPtr m_transformationmatrix,
DNekMatSharedPtr  m_inversetransformationmatrix,
DNekMatSharedPtr  m_inversetransposedtransformationmatrix 
)
protected

◆ v_AlignVectorToCollapsedDir()

void Nektar::LocalRegions::TetExp::v_AlignVectorToCollapsedDir ( const int  dir,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
overrideprotectedvirtual

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 386 of file TetExp.cpp.

389{
390 int i, j;
391
392 const int nquad0 = m_base[0]->GetNumPoints();
393 const int nquad1 = m_base[1]->GetNumPoints();
394 const int nquad2 = m_base[2]->GetNumPoints();
395 const int nqtot = nquad0 * nquad1 * nquad2;
396
397 const Array<OneD, const NekDouble> &z0 = m_base[0]->GetZ();
398 const Array<OneD, const NekDouble> &z1 = m_base[1]->GetZ();
399 const Array<OneD, const NekDouble> &z2 = m_base[2]->GetZ();
400
401 Array<OneD, NekDouble> tmp2(nqtot);
402 Array<OneD, NekDouble> tmp3(nqtot);
403
404 const Array<TwoD, const NekDouble> &df =
405 m_metricinfo->GetDerivFactors(GetPointsKeys());
406
407 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
408 {
409 Vmath::Vmul(nqtot, &df[3 * dir][0], 1, inarray.data(), 1, tmp2.data(),
410 1);
411 Vmath::Vmul(nqtot, &df[3 * dir + 1][0], 1, inarray.data(), 1,
412 tmp3.data(), 1);
413 Vmath::Vmul(nqtot, &df[3 * dir + 2][0], 1, inarray.data(), 1,
414 outarray[2].data(), 1);
415 }
416 else
417 {
418 Vmath::Smul(nqtot, df[3 * dir][0], inarray.data(), 1, tmp2.data(), 1);
419 Vmath::Smul(nqtot, df[3 * dir + 1][0], inarray.data(), 1, tmp3.data(),
420 1);
421 Vmath::Smul(nqtot, df[3 * dir + 2][0], inarray.data(), 1,
422 outarray[2].data(), 1);
423 }
424
425 NekDouble g0, g1, g1a, g2, g3;
426 int k, cnt;
427
428 for (cnt = 0, k = 0; k < nquad2; ++k)
429 {
430 for (j = 0; j < nquad1; ++j)
431 {
432 g2 = 2.0 / (1.0 - z2[k]);
433 g1 = g2 / (1.0 - z1[j]);
434 g0 = 2.0 * g1;
435 g3 = (1.0 + z1[j]) * g2 * 0.5;
436
437 for (i = 0; i < nquad0; ++i, ++cnt)
438 {
439 g1a = g1 * (1 + z0[i]);
440
441 outarray[0][cnt] =
442 g0 * tmp2[cnt] + g1a * (tmp3[cnt] + outarray[2][cnt]);
443
444 outarray[1][cnt] = g2 * tmp3[cnt] + g3 * outarray[2][cnt];
445 }
446 }
447 }
448}
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
Definition: Expansion.h:274
const LibUtilities::PointsKeyVector GetPointsKeys() const
Array< OneD, LibUtilities::BasisSharedPtr > m_base
@ eDeformed
Geometry is curved or has non-constant factors.
double NekDouble
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.hpp:72
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.hpp:100

References Nektar::SpatialDomains::eDeformed, Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::StdRegions::StdExpansion::m_base, Nektar::LocalRegions::Expansion::m_metricinfo, Vmath::Smul(), and Vmath::Vmul().

Referenced by v_IProductWRTDerivBase().

◆ v_ComputeLaplacianMetric()

void Nektar::LocalRegions::TetExp::v_ComputeLaplacianMetric ( )
overrideprotectedvirtual

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 1171 of file TetExp.cpp.

1172{
1173 if (m_metrics.count(eMetricQuadrature) == 0)
1174 {
1176 }
1177
1178 int i, j;
1179 const unsigned int nqtot = GetTotPoints();
1180 const unsigned int dim = 3;
1181 const MetricType m[3][3] = {
1185
1186 for (unsigned int i = 0; i < dim; ++i)
1187 {
1188 for (unsigned int j = i; j < dim; ++j)
1189 {
1190 m_metrics[m[i][j]] = Array<OneD, NekDouble>(nqtot);
1191 }
1192 }
1193
1194 // Define shorthand synonyms for m_metrics storage
1195 Array<OneD, NekDouble> g0(m_metrics[m[0][0]]);
1196 Array<OneD, NekDouble> g1(m_metrics[m[1][1]]);
1197 Array<OneD, NekDouble> g2(m_metrics[m[2][2]]);
1198 Array<OneD, NekDouble> g3(m_metrics[m[0][1]]);
1199 Array<OneD, NekDouble> g4(m_metrics[m[0][2]]);
1200 Array<OneD, NekDouble> g5(m_metrics[m[1][2]]);
1201
1202 // Allocate temporary storage
1203 Array<OneD, NekDouble> alloc(7 * nqtot, 0.0);
1204 Array<OneD, NekDouble> h0(alloc); // h0
1205 Array<OneD, NekDouble> h1(alloc + 1 * nqtot); // h1
1206 Array<OneD, NekDouble> h2(alloc + 2 * nqtot); // h2
1207 Array<OneD, NekDouble> h3(alloc + 3 * nqtot); // h3
1208 Array<OneD, NekDouble> wsp4(alloc + 4 * nqtot); // wsp4
1209 Array<OneD, NekDouble> wsp5(alloc + 5 * nqtot); // wsp5
1210 Array<OneD, NekDouble> wsp6(alloc + 6 * nqtot); // wsp6
1211 // Reuse some of the storage as workspace
1212 Array<OneD, NekDouble> wsp7(alloc); // wsp7
1213 Array<OneD, NekDouble> wsp8(alloc + 1 * nqtot); // wsp8
1214 Array<OneD, NekDouble> wsp9(alloc + 2 * nqtot); // wsp9
1215
1216 const Array<TwoD, const NekDouble> &df =
1217 m_metricinfo->GetDerivFactors(GetPointsKeys());
1218 const Array<OneD, const NekDouble> &z0 = m_base[0]->GetZ();
1219 const Array<OneD, const NekDouble> &z1 = m_base[1]->GetZ();
1220 const Array<OneD, const NekDouble> &z2 = m_base[2]->GetZ();
1221 const unsigned int nquad0 = m_base[0]->GetNumPoints();
1222 const unsigned int nquad1 = m_base[1]->GetNumPoints();
1223 const unsigned int nquad2 = m_base[2]->GetNumPoints();
1224
1225 for (j = 0; j < nquad2; ++j)
1226 {
1227 for (i = 0; i < nquad1; ++i)
1228 {
1229 Vmath::Fill(nquad0, 4.0 / (1.0 - z1[i]) / (1.0 - z2[j]),
1230 &h0[0] + i * nquad0 + j * nquad0 * nquad1, 1);
1231 Vmath::Fill(nquad0, 2.0 / (1.0 - z1[i]) / (1.0 - z2[j]),
1232 &h1[0] + i * nquad0 + j * nquad0 * nquad1, 1);
1233 Vmath::Fill(nquad0, 2.0 / (1.0 - z2[j]),
1234 &h2[0] + i * nquad0 + j * nquad0 * nquad1, 1);
1235 Vmath::Fill(nquad0, (1.0 + z1[i]) / (1.0 - z2[j]),
1236 &h3[0] + i * nquad0 + j * nquad0 * nquad1, 1);
1237 }
1238 }
1239 for (i = 0; i < nquad0; i++)
1240 {
1241 Blas::Dscal(nquad1 * nquad2, 1 + z0[i], &h1[0] + i, nquad0);
1242 }
1243
1244 // Step 3. Construct combined metric terms for physical space to
1245 // collapsed coordinate system.
1246 // Order of construction optimised to minimise temporary storage
1247 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
1248 {
1249 // wsp4
1250 Vmath::Vadd(nqtot, &df[1][0], 1, &df[2][0], 1, &wsp4[0], 1);
1251 Vmath::Vvtvvtp(nqtot, &df[0][0], 1, &h0[0], 1, &wsp4[0], 1, &h1[0], 1,
1252 &wsp4[0], 1);
1253 // wsp5
1254 Vmath::Vadd(nqtot, &df[4][0], 1, &df[5][0], 1, &wsp5[0], 1);
1255 Vmath::Vvtvvtp(nqtot, &df[3][0], 1, &h0[0], 1, &wsp5[0], 1, &h1[0], 1,
1256 &wsp5[0], 1);
1257 // wsp6
1258 Vmath::Vadd(nqtot, &df[7][0], 1, &df[8][0], 1, &wsp6[0], 1);
1259 Vmath::Vvtvvtp(nqtot, &df[6][0], 1, &h0[0], 1, &wsp6[0], 1, &h1[0], 1,
1260 &wsp6[0], 1);
1261
1262 // g0
1263 Vmath::Vvtvvtp(nqtot, &wsp4[0], 1, &wsp4[0], 1, &wsp5[0], 1, &wsp5[0],
1264 1, &g0[0], 1);
1265 Vmath::Vvtvp(nqtot, &wsp6[0], 1, &wsp6[0], 1, &g0[0], 1, &g0[0], 1);
1266
1267 // g4
1268 Vmath::Vvtvvtp(nqtot, &df[2][0], 1, &wsp4[0], 1, &df[5][0], 1, &wsp5[0],
1269 1, &g4[0], 1);
1270 Vmath::Vvtvp(nqtot, &df[8][0], 1, &wsp6[0], 1, &g4[0], 1, &g4[0], 1);
1271
1272 // overwrite h0, h1, h2
1273 // wsp7 (h2f1 + h3f2)
1274 Vmath::Vvtvvtp(nqtot, &df[1][0], 1, &h2[0], 1, &df[2][0], 1, &h3[0], 1,
1275 &wsp7[0], 1);
1276 // wsp8 (h2f4 + h3f5)
1277 Vmath::Vvtvvtp(nqtot, &df[4][0], 1, &h2[0], 1, &df[5][0], 1, &h3[0], 1,
1278 &wsp8[0], 1);
1279 // wsp9 (h2f7 + h3f8)
1280 Vmath::Vvtvvtp(nqtot, &df[7][0], 1, &h2[0], 1, &df[8][0], 1, &h3[0], 1,
1281 &wsp9[0], 1);
1282
1283 // g3
1284 Vmath::Vvtvvtp(nqtot, &wsp4[0], 1, &wsp7[0], 1, &wsp5[0], 1, &wsp8[0],
1285 1, &g3[0], 1);
1286 Vmath::Vvtvp(nqtot, &wsp6[0], 1, &wsp9[0], 1, &g3[0], 1, &g3[0], 1);
1287
1288 // overwrite wsp4, wsp5, wsp6
1289 // g1
1290 Vmath::Vvtvvtp(nqtot, &wsp7[0], 1, &wsp7[0], 1, &wsp8[0], 1, &wsp8[0],
1291 1, &g1[0], 1);
1292 Vmath::Vvtvp(nqtot, &wsp9[0], 1, &wsp9[0], 1, &g1[0], 1, &g1[0], 1);
1293
1294 // g5
1295 Vmath::Vvtvvtp(nqtot, &df[2][0], 1, &wsp7[0], 1, &df[5][0], 1, &wsp8[0],
1296 1, &g5[0], 1);
1297 Vmath::Vvtvp(nqtot, &df[8][0], 1, &wsp9[0], 1, &g5[0], 1, &g5[0], 1);
1298
1299 // g2
1300 Vmath::Vvtvvtp(nqtot, &df[2][0], 1, &df[2][0], 1, &df[5][0], 1,
1301 &df[5][0], 1, &g2[0], 1);
1302 Vmath::Vvtvp(nqtot, &df[8][0], 1, &df[8][0], 1, &g2[0], 1, &g2[0], 1);
1303 }
1304 else
1305 {
1306 // wsp4
1307 Vmath::Svtsvtp(nqtot, df[0][0], &h0[0], 1, df[1][0] + df[2][0], &h1[0],
1308 1, &wsp4[0], 1);
1309 // wsp5
1310 Vmath::Svtsvtp(nqtot, df[3][0], &h0[0], 1, df[4][0] + df[5][0], &h1[0],
1311 1, &wsp5[0], 1);
1312 // wsp6
1313 Vmath::Svtsvtp(nqtot, df[6][0], &h0[0], 1, df[7][0] + df[8][0], &h1[0],
1314 1, &wsp6[0], 1);
1315
1316 // g0
1317 Vmath::Vvtvvtp(nqtot, &wsp4[0], 1, &wsp4[0], 1, &wsp5[0], 1, &wsp5[0],
1318 1, &g0[0], 1);
1319 Vmath::Vvtvp(nqtot, &wsp6[0], 1, &wsp6[0], 1, &g0[0], 1, &g0[0], 1);
1320
1321 // g4
1322 Vmath::Svtsvtp(nqtot, df[2][0], &wsp4[0], 1, df[5][0], &wsp5[0], 1,
1323 &g4[0], 1);
1324 Vmath::Svtvp(nqtot, df[8][0], &wsp6[0], 1, &g4[0], 1, &g4[0], 1);
1325
1326 // overwrite h0, h1, h2
1327 // wsp7 (h2f1 + h3f2)
1328 Vmath::Svtsvtp(nqtot, df[1][0], &h2[0], 1, df[2][0], &h3[0], 1,
1329 &wsp7[0], 1);
1330 // wsp8 (h2f4 + h3f5)
1331 Vmath::Svtsvtp(nqtot, df[4][0], &h2[0], 1, df[5][0], &h3[0], 1,
1332 &wsp8[0], 1);
1333 // wsp9 (h2f7 + h3f8)
1334 Vmath::Svtsvtp(nqtot, df[7][0], &h2[0], 1, df[8][0], &h3[0], 1,
1335 &wsp9[0], 1);
1336
1337 // g3
1338 Vmath::Vvtvvtp(nqtot, &wsp4[0], 1, &wsp7[0], 1, &wsp5[0], 1, &wsp8[0],
1339 1, &g3[0], 1);
1340 Vmath::Vvtvp(nqtot, &wsp6[0], 1, &wsp9[0], 1, &g3[0], 1, &g3[0], 1);
1341
1342 // overwrite wsp4, wsp5, wsp6
1343 // g1
1344 Vmath::Vvtvvtp(nqtot, &wsp7[0], 1, &wsp7[0], 1, &wsp8[0], 1, &wsp8[0],
1345 1, &g1[0], 1);
1346 Vmath::Vvtvp(nqtot, &wsp9[0], 1, &wsp9[0], 1, &g1[0], 1, &g1[0], 1);
1347
1348 // g5
1349 Vmath::Svtsvtp(nqtot, df[2][0], &wsp7[0], 1, df[5][0], &wsp8[0], 1,
1350 &g5[0], 1);
1351 Vmath::Svtvp(nqtot, df[8][0], &wsp9[0], 1, &g5[0], 1, &g5[0], 1);
1352
1353 // g2
1354 Vmath::Fill(nqtot,
1355 df[2][0] * df[2][0] + df[5][0] * df[5][0] +
1356 df[8][0] * df[8][0],
1357 &g2[0], 1);
1358 }
1359
1360 for (unsigned int i = 0; i < dim; ++i)
1361 {
1362 for (unsigned int j = i; j < dim; ++j)
1363 {
1365 }
1366 }
1367}
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
Definition: StdExpansion.h:134
void MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdExpansion.h:732
static void Dscal(const int &n, const double &alpha, double *x, const int &incx)
BLAS level 1: x = alpha x.
Definition: Blas.hpp:149
void Svtsvtp(int n, const T alpha, const T *x, int incx, const T beta, const T *y, int incy, T *z, int incz)
Svtsvtp (scalar times vector plus scalar times vector):
Definition: Vmath.hpp:473
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Svtvp (scalar times vector plus vector): z = alpha*x + y.
Definition: Vmath.hpp:396
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.hpp:366
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.hpp:180
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.hpp:54
void Vvtvvtp(int n, const T *v, int incv, const T *w, int incw, const T *x, int incx, const T *y, int incy, T *z, int incz)
vvtvvtp (vector times vector plus vector times vector):
Definition: Vmath.hpp:439

◆ v_ComputeTraceNormal()

void Nektar::LocalRegions::TetExp::v_ComputeTraceNormal ( const int  face)
overrideprotectedvirtual

Compute the normal of a triangular face.

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 703 of file TetExp.cpp.

704{
705 int i;
706 const SpatialDomains::GeomFactorsSharedPtr &geomFactors =
707 GetGeom()->GetMetricInfo();
708
710 for (int i = 0; i < ptsKeys.size(); ++i)
711 {
712 // Need at least 2 points for computing normals
713 if (ptsKeys[i].GetNumPoints() == 1)
714 {
715 LibUtilities::PointsKey pKey(2, ptsKeys[i].GetPointsType());
716 ptsKeys[i] = pKey;
717 }
718 }
719
720 SpatialDomains::GeomType type = geomFactors->GetGtype();
721 const Array<TwoD, const NekDouble> &df =
722 geomFactors->GetDerivFactors(ptsKeys);
723 const Array<OneD, const NekDouble> &jac = geomFactors->GetJac(ptsKeys);
724
725 LibUtilities::BasisKey tobasis0 = GetTraceBasisKey(face, 0);
726 LibUtilities::BasisKey tobasis1 = GetTraceBasisKey(face, 1);
727
728 // number of face quadrature points
729 int nq_face = tobasis0.GetNumPoints() * tobasis1.GetNumPoints();
730
731 int vCoordDim = GetCoordim();
732
733 m_traceNormals[face] = Array<OneD, Array<OneD, NekDouble>>(vCoordDim);
734 Array<OneD, Array<OneD, NekDouble>> &normal = m_traceNormals[face];
735 for (i = 0; i < vCoordDim; ++i)
736 {
737 normal[i] = Array<OneD, NekDouble>(nq_face);
738 }
739
740 size_t nqb = nq_face;
741 size_t nbnd = face;
742 m_elmtBndNormDirElmtLen[nbnd] = Array<OneD, NekDouble>{nqb, 0.0};
743 Array<OneD, NekDouble> &length = m_elmtBndNormDirElmtLen[nbnd];
744
745 // Regular geometry case
746 if (type == SpatialDomains::eRegular ||
748 {
749 NekDouble fac;
750
751 // Set up normals
752 switch (face)
753 {
754 case 0:
755 {
756 for (i = 0; i < vCoordDim; ++i)
757 {
758 normal[i][0] = -df[3 * i + 2][0];
759 }
760
761 break;
762 }
763 case 1:
764 {
765 for (i = 0; i < vCoordDim; ++i)
766 {
767 normal[i][0] = -df[3 * i + 1][0];
768 }
769
770 break;
771 }
772 case 2:
773 {
774 for (i = 0; i < vCoordDim; ++i)
775 {
776 normal[i][0] =
777 df[3 * i][0] + df[3 * i + 1][0] + df[3 * i + 2][0];
778 }
779
780 break;
781 }
782 case 3:
783 {
784 for (i = 0; i < vCoordDim; ++i)
785 {
786 normal[i][0] = -df[3 * i][0];
787 }
788 break;
789 }
790 default:
791 ASSERTL0(false, "face is out of range (edge < 3)");
792 }
793
794 // normalise
795 fac = 0.0;
796 for (i = 0; i < vCoordDim; ++i)
797 {
798 fac += normal[i][0] * normal[i][0];
799 }
800 fac = 1.0 / sqrt(fac);
801 Vmath::Fill(nqb, fac, length, 1);
802
803 for (i = 0; i < vCoordDim; ++i)
804 {
805 Vmath::Fill(nq_face, fac * normal[i][0], normal[i], 1);
806 }
807 }
808 else
809 {
810 // Set up deformed normals
811 int j, k;
812
813 int nq0 = ptsKeys[0].GetNumPoints();
814 int nq1 = ptsKeys[1].GetNumPoints();
815 int nq2 = ptsKeys[2].GetNumPoints();
816 int nqtot;
817 int nq01 = nq0 * nq1;
818
819 // number of elemental quad points
820 if (face == 0)
821 {
822 nqtot = nq01;
823 }
824 else if (face == 1)
825 {
826 nqtot = nq0 * nq2;
827 }
828 else
829 {
830 nqtot = nq1 * nq2;
831 }
832
833 LibUtilities::PointsKey points0;
834 LibUtilities::PointsKey points1;
835
836 Array<OneD, NekDouble> faceJac(nqtot);
837 Array<OneD, NekDouble> normals(vCoordDim * nqtot, 0.0);
838
839 // Extract Jacobian along face and recover local derivates
840 // (dx/dr) for polynomial interpolation by multiplying m_gmat by
841 // jacobian
842 switch (face)
843 {
844 case 0:
845 {
846 for (j = 0; j < nq01; ++j)
847 {
848 normals[j] = -df[2][j] * jac[j];
849 normals[nqtot + j] = -df[5][j] * jac[j];
850 normals[2 * nqtot + j] = -df[8][j] * jac[j];
851 faceJac[j] = jac[j];
852 }
853
854 points0 = ptsKeys[0];
855 points1 = ptsKeys[1];
856 break;
857 }
858
859 case 1:
860 {
861 for (j = 0; j < nq0; ++j)
862 {
863 for (k = 0; k < nq2; ++k)
864 {
865 int tmp = j + nq01 * k;
866 normals[j + k * nq0] = -df[1][tmp] * jac[tmp];
867 normals[nqtot + j + k * nq0] = -df[4][tmp] * jac[tmp];
868 normals[2 * nqtot + j + k * nq0] =
869 -df[7][tmp] * jac[tmp];
870 faceJac[j + k * nq0] = jac[tmp];
871 }
872 }
873
874 points0 = ptsKeys[0];
875 points1 = ptsKeys[2];
876 break;
877 }
878
879 case 2:
880 {
881 for (j = 0; j < nq1; ++j)
882 {
883 for (k = 0; k < nq2; ++k)
884 {
885 int tmp = nq0 - 1 + nq0 * j + nq01 * k;
886 normals[j + k * nq1] =
887 (df[0][tmp] + df[1][tmp] + df[2][tmp]) * jac[tmp];
888 normals[nqtot + j + k * nq1] =
889 (df[3][tmp] + df[4][tmp] + df[5][tmp]) * jac[tmp];
890 normals[2 * nqtot + j + k * nq1] =
891 (df[6][tmp] + df[7][tmp] + df[8][tmp]) * jac[tmp];
892 faceJac[j + k * nq1] = jac[tmp];
893 }
894 }
895
896 points0 = ptsKeys[1];
897 points1 = ptsKeys[2];
898 break;
899 }
900
901 case 3:
902 {
903 for (j = 0; j < nq1; ++j)
904 {
905 for (k = 0; k < nq2; ++k)
906 {
907 int tmp = j * nq0 + nq01 * k;
908 normals[j + k * nq1] = -df[0][tmp] * jac[tmp];
909 normals[nqtot + j + k * nq1] = -df[3][tmp] * jac[tmp];
910 normals[2 * nqtot + j + k * nq1] =
911 -df[6][tmp] * jac[tmp];
912 faceJac[j + k * nq1] = jac[tmp];
913 }
914 }
915
916 points0 = ptsKeys[1];
917 points1 = ptsKeys[2];
918 break;
919 }
920
921 default:
922 ASSERTL0(false, "face is out of range (face < 3)");
923 }
924
925 Array<OneD, NekDouble> work(nq_face, 0.0);
926 // Interpolate Jacobian and invert
927 LibUtilities::Interp2D(points0, points1, faceJac,
928 tobasis0.GetPointsKey(), tobasis1.GetPointsKey(),
929 work);
930 Vmath::Sdiv(nq_face, 1.0, &work[0], 1, &work[0], 1);
931
932 // Interpolate normal and multiply by inverse Jacobian.
933 for (i = 0; i < vCoordDim; ++i)
934 {
935 LibUtilities::Interp2D(points0, points1, &normals[i * nqtot],
936 tobasis0.GetPointsKey(),
937 tobasis1.GetPointsKey(), &normal[i][0]);
938 Vmath::Vmul(nq_face, work, 1, normal[i], 1, normal[i], 1);
939 }
940
941 // Normalise to obtain unit normals.
942 Vmath::Zero(nq_face, work, 1);
943 for (i = 0; i < GetCoordim(); ++i)
944 {
945 Vmath::Vvtvp(nq_face, normal[i], 1, normal[i], 1, work, 1, work, 1);
946 }
947
948 Vmath::Vsqrt(nq_face, work, 1, work, 1);
949 Vmath::Sdiv(nq_face, 1.0, work, 1, work, 1);
950
951 Vmath::Vcopy(nqb, work, 1, length, 1);
952
953 for (i = 0; i < GetCoordim(); ++i)
954 {
955 Vmath::Vmul(nq_face, normal[i], 1, work, 1, normal[i], 1);
956 }
957 }
958}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
std::map< int, NormalVector > m_traceNormals
Definition: Expansion.h:276
std::map< int, Array< OneD, NekDouble > > m_elmtBndNormDirElmtLen
the element length in each element boundary(Vertex, edge or face) normal direction calculated based o...
Definition: Expansion.h:286
SpatialDomains::GeometrySharedPtr GetGeom() const
Definition: Expansion.cpp:167
LibUtilities::PointsType GetPointsType(const int dir) const
This function returns the type of quadrature points used in the dir direction.
Definition: StdExpansion.h:205
int GetNumPoints(const int dir) const
This function returns the number of quadrature points in the dir direction.
Definition: StdExpansion.h:218
const LibUtilities::BasisKey GetTraceBasisKey(const int i, int k=-1, bool UseGLL=false) const
This function returns the basis key belonging to the i-th trace.
Definition: StdExpansion.h:301
void Interp2D(const BasisKey &fbasis0, const BasisKey &fbasis1, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, Array< OneD, NekDouble > &to)
this function interpolates a 2D function evaluated at the quadrature points of the 2D basis,...
Definition: Interp.cpp:101
std::vector< PointsKey > PointsKeyVector
Definition: Points.h:231
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
Definition: GeomFactors.h:60
GeomType
Indicates the type of element geometry.
@ eRegular
Geometry is straight-sided with constant geometric factors.
@ eMovingRegular
Currently unused.
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.hpp:340
void Sdiv(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha/x.
Definition: Vmath.hpp:154
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References ASSERTL0, Nektar::SpatialDomains::eMovingRegular, Nektar::SpatialDomains::eRegular, Vmath::Fill(), Nektar::StdRegions::StdExpansion::GetCoordim(), Nektar::LocalRegions::Expansion::GetGeom(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::StdRegions::StdExpansion::GetNumPoints(), Nektar::LibUtilities::BasisKey::GetPointsKey(), Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::StdRegions::StdExpansion::GetPointsType(), Nektar::StdRegions::StdExpansion::GetTraceBasisKey(), Nektar::LibUtilities::Interp2D(), Nektar::LocalRegions::Expansion::m_elmtBndNormDirElmtLen, Nektar::LocalRegions::Expansion::m_traceNormals, Vmath::Sdiv(), tinysimd::sqrt(), Vmath::Vcopy(), Vmath::Vmul(), Vmath::Vsqrt(), Vmath::Vvtvp(), and Vmath::Zero().

◆ v_CreateStdMatrix()

DNekMatSharedPtr Nektar::LocalRegions::TetExp::v_CreateStdMatrix ( const StdRegions::StdMatrixKey mkey)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 1037 of file TetExp.cpp.

1038{
1039 LibUtilities::BasisKey bkey0 = m_base[0]->GetBasisKey();
1040 LibUtilities::BasisKey bkey1 = m_base[1]->GetBasisKey();
1041 LibUtilities::BasisKey bkey2 = m_base[2]->GetBasisKey();
1044
1045 return tmp->GetStdMatrix(mkey);
1046}
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< StdTetExp > StdTetExpSharedPtr
Definition: StdTetExp.h:227

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and Nektar::StdRegions::StdExpansion::m_base.

◆ v_DetShapeType()

LibUtilities::ShapeType Nektar::LocalRegions::TetExp::v_DetShapeType ( ) const
overrideprotectedvirtual

Return Shape of region, using ShapeType enum list.

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 530 of file TetExp.cpp.

531{
533}

References Nektar::LibUtilities::eTetrahedron.

◆ v_DropLocMatrix()

void Nektar::LocalRegions::TetExp::v_DropLocMatrix ( const MatrixKey mkey)
overrideprotectedvirtual

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 1053 of file TetExp.cpp.

1054{
1055 m_matrixManager.DeleteObject(mkey);
1056}

References m_matrixManager.

◆ v_DropLocStaticCondMatrix()

void Nektar::LocalRegions::TetExp::v_DropLocStaticCondMatrix ( const MatrixKey mkey)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1063 of file TetExp.cpp.

1064{
1065 m_staticCondMatrixManager.DeleteObject(mkey);
1066}

References m_staticCondMatrixManager.

◆ v_ExtractDataToCoeffs()

void Nektar::LocalRegions::TetExp::v_ExtractDataToCoeffs ( const NekDouble data,
const std::vector< unsigned int > &  nummodes,
const int  mode_offset,
NekDouble coeffs,
std::vector< LibUtilities::BasisType > &  fromType 
)
overrideprotectedvirtual

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 555 of file TetExp.cpp.

559{
560 int data_order0 = nummodes[mode_offset];
561 int fillorder0 = min(m_base[0]->GetNumModes(), data_order0);
562 int data_order1 = nummodes[mode_offset + 1];
563 int order1 = m_base[1]->GetNumModes();
564 int fillorder1 = min(order1, data_order1);
565 int data_order2 = nummodes[mode_offset + 2];
566 int order2 = m_base[2]->GetNumModes();
567 int fillorder2 = min(order2, data_order2);
568
569 switch (m_base[0]->GetBasisType())
570 {
572 {
573 int i, j;
574 int cnt = 0;
575 int cnt1 = 0;
576
578 "Extraction routine not set up for this basis");
580 "Extraction routine not set up for this basis");
581
582 Vmath::Zero(m_ncoeffs, coeffs, 1);
583 for (j = 0; j < fillorder0; ++j)
584 {
585 for (i = 0; i < fillorder1 - j; ++i)
586 {
587 Vmath::Vcopy(fillorder2 - j - i, &data[cnt], 1,
588 &coeffs[cnt1], 1);
589 cnt += data_order2 - j - i;
590 cnt1 += order2 - j - i;
591 }
592
593 // count out data for j iteration
594 for (i = fillorder1 - j; i < data_order1 - j; ++i)
595 {
596 cnt += data_order2 - j - i;
597 }
598
599 for (i = fillorder1 - j; i < order1 - j; ++i)
600 {
601 cnt1 += order2 - j - i;
602 }
603 }
604 }
605 break;
606 default:
607 ASSERTL0(false, "basis is either not set up or not "
608 "hierarchicial");
609 }
610}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
Definition: StdExpansion.h:156
@ eModified_B
Principle Modified Functions .
Definition: BasisType.h:49
@ eModified_C
Principle Modified Functions .
Definition: BasisType.h:50
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:48

References ASSERTL0, ASSERTL1, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModified_B, Nektar::LibUtilities::eModified_C, Nektar::StdRegions::StdExpansion::GetBasisType(), Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::m_ncoeffs, Vmath::Vcopy(), and Vmath::Zero().

◆ v_FwdTrans()

void Nektar::LocalRegions::TetExp::v_FwdTrans ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
overrideprotectedvirtual

Forward transform from physical quadrature space stored in inarray and evaluate the expansion coefficients and store in (this)->_coeffs.

Parameters
inarrayArray of physical quadrature points to be transformed.
outarrayArray of coefficients to update.

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 223 of file TetExp.cpp.

225{
226 if ((m_base[0]->Collocation()) && (m_base[1]->Collocation()) &&
227 (m_base[2]->Collocation()))
228 {
229 Vmath::Vcopy(GetNcoeffs(), &inarray[0], 1, &outarray[0], 1);
230 }
231 else
232 {
233 IProductWRTBase(inarray, outarray);
234
235 // get Mass matrix inverse
236 MatrixKey masskey(StdRegions::eInvMass, DetShapeType(), *this);
237 DNekScalMatSharedPtr matsys = m_matrixManager[masskey];
238
239 // copy inarray in case inarray == outarray
240 DNekVec in(m_ncoeffs, outarray);
241 DNekVec out(m_ncoeffs, outarray, eWrapper);
242
243 out = (*matsys) * in;
244 }
245}
int GetNcoeffs(void) const
This function returns the total number of coefficients used in the expansion.
Definition: StdExpansion.h:124
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 expa...
Definition: StdExpansion.h:537
LibUtilities::ShapeType DetShapeType() const
Definition: StdTetExp.h:56
NekVector< NekDouble > DNekVec
Definition: NekTypeDefs.hpp:48

References Nektar::StdRegions::StdTetExp::DetShapeType(), Nektar::StdRegions::eInvMass, Nektar::eWrapper, Nektar::StdRegions::StdExpansion::GetNcoeffs(), Nektar::StdRegions::StdExpansion::IProductWRTBase(), Nektar::StdRegions::StdExpansion::m_base, m_matrixManager, Nektar::StdRegions::StdExpansion::m_ncoeffs, and Vmath::Vcopy().

◆ v_GenMatrix()

DNekMatSharedPtr Nektar::LocalRegions::TetExp::v_GenMatrix ( const StdRegions::StdMatrixKey mkey)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 1015 of file TetExp.cpp.

1016{
1017 DNekMatSharedPtr returnval;
1018
1019 switch (mkey.GetMatrixType())
1020 {
1028 returnval = Expansion3D::v_GenMatrix(mkey);
1029 break;
1030 default:
1031 returnval = StdTetExp::v_GenMatrix(mkey);
1032 }
1033
1034 return returnval;
1035}
DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey) override
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75

References Nektar::StdRegions::eHybridDGHelmBndLam, Nektar::StdRegions::eHybridDGHelmholtz, Nektar::StdRegions::eHybridDGLamToQ0, Nektar::StdRegions::eHybridDGLamToQ1, Nektar::StdRegions::eHybridDGLamToQ2, Nektar::StdRegions::eHybridDGLamToU, Nektar::StdRegions::eInvLaplacianWithUnityMean, Nektar::StdRegions::StdMatrixKey::GetMatrixType(), and Nektar::LocalRegions::Expansion3D::v_GenMatrix().

◆ v_GetCoord()

void Nektar::LocalRegions::TetExp::v_GetCoord ( const Array< OneD, const NekDouble > &  Lcoords,
Array< OneD, NekDouble > &  coords 
)
overrideprotectedvirtual

Get the coordinates "coords" at the local coordinates "Lcoords".

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 499 of file TetExp.cpp.

501{
502 int i;
503
504 ASSERTL1(Lcoords[0] <= -1.0 && Lcoords[0] >= 1.0 && Lcoords[1] <= -1.0 &&
505 Lcoords[1] >= 1.0 && Lcoords[2] <= -1.0 && Lcoords[2] >= 1.0,
506 "Local coordinates are not in region [-1,1]");
507
508 // m_geom->FillGeom(); // TODO: implement FillGeom()
509
510 for (i = 0; i < m_geom->GetCoordim(); ++i)
511 {
512 coords[i] = m_geom->GetCoord(i, Lcoords);
513 }
514}
SpatialDomains::GeometrySharedPtr m_geom
Definition: Expansion.h:273

References ASSERTL1, and Nektar::LocalRegions::Expansion::m_geom.

◆ v_GetCoords()

void Nektar::LocalRegions::TetExp::v_GetCoords ( Array< OneD, NekDouble > &  coords_1,
Array< OneD, NekDouble > &  coords_2,
Array< OneD, NekDouble > &  coords_3 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 516 of file TetExp.cpp.

519{
520 Expansion::v_GetCoords(coords_0, coords_1, coords_2);
521}
void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
Definition: Expansion.cpp:534

References Nektar::LocalRegions::Expansion::v_GetCoords().

◆ v_GetLinStdExp()

StdRegions::StdExpansionSharedPtr Nektar::LocalRegions::TetExp::v_GetLinStdExp ( void  ) const
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 542 of file TetExp.cpp.

543{
544 LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(), 2,
545 m_base[0]->GetPointsKey());
546 LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(), 2,
547 m_base[1]->GetPointsKey());
548 LibUtilities::BasisKey bkey2(m_base[2]->GetBasisType(), 2,
549 m_base[2]->GetPointsKey());
550
552 bkey2);
553}

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::StdRegions::StdExpansion::GetBasisType(), and Nektar::StdRegions::StdExpansion::m_base.

◆ v_GetLocMatrix()

DNekScalMatSharedPtr Nektar::LocalRegions::TetExp::v_GetLocMatrix ( const MatrixKey mkey)
overrideprotectedvirtual

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 1048 of file TetExp.cpp.

1049{
1050 return m_matrixManager[mkey];
1051}

References m_matrixManager.

◆ v_GetLocStaticCondMatrix()

DNekScalBlkMatSharedPtr Nektar::LocalRegions::TetExp::v_GetLocStaticCondMatrix ( const MatrixKey mkey)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1058 of file TetExp.cpp.

1059{
1060 return m_staticCondMatrixManager[mkey];
1061}

References m_staticCondMatrixManager.

◆ v_GetStdExp()

StdRegions::StdExpansionSharedPtr Nektar::LocalRegions::TetExp::v_GetStdExp ( void  ) const
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 535 of file TetExp.cpp.

536{
538 m_base[0]->GetBasisKey(), m_base[1]->GetBasisKey(),
539 m_base[2]->GetBasisKey());
540}

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and Nektar::StdRegions::StdExpansion::m_base.

◆ v_GetTracePhysMap()

void Nektar::LocalRegions::TetExp::v_GetTracePhysMap ( const int  face,
Array< OneD, int > &  outarray 
)
overrideprotectedvirtual

Returns the physical values at the quadrature points of a face.

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 615 of file TetExp.cpp.

616{
617 int nquad0 = m_base[0]->GetNumPoints();
618 int nquad1 = m_base[1]->GetNumPoints();
619 int nquad2 = m_base[2]->GetNumPoints();
620
621 int nq0 = 0;
622 int nq1 = 0;
623
624 // get forward aligned faces.
625 switch (face)
626 {
627 case 0:
628 {
629 nq0 = nquad0;
630 nq1 = nquad1;
631 if (outarray.size() != nq0 * nq1)
632 {
633 outarray = Array<OneD, int>(nq0 * nq1);
634 }
635
636 for (int i = 0; i < nquad0 * nquad1; ++i)
637 {
638 outarray[i] = i;
639 }
640
641 break;
642 }
643 case 1:
644 {
645 nq0 = nquad0;
646 nq1 = nquad2;
647 if (outarray.size() != nq0 * nq1)
648 {
649 outarray = Array<OneD, int>(nq0 * nq1);
650 }
651
652 // Direction A and B positive
653 for (int k = 0; k < nquad2; k++)
654 {
655 for (int i = 0; i < nquad0; ++i)
656 {
657 outarray[k * nquad0 + i] = (nquad0 * nquad1 * k) + i;
658 }
659 }
660 break;
661 }
662 case 2:
663 {
664 nq0 = nquad1;
665 nq1 = nquad2;
666 if (outarray.size() != nq0 * nq1)
667 {
668 outarray = Array<OneD, int>(nq0 * nq1);
669 }
670
671 // Directions A and B positive
672 for (int j = 0; j < nquad1 * nquad2; ++j)
673 {
674 outarray[j] = nquad0 - 1 + j * nquad0;
675 }
676 break;
677 }
678 case 3:
679 {
680 nq0 = nquad1;
681 nq1 = nquad2;
682 if (outarray.size() != nq0 * nq1)
683 {
684 outarray = Array<OneD, int>(nq0 * nq1);
685 }
686
687 // Directions A and B positive
688 for (int j = 0; j < nquad1 * nquad2; ++j)
689 {
690 outarray[j] = j * nquad0;
691 }
692 }
693 break;
694 default:
695 ASSERTL0(false, "face value (> 3) is out of range");
696 break;
697 }
698}

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

◆ v_HelmholtzMatrixOp()

void Nektar::LocalRegions::TetExp::v_HelmholtzMatrixOp ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
const StdRegions::StdMatrixKey mkey 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 963 of file TetExp.cpp.

966{
967 TetExp::v_HelmholtzMatrixOp_MatFree(inarray, outarray, mkey);
968}
void v_HelmholtzMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override

References Nektar::StdRegions::StdExpansion3D::v_HelmholtzMatrixOp_MatFree().

◆ v_Integral()

NekDouble Nektar::LocalRegions::TetExp::v_Integral ( const Array< OneD, const NekDouble > &  inarray)
overrideprotectedvirtual

Integrate the physical point list inarray over region.

Parameters
inarrayDefinition of function to be returned at quadrature point of expansion.
Returns
\(\int^1_{-1}\int^1_{-1} \int^1_{-1} u(\eta_1, \eta_2, \eta_3) J[i,j,k] d \eta_1 d \eta_2 d \eta_3 \) where \(inarray[i,j,k] = u(\eta_{1i},\eta_{2j},\eta_{3k}) \) and \( J[i,j,k] \) is the Jacobian evaluated at the quadrature point.

Reimplemented from Nektar::StdRegions::StdExpansion3D.

Definition at line 102 of file TetExp.cpp.

103{
104 int nquad0 = m_base[0]->GetNumPoints();
105 int nquad1 = m_base[1]->GetNumPoints();
106 int nquad2 = m_base[2]->GetNumPoints();
107 Array<OneD, const NekDouble> jac = m_metricinfo->GetJac(GetPointsKeys());
108 NekDouble retrunVal;
109 Array<OneD, NekDouble> tmp(nquad0 * nquad1 * nquad2);
110
111 // multiply inarray with Jacobian
112 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
113 {
114 Vmath::Vmul(nquad0 * nquad1 * nquad2, &jac[0], 1,
115 (NekDouble *)&inarray[0], 1, &tmp[0], 1);
116 }
117 else
118 {
119 Vmath::Smul(nquad0 * nquad1 * nquad2, (NekDouble)jac[0],
120 (NekDouble *)&inarray[0], 1, &tmp[0], 1);
121 }
122
123 // call StdTetExp version;
124 retrunVal = StdTetExp::v_Integral(tmp);
125
126 return retrunVal;
127}

References Nektar::SpatialDomains::eDeformed, Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::StdRegions::StdExpansion::m_base, Nektar::LocalRegions::Expansion::m_metricinfo, Vmath::Smul(), and Vmath::Vmul().

◆ v_IProductWRTBase()

void Nektar::LocalRegions::TetExp::v_IProductWRTBase ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
overrideprotectedvirtual

Calculate the inner product of inarray with respect to the basis B=m_base0*m_base1*m_base2 and put into outarray:

\( \begin{array}{rcl} I_{pqr} = (\phi_{pqr}, u)_{\delta} & = & \sum_{i=0}^{nq_0} \sum_{j=0}^{nq_1} \sum_{k=0}^{nq_2} \psi_{p}^{a} (\eta_{1i}) \psi_{pq}^{b} (\eta_{2j}) \psi_{pqr}^{c} (\eta_{3k}) w_i w_j w_k u(\eta_{1,i} \eta_{2,j} \eta_{3,k}) J_{i,j,k}\\ & = & \sum_{i=0}^{nq_0} \psi_p^a(\eta_{1,i}) \sum_{j=0}^{nq_1} \psi_{pq}^b(\eta_{2,j}) \sum_{k=0}^{nq_2} \psi_{pqr}^c u(\eta_{1i},\eta_{2j},\eta_{3k}) J_{i,j,k} \end{array} \)
where \( \phi_{pqr} (\xi_1 , \xi_2 , \xi_3) = \psi_p^a (\eta_1) \psi_{pq}^b (\eta_2) \psi_{pqr}^c (\eta_3) \) which can be implemented as
\(f_{pqr} (\xi_{3k}) = \sum_{k=0}^{nq_3} \psi_{pqr}^c u(\eta_{1i},\eta_{2j},\eta_{3k}) J_{i,j,k} = {\bf B_3 U} \)
\( g_{pq} (\xi_{3k}) = \sum_{j=0}^{nq_1} \psi_{pq}^b (\xi_{2j}) f_{pqr} (\xi_{3k}) = {\bf B_2 F} \)
\( (\phi_{pqr}, u)_{\delta} = \sum_{k=0}^{nq_0} \psi_{p}^a (\xi_{3k}) g_{pq} (\xi_{3k}) = {\bf B_1 G} \)

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 276 of file TetExp.cpp.

278{
279 v_IProductWRTBase_SumFac(inarray, outarray);
280}
void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
Definition: TetExp.cpp:282

References v_IProductWRTBase_SumFac().

◆ v_IProductWRTBase_SumFac()

void Nektar::LocalRegions::TetExp::v_IProductWRTBase_SumFac ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
bool  multiplybyweights = true 
)
overrideprotectedvirtual
Parameters
inarrayFunction evaluated at physical collocation points.
outarrayInner product with respect to each basis function over the element.

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 282 of file TetExp.cpp.

285{
286 const int nquad0 = m_base[0]->GetNumPoints();
287 const int nquad1 = m_base[1]->GetNumPoints();
288 const int nquad2 = m_base[2]->GetNumPoints();
289 const int order0 = m_base[0]->GetNumModes();
290 const int order1 = m_base[1]->GetNumModes();
291 Array<OneD, NekDouble> wsp(nquad1 * nquad2 * order0 +
292 nquad2 * order0 * (order1 + 1) / 2);
293
294 if (multiplybyweights)
295 {
296 Array<OneD, NekDouble> tmp(nquad0 * nquad1 * nquad2);
297
298 MultiplyByQuadratureMetric(inarray, tmp);
300 m_base[0]->GetBdata(), m_base[1]->GetBdata(), m_base[2]->GetBdata(),
301 tmp, outarray, wsp, true, true, true);
302 }
303 else
304 {
306 m_base[0]->GetBdata(), m_base[1]->GetBdata(), m_base[2]->GetBdata(),
307 inarray, outarray, wsp, true, true, true);
308 }
309}
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)

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

Referenced by v_IProductWRTBase().

◆ v_IProductWRTDerivBase()

void Nektar::LocalRegions::TetExp::v_IProductWRTDerivBase ( const int  dir,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
overrideprotectedvirtual

Calculates the inner product \( I_{pqr} = (u, \partial_{x_i} \phi_{pqr}) \).

The derivative of the basis functions is performed using the chain rule in order to incorporate the geometric factors. Assuming that the basis functions are a tensor product \(\phi_{pqr}(\eta_1,\eta_2,\eta_3) = \phi_1(\eta_1)\phi_2(\eta_2)\phi_3(\eta_3)\), this yields the result

\[ I_{pqr} = \sum_{j=1}^3 \left(u, \frac{\partial u}{\partial \eta_j} \frac{\partial \eta_j}{\partial x_i}\right) \]

In the prismatic element, we must also incorporate a second set of geometric factors which incorporate the collapsed co-ordinate system, so that

\[ \frac{\partial\eta_j}{\partial x_i} = \sum_{k=1}^3 \frac{\partial\eta_j}{\partial\xi_k}\frac{\partial\xi_k}{\partial x_i} \]

These derivatives can be found on p152 of Sherwin & Karniadakis.

Parameters
dirDirection in which to take the derivative.
inarrayThe function \( u \).
outarrayValue of the inner product.

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 341 of file TetExp.cpp.

344{
345 const int nquad0 = m_base[0]->GetNumPoints();
346 const int nquad1 = m_base[1]->GetNumPoints();
347 const int nquad2 = m_base[2]->GetNumPoints();
348 const int order0 = m_base[0]->GetNumModes();
349 const int order1 = m_base[1]->GetNumModes();
350 const int nqtot = nquad0 * nquad1 * nquad2;
351
352 Array<OneD, NekDouble> tmp1(nqtot);
353 Array<OneD, NekDouble> tmp2(nqtot);
354 Array<OneD, NekDouble> tmp3(nqtot);
355 Array<OneD, NekDouble> tmp4(nqtot);
356 Array<OneD, NekDouble> tmp6(m_ncoeffs);
357 Array<OneD, NekDouble> wsp(nquad1 * nquad2 * order0 +
358 nquad2 * order0 * (order1 + 1) / 2);
359
360 MultiplyByQuadratureMetric(inarray, tmp1);
361
362 Array<OneD, Array<OneD, NekDouble>> tmp2D{3};
363 tmp2D[0] = tmp2;
364 tmp2D[1] = tmp3;
365 tmp2D[2] = tmp4;
366
367 TetExp::v_AlignVectorToCollapsedDir(dir, tmp1, tmp2D);
368
369 IProductWRTBase_SumFacKernel(m_base[0]->GetDbdata(), m_base[1]->GetBdata(),
370 m_base[2]->GetBdata(), tmp2, outarray, wsp,
371 false, true, true);
372
373 IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(), m_base[1]->GetDbdata(),
374 m_base[2]->GetBdata(), tmp3, tmp6, wsp, true,
375 false, true);
376
377 Vmath::Vadd(m_ncoeffs, tmp6, 1, outarray, 1, outarray, 1);
378
379 IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(), m_base[1]->GetBdata(),
380 m_base[2]->GetDbdata(), tmp4, tmp6, wsp, true,
381 true, false);
382
383 Vmath::Vadd(m_ncoeffs, tmp6, 1, outarray, 1, outarray, 1);
384}
void v_AlignVectorToCollapsedDir(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
Definition: TetExp.cpp:386

References Nektar::StdRegions::StdExpansion3D::IProductWRTBase_SumFacKernel(), Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::m_ncoeffs, Nektar::StdRegions::StdExpansion::MultiplyByQuadratureMetric(), v_AlignVectorToCollapsedDir(), and Vmath::Vadd().

◆ v_LaplacianMatrixOp() [1/2]

void Nektar::LocalRegions::TetExp::v_LaplacianMatrixOp ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
const StdRegions::StdMatrixKey mkey 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 970 of file TetExp.cpp.

973{
974 TetExp::v_LaplacianMatrixOp_MatFree(inarray, outarray, mkey);
975}
void v_LaplacianMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override

References Nektar::StdRegions::StdExpansion3D::v_LaplacianMatrixOp_MatFree().

◆ v_LaplacianMatrixOp() [2/2]

void Nektar::LocalRegions::TetExp::v_LaplacianMatrixOp ( const int  k1,
const int  k2,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
const StdRegions::StdMatrixKey mkey 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 977 of file TetExp.cpp.

981{
982 StdExpansion::LaplacianMatrixOp_MatFree(k1, k2, inarray, outarray, mkey);
983}

◆ v_LaplacianMatrixOp_MatFree_Kernel()

void Nektar::LocalRegions::TetExp::v_LaplacianMatrixOp_MatFree_Kernel ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
Array< OneD, NekDouble > &  wsp 
)
overrideprivatevirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1091 of file TetExp.cpp.

1094{
1095 // This implementation is only valid when there are no
1096 // coefficients associated to the Laplacian operator
1097 if (m_metrics.count(eMetricLaplacian00) == 0)
1098 {
1100 }
1101
1102 int nquad0 = m_base[0]->GetNumPoints();
1103 int nquad1 = m_base[1]->GetNumPoints();
1104 int nquad2 = m_base[2]->GetNumPoints();
1105 int nqtot = nquad0 * nquad1 * nquad2;
1106
1107 ASSERTL1(wsp.size() >= 6 * nqtot, "Insufficient workspace size.");
1108 ASSERTL1(m_ncoeffs <= nqtot, "Workspace not set up for ncoeffs > nqtot");
1109
1110 const Array<OneD, const NekDouble> &base0 = m_base[0]->GetBdata();
1111 const Array<OneD, const NekDouble> &base1 = m_base[1]->GetBdata();
1112 const Array<OneD, const NekDouble> &base2 = m_base[2]->GetBdata();
1113 const Array<OneD, const NekDouble> &dbase0 = m_base[0]->GetDbdata();
1114 const Array<OneD, const NekDouble> &dbase1 = m_base[1]->GetDbdata();
1115 const Array<OneD, const NekDouble> &dbase2 = m_base[2]->GetDbdata();
1116 const Array<OneD, const NekDouble> &metric00 =
1117 m_metrics[eMetricLaplacian00];
1118 const Array<OneD, const NekDouble> &metric01 =
1119 m_metrics[eMetricLaplacian01];
1120 const Array<OneD, const NekDouble> &metric02 =
1121 m_metrics[eMetricLaplacian02];
1122 const Array<OneD, const NekDouble> &metric11 =
1123 m_metrics[eMetricLaplacian11];
1124 const Array<OneD, const NekDouble> &metric12 =
1125 m_metrics[eMetricLaplacian12];
1126 const Array<OneD, const NekDouble> &metric22 =
1127 m_metrics[eMetricLaplacian22];
1128
1129 // Allocate temporary storage
1130 Array<OneD, NekDouble> wsp0(2 * nqtot, wsp);
1131 Array<OneD, NekDouble> wsp1(nqtot, wsp + 1 * nqtot);
1132 Array<OneD, NekDouble> wsp2(nqtot, wsp + 2 * nqtot);
1133 Array<OneD, NekDouble> wsp3(nqtot, wsp + 3 * nqtot);
1134 Array<OneD, NekDouble> wsp4(nqtot, wsp + 4 * nqtot);
1135 Array<OneD, NekDouble> wsp5(nqtot, wsp + 5 * nqtot);
1136
1137 // LAPLACIAN MATRIX OPERATION
1138 // wsp1 = du_dxi1 = D_xi1 * inarray = D_xi1 * u
1139 // wsp2 = du_dxi2 = D_xi2 * inarray = D_xi2 * u
1140 // wsp2 = du_dxi3 = D_xi3 * inarray = D_xi3 * u
1141 StdExpansion3D::PhysTensorDeriv(inarray, wsp0, wsp1, wsp2);
1142
1143 // wsp0 = k = g0 * wsp1 + g1 * wsp2 = g0 * du_dxi1 + g1 * du_dxi2
1144 // wsp2 = l = g1 * wsp1 + g2 * wsp2 = g0 * du_dxi1 + g1 * du_dxi2
1145 // where g0, g1 and g2 are the metric terms set up in the GeomFactors class
1146 // especially for this purpose
1147 Vmath::Vvtvvtp(nqtot, &metric00[0], 1, &wsp0[0], 1, &metric01[0], 1,
1148 &wsp1[0], 1, &wsp3[0], 1);
1149 Vmath::Vvtvp(nqtot, &metric02[0], 1, &wsp2[0], 1, &wsp3[0], 1, &wsp3[0], 1);
1150 Vmath::Vvtvvtp(nqtot, &metric01[0], 1, &wsp0[0], 1, &metric11[0], 1,
1151 &wsp1[0], 1, &wsp4[0], 1);
1152 Vmath::Vvtvp(nqtot, &metric12[0], 1, &wsp2[0], 1, &wsp4[0], 1, &wsp4[0], 1);
1153 Vmath::Vvtvvtp(nqtot, &metric02[0], 1, &wsp0[0], 1, &metric12[0], 1,
1154 &wsp1[0], 1, &wsp5[0], 1);
1155 Vmath::Vvtvp(nqtot, &metric22[0], 1, &wsp2[0], 1, &wsp5[0], 1, &wsp5[0], 1);
1156
1157 // outarray = m = (D_xi1 * B)^T * k
1158 // wsp1 = n = (D_xi2 * B)^T * l
1159 IProductWRTBase_SumFacKernel(dbase0, base1, base2, wsp3, outarray, wsp0,
1160 false, true, true);
1161 IProductWRTBase_SumFacKernel(base0, dbase1, base2, wsp4, wsp2, wsp0, true,
1162 false, true);
1163 Vmath::Vadd(m_ncoeffs, wsp2.data(), 1, outarray.data(), 1, outarray.data(),
1164 1);
1165 IProductWRTBase_SumFacKernel(base0, base1, dbase2, wsp5, wsp2, wsp0, true,
1166 true, false);
1167 Vmath::Vadd(m_ncoeffs, wsp2.data(), 1, outarray.data(), 1, outarray.data(),
1168 1);
1169}

References ASSERTL1, Nektar::LocalRegions::Expansion::ComputeLaplacianMetric(), Nektar::LocalRegions::eMetricLaplacian00, Nektar::StdRegions::StdExpansion::m_base, and Nektar::LocalRegions::Expansion::m_metrics.

◆ v_NormalTraceDerivFactors()

void Nektar::LocalRegions::TetExp::v_NormalTraceDerivFactors ( Array< OneD, Array< OneD, NekDouble > > &  d0factors,
Array< OneD, Array< OneD, NekDouble > > &  d1factors,
Array< OneD, Array< OneD, NekDouble > > &  d2factors 
)
overrideprotectedvirtual

: This method gets all of the factors which are required as part of the Gradient Jump Penalty stabilisation and involves the product of the normal and geometric factors along the element trace.

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 1374 of file TetExp.cpp.

1378{
1379 int nquad0 = GetNumPoints(0);
1380 int nquad1 = GetNumPoints(1);
1381 int nquad2 = GetNumPoints(2);
1382
1383 const Array<TwoD, const NekDouble> &df =
1384 m_metricinfo->GetDerivFactors(GetPointsKeys());
1385
1386 if (d0factors.size() != 4)
1387 {
1388 d0factors = Array<OneD, Array<OneD, NekDouble>>(4);
1389 d1factors = Array<OneD, Array<OneD, NekDouble>>(4);
1390 d2factors = Array<OneD, Array<OneD, NekDouble>>(4);
1391 }
1392
1393 if (d0factors[0].size() != nquad0 * nquad1)
1394 {
1395 d0factors[0] = Array<OneD, NekDouble>(nquad0 * nquad1);
1396 d1factors[0] = Array<OneD, NekDouble>(nquad0 * nquad1);
1397 d2factors[0] = Array<OneD, NekDouble>(nquad0 * nquad1);
1398 }
1399
1400 if (d0factors[1].size() != nquad0 * nquad2)
1401 {
1402 d0factors[1] = Array<OneD, NekDouble>(nquad0 * nquad2);
1403 d1factors[1] = Array<OneD, NekDouble>(nquad0 * nquad2);
1404 d2factors[1] = Array<OneD, NekDouble>(nquad0 * nquad2);
1405 }
1406
1407 if (d0factors[2].size() != nquad1 * nquad2)
1408 {
1409 d0factors[2] = Array<OneD, NekDouble>(nquad1 * nquad2);
1410 d0factors[3] = Array<OneD, NekDouble>(nquad1 * nquad2);
1411 d1factors[2] = Array<OneD, NekDouble>(nquad1 * nquad2);
1412 d1factors[3] = Array<OneD, NekDouble>(nquad1 * nquad2);
1413 d2factors[2] = Array<OneD, NekDouble>(nquad1 * nquad2);
1414 d2factors[3] = Array<OneD, NekDouble>(nquad1 * nquad2);
1415 }
1416
1417 // Outwards normals
1418 const Array<OneD, const Array<OneD, NekDouble>> &normal_0 =
1419 GetTraceNormal(0);
1420 const Array<OneD, const Array<OneD, NekDouble>> &normal_1 =
1421 GetTraceNormal(1);
1422 const Array<OneD, const Array<OneD, NekDouble>> &normal_2 =
1423 GetTraceNormal(2);
1424 const Array<OneD, const Array<OneD, NekDouble>> &normal_3 =
1425 GetTraceNormal(3);
1426
1427 int ncoords = normal_0.size();
1428
1429 // first gather together standard cartesian inner products
1430 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
1431 {
1432 // face 0
1433 for (int i = 0; i < nquad0 * nquad1; ++i)
1434 {
1435 d0factors[0][i] = df[0][i] * normal_0[0][i];
1436 d1factors[0][i] = df[1][i] * normal_0[0][i];
1437 d2factors[0][i] = df[2][i] * normal_0[0][i];
1438 }
1439
1440 for (int n = 1; n < ncoords; ++n)
1441 {
1442 for (int i = 0; i < nquad0 * nquad1; ++i)
1443 {
1444 d0factors[0][i] += df[3 * n][i] * normal_0[n][i];
1445 d1factors[0][i] += df[3 * n + 1][i] * normal_0[n][i];
1446 d2factors[0][i] += df[3 * n + 2][i] * normal_0[n][i];
1447 }
1448 }
1449
1450 // face 1
1451 for (int j = 0; j < nquad2; ++j)
1452 {
1453 for (int i = 0; i < nquad0; ++i)
1454 {
1455 d0factors[1][j * nquad0 + i] = df[0][j * nquad0 * nquad1 + i] *
1456 normal_1[0][j * nquad0 + i];
1457 d1factors[1][j * nquad0 + i] = df[1][j * nquad0 * nquad1 + i] *
1458 normal_1[0][j * nquad0 + i];
1459 d2factors[1][j * nquad0 + i] = df[2][j * nquad0 * nquad1 + i] *
1460 normal_1[0][j * nquad0 + i];
1461 }
1462 }
1463
1464 for (int n = 1; n < ncoords; ++n)
1465 {
1466 for (int j = 0; j < nquad2; ++j)
1467 {
1468 for (int i = 0; i < nquad0; ++i)
1469 {
1470 d0factors[1][j * nquad0 + i] +=
1471 df[3 * n][j * nquad0 * nquad1 + i] *
1472 normal_1[0][j * nquad0 + i];
1473 d1factors[1][j * nquad0 + i] +=
1474 df[3 * n + 1][j * nquad0 * nquad1 + i] *
1475 normal_1[0][j * nquad0 + i];
1476 d2factors[1][j * nquad0 + i] +=
1477 df[3 * n + 2][j * nquad0 * nquad1 + i] *
1478 normal_1[0][j * nquad0 + i];
1479 }
1480 }
1481 }
1482
1483 // faces 2 and 3
1484 for (int j = 0; j < nquad2; ++j)
1485 {
1486 for (int i = 0; i < nquad1; ++i)
1487 {
1488 d0factors[2][j * nquad1 + i] =
1489 df[0][(j * nquad1 + i + 1) * nquad0 - 1] *
1490 normal_2[0][j * nquad1 + i];
1491 d1factors[2][j * nquad1 + i] =
1492 df[1][(j * nquad1 + i + 1) * nquad0 - 1] *
1493 normal_2[0][j * nquad1 + i];
1494 d2factors[2][j * nquad1 + i] =
1495 df[2][(j * nquad1 + i + 1) * nquad0 - 1] *
1496 normal_2[0][j * nquad1 + i];
1497
1498 d0factors[3][j * nquad1 + i] =
1499 df[0][(j * nquad1 + i) * nquad0] *
1500 normal_3[0][j * nquad1 + i];
1501 d1factors[3][j * nquad1 + i] =
1502 df[1][(j * nquad1 + i) * nquad0] *
1503 normal_3[0][j * nquad1 + i];
1504 d2factors[3][j * nquad1 + i] =
1505 df[2][(j * nquad1 + i) * nquad0] *
1506 normal_3[0][j * nquad1 + i];
1507 }
1508 }
1509
1510 for (int n = 1; n < ncoords; ++n)
1511 {
1512 for (int j = 0; j < nquad2; ++j)
1513 {
1514 for (int i = 0; i < nquad1; ++i)
1515 {
1516 d0factors[2][j * nquad1 + i] +=
1517 df[3 * n][(j * nquad1 + i + 1) * nquad0 - 1] *
1518 normal_2[n][j * nquad1 + i];
1519 d1factors[2][j * nquad1 + i] +=
1520 df[3 * n + 1][(j * nquad1 + i + 1) * nquad0 - 1] *
1521 normal_2[n][j * nquad1 + i];
1522 d2factors[2][j * nquad1 + i] +=
1523 df[3 * n + 2][(j * nquad1 + i + 1) * nquad0 - 1] *
1524 normal_2[n][j * nquad1 + i];
1525
1526 d0factors[3][j * nquad1 + i] +=
1527 df[3 * n][(j * nquad1 + i) * nquad0] *
1528 normal_3[n][j * nquad1 + i];
1529 d1factors[3][j * nquad1 + i] +=
1530 df[3 * n + 1][(j * nquad1 + i) * nquad0] *
1531 normal_3[n][j * nquad1 + i];
1532 d2factors[3][j * nquad1 + i] +=
1533 df[3 * n + 2][(j * nquad1 + i) * nquad0] *
1534 normal_3[n][j * nquad1 + i];
1535 }
1536 }
1537 }
1538 }
1539 else
1540 {
1541 // Face 0
1542 for (int i = 0; i < nquad0 * nquad1; ++i)
1543 {
1544 d0factors[0][i] = df[0][0] * normal_0[0][i];
1545 d1factors[0][i] = df[1][0] * normal_0[0][i];
1546 d2factors[0][i] = df[2][0] * normal_0[0][i];
1547 }
1548
1549 for (int n = 1; n < ncoords; ++n)
1550 {
1551 for (int i = 0; i < nquad0 * nquad1; ++i)
1552 {
1553 d0factors[0][i] += df[3 * n][0] * normal_0[n][i];
1554 d1factors[0][i] += df[3 * n + 1][0] * normal_0[n][i];
1555 d2factors[0][i] += df[3 * n + 2][0] * normal_0[n][i];
1556 }
1557 }
1558
1559 // face 1
1560 for (int i = 0; i < nquad0 * nquad2; ++i)
1561 {
1562 d0factors[1][i] = df[0][0] * normal_1[0][i];
1563 d1factors[1][i] = df[1][0] * normal_1[0][i];
1564 d2factors[1][i] = df[2][0] * normal_1[0][i];
1565 }
1566
1567 for (int n = 1; n < ncoords; ++n)
1568 {
1569 for (int i = 0; i < nquad0 * nquad2; ++i)
1570 {
1571 d0factors[1][i] += df[3 * n][0] * normal_1[n][i];
1572 d1factors[1][i] += df[3 * n + 1][0] * normal_1[n][i];
1573 d2factors[1][i] += df[3 * n + 2][0] * normal_1[n][i];
1574 }
1575 }
1576
1577 // faces 2 and 3
1578 for (int i = 0; i < nquad1 * nquad2; ++i)
1579 {
1580 d0factors[2][i] = df[0][0] * normal_2[0][i];
1581 d0factors[3][i] = df[0][0] * normal_3[0][i];
1582
1583 d1factors[2][i] = df[1][0] * normal_2[0][i];
1584 d1factors[3][i] = df[1][0] * normal_3[0][i];
1585
1586 d2factors[2][i] = df[2][0] * normal_2[0][i];
1587 d2factors[3][i] = df[2][0] * normal_3[0][i];
1588 }
1589
1590 for (int n = 1; n < ncoords; ++n)
1591 {
1592 for (int i = 0; i < nquad1 * nquad2; ++i)
1593 {
1594 d0factors[2][i] += df[3 * n][0] * normal_2[n][i];
1595 d0factors[3][i] += df[3 * n][0] * normal_3[n][i];
1596
1597 d1factors[2][i] += df[3 * n + 1][0] * normal_2[n][i];
1598 d1factors[3][i] += df[3 * n + 1][0] * normal_3[n][i];
1599
1600 d2factors[2][i] += df[3 * n + 2][0] * normal_2[n][i];
1601 d2factors[3][i] += df[3 * n + 2][0] * normal_3[n][i];
1602 }
1603 }
1604 }
1605}
const NormalVector & GetTraceNormal(const int id)
Definition: Expansion.cpp:251

◆ v_PhysDeriv()

void Nektar::LocalRegions::TetExp::v_PhysDeriv ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  out_d0,
Array< OneD, NekDouble > &  out_d1,
Array< OneD, NekDouble > &  out_d2 
)
overrideprotectedvirtual

Differentiate inarray in the three coordinate directions.

Parameters
inarrayInput array of values at quadrature points to be differentiated.
out_d0Derivative in first coordinate direction.
out_d1Derivative in second coordinate direction.
out_d2Derivative in third coordinate direction.

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 141 of file TetExp.cpp.

145{
146 int TotPts = m_base[0]->GetNumPoints() * m_base[1]->GetNumPoints() *
147 m_base[2]->GetNumPoints();
148
149 Array<TwoD, const NekDouble> df =
150 m_metricinfo->GetDerivFactors(GetPointsKeys());
151 Array<OneD, NekDouble> Diff0 = Array<OneD, NekDouble>(3 * TotPts);
152 Array<OneD, NekDouble> Diff1 = Diff0 + TotPts;
153 Array<OneD, NekDouble> Diff2 = Diff1 + TotPts;
154
155 StdTetExp::v_PhysDeriv(inarray, Diff0, Diff1, Diff2);
156
157 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
158 {
159 if (out_d0.size())
160 {
161 Vmath::Vmul(TotPts, &df[0][0], 1, &Diff0[0], 1, &out_d0[0], 1);
162 Vmath::Vvtvp(TotPts, &df[1][0], 1, &Diff1[0], 1, &out_d0[0], 1,
163 &out_d0[0], 1);
164 Vmath::Vvtvp(TotPts, &df[2][0], 1, &Diff2[0], 1, &out_d0[0], 1,
165 &out_d0[0], 1);
166 }
167
168 if (out_d1.size())
169 {
170 Vmath::Vmul(TotPts, &df[3][0], 1, &Diff0[0], 1, &out_d1[0], 1);
171 Vmath::Vvtvp(TotPts, &df[4][0], 1, &Diff1[0], 1, &out_d1[0], 1,
172 &out_d1[0], 1);
173 Vmath::Vvtvp(TotPts, &df[5][0], 1, &Diff2[0], 1, &out_d1[0], 1,
174 &out_d1[0], 1);
175 }
176
177 if (out_d2.size())
178 {
179 Vmath::Vmul(TotPts, &df[6][0], 1, &Diff0[0], 1, &out_d2[0], 1);
180 Vmath::Vvtvp(TotPts, &df[7][0], 1, &Diff1[0], 1, &out_d2[0], 1,
181 &out_d2[0], 1);
182 Vmath::Vvtvp(TotPts, &df[8][0], 1, &Diff2[0], 1, &out_d2[0], 1,
183 &out_d2[0], 1);
184 }
185 }
186 else // regular geometry
187 {
188 if (out_d0.size())
189 {
190 Vmath::Smul(TotPts, df[0][0], &Diff0[0], 1, &out_d0[0], 1);
191 Blas::Daxpy(TotPts, df[1][0], &Diff1[0], 1, &out_d0[0], 1);
192 Blas::Daxpy(TotPts, df[2][0], &Diff2[0], 1, &out_d0[0], 1);
193 }
194
195 if (out_d1.size())
196 {
197 Vmath::Smul(TotPts, df[3][0], &Diff0[0], 1, &out_d1[0], 1);
198 Blas::Daxpy(TotPts, df[4][0], &Diff1[0], 1, &out_d1[0], 1);
199 Blas::Daxpy(TotPts, df[5][0], &Diff2[0], 1, &out_d1[0], 1);
200 }
201
202 if (out_d2.size())
203 {
204 Vmath::Smul(TotPts, df[6][0], &Diff0[0], 1, &out_d2[0], 1);
205 Blas::Daxpy(TotPts, df[7][0], &Diff1[0], 1, &out_d2[0], 1);
206 Blas::Daxpy(TotPts, df[8][0], &Diff2[0], 1, &out_d2[0], 1);
207 }
208 }
209}
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:135

References Blas::Daxpy(), Nektar::SpatialDomains::eDeformed, Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::StdRegions::StdExpansion::m_base, Nektar::LocalRegions::Expansion::m_metricinfo, Vmath::Smul(), Vmath::Vmul(), and Vmath::Vvtvp().

◆ v_PhysEvalFirstDeriv()

NekDouble Nektar::LocalRegions::TetExp::v_PhysEvalFirstDeriv ( const Array< OneD, NekDouble > &  coord,
const Array< OneD, const NekDouble > &  inarray,
std::array< NekDouble, 3 > &  firstOrderDerivs 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 485 of file TetExp.cpp.

489{
490 Array<OneD, NekDouble> Lcoord(3);
491 ASSERTL0(m_geom, "m_geom not defined");
492 m_geom->GetLocCoords(coord, Lcoord);
493 return StdTetExp::v_PhysEvalFirstDeriv(Lcoord, inarray, firstOrderDerivs);
494}

References ASSERTL0, and Nektar::LocalRegions::Expansion::m_geom.

◆ v_PhysEvaluate()

NekDouble Nektar::LocalRegions::TetExp::v_PhysEvaluate ( const Array< OneD, const NekDouble > &  coord,
const Array< OneD, const NekDouble > &  physvals 
)
overrideprotectedvirtual
Parameters
coordPhysical space coordinate
Returns
Evaluation of expansion at given coordinate.

Reimplemented from Nektar::StdRegions::StdExpansion3D.

Definition at line 471 of file TetExp.cpp.

473{
474 ASSERTL0(m_geom, "m_geom not defined");
475
476 Array<OneD, NekDouble> Lcoord = Array<OneD, NekDouble>(3);
477
478 // Get the local (eta) coordinates of the point
479 m_geom->GetLocCoords(coord, Lcoord);
480
481 // Evaluate point in local (eta) coordinates.
482 return StdExpansion3D::v_PhysEvaluate(Lcoord, physvals);
483}

References ASSERTL0, and Nektar::LocalRegions::Expansion::m_geom.

◆ v_StdPhysEvaluate()

NekDouble Nektar::LocalRegions::TetExp::v_StdPhysEvaluate ( const Array< OneD, const NekDouble > &  Lcoord,
const Array< OneD, const NekDouble > &  physvals 
)
overrideprotectedvirtual

Given the local cartesian coordinate Lcoord evaluate the value of physvals at this point by calling through to the StdExpansion method

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 459 of file TetExp.cpp.

462{
463 // Evaluate point in local (eta) coordinates.
464 return StdExpansion3D::v_PhysEvaluate(Lcoord, physvals);
465}

◆ v_SVVLaplacianFilter()

void Nektar::LocalRegions::TetExp::v_SVVLaplacianFilter ( Array< OneD, NekDouble > &  array,
const StdRegions::StdMatrixKey mkey 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdTetExp.

Definition at line 985 of file TetExp.cpp.

987{
988 int nq = GetTotPoints();
989
990 // Calculate sqrt of the Jacobian
991 Array<OneD, const NekDouble> jac = m_metricinfo->GetJac(GetPointsKeys());
992 Array<OneD, NekDouble> sqrt_jac(nq);
993 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
994 {
995 Vmath::Vsqrt(nq, jac, 1, sqrt_jac, 1);
996 }
997 else
998 {
999 Vmath::Fill(nq, sqrt(jac[0]), sqrt_jac, 1);
1000 }
1001
1002 // Multiply array by sqrt(Jac)
1003 Vmath::Vmul(nq, sqrt_jac, 1, array, 1, array, 1);
1004
1005 // Apply std region filter
1006 StdTetExp::v_SVVLaplacianFilter(array, mkey);
1007
1008 // Divide by sqrt(Jac)
1009 Vmath::Vdiv(nq, array, 1, sqrt_jac, 1, array, 1);
1010}
void Vdiv(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.hpp:126

References Nektar::SpatialDomains::eDeformed, Vmath::Fill(), Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::StdRegions::StdExpansion::GetTotPoints(), Nektar::LocalRegions::Expansion::m_metricinfo, tinysimd::sqrt(), Vmath::Vdiv(), Vmath::Vmul(), and Vmath::Vsqrt().

Member Data Documentation

◆ m_matrixManager

LibUtilities::NekManager<MatrixKey, DNekScalMat, MatrixKey::opLess> Nektar::LocalRegions::TetExp::m_matrixManager
private

Definition at line 200 of file TetExp.h.

Referenced by v_DropLocMatrix(), v_FwdTrans(), and v_GetLocMatrix().

◆ m_staticCondMatrixManager

LibUtilities::NekManager<MatrixKey, DNekScalBlkMat, MatrixKey::opLess> Nektar::LocalRegions::TetExp::m_staticCondMatrixManager
private

Definition at line 202 of file TetExp.h.

Referenced by v_DropLocStaticCondMatrix(), and v_GetLocStaticCondMatrix().