Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
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...
 
virtual ~TetExp () override=default
 
- Public Member Functions inherited from Nektar::StdRegions::StdTetExp
 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)
 
virtual ~StdTetExp () override
 
LibUtilities::ShapeType DetShapeType () const
 
NekDouble PhysEvaluate3D (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
 Single Point Evaluation. More...
 
- Public Member Functions inherited from Nektar::StdRegions::StdExpansion3D
 StdExpansion3D ()
 
 StdExpansion3D (int numcoeffs, const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc)
 
 StdExpansion3D (const StdExpansion3D &T)
 
virtual ~StdExpansion3D () override
 
void PhysTensorDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray_d1, Array< OneD, NekDouble > &outarray_d2, Array< OneD, NekDouble > &outarray_d3)
 Calculate the 3D derivative in the local tensor/collapsed coordinate at the physical points. More...
 
void BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
 
void IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
 
int GetNedges () const
 return the number of edges in 3D expansion More...
 
int GetEdgeNcoeffs (const int i) const
 This function returns the number of expansion coefficients belonging to the i-th edge. More...
 
void GetEdgeInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards)
 
- Public Member Functions inherited from Nektar::StdRegions::StdExpansion
 StdExpansion ()
 Default Constructor. More...
 
 StdExpansion (const int numcoeffs, const int numbases, const LibUtilities::BasisKey &Ba=LibUtilities::NullBasisKey, const LibUtilities::BasisKey &Bb=LibUtilities::NullBasisKey, const LibUtilities::BasisKey &Bc=LibUtilities::NullBasisKey)
 Constructor. More...
 
 StdExpansion (const StdExpansion &T)
 Copy Constructor. More...
 
virtual ~StdExpansion ()
 Destructor. More...
 
int GetNumBases () const
 This function returns the number of 1D bases used in the expansion. More...
 
const Array< OneD, const LibUtilities::BasisSharedPtr > & GetBase () const
 This function gets the shared point to basis. More...
 
const LibUtilities::BasisSharedPtrGetBasis (int dir) const
 This function gets the shared point to basis in the dir direction. More...
 
int GetNcoeffs (void) const
 This function returns the total number of coefficients used in the expansion. More...
 
int GetTotPoints () const
 This function returns the total number of quadrature points used in the element. More...
 
LibUtilities::BasisType GetBasisType (const int dir) const
 This function returns the type of basis used in the dir direction. More...
 
int GetBasisNumModes (const int dir) const
 This function returns the number of expansion modes in the dir direction. More...
 
int EvalBasisNumModesMax (void) const
 This function returns the maximum number of expansion modes over all local directions. More...
 
LibUtilities::PointsType GetPointsType (const int dir) const
 This function returns the type of quadrature points used in the dir direction. More...
 
int GetNumPoints (const int dir) const
 This function returns the number of quadrature points in the dir direction. More...
 
const Array< OneD, const NekDouble > & GetPoints (const int dir) const
 This function returns a pointer to the array containing the quadrature points in dir direction. More...
 
int GetNverts () const
 This function returns the number of vertices of the expansion domain. More...
 
int GetTraceNcoeffs (const int i) const
 This function returns the number of expansion coefficients belonging to the i-th trace. More...
 
int GetTraceIntNcoeffs (const int i) const
 
int GetTraceNumPoints (const int i) const
 This function returns the number of quadrature points belonging to the i-th trace. More...
 
const LibUtilities::BasisKey GetTraceBasisKey (const int i, int k=-1) const
 This function returns the basis key belonging to the i-th trace. More...
 
LibUtilities::PointsKey GetTracePointsKey (const int i, int k=-1) const
 This function returns the basis key belonging to the i-th trace. More...
 
int NumBndryCoeffs (void) const
 
int NumDGBndryCoeffs (void) const
 
const LibUtilities::PointsKey GetNodalPointsKey () const
 This function returns the type of expansion Nodal point type if defined. More...
 
int GetNtraces () const
 Returns the number of trace elements connected to this element. More...
 
LibUtilities::ShapeType DetShapeType () const
 This function returns the shape of the expansion domain. More...
 
std::shared_ptr< StdExpansionGetStdExp () const
 
std::shared_ptr< StdExpansionGetLinStdExp (void) const
 
int GetShapeDimension () const
 
bool IsBoundaryInteriorExpansion () const
 
bool IsNodalNonTensorialExp ()
 
void BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs the Backward transformation from coefficient space to physical space. More...
 
void FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs the Forward transformation from physical space to coefficient space. More...
 
void 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...
 
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)
 
virtual ~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)
 
virtual ~Expansion ()
 
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

virtual NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray) override
 Integrate the physical point list inarray over region. More...
 
virtual void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2) override
 Differentiate inarray in the three coordinate directions. More...
 
virtual 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...
 
virtual 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...
 
virtual void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
virtual 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...
 
virtual void v_AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
virtual NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals) override
 
virtual NekDouble v_PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals) override
 
virtual NekDouble v_PhysEvaluate (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) override
 
virtual void v_GetCoord (const Array< OneD, const NekDouble > &Lcoords, Array< OneD, NekDouble > &coords) override
 Get the coordinates "coords" at the local coordinates "Lcoords". More...
 
virtual void v_GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
 
virtual LibUtilities::ShapeType v_DetShapeType () const override
 Return Shape of region, using ShapeType enum list. More...
 
virtual StdRegions::StdExpansionSharedPtr v_GetStdExp (void) const override
 
virtual StdRegions::StdExpansionSharedPtr v_GetLinStdExp (void) const override
 
virtual void v_ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int mode_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType) override
 
virtual 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...
 
virtual void v_HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
virtual void v_LaplacianMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
virtual void v_LaplacianMatrixOp (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
virtual void v_SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdRegions::StdMatrixKey &mkey) override
 
virtual DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
virtual DNekMatSharedPtr v_CreateStdMatrix (const StdRegions::StdMatrixKey &mkey) override
 
virtual DNekScalMatSharedPtr v_GetLocMatrix (const MatrixKey &mkey) override
 
void v_DropLocMatrix (const MatrixKey &mkey) override
 
virtual 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)
 
virtual void v_ComputeLaplacianMetric () override
 
virtual 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...
 
virtual void v_PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_StdPhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2) override
 
virtual 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
 
virtual void v_BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) 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) override
 
virtual 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
 
virtual void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
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) override
 
virtual void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_LocCoordToLocCollapsed (const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta) override
 
virtual void v_LocCollapsedToLocCoord (const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi) override
 
virtual void v_GetCoords (Array< OneD, NekDouble > &coords_x, Array< OneD, NekDouble > &coords_y, Array< OneD, NekDouble > &coords_z) override
 
virtual void v_FillMode (const int mode, Array< OneD, NekDouble > &outarray) override
 
NekDouble v_PhysEvaluateBasis (const Array< OneD, const NekDouble > &coords, int mode) final override
 
NekDouble v_PhysEvaluate (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
 
virtual int v_GetNverts () const override
 
virtual int v_GetNedges () const override
 
virtual int v_GetNtraces () const override
 
virtual LibUtilities::ShapeType v_DetShapeType () const override
 
virtual int v_NumBndryCoeffs () const override
 
virtual int v_NumDGBndryCoeffs () const override
 
virtual int v_GetTraceNcoeffs (const int i) const override
 
virtual int v_GetTraceIntNcoeffs (const int i) const override
 
virtual int v_GetTraceNumPoints (const int i) const override
 
virtual int v_GetEdgeNcoeffs (const int i) const override
 
virtual LibUtilities::PointsKey v_GetTracePointsKey (const int i, const int j) const override
 
virtual int v_CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset) override
 
virtual const LibUtilities::BasisKey v_GetTraceBasisKey (const int i, const int k) const override
 
virtual bool v_IsBoundaryInteriorExpansion () const override
 
virtual int v_GetVertexMap (int localVertexId, bool useCoeffPacking=false) override
 
virtual void v_GetInteriorMap (Array< OneD, unsigned int > &outarray) override
 
virtual void v_GetBoundaryMap (Array< OneD, unsigned int > &outarray) override
 
virtual 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
 
virtual DNekMatSharedPtr v_GenMatrix (const StdMatrixKey &mkey) override
 
virtual DNekMatSharedPtr v_CreateStdMatrix (const StdMatrixKey &mkey) override
 
void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual 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
 
virtual void v_GetSimplexEquiSpacedConnectivity (Array< OneD, int > &conn, bool standard=true) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion3D
virtual 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...
 
virtual NekDouble v_PhysEvaluate (const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals) override
 
virtual NekDouble v_PhysEvaluate (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) 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
 
virtual void v_LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
virtual void v_HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
virtual 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)
 
virtual void v_GetTraceToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient, int P, int Q) override
 
virtual void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat) 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 void v_SetCoeffsToOrientation (Array< OneD, NekDouble > &coeffs, StdRegions::Orientation dir)
 
virtual NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
 
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
virtual 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...
 
virtual DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
virtual void v_AddFaceNormBoundaryInt (const int face, const ExpansionSharedPtr &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray) override
 
virtual void v_AddRobinMassMatrix (const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat) override
 
virtual StdRegions::Orientation v_GetTraceOrient (int face) override
 
virtual 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...
 
virtual 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)
 
virtual DNekMatSharedPtr v_BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType) override
 
virtual DNekMatSharedPtr v_BuildInverseTransformationMatrix (const DNekScalMatSharedPtr &transformationmatrix) override
 Build inverse and inverse transposed transformation matrix: \(\mathbf{R^{-1}}\) and \(\mathbf{R^{-T}}\). More...
 
virtual DNekMatSharedPtr v_BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd) override
 
virtual 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)
 
virtual 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 ()
 
virtual int v_GetCoordim () const override
 
virtual 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)
 
virtual 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)
 
virtual 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 50 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 61 of file TetExp.cpp.

66 Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes()),
67 3, Ba, Bb, Bc),
69 Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes()),
70 Ba, Bb, Bc),
71 StdTetExp(Ba, Bb, Bc), Expansion(geom), Expansion3D(geom),
73 std::bind(&Expansion3D::CreateMatrix, this, std::placeholders::_1),
74 std::string("TetExpMatrix")),
76 this, std::placeholders::_1),
77 std::string("TetExpStaticCondMatrix"))
78{
79}
Expansion3D(SpatialDomains::Geometry3DSharedPtr pGeom)
Definition: Expansion3D.h:61
DNekScalMatSharedPtr CreateMatrix(const MatrixKey &mkey)
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
Definition: Expansion.cpp:277
Expansion(SpatialDomains::GeometrySharedPtr pGeom)
Definition: Expansion.cpp:47
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLess > m_staticCondMatrixManager
Definition: TetExp.h:204
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager
Definition: TetExp.h:202
StdExpansion()
Default Constructor.
int getNumberOfCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:192

◆ TetExp() [2/2]

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

Copy Constructor.

Definition at line 84 of file TetExp.cpp.

85 : StdRegions::StdExpansion(T), StdRegions::StdExpansion3D(T),
86 StdRegions::StdTetExp(T), Expansion(T), Expansion3D(T),
87 m_matrixManager(T.m_matrixManager),
88 m_staticCondMatrixManager(T.m_staticCondMatrixManager)
89{
90}

◆ ~TetExp()

virtual Nektar::LocalRegions::TetExp::~TetExp ( )
overridevirtualdefault

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 1071 of file TetExp.cpp.

1074{
1076
1077 if (inarray.get() == outarray.get())
1078 {
1079 Array<OneD, NekDouble> tmp(m_ncoeffs);
1080 Vmath::Vcopy(m_ncoeffs, inarray.get(), 1, tmp.get(), 1);
1081
1082 Blas::Dgemv('N', m_ncoeffs, m_ncoeffs, mat->Scale(),
1083 (mat->GetOwnedMatrix())->GetPtr().get(), m_ncoeffs,
1084 tmp.get(), 1, 0.0, outarray.get(), 1);
1085 }
1086 else
1087 {
1088 Blas::Dgemv('N', m_ncoeffs, m_ncoeffs, mat->Scale(),
1089 (mat->GetOwnedMatrix())->GetPtr().get(), m_ncoeffs,
1090 inarray.get(), 1, 0.0, outarray.get(), 1);
1091 }
1092}
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:88
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:213
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1191

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 390 of file TetExp.cpp.

393{
394 int i, j;
395
396 const int nquad0 = m_base[0]->GetNumPoints();
397 const int nquad1 = m_base[1]->GetNumPoints();
398 const int nquad2 = m_base[2]->GetNumPoints();
399 const int nqtot = nquad0 * nquad1 * nquad2;
400
401 const Array<OneD, const NekDouble> &z0 = m_base[0]->GetZ();
402 const Array<OneD, const NekDouble> &z1 = m_base[1]->GetZ();
403 const Array<OneD, const NekDouble> &z2 = m_base[2]->GetZ();
404
405 Array<OneD, NekDouble> tmp2(nqtot);
406 Array<OneD, NekDouble> tmp3(nqtot);
407
408 const Array<TwoD, const NekDouble> &df =
409 m_metricinfo->GetDerivFactors(GetPointsKeys());
410
411 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
412 {
413 Vmath::Vmul(nqtot, &df[3 * dir][0], 1, inarray.get(), 1, tmp2.get(), 1);
414 Vmath::Vmul(nqtot, &df[3 * dir + 1][0], 1, inarray.get(), 1, tmp3.get(),
415 1);
416 Vmath::Vmul(nqtot, &df[3 * dir + 2][0], 1, inarray.get(), 1,
417 outarray[2].get(), 1);
418 }
419 else
420 {
421 Vmath::Smul(nqtot, df[3 * dir][0], inarray.get(), 1, tmp2.get(), 1);
422 Vmath::Smul(nqtot, df[3 * dir + 1][0], inarray.get(), 1, tmp3.get(), 1);
423 Vmath::Smul(nqtot, df[3 * dir + 2][0], inarray.get(), 1,
424 outarray[2].get(), 1);
425 }
426
427 NekDouble g0, g1, g1a, g2, g3;
428 int k, cnt;
429
430 for (cnt = 0, k = 0; k < nquad2; ++k)
431 {
432 for (j = 0; j < nquad1; ++j)
433 {
434 g2 = 2.0 / (1.0 - z2[k]);
435 g1 = g2 / (1.0 - z1[j]);
436 g0 = 2.0 * g1;
437 g3 = (1.0 + z1[j]) * g2 * 0.5;
438
439 for (i = 0; i < nquad0; ++i, ++cnt)
440 {
441 g1a = g1 * (1 + z0[i]);
442
443 outarray[0][cnt] =
444 g0 * tmp2[cnt] + g1a * (tmp3[cnt] + outarray[2][cnt]);
445
446 outarray[1][cnt] = g2 * tmp3[cnt] + g3 * outarray[2][cnt];
447 }
448 }
449 }
450}
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
Definition: Expansion.h:276
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.cpp:207
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.cpp:245

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 1172 of file TetExp.cpp.

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

◆ 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 706 of file TetExp.cpp.

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

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 1040 of file TetExp.cpp.

1041{
1042 LibUtilities::BasisKey bkey0 = m_base[0]->GetBasisKey();
1043 LibUtilities::BasisKey bkey1 = m_base[1]->GetBasisKey();
1044 LibUtilities::BasisKey bkey2 = m_base[2]->GetBasisKey();
1047
1048 return tmp->GetStdMatrix(mkey);
1049}
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:255

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 531 of file TetExp.cpp.

532{
534}

References Nektar::LibUtilities::eTetrahedron.

◆ v_DropLocMatrix()

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

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 1056 of file TetExp.cpp.

1057{
1058 m_matrixManager.DeleteObject(mkey);
1059}

References m_matrixManager.

◆ v_DropLocStaticCondMatrix()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1066 of file TetExp.cpp.

1067{
1068 m_staticCondMatrixManager.DeleteObject(mkey);
1069}

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 556 of file TetExp.cpp.

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

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 227 of file TetExp.cpp.

229{
230 if ((m_base[0]->Collocation()) && (m_base[1]->Collocation()) &&
231 (m_base[2]->Collocation()))
232 {
233 Vmath::Vcopy(GetNcoeffs(), &inarray[0], 1, &outarray[0], 1);
234 }
235 else
236 {
237 IProductWRTBase(inarray, outarray);
238
239 // get Mass matrix inverse
240 MatrixKey masskey(StdRegions::eInvMass, DetShapeType(), *this);
241 DNekScalMatSharedPtr matsys = m_matrixManager[masskey];
242
243 // copy inarray in case inarray == outarray
244 DNekVec in(m_ncoeffs, outarray);
245 DNekVec out(m_ncoeffs, outarray, eWrapper);
246
247 out = (*matsys) * in;
248 }
249}
int GetNcoeffs(void) const
This function returns the total number of coefficients used in the expansion.
Definition: StdExpansion.h:130
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:534
LibUtilities::ShapeType DetShapeType() const
Definition: StdTetExp.h:64
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 1018 of file TetExp.cpp.

1019{
1020 DNekMatSharedPtr returnval;
1021
1022 switch (mkey.GetMatrixType())
1023 {
1031 returnval = Expansion3D::v_GenMatrix(mkey);
1032 break;
1033 default:
1034 returnval = StdTetExp::v_GenMatrix(mkey);
1035 }
1036
1037 return returnval;
1038}
virtual 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 500 of file TetExp.cpp.

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

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 517 of file TetExp.cpp.

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

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 543 of file TetExp.cpp.

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

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 1051 of file TetExp.cpp.

1052{
1053 return m_matrixManager[mkey];
1054}

References m_matrixManager.

◆ v_GetLocStaticCondMatrix()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1061 of file TetExp.cpp.

1062{
1063 return m_staticCondMatrixManager[mkey];
1064}

References m_staticCondMatrixManager.

◆ v_GetStdExp()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 536 of file TetExp.cpp.

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

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 618 of file TetExp.cpp.

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

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 966 of file TetExp.cpp.

969{
970 TetExp::v_HelmholtzMatrixOp_MatFree(inarray, outarray, mkey);
971}
virtual 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 106 of file TetExp.cpp.

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

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 280 of file TetExp.cpp.

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

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 286 of file TetExp.cpp.

289{
290 const int nquad0 = m_base[0]->GetNumPoints();
291 const int nquad1 = m_base[1]->GetNumPoints();
292 const int nquad2 = m_base[2]->GetNumPoints();
293 const int order0 = m_base[0]->GetNumModes();
294 const int order1 = m_base[1]->GetNumModes();
295 Array<OneD, NekDouble> wsp(nquad1 * nquad2 * order0 +
296 nquad2 * order0 * (order1 + 1) / 2);
297
298 if (multiplybyweights)
299 {
300 Array<OneD, NekDouble> tmp(nquad0 * nquad1 * nquad2);
301
302 MultiplyByQuadratureMetric(inarray, tmp);
304 m_base[0]->GetBdata(), m_base[1]->GetBdata(), m_base[2]->GetBdata(),
305 tmp, outarray, wsp, true, true, true);
306 }
307 else
308 {
310 m_base[0]->GetBdata(), m_base[1]->GetBdata(), m_base[2]->GetBdata(),
311 inarray, outarray, wsp, true, true, true);
312 }
313}
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 345 of file TetExp.cpp.

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

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 973 of file TetExp.cpp.

976{
977 TetExp::v_LaplacianMatrixOp_MatFree(inarray, outarray, mkey);
978}
virtual 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 980 of file TetExp.cpp.

984{
985 StdExpansion::LaplacianMatrixOp_MatFree(k1, k2, inarray, outarray, mkey);
986}

◆ 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 1094 of file TetExp.cpp.

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

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 1375 of file TetExp.cpp.

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

◆ 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 145 of file TetExp.cpp.

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

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_PhysEvaluate() [1/2]

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 473 of file TetExp.cpp.

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

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

◆ v_PhysEvaluate() [2/2]

NekDouble Nektar::LocalRegions::TetExp::v_PhysEvaluate ( 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 487 of file TetExp.cpp.

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

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 461 of file TetExp.cpp.

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

◆ 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 988 of file TetExp.cpp.

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

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 202 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 204 of file TetExp.h.

Referenced by v_DropLocStaticCondMatrix(), and v_GetLocStaticCondMatrix().