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

Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expansions. More...

#include <ExpList2D.h>

Inheritance diagram for Nektar::MultiRegions::ExpList2D:
Inheritance graph
[legend]
Collaboration diagram for Nektar::MultiRegions::ExpList2D:
Collaboration graph
[legend]

Public Member Functions

 ExpList2D ()
 Default constructor. More...
 
 ExpList2D (const ExpList2D &In, const bool DeclareCoeffPhysArrays=true)
 Copy constructor. More...
 
 ExpList2D (const ExpList2D &In, const std::vector< unsigned int > &eIDs, const bool DeclareCoeffPhysArrays=true)
 Constructor copying only elements defined in eIds. More...
 
 ExpList2D (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &graph2D, const bool DelcareCoeffPhysArrays=true, const std::string &var="DefaultVar")
 Sets up a list of local expansions based on an input mesh. More...
 
 ExpList2D (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::ExpansionMap &expansions, const bool DeclareCoeffPhysArrays=true)
 Sets up a list of local expansions based on an expansion Map. More...
 
 ExpList2D (const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::BasisKey &TriBa, const LibUtilities::BasisKey &TriBb, const LibUtilities::BasisKey &QuadBa, const LibUtilities::BasisKey &QuadBb, const SpatialDomains::MeshGraphSharedPtr &graph2D, const LibUtilities::PointsType TriNb=LibUtilities::SIZE_PointsType)
 Sets up a list of local expansions based on an input mesh and separately defined basiskeys. More...
 
 ExpList2D (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, const ExpListSharedPtr > &bndConstraint, const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > &bndCond, const LocalRegions::ExpansionVector &locexp, const SpatialDomains::MeshGraphSharedPtr &graph3D, const PeriodicMap &periodicFaces, const bool DeclareCoeffPhysArrays=true, const std::string variable="DefaultVar")
 
 ExpList2D (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::CompositeMap &domain, const SpatialDomains::MeshGraphSharedPtr &graph3D, const std::string variable="DefaultVar")
 Specialised constructor for Neumann boundary conditions in DisContField3D and ContField3D. More...
 
virtual ~ExpList2D ()
 Destructor. More...
 
- Public Member Functions inherited from Nektar::MultiRegions::ExpList
 ExpList ()
 The default constructor. More...
 
 ExpList (const LibUtilities::SessionReaderSharedPtr &pSession)
 The default constructor. More...
 
 ExpList (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 The default constructor. More...
 
 ExpList (const ExpList &in, const std::vector< unsigned int > &eIDs, const bool DeclareCoeffPhysArrays=true)
 Constructor copying only elements defined in eIds. More...
 
 ExpList (const ExpList &in, const bool DeclareCoeffPhysArrays=true)
 The copy constructor. More...
 
virtual ~ExpList ()
 The default destructor. More...
 
int GetNcoeffs (void) const
 Returns the total number of local degrees of freedom $N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_m$. More...
 
int GetNcoeffs (const int eid) const
 Returns the total number of local degrees of freedom for element eid. More...
 
ExpansionType GetExpType (void)
 Returns the type of the expansion. More...
 
void SetExpType (ExpansionType Type)
 Returns the type of the expansion. More...
 
int EvalBasisNumModesMax (void) const
 Evaulates the maximum number of modes in the elemental basis order over all elements. More...
 
const Array< OneD, int > EvalBasisNumModesMaxPerExp (void) const
 Returns the vector of the number of modes in the elemental basis order over all elements. More...
 
int GetTotPoints (void) const
 Returns the total number of quadrature points m_npoints $=Q_{\mathrm{tot}}$. More...
 
int GetTotPoints (const int eid) const
 Returns the total number of quadrature points for eid's element $=Q_{\mathrm{tot}}$. More...
 
int GetNpoints (void) const
 Returns the total number of quadrature points m_npoints $=Q_{\mathrm{tot}}$. More...
 
int Get1DScaledTotPoints (const NekDouble scale) const
 Returns the total number of qudature points scaled by the factor scale on each 1D direction. More...
 
void SetWaveSpace (const bool wavespace)
 Sets the wave space to the one of the possible configuration true or false. More...
 
void SetModifiedBasis (const bool modbasis)
 Set Modified Basis for the stability analysis. More...
 
void SetPhys (int i, NekDouble val)
 Set the i th value of m_phys to value val. More...
 
bool GetWaveSpace (void) const
 This function returns the third direction expansion condition, which can be in wave space (coefficient) or not It is stored in the variable m_WaveSpace. More...
 
void SetPhys (const Array< OneD, const NekDouble > &inarray)
 Fills the array m_phys. More...
 
void SetPhysArray (Array< OneD, NekDouble > &inarray)
 Sets the array m_phys. More...
 
void SetPhysState (const bool physState)
 This function manually sets whether the array of physical values $\boldsymbol{u}_l$ (implemented as m_phys) is filled or not. More...
 
bool GetPhysState (void) const
 This function indicates whether the array of physical values $\boldsymbol{u}_l$ (implemented as m_phys) is filled or not. More...
 
NekDouble PhysIntegral (void)
 This function integrates a function $f(\boldsymbol{x})$ over the domain consisting of all the elements of the expansion. More...
 
NekDouble PhysIntegral (const Array< OneD, const NekDouble > &inarray)
 This function integrates a function $f(\boldsymbol{x})$ over the domain consisting of all the elements of the expansion. More...
 
void IProductWRTBase_IterPerExp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function calculates the inner product of a function $f(\boldsymbol{x})$ with respect to all {local} expansion modes $\phi_n^e(\boldsymbol{x})$. More...
 
void IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
 
void IProductWRTDerivBase (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function calculates the inner product of a function $f(\boldsymbol{x})$ with respect to the derivative (in direction. More...
 
void IProductWRTDerivBase (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &outarray)
 This function calculates the inner product of a function $f(\boldsymbol{x})$ with respect to the derivative (in direction. More...
 
void FwdTrans_IterPerExp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function elementally evaluates the forward transformation of a function $u(\boldsymbol{x})$ onto the global spectral/hp expansion. More...
 
void FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
 
void MultiplyByElmtInvMass (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function elementally mulplies the coefficient space of Sin my the elemental inverse of the mass matrix. More...
 
void MultiplyByInvMassMatrix (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
 
void SmoothField (Array< OneD, NekDouble > &field)
 Smooth a field across elements. More...
 
void HelmSolve (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const FlagList &flags, const StdRegions::ConstFactorMap &factors, const StdRegions::VarCoeffMap &varcoeff=StdRegions::NullVarCoeffMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
 Solve helmholtz problem. More...
 
void LinearAdvectionDiffusionReactionSolve (const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
 Solve Advection Diffusion Reaction. More...
 
void LinearAdvectionReactionSolve (const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
 Solve Advection Diffusion Reaction. More...
 
void FwdTrans_BndConstrained (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void BwdTrans_IterPerExp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function elementally evaluates the backward transformation of the global spectral/hp element expansion. More...
 
void BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
 
void GetCoords (Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1=NullNekDouble1DArray, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
 This function calculates the coordinates of all the elemental quadrature points $\boldsymbol{x}_i$. More...
 
void HomogeneousFwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
 
void HomogeneousBwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
 
void DealiasedProd (const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
 
void GetBCValues (Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
 
void NormVectorIProductWRTBase (Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
 
void NormVectorIProductWRTBase (Array< OneD, Array< OneD, NekDouble > > &V, Array< OneD, NekDouble > &outarray)
 
void ApplyGeomInfo ()
 Apply geometry information to each expansion. More...
 
void Reset ()
 Reset geometry information and reset matrices. More...
 
void WriteTecplotHeader (std::ostream &outfile, std::string var="")
 
void WriteTecplotZone (std::ostream &outfile, int expansion=-1)
 
void WriteTecplotField (std::ostream &outfile, int expansion=-1)
 
void WriteTecplotConnectivity (std::ostream &outfile, int expansion=-1)
 
void WriteVtkHeader (std::ostream &outfile)
 
void WriteVtkFooter (std::ostream &outfile)
 
void WriteVtkPieceHeader (std::ostream &outfile, int expansion, int istrip=0)
 
void WriteVtkPieceFooter (std::ostream &outfile, int expansion)
 
void WriteVtkPieceData (std::ostream &outfile, int expansion, std::string var="v")
 
int GetCoordim (int eid)
 This function returns the dimension of the coordinates of the element eid. More...
 
void SetCoeff (int i, NekDouble val)
 Set the i th coefficiient in m_coeffs to value val. More...
 
void SetCoeffs (int i, NekDouble val)
 Set the i th coefficiient in m_coeffs to value val. More...
 
void SetCoeffsArray (Array< OneD, NekDouble > &inarray)
 Set the m_coeffs array to inarray. More...
 
const Array< OneD, const
NekDouble > & 
GetCoeffs () const
 This function returns (a reference to) the array $\boldsymbol{\hat{u}}_l$ (implemented as m_coeffs) containing all local expansion coefficients. More...
 
void ImposeDirichletConditions (Array< OneD, NekDouble > &outarray)
 Impose Dirichlet Boundary Conditions onto Array. More...
 
void FillBndCondFromField (void)
 Fill Bnd Condition expansion from the values stored in expansion. More...
 
void LocalToGlobal (void)
 Put the coefficients into global ordering using m_coeffs. More...
 
void GlobalToLocal (void)
 Put the coefficients into local ordering and place in m_coeffs. More...
 
NekDouble GetCoeff (int i)
 Get the i th value (coefficient) of m_coeffs. More...
 
NekDouble GetCoeffs (int i)
 Get the i th value (coefficient) of m_coeffs. More...
 
const Array< OneD, const
NekDouble > & 
GetPhys () const
 This function returns (a reference to) the array $\boldsymbol{u}_l$ (implemented as m_phys) containing the function $u^{\delta}(\boldsymbol{x})$ evaluated at the quadrature points. More...
 
NekDouble Linf (const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
 This function calculates the $L_\infty$ error of the global spectral/hp element approximation. More...
 
NekDouble L2 (const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
 This function calculates the $L_2$ error with respect to soln of the global spectral/hp element approximation. More...
 
NekDouble H1 (const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
 Calculates the $H^1$ error of the global spectral/hp element approximation. More...
 
NekDouble Integral (const Array< OneD, const NekDouble > &inarray)
 
Array< OneD, const NekDoubleHomogeneousEnergy (void)
 This function calculates the energy associated with each one of the modesof a 3D homogeneous nD expansion. More...
 
void SetHomo1DSpecVanVisc (Array< OneD, NekDouble > visc)
 This function sets the Spectral Vanishing Viscosity in homogeneous1D expansion. More...
 
Array< OneD, const unsigned int > GetZIDs (void)
 This function returns a vector containing the wave numbers in z-direction associated with the 3D homogenous expansion. Required if a parellelisation is applied in the Fourier direction. More...
 
LibUtilities::TranspositionSharedPtr GetTransposition (void)
 This function returns the transposition class associaed with the homogeneous expansion. More...
 
NekDouble GetHomoLen (void)
 This function returns the Width of homogeneous direction associaed with the homogeneous expansion. More...
 
Array< OneD, const unsigned int > GetYIDs (void)
 This function returns a vector containing the wave numbers in y-direction associated with the 3D homogenous expansion. Required if a parellelisation is applied in the Fourier direction. More...
 
void PhysInterp1DScaled (const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function interpolates the physical space points in inarray to outarray using the same points defined in the expansion but where the number of points are rescaled by 1DScale. More...
 
void PhysGalerkinProjection1DScaled (const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 This function Galerkin projects the physical space points in inarray to outarray where inarray is assumed to be defined in the expansion but where the number of points are rescaled by 1DScale. More...
 
int GetExpSize (void)
 This function returns the number of elements in the expansion. More...
 
int GetNumElmts (void)
 This function returns the number of elements in the expansion which may be different for a homogeoenous extended expansionp. More...
 
const boost::shared_ptr
< LocalRegions::ExpansionVector
GetExp () const
 This function returns the vector of elements in the expansion. More...
 
LocalRegions::ExpansionSharedPtrGetExp (int n) const
 This function returns (a shared pointer to) the local elemental expansion of the $n^{\mathrm{th}}$ element. More...
 
LocalRegions::ExpansionSharedPtrGetExp (const Array< OneD, const NekDouble > &gloCoord)
 This function returns (a shared pointer to) the local elemental expansion containing the arbitrary point given by gloCoord. More...
 
int GetExpIndex (const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0, bool returnNearestElmt=false)
 
int GetExpIndex (const Array< OneD, const NekDouble > &gloCoords, Array< OneD, NekDouble > &locCoords, NekDouble tol=0.0, bool returnNearestElmt=false)
 
int GetCoeff_Offset (int n) const
 Get the start offset position for a global list of m_coeffs correspoinding to element n. More...
 
int GetPhys_Offset (int n) const
 Get the start offset position for a global list of m_phys correspoinding to element n. More...
 
int GetOffset_Elmt_Id (int n) const
 Get the element id associated with the n th consecutive block of data in m_phys and m_coeffs. More...
 
Array< OneD, NekDouble > & UpdateCoeffs ()
 This function returns (a reference to) the array $\boldsymbol{\hat{u}}_l$ (implemented as m_coeffs) containing all local expansion coefficients. More...
 
Array< OneD, NekDouble > & UpdatePhys ()
 This function returns (a reference to) the array $\boldsymbol{u}_l$ (implemented as m_phys) containing the function $u^{\delta}(\boldsymbol{x})$ evaluated at the quadrature points. More...
 
void PhysDeriv (Direction edir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d)
 
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)
 This function discretely evaluates the derivative of a function $f(\boldsymbol{x})$ on the domain consisting of all elements of the expansion. More...
 
void PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d)
 
const Array< OneD, const
boost::shared_ptr< ExpList > > & 
GetBndCondExpansions ()
 
boost::shared_ptr< ExpList > & UpdateBndCondExpansion (int i)
 
void Upwind (const Array< OneD, const Array< OneD, NekDouble > > &Vec, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &Upwind)
 
void Upwind (const Array< OneD, const NekDouble > &Vn, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &Upwind)
 
boost::shared_ptr< ExpList > & GetTrace ()
 
boost::shared_ptr
< AssemblyMapDG > & 
GetTraceMap (void)
 
const Array< OneD, const int > & GetTraceBndMap (void)
 
void GetNormals (Array< OneD, Array< OneD, NekDouble > > &normals)
 
void AddTraceIntegral (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
void AddTraceIntegral (const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
void AddFwdBwdTraceIntegral (const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &outarray)
 
void GetFwdBwdTracePhys (Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
 
void GetFwdBwdTracePhys (const Array< OneD, const NekDouble > &field, Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
 
const vector< bool > & GetLeftAdjacentFaces (void) const
 
void ExtractTracePhys (Array< OneD, NekDouble > &outarray)
 
void ExtractTracePhys (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
const Array< OneD, const
SpatialDomains::BoundaryConditionShPtr > & 
GetBndConditions ()
 
Array< OneD,
SpatialDomains::BoundaryConditionShPtr > & 
UpdateBndConditions ()
 
void EvaluateBoundaryConditions (const NekDouble time=0.0, const std::string varName="", const NekDouble=NekConstants::kNekUnsetDouble, const NekDouble=NekConstants::kNekUnsetDouble)
 
void GeneralMatrixOp (const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
 This function calculates the result of the multiplication of a matrix of type specified by mkey with a vector given by inarray. More...
 
void GeneralMatrixOp_IterPerExp (const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
void SetUpPhysNormals ()
 
void GetBoundaryToElmtMap (Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
 
void GetBndElmtExpansion (int i, boost::shared_ptr< ExpList > &result)
 
void ExtractElmtToBndPhys (int i, Array< OneD, NekDouble > &elmt, Array< OneD, NekDouble > &boundary)
 
void ExtractPhysToBndElmt (int i, const Array< OneD, const NekDouble > &phys, Array< OneD, NekDouble > &bndElmt)
 
void GetBoundaryNormals (int i, Array< OneD, Array< OneD, NekDouble > > &normals)
 
void GeneralGetFieldDefinitions (std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef, int NumHomoDir=0, Array< OneD, LibUtilities::BasisSharedPtr > &HomoBasis=LibUtilities::NullBasisSharedPtr1DArray, std::vector< NekDouble > &HomoLen=LibUtilities::NullNekDoubleVector, bool homoStrips=false, std::vector< unsigned int > &HomoSIDs=LibUtilities::NullUnsignedIntVector, std::vector< unsigned int > &HomoZIDs=LibUtilities::NullUnsignedIntVector, std::vector< unsigned int > &HomoYIDs=LibUtilities::NullUnsignedIntVector)
 
const
NekOptimize::GlobalOptParamSharedPtr
GetGlobalOptParam (void)
 
map< int, RobinBCInfoSharedPtrGetRobinBCInfo ()
 
void GetPeriodicEntities (PeriodicMap &periodicVerts, PeriodicMap &periodicEdges, PeriodicMap &periodicFaces=NullPeriodicMap)
 
std::vector
< LibUtilities::FieldDefinitionsSharedPtr
GetFieldDefinitions ()
 
void GetFieldDefinitions (std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef)
 
void AppendFieldData (LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
 Append the element data listed in elements fielddef->m_ElementIDs onto fielddata. More...
 
void AppendFieldData (LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, Array< OneD, NekDouble > &coeffs)
 Append the data in coeffs listed in elements fielddef->m_ElementIDs onto fielddata. More...
 
void ExtractElmtDataToCoeffs (LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
 Extract the data in fielddata into the coeffs using the basic ExpList Elemental expansions rather than planes in homogeneous case. More...
 
void ExtractCoeffsToCoeffs (const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
 Extract the data from fromField using fromExpList the coeffs using the basic ExpList Elemental expansions rather than planes in homogeneous case. More...
 
void ExtractDataToCoeffs (LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
 Extract the data in fielddata into the coeffs. More...
 
boost::shared_ptr< ExpListGetSharedThisPtr ()
 Returns a shared pointer to the current object. More...
 
boost::shared_ptr
< LibUtilities::SessionReader
GetSession ()
 Returns the session object. More...
 
boost::shared_ptr
< LibUtilities::Comm
GetComm ()
 Returns the comm object. More...
 
SpatialDomains::MeshGraphSharedPtr GetGraph ()
 
LibUtilities::BasisSharedPtr GetHomogeneousBasis (void)
 
boost::shared_ptr< ExpList > & GetPlane (int n)
 
void CreateCollections (Collections::ImplementationType ImpType=Collections::eNoImpType)
 Construct collections of elements containing a single element type and polynomial order from the list of expansions. More...
 
void ClearGlobalLinSysManager (void)
 

Protected Member Functions

void v_Upwind (const Array< OneD, const NekDouble > &Vn, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &Upwind)
 Upwind the Fwd and Bwd states based on the one- dimensional normal velocity field given by Vn. More...
 
void v_GetNormals (Array< OneD, Array< OneD, NekDouble > > &normals)
 For each local element, copy the normals stored in the element list into the array normals. More...
 
- Protected Member Functions inherited from Nektar::MultiRegions::ExpList
boost::shared_ptr< DNekMatGenGlobalMatrixFull (const GlobalLinSysKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
 
const DNekScalBlkMatSharedPtr GenBlockMatrix (const GlobalMatrixKey &gkey)
 This function assembles the block diagonal matrix of local matrices of the type mtype. More...
 
const DNekScalBlkMatSharedPtrGetBlockMatrix (const GlobalMatrixKey &gkey)
 
void MultiplyByBlockMatrix (const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
boost::shared_ptr< GlobalMatrixGenGlobalMatrix (const GlobalMatrixKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
 Generates a global matrix from the given key and map. More...
 
void GlobalEigenSystem (const boost::shared_ptr< DNekMat > &Gmat, Array< OneD, NekDouble > &EigValsReal, Array< OneD, NekDouble > &EigValsImag, Array< OneD, NekDouble > &EigVecs=NullNekDouble1DArray)
 
boost::shared_ptr< GlobalLinSysGenGlobalLinSys (const GlobalLinSysKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
 This operation constructs the global linear system of type mkey. More...
 
boost::shared_ptr< GlobalLinSysGenGlobalBndLinSys (const GlobalLinSysKey &mkey, const AssemblyMapSharedPtr &locToGloMap)
 Generate a GlobalLinSys from information provided by the key "mkey" and the mapping provided in LocToGloBaseMap. More...
 
void ReadGlobalOptimizationParameters ()
 
virtual int v_GetNumElmts (void)
 
virtual const Array< OneD,
const boost::shared_ptr
< ExpList > > & 
v_GetBndCondExpansions (void)
 
virtual boost::shared_ptr
< ExpList > & 
v_UpdateBndCondExpansion (int i)
 
virtual void v_Upwind (const Array< OneD, const Array< OneD, NekDouble > > &Vec, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &Upwind)
 
virtual boost::shared_ptr
< ExpList > & 
v_GetTrace ()
 
virtual boost::shared_ptr
< AssemblyMapDG > & 
v_GetTraceMap ()
 
virtual const Array< OneD,
const int > & 
v_GetTraceBndMap ()
 
virtual void v_AddTraceIntegral (const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddTraceIntegral (const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
 
virtual void v_AddFwdBwdTraceIntegral (const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &outarray)
 
virtual void v_GetFwdBwdTracePhys (Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
 
virtual void v_GetFwdBwdTracePhys (const Array< OneD, const NekDouble > &field, Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
 
virtual const vector< bool > & v_GetLeftAdjacentFaces (void) const
 
virtual void v_ExtractTracePhys (Array< OneD, NekDouble > &outarray)
 
virtual void v_ExtractTracePhys (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_MultiplyByInvMassMatrix (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
 
virtual void v_HelmSolve (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const FlagList &flags, const StdRegions::ConstFactorMap &factors, const StdRegions::VarCoeffMap &varcoeff, const Array< OneD, const NekDouble > &dirForcing)
 
virtual void v_LinearAdvectionDiffusionReactionSolve (const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
 
virtual void v_LinearAdvectionReactionSolve (const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
 
virtual void v_ImposeDirichletConditions (Array< OneD, NekDouble > &outarray)
 
virtual void v_FillBndCondFromField ()
 
virtual void v_Reset ()
 Reset geometry information, metrics, matrix managers and geometry information. More...
 
virtual void v_LocalToGlobal (void)
 
virtual void v_GlobalToLocal (void)
 
virtual void v_BwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
 
virtual void v_BwdTrans_IterPerExp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_FwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
 
virtual void v_FwdTrans_IterPerExp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_SmoothField (Array< OneD, NekDouble > &field)
 
virtual void v_IProductWRTBase (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
 
virtual void v_IProductWRTBase_IterPerExp (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_GeneralMatrixOp (const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
 
virtual void v_GetCoords (Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
 
virtual void v_PhysDeriv (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
 
virtual void v_PhysDeriv (const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d)
 
virtual void v_PhysDeriv (Direction edir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d)
 
virtual void v_HomogeneousFwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
 
virtual void v_HomogeneousBwdTrans (const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
 
virtual void v_DealiasedProd (const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
 
virtual void v_GetBCValues (Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
 
virtual void v_NormVectorIProductWRTBase (Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
 
virtual void v_NormVectorIProductWRTBase (Array< OneD, Array< OneD, NekDouble > > &V, Array< OneD, NekDouble > &outarray)
 
virtual void v_GetBoundaryToElmtMap (Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
 
virtual void v_GetBndElmtExpansion (int i, boost::shared_ptr< ExpList > &result)
 
virtual void v_ExtractElmtToBndPhys (int i, Array< OneD, NekDouble > &elmt, Array< OneD, NekDouble > &boundary)
 
virtual void v_ExtractPhysToBndElmt (int i, const Array< OneD, const NekDouble > &phys, Array< OneD, NekDouble > &bndElmt)
 
virtual void v_GetBoundaryNormals (int i, Array< OneD, Array< OneD, NekDouble > > &normals)
 
virtual std::vector
< LibUtilities::FieldDefinitionsSharedPtr
v_GetFieldDefinitions (void)
 
virtual void v_GetFieldDefinitions (std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef)
 
virtual void v_AppendFieldData (LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
 
virtual void v_AppendFieldData (LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, Array< OneD, NekDouble > &coeffs)
 
virtual void v_ExtractDataToCoeffs (LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
 Extract data from raw field data into expansion list. More...
 
virtual void v_ExtractCoeffsToCoeffs (const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
 
virtual void v_WriteTecplotHeader (std::ostream &outfile, std::string var="")
 
virtual void v_WriteTecplotZone (std::ostream &outfile, int expansion)
 
virtual void v_WriteTecplotField (std::ostream &outfile, int expansion)
 
virtual void v_WriteTecplotConnectivity (std::ostream &outfile, int expansion)
 
virtual void v_WriteVtkPieceData (std::ostream &outfile, int expansion, std::string var)
 
virtual NekDouble v_L2 (const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
 
virtual NekDouble v_Integral (const Array< OneD, const NekDouble > &inarray)
 
virtual Array< OneD, const
NekDouble
v_HomogeneousEnergy (void)
 
virtual
LibUtilities::TranspositionSharedPtr 
v_GetTransposition (void)
 
virtual NekDouble v_GetHomoLen (void)
 
virtual Array< OneD, const
unsigned int > 
v_GetZIDs (void)
 
virtual Array< OneD, const
unsigned int > 
v_GetYIDs (void)
 
virtual void v_ClearGlobalLinSysManager (void)
 
void ExtractFileBCs (const std::string &fileName, const std::string &varName, const boost::shared_ptr< ExpList > locExpList)
 

Private Member Functions

void SetCoeffPhysOffsets (void)
 Definition of the total number of degrees of freedom and quadrature points and offsets to access datax. More...
 
virtual void v_SetUpPhysNormals ()
 Set up the normals on each expansion. More...
 
virtual void v_ReadGlobalOptimizationParameters ()
 
virtual void v_WriteVtkPieceHeader (std::ostream &outfile, int expansion, int istrip)
 
virtual void v_PhysInterp1DScaled (const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 
virtual void v_PhysGalerkinProjection1DScaled (const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
 

Additional Inherited Members

- Public Attributes inherited from Nektar::MultiRegions::ExpList
ExpansionType m_expType
 
- Static Protected Member Functions inherited from Nektar::MultiRegions::ExpList
static
SpatialDomains::BoundaryConditionShPtr 
GetBoundaryCondition (const SpatialDomains::BoundaryConditionCollection &collection, unsigned int index, const std::string &variable)
 
- Protected Attributes inherited from Nektar::MultiRegions::ExpList
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
LibUtilities::SessionReaderSharedPtr m_session
 Session. More...
 
SpatialDomains::MeshGraphSharedPtr m_graph
 Mesh associated with this expansion list. More...
 
int m_ncoeffs
 The total number of local degrees of freedom. m_ncoeffs $=N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_l$. More...
 
int m_npoints
 
Array< OneD, NekDoublem_coeffs
 Concatenation of all local expansion coefficients. More...
 
Array< OneD, NekDoublem_phys
 The global expansion evaluated at the quadrature points. More...
 
bool m_physState
 The state of the array m_phys. More...
 
boost::shared_ptr
< LocalRegions::ExpansionVector
m_exp
 The list of local expansions. More...
 
Collections::CollectionVector m_collections
 
std::vector< int > m_coll_coeff_offset
 Offset of elemental data into the array m_coeffs. More...
 
std::vector< int > m_coll_phys_offset
 Offset of elemental data into the array m_phys. More...
 
Array< OneD, int > m_coeff_offset
 Offset of elemental data into the array m_coeffs. More...
 
Array< OneD, int > m_phys_offset
 Offset of elemental data into the array m_phys. More...
 
Array< OneD, int > m_offset_elmt_id
 Array containing the element id m_offset_elmt_id[n] that the n^th consecutive block of data in m_coeffs and m_phys is associated, i.e. for an array of constant expansion size and single shape elements m_phys[n*m_npoints] is the data related to m_exp[m_offset_elmt_id[n]];. More...
 
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
 
BlockMatrixMapShPtr m_blockMat
 
bool m_WaveSpace
 

Detailed Description

Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expansions.

This multi-elemental expansion, which does not exhibit any coupling between the expansion on the separate elements, can be formulated as,

\[u^{\delta}(\boldsymbol{x}_i)=\sum_{e=1}^{{N_{\mathrm{el}}}} \sum_{n=0}^{N^{e}_m-1}\hat{u}_n^e\phi_n^e(\boldsymbol{x}_i).\]

where ${N_{\mathrm{el}}}$ is the number of elements and $N^{e}_m$ is the local elemental number of expansion modes. This class inherits all its variables and member functions from the base class ExpList.

Definition at line 60 of file ExpList2D.h.

Constructor & Destructor Documentation

Nektar::MultiRegions::ExpList2D::ExpList2D ( )

Default constructor.

Definition at line 68 of file ExpList2D.cpp.

References Nektar::MultiRegions::e2D, and Nektar::MultiRegions::ExpList::SetExpType().

68  :
69  ExpList()
70  {
71  SetExpType(e2D);
72  }
ExpList()
The default constructor.
Definition: ExpList.cpp:93
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
Nektar::MultiRegions::ExpList2D::ExpList2D ( const ExpList2D In,
const bool  DeclareCoeffPhysArrays = true 
)

Copy constructor.

Parameters
InExpList2D object to copy.

Definition at line 86 of file ExpList2D.cpp.

References Nektar::MultiRegions::e2D, and Nektar::MultiRegions::ExpList::SetExpType().

88  :
89  ExpList(In,DeclareCoeffPhysArrays)
90  {
91  SetExpType(e2D);
92  }
ExpList()
The default constructor.
Definition: ExpList.cpp:93
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
Nektar::MultiRegions::ExpList2D::ExpList2D ( const ExpList2D In,
const std::vector< unsigned int > &  eIDs,
const bool  DeclareCoeffPhysArrays = true 
)

Constructor copying only elements defined in eIds.

Parameters
InExpList2D object to copy.
eIDsId of elements that should be copied.

Definition at line 98 of file ExpList2D.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::MultiRegions::ExpList::CreateCollections(), Nektar::MultiRegions::e2D, Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::MultiRegions::ExpList::m_globalOptParam, Nektar::MultiRegions::ExpList::ReadGlobalOptimizationParameters(), SetCoeffPhysOffsets(), and Nektar::MultiRegions::ExpList::SetExpType().

101  :
102  ExpList(In,eIDs,DeclareCoeffPhysArrays)
103  {
104  SetExpType(e2D);
105 
106  // Setup Default optimisation information.
107  int nel = GetExpSize();
110 
111  // set up offset arrays.
113 
116  }
ExpList()
The default constructor.
Definition: ExpList.cpp:93
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
void SetCoeffPhysOffsets(void)
Definition of the total number of degrees of freedom and quadrature points and offsets to access data...
Definition: ExpList2D.cpp:1120
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2950
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
Nektar::MultiRegions::ExpList2D::ExpList2D ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr graph2D,
const bool  DeclareCoeffPhysArrays = true,
const std::string &  var = "DefaultVar" 
)

Sets up a list of local expansions based on an input mesh.

Given a mesh graph2D, containing information about the domain and the spectral/hp element expansion, this constructor fills the list of local expansions {m_exp} with the proper expansions, calculates the total number of quadrature points $\boldsymbol{x}_i$ and local expansion coefficients $\hat{u}^e_n$ and allocates memory for the arrays m_coeffs and m_phys.

Parameters
graph2DA mesh, containing information about the domain and the spectral/hp element expansion.

Definition at line 131 of file ExpList2D.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::MultiRegions::ExpList::CreateCollections(), Nektar::MultiRegions::e2D, Nektar::LibUtilities::eGLL_Lagrange, Nektar::LibUtilities::eNodalTriElec, Nektar::LibUtilities::eOrtho_A, Nektar::LibUtilities::BasisKey::GetBasisType(), Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::LibUtilities::BasisKey::GetNumModes(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::LibUtilities::BasisKey::GetPointsKey(), Nektar::MultiRegions::ExpList::m_coeffs, Nektar::MultiRegions::ExpList::m_globalOptParam, Nektar::MultiRegions::ExpList::m_ncoeffs, Nektar::MultiRegions::ExpList::m_npoints, Nektar::MultiRegions::ExpList::m_phys, Nektar::MultiRegions::ExpList::ReadGlobalOptimizationParameters(), SetCoeffPhysOffsets(), and Nektar::MultiRegions::ExpList::SetExpType().

135  :
136  ExpList(pSession,graph2D)
137  {
138  SetExpType(e2D);
139 
140  int elmtid=0;
146 
147  const SpatialDomains::ExpansionMap &expansions
148  = graph2D->GetExpansions(var);
149 
150  SpatialDomains::ExpansionMap::const_iterator expIt;
151  for (expIt = expansions.begin(); expIt != expansions.end(); ++expIt)
152  {
154  SpatialDomains::QuadGeomSharedPtr QuadrilateralGeom;
155 
156  if ((TriangleGeom = boost::dynamic_pointer_cast<SpatialDomains
157  ::TriGeom>(expIt->second->m_geomShPtr)))
158  {
160  = expIt->second->m_basisKeyVector[0];
162  = expIt->second->m_basisKeyVector[1];
163 
164  // This is not elegantly implemented needs re-thinking.
166  {
168  TriBa.GetNumModes(),
169  TriBa.GetPointsKey());
170 
173  ::AllocateSharedPtr(newBa,TriBb,TriNb,
174  TriangleGeom);
175  Ntri->SetElmtId(elmtid++);
176  (*m_exp).push_back(Ntri);
177  }
178  else
179  {
181  ::AllocateSharedPtr(TriBa,TriBb,
182  TriangleGeom);
183  tri->SetElmtId(elmtid++);
184  (*m_exp).push_back(tri);
185  }
186  m_ncoeffs += (TriBa.GetNumModes()*(TriBa.GetNumModes()+1))/2
187  + TriBa.GetNumModes()*(TriBb.GetNumModes()
188  -TriBa.GetNumModes());
189  m_npoints += TriBa.GetNumPoints()*TriBb.GetNumPoints();
190  }
191  else if ((QuadrilateralGeom = boost::dynamic_pointer_cast<
192  SpatialDomains::QuadGeom>(expIt->second->m_geomShPtr)))
193  {
195  = expIt->second->m_basisKeyVector[0];
197  = expIt->second->m_basisKeyVector[1];
198 
200  ::AllocateSharedPtr(QuadBa,QuadBb,
201  QuadrilateralGeom);
202  quad->SetElmtId(elmtid++);
203  (*m_exp).push_back(quad);
204 
205  m_ncoeffs += QuadBa.GetNumModes()*QuadBb.GetNumModes();
206  m_npoints += QuadBa.GetNumPoints()*QuadBb.GetNumPoints();
207  }
208  else
209  {
210  ASSERTL0(false, "dynamic cast to a proper Geometry2D "
211  "failed");
212  }
213 
214  }
215 
216  // set up element numbering
217  for(int i = 0; i < (*m_exp).size(); ++i)
218  {
219  (*m_exp)[i]->SetElmtId(i);
220  }
221 
222  // Setup Default optimisation information.
223  int nel = GetExpSize();
226 
227 
228  // set up offset arrays.
230 
231  if (DeclareCoeffPhysArrays)
232  {
233  // Set up m_coeffs, m_phys.
236  }
237 
240  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
ExpList()
The default constructor.
Definition: ExpList.cpp:93
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
BasisType GetBasisType() const
Return type of expansion basis.
Definition: Basis.h:139
boost::shared_ptr< QuadGeom > QuadGeomSharedPtr
Definition: QuadGeom.h:54
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:956
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:939
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:917
int GetNumPoints() const
Return points order at which basis is defined.
Definition: Basis.h:128
Principle Orthogonal Functions .
Definition: BasisType.h:46
boost::shared_ptr< QuadExp > QuadExpSharedPtr
Definition: QuadExp.h:293
void SetCoeffPhysOffsets(void)
Definition of the total number of degrees of freedom and quadrature points and offsets to access data...
Definition: ExpList2D.cpp:1120
boost::shared_ptr< GeometryVector > Composite
Definition: MeshGraph.h:114
PointsKey GetPointsKey() const
Return distribution of points.
Definition: Basis.h:145
boost::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:58
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2950
int GetNumModes() const
Returns the order of the basis.
Definition: Basis.h:84
Lagrange for SEM basis .
Definition: BasisType.h:53
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
Describes the specification for a Basis.
Definition: Basis.h:50
2D Nodal Electrostatic Points on a Triangle
Definition: PointsType.h:68
std::map< int, ExpansionShPtr > ExpansionMap
Definition: MeshGraph.h:174
boost::shared_ptr< TriExp > TriExpSharedPtr
Definition: TriExp.h:291
boost::shared_ptr< NodalTriExp > NodalTriExpSharedPtr
Definition: NodalTriExp.h:379
Nektar::MultiRegions::ExpList2D::ExpList2D ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::ExpansionMap expansions,
const bool  DeclareCoeffPhysArrays = true 
)

Sets up a list of local expansions based on an expansion Map.

Given an expansion vector expansions, containing information about the domain and the spectral/hp element expansion, this constructor fills the list of local expansions {m_exp} with the proper expansions, calculates the total number of quadrature points $\boldsymbol{x}_i$ and local expansion coefficients $\hat{u}^e_n$ and allocates memory for the arrays m_coeffs and m_phys.

Parameters
expansionsA vector containing information about the domain and the spectral/hp element expansion.

Definition at line 257 of file ExpList2D.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::MultiRegions::ExpList::CreateCollections(), Nektar::MultiRegions::e2D, Nektar::LibUtilities::eGLL_Lagrange, Nektar::LibUtilities::eNodalTriElec, Nektar::LibUtilities::eOrtho_A, Nektar::LibUtilities::BasisKey::GetBasisType(), Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::LibUtilities::BasisKey::GetNumModes(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::LibUtilities::BasisKey::GetPointsKey(), Nektar::MultiRegions::ExpList::m_coeffs, Nektar::MultiRegions::ExpList::m_globalOptParam, Nektar::MultiRegions::ExpList::m_ncoeffs, Nektar::MultiRegions::ExpList::m_npoints, Nektar::MultiRegions::ExpList::m_phys, Nektar::MultiRegions::ExpList::ReadGlobalOptimizationParameters(), SetCoeffPhysOffsets(), and Nektar::MultiRegions::ExpList::SetExpType().

260  :ExpList(pSession)
261  {
262  SetExpType(e2D);
263 
264  int elmtid=0;
270 
272  for (expIt = expansions.begin(); expIt != expansions.end(); ++expIt)
273  {
275  SpatialDomains::QuadGeomSharedPtr QuadrilateralGeom;
276 
277  if ((TriangleGeom = boost::dynamic_pointer_cast<SpatialDomains
278  ::TriGeom>(expIt->second->m_geomShPtr)))
279  {
281  = expIt->second->m_basisKeyVector[0];
283  = expIt->second->m_basisKeyVector[1];
284 
285  // This is not elegantly implemented needs re-thinking.
287  {
289  TriBa.GetNumModes(),
290  TriBa.GetPointsKey());
291 
294  ::AllocateSharedPtr(newBa,TriBb,TriNb,
295  TriangleGeom);
296  Ntri->SetElmtId(elmtid++);
297  (*m_exp).push_back(Ntri);
298  }
299  else
300  {
302  ::AllocateSharedPtr(TriBa,TriBb,
303  TriangleGeom);
304  tri->SetElmtId(elmtid++);
305  (*m_exp).push_back(tri);
306  }
307  m_ncoeffs += (TriBa.GetNumModes()*(TriBa.GetNumModes()+1))/2
308  + TriBa.GetNumModes()*(TriBb.GetNumModes()
309  -TriBa.GetNumModes());
310  m_npoints += TriBa.GetNumPoints()*TriBb.GetNumPoints();
311  }
312  else if ((QuadrilateralGeom = boost::dynamic_pointer_cast<
313  SpatialDomains::QuadGeom>(expIt->second->m_geomShPtr)))
314  {
316  = expIt->second->m_basisKeyVector[0];
318  = expIt->second->m_basisKeyVector[1];
319 
321  ::AllocateSharedPtr(QuadBa,QuadBb,
322  QuadrilateralGeom);
323  quad->SetElmtId(elmtid++);
324  (*m_exp).push_back(quad);
325 
326  m_ncoeffs += QuadBa.GetNumModes()*QuadBb.GetNumModes();
327  m_npoints += QuadBa.GetNumPoints()*QuadBb.GetNumPoints();
328  }
329  else
330  {
331  ASSERTL0(false, "dynamic cast to a proper Geometry2D "
332  "failed");
333  }
334 
335  }
336 
337  // Setup Default optimisation information.
338  int nel = GetExpSize();
341 
342 
343  // set up offset arrays.
345 
346  if (DeclareCoeffPhysArrays)
347  {
348  // Set up m_coeffs, m_phys.
351  }
352 
355  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::map< int, ExpansionShPtr >::const_iterator ExpansionMapConstIter
Definition: MeshGraph.h:176
ExpList()
The default constructor.
Definition: ExpList.cpp:93
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
BasisType GetBasisType() const
Return type of expansion basis.
Definition: Basis.h:139
boost::shared_ptr< QuadGeom > QuadGeomSharedPtr
Definition: QuadGeom.h:54
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:956
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:939
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:917
int GetNumPoints() const
Return points order at which basis is defined.
Definition: Basis.h:128
Principle Orthogonal Functions .
Definition: BasisType.h:46
boost::shared_ptr< QuadExp > QuadExpSharedPtr
Definition: QuadExp.h:293
void SetCoeffPhysOffsets(void)
Definition of the total number of degrees of freedom and quadrature points and offsets to access data...
Definition: ExpList2D.cpp:1120
boost::shared_ptr< GeometryVector > Composite
Definition: MeshGraph.h:114
PointsKey GetPointsKey() const
Return distribution of points.
Definition: Basis.h:145
boost::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:58
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2950
int GetNumModes() const
Returns the order of the basis.
Definition: Basis.h:84
Lagrange for SEM basis .
Definition: BasisType.h:53
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
Describes the specification for a Basis.
Definition: Basis.h:50
2D Nodal Electrostatic Points on a Triangle
Definition: PointsType.h:68
boost::shared_ptr< TriExp > TriExpSharedPtr
Definition: TriExp.h:291
boost::shared_ptr< NodalTriExp > NodalTriExpSharedPtr
Definition: NodalTriExp.h:379
Nektar::MultiRegions::ExpList2D::ExpList2D ( const LibUtilities::SessionReaderSharedPtr pSession,
const LibUtilities::BasisKey TriBa,
const LibUtilities::BasisKey TriBb,
const LibUtilities::BasisKey QuadBa,
const LibUtilities::BasisKey QuadBb,
const SpatialDomains::MeshGraphSharedPtr graph2D,
const LibUtilities::PointsType  TriNb = LibUtilities::SIZE_PointsType 
)

Sets up a list of local expansions based on an input mesh and separately defined basiskeys.

Given a mesh graph2D, containing information about the domain and the a list of basiskeys, this constructor fills the list of local expansions {m_exp} with the proper expansions, calculates the total number of quadrature points $\boldsymbol{x}_i$ and local expansion coefficients $\hat{u}^e_n$ and allocates memory for the arrays m_coeffs and m_phys.

Parameters
TriBaA BasisKey, containing the definition of the basis in the first coordinate direction for triangular elements
TriBbA BasisKey, containing the definition of the basis in the second coordinate direction for triangular elements
QuadBaA BasisKey, containing the definition of the basis in the first coordinate direction for quadrilateral elements
QuadBbA BasisKey, containing the definition of the basis in the second coordinate direction for quadrilateral elements
graph2DA mesh, containing information about the domain and the spectral/hp element expansion.
TriNbThe PointsType of possible nodal points

Definition at line 383 of file ExpList2D.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::MultiRegions::ExpList::CreateCollections(), Nektar::MultiRegions::e2D, Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::LibUtilities::BasisKey::GetNumModes(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::MultiRegions::ExpList::m_coeffs, Nektar::MultiRegions::ExpList::m_globalOptParam, Nektar::MultiRegions::ExpList::m_ncoeffs, Nektar::MultiRegions::ExpList::m_npoints, Nektar::MultiRegions::ExpList::m_phys, Nektar::MultiRegions::ExpList::m_physState, Nektar::MultiRegions::ExpList::ReadGlobalOptimizationParameters(), SetCoeffPhysOffsets(), Nektar::MultiRegions::ExpList::SetExpType(), and Nektar::LibUtilities::SIZE_PointsType.

390  :ExpList(pSession, graph2D)
391  {
392  SetExpType(e2D);
393 
394  int elmtid=0;
399 
400  const SpatialDomains::ExpansionMap &expansions =
401  graph2D->GetExpansions();
402  m_ncoeffs = 0;
403  m_npoints = 0;
404 
405  m_physState = false;
406 
407  SpatialDomains::ExpansionMap::const_iterator expIt;
408  for (expIt = expansions.begin(); expIt != expansions.end(); ++expIt)
409  {
411  SpatialDomains::QuadGeomSharedPtr QuadrilateralGeom;
412 
413  if ((TriangleGeom = boost::dynamic_pointer_cast<SpatialDomains::
414  TriGeom>(expIt->second->m_geomShPtr)))
415  {
416  if (TriNb < LibUtilities::SIZE_PointsType)
417  {
419  AllocateSharedPtr(TriBa, TriBb, TriNb,
420  TriangleGeom);
421  Ntri->SetElmtId(elmtid++);
422  (*m_exp).push_back(Ntri);
423  }
424  else
425  {
427  AllocateSharedPtr(TriBa, TriBb, TriangleGeom);
428  tri->SetElmtId(elmtid++);
429  (*m_exp).push_back(tri);
430  }
431 
432  m_ncoeffs += (TriBa.GetNumModes()*(TriBa.GetNumModes()+1))/2
433  + TriBa.GetNumModes() * (TriBb.GetNumModes() -
434  TriBa.GetNumModes());
435  m_npoints += TriBa.GetNumPoints()*TriBb.GetNumPoints();
436  }
437  else if ((QuadrilateralGeom = boost::dynamic_pointer_cast<
438  SpatialDomains::QuadGeom>(expIt->second->m_geomShPtr)))
439  {
441  AllocateSharedPtr(QuadBa, QuadBb, QuadrilateralGeom);
442  quad->SetElmtId(elmtid++);
443  (*m_exp).push_back(quad);
444 
445  m_ncoeffs += QuadBa.GetNumModes()*QuadBb.GetNumModes();
446  m_npoints += QuadBa.GetNumPoints()*QuadBb.GetNumPoints();
447  }
448  else
449  {
450  ASSERTL0(false,
451  "dynamic cast to a proper Geometry2D failed");
452  }
453 
454  }
455 
456  // Setup Default optimisation information.
457  int nel = GetExpSize();
460 
461  // Set up m_coeffs, m_phys and offset arrays.
465 
468  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
ExpList()
The default constructor.
Definition: ExpList.cpp:93
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
boost::shared_ptr< QuadGeom > QuadGeomSharedPtr
Definition: QuadGeom.h:54
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:956
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:939
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
bool m_physState
The state of the array m_phys.
Definition: ExpList.h:965
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:917
int GetNumPoints() const
Return points order at which basis is defined.
Definition: Basis.h:128
boost::shared_ptr< QuadExp > QuadExpSharedPtr
Definition: QuadExp.h:293
void SetCoeffPhysOffsets(void)
Definition of the total number of degrees of freedom and quadrature points and offsets to access data...
Definition: ExpList2D.cpp:1120
boost::shared_ptr< GeometryVector > Composite
Definition: MeshGraph.h:114
boost::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:58
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2950
int GetNumModes() const
Returns the order of the basis.
Definition: Basis.h:84
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
std::map< int, ExpansionShPtr > ExpansionMap
Definition: MeshGraph.h:174
boost::shared_ptr< TriExp > TriExpSharedPtr
Definition: TriExp.h:291
boost::shared_ptr< NodalTriExp > NodalTriExpSharedPtr
Definition: NodalTriExp.h:379
Nektar::MultiRegions::ExpList2D::ExpList2D ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, const ExpListSharedPtr > &  bndConstraint,
const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > &  bndCond,
const LocalRegions::ExpansionVector locexp,
const SpatialDomains::MeshGraphSharedPtr graph3D,
const PeriodicMap periodicFaces,
const bool  DeclareCoeffPhysArrays = true,
const std::string  variable = "DefaultVar" 
)

Specialized constructor for trace expansions. Store expansions for the trace space used in DisContField3D

Parameters
bndConstraintArray of ExpList2D objects each containing a 2D spectral/hp element expansion on a single boundary region.
bndCondArray of BoundaryCondition objects which contain information about the boundary conditions on the different boundary regions.
locexpComplete domain expansion list.
graph3D3D mesh corresponding to the expansion list.
periodicFacesList of periodic faces.
DeclareCoeffPhysArraysIf true, set up m_coeffs, m_phys arrays

Definition at line 486 of file ExpList2D.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::MultiRegions::ExpList::CreateCollections(), Nektar::MultiRegions::e2D, Nektar::SpatialDomains::eDirichlet, Nektar::LibUtilities::BasisKey::GetBasisType(), Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::LocalRegions::Expansion2D::GetGeom2D(), Nektar::LocalRegions::Expansion3D::GetGeom3D(), Nektar::StdRegions::StdExpansion::GetNfaces(), Nektar::LibUtilities::BasisKey::GetNumModes(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::LibUtilities::BasisKey::GetPointsType(), Nektar::iterator, Nektar::MultiRegions::ExpList::m_coeffs, Nektar::MultiRegions::ExpList::m_globalOptParam, Nektar::MultiRegions::ExpList::m_ncoeffs, Nektar::MultiRegions::ExpList::m_npoints, Nektar::MultiRegions::ExpList::m_phys, Nektar::LibUtilities::ReduceSum, SetCoeffPhysOffsets(), Nektar::MultiRegions::ExpList::SetExpType(), and Vmath::Vsum().

494  :
495  ExpList(pSession, graph3D)
496  {
497  SetExpType(e2D);
498 
499  int i, j, id, elmtid=0;
500  set<int> facesDone;
501 
505  LocalRegions::QuadExpSharedPtr FaceQuadExp;
509 
510  // First loop over boundary conditions to renumber
511  // Dirichlet boundaries
512  for (i = 0; i < bndCond.num_elements(); ++i)
513  {
514  if (bndCond[i]->GetBoundaryConditionType()
516  {
517  for (j = 0; j < bndConstraint[i]->GetExpSize(); ++j)
518  {
519  LibUtilities::BasisKey bkey0 = bndConstraint[i]
520  ->GetExp(j)->GetBasis(0)->GetBasisKey();
521  LibUtilities::BasisKey bkey1 = bndConstraint[i]
522  ->GetExp(j)->GetBasis(1)->GetBasisKey();
523  exp2D = bndConstraint[i]->GetExp(j)
525  FaceGeom = exp2D->GetGeom2D();
526 
527  //if face is a quad
528  if((FaceQuadGeom = boost::dynamic_pointer_cast<
529  SpatialDomains::QuadGeom>(FaceGeom)))
530  {
532  ::AllocateSharedPtr(bkey0, bkey1, FaceQuadGeom);
533  facesDone.insert(FaceQuadGeom->GetFid());
534  FaceQuadExp->SetElmtId(elmtid++);
535  (*m_exp).push_back(FaceQuadExp);
536  }
537  //if face is a triangle
538  else if((FaceTriGeom = boost::dynamic_pointer_cast<
539  SpatialDomains::TriGeom>(FaceGeom)))
540  {
542  ::AllocateSharedPtr(bkey0, bkey1, FaceTriGeom);
543  facesDone.insert(FaceTriGeom->GetFid());
544  FaceTriExp->SetElmtId(elmtid++);
545  (*m_exp).push_back(FaceTriExp);
546  }
547  else
548  {
549  ASSERTL0(false,"dynamic cast to a proper face geometry failed");
550  }
551  }
552  }
553  }
554 
557  LibUtilities::BasisKey> > > faceOrders;
559  pair<LibUtilities::BasisKey,
560  LibUtilities::BasisKey> > >::iterator it;
561 
562  for(i = 0; i < locexp.size(); ++i)
563  {
564  exp3D = locexp[i]->as<LocalRegions::Expansion3D>();
565  for (j = 0; j < exp3D->GetNfaces(); ++j)
566  {
567  FaceGeom = exp3D->GetGeom3D()->GetFace(j);
568  id = FaceGeom->GetFid();
569 
570  if(facesDone.count(id) != 0)
571  {
572  continue;
573  }
574  it = faceOrders.find(id);
575 
576  if (it == faceOrders.end())
577  {
578  LibUtilities::BasisKey face_dir0
579  = locexp[i]->DetFaceBasisKey(j,0);
580  LibUtilities::BasisKey face_dir1
581  = locexp[i]->DetFaceBasisKey(j,1);
582 
583  faceOrders.insert(
584  std::make_pair(
585  id, std::make_pair(
586  FaceGeom,
587  std::make_pair(face_dir0, face_dir1))));
588  }
589  else // variable modes/points
590  {
591  LibUtilities::BasisKey face0 =
592  locexp[i]->DetFaceBasisKey(j,0);
593  LibUtilities::BasisKey face1 =
594  locexp[i]->DetFaceBasisKey(j,1);
595  LibUtilities::BasisKey existing0 =
596  it->second.second.first;
597  LibUtilities::BasisKey existing1 =
598  it->second.second.second;
599 
600  int np11 = face0 .GetNumPoints();
601  int np12 = face1 .GetNumPoints();
602  int np21 = existing0.GetNumPoints();
603  int np22 = existing1.GetNumPoints();
604  int nm11 = face0 .GetNumModes ();
605  int nm12 = face1 .GetNumModes ();
606  int nm21 = existing0.GetNumModes ();
607  int nm22 = existing1.GetNumModes ();
608 
609  if ((np22 >= np12 || np21 >= np11) &&
610  (nm22 >= nm12 || nm21 >= nm11))
611  {
612  continue;
613  }
614  else if((np22 < np12 || np21 < np11) &&
615  (nm22 < nm12 || nm21 < nm11))
616  {
617  it->second.second.first = face0;
618  it->second.second.second = face1;
619  }
620  else
621  {
622  ASSERTL0(false,
623  "inappropriate number of points/modes (max "
624  "num of points is not set with max order)");
625  }
626  }
627  }
628  }
629 
630  LibUtilities::CommSharedPtr vComm = pSession->GetComm();
631  int nproc = vComm->GetSize(); // number of processors
632  int facepr = vComm->GetRank(); // ID processor
633 
634  if (nproc > 1)
635  {
636  int fCnt = 0;
637 
638  // Count the number of faces on each partition
639  for(i = 0; i < locexp.size(); ++i)
640  {
641  fCnt += locexp[i]->GetNfaces();
642  }
643 
644  // Set up the offset and the array that will contain the list of
645  // face IDs, then reduce this across processors.
646  Array<OneD, int> faceCnt(nproc,0);
647  faceCnt[facepr] = fCnt;
648  vComm->AllReduce(faceCnt, LibUtilities::ReduceSum);
649 
650  int totFaceCnt = Vmath::Vsum(nproc, faceCnt, 1);
651  Array<OneD, int> fTotOffsets(nproc,0);
652 
653  for (i = 1; i < nproc; ++i)
654  {
655  fTotOffsets[i] = fTotOffsets[i-1] + faceCnt[i-1];
656  }
657 
658  // Local list of the edges per element
659 
660  Array<OneD, int> FacesTotID (totFaceCnt, 0);
661  Array<OneD, int> FacesTotNm0 (totFaceCnt, 0);
662  Array<OneD, int> FacesTotNm1 (totFaceCnt, 0);
663  Array<OneD, int> FacesTotPnts0(totFaceCnt, 0);
664  Array<OneD, int> FacesTotPnts1(totFaceCnt, 0);
665 
666  int cntr = fTotOffsets[facepr];
667 
668  for(i = 0; i < locexp.size(); ++i)
669  {
670  exp3D = locexp[i]->as<LocalRegions::Expansion3D>();
671 
672  int nfaces = locexp[i]->GetNfaces();
673 
674  for(j = 0; j < nfaces; ++j, ++cntr)
675  {
676  LibUtilities::BasisKey face_dir0
677  = locexp[i]->DetFaceBasisKey(j,0);
678  LibUtilities::BasisKey face_dir1
679  = locexp[i]->DetFaceBasisKey(j,1);
680 
681  FacesTotID[cntr] = exp3D->GetGeom3D()->GetFid(j);
682  FacesTotNm0[cntr] = face_dir0.GetNumModes ();
683  FacesTotNm1[cntr] = face_dir1.GetNumModes ();
684  FacesTotPnts0[cntr] = face_dir0.GetNumPoints();
685  FacesTotPnts1[cntr] = face_dir1.GetNumPoints();
686  }
687  }
688 
689  vComm->AllReduce(FacesTotID, LibUtilities::ReduceSum);
690  vComm->AllReduce(FacesTotNm0, LibUtilities::ReduceSum);
691  vComm->AllReduce(FacesTotNm1, LibUtilities::ReduceSum);
692  vComm->AllReduce(FacesTotPnts0, LibUtilities::ReduceSum);
693  vComm->AllReduce(FacesTotPnts1, LibUtilities::ReduceSum);
694 
695  for (i = 0; i < totFaceCnt; ++i)
696  {
697  it = faceOrders.find(FacesTotID[i]);
698 
699  if (it == faceOrders.end())
700  {
701  continue;
702  }
703 
704  LibUtilities::BasisKey existing0 =
705  it->second.second.first;
706  LibUtilities::BasisKey existing1 =
707  it->second.second.second;
708  LibUtilities::BasisKey face0(
709  existing0.GetBasisType(), FacesTotNm0[i],
710  LibUtilities::PointsKey(FacesTotPnts0[i],
711  existing0.GetPointsType()));
712  LibUtilities::BasisKey face1(
713  existing1.GetBasisType(), FacesTotNm1[i],
714  LibUtilities::PointsKey(FacesTotPnts1[i],
715  existing1.GetPointsType()));
716 
717  int np11 = face0 .GetNumPoints();
718  int np12 = face1 .GetNumPoints();
719  int np21 = existing0.GetNumPoints();
720  int np22 = existing1.GetNumPoints();
721  int nm11 = face0 .GetNumModes ();
722  int nm12 = face1 .GetNumModes ();
723  int nm21 = existing0.GetNumModes ();
724  int nm22 = existing1.GetNumModes ();
725 
726  if ((np22 >= np12 || np21 >= np11) &&
727  (nm22 >= nm12 || nm21 >= nm11))
728  {
729  continue;
730  }
731  else if((np22 < np12 || np21 < np11) &&
732  (nm22 < nm12 || nm21 < nm11))
733  {
734  it->second.second.first = face0;
735  it->second.second.second = face1;
736  }
737  else
738  {
739  ASSERTL0(false,
740  "inappropriate number of points/modes (max "
741  "num of points is not set with max order)");
742  }
743  }
744  }
745 
746  for (it = faceOrders.begin(); it != faceOrders.end(); ++it)
747  {
748  FaceGeom = it->second.first;
749 
750  if ((FaceQuadGeom = boost::dynamic_pointer_cast<
751  SpatialDomains::QuadGeom>(FaceGeom)))
752  {
754  ::AllocateSharedPtr(it->second.second.first,
755  it->second.second.second,
756  FaceQuadGeom);
757  FaceQuadExp->SetElmtId(elmtid++);
758  (*m_exp).push_back(FaceQuadExp);
759  }
760  else if ((FaceTriGeom = boost::dynamic_pointer_cast<
761  SpatialDomains::TriGeom>(FaceGeom)))
762  {
764  ::AllocateSharedPtr(it->second.second.first,
765  it->second.second.second,
766  FaceTriGeom);
767  FaceTriExp->SetElmtId(elmtid++);
768  (*m_exp).push_back(FaceTriExp);
769  }
770  }
771 
772  // Setup Default optimisation information.
773  int nel = GetExpSize();
774 
777 
778  // Set up offset information and array sizes
780 
781  // Set up m_coeffs, m_phys.
782  if(DeclareCoeffPhysArrays)
783  {
786  }
787 
789  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
ExpList()
The default constructor.
Definition: ExpList.cpp:93
int GetNfaces() const
This function returns the number of faces of the expansion domain.
Definition: StdExpansion.h:438
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
BasisType GetBasisType() const
Return type of expansion basis.
Definition: Basis.h:139
boost::shared_ptr< QuadGeom > QuadGeomSharedPtr
Definition: QuadGeom.h:54
PointsType GetPointsType() const
Return type of quadrature.
Definition: Basis.h:151
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:956
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:939
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
SpatialDomains::Geometry2DSharedPtr GetGeom2D() const
Definition: Expansion2D.h:269
boost::shared_ptr< Expansion3D > Expansion3DSharedPtr
Definition: Expansion2D.h:48
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:917
int GetNumPoints() const
Return points order at which basis is defined.
Definition: Basis.h:128
Defines a specification for a set of points.
Definition: Points.h:58
boost::shared_ptr< QuadExp > QuadExpSharedPtr
Definition: QuadExp.h:293
void SetCoeffPhysOffsets(void)
Definition of the total number of degrees of freedom and quadrature points and offsets to access data...
Definition: ExpList2D.cpp:1120
boost::shared_ptr< Geometry2D > Geometry2DSharedPtr
Definition: Geometry2D.h:59
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
boost::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:58
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2950
const LibUtilities::BasisKey DetFaceBasisKey(const int i, const int k) const
Definition: StdExpansion.h:324
int GetNumModes() const
Returns the order of the basis.
Definition: Basis.h:84
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.cpp:723
Describes the specification for a Basis.
Definition: Basis.h:50
boost::shared_ptr< Expansion2D > Expansion2DSharedPtr
Definition: Expansion1D.h:49
boost::shared_ptr< TriExp > TriExpSharedPtr
Definition: TriExp.h:291
Nektar::MultiRegions::ExpList2D::ExpList2D ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::CompositeMap domain,
const SpatialDomains::MeshGraphSharedPtr graph3D,
const std::string  variable = "DefaultVar" 
)

Specialised constructor for Neumann boundary conditions in DisContField3D and ContField3D.

Fills the list of local expansions with the segments from the 3D mesh specified by domain. This CompositeMap contains a list of Composites which define the Neumann boundary.

See also
ExpList2D::ExpList2D(SpatialDomains::MeshGraph2D&) for details.
Parameters
domainA domain, comprising of one or more composite regions.
graph3DA mesh, containing information about the domain and the spectral/hp element expansions.

Definition at line 802 of file ExpList2D.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::MultiRegions::ExpList::CreateCollections(), Nektar::MultiRegions::e2D, Nektar::LibUtilities::eGLL_Lagrange, Nektar::LibUtilities::eNodalTriElec, Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::LibUtilities::BasisKey::GetNumModes(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::MultiRegions::ExpList::m_coeffs, Nektar::MultiRegions::ExpList::m_globalOptParam, Nektar::MultiRegions::ExpList::m_ncoeffs, Nektar::MultiRegions::ExpList::m_npoints, Nektar::MultiRegions::ExpList::m_phys, Nektar::MultiRegions::ExpList::ReadGlobalOptimizationParameters(), SetCoeffPhysOffsets(), and Nektar::MultiRegions::ExpList::SetExpType().

806  :ExpList(pSession, graph3D)
807  {
808 
809  SetExpType(e2D);
810 
811  ASSERTL0(boost::dynamic_pointer_cast<
812  SpatialDomains::MeshGraph3D>(graph3D),
813  "Expected a MeshGraph3D object.");
814 
815  int j, elmtid=0;
816  int nel = 0;
817 
820  SpatialDomains::QuadGeomSharedPtr QuadrilateralGeom;
821 
826 
827  SpatialDomains::CompositeMap::const_iterator compIt;
828  for (compIt = domain.begin(); compIt != domain.end(); ++compIt)
829  {
830  nel += (compIt->second)->size();
831  }
832 
833  for (compIt = domain.begin(); compIt != domain.end(); ++compIt)
834  {
835  for (j = 0; j < compIt->second->size(); ++j)
836  {
837  if ((TriangleGeom = boost::dynamic_pointer_cast<
838  SpatialDomains::TriGeom>((*compIt->second)[j])))
839  {
841  = boost::dynamic_pointer_cast<
842  SpatialDomains::MeshGraph3D>(graph3D)->
843  GetFaceBasisKey(TriangleGeom, 0, variable);
845  = boost::dynamic_pointer_cast<
846  SpatialDomains::MeshGraph3D>(graph3D)->
847  GetFaceBasisKey(TriangleGeom,1,variable);
848 
849  if (graph3D->GetExpansions().begin()->second->
850  m_basisKeyVector[0].GetBasisType() ==
852  {
853  ASSERTL0(false,"This method needs sorting");
855 
857  ::AllocateSharedPtr(TriBa,TriBb,TriNb,
858  TriangleGeom);
859  Ntri->SetElmtId(elmtid++);
860  (*m_exp).push_back(Ntri);
861  }
862  else
863  {
865  ::AllocateSharedPtr(TriBa, TriBb,
866  TriangleGeom);
867  tri->SetElmtId(elmtid++);
868  (*m_exp).push_back(tri);
869  }
870 
871  m_ncoeffs
872  += (TriBa.GetNumModes()*(TriBa.GetNumModes()+1))/2
873  + TriBa.GetNumModes()*(TriBb.GetNumModes()
874  -TriBa.GetNumModes());
875  m_npoints += TriBa.GetNumPoints()*TriBb.GetNumPoints();
876  }
877  else if ((QuadrilateralGeom = boost::dynamic_pointer_cast<
878  SpatialDomains::QuadGeom>((*compIt->second)[j])))
879  {
881  = boost::dynamic_pointer_cast<
882  SpatialDomains::MeshGraph3D>(graph3D)->
883  GetFaceBasisKey(QuadrilateralGeom, 0,
884  variable);
886  = boost::dynamic_pointer_cast<
887  SpatialDomains::MeshGraph3D>(graph3D)->
888  GetFaceBasisKey(QuadrilateralGeom, 1,
889  variable);
890 
892  ::AllocateSharedPtr(QuadBa, QuadBb,
893  QuadrilateralGeom);
894  quad->SetElmtId(elmtid++);
895  (*m_exp).push_back(quad);
896 
897  m_ncoeffs += QuadBa.GetNumModes()*QuadBb.GetNumModes();
898  m_npoints += QuadBa.GetNumPoints()
899  * QuadBb.GetNumPoints();
900  }
901  else
902  {
903  ASSERTL0(false,
904  "dynamic cast to a proper Geometry2D failed");
905  }
906  }
907 
908  }
909 
910  // Setup Default optimisation information.
911  nel = GetExpSize();
914 
915  // Set up m_coeffs, m_phys and offset arrays.
919 
922  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
ExpList()
The default constructor.
Definition: ExpList.cpp:93
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
boost::shared_ptr< QuadGeom > QuadGeomSharedPtr
Definition: QuadGeom.h:54
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:956
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:939
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:917
int GetNumPoints() const
Return points order at which basis is defined.
Definition: Basis.h:128
boost::shared_ptr< QuadExp > QuadExpSharedPtr
Definition: QuadExp.h:293
void SetCoeffPhysOffsets(void)
Definition of the total number of degrees of freedom and quadrature points and offsets to access data...
Definition: ExpList2D.cpp:1120
boost::shared_ptr< GeometryVector > Composite
Definition: MeshGraph.h:114
boost::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:58
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2950
int GetNumModes() const
Returns the order of the basis.
Definition: Basis.h:84
Lagrange for SEM basis .
Definition: BasisType.h:53
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
Describes the specification for a Basis.
Definition: Basis.h:50
2D Nodal Electrostatic Points on a Triangle
Definition: PointsType.h:68
boost::shared_ptr< TriExp > TriExpSharedPtr
Definition: TriExp.h:291
boost::shared_ptr< NodalTriExp > NodalTriExpSharedPtr
Definition: NodalTriExp.h:379
Nektar::MultiRegions::ExpList2D::~ExpList2D ( )
virtual

Destructor.

Definition at line 78 of file ExpList2D.cpp.

79  {
80  }

Member Function Documentation

void Nektar::MultiRegions::ExpList2D::SetCoeffPhysOffsets ( void  )
private

Definition of the total number of degrees of freedom and quadrature points and offsets to access datax.

Each expansion (local element) is processed in turn to determine the number of coefficients and physical data points it contributes to the domain. Three arrays, m_coeff_offset, m_phys_offset and m_offset_elmt_id, are initialised and updated to store the data offsets of each element in the m_coeffs and m_phys arrays, and the element id that each consecutive block is associated respectively.

Definition at line 1120 of file ExpList2D.cpp.

References Nektar::LibUtilities::eQuadrilateral, Nektar::LibUtilities::eTriangle, Nektar::MultiRegions::ExpList::m_coeff_offset, Nektar::MultiRegions::ExpList::m_exp, Nektar::MultiRegions::ExpList::m_ncoeffs, Nektar::MultiRegions::ExpList::m_npoints, Nektar::MultiRegions::ExpList::m_offset_elmt_id, and Nektar::MultiRegions::ExpList::m_phys_offset.

Referenced by ExpList2D().

1121  {
1122  int i;
1123 
1124  // Set up offset information and array sizes
1126  m_phys_offset = Array<OneD,int>(m_exp->size());
1128 
1129  m_ncoeffs = m_npoints = 0;
1130 
1131  int cnt = 0;
1132  for(i = 0; i < m_exp->size(); ++i)
1133  {
1134  if((*m_exp)[i]->DetShapeType() == LibUtilities::eTriangle)
1135  {
1137  m_phys_offset [i] = m_npoints;
1138  m_offset_elmt_id[cnt++] = i;
1139  m_ncoeffs += (*m_exp)[i]->GetNcoeffs();
1140  m_npoints += (*m_exp)[i]->GetTotPoints();
1141  }
1142  }
1143 
1144  for(i = 0; i < m_exp->size(); ++i)
1145  {
1146  if((*m_exp)[i]->DetShapeType() == LibUtilities::eQuadrilateral)
1147  {
1149  m_phys_offset [i] = m_npoints;
1150  m_offset_elmt_id[cnt++] = i;
1151  m_ncoeffs += (*m_exp)[i]->GetNcoeffs();
1152  m_npoints += (*m_exp)[i]->GetTotPoints();
1153  }
1154  }
1155  }
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:988
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:991
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:917
Array< OneD, int > m_offset_elmt_id
Array containing the element id m_offset_elmt_id[n] that the n^th consecutive block of data in m_coef...
Definition: ExpList.h:999
void Nektar::MultiRegions::ExpList2D::v_GetNormals ( Array< OneD, Array< OneD, NekDouble > > &  normals)
protectedvirtual

For each local element, copy the normals stored in the element list into the array normals.

Parameters
normalsMultidimensional array in which to copy normals to. Must have dimension equal to or larger than the spatial dimension of the elements.

Reimplemented from Nektar::MultiRegions::ExpList.

Definition at line 1038 of file ExpList2D.cpp.

References Nektar::MultiRegions::AlignFace(), ASSERTL1, Nektar::StdRegions::eDir1FwdDir2_Dir2FwdDir1, Nektar::MultiRegions::ExpList::GetCoordim(), Nektar::LocalRegions::Expansion2D::GetLeftAdjacentElementExp(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::LibUtilities::BasisKey::GetPointsKey(), Nektar::LibUtilities::Interp2D(), Nektar::MultiRegions::ExpList::m_exp, and Nektar::MultiRegions::ExpList::m_phys_offset.

1040  {
1042  int i, j;
1043  const int coordim = GetCoordim(0);
1044 
1045  ASSERTL1(normals.num_elements() >= coordim,
1046  "Output vector does not have sufficient dimensions to "
1047  "match coordim");
1048 
1049  // Process each expansion.
1050  for (i = 0; i < m_exp->size(); ++i)
1051  {
1052  LocalRegions::Expansion2DSharedPtr traceExp = (*m_exp)[i]->as<
1055  traceExp->GetLeftAdjacentElementExp();
1056 
1057  // Get the number of points and normals for this expansion.
1058  int faceNum = traceExp->GetLeftAdjacentElementFace();
1059  int offset = m_phys_offset[i];
1060 
1061  const Array<OneD, const Array<OneD, NekDouble> > &locNormals
1062  = exp3D->GetFaceNormal(faceNum);
1063 
1064  // Project normals from 3D element onto the same orientation as
1065  // the trace expansion.
1066  StdRegions::Orientation orient = exp3D->GetForient(faceNum);
1067 
1068 
1069  int fromid0,fromid1;
1070 
1072  {
1073  fromid0 = 0;
1074  fromid1 = 1;
1075  }
1076  else
1077  {
1078  fromid0 = 1;
1079  fromid1 = 0;
1080  }
1081 
1082  LibUtilities::BasisKey faceBasis0
1083  = exp3D->DetFaceBasisKey(faceNum, fromid0);
1084  LibUtilities::BasisKey faceBasis1
1085  = exp3D->DetFaceBasisKey(faceNum, fromid1);
1086  LibUtilities::BasisKey traceBasis0
1087  = traceExp->GetBasis(0)->GetBasisKey();
1088  LibUtilities::BasisKey traceBasis1
1089  = traceExp->GetBasis(1)->GetBasisKey();
1090 
1091  const int faceNq0 = faceBasis0.GetNumPoints();
1092  const int faceNq1 = faceBasis1.GetNumPoints();
1093 
1094  for (j = 0; j < coordim; ++j)
1095  {
1096  Array<OneD, NekDouble> traceNormals(faceNq0 * faceNq1);
1097  AlignFace(orient, faceNq0, faceNq1,
1098  locNormals[j], traceNormals);
1100  faceBasis0.GetPointsKey(),
1101  faceBasis1.GetPointsKey(),
1102  traceNormals,
1103  traceBasis0.GetPointsKey(),
1104  traceBasis1.GetPointsKey(),
1105  tmp = normals[j]+offset);
1106  }
1107  }
1108  }
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:116
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:991
boost::shared_ptr< Expansion3D > Expansion3DSharedPtr
Definition: Expansion2D.h:48
int GetNumPoints() const
Return points order at which basis is defined.
Definition: Basis.h:128
Expansion3DSharedPtr GetLeftAdjacentElementExp() const
Definition: Expansion2D.h:223
PointsKey GetPointsKey() const
Return distribution of points.
Definition: Basis.h:145
int GetCoordim(int eid)
This function returns the dimension of the coordinates of the element eid.
Definition: ExpList.h:1794
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
Describes the specification for a Basis.
Definition: Basis.h:50
boost::shared_ptr< Expansion2D > Expansion2DSharedPtr
Definition: Expansion1D.h:49
void AlignFace(const StdRegions::Orientation orient, const int nquad1, const int nquad2, const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out)
Helper function to re-align face to a given orientation.
Definition: ExpList2D.cpp:966
void Nektar::MultiRegions::ExpList2D::v_PhysGalerkinProjection1DScaled ( const NekDouble  scale,
const Array< OneD, NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
privatevirtual

Reimplemented from Nektar::MultiRegions::ExpList.

Definition at line 1298 of file ExpList2D.cpp.

References Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::MultiRegions::ExpList::m_exp, and Nektar::LibUtilities::PhysGalerkinProject2D().

1302  {
1303  int cnt,cnt1;
1304 
1305  cnt = cnt1 = 0;
1306  for(int i = 0; i < GetExpSize(); ++i)
1307  {
1308  // get new points key
1309  int pt0 = (*m_exp)[i]->GetNumPoints(0);
1310  int pt1 = (*m_exp)[i]->GetNumPoints(1);
1311  int npt0 = (int) pt0*scale;
1312  int npt1 = (int) pt1*scale;
1313 
1314  LibUtilities::PointsKey newPointsKey0(npt0,
1315  (*m_exp)[i]->GetPointsType(0));
1316  LibUtilities::PointsKey newPointsKey1(npt1,
1317  (*m_exp)[i]->GetPointsType(1));
1318 
1319  // Project points;
1321  newPointsKey0,
1322  newPointsKey1,
1323  &inarray[cnt],
1324  (*m_exp)[i]->GetBasis(0)->GetPointsKey(),
1325  (*m_exp)[i]->GetBasis(1)->GetPointsKey(),
1326  &outarray[cnt1]);
1327 
1328  cnt += npt0*npt1;
1329  cnt1 += pt0*pt1;
1330  }
1331 
1332  }
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
Defines a specification for a set of points.
Definition: Points.h:58
void PhysGalerkinProject2D(const BasisKey &fbasis0, const BasisKey &fbasis1, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, Array< OneD, NekDouble > &to)
void Nektar::MultiRegions::ExpList2D::v_PhysInterp1DScaled ( const NekDouble  scale,
const Array< OneD, NekDouble > &  inarray,
Array< OneD, NekDouble > &  outarray 
)
privatevirtual

Reimplemented from Nektar::MultiRegions::ExpList.

Definition at line 1266 of file ExpList2D.cpp.

References Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::LibUtilities::Interp2D(), and Nektar::MultiRegions::ExpList::m_exp.

1270  {
1271  int cnt,cnt1;
1272 
1273  cnt = cnt1 = 0;
1274  for(int i = 0; i < GetExpSize(); ++i)
1275  {
1276  // get new points key
1277  int pt0 = (*m_exp)[i]->GetNumPoints(0);
1278  int pt1 = (*m_exp)[i]->GetNumPoints(1);
1279  int npt0 = (int) pt0*scale;
1280  int npt1 = (int) pt1*scale;
1281 
1282  LibUtilities::PointsKey newPointsKey0(npt0,
1283  (*m_exp)[i]->GetPointsType(0));
1284  LibUtilities::PointsKey newPointsKey1(npt1,
1285  (*m_exp)[i]->GetPointsType(1));
1286 
1287  // Interpolate points;
1288  LibUtilities::Interp2D((*m_exp)[i]->GetBasis(0)->GetPointsKey(),
1289  (*m_exp)[i]->GetBasis(1)->GetPointsKey(),
1290  &inarray[cnt],newPointsKey0,
1291  newPointsKey1,&outarray[cnt1]);
1292 
1293  cnt += pt0*pt1;
1294  cnt1 += npt0*npt1;
1295  }
1296  }
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
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:116
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
Defines a specification for a set of points.
Definition: Points.h:58
void Nektar::MultiRegions::ExpList2D::v_ReadGlobalOptimizationParameters ( )
privatevirtual

Reimplemented from Nektar::MultiRegions::ExpList.

Definition at line 1174 of file ExpList2D.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::eTriangle, Nektar::MultiRegions::ExpList::GetExpSize(), Nektar::MultiRegions::ExpList::m_exp, Nektar::MultiRegions::ExpList::m_globalOptParam, and Nektar::MultiRegions::ExpList::m_session.

1175  {
1176  Array<OneD, int> NumShape(2,0);
1177 
1178  for(int i = 0; i < GetExpSize(); ++i)
1179  {
1180  if((*m_exp)[i]->DetShapeType() == LibUtilities::eTriangle)
1181  {
1182  NumShape[0] += 1;
1183  }
1184  else // Quadrilateral element
1185  {
1186  NumShape[1] += 1;
1187  }
1188  }
1189 
1191  ::AllocateSharedPtr(m_session,2,NumShape);
1192  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:910
void Nektar::MultiRegions::ExpList2D::v_SetUpPhysNormals ( )
privatevirtual

Set up the normals on each expansion.

Reimplemented from Nektar::MultiRegions::ExpList.

Definition at line 1161 of file ExpList2D.cpp.

References Nektar::MultiRegions::ExpList::m_exp.

1162  {
1163  int i, j;
1164  for (i = 0; i < m_exp->size(); ++i)
1165  {
1166  for (j = 0; j < (*m_exp)[i]->GetNedges(); ++j)
1167  {
1168  (*m_exp)[i]->ComputeEdgeNormal(j);
1169  }
1170  }
1171  }
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
void Nektar::MultiRegions::ExpList2D::v_Upwind ( const Array< OneD, const NekDouble > &  Vn,
const Array< OneD, const NekDouble > &  Fwd,
const Array< OneD, const NekDouble > &  Bwd,
Array< OneD, NekDouble > &  Upwind 
)
protectedvirtual

Upwind the Fwd and Bwd states based on the one- dimensional normal velocity field given by Vn.

One-dimensional upwind.

Parameters
VnVelocity field.
FwdLeft state.
BwdRight state.
UpwindOutput vector.

Reimplemented from Nektar::MultiRegions::ExpList.

Definition at line 931 of file ExpList2D.cpp.

References Nektar::MultiRegions::ExpList::m_exp, and Nektar::MultiRegions::ExpList::m_phys_offset.

936  {
937  int i,j,f_npoints,offset;
938 
939  // Process each expansion.
940  for(i = 0; i < m_exp->size(); ++i)
941  {
942  // Get the number of points and the data offset.
943  f_npoints = (*m_exp)[i]->GetNumPoints(0)*
944  (*m_exp)[i]->GetNumPoints(1);
945  offset = m_phys_offset[i];
946 
947  // Process each point in the expansion.
948  for(j = 0; j < f_npoints; ++j)
949  {
950  // Upwind based on one-dimensional velocity.
951  if(Vn[offset + j] > 0.0)
952  {
953  Upwind[offset + j] = Fwd[offset + j];
954  }
955  else
956  {
957  Upwind[offset + j] = Bwd[offset + j];
958  }
959  }
960  }
961  }
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:991
void Nektar::MultiRegions::ExpList2D::v_WriteVtkPieceHeader ( std::ostream &  outfile,
int  expansion,
int  istrip 
)
privatevirtual

Reimplemented from Nektar::MultiRegions::ExpList.

Definition at line 1194 of file ExpList2D.cpp.

1198  {
1199  int i,j;
1200  int nquad0 = (*m_exp)[expansion]->GetNumPoints(0);
1201  int nquad1 = (*m_exp)[expansion]->GetNumPoints(1);
1202  int ntot = nquad0*nquad1;
1203  int ntotminus = (nquad0-1)*(nquad1-1);
1204 
1205  Array<OneD,NekDouble> coords[3];
1206  coords[0] = Array<OneD,NekDouble>(ntot,0.0);
1207  coords[1] = Array<OneD,NekDouble>(ntot,0.0);
1208  coords[2] = Array<OneD,NekDouble>(ntot,0.0);
1209  (*m_exp)[expansion]->GetCoords(coords[0],coords[1],coords[2]);
1210 
1211  outfile << " <Piece NumberOfPoints=\""
1212  << ntot << "\" NumberOfCells=\""
1213  << ntotminus << "\">" << endl;
1214  outfile << " <Points>" << endl;
1215  outfile << " <DataArray type=\"Float64\" "
1216  << "NumberOfComponents=\"3\" format=\"ascii\">" << endl;
1217  outfile << " ";
1218  for (i = 0; i < ntot; ++i)
1219  {
1220  for (j = 0; j < 3; ++j)
1221  {
1222  outfile << setprecision(8) << scientific
1223  << (float)coords[j][i] << " ";
1224  }
1225  outfile << endl;
1226  }
1227  outfile << endl;
1228  outfile << " </DataArray>" << endl;
1229  outfile << " </Points>" << endl;
1230  outfile << " <Cells>" << endl;
1231  outfile << " <DataArray type=\"Int32\" "
1232  << "Name=\"connectivity\" format=\"ascii\">" << endl;
1233  for (i = 0; i < nquad0-1; ++i)
1234  {
1235  for (j = 0; j < nquad1-1; ++j)
1236  {
1237  outfile << j*nquad0 + i << " "
1238  << j*nquad0 + i + 1 << " "
1239  << (j+1)*nquad0 + i + 1 << " "
1240  << (j+1)*nquad0 + i << endl;
1241  }
1242  }
1243  outfile << endl;
1244  outfile << " </DataArray>" << endl;
1245  outfile << " <DataArray type=\"Int32\" "
1246  << "Name=\"offsets\" format=\"ascii\">" << endl;
1247  for (i = 0; i < ntotminus; ++i)
1248  {
1249  outfile << i*4+4 << " ";
1250  }
1251  outfile << endl;
1252  outfile << " </DataArray>" << endl;
1253  outfile << " <DataArray type=\"UInt8\" "
1254  << "Name=\"types\" format=\"ascii\">" << endl;
1255  for (i = 0; i < ntotminus; ++i)
1256  {
1257  outfile << "9 ";
1258  }
1259  outfile << endl;
1260  outfile << " </DataArray>" << endl;
1261  outfile << " </Cells>" << endl;
1262  outfile << " <PointData>" << endl;
1263  }