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

#include <NodalTriExp.h>

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

Public Member Functions

 NodalTriExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::PointsType Ntype, const SpatialDomains::TriGeomSharedPtr &geom)
 Constructor using BasisKey class for quadrature points and order definition. More...
 
 NodalTriExp (const NodalTriExp &T)
 
 ~NodalTriExp () override=default
 
- Public Member Functions inherited from Nektar::StdRegions::StdNodalTriExp
 StdNodalTriExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::PointsType Ntype)
 
 StdNodalTriExp ()=default
 
 StdNodalTriExp (const StdNodalTriExp &T)=default
 
 ~StdNodalTriExp () override=default
 
void v_NodalToModal (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void NodalToModalTranspose (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void ModalToNodal (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void GetNodalPoints (Array< OneD, const NekDouble > &x, Array< OneD, const NekDouble > &y)
 
DNekMatSharedPtr GenNBasisTransMatrix ()
 
- Public Member Functions inherited from Nektar::StdRegions::StdTriExp
 StdTriExp (const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb)
 
 StdTriExp ()=default
 
 StdTriExp (const StdTriExp &T)=default
 
 ~StdTriExp () override=default
 
- Public Member Functions inherited from Nektar::StdRegions::StdExpansion2D
 StdExpansion2D (int numcoeffs, const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb)
 
 StdExpansion2D ()=default
 
 StdExpansion2D (const StdExpansion2D &T)=default
 
 ~StdExpansion2D () override=default
 
void PhysTensorDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray_d0, Array< OneD, NekDouble > &outarray_d1)
 Calculate the 2D derivative in the local tensor/collapsed coordinate at the physical points. More...
 
NekDouble Integral (const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &w0, const Array< OneD, const NekDouble > &w1)
 
NekDouble BaryTensorDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs)
 
void BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0=true, bool doCheckCollDir1=true)
 
void IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0=true, bool doCheckCollDir1=true)
 
- Public Member Functions inherited from Nektar::StdRegions::StdExpansion
 StdExpansion ()
 Default Constructor. More...
 
 StdExpansion (const int numcoeffs, const int numbases, const LibUtilities::BasisKey &Ba=LibUtilities::NullBasisKey, const LibUtilities::BasisKey &Bb=LibUtilities::NullBasisKey, const LibUtilities::BasisKey &Bc=LibUtilities::NullBasisKey)
 Constructor. More...
 
 StdExpansion (const StdExpansion &T)
 Copy Constructor. More...
 
virtual ~StdExpansion ()
 Destructor. More...
 
int GetNumBases () const
 This function returns the number of 1D bases used in the expansion. More...
 
const Array< OneD, const LibUtilities::BasisSharedPtr > & GetBase () const
 This function gets the shared point to basis. More...
 
const LibUtilities::BasisSharedPtrGetBasis (int dir) const
 This function gets the shared point to basis in the dir direction. More...
 
int GetNcoeffs (void) const
 This function returns the total number of coefficients used in the expansion. More...
 
int GetTotPoints () const
 This function returns the total number of quadrature points used in the element. More...
 
LibUtilities::BasisType GetBasisType (const int dir) const
 This function returns the type of basis used in the dir direction. More...
 
int GetBasisNumModes (const int dir) const
 This function returns the number of expansion modes in the dir direction. More...
 
int EvalBasisNumModesMax (void) const
 This function returns the maximum number of expansion modes over all local directions. More...
 
LibUtilities::PointsType GetPointsType (const int dir) const
 This function returns the type of quadrature points used in the dir direction. More...
 
int GetNumPoints (const int dir) const
 This function returns the number of quadrature points in the dir direction. More...
 
const Array< OneD, const NekDouble > & GetPoints (const int dir) const
 This function returns a pointer to the array containing the quadrature points in dir direction. More...
 
int GetNverts () const
 This function returns the number of vertices of the expansion domain. More...
 
int GetTraceNcoeffs (const int i) const
 This function returns the number of expansion coefficients belonging to the i-th trace. More...
 
int GetTraceIntNcoeffs (const int i) const
 
int GetTraceNumPoints (const int i) const
 This function returns the number of quadrature points belonging to the i-th trace. More...
 
const LibUtilities::BasisKey GetTraceBasisKey (const int i, int k=-1, bool UseGLL=false) const
 This function returns the basis key belonging to the i-th trace. More...
 
LibUtilities::PointsKey GetTracePointsKey (const int i, int k=-1) const
 This function returns the basis key belonging to the i-th trace. More...
 
int NumBndryCoeffs (void) const
 
int NumDGBndryCoeffs (void) const
 
const LibUtilities::PointsKey GetNodalPointsKey () const
 This function returns the type of expansion Nodal point type if defined. More...
 
int GetNtraces () const
 Returns the number of trace elements connected to this element. More...
 
LibUtilities::ShapeType DetShapeType () const
 This function returns the shape of the expansion domain. More...
 
std::shared_ptr< StdExpansionGetStdExp () const
 
std::shared_ptr< StdExpansionGetLinStdExp (void) const
 
int GetShapeDimension () const
 
bool IsBoundaryInteriorExpansion () const
 
bool IsNodalNonTensorialExp ()
 
void NodalToModal (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs the Backward transformation from coefficient space to physical space. More...
 
void FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs the Forward transformation from physical space to coefficient space. More...
 
void FwdTransBndConstrained (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
NekDouble Integral (const Array< OneD, const NekDouble > &inarray)
 This function integrates the specified function over the domain. More...
 
void FillMode (const int mode, Array< OneD, NekDouble > &outarray)
 This function fills the array outarray with the mode-th mode of the expansion. More...
 
void IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 this function calculates the inner product of a given function f with the different modes of the expansion More...
 
void IProductWRTBase (const Array< OneD, const NekDouble > &base, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, int coll_check)
 
void IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDirectionalDerivBase (const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
int GetElmtId ()
 Get the element id of this expansion when used in a list by returning value of m_elmt_id. More...
 
void SetElmtId (const int id)
 Set the element id of this expansion when used in a list by returning value of m_elmt_id. More...
 
void GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2=NullNekDouble1DArray, Array< OneD, NekDouble > &coords_3=NullNekDouble1DArray)
 this function returns the physical coordinates of the quadrature points of the expansion More...
 
void GetCoord (const Array< OneD, const NekDouble > &Lcoord, Array< OneD, NekDouble > &coord)
 given the coordinates of a point of the element in the local collapsed coordinate system, this function calculates the physical coordinates of the point More...
 
DNekMatSharedPtr GetStdMatrix (const StdMatrixKey &mkey)
 
DNekBlkMatSharedPtr GetStdStaticCondMatrix (const StdMatrixKey &mkey)
 
void NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
 
void NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
void NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, const Array< OneD, const NekDouble > &Fz, Array< OneD, NekDouble > &outarray)
 
void NormVectorIProductWRTBase (const Array< OneD, const Array< OneD, NekDouble > > &Fvec, Array< OneD, NekDouble > &outarray)
 
DNekScalBlkMatSharedPtr GetLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
void DropLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
int CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset)
 
NekDouble StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
 
int GetCoordim ()
 
void GetBoundaryMap (Array< OneD, unsigned int > &outarray)
 
void GetInteriorMap (Array< OneD, unsigned int > &outarray)
 
int GetVertexMap (const int localVertexId, bool useCoeffPacking=false)
 
void GetTraceToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
 
void GetTraceCoeffMap (const unsigned int traceid, Array< OneD, unsigned int > &maparray)
 
void GetElmtTraceToTraceMap (const unsigned int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards, int P=-1, int Q=-1)
 
void GetTraceInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eForwards)
 
void GetTraceNumModes (const int tid, int &numModes0, int &numModes1, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2)
 
void MultiplyByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
DNekMatSharedPtr CreateGeneralMatrix (const StdMatrixKey &mkey)
 this function generates the mass matrix \(\mathbf{M}[i][j] = \int \phi_i(\mathbf{x}) \phi_j(\mathbf{x}) d\mathbf{x}\) More...
 
void GeneralMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void ReduceOrderCoeffs (int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdMatrixKey &mkey)
 
void ExponentialFilter (Array< OneD, NekDouble > &array, const NekDouble alpha, const NekDouble exponent, const NekDouble cutoff)
 
void LaplacianMatrixOp (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDerivMatrixOp (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDirectionalDerivMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassLevelCurvatureMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionDiffusionReactionMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
 
void HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
DNekMatSharedPtr GenMatrix (const StdMatrixKey &mkey)
 
void PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
 
void PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void PhysDeriv_s (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_ds)
 
void PhysDeriv_n (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_dn)
 
void PhysDirectionalDeriv (const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &direction, Array< OneD, NekDouble > &outarray)
 
void StdPhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
 
void StdPhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
NekDouble PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble PhysEvaluate (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs)
 This function evaluates the first derivative of the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble PhysEvaluate (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs, std::array< NekDouble, 6 > &secondOrderDerivs)
 
NekDouble PhysEvaluate (const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals)
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble PhysEvaluateBasis (const Array< OneD, const NekDouble > &coords, int mode)
 This function evaluates the basis function mode mode at a point coords of the domain. More...
 
void LocCoordToLocCollapsed (const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta)
 Convert local cartesian coordinate xi into local collapsed coordinates eta. More...
 
void LocCollapsedToLocCoord (const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi)
 Convert local collapsed coordinates eta into local cartesian coordinate xi. More...
 
void PhysInterp (std::shared_ptr< StdExpansion > fromExp, const Array< OneD, const NekDouble > &fromData, Array< OneD, NekDouble > &toData)
 interpolate from one set of quadrature points available from FromExp to the set of quadrature points in the current expansion. If the points are the same this routine will just copy the data More...
 
virtual int v_CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset)
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, const Array< OneD, const NekDouble > &Fz, Array< OneD, NekDouble > &outarray)
 
virtual void v_NormVectorIProductWRTBase (const Array< OneD, const Array< OneD, NekDouble > > &Fvec, Array< OneD, NekDouble > &outarray)
 
virtual DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual void v_DropLocStaticCondMatrix (const LocalRegions::MatrixKey &mkey)
 
NekDouble Linf (const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
 Function to evaluate the discrete \( L_\infty\) error \( |\epsilon|_\infty = \max |u - u_{exact}|\) where \( u_{exact}\) is given by the array sol. More...
 
NekDouble L2 (const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
 Function to evaluate the discrete \( L_2\) error, \( | \epsilon |_{2} = \left [ \int^1_{-1} [u - u_{exact}]^2 dx \right]^{1/2} d\xi_1 \) where \( u_{exact}\) is given by the array sol. More...
 
NekDouble H1 (const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &sol=NullNekDouble1DArray)
 Function to evaluate the discrete \( H^1\) error, \( | \epsilon |^1_{2} = \left [ \int^1_{-1} [u - u_{exact}]^2 + \nabla(u - u_{exact})\cdot\nabla(u - u_{exact})\cdot dx \right]^{1/2} d\xi_1 \) where \( u_{exact}\) is given by the array sol. More...
 
const LibUtilities::PointsKeyVector GetPointsKeys () const
 
DNekMatSharedPtr BuildInverseTransformationMatrix (const DNekScalMatSharedPtr &m_transformationmatrix)
 
void PhysInterpToSimplexEquiSpaced (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, int npset=-1)
 This function performs an interpolation from the physical space points provided at input into an array of equispaced points which are not the collapsed coordinate. So for a tetrahedron you will only get a tetrahedral number of values. More...
 
void GetSimplexEquiSpacedConnectivity (Array< OneD, int > &conn, bool standard=true)
 This function provides the connectivity of local simplices (triangles or tets) to connect the equispaced data points provided by PhysInterpToSimplexEquiSpaced. More...
 
void EquiSpacedToCoeffs (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function performs a projection/interpolation from the equispaced points sometimes used in post-processing onto the coefficient space. More...
 
template<class T >
std::shared_ptr< T > as ()
 
void IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
 
void GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat)
 
- Public Member Functions inherited from Nektar::LocalRegions::Expansion2D
 Expansion2D (SpatialDomains::Geometry2DSharedPtr pGeom)
 
 ~Expansion2D () override=default
 
DNekScalMatSharedPtr CreateMatrix (const MatrixKey &mkey)
 
void SetTraceToGeomOrientation (Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, NekDouble > &inout)
 
Array< OneD, unsigned int > GetTraceInverseBoundaryMap (int eid)
 
void AddNormTraceInt (const int dir, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &edgeCoeffs, Array< OneD, NekDouble > &outarray)
 
void AddNormTraceInt (const int dir, Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, NekDouble > &outarray, const StdRegions::VarCoeffMap &varcoeffs)
 
void AddEdgeBoundaryInt (const int edge, ExpansionSharedPtr &EdgeExp, Array< OneD, NekDouble > &edgePhys, Array< OneD, NekDouble > &outarray, const StdRegions::VarCoeffMap &varcoeffs=StdRegions::NullVarCoeffMap)
 
void AddHDGHelmholtzEdgeTerms (const NekDouble tau, const int edge, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, NekDouble > &edgePhys, const StdRegions::VarCoeffMap &dirForcing, Array< OneD, NekDouble > &outarray)
 
void AddHDGHelmholtzTraceTerms (const NekDouble tau, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, const StdRegions::VarCoeffMap &dirForcing, Array< OneD, NekDouble > &outarray)
 
SpatialDomains::Geometry2DSharedPtr GetGeom2D () const
 
void ReOrientEdgePhysMap (const int nvert, const StdRegions::Orientation orient, const int nq0, Array< OneD, int > &idmap)
 
DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
void v_GenTraceExp (const int traceid, ExpansionSharedPtr &exp) override
 
- Public Member Functions inherited from Nektar::LocalRegions::Expansion
 Expansion (SpatialDomains::GeometrySharedPtr pGeom)
 
 Expansion (const Expansion &pSrc)
 
 ~Expansion () override
 
void SetTraceExp (const int traceid, ExpansionSharedPtr &f)
 
ExpansionSharedPtr GetTraceExp (const int traceid)
 
DNekScalMatSharedPtr GetLocMatrix (const LocalRegions::MatrixKey &mkey)
 
void DropLocMatrix (const LocalRegions::MatrixKey &mkey)
 
DNekScalMatSharedPtr GetLocMatrix (const StdRegions::MatrixType mtype, const StdRegions::ConstFactorMap &factors=StdRegions::NullConstFactorMap, const StdRegions::VarCoeffMap &varcoeffs=StdRegions::NullVarCoeffMap)
 
SpatialDomains::GeometrySharedPtr GetGeom () const
 
void Reset ()
 
IndexMapValuesSharedPtr CreateIndexMap (const IndexMapKey &ikey)
 
DNekScalBlkMatSharedPtr CreateStaticCondMatrix (const MatrixKey &mkey)
 
const SpatialDomains::GeomFactorsSharedPtrGetMetricInfo () const
 
DNekMatSharedPtr BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
 
DNekMatSharedPtr BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd)
 
void ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
 
void AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
void AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
void AddFaceNormBoundaryInt (const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
void DGDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &coeffs, Array< OneD, NekDouble > &outarray)
 
NekDouble VectorFlux (const Array< OneD, Array< OneD, NekDouble > > &vec)
 
void NormalTraceDerivFactors (Array< OneD, Array< OneD, NekDouble > > &factors, Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors)
 
IndexMapValuesSharedPtr GetIndexMap (const IndexMapKey &ikey)
 
void AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
ExpansionSharedPtr GetLeftAdjacentElementExp () const
 
ExpansionSharedPtr GetRightAdjacentElementExp () const
 
int GetLeftAdjacentElementTrace () const
 
int GetRightAdjacentElementTrace () const
 
void SetAdjacentElementExp (int traceid, ExpansionSharedPtr &e)
 
StdRegions::Orientation GetTraceOrient (int trace)
 
void SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void DivideByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 Divided by the metric jacobi and quadrature weights. More...
 
void GetTraceQFactors (const int trace, Array< OneD, NekDouble > &outarray)
 Extract the metric factors to compute the contravariant fluxes along edge edge and stores them into outarray following the local edge orientation (i.e. anticlockwise convention). More...
 
void GetTracePhysVals (const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient=StdRegions::eNoOrientation)
 
void GetTracePhysMap (const int edge, Array< OneD, int > &outarray)
 
void ReOrientTracePhysMap (const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1)
 
const NormalVectorGetTraceNormal (const int id)
 
void ComputeTraceNormal (const int id)
 
const Array< OneD, const NekDouble > & GetPhysNormals (void)
 
void SetPhysNormals (Array< OneD, const NekDouble > &normal)
 
void SetUpPhysNormals (const int trace)
 
void AddRobinMassMatrix (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
 
void TraceNormLen (const int traceid, NekDouble &h, NekDouble &p)
 
void AddRobinTraceContribution (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
 
const Array< OneD, const NekDouble > & GetElmtBndNormDirElmtLen (const int nbnd) const
 
void StdDerivBaseOnTraceMat (Array< OneD, DNekMatSharedPtr > &DerivMat)
 

Protected Member Functions

NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray) override
 Integrates the specified function over the domain. More...
 
void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray) override
 Calculate the derivative of the physical points. More...
 
void v_PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculate the derivative of the physical points in a given direction. More...
 
void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Transform a given function from physical quadrature space to coefficient space. More...
 
void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculates the inner product of a given function f with the different modes of the expansion. More...
 
void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
void v_IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray) override
 
StdRegions::StdExpansionSharedPtr v_GetStdExp (void) const override
 
StdRegions::StdExpansionSharedPtr v_GetLinStdExp (void) const override
 
void v_GetCoords (Array< OneD, NekDouble > &coords_0, Array< OneD, NekDouble > &coords_1=NullNekDouble1DArray, Array< OneD, NekDouble > &coords_2=NullNekDouble1DArray) override
 
void v_GetCoord (const Array< OneD, const NekDouble > &lcoord, Array< OneD, NekDouble > &coord) override
 
NekDouble v_PhysEvaluate (const Array< OneD, const NekDouble > &coord, const Array< OneD, const NekDouble > &physvals) override
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
void v_ComputeTraceNormal (const int edge) override
 
void v_ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int mode_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType) override
 
void v_GetTracePhysVals (const int edge, const StdRegions::StdExpansionSharedPtr &EdgeExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient) override
 
DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_CreateStdMatrix (const StdRegions::StdMatrixKey &mkey) override
 
DNekScalMatSharedPtr v_GetLocMatrix (const MatrixKey &mkey) override
 
DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix (const MatrixKey &mkey) override
 
void v_DropLocMatrix (const MatrixKey &mkey) override
 
void v_MassMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_LaplacianMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_LaplacianMatrixOp (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_WeakDerivMatrixOp (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdNodalTriExp
const LibUtilities::PointsKey v_GetNodalPointsKey () const override
 
bool v_IsNodalNonTensorialExp () override
 
void v_BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Transform a given function from physical quadrature space to coefficient space. More...
 
void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculates the inner product of a given function f with the different modes of the expansion. More...
 
void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_FillMode (const int mode, Array< OneD, NekDouble > &outarray) override
 
int v_NumBndryCoeffs () const override
 
int v_GetVertexMap (int localVertexId, bool useCoeffPacking=false) override
 
void v_GetTraceToElementMap (const int eid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation edgeOrient=eForwards, int P=-1, int Q=-1) override
 
void v_GetTraceInteriorToElementMap (const int eid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation edgeOrient=eForwards) override
 
void v_GetInteriorMap (Array< OneD, unsigned int > &outarray) override
 
void v_GetBoundaryMap (Array< OneD, unsigned int > &outarray) override
 
void v_MassMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_LaplacianMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_LaplacianMatrixOp (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_WeakDerivMatrixOp (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_GenMatrix (const StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_CreateStdMatrix (const StdMatrixKey &mkey) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdTriExp
NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray) override
 Integrates the specified function over the domain. More...
 
void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray) override
 Calculate the derivative of the physical points. More...
 
void v_PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculate the derivative of the physical points in a given direction. More...
 
void v_StdPhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray) override
 
void v_StdPhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Backward tranform for triangular elements. More...
 
void v_BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1) override
 
void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Transform a given function from physical quadrature space to coefficient space. More...
 
void v_FwdTransBndConstrained (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
 Calculate the inner product of inarray with respect to the basis B=base0[p]*base1[pq] and put into outarray. More...
 
void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
void v_IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1) override
 
void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_LocCoordToLocCollapsed (const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta) override
 
void v_LocCollapsedToLocCoord (const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi) override
 
void v_FillMode (const int mode, Array< OneD, NekDouble > &outarray) override
 
NekDouble v_PhysEvaluateBasis (const Array< OneD, const NekDouble > &coords, int mode) final
 
NekDouble v_PhysEvalFirstDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) override
 
int v_GetNverts () const final
 
int v_GetNtraces () const final
 
LibUtilities::ShapeType v_DetShapeType () const final
 
int v_NumBndryCoeffs () const override
 
int v_NumDGBndryCoeffs () const override
 
int v_GetTraceNcoeffs (const int i) const override
 
int v_GetTraceIntNcoeffs (const int i) const override
 
int v_GetTraceNumPoints (const int i) const override
 
int v_CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset) override
 
void v_GetCoords (Array< OneD, NekDouble > &coords_x, Array< OneD, NekDouble > &coords_y, Array< OneD, NekDouble > &coords_z) override
 
bool v_IsBoundaryInteriorExpansion () const override
 
const LibUtilities::BasisKey v_GetTraceBasisKey (const int i, const int j, bool UseGLL=false) const override
 
int v_GetVertexMap (int localVertexId, bool useCoeffPacking=false) override
 
void v_GetInteriorMap (Array< OneD, unsigned int > &outarray) override
 
void v_GetBoundaryMap (Array< OneD, unsigned int > &outarray) override
 
void v_GetTraceCoeffMap (const unsigned int traceid, Array< OneD, unsigned int > &maparray) override
 
void v_GetTraceInteriorToElementMap (const int eid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation edgeOrient=eForwards) override
 
DNekMatSharedPtr v_GenMatrix (const StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_CreateStdMatrix (const StdMatrixKey &mkey) override
 
void v_MassMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_LaplacianMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdMatrixKey &mkey) override
 
void v_ReduceOrderCoeffs (int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_LaplacianMatrixOp (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_WeakDerivMatrixOp (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey) override
 
void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_GetSimplexEquiSpacedConnectivity (Array< OneD, int > &conn, bool standard=true) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion2D
NekDouble v_PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals) override
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble v_PhysEvaluateInterp (const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals) override
 
virtual void v_BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1)=0
 
virtual void v_IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1)=0
 
void v_LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_GetTraceCoeffMap (const unsigned int traceid, Array< OneD, unsigned int > &maparray) override
 
void v_GetElmtTraceToTraceMap (const unsigned int eid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation edgeOrient, int P, int Q) override
 Determine the mapping to re-orientate the coefficients along the element trace (assumed to align with the standard element) into the orientation of the local trace given by edgeOrient. More...
 
void v_GetTraceToElementMap (const int eid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation edgeOrient=eForwards, int P=-1, int Q=-1) override
 
void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat) override
 
void v_PhysInterp (std::shared_ptr< StdExpansion > fromExp, const Array< OneD, const NekDouble > &fromData, Array< OneD, NekDouble > &toData) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion
DNekMatSharedPtr CreateStdMatrix (const StdMatrixKey &mkey)
 
DNekBlkMatSharedPtr CreateStdStaticCondMatrix (const StdMatrixKey &mkey)
 Create the static condensation of a matrix when using a boundary interior decomposition. More...
 
void BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDirectionalDerivBase_SumFac (const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void GeneralMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree_Kernel (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp)
 
void LaplacianMatrixOp_MatFree_GenericImpl (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDerivMatrixOp_MatFree (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDirectionalDerivMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassLevelCurvatureMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionDiffusionReactionMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
 
void HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void HelmholtzMatrixOp_MatFree_GenericImpl (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
virtual void v_SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
 
virtual void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat)
 
virtual void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals, NekDouble &deriv, NekDouble &deriv2)
 This function performs the barycentric interpolation of the polynomial stored in coord at a point physvals using barycentric interpolation weights in direction. More...
 
template<int DIR>
NekDouble BaryEvaluateBasis (const NekDouble &coord, const int &mode)
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals)
 Helper function to pass an unused value by reference into BaryEvaluate. More...
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals, NekDouble &deriv)
 
- Protected Member Functions inherited from Nektar::LocalRegions::Expansion2D
void v_DGDeriv (const int dir, const Array< OneD, const NekDouble > &incoeffs, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &edgeCoeffs, Array< OneD, NekDouble > &out_d) override
 
void v_AddEdgeNormBoundaryInt (const int edge, const ExpansionSharedPtr &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray) override
 
void v_AddEdgeNormBoundaryInt (const int edge, const ExpansionSharedPtr &EdgeExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray) override
 
void v_AddRobinMassMatrix (const int edgeid, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat) override
 
void v_AddRobinTraceContribution (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs) override
 
DNekMatSharedPtr v_BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd) override
 
void v_ReOrientTracePhysMap (const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1) override
 
void v_SetUpPhysNormals (const int edge) override
 
NekDouble v_VectorFlux (const Array< OneD, Array< OneD, NekDouble > > &vec) override
 
void v_TraceNormLen (const int traceid, NekDouble &h, NekDouble &p) override
 
- Protected Member Functions inherited from Nektar::LocalRegions::Expansion
void ComputeLaplacianMetric ()
 
void ComputeQuadratureMetric ()
 
void ComputeGmatcdotMF (const Array< TwoD, const NekDouble > &df, const Array< OneD, const NekDouble > &direction, Array< OneD, Array< OneD, NekDouble > > &dfdir)
 
Array< OneD, NekDoubleGetMF (const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
 
Array< OneD, NekDoubleGetMFDiv (const int dir, const StdRegions::VarCoeffMap &varcoeffs)
 
Array< OneD, NekDoubleGetMFMag (const int dir, const StdRegions::VarCoeffMap &varcoeffs)
 
void v_MultiplyByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_DivideByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_ComputeLaplacianMetric ()
 
int v_GetCoordim () const override
 
void v_GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
 
virtual DNekScalMatSharedPtr v_GetLocMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual void v_DropLocMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual DNekMatSharedPtr v_BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
 
virtual DNekMatSharedPtr v_BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd)
 
virtual void v_ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
 
virtual void v_AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddFaceNormBoundaryInt (const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
virtual void v_DGDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &coeffs, Array< OneD, NekDouble > &outarray)
 
virtual NekDouble v_VectorFlux (const Array< OneD, Array< OneD, NekDouble > > &vec)
 
virtual void v_NormalTraceDerivFactors (Array< OneD, Array< OneD, NekDouble > > &factors, Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors)
 
virtual void v_AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual StdRegions::Orientation v_GetTraceOrient (int trace)
 
void v_SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_GetTraceQFactors (const int trace, Array< OneD, NekDouble > &outarray)
 
virtual void v_GetTracePhysVals (const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient)
 
virtual void v_GetTracePhysMap (const int edge, Array< OneD, int > &outarray)
 
virtual void v_ReOrientTracePhysMap (const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1=-1)
 
virtual void v_ComputeTraceNormal (const int id)
 
virtual const Array< OneD, const NekDouble > & v_GetPhysNormals ()
 
virtual void v_SetPhysNormals (Array< OneD, const NekDouble > &normal)
 
virtual void v_SetUpPhysNormals (const int id)
 
virtual void v_AddRobinMassMatrix (const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
 
virtual void v_AddRobinTraceContribution (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
 
virtual void v_TraceNormLen (const int traceid, NekDouble &h, NekDouble &p)
 
virtual void v_GenTraceExp (const int traceid, ExpansionSharedPtr &exp)
 

Private Attributes

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

Additional Inherited Members

- Protected Attributes inherited from Nektar::StdRegions::StdNodalTriExp
LibUtilities::PointsKey m_nodalPointsKey
 
- 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::Expansion2D
std::vector< bool > m_requireNeg
 
- 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

Definition at line 48 of file NodalTriExp.h.

Constructor & Destructor Documentation

◆ NodalTriExp() [1/2]

Nektar::LocalRegions::NodalTriExp::NodalTriExp ( const LibUtilities::BasisKey Ba,
const LibUtilities::BasisKey Bb,
const LibUtilities::PointsType  Ntype,
const SpatialDomains::TriGeomSharedPtr geom 
)

Constructor using BasisKey class for quadrature points and order definition.

Definition at line 42 of file NodalTriExp.cpp.

47 Ba.GetNumModes(), (Bb.GetNumModes())),
48 2, Ba, Bb),
50 Ba.GetNumModes(), (Bb.GetNumModes())),
51 Ba, Bb),
52 StdNodalTriExp(Ba, Bb, Ntype), Expansion(geom), Expansion2D(geom),
54 std::bind(&Expansion2D::CreateMatrix, this, std::placeholders::_1),
55 std::string("NodalTriExpMatrix")),
57 this, std::placeholders::_1),
58 std::string("NodalTriExpStaticCondMatrix"))
59{
60}
DNekScalMatSharedPtr CreateMatrix(const MatrixKey &mkey)
Definition: Expansion2D.cpp:55
Expansion2D(SpatialDomains::Geometry2DSharedPtr pGeom)
Definition: Expansion2D.cpp:50
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
Definition: Expansion.cpp:272
Expansion(SpatialDomains::GeometrySharedPtr pGeom)
Definition: Expansion.cpp:43
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager
Definition: NodalTriExp.h:177
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLess > m_staticCondMatrixManager
Definition: NodalTriExp.h:179
StdExpansion()
Default Constructor.
int getNumberOfCoefficients(int Na, int Nb)
Definition: ShapeType.hpp:109

◆ NodalTriExp() [2/2]

Nektar::LocalRegions::NodalTriExp::NodalTriExp ( const NodalTriExp T)

Definition at line 62 of file NodalTriExp.cpp.

63 : StdExpansion(T), StdExpansion2D(T), StdRegions::StdTriExp(T),
64 StdRegions::StdNodalTriExp(T), Expansion(T), Expansion2D(T),
65 m_matrixManager(T.m_matrixManager),
66 m_staticCondMatrixManager(T.m_staticCondMatrixManager)
67{
68}

◆ ~NodalTriExp()

Nektar::LocalRegions::NodalTriExp::~NodalTriExp ( )
overridedefault

Member Function Documentation

◆ v_AlignVectorToCollapsedDir()

void Nektar::LocalRegions::NodalTriExp::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 269 of file NodalTriExp.cpp.

272{
273 ASSERTL1((dir == 0) || (dir == 1) || (dir == 2), "Invalid direction.");
274 ASSERTL1((dir == 2) ? (m_geom->GetCoordim() == 3) : true,
275 "Invalid direction.");
276
277 int nquad0 = m_base[0]->GetNumPoints();
278 int nquad1 = m_base[1]->GetNumPoints();
279 int nqtot = nquad0 * nquad1;
280 int wspsize = max(nqtot, m_ncoeffs);
281
282 const Array<TwoD, const NekDouble> &df =
283 m_metricinfo->GetDerivFactors(GetPointsKeys());
284
285 Array<OneD, NekDouble> tmp0(4 * wspsize);
286 Array<OneD, NekDouble> tmp3(tmp0 + wspsize);
287 Array<OneD, NekDouble> gfac0(tmp0 + 2 * wspsize);
288 Array<OneD, NekDouble> gfac1(tmp0 + 3 * wspsize);
289
290 Array<OneD, NekDouble> tmp1 = outarray[0];
291 Array<OneD, NekDouble> tmp2 = outarray[1];
292
293 const Array<OneD, const NekDouble> &z0 = m_base[0]->GetZ();
294 const Array<OneD, const NekDouble> &z1 = m_base[1]->GetZ();
295
296 // set up geometric factor: 2/(1-z1)
297 for (int i = 0; i < nquad1; ++i)
298 {
299 gfac0[i] = 2.0 / (1 - z1[i]);
300 }
301 for (int i = 0; i < nquad0; ++i)
302 {
303 gfac1[i] = 0.5 * (1 + z0[i]);
304 }
305
306 for (int i = 0; i < nquad1; ++i)
307 {
308 Vmath::Smul(nquad0, gfac0[i], &inarray[0] + i * nquad0, 1,
309 &tmp0[0] + i * nquad0, 1);
310 }
311
312 for (int i = 0; i < nquad1; ++i)
313 {
314 Vmath::Vmul(nquad0, &gfac1[0], 1, &tmp0[0] + i * nquad0, 1,
315 &tmp1[0] + i * nquad0, 1);
316 }
317
318 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
319 {
320 Vmath::Vmul(nqtot, &df[2 * dir][0], 1, &tmp0[0], 1, &tmp0[0], 1);
321 Vmath::Vmul(nqtot, &df[2 * dir + 1][0], 1, &tmp1[0], 1, &tmp1[0], 1);
322 Vmath::Vmul(nqtot, &df[2 * dir + 1][0], 1, &inarray[0], 1, &tmp2[0], 1);
323 }
324 else
325 {
326 Vmath::Smul(nqtot, df[2 * dir][0], tmp0, 1, tmp0, 1);
327 Vmath::Smul(nqtot, df[2 * dir + 1][0], tmp1, 1, tmp1, 1);
328 Vmath::Smul(nqtot, df[2 * dir + 1][0], inarray, 1, tmp2, 1);
329 }
330 Vmath::Vadd(nqtot, tmp0, 1, tmp1, 1, tmp1, 1);
331}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242
SpatialDomains::GeometrySharedPtr m_geom
Definition: Expansion.h:273
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
Definition: Expansion.h:274
const LibUtilities::PointsKeyVector GetPointsKeys() const
Array< OneD, LibUtilities::BasisSharedPtr > m_base
@ eDeformed
Geometry is curved or has non-constant factors.
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.hpp:72
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.hpp:180
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100

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

◆ v_ComputeTraceNormal()

void Nektar::LocalRegions::NodalTriExp::v_ComputeTraceNormal ( const int  edge)
overrideprotectedvirtual

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 389 of file NodalTriExp.cpp.

390{
391 int i;
392 const SpatialDomains::GeomFactorsSharedPtr &geomFactors =
393 GetGeom()->GetMetricInfo();
394 const SpatialDomains::GeomType type = geomFactors->GetGtype();
395
397 const Array<TwoD, const NekDouble> &df =
398 geomFactors->GetDerivFactors(ptsKeys);
399 const Array<OneD, const NekDouble> &jac = geomFactors->GetJac(ptsKeys);
400 int nqe = m_base[0]->GetNumPoints();
401 int dim = GetCoordim();
402
403 m_traceNormals[edge] = Array<OneD, Array<OneD, NekDouble>>(dim);
404 Array<OneD, Array<OneD, NekDouble>> &normal = m_traceNormals[edge];
405 for (i = 0; i < dim; ++i)
406 {
407 normal[i] = Array<OneD, NekDouble>(nqe);
408 }
409
410 size_t nqb = nqe;
411 size_t nbnd = edge;
412 m_elmtBndNormDirElmtLen[nbnd] = Array<OneD, NekDouble>{nqb, 0.0};
413 Array<OneD, NekDouble> &length = m_elmtBndNormDirElmtLen[nbnd];
414
415 // Regular geometry case
416 if ((type == SpatialDomains::eRegular) ||
418 {
419 NekDouble fac;
420 // Set up normals
421 switch (edge)
422 {
423 case 0:
424 for (i = 0; i < GetCoordim(); ++i)
425 {
426 Vmath::Fill(nqe, -df[2 * i + 1][0], normal[i], 1);
427 }
428 break;
429 case 1:
430 for (i = 0; i < GetCoordim(); ++i)
431 {
432 Vmath::Fill(nqe, df[2 * i + 1][0] + df[2 * i][0], normal[i],
433 1);
434 }
435 break;
436 case 2:
437 for (i = 0; i < GetCoordim(); ++i)
438 {
439 Vmath::Fill(nqe, -df[2 * i][0], normal[i], 1);
440 }
441 break;
442 default:
443 ASSERTL0(false, "Edge is out of range (edge < 3)");
444 }
445
446 // normalise
447 fac = 0.0;
448 for (i = 0; i < GetCoordim(); ++i)
449 {
450 fac += normal[i][0] * normal[i][0];
451 }
452 fac = 1.0 / sqrt(fac);
453
454 Vmath::Fill(nqb, fac, length, 1);
455
456 for (i = 0; i < GetCoordim(); ++i)
457 {
458 Vmath::Smul(nqe, fac, normal[i], 1, normal[i], 1);
459 }
460 }
461 else // Set up deformed normals
462 {
463 int j;
464
465 int nquad0 = ptsKeys[0].GetNumPoints();
466 int nquad1 = ptsKeys[1].GetNumPoints();
467
468 LibUtilities::PointsKey from_key;
469
470 Array<OneD, NekDouble> normals(GetCoordim() * max(nquad0, nquad1), 0.0);
471 Array<OneD, NekDouble> edgejac(GetCoordim() * max(nquad0, nquad1), 0.0);
472
473 // Extract Jacobian along edges and recover local
474 // derivates (dx/dr) for polynomial interpolation by
475 // multiplying m_gmat by jacobian
476 switch (edge)
477 {
478 case 0:
479 for (j = 0; j < nquad0; ++j)
480 {
481 edgejac[j] = jac[j];
482 for (i = 0; i < GetCoordim(); ++i)
483 {
484 normals[i * nquad0 + j] =
485 -df[2 * i + 1][j] * edgejac[j];
486 }
487 }
488 from_key = ptsKeys[0];
489 break;
490 case 1:
491 for (j = 0; j < nquad1; ++j)
492 {
493 edgejac[j] = jac[nquad0 * j + nquad0 - 1];
494 for (i = 0; i < GetCoordim(); ++i)
495 {
496 normals[i * nquad1 + j] =
497 (df[2 * i][nquad0 * j + nquad0 - 1] +
498 df[2 * i + 1][nquad0 * j + nquad0 - 1]) *
499 edgejac[j];
500 }
501 }
502 from_key = ptsKeys[1];
503 break;
504 case 2:
505 for (j = 0; j < nquad1; ++j)
506 {
507 edgejac[j] = jac[nquad0 * j];
508 for (i = 0; i < GetCoordim(); ++i)
509 {
510 normals[i * nquad1 + j] =
511 -df[2 * i][nquad0 * j] * edgejac[j];
512 }
513 }
514 from_key = ptsKeys[1];
515 break;
516 default:
517 ASSERTL0(false, "edge is out of range (edge < 3)");
518 }
519
520 int nq = from_key.GetNumPoints();
521 Array<OneD, NekDouble> work(nqe, 0.0);
522
523 // interpolate Jacobian and invert
524 LibUtilities::Interp1D(from_key, jac, m_base[0]->GetPointsKey(), work);
525 Vmath::Sdiv(nq, 1.0, &work[0], 1, &work[0], 1);
526
527 // interpolate
528 for (i = 0; i < GetCoordim(); ++i)
529 {
530 LibUtilities::Interp1D(from_key, &normals[i * nq],
531 m_base[0]->GetPointsKey(), &normal[i][0]);
532 Vmath::Vmul(nqe, work, 1, normal[i], 1, normal[i], 1);
533 }
534
535 // normalise normal vectors
536 Vmath::Zero(nqe, work, 1);
537 for (i = 0; i < GetCoordim(); ++i)
538 {
539 Vmath::Vvtvp(nqe, normal[i], 1, normal[i], 1, work, 1, work, 1);
540 }
541
542 Vmath::Vsqrt(nqe, work, 1, work, 1);
543 Vmath::Sdiv(nqe, 1.0, work, 1, work, 1);
544
545 Vmath::Vcopy(nqb, work, 1, length, 1);
546
547 for (i = 0; i < GetCoordim(); ++i)
548 {
549 Vmath::Vmul(nqe, normal[i], 1, work, 1, normal[i], 1);
550 }
551
552 // Reverse direction so that points are in
553 // anticlockwise direction if edge >=2
554 if (edge >= 2)
555 {
556 for (i = 0; i < GetCoordim(); ++i)
557 {
558 Vmath::Reverse(nqe, normal[i], 1, normal[i], 1);
559 }
560 }
561 }
562}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
std::map< int, NormalVector > m_traceNormals
Definition: Expansion.h:276
std::map< int, Array< OneD, NekDouble > > m_elmtBndNormDirElmtLen
the element length in each element boundary(Vertex, edge or face) normal direction calculated based o...
Definition: Expansion.h:286
SpatialDomains::GeometrySharedPtr GetGeom() const
Definition: Expansion.cpp:167
void Interp1D(const BasisKey &fbasis0, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, Array< OneD, NekDouble > &to)
this function interpolates a 1D function evaluated at the quadrature points of the basis fbasis0 to ...
Definition: Interp.cpp:47
std::vector< PointsKey > PointsKeyVector
Definition: Points.h:231
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
Definition: GeomFactors.h:60
GeomType
Indicates the type of element geometry.
@ eRegular
Geometry is straight-sided with constant geometric factors.
@ eMovingRegular
Currently unused.
double NekDouble
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.hpp:340
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.hpp:366
void Sdiv(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha/x.
Definition: Vmath.hpp:154
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.hpp:54
void Reverse(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:844
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References ASSERTL0, Nektar::SpatialDomains::eMovingRegular, Nektar::SpatialDomains::eRegular, Vmath::Fill(), Nektar::StdRegions::StdExpansion::GetCoordim(), Nektar::LocalRegions::Expansion::GetGeom(), Nektar::LibUtilities::PointsKey::GetNumPoints(), Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::LibUtilities::Interp1D(), Nektar::StdRegions::StdExpansion::m_base, Nektar::LocalRegions::Expansion::m_elmtBndNormDirElmtLen, Nektar::LocalRegions::Expansion::m_traceNormals, Vmath::Reverse(), Vmath::Sdiv(), Vmath::Smul(), tinysimd::sqrt(), Vmath::Vcopy(), Vmath::Vmul(), Vmath::Vsqrt(), Vmath::Vvtvp(), and Vmath::Zero().

◆ v_CreateStdMatrix()

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

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 657 of file NodalTriExp.cpp.

659{
660 LibUtilities::BasisKey bkey0 = m_base[0]->GetBasisKey();
661 LibUtilities::BasisKey bkey1 = m_base[1]->GetBasisKey();
665
666 return tmp->GetStdMatrix(mkey);
667}
PointsType GetPointsType() const
Definition: Points.h:90
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
LibUtilities::PointsKey m_nodalPointsKey
std::shared_ptr< StdNodalTriExp > StdNodalTriExpSharedPtr

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::PointsKey::GetPointsType(), Nektar::StdRegions::StdExpansion::m_base, and Nektar::StdRegions::StdNodalTriExp::m_nodalPointsKey.

◆ v_DropLocMatrix()

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

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 678 of file NodalTriExp.cpp.

679{
680 m_matrixManager.DeleteObject(mkey);
681}

References m_matrixManager.

◆ v_ExtractDataToCoeffs()

void Nektar::LocalRegions::NodalTriExp::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 564 of file NodalTriExp.cpp.

568{
569 Array<OneD, NekDouble> modes(m_ncoeffs);
570 Expansion::ExtractDataToCoeffs(data, nummodes, mode_offset, &modes[0],
571 fromType);
572
573 Array<OneD, NekDouble> nodes(m_ncoeffs, coeffs, eArrayWrapper);
574 ModalToNodal(modes, nodes);
575}
void ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
Definition: Expansion.cpp:106
void ModalToNodal(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)

References Nektar::eArrayWrapper, Nektar::LocalRegions::Expansion::ExtractDataToCoeffs(), Nektar::StdRegions::StdExpansion::m_ncoeffs, Nektar::StdRegions::StdNodalTriExp::ModalToNodal(), and CG_Iterations::modes.

◆ v_FwdTrans()

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

Transform a given function from physical quadrature space to coefficient space.

See also
StdExpansion::FwdTrans

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 176 of file NodalTriExp.cpp.

178{
179 IProductWRTBase(inarray, outarray);
180
181 // get Mass matrix inverse
182 MatrixKey masskey(StdRegions::eInvMass, DetShapeType(), *this,
186 DNekScalMatSharedPtr matsys = m_matrixManager[masskey];
187
188 // copy inarray in case inarray == outarray
189 NekVector<NekDouble> in(m_ncoeffs, outarray, eCopy);
190 NekVector<NekDouble> out(m_ncoeffs, outarray, eWrapper);
191
192 out = (*matsys) * in;
193}
void IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
this function calculates the inner product of a given function f with the different modes of the expa...
Definition: StdExpansion.h:537
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
Definition: StdExpansion.h:370
static ConstFactorMap NullConstFactorMap
Definition: StdRegions.hpp:431
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:376
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr

References Nektar::StdRegions::StdExpansion::DetShapeType(), Nektar::eCopy, Nektar::StdRegions::eInvMass, Nektar::eWrapper, Nektar::LibUtilities::PointsKey::GetPointsType(), Nektar::StdRegions::StdExpansion::IProductWRTBase(), m_matrixManager, Nektar::StdRegions::StdExpansion::m_ncoeffs, Nektar::StdRegions::StdNodalTriExp::m_nodalPointsKey, Nektar::StdRegions::NullConstFactorMap, and Nektar::StdRegions::NullVarCoeffMap.

◆ v_GenMatrix()

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

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 636 of file NodalTriExp.cpp.

637{
638 DNekMatSharedPtr returnval;
639
640 switch (mkey.GetMatrixType())
641 {
648 returnval = Expansion2D::v_GenMatrix(mkey);
649 break;
650 default:
651 returnval = StdNodalTriExp::v_GenMatrix(mkey);
652 break;
653 }
654 return returnval;
655}
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::StdMatrixKey::GetMatrixType(), and Nektar::LocalRegions::Expansion2D::v_GenMatrix().

◆ v_GetCoord()

void Nektar::LocalRegions::NodalTriExp::v_GetCoord ( const Array< OneD, const NekDouble > &  lcoord,
Array< OneD, NekDouble > &  coord 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 359 of file NodalTriExp.cpp.

361{
362 int i;
363
364 ASSERTL1(Lcoords[0] >= -1.0 && Lcoords[1] <= 1.0 && Lcoords[1] >= -1.0 &&
365 Lcoords[1] <= 1.0,
366 "Local coordinates are not in region [-1,1]");
367
368 m_geom->FillGeom();
369
370 for (i = 0; i < m_geom->GetCoordim(); ++i)
371 {
372 coords[i] = m_geom->GetCoord(i, Lcoords);
373 }
374}

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

◆ v_GetCoords()

void Nektar::LocalRegions::NodalTriExp::v_GetCoords ( Array< OneD, NekDouble > &  coords_0,
Array< OneD, NekDouble > &  coords_1 = NullNekDouble1DArray,
Array< OneD, NekDouble > &  coords_2 = NullNekDouble1DArray 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdTriExp.

Definition at line 352 of file NodalTriExp.cpp.

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

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

◆ v_GetLinStdExp()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 341 of file NodalTriExp.cpp.

342{
343 LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(), 2,
344 m_base[0]->GetPointsKey());
345 LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(), 2,
346 m_base[1]->GetPointsKey());
347
349 bkey0, bkey1, m_nodalPointsKey.GetPointsType());
350}
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
Definition: StdExpansion.h:156

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::StdRegions::StdExpansion::GetBasisType(), Nektar::LibUtilities::PointsKey::GetPointsType(), Nektar::StdRegions::StdExpansion::m_base, and Nektar::StdRegions::StdNodalTriExp::m_nodalPointsKey.

◆ v_GetLocMatrix()

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

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 669 of file NodalTriExp.cpp.

670{
671 return m_matrixManager[mkey];
672}

References m_matrixManager.

◆ v_GetLocStaticCondMatrix()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 673 of file NodalTriExp.cpp.

675{
676 return m_staticCondMatrixManager[mkey];
677}

References m_staticCondMatrixManager.

◆ v_GetStdExp()

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

◆ v_GetTracePhysVals()

void Nektar::LocalRegions::NodalTriExp::v_GetTracePhysVals ( const int  edge,
const StdRegions::StdExpansionSharedPtr EdgeExp,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
StdRegions::Orientation  orient 
)
overrideprotectedvirtual

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 577 of file NodalTriExp.cpp.

581{
582 int nquad0 = m_base[0]->GetNumPoints();
583 int nquad1 = m_base[1]->GetNumPoints();
584 int nt = 0;
585 // Extract in Cartesian direction because we have to deal with
586 // e.g. Gauss-Radau points.
587 switch (edge)
588 {
589 case 0:
590 Vmath::Vcopy(nquad0, &(inarray[0]), 1, &(outarray[0]), 1);
591 nt = nquad0;
592 break;
593 case 1:
594 Vmath::Vcopy(nquad1, &(inarray[0]) + (nquad0 - 1), nquad0,
595 &(outarray[0]), 1);
596 nt = nquad1;
597 break;
598 case 2:
599 Vmath::Vcopy(nquad1, &(inarray[0]), nquad0, &(outarray[0]), 1);
600 nt = nquad1;
601 break;
602 default:
603 ASSERTL0(false, "edge value (< 3) is out of range");
604 break;
605 }
606
607 ASSERTL1(EdgeExp->GetBasis(0)->GetPointsType() ==
609 "Edge expansion should be GLL");
610
611 // Interpolate if required
612 if (m_base[edge ? 1 : 0]->GetPointsKey() !=
613 EdgeExp->GetBasis(0)->GetPointsKey())
614 {
615 Array<OneD, NekDouble> outtmp(max(nquad0, nquad1));
616
617 Vmath::Vcopy(nt, outarray, 1, outtmp, 1);
618
619 LibUtilities::Interp1D(m_base[edge ? 1 : 0]->GetPointsKey(), outtmp,
620 EdgeExp->GetBasis(0)->GetPointsKey(), outarray);
621 }
622
623 if (orient == StdRegions::eNoOrientation)
624 {
625 orient = GetTraceOrient(edge);
626 }
627
628 // Reverse data if necessary
629 if (orient == StdRegions::eBackwards)
630 {
631 Vmath::Reverse(EdgeExp->GetNumPoints(0), &outarray[0], 1, &outarray[0],
632 1);
633 }
634}
StdRegions::Orientation GetTraceOrient(int trace)
Definition: Expansion.h:168
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:51

References ASSERTL0, ASSERTL1, Nektar::StdRegions::eBackwards, Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::StdRegions::eNoOrientation, Nektar::LocalRegions::Expansion::GetTraceOrient(), Nektar::LibUtilities::Interp1D(), Nektar::StdRegions::StdExpansion::m_base, Vmath::Reverse(), and Vmath::Vcopy().

◆ v_HelmholtzMatrixOp()

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

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 712 of file NodalTriExp.cpp.

715{
716 StdExpansion::HelmholtzMatrixOp_MatFree_GenericImpl(inarray, outarray,
717 mkey);
718}

◆ v_Integral()

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

Integrates the specified function over the domain.

See also
StdRegions::StdExpansion::Integral.

Reimplemented from Nektar::StdRegions::StdTriExp.

Definition at line 70 of file NodalTriExp.cpp.

71{
72 int nquad0 = m_base[0]->GetNumPoints();
73 int nquad1 = m_base[1]->GetNumPoints();
74 Array<OneD, const NekDouble> jac = m_metricinfo->GetJac(GetPointsKeys());
75 NekDouble ival;
76 Array<OneD, NekDouble> tmp(nquad0 * nquad1);
77
78 // multiply inarray with Jacobian
79 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
80 {
81 Vmath::Vmul(nquad0 * nquad1, jac, 1, inarray, 1, tmp, 1);
82 }
83 else
84 {
85 Vmath::Smul(nquad0 * nquad1, jac[0], inarray, 1, tmp, 1);
86 }
87
88 // call StdQuadExp version;
89 ival = StdNodalTriExp::v_Integral(tmp);
90 return ival;
91}

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::NodalTriExp::v_IProductWRTBase ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
overrideprotectedvirtual

Calculates the inner product of a given function f with the different modes of the expansion.

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 195 of file NodalTriExp.cpp.

197{
198 v_IProductWRTBase_SumFac(inarray, outarray);
199}
void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override

References v_IProductWRTBase_SumFac().

◆ v_IProductWRTBase_SumFac()

void Nektar::LocalRegions::NodalTriExp::v_IProductWRTBase_SumFac ( const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
bool  multiplybyweights = true 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 208 of file NodalTriExp.cpp.

211{
212 int nquad0 = m_base[0]->GetNumPoints();
213 int nquad1 = m_base[1]->GetNumPoints();
214 int order1 = m_base[1]->GetNumModes();
215
216 if (multiplybyweights)
217 {
218 Array<OneD, NekDouble> tmp(nquad0 * nquad1 + nquad0 * order1);
219 Array<OneD, NekDouble> wsp(tmp + nquad0 * nquad1);
220
221 MultiplyByQuadratureMetric(inarray, tmp);
222 StdTriExp::IProductWRTBase_SumFacKernel(
223 m_base[0]->GetBdata(), m_base[1]->GetBdata(), tmp, outarray, wsp);
224 NodalToModalTranspose(outarray, outarray);
225 }
226 else
227 {
228 Array<OneD, NekDouble> wsp(nquad0 * order1);
229
230 StdTriExp::IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(),
231 m_base[1]->GetBdata(), inarray,
232 outarray, wsp);
233 NodalToModalTranspose(outarray, outarray);
234 }
235}
void MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdExpansion.h:732
void NodalToModalTranspose(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)

References Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::MultiplyByQuadratureMetric(), and Nektar::StdRegions::StdNodalTriExp::NodalToModalTranspose().

Referenced by v_IProductWRTBase().

◆ v_IProductWRTDerivBase()

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

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 201 of file NodalTriExp.cpp.

204{
205 v_IProductWRTDerivBase_SumFac(dir, inarray, outarray);
206}
void v_IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override

References v_IProductWRTDerivBase_SumFac().

◆ v_IProductWRTDerivBase_SumFac()

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

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 237 of file NodalTriExp.cpp.

240{
241 int nquad0 = m_base[0]->GetNumPoints();
242 int nquad1 = m_base[1]->GetNumPoints();
243 int nqtot = nquad0 * nquad1;
244 int wspsize = max(nqtot, m_ncoeffs);
245
246 Array<OneD, NekDouble> tmp0(4 * wspsize);
247 Array<OneD, NekDouble> tmp1(tmp0 + wspsize);
248 Array<OneD, NekDouble> tmp2(tmp0 + 2 * wspsize);
249 Array<OneD, NekDouble> tmp3(tmp0 + 3 * wspsize);
250
251 Array<OneD, Array<OneD, NekDouble>> tmp2D{2};
252 tmp2D[0] = tmp1;
253 tmp2D[1] = tmp2;
254
255 AlignVectorToCollapsedDir(dir, inarray, tmp2D);
256
257 MultiplyByQuadratureMetric(tmp1, tmp1);
258 MultiplyByQuadratureMetric(tmp2, tmp2);
259
260 IProductWRTBase_SumFacKernel(m_base[0]->GetDbdata(), m_base[1]->GetBdata(),
261 tmp1, tmp3, tmp0);
262 IProductWRTBase_SumFacKernel(m_base[0]->GetBdata(), m_base[1]->GetDbdata(),
263 tmp2, outarray, tmp0);
264 Vmath::Vadd(m_ncoeffs, tmp3, 1, outarray, 1, outarray, 1);
265
266 NodalToModalTranspose(outarray, outarray);
267}
void AlignVectorToCollapsedDir(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: Expansion.h:151
void IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0=true, bool doCheckCollDir1=true)

References Nektar::LocalRegions::Expansion::AlignVectorToCollapsedDir(), Nektar::StdRegions::StdExpansion2D::IProductWRTBase_SumFacKernel(), Nektar::StdRegions::StdExpansion::m_base, Nektar::StdRegions::StdExpansion::m_ncoeffs, Nektar::StdRegions::StdExpansion::MultiplyByQuadratureMetric(), Nektar::StdRegions::StdNodalTriExp::NodalToModalTranspose(), and Vmath::Vadd().

Referenced by v_IProductWRTDerivBase().

◆ v_LaplacianMatrixOp() [1/2]

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

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 690 of file NodalTriExp.cpp.

693{
694 StdExpansion::LaplacianMatrixOp_MatFree_GenericImpl(inarray, outarray,
695 mkey);
696}

◆ v_LaplacianMatrixOp() [2/2]

void Nektar::LocalRegions::NodalTriExp::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::StdNodalTriExp.

Definition at line 698 of file NodalTriExp.cpp.

701{
702 StdExpansion::LaplacianMatrixOp_MatFree(k1, k2, inarray, outarray, mkey);
703}

◆ v_MassMatrixOp()

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

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 683 of file NodalTriExp.cpp.

686{
687 StdExpansion::MassMatrixOp_MatFree(inarray, outarray, mkey);
688}

◆ v_PhysDeriv() [1/2]

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

Calculate the derivative of the physical points.

\( \frac{\partial u}{\partial x_1} = \left . \frac{2.0}{1-\eta_2} \frac{\partial u}{\partial d\eta_1} \right |_{\eta_2}\)

\( \frac{\partial u}{\partial x_2} = \left . \frac{1+\eta_1}{1-\eta_2} \frac{\partial u}{\partial d\eta_1} \right |_{\eta_2} + \left . \frac{\partial u}{\partial d\eta_2} \right |_{\eta_1} \)

Reimplemented from Nektar::StdRegions::StdTriExp.

Definition at line 93 of file NodalTriExp.cpp.

97{
98 int nquad0 = m_base[0]->GetNumPoints();
99 int nquad1 = m_base[1]->GetNumPoints();
100 int nqtot = nquad0 * nquad1;
101 const Array<TwoD, const NekDouble> &df =
102 m_metricinfo->GetDerivFactors(GetPointsKeys());
103
104 Array<OneD, NekDouble> diff0(2 * nqtot);
105 Array<OneD, NekDouble> diff1(diff0 + nqtot);
106
107 StdNodalTriExp::v_PhysDeriv(inarray, diff0, diff1);
108
109 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
110 {
111 if (out_d0.size())
112 {
113 Vmath::Vmul(nqtot, df[0], 1, diff0, 1, out_d0, 1);
114 Vmath::Vvtvp(nqtot, df[1], 1, diff1, 1, out_d0, 1, out_d0, 1);
115 }
116
117 if (out_d1.size())
118 {
119 Vmath::Vmul(nqtot, df[2], 1, diff0, 1, out_d1, 1);
120 Vmath::Vvtvp(nqtot, df[3], 1, diff1, 1, out_d1, 1, out_d1, 1);
121 }
122
123 if (out_d2.size())
124 {
125 Vmath::Vmul(nqtot, df[4], 1, diff0, 1, out_d2, 1);
126 Vmath::Vvtvp(nqtot, df[5], 1, diff1, 1, out_d2, 1, out_d2, 1);
127 }
128 }
129 else // regular geometry
130 {
131 if (out_d0.size())
132 {
133 Vmath::Smul(nqtot, df[0][0], diff0, 1, out_d0, 1);
134 Blas::Daxpy(nqtot, df[1][0], diff1, 1, out_d0, 1);
135 }
136
137 if (out_d1.size())
138 {
139 Vmath::Smul(nqtot, df[2][0], diff0, 1, out_d1, 1);
140 Blas::Daxpy(nqtot, df[3][0], diff1, 1, out_d1, 1);
141 }
142
143 if (out_d2.size())
144 {
145 Vmath::Smul(nqtot, df[4][0], diff0, 1, out_d2, 1);
146 Blas::Daxpy(nqtot, df[5][0], diff1, 1, out_d2, 1);
147 }
148 }
149}
static void Daxpy(const int &n, const double &alpha, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: y = alpha x plus y.
Definition: Blas.hpp:135

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

◆ v_PhysDeriv() [2/2]

void Nektar::LocalRegions::NodalTriExp::v_PhysDeriv ( const int  dir,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  out_d0 
)
overrideprotectedvirtual

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

See also
StdRegions::StdExpansion::PhysDeriv

Reimplemented from Nektar::StdRegions::StdTriExp.

Definition at line 151 of file NodalTriExp.cpp.

154{
155 Array<OneD, NekDouble> tmp;
156 switch (dir)
157 {
158 case 0:
159 {
160 PhysDeriv(inarray, outarray, tmp);
161 }
162 break;
163 case 1:
164 {
165 PhysDeriv(inarray, tmp, outarray);
166 }
167 break;
168 default:
169 {
170 ASSERTL1(dir >= 0 && dir < 2, "input dir is out of range");
171 }
172 break;
173 }
174}
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)
Definition: StdExpansion.h:858

References ASSERTL1, and Nektar::StdRegions::StdExpansion::PhysDeriv().

◆ v_PhysEvaluate()

NekDouble Nektar::LocalRegions::NodalTriExp::v_PhysEvaluate ( const Array< OneD, const NekDouble > &  coords,
const Array< OneD, const NekDouble > &  physvals 
)
overrideprotectedvirtual

This function evaluates the expansion at a single (arbitrary) point of the domain.

This function is a wrapper around the virtual function v_PhysEvaluate()

Based on the value of the expansion at the quadrature points, this function calculates the value of the expansion at an arbitrary single points (with coordinates \( \mathbf{x_c}\) given by the pointer coords). This operation, equivalent to

\[ u(\mathbf{x_c}) = \sum_p \phi_p(\mathbf{x_c}) \hat{u}_p \]

is evaluated using Lagrangian interpolants through the quadrature points:

\[ u(\mathbf{x_c}) = \sum_p h_p(\mathbf{x_c}) u_p\]

This function requires that the physical value array \(\mathbf{u}\) (implemented as the attribute #m_phys) is set.

Parameters
coordsthe coordinates of the single point
Returns
returns the value of the expansion at the single point

Reimplemented from Nektar::StdRegions::StdExpansion2D.

Definition at line 376 of file NodalTriExp.cpp.

380{
381 Array<OneD, NekDouble> Lcoord = Array<OneD, NekDouble>(2);
382
383 ASSERTL0(m_geom, "m_geom not defined");
384 m_geom->GetLocCoords(coord, Lcoord);
385
386 return StdExpansion2D::v_PhysEvaluate(Lcoord, physvals);
387}

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

◆ v_WeakDerivMatrixOp()

void Nektar::LocalRegions::NodalTriExp::v_WeakDerivMatrixOp ( const int  i,
const Array< OneD, const NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray,
const StdRegions::StdMatrixKey mkey 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdNodalTriExp.

Definition at line 705 of file NodalTriExp.cpp.

708{
709 StdExpansion::WeakDerivMatrixOp_MatFree(i, inarray, outarray, mkey);
710}

Member Data Documentation

◆ m_matrixManager

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

Definition at line 177 of file NodalTriExp.h.

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

◆ m_staticCondMatrixManager

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

Definition at line 179 of file NodalTriExp.h.

Referenced by v_GetLocStaticCondMatrix().