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

#include <PrismExp.h>

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

Public Member Functions

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

Protected Member Functions

NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray) override
 Integrate the physical point list inarray over prismatic region and return the value. More...
 
void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2) override
 Calculate the derivative of the physical points. More...
 
void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Forward transform from physical quadrature space stored in inarray and evaluate the expansion coefficients and store in (this)->m_coeffs. More...
 
void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculate the inner product of inarray with respect to the basis B=base0*base1*base2 and put into outarray: More...
 
void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculates the inner product \( I_{pqr} = (u, \partial_{x_i} \phi_{pqr}) \). More...
 
void v_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
 
void v_GetCoord (const Array< OneD, const NekDouble > &Lcoords, Array< OneD, NekDouble > &coords) override
 Get the coordinates #coords at the local coordinates #Lcoords. More...
 
void v_GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
 
NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals) 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...
 
NekDouble v_PhysEvalFirstDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) override
 
StdRegions::StdExpansionSharedPtr v_GetStdExp (void) const override
 
StdRegions::StdExpansionSharedPtr v_GetLinStdExp (void) const override
 
void v_ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int mode_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType) override
 
void v_GetTracePhysMap (const int face, Array< OneD, int > &outarray) override
 
void v_ComputeTraceNormal (const int face) override
 Get the normals along specficied face Get the face normals interplated to a points0 x points 0 type distribution. More...
 
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_HelmholtzMatrixOp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdRegions::StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_CreateStdMatrix (const StdRegions::StdMatrixKey &mkey) override
 
DNekScalMatSharedPtr v_GetLocMatrix (const MatrixKey &mkey) override
 
DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix (const MatrixKey &mkey) override
 
void v_DropLocMatrix (const MatrixKey &mkey) override
 
void v_DropLocStaticCondMatrix (const MatrixKey &mkey) override
 
void v_GetSimplexEquiSpacedConnectivity (Array< OneD, int > &conn, bool standard=true) override
 
void v_NormalTraceDerivFactors (Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors, Array< OneD, Array< OneD, NekDouble > > &d2factors) override
 : This method gets all of the factors which are required as part of the Gradient Jump Penalty stabilisation and involves the product of the normal and geometric factors along the element trace. More...
 
- Protected Member Functions inherited from Nektar::StdRegions::StdPrismExp
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
 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) override
 
void v_StdPhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2) override
 
void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Forward transform from physical quadrature space stored in inarray and evaluate the expansion coefficients and store in outarray. More...
 
void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Calculate the inner product of inarray with respect to the basis B=base0*base1*base2 and put into outarray: More...
 
void v_IProductWRTBase_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
 
void v_IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2) override
 
void v_IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 Inner product of inarray over region with respect to the object's default expansion basis; output in outarray. More...
 
void v_IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_LocCoordToLocCollapsed (const Array< OneD, const NekDouble > &xi, Array< OneD, NekDouble > &eta) override
 
void v_LocCollapsedToLocCoord (const Array< OneD, const NekDouble > &eta, Array< OneD, NekDouble > &xi) override
 
void v_GetCoords (Array< OneD, NekDouble > &xi_x, Array< OneD, NekDouble > &xi_y, Array< OneD, NekDouble > &xi_z) override
 
void v_FillMode (const int mode, Array< OneD, NekDouble > &outarray) override
 
NekDouble v_PhysEvaluateBasis (const Array< OneD, const NekDouble > &coords, int mode) final
 
void v_GetTraceNumModes (const int fid, int &numModes0, int &numModes1, Orientation faceOrient=eDir1FwdDir1_Dir2FwdDir2) override
 
NekDouble v_PhysEvalFirstDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) override
 
int v_GetNverts () const override
 
int v_GetNedges () const override
 
int v_GetNtraces () const override
 
LibUtilities::ShapeType v_DetShapeType () const override
 Return Shape of region, using ShapeType enum list; i.e. prism. More...
 
int v_NumBndryCoeffs () const override
 
int v_NumDGBndryCoeffs () const override
 
int v_GetTraceNcoeffs (const int i) const override
 
int v_GetTraceIntNcoeffs (const int i) const override
 
int v_GetTraceNumPoints (const int i) const override
 
int v_GetEdgeNcoeffs (const int i) const override
 
const LibUtilities::BasisKey v_GetTraceBasisKey (const int i, const int k, bool UseGLL=false) const override
 
LibUtilities::PointsKey v_GetTracePointsKey (const int i, const int j) const override
 
int v_CalcNumberOfCoefficients (const std::vector< unsigned int > &nummodes, int &modes_offset) override
 
bool v_IsBoundaryInteriorExpansion () const override
 
int v_GetVertexMap (int localVertexId, bool useCoeffPacking=false) override
 
void v_GetInteriorMap (Array< OneD, unsigned int > &outarray) override
 
void v_GetBoundaryMap (Array< OneD, unsigned int > &outarray) override
 
void v_GetTraceCoeffMap (const unsigned int fid, Array< OneD, unsigned int > &maparray) override
 
void v_GetElmtTraceToTraceMap (const unsigned int fid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation faceOrient, int P, int Q) override
 
void v_GetEdgeInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2) override
 
void v_GetTraceInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, const Orientation traceOrient=eDir1FwdDir1_Dir2FwdDir2) override
 
DNekMatSharedPtr v_GenMatrix (const StdMatrixKey &mkey) override
 
DNekMatSharedPtr v_CreateStdMatrix (const StdMatrixKey &mkey) override
 
void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
void v_SVVLaplacianFilter (Array< OneD, NekDouble > &array, const StdMatrixKey &mkey) override
 
void v_ReduceOrderCoeffs (int numMin, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion3D
NekDouble v_PhysEvaluate (const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals) override
 This function evaluates the expansion at a single (arbitrary) point of the domain. More...
 
NekDouble v_PhysEvaluateInterp (const Array< OneD, DNekMatSharedPtr > &I, const Array< OneD, const NekDouble > &physvals) override
 
virtual void v_BwdTrans_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)=0
 
virtual void v_IProductWRTBase_SumFacKernel (const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)=0
 
void v_LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
void v_HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
 
NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray) override
 Integrates the specified function over the domain. More...
 
virtual int v_GetNedges (void) const
 
virtual int v_GetEdgeNcoeffs (const int i) const
 
NekDouble BaryTensorDeriv (const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs)
 
virtual void v_GetEdgeInteriorToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient=eForwards)
 
void v_GetTraceToElementMap (const int tid, Array< OneD, unsigned int > &maparray, Array< OneD, int > &signarray, Orientation traceOrient, int P, int Q) override
 
void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat) override
 
void v_PhysInterp (std::shared_ptr< StdExpansion > fromExp, const Array< OneD, const NekDouble > &fromData, Array< OneD, NekDouble > &toData) override
 
- Protected Member Functions inherited from Nektar::StdRegions::StdExpansion
DNekMatSharedPtr CreateStdMatrix (const StdMatrixKey &mkey)
 
DNekBlkMatSharedPtr CreateStdStaticCondMatrix (const StdMatrixKey &mkey)
 Create the static condensation of a matrix when using a boundary interior decomposition. More...
 
void BwdTrans_SumFac (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDerivBase_SumFac (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void IProductWRTDirectionalDerivBase_SumFac (const Array< OneD, const NekDouble > &direction, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void GeneralMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree_Kernel (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp)
 
void LaplacianMatrixOp_MatFree_GenericImpl (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LaplacianMatrixOp_MatFree (const int k1, const int k2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDerivMatrixOp_MatFree (const int i, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void WeakDirectionalDerivMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void MassLevelCurvatureMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void LinearAdvectionDiffusionReactionMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey, bool addDiffusionTerm=true)
 
void HelmholtzMatrixOp_MatFree (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
void HelmholtzMatrixOp_MatFree_GenericImpl (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
 
virtual void v_SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual NekDouble v_StdPhysEvaluate (const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
 
virtual void v_GenStdMatBwdDeriv (const int dir, DNekMatSharedPtr &mat)
 
virtual void v_MultiplyByStdQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals, NekDouble &deriv, NekDouble &deriv2)
 This function performs the barycentric interpolation of the polynomial stored in coord at a point physvals using barycentric interpolation weights in direction. More...
 
template<int DIR>
NekDouble BaryEvaluateBasis (const NekDouble &coord, const int &mode)
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals)
 Helper function to pass an unused value by reference into BaryEvaluate. More...
 
template<int DIR, bool DERIV = false, bool DERIV2 = false>
NekDouble BaryEvaluate (const NekDouble &coord, const NekDouble *physvals, NekDouble &deriv)
 
- Protected Member Functions inherited from Nektar::LocalRegions::Expansion3D
void v_DGDeriv (const int dir, const Array< OneD, const NekDouble > &incoeffs, Array< OneD, ExpansionSharedPtr > &FaceExp, Array< OneD, Array< OneD, NekDouble > > &faceCoeffs, Array< OneD, NekDouble > &out_d) override
 Evaluate coefficients of weak deriviative in the direction dir given the input coefficicents incoeffs and the imposed boundary values in EdgeExp (which will have its phys space updated). More...
 
DNekMatSharedPtr v_GenMatrix (const StdRegions::StdMatrixKey &mkey) override
 
void v_AddFaceNormBoundaryInt (const int face, const ExpansionSharedPtr &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray) override
 
void v_AddRobinMassMatrix (const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat) override
 
StdRegions::Orientation v_GetTraceOrient (int face) override
 
void v_GetTracePhysVals (const int face, const StdRegions::StdExpansionSharedPtr &FaceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient) override
 Extract the physical values along face face from inarray into outarray following the local face orientation and point distribution defined by defined in FaceExp. More...
 
void v_GenTraceExp (const int traceid, ExpansionSharedPtr &exp) override
 
void GetPhysFaceVarCoeffsFromElement (const int face, ExpansionSharedPtr &FaceExp, const Array< OneD, const NekDouble > &varcoeff, Array< OneD, NekDouble > &outarray)
 
DNekMatSharedPtr v_BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType) override
 
DNekMatSharedPtr v_BuildInverseTransformationMatrix (const DNekScalMatSharedPtr &transformationmatrix) override
 Build inverse and inverse transposed transformation matrix: \(\mathbf{R^{-1}}\) and \(\mathbf{R^{-T}}\). More...
 
DNekMatSharedPtr v_BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd) override
 
void v_TraceNormLen (const int traceid, NekDouble &h, NekDouble &p) override
 
- Protected Member Functions inherited from Nektar::LocalRegions::Expansion
void ComputeLaplacianMetric ()
 
void ComputeQuadratureMetric ()
 
void ComputeGmatcdotMF (const Array< TwoD, const NekDouble > &df, const Array< OneD, const NekDouble > &direction, Array< OneD, Array< OneD, NekDouble > > &dfdir)
 
Array< OneD, NekDoubleGetMF (const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
 
Array< OneD, NekDoubleGetMFDiv (const int dir, const StdRegions::VarCoeffMap &varcoeffs)
 
Array< OneD, NekDoubleGetMFMag (const int dir, const StdRegions::VarCoeffMap &varcoeffs)
 
void v_MultiplyByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_DivideByQuadratureMetric (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_ComputeLaplacianMetric ()
 
int v_GetCoordim () const override
 
void v_GetCoords (Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
 
virtual DNekScalMatSharedPtr v_GetLocMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual void v_DropLocMatrix (const LocalRegions::MatrixKey &mkey)
 
virtual DNekMatSharedPtr v_BuildTransformationMatrix (const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
 
virtual DNekMatSharedPtr v_BuildVertexMatrix (const DNekScalMatSharedPtr &r_bnd)
 
virtual void v_ExtractDataToCoeffs (const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
 
virtual void v_AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddEdgeNormBoundaryInt (const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddFaceNormBoundaryInt (const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
virtual void v_DGDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble > > &coeffs, Array< OneD, NekDouble > &outarray)
 
virtual NekDouble v_VectorFlux (const Array< OneD, Array< OneD, NekDouble > > &vec)
 
virtual void v_NormalTraceDerivFactors (Array< OneD, Array< OneD, NekDouble > > &factors, Array< OneD, Array< OneD, NekDouble > > &d0factors, Array< OneD, Array< OneD, NekDouble > > &d1factors)
 
virtual void v_AlignVectorToCollapsedDir (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual StdRegions::Orientation v_GetTraceOrient (int trace)
 
void v_SetCoeffsToOrientation (StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
 
virtual void v_GetTraceQFactors (const int trace, Array< OneD, NekDouble > &outarray)
 
virtual void v_GetTracePhysVals (const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient)
 
virtual void v_GetTracePhysMap (const int edge, Array< OneD, int > &outarray)
 
virtual void v_ReOrientTracePhysMap (const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1=-1)
 
virtual void v_ComputeTraceNormal (const int id)
 
virtual const Array< OneD, const NekDouble > & v_GetPhysNormals ()
 
virtual void v_SetPhysNormals (Array< OneD, const NekDouble > &normal)
 
virtual void v_SetUpPhysNormals (const int id)
 
virtual void v_AddRobinMassMatrix (const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
 
virtual void v_AddRobinTraceContribution (const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
 
virtual void v_TraceNormLen (const int traceid, NekDouble &h, NekDouble &p)
 
virtual void v_GenTraceExp (const int traceid, ExpansionSharedPtr &exp)
 

Private Member Functions

void v_LaplacianMatrixOp_MatFree_Kernel (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp) override
 Calculate the Laplacian multiplication in a matrix-free manner. More...
 

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

Definition at line 47 of file PrismExp.h.

Constructor & Destructor Documentation

◆ PrismExp() [1/2]

Nektar::LocalRegions::PrismExp::PrismExp ( const LibUtilities::BasisKey Ba,
const LibUtilities::BasisKey Bb,
const LibUtilities::BasisKey Bc,
const SpatialDomains::PrismGeomSharedPtr geom 
)

Constructor using BasisKey class for quadrature points and order definition.

Definition at line 45 of file PrismExp.cpp.

50 Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes()),
51 3, Ba, Bb, Bc),
53 Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes()),
54 Ba, Bb, Bc),
55 StdPrismExp(Ba, Bb, Bc), Expansion(geom), Expansion3D(geom),
57 std::bind(&Expansion3D::CreateMatrix, this, std::placeholders::_1),
58 std::string("PrismExpMatrix")),
60 this, std::placeholders::_1),
61 std::string("PrismExpStaticCondMatrix"))
62{
63}
Expansion3D(SpatialDomains::Geometry3DSharedPtr pGeom)
Definition: Expansion3D.h:59
DNekScalMatSharedPtr CreateMatrix(const MatrixKey &mkey)
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
Definition: Expansion.cpp:272
Expansion(SpatialDomains::GeometrySharedPtr pGeom)
Definition: Expansion.cpp:43
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLess > m_staticCondMatrixManager
Definition: PrismExp.h:200
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager
Definition: PrismExp.h:198
StdExpansion()
Default Constructor.
StdPrismExp(const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::BasisKey &Bc)
Definition: StdPrismExp.cpp:43
int getNumberOfCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:279

◆ PrismExp() [2/2]

Nektar::LocalRegions::PrismExp::PrismExp ( const PrismExp T)

Definition at line 65 of file PrismExp.cpp.

67 Expansion3D(T), m_matrixManager(T.m_matrixManager),
68 m_staticCondMatrixManager(T.m_staticCondMatrixManager)
69{
70}

◆ ~PrismExp()

Nektar::LocalRegions::PrismExp::~PrismExp ( )
overridedefault

Member Function Documentation

◆ v_AlignVectorToCollapsedDir()

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

381{
382 const int nquad0 = m_base[0]->GetNumPoints();
383 const int nquad1 = m_base[1]->GetNumPoints();
384 const int nquad2 = m_base[2]->GetNumPoints();
385 const int order0 = m_base[0]->GetNumModes();
386 const int order1 = m_base[1]->GetNumModes();
387 const int nqtot = nquad0 * nquad1 * nquad2;
388
389 const Array<OneD, const NekDouble> &z0 = m_base[0]->GetZ();
390 const Array<OneD, const NekDouble> &z2 = m_base[2]->GetZ();
391
392 Array<OneD, NekDouble> gfac0(nquad0);
393 Array<OneD, NekDouble> gfac2(nquad2);
394 Array<OneD, NekDouble> tmp1(nqtot);
395 Array<OneD, NekDouble> tmp5(nqtot);
396 Array<OneD, NekDouble> tmp6(m_ncoeffs);
397 Array<OneD, NekDouble> wsp(order0 * nquad2 * (nquad1 + order1));
398
399 Array<OneD, NekDouble> tmp2 = outarray[0];
400 Array<OneD, NekDouble> tmp3 = outarray[1];
401 Array<OneD, NekDouble> tmp4 = outarray[2];
402
403 const Array<TwoD, const NekDouble> &df =
404 m_metricinfo->GetDerivFactors(GetPointsKeys());
405
406 Vmath::Vcopy(nqtot, inarray, 1, tmp1, 1); // Dir3 metric
407
408 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
409 {
410 Vmath::Vmul(nqtot, &df[3 * dir][0], 1, tmp1.data(), 1, tmp2.data(), 1);
411 Vmath::Vmul(nqtot, &df[3 * dir + 1][0], 1, tmp1.data(), 1, tmp3.data(),
412 1);
413 Vmath::Vmul(nqtot, &df[3 * dir + 2][0], 1, tmp1.data(), 1, tmp4.data(),
414 1);
415 }
416 else
417 {
418 Vmath::Smul(nqtot, df[3 * dir][0], tmp1.data(), 1, tmp2.data(), 1);
419 Vmath::Smul(nqtot, df[3 * dir + 1][0], tmp1.data(), 1, tmp3.data(), 1);
420 Vmath::Smul(nqtot, df[3 * dir + 2][0], tmp1.data(), 1, tmp4.data(), 1);
421 }
422
423 // set up geometric factor: (1+z0)/2
424 for (int i = 0; i < nquad0; ++i)
425 {
426 gfac0[i] = 0.5 * (1 + z0[i]);
427 }
428
429 // Set up geometric factor: 2/(1-z2)
430 for (int i = 0; i < nquad2; ++i)
431 {
432 gfac2[i] = 2.0 / (1 - z2[i]);
433 }
434
435 const int nq01 = nquad0 * nquad1;
436
437 for (int i = 0; i < nquad2; ++i)
438 {
439 Vmath::Smul(nq01, gfac2[i], &tmp2[0] + i * nq01, 1, &tmp2[0] + i * nq01,
440 1);
441 Vmath::Smul(nq01, gfac2[i], &tmp4[0] + i * nq01, 1, &tmp5[0] + i * nq01,
442 1);
443 }
444
445 for (int i = 0; i < nquad1 * nquad2; ++i)
446 {
447 Vmath::Vmul(nquad0, &gfac0[0], 1, &tmp5[0] + i * nquad0, 1,
448 &tmp5[0] + i * nquad0, 1);
449 }
450
451 Vmath::Vadd(nqtot, &tmp2[0], 1, &tmp5[0], 1, &tmp2[0], 1);
452}
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
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

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

Referenced by v_IProductWRTDerivBase_SumFac().

◆ v_ComputeTraceNormal()

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

Get the normals along specficied face Get the face normals interplated to a points0 x points 0 type distribution.

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 704 of file PrismExp.cpp.

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

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

◆ v_CreateStdMatrix()

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 1069 of file PrismExp.cpp.

1071{
1072 LibUtilities::BasisKey bkey0 = m_base[0]->GetBasisKey();
1073 LibUtilities::BasisKey bkey1 = m_base[1]->GetBasisKey();
1074 LibUtilities::BasisKey bkey2 = m_base[2]->GetBasisKey();
1077
1078 return tmp->GetStdMatrix(mkey);
1079}
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< StdPrismExp > StdPrismExpSharedPtr
Definition: StdPrismExp.h:218

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

◆ v_DropLocMatrix()

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

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 1086 of file PrismExp.cpp.

1087{
1088 m_matrixManager.DeleteObject(mkey);
1089}

References m_matrixManager.

◆ v_DropLocStaticCondMatrix()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1097 of file PrismExp.cpp.

1098{
1099 m_staticCondMatrixManager.DeleteObject(mkey);
1100}

References m_staticCondMatrixManager.

◆ v_ExtractDataToCoeffs()

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

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

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

◆ v_FwdTrans()

void Nektar::LocalRegions::PrismExp::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)->m_coeffs.

Inputs:

  • inarray: array of physical quadrature points to be transformed

Outputs:

  • (this)->_coeffs: updated array of expansion coefficients.

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 209 of file PrismExp.cpp.

211{
212 if (m_base[0]->Collocation() && m_base[1]->Collocation() &&
213 m_base[2]->Collocation())
214 {
215 Vmath::Vcopy(GetNcoeffs(), &inarray[0], 1, &outarray[0], 1);
216 }
217 else
218 {
219 v_IProductWRTBase(inarray, outarray);
220
221 // get Mass matrix inverse
222 MatrixKey masskey(StdRegions::eInvMass, DetShapeType(), *this);
223 DNekScalMatSharedPtr matsys = m_matrixManager[masskey];
224
225 // copy inarray in case inarray == outarray
226 DNekVec in(m_ncoeffs, outarray);
227 DNekVec out(m_ncoeffs, outarray, eWrapper);
228
229 out = (*matsys) * in;
230 }
231}
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*base1*base2 and put into out...
Definition: PrismExp.cpp:261
int GetNcoeffs(void) const
This function returns the total number of coefficients used in the expansion.
Definition: StdExpansion.h:124
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
Definition: StdExpansion.h:370
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
NekVector< NekDouble > DNekVec
Definition: NekTypeDefs.hpp:48

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

◆ v_GenMatrix()

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 1046 of file PrismExp.cpp.

1047{
1048 DNekMatSharedPtr returnval;
1049
1050 switch (mkey.GetMatrixType())
1051 {
1059 returnval = Expansion3D::v_GenMatrix(mkey);
1060 break;
1061 default:
1062 returnval = StdPrismExp::v_GenMatrix(mkey);
1063 break;
1064 }
1065
1066 return returnval;
1067}
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::PrismExp::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 482 of file PrismExp.cpp.

484{
485 int i;
486
487 ASSERTL1(Lcoords[0] <= -1.0 && Lcoords[0] >= 1.0 && Lcoords[1] <= -1.0 &&
488 Lcoords[1] >= 1.0 && Lcoords[2] <= -1.0 && Lcoords[2] >= 1.0,
489 "Local coordinates are not in region [-1,1]");
490
491 m_geom->FillGeom();
492
493 for (i = 0; i < m_geom->GetCoordim(); ++i)
494 {
495 coords[i] = m_geom->GetCoord(i, Lcoords);
496 }
497}
SpatialDomains::GeometrySharedPtr m_geom
Definition: Expansion.h:273

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

◆ v_GetCoords()

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 499 of file PrismExp.cpp.

502{
503 Expansion::v_GetCoords(coords_0, coords_1, coords_2);
504}
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::PrismExp::v_GetLinStdExp ( void  ) const
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 465 of file PrismExp.cpp.

466{
467 LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(), 2,
468 m_base[0]->GetPointsKey());
469 LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(), 2,
470 m_base[1]->GetPointsKey());
471 LibUtilities::BasisKey bkey2(m_base[2]->GetBasisType(), 2,
472 m_base[2]->GetPointsKey());
473
475 bkey0, bkey1, bkey2);
476}

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

◆ v_GetLocMatrix()

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

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 1081 of file PrismExp.cpp.

1082{
1083 return m_matrixManager[mkey];
1084}

References m_matrixManager.

◆ v_GetLocStaticCondMatrix()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1091 of file PrismExp.cpp.

1093{
1094 return m_staticCondMatrixManager[mkey];
1095}

References m_staticCondMatrixManager.

◆ v_GetSimplexEquiSpacedConnectivity()

void Nektar::LocalRegions::PrismExp::v_GetSimplexEquiSpacedConnectivity ( Array< OneD, int > &  conn,
bool  standard = true 
)
overrideprotectedvirtual

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1308 of file PrismExp.cpp.

1310{
1311 int np0 = m_base[0]->GetNumPoints();
1312 int np1 = m_base[1]->GetNumPoints();
1313 int np2 = m_base[2]->GetNumPoints();
1314 int np = max(np0, max(np1, np2));
1315 Array<OneD, int> prismpt(6);
1316 bool standard = true;
1317
1318 int vid0 = m_geom->GetVid(0);
1319 int vid1 = m_geom->GetVid(1);
1320 int vid2 = m_geom->GetVid(4);
1321 int rotate = 0;
1322
1323 // sort out prism rotation according to
1324 if ((vid2 < vid1) && (vid2 < vid0)) // top triangle vertex is lowest id
1325 {
1326 rotate = 0;
1327 if (vid0 > vid1)
1328 {
1329 standard = false; // reverse base direction
1330 }
1331 }
1332 else if ((vid1 < vid2) && (vid1 < vid0))
1333 {
1334 rotate = 1;
1335 if (vid2 > vid0)
1336 {
1337 standard = false; // reverse base direction
1338 }
1339 }
1340 else if ((vid0 < vid2) && (vid0 < vid1))
1341 {
1342 rotate = 2;
1343 if (vid1 > vid2)
1344 {
1345 standard = false; // reverse base direction
1346 }
1347 }
1348
1349 conn = Array<OneD, int>(12 * (np - 1) * (np - 1) * (np - 1));
1350
1351 int row = 0;
1352 int rowp1 = 0;
1353 int plane = 0;
1354 int row1 = 0;
1355 int row1p1 = 0;
1356 int planep1 = 0;
1357 int cnt = 0;
1358
1359 Array<OneD, int> rot(3);
1360
1361 rot[0] = (0 + rotate) % 3;
1362 rot[1] = (1 + rotate) % 3;
1363 rot[2] = (2 + rotate) % 3;
1364
1365 // lower diagonal along 1-3 on base
1366 for (int i = 0; i < np - 1; ++i)
1367 {
1368 planep1 += (np - i) * np;
1369 row = 0; // current plane row offset
1370 rowp1 = 0; // current plane row plus one offset
1371 row1 = 0; // next plane row offset
1372 row1p1 = 0; // nex plane row plus one offset
1373 if (standard == false)
1374 {
1375 for (int j = 0; j < np - 1; ++j)
1376 {
1377 rowp1 += np - i;
1378 row1p1 += np - i - 1;
1379 for (int k = 0; k < np - i - 2; ++k)
1380 {
1381 // bottom prism block
1382 prismpt[rot[0]] = plane + row + k;
1383 prismpt[rot[1]] = plane + row + k + 1;
1384 prismpt[rot[2]] = planep1 + row1 + k;
1385
1386 prismpt[3 + rot[0]] = plane + rowp1 + k;
1387 prismpt[3 + rot[1]] = plane + rowp1 + k + 1;
1388 prismpt[3 + rot[2]] = planep1 + row1p1 + k;
1389
1390 conn[cnt++] = prismpt[0];
1391 conn[cnt++] = prismpt[1];
1392 conn[cnt++] = prismpt[3];
1393 conn[cnt++] = prismpt[2];
1394
1395 conn[cnt++] = prismpt[5];
1396 conn[cnt++] = prismpt[2];
1397 conn[cnt++] = prismpt[3];
1398 conn[cnt++] = prismpt[4];
1399
1400 conn[cnt++] = prismpt[3];
1401 conn[cnt++] = prismpt[1];
1402 conn[cnt++] = prismpt[4];
1403 conn[cnt++] = prismpt[2];
1404
1405 // upper prism block.
1406 prismpt[rot[0]] = planep1 + row1 + k + 1;
1407 prismpt[rot[1]] = planep1 + row1 + k;
1408 prismpt[rot[2]] = plane + row + k + 1;
1409
1410 prismpt[3 + rot[0]] = planep1 + row1p1 + k + 1;
1411 prismpt[3 + rot[1]] = planep1 + row1p1 + k;
1412 prismpt[3 + rot[2]] = plane + rowp1 + k + 1;
1413
1414 conn[cnt++] = prismpt[0];
1415 conn[cnt++] = prismpt[1];
1416 conn[cnt++] = prismpt[2];
1417 conn[cnt++] = prismpt[5];
1418
1419 conn[cnt++] = prismpt[5];
1420 conn[cnt++] = prismpt[0];
1421 conn[cnt++] = prismpt[4];
1422 conn[cnt++] = prismpt[1];
1423
1424 conn[cnt++] = prismpt[3];
1425 conn[cnt++] = prismpt[4];
1426 conn[cnt++] = prismpt[0];
1427 conn[cnt++] = prismpt[5];
1428 }
1429
1430 // bottom prism block
1431 prismpt[rot[0]] = plane + row + np - i - 2;
1432 prismpt[rot[1]] = plane + row + np - i - 1;
1433 prismpt[rot[2]] = planep1 + row1 + np - i - 2;
1434
1435 prismpt[3 + rot[0]] = plane + rowp1 + np - i - 2;
1436 prismpt[3 + rot[1]] = plane + rowp1 + np - i - 1;
1437 prismpt[3 + rot[2]] = planep1 + row1p1 + np - i - 2;
1438
1439 conn[cnt++] = prismpt[0];
1440 conn[cnt++] = prismpt[1];
1441 conn[cnt++] = prismpt[3];
1442 conn[cnt++] = prismpt[2];
1443
1444 conn[cnt++] = prismpt[5];
1445 conn[cnt++] = prismpt[2];
1446 conn[cnt++] = prismpt[3];
1447 conn[cnt++] = prismpt[4];
1448
1449 conn[cnt++] = prismpt[3];
1450 conn[cnt++] = prismpt[1];
1451 conn[cnt++] = prismpt[4];
1452 conn[cnt++] = prismpt[2];
1453
1454 row += np - i;
1455 row1 += np - i - 1;
1456 }
1457 }
1458 else
1459 { // lower diagonal along 0-4 on base
1460 for (int j = 0; j < np - 1; ++j)
1461 {
1462 rowp1 += np - i;
1463 row1p1 += np - i - 1;
1464 for (int k = 0; k < np - i - 2; ++k)
1465 {
1466 // bottom prism block
1467 prismpt[rot[0]] = plane + row + k;
1468 prismpt[rot[1]] = plane + row + k + 1;
1469 prismpt[rot[2]] = planep1 + row1 + k;
1470
1471 prismpt[3 + rot[0]] = plane + rowp1 + k;
1472 prismpt[3 + rot[1]] = plane + rowp1 + k + 1;
1473 prismpt[3 + rot[2]] = planep1 + row1p1 + k;
1474
1475 conn[cnt++] = prismpt[0];
1476 conn[cnt++] = prismpt[1];
1477 conn[cnt++] = prismpt[4];
1478 conn[cnt++] = prismpt[2];
1479
1480 conn[cnt++] = prismpt[4];
1481 conn[cnt++] = prismpt[3];
1482 conn[cnt++] = prismpt[0];
1483 conn[cnt++] = prismpt[2];
1484
1485 conn[cnt++] = prismpt[3];
1486 conn[cnt++] = prismpt[4];
1487 conn[cnt++] = prismpt[5];
1488 conn[cnt++] = prismpt[2];
1489
1490 // upper prism block.
1491 prismpt[rot[0]] = planep1 + row1 + k + 1;
1492 prismpt[rot[1]] = planep1 + row1 + k;
1493 prismpt[rot[2]] = plane + row + k + 1;
1494
1495 prismpt[3 + rot[0]] = planep1 + row1p1 + k + 1;
1496 prismpt[3 + rot[1]] = planep1 + row1p1 + k;
1497 prismpt[3 + rot[2]] = plane + rowp1 + k + 1;
1498
1499 conn[cnt++] = prismpt[0];
1500 conn[cnt++] = prismpt[2];
1501 conn[cnt++] = prismpt[1];
1502 conn[cnt++] = prismpt[5];
1503
1504 conn[cnt++] = prismpt[3];
1505 conn[cnt++] = prismpt[5];
1506 conn[cnt++] = prismpt[0];
1507 conn[cnt++] = prismpt[1];
1508
1509 conn[cnt++] = prismpt[5];
1510 conn[cnt++] = prismpt[3];
1511 conn[cnt++] = prismpt[4];
1512 conn[cnt++] = prismpt[1];
1513 }
1514
1515 // bottom prism block
1516 prismpt[rot[0]] = plane + row + np - i - 2;
1517 prismpt[rot[1]] = plane + row + np - i - 1;
1518 prismpt[rot[2]] = planep1 + row1 + np - i - 2;
1519
1520 prismpt[3 + rot[0]] = plane + rowp1 + np - i - 2;
1521 prismpt[3 + rot[1]] = plane + rowp1 + np - i - 1;
1522 prismpt[3 + rot[2]] = planep1 + row1p1 + np - i - 2;
1523
1524 conn[cnt++] = prismpt[0];
1525 conn[cnt++] = prismpt[1];
1526 conn[cnt++] = prismpt[4];
1527 conn[cnt++] = prismpt[2];
1528
1529 conn[cnt++] = prismpt[4];
1530 conn[cnt++] = prismpt[3];
1531 conn[cnt++] = prismpt[0];
1532 conn[cnt++] = prismpt[2];
1533
1534 conn[cnt++] = prismpt[3];
1535 conn[cnt++] = prismpt[4];
1536 conn[cnt++] = prismpt[5];
1537 conn[cnt++] = prismpt[2];
1538
1539 row += np - i;
1540 row1 += np - i - 1;
1541 }
1542 }
1543 plane += (np - i) * np;
1544 }
1545}

References Nektar::StdRegions::StdExpansion::m_base, and Nektar::LocalRegions::Expansion::m_geom.

◆ v_GetStdExp()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 458 of file PrismExp.cpp.

459{
461 m_base[0]->GetBasisKey(), m_base[1]->GetBasisKey(),
462 m_base[2]->GetBasisKey());
463}

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

◆ v_GetTracePhysMap()

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

Reimplemented from Nektar::LocalRegions::Expansion.

Definition at line 603 of file PrismExp.cpp.

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

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

◆ v_HelmholtzMatrixOp()

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

1011{
1012 PrismExp::v_HelmholtzMatrixOp_MatFree(inarray, outarray, mkey);
1013}
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::PrismExp::v_Integral ( const Array< OneD, const NekDouble > &  inarray)
overrideprotectedvirtual

Integrate the physical point list inarray over prismatic region and return the value.

Inputs:

  • inarray: definition of function to be returned at quadrature point of expansion.

Outputs:

  • returns \(\int^1_{-1}\int^1_{-1}\int^1_{-1} u(\bar \eta_1, \xi_2, \xi_3) J[i,j,k] d \bar \eta_1 d \xi_2 d \xi_3 \)
    \( = \sum_{i=0}^{Q_1 - 1} \sum_{j=0}^{Q_2 - 1} \sum_{k=0}^{Q_3 - 1} u(\bar \eta_{1i}^{0,0}, \xi_{2j}^{0,0},\xi_{3k}^{1,0})w_{i}^{0,0} w_{j}^{0,0} \hat w_{k}^{1,0} \)
    where \( inarray[i,j, k] = u(\bar \eta_{1i}^{0,0}, \xi_{2j}^{0,0},\xi_{3k}^{1,0}) \),
    \(\hat w_{i}^{1,0} = \frac {w_{j}^{1,0}} {2} \)
    and \( J[i,j,k] \) is the Jacobian evaluated at the quadrature point.

Reimplemented from Nektar::StdRegions::StdExpansion3D.

Definition at line 96 of file PrismExp.cpp.

97{
98 int nquad0 = m_base[0]->GetNumPoints();
99 int nquad1 = m_base[1]->GetNumPoints();
100 int nquad2 = m_base[2]->GetNumPoints();
101 Array<OneD, const NekDouble> jac = m_metricinfo->GetJac(GetPointsKeys());
102 Array<OneD, NekDouble> tmp(nquad0 * nquad1 * nquad2);
103
104 // Multiply inarray with Jacobian
105 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
106 {
107 Vmath::Vmul(nquad0 * nquad1 * nquad2, &jac[0], 1,
108 (NekDouble *)&inarray[0], 1, &tmp[0], 1);
109 }
110 else
111 {
112 Vmath::Smul(nquad0 * nquad1 * nquad2, (NekDouble)jac[0],
113 (NekDouble *)&inarray[0], 1, &tmp[0], 1);
114 }
115
116 // Call StdPrismExp version.
117 return StdPrismExp::v_Integral(tmp);
118}

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::PrismExp::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=base0*base1*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} (\bar \eta_{1i}) \psi_{q}^{a} (\xi_{2j}) \psi_{pr}^{b} (\xi_{3k}) w_i w_j w_k u(\bar \eta_{1,i} \xi_{2,j} \xi_{3,k}) J_{i,j,k}\\ & = & \sum_{i=0}^{nq_0} \psi_p^a(\bar \eta_{1,i}) \sum_{j=0}^{nq_1} \psi_{q}^a(\xi_{2,j}) \sum_{k=0}^{nq_2} \psi_{pr}^b u(\bar \eta_{1i},\xi_{2j},\xi_{3k}) J_{i,j,k} \end{array} \)
where

\( \phi_{pqr} (\xi_1 , \xi_2 , \xi_3) = \psi_p^a (\bar \eta_1) \psi_{q}^a (\xi_2) \psi_{pr}^b (\xi_3) \)
which can be implemented as
\(f_{pr} (\xi_{3k}) = \sum_{k=0}^{nq_3} \psi_{pr}^b u(\bar \eta_{1i},\xi_{2j},\xi_{3k}) J_{i,j,k} = {\bf B_3 U} \)
\( g_{q} (\xi_{3k}) = \sum_{j=0}^{nq_1} \psi_{q}^a (\xi_{2j}) f_{pr} (\xi_{3k}) = {\bf B_2 F} \)
\( (\phi_{pqr}, u)_{\delta} = \sum_{k=0}^{nq_0} \psi_{p}^a (\xi_{3k}) g_{q} (\xi_{3k}) = {\bf B_1 G} \)

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 261 of file PrismExp.cpp.

263{
264 v_IProductWRTBase_SumFac(inarray, outarray);
265}
void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true) override
Definition: PrismExp.cpp:267

References v_IProductWRTBase_SumFac().

Referenced by v_FwdTrans().

◆ v_IProductWRTBase_SumFac()

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 267 of file PrismExp.cpp.

270{
271 const int nquad0 = m_base[0]->GetNumPoints();
272 const int nquad1 = m_base[1]->GetNumPoints();
273 const int nquad2 = m_base[2]->GetNumPoints();
274 const int order0 = m_base[0]->GetNumModes();
275 const int order1 = m_base[1]->GetNumModes();
276
277 Array<OneD, NekDouble> wsp(order0 * nquad2 * (nquad1 + order1));
278
279 if (multiplybyweights)
280 {
281 Array<OneD, NekDouble> tmp(nquad0 * nquad1 * nquad2);
282
283 MultiplyByQuadratureMetric(inarray, tmp);
284
286 m_base[0]->GetBdata(), m_base[1]->GetBdata(), m_base[2]->GetBdata(),
287 tmp, outarray, wsp, true, true, true);
288 }
289 else
290 {
292 m_base[0]->GetBdata(), m_base[1]->GetBdata(), m_base[2]->GetBdata(),
293 inarray, outarray, wsp, true, true, true);
294 }
295}
void IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &base2, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0, bool doCheckCollDir1, bool doCheckCollDir2)
void MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdExpansion.h:732

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::PrismExp::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 tetrahedral 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::StdPrismExp.

Definition at line 327 of file PrismExp.cpp.

330{
331 v_IProductWRTDerivBase_SumFac(dir, inarray, outarray);
332}
void v_IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Definition: PrismExp.cpp:334

References v_IProductWRTDerivBase_SumFac().

◆ v_IProductWRTDerivBase_SumFac()

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 334 of file PrismExp.cpp.

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

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().

Referenced by v_IProductWRTDerivBase().

◆ v_LaplacianMatrixOp() [1/2]

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

996{
997 PrismExp::LaplacianMatrixOp_MatFree(inarray, outarray, mkey);
998}
void LaplacianMatrixOp_MatFree(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)

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

◆ v_LaplacianMatrixOp() [2/2]

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

1004{
1005 StdExpansion::LaplacianMatrixOp_MatFree(k1, k2, inarray, outarray, mkey);
1006}

◆ v_LaplacianMatrixOp_MatFree_Kernel()

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

Calculate the Laplacian multiplication in a matrix-free manner.

This function is the kernel of the Laplacian matrix-free operator, and is used in v_HelmholtzMatrixOp_MatFree to determine the effect of the Helmholtz operator in a similar fashion.

The majority of the calculation is precisely the same as in the hexahedral expansion; however the collapsed co-ordinate system must be taken into account when constructing the geometric factors. How this is done is detailed more exactly in the tetrahedral expansion. On entry to this function, the input #inarray must be in its backwards-transformed state (i.e. \(\mathbf{u} = \mathbf{B}\hat{\mathbf{u}}\)). The output is in coefficient space.

See also
TetExp::v_HelmholtzMatrixOp_MatFree

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 1120 of file PrismExp.cpp.

1123{
1124 int nquad0 = m_base[0]->GetNumPoints();
1125 int nquad1 = m_base[1]->GetNumPoints();
1126 int nquad2 = m_base[2]->GetNumPoints();
1127 int nqtot = nquad0 * nquad1 * nquad2;
1128 int i;
1129
1130 // Set up temporary storage.
1131 Array<OneD, NekDouble> alloc(11 * nqtot, 0.0);
1132 Array<OneD, NekDouble> wsp1(alloc); // TensorDeriv 1
1133 Array<OneD, NekDouble> wsp2(alloc + 1 * nqtot); // TensorDeriv 2
1134 Array<OneD, NekDouble> wsp3(alloc + 2 * nqtot); // TensorDeriv 3
1135 Array<OneD, NekDouble> g0(alloc + 3 * nqtot); // g0
1136 Array<OneD, NekDouble> g1(alloc + 4 * nqtot); // g1
1137 Array<OneD, NekDouble> g2(alloc + 5 * nqtot); // g2
1138 Array<OneD, NekDouble> g3(alloc + 6 * nqtot); // g3
1139 Array<OneD, NekDouble> g4(alloc + 7 * nqtot); // g4
1140 Array<OneD, NekDouble> g5(alloc + 8 * nqtot); // g5
1141 Array<OneD, NekDouble> h0(alloc + 3 * nqtot); // h0 == g0
1142 Array<OneD, NekDouble> h1(alloc + 6 * nqtot); // h1 == g3
1143 Array<OneD, NekDouble> wsp4(alloc + 4 * nqtot); // wsp4 == g1
1144 Array<OneD, NekDouble> wsp5(alloc + 5 * nqtot); // wsp5 == g2
1145 Array<OneD, NekDouble> wsp6(alloc + 8 * nqtot); // wsp6 == g5
1146 Array<OneD, NekDouble> wsp7(alloc + 3 * nqtot); // wsp7 == g0
1147 Array<OneD, NekDouble> wsp8(alloc + 9 * nqtot); // wsp8
1148 Array<OneD, NekDouble> wsp9(alloc + 10 * nqtot); // wsp9
1149
1150 const Array<OneD, const NekDouble> &base0 = m_base[0]->GetBdata();
1151 const Array<OneD, const NekDouble> &base1 = m_base[1]->GetBdata();
1152 const Array<OneD, const NekDouble> &base2 = m_base[2]->GetBdata();
1153 const Array<OneD, const NekDouble> &dbase0 = m_base[0]->GetDbdata();
1154 const Array<OneD, const NekDouble> &dbase1 = m_base[1]->GetDbdata();
1155 const Array<OneD, const NekDouble> &dbase2 = m_base[2]->GetDbdata();
1156
1157 // Step 1. LAPLACIAN MATRIX OPERATION
1158 // wsp1 = du_dxi1 = D_xi1 * wsp0 = D_xi1 * u
1159 // wsp2 = du_dxi2 = D_xi2 * wsp0 = D_xi2 * u
1160 // wsp3 = du_dxi3 = D_xi3 * wsp0 = D_xi3 * u
1161 StdExpansion3D::PhysTensorDeriv(inarray, wsp1, wsp2, wsp3);
1162
1163 const Array<TwoD, const NekDouble> &df =
1164 m_metricinfo->GetDerivFactors(GetPointsKeys());
1165 const Array<OneD, const NekDouble> &z0 = m_base[0]->GetZ();
1166 const Array<OneD, const NekDouble> &z2 = m_base[2]->GetZ();
1167
1168 // Step 2. Calculate the metric terms of the collapsed
1169 // coordinate transformation (Spencer's book P152)
1170 for (i = 0; i < nquad2; ++i)
1171 {
1172 Vmath::Fill(nquad0 * nquad1, 2.0 / (1.0 - z2[i]),
1173 &h0[0] + i * nquad0 * nquad1, 1);
1174 Vmath::Fill(nquad0 * nquad1, 2.0 / (1.0 - z2[i]),
1175 &h1[0] + i * nquad0 * nquad1, 1);
1176 }
1177 for (i = 0; i < nquad0; i++)
1178 {
1179 Blas::Dscal(nquad1 * nquad2, 0.5 * (1 + z0[i]), &h1[0] + i, nquad0);
1180 }
1181
1182 // Step 3. Construct combined metric terms for physical space to
1183 // collapsed coordinate system. Order of construction optimised
1184 // to minimise temporary storage
1185 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
1186 {
1187 // wsp4 = d eta_1/d x_1
1188 Vmath::Vvtvvtp(nqtot, &df[0][0], 1, &h0[0], 1, &df[2][0], 1, &h1[0], 1,
1189 &wsp4[0], 1);
1190 // wsp5 = d eta_2/d x_1
1191 Vmath::Vvtvvtp(nqtot, &df[3][0], 1, &h0[0], 1, &df[5][0], 1, &h1[0], 1,
1192 &wsp5[0], 1);
1193 // wsp6 = d eta_3/d x_1d
1194 Vmath::Vvtvvtp(nqtot, &df[6][0], 1, &h0[0], 1, &df[8][0], 1, &h1[0], 1,
1195 &wsp6[0], 1);
1196
1197 // g0 (overwrites h0)
1198 Vmath::Vvtvvtp(nqtot, &wsp4[0], 1, &wsp4[0], 1, &wsp5[0], 1, &wsp5[0],
1199 1, &g0[0], 1);
1200 Vmath::Vvtvp(nqtot, &wsp6[0], 1, &wsp6[0], 1, &g0[0], 1, &g0[0], 1);
1201
1202 // g3 (overwrites h1)
1203 Vmath::Vvtvvtp(nqtot, &df[1][0], 1, &wsp4[0], 1, &df[4][0], 1, &wsp5[0],
1204 1, &g3[0], 1);
1205 Vmath::Vvtvp(nqtot, &df[7][0], 1, &wsp6[0], 1, &g3[0], 1, &g3[0], 1);
1206
1207 // g4
1208 Vmath::Vvtvvtp(nqtot, &df[2][0], 1, &wsp4[0], 1, &df[5][0], 1, &wsp5[0],
1209 1, &g4[0], 1);
1210 Vmath::Vvtvp(nqtot, &df[8][0], 1, &wsp6[0], 1, &g4[0], 1, &g4[0], 1);
1211
1212 // Overwrite wsp4/5/6 with g1/2/5
1213 // g1
1214 Vmath::Vvtvvtp(nqtot, &df[1][0], 1, &df[1][0], 1, &df[4][0], 1,
1215 &df[4][0], 1, &g1[0], 1);
1216 Vmath::Vvtvp(nqtot, &df[7][0], 1, &df[7][0], 1, &g1[0], 1, &g1[0], 1);
1217
1218 // g2
1219 Vmath::Vvtvvtp(nqtot, &df[2][0], 1, &df[2][0], 1, &df[5][0], 1,
1220 &df[5][0], 1, &g2[0], 1);
1221 Vmath::Vvtvp(nqtot, &df[8][0], 1, &df[8][0], 1, &g2[0], 1, &g2[0], 1);
1222
1223 // g5
1224 Vmath::Vvtvvtp(nqtot, &df[1][0], 1, &df[2][0], 1, &df[4][0], 1,
1225 &df[5][0], 1, &g5[0], 1);
1226 Vmath::Vvtvp(nqtot, &df[7][0], 1, &df[8][0], 1, &g5[0], 1, &g5[0], 1);
1227 }
1228 else
1229 {
1230 // wsp4 = d eta_1/d x_1
1231 Vmath::Svtsvtp(nqtot, df[0][0], &h0[0], 1, df[2][0], &h1[0], 1,
1232 &wsp4[0], 1);
1233 // wsp5 = d eta_2/d x_1
1234 Vmath::Svtsvtp(nqtot, df[3][0], &h0[0], 1, df[5][0], &h1[0], 1,
1235 &wsp5[0], 1);
1236 // wsp6 = d eta_3/d x_1
1237 Vmath::Svtsvtp(nqtot, df[6][0], &h0[0], 1, df[8][0], &h1[0], 1,
1238 &wsp6[0], 1);
1239
1240 // g0 (overwrites h0)
1241 Vmath::Vvtvvtp(nqtot, &wsp4[0], 1, &wsp4[0], 1, &wsp5[0], 1, &wsp5[0],
1242 1, &g0[0], 1);
1243 Vmath::Vvtvp(nqtot, &wsp6[0], 1, &wsp6[0], 1, &g0[0], 1, &g0[0], 1);
1244
1245 // g3 (overwrites h1)
1246 Vmath::Svtsvtp(nqtot, df[1][0], &wsp4[0], 1, df[4][0], &wsp5[0], 1,
1247 &g3[0], 1);
1248 Vmath::Svtvp(nqtot, df[7][0], &wsp6[0], 1, &g3[0], 1, &g3[0], 1);
1249
1250 // g4
1251 Vmath::Svtsvtp(nqtot, df[2][0], &wsp4[0], 1, df[5][0], &wsp5[0], 1,
1252 &g4[0], 1);
1253 Vmath::Svtvp(nqtot, df[8][0], &wsp6[0], 1, &g4[0], 1, &g4[0], 1);
1254
1255 // Overwrite wsp4/5/6 with g1/2/5
1256 // g1
1257 Vmath::Fill(nqtot,
1258 df[1][0] * df[1][0] + df[4][0] * df[4][0] +
1259 df[7][0] * df[7][0],
1260 &g1[0], 1);
1261
1262 // g2
1263 Vmath::Fill(nqtot,
1264 df[2][0] * df[2][0] + df[5][0] * df[5][0] +
1265 df[8][0] * df[8][0],
1266 &g2[0], 1);
1267
1268 // g5
1269 Vmath::Fill(nqtot,
1270 df[1][0] * df[2][0] + df[4][0] * df[5][0] +
1271 df[7][0] * df[8][0],
1272 &g5[0], 1);
1273 }
1274 // Compute component derivatives into wsp7, 8, 9 (wsp7 overwrites
1275 // g0).
1276 Vmath::Vvtvvtp(nqtot, &g0[0], 1, &wsp1[0], 1, &g3[0], 1, &wsp2[0], 1,
1277 &wsp7[0], 1);
1278 Vmath::Vvtvp(nqtot, &g4[0], 1, &wsp3[0], 1, &wsp7[0], 1, &wsp7[0], 1);
1279 Vmath::Vvtvvtp(nqtot, &g1[0], 1, &wsp2[0], 1, &g3[0], 1, &wsp1[0], 1,
1280 &wsp8[0], 1);
1281 Vmath::Vvtvp(nqtot, &g5[0], 1, &wsp3[0], 1, &wsp8[0], 1, &wsp8[0], 1);
1282 Vmath::Vvtvvtp(nqtot, &g2[0], 1, &wsp3[0], 1, &g4[0], 1, &wsp1[0], 1,
1283 &wsp9[0], 1);
1284 Vmath::Vvtvp(nqtot, &g5[0], 1, &wsp2[0], 1, &wsp9[0], 1, &wsp9[0], 1);
1285
1286 // Step 4.
1287 // Multiply by quadrature metric
1288 MultiplyByQuadratureMetric(wsp7, wsp7);
1289 MultiplyByQuadratureMetric(wsp8, wsp8);
1290 MultiplyByQuadratureMetric(wsp9, wsp9);
1291
1292 // Perform inner product w.r.t derivative bases.
1293 IProductWRTBase_SumFacKernel(dbase0, base1, base2, wsp7, wsp1, wsp, false,
1294 true, true);
1295 IProductWRTBase_SumFacKernel(base0, dbase1, base2, wsp8, wsp2, wsp, true,
1296 false, true);
1297 IProductWRTBase_SumFacKernel(base0, base1, dbase2, wsp9, outarray, wsp,
1298 true, true, false);
1299
1300 // Step 5.
1301 // Sum contributions from wsp1, wsp2 and outarray.
1302 Vmath::Vadd(m_ncoeffs, wsp1.data(), 1, outarray.data(), 1, outarray.data(),
1303 1);
1304 Vmath::Vadd(m_ncoeffs, wsp2.data(), 1, outarray.data(), 1, outarray.data(),
1305 1);
1306}
static void Dscal(const int &n, const double &alpha, double *x, const int &incx)
BLAS level 1: x = alpha x.
Definition: Blas.hpp:149
void Svtsvtp(int n, const T alpha, const T *x, int incx, const T beta, const T *y, int incy, T *z, int incz)
Svtsvtp (scalar times vector plus scalar times vector):
Definition: Vmath.hpp:473
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Svtvp (scalar times vector plus vector): z = alpha*x + y.
Definition: Vmath.hpp:396
void Vvtvvtp(int n, const T *v, int incv, const T *w, int incw, const T *x, int incx, const T *y, int incy, T *z, int incz)
vvtvvtp (vector times vector plus vector times vector):
Definition: Vmath.hpp:439

References Blas::Dscal(), Nektar::SpatialDomains::eDeformed, Vmath::Fill(), Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::StdRegions::StdExpansion3D::IProductWRTBase_SumFacKernel(), Nektar::StdRegions::StdExpansion::m_base, Nektar::LocalRegions::Expansion::m_metricinfo, Nektar::StdRegions::StdExpansion::m_ncoeffs, Nektar::StdRegions::StdExpansion::MultiplyByQuadratureMetric(), Vmath::Svtsvtp(), Vmath::Svtvp(), Vmath::Vadd(), Vmath::Vvtvp(), and Vmath::Vvtvvtp().

◆ v_MassMatrixOp()

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

Reimplemented from Nektar::StdRegions::StdExpansion.

Definition at line 986 of file PrismExp.cpp.

989{
990 StdExpansion::MassMatrixOp_MatFree(inarray, outarray, mkey);
991}

◆ v_NormalTraceDerivFactors()

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

1556{
1557 int nquad0 = GetNumPoints(0);
1558 int nquad1 = GetNumPoints(1);
1559 int nquad2 = GetNumPoints(2);
1560
1561 const Array<TwoD, const NekDouble> &df =
1562 m_metricinfo->GetDerivFactors(GetPointsKeys());
1563
1564 if (d0factors.size() != 5)
1565 {
1566 d0factors = Array<OneD, Array<OneD, NekDouble>>(5);
1567 d1factors = Array<OneD, Array<OneD, NekDouble>>(5);
1568 d2factors = Array<OneD, Array<OneD, NekDouble>>(5);
1569 }
1570
1571 if (d0factors[0].size() != nquad0 * nquad1)
1572 {
1573 d0factors[0] = Array<OneD, NekDouble>(nquad0 * nquad1);
1574 d1factors[0] = Array<OneD, NekDouble>(nquad0 * nquad1);
1575 d2factors[0] = Array<OneD, NekDouble>(nquad0 * nquad1);
1576 }
1577
1578 if (d0factors[1].size() != nquad0 * nquad2)
1579 {
1580 d0factors[1] = Array<OneD, NekDouble>(nquad0 * nquad2);
1581 d0factors[3] = Array<OneD, NekDouble>(nquad0 * nquad2);
1582 d1factors[1] = Array<OneD, NekDouble>(nquad0 * nquad2);
1583 d1factors[3] = Array<OneD, NekDouble>(nquad0 * nquad2);
1584 d2factors[1] = Array<OneD, NekDouble>(nquad0 * nquad2);
1585 d2factors[3] = Array<OneD, NekDouble>(nquad0 * nquad2);
1586 }
1587
1588 if (d0factors[2].size() != nquad1 * nquad2)
1589 {
1590 d0factors[2] = Array<OneD, NekDouble>(nquad1 * nquad2);
1591 d0factors[4] = Array<OneD, NekDouble>(nquad1 * nquad2);
1592 d1factors[2] = Array<OneD, NekDouble>(nquad1 * nquad2);
1593 d1factors[4] = Array<OneD, NekDouble>(nquad1 * nquad2);
1594 d2factors[2] = Array<OneD, NekDouble>(nquad1 * nquad2);
1595 d2factors[4] = Array<OneD, NekDouble>(nquad1 * nquad2);
1596 }
1597
1598 // Outwards normals
1599 const Array<OneD, const Array<OneD, NekDouble>> &normal_0 =
1600 GetTraceNormal(0);
1601 const Array<OneD, const Array<OneD, NekDouble>> &normal_1 =
1602 GetTraceNormal(1);
1603 const Array<OneD, const Array<OneD, NekDouble>> &normal_2 =
1604 GetTraceNormal(2);
1605 const Array<OneD, const Array<OneD, NekDouble>> &normal_3 =
1606 GetTraceNormal(3);
1607 const Array<OneD, const Array<OneD, NekDouble>> &normal_4 =
1608 GetTraceNormal(4);
1609
1610 int ncoords = normal_0.size();
1611
1612 // first gather together standard cartesian inner products
1613 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
1614 {
1615 // face 0
1616 for (int i = 0; i < nquad0 * nquad1; ++i)
1617 {
1618 d0factors[0][i] = df[0][i] * normal_0[0][i];
1619 d1factors[0][i] = df[1][i] * normal_0[0][i];
1620 d2factors[0][i] = df[2][i] * normal_0[0][i];
1621 }
1622
1623 for (int n = 1; n < ncoords; ++n)
1624 {
1625 for (int i = 0; i < nquad0 * nquad1; ++i)
1626 {
1627 d0factors[0][i] += df[3 * n][i] * normal_0[n][i];
1628 d1factors[0][i] += df[3 * n + 1][i] * normal_0[n][i];
1629 d2factors[0][i] += df[3 * n + 2][i] * normal_0[n][i];
1630 }
1631 }
1632
1633 // faces 1 and 3
1634 for (int j = 0; j < nquad2; ++j)
1635 {
1636 for (int i = 0; i < nquad0; ++i)
1637 {
1638 d0factors[1][j * nquad0 + i] = df[0][j * nquad0 * nquad1 + i] *
1639 normal_1[0][j * nquad0 + i];
1640 d1factors[1][j * nquad0 + i] = df[1][j * nquad0 * nquad1 + i] *
1641 normal_1[0][j * nquad0 + i];
1642 d2factors[1][j * nquad0 + i] = df[2][j * nquad0 * nquad1 + i] *
1643 normal_1[0][j * nquad0 + i];
1644
1645 d0factors[3][j * nquad0 + i] =
1646 df[0][(j + 1) * nquad0 * nquad1 - nquad0 + i] *
1647 normal_3[0][j * nquad0 + i];
1648 d1factors[3][j * nquad0 + i] =
1649 df[1][(j + 1) * nquad0 * nquad1 - nquad0 + i] *
1650 normal_3[0][j * nquad0 + i];
1651 d2factors[3][j * nquad0 + i] =
1652 df[2][(j + 1) * nquad0 * nquad1 - nquad0 + i] *
1653 normal_3[0][j * nquad0 + i];
1654 }
1655 }
1656
1657 for (int n = 1; n < ncoords; ++n)
1658 {
1659 for (int j = 0; j < nquad2; ++j)
1660 {
1661 for (int i = 0; i < nquad0; ++i)
1662 {
1663 d0factors[1][j * nquad0 + i] +=
1664 df[3 * n][j * nquad0 * nquad1 + i] *
1665 normal_1[n][j * nquad0 + i];
1666 d1factors[1][j * nquad0 + i] +=
1667 df[3 * n + 1][j * nquad0 * nquad1 + i] *
1668 normal_1[n][j * nquad0 + i];
1669 d2factors[1][j * nquad0 + i] +=
1670 df[3 * n + 2][j * nquad0 * nquad1 + i] *
1671 normal_1[n][j * nquad0 + i];
1672
1673 d0factors[3][j * nquad0 + i] +=
1674 df[3 * n][(j + 1) * nquad0 * nquad1 - nquad0 + i] *
1675 normal_3[n][j * nquad0 + i];
1676 d1factors[3][j * nquad0 + i] +=
1677 df[3 * n + 1][(j + 1) * nquad0 * nquad1 - nquad0 + i] *
1678 normal_3[n][j * nquad0 + i];
1679 d2factors[3][j * nquad0 + i] +=
1680 df[3 * n + 2][(j + 1) * nquad0 * nquad1 - nquad0 + i] *
1681 normal_3[n][j * nquad0 + i];
1682 }
1683 }
1684 }
1685
1686 // faces 2 and 4
1687 for (int j = 0; j < nquad2; ++j)
1688 {
1689 for (int i = 0; i < nquad1; ++i)
1690 {
1691 d0factors[2][j * nquad1 + i] =
1692 df[0][j * nquad0 * nquad1 + (i + 1) * nquad0 - 1] *
1693 normal_2[0][j * nquad1 + i];
1694 d1factors[2][j * nquad1 + i] =
1695 df[1][j * nquad0 * nquad1 + (i + 1) * nquad0 - 1] *
1696 normal_2[0][j * nquad1 + i];
1697 d2factors[2][j * nquad1 + i] =
1698 df[2][j * nquad0 * nquad1 + (i + 1) * nquad0 - 1] *
1699 normal_2[0][j * nquad1 + i];
1700
1701 d0factors[4][j * nquad1 + i] =
1702 df[0][j * nquad0 * nquad1 + i * nquad0] *
1703 normal_4[0][j * nquad1 + i];
1704 d1factors[4][j * nquad1 + i] =
1705 df[1][j * nquad0 * nquad1 + i * nquad0] *
1706 normal_4[0][j * nquad1 + i];
1707 d2factors[4][j * nquad1 + i] =
1708 df[2][j * nquad0 * nquad1 + i * nquad0] *
1709 normal_4[0][j * nquad1 + i];
1710 }
1711 }
1712
1713 for (int n = 1; n < ncoords; ++n)
1714 {
1715 for (int j = 0; j < nquad2; ++j)
1716 {
1717 for (int i = 0; i < nquad1; ++i)
1718 {
1719 d0factors[2][j * nquad1 + i] +=
1720 df[3 * n][j * nquad0 * nquad1 + (i + 1) * nquad0 - 1] *
1721 normal_2[n][j * nquad1 + i];
1722 d1factors[2][j * nquad1 + i] +=
1723 df[3 * n + 1]
1724 [j * nquad0 * nquad1 + (i + 1) * nquad0 - 1] *
1725 normal_2[n][j * nquad1 + i];
1726 d2factors[2][j * nquad1 + i] +=
1727 df[3 * n + 2]
1728 [j * nquad0 * nquad1 + (i + 1) * nquad0 - 1] *
1729 normal_2[n][j * nquad1 + i];
1730
1731 d0factors[4][j * nquad1 + i] +=
1732 df[3 * n][j * nquad0 * nquad1 + i * nquad0] *
1733 normal_4[n][j * nquad1 + i];
1734 d1factors[4][j * nquad1 + i] +=
1735 df[3 * n + 1][j * nquad0 * nquad1 + i * nquad0] *
1736 normal_4[n][j * nquad1 + i];
1737 d2factors[4][j * nquad1 + i] +=
1738 df[3 * n + 2][j * nquad0 * nquad1 + i * nquad0] *
1739 normal_4[n][j * nquad1 + i];
1740 }
1741 }
1742 }
1743 }
1744 else
1745 {
1746 // Face 0
1747 for (int i = 0; i < nquad0 * nquad1; ++i)
1748 {
1749 d0factors[0][i] = df[0][0] * normal_0[0][i];
1750 d1factors[0][i] = df[1][0] * normal_0[0][i];
1751 d2factors[0][i] = df[2][0] * normal_0[0][i];
1752 }
1753
1754 for (int n = 1; n < ncoords; ++n)
1755 {
1756 for (int i = 0; i < nquad0 * nquad1; ++i)
1757 {
1758 d0factors[0][i] += df[3 * n][0] * normal_0[n][i];
1759 d1factors[0][i] += df[3 * n + 1][0] * normal_0[n][i];
1760 d2factors[0][i] += df[3 * n + 2][0] * normal_0[n][i];
1761 }
1762 }
1763
1764 // faces 1 and 3
1765 for (int i = 0; i < nquad0 * nquad2; ++i)
1766 {
1767 d0factors[1][i] = df[0][0] * normal_1[0][i];
1768 d0factors[3][i] = df[0][0] * normal_3[0][i];
1769
1770 d1factors[1][i] = df[1][0] * normal_1[0][i];
1771 d1factors[3][i] = df[1][0] * normal_3[0][i];
1772
1773 d2factors[1][i] = df[2][0] * normal_1[0][i];
1774 d2factors[3][i] = df[2][0] * normal_3[0][i];
1775 }
1776
1777 for (int n = 1; n < ncoords; ++n)
1778 {
1779 for (int i = 0; i < nquad0 * nquad2; ++i)
1780 {
1781 d0factors[1][i] += df[3 * n][0] * normal_1[n][i];
1782 d0factors[3][i] += df[3 * n][0] * normal_3[n][i];
1783
1784 d1factors[1][i] += df[3 * n + 1][0] * normal_1[n][i];
1785 d1factors[3][i] += df[3 * n + 1][0] * normal_3[n][i];
1786
1787 d2factors[1][i] += df[3 * n + 2][0] * normal_1[n][i];
1788 d2factors[3][i] += df[3 * n + 2][0] * normal_3[n][i];
1789 }
1790 }
1791
1792 // faces 2 and 4
1793 for (int i = 0; i < nquad1 * nquad2; ++i)
1794 {
1795 d0factors[2][i] = df[0][0] * normal_2[0][i];
1796 d0factors[4][i] = df[0][0] * normal_4[0][i];
1797
1798 d1factors[2][i] = df[1][0] * normal_2[0][i];
1799 d1factors[4][i] = df[1][0] * normal_4[0][i];
1800
1801 d2factors[2][i] = df[2][0] * normal_2[0][i];
1802 d2factors[4][i] = df[2][0] * normal_4[0][i];
1803 }
1804
1805 for (int n = 1; n < ncoords; ++n)
1806 {
1807 for (int i = 0; i < nquad1 * nquad2; ++i)
1808 {
1809 d0factors[2][i] += df[3 * n][0] * normal_2[n][i];
1810 d0factors[4][i] += df[3 * n][0] * normal_4[n][i];
1811
1812 d1factors[2][i] += df[3 * n + 1][0] * normal_2[n][i];
1813 d1factors[4][i] += df[3 * n + 1][0] * normal_4[n][i];
1814
1815 d2factors[2][i] += df[3 * n + 2][0] * normal_2[n][i];
1816 d2factors[4][i] += df[3 * n + 2][0] * normal_4[n][i];
1817 }
1818 }
1819 }
1820}
const NormalVector & GetTraceNormal(const int id)
Definition: Expansion.cpp:251

References Nektar::SpatialDomains::eDeformed, Nektar::StdRegions::StdExpansion::GetNumPoints(), Nektar::StdRegions::StdExpansion::GetPointsKeys(), Nektar::LocalRegions::Expansion::GetTraceNormal(), and Nektar::LocalRegions::Expansion::m_metricinfo.

◆ v_PhysDeriv()

void Nektar::LocalRegions::PrismExp::v_PhysDeriv ( const Array< OneD, const NekDouble > &  u_physical,
Array< OneD, NekDouble > &  out_dxi1,
Array< OneD, NekDouble > &  out_dxi2,
Array< OneD, NekDouble > &  out_dxi3 
)
overrideprotectedvirtual

Calculate the derivative of the physical points.

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

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 123 of file PrismExp.cpp.

127{
128 int nqtot = GetTotPoints();
129
130 Array<TwoD, const NekDouble> df =
131 m_metricinfo->GetDerivFactors(GetPointsKeys());
132 Array<OneD, NekDouble> diff0(nqtot);
133 Array<OneD, NekDouble> diff1(nqtot);
134 Array<OneD, NekDouble> diff2(nqtot);
135
136 StdPrismExp::v_PhysDeriv(inarray, diff0, diff1, diff2);
137
138 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
139 {
140 if (out_d0.size())
141 {
142 Vmath::Vmul(nqtot, &df[0][0], 1, &diff0[0], 1, &out_d0[0], 1);
143 Vmath::Vvtvp(nqtot, &df[1][0], 1, &diff1[0], 1, &out_d0[0], 1,
144 &out_d0[0], 1);
145 Vmath::Vvtvp(nqtot, &df[2][0], 1, &diff2[0], 1, &out_d0[0], 1,
146 &out_d0[0], 1);
147 }
148
149 if (out_d1.size())
150 {
151 Vmath::Vmul(nqtot, &df[3][0], 1, &diff0[0], 1, &out_d1[0], 1);
152 Vmath::Vvtvp(nqtot, &df[4][0], 1, &diff1[0], 1, &out_d1[0], 1,
153 &out_d1[0], 1);
154 Vmath::Vvtvp(nqtot, &df[5][0], 1, &diff2[0], 1, &out_d1[0], 1,
155 &out_d1[0], 1);
156 }
157
158 if (out_d2.size())
159 {
160 Vmath::Vmul(nqtot, &df[6][0], 1, &diff0[0], 1, &out_d2[0], 1);
161 Vmath::Vvtvp(nqtot, &df[7][0], 1, &diff1[0], 1, &out_d2[0], 1,
162 &out_d2[0], 1);
163 Vmath::Vvtvp(nqtot, &df[8][0], 1, &diff2[0], 1, &out_d2[0], 1,
164 &out_d2[0], 1);
165 }
166 }
167 else // regular geometry
168 {
169 if (out_d0.size())
170 {
171 Vmath::Smul(nqtot, df[0][0], &diff0[0], 1, &out_d0[0], 1);
172 Blas::Daxpy(nqtot, df[1][0], &diff1[0], 1, &out_d0[0], 1);
173 Blas::Daxpy(nqtot, df[2][0], &diff2[0], 1, &out_d0[0], 1);
174 }
175
176 if (out_d1.size())
177 {
178 Vmath::Smul(nqtot, df[3][0], &diff0[0], 1, &out_d1[0], 1);
179 Blas::Daxpy(nqtot, df[4][0], &diff1[0], 1, &out_d1[0], 1);
180 Blas::Daxpy(nqtot, df[5][0], &diff2[0], 1, &out_d1[0], 1);
181 }
182
183 if (out_d2.size())
184 {
185 Vmath::Smul(nqtot, df[6][0], &diff0[0], 1, &out_d2[0], 1);
186 Blas::Daxpy(nqtot, df[7][0], &diff1[0], 1, &out_d2[0], 1);
187 Blas::Daxpy(nqtot, df[8][0], &diff2[0], 1, &out_d2[0], 1);
188 }
189 }
190}
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
Definition: StdExpansion.h:134
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::GetTotPoints(), Nektar::LocalRegions::Expansion::m_metricinfo, Vmath::Smul(), Vmath::Vmul(), and Vmath::Vvtvp().

◆ v_PhysEvalFirstDeriv()

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 531 of file PrismExp.cpp.

535{
536 Array<OneD, NekDouble> Lcoord(3);
537 ASSERTL0(m_geom, "m_geom not defined");
538 m_geom->GetLocCoords(coord, Lcoord);
539 return StdPrismExp::v_PhysEvalFirstDeriv(Lcoord, inarray, firstOrderDerivs);
540}

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

◆ v_PhysEvaluate()

NekDouble Nektar::LocalRegions::PrismExp::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.

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 #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::StdExpansion3D.

Definition at line 519 of file PrismExp.cpp.

521{
522 Array<OneD, NekDouble> Lcoord(3);
523
524 ASSERTL0(m_geom, "m_geom not defined");
525
526 m_geom->GetLocCoords(coord, Lcoord);
527
528 return StdExpansion3D::v_PhysEvaluate(Lcoord, physvals);
529}

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

◆ v_StdPhysEvaluate()

NekDouble Nektar::LocalRegions::PrismExp::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 511 of file PrismExp.cpp.

514{
515 // Evaluate point in local (eta) coordinates.
516 return StdExpansion3D::v_PhysEvaluate(Lcoord, physvals);
517}

◆ v_SVVLaplacianFilter()

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

Reimplemented from Nektar::StdRegions::StdPrismExp.

Definition at line 1015 of file PrismExp.cpp.

1017{
1018 int nq = GetTotPoints();
1019
1020 // Calculate sqrt of the Jacobian
1021 Array<OneD, const NekDouble> jac = m_metricinfo->GetJac(GetPointsKeys());
1022 Array<OneD, NekDouble> sqrt_jac(nq);
1023 if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
1024 {
1025 Vmath::Vsqrt(nq, jac, 1, sqrt_jac, 1);
1026 }
1027 else
1028 {
1029 Vmath::Fill(nq, sqrt(jac[0]), sqrt_jac, 1);
1030 }
1031
1032 // Multiply array by sqrt(Jac)
1033 Vmath::Vmul(nq, sqrt_jac, 1, array, 1, array, 1);
1034
1035 // Apply std region filter
1036 StdPrismExp::v_SVVLaplacianFilter(array, mkey);
1037
1038 // Divide by sqrt(Jac)
1039 Vmath::Vdiv(nq, array, 1, sqrt_jac, 1, array, 1);
1040}
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
Definition: Vmath.hpp:126

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

Member Data Documentation

◆ m_matrixManager

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

Definition at line 198 of file PrismExp.h.

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

◆ m_staticCondMatrixManager

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

Definition at line 200 of file PrismExp.h.

Referenced by v_DropLocStaticCondMatrix(), and v_GetLocStaticCondMatrix().