Nektar++
Public Member Functions | List of all members
Nektar::CoupledLocalToGlobalC0ContMap Class Reference

#include <CoupledLocalToGlobalC0ContMap.h>

Inheritance diagram for Nektar::CoupledLocalToGlobalC0ContMap:
[legend]

Public Member Functions

 CoupledLocalToGlobalC0ContMap (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &graph, const SpatialDomains::BoundaryConditionsSharedPtr &boundaryConditions, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const MultiRegions::ExpListSharedPtr &pressure, const size_t nz_loc, const bool CheeckForSingularSys=true)
 
void FindEdgeIdToAddMeanPressure (std::vector< std::map< int, int > > &ReorderedGraphVertId, size_t &nel, const LocalRegions::ExpansionVector &locExpVector, size_t &edgeId, size_t &vertId, int &firstNonDirGraphVertId, std::map< int, int > &IsDirEdgeDof, MultiRegions::BottomUpSubStructuredGraphSharedPtr &bottomUpGraph, Array< OneD, int > &AddMeanPressureToEdgeId)
 
- Public Member Functions inherited from Nektar::MultiRegions::AssemblyMapCG
 AssemblyMapCG (const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &comm, const std::string variable="DefaultVar")
 Default constructor. More...
 
 AssemblyMapCG (const LibUtilities::SessionReaderSharedPtr &pSession, const int numLocalCoeffs, const ExpList &locExp, const BndCondExp &bndCondExp=NullExpListSharedPtrArray, const BndCond &bndConditions=SpatialDomains::NullBoundaryConditionShPtrArray, const bool checkIfSingular=false, const std::string variable="defaultVar", const PeriodicMap &periodicVerts=NullPeriodicMap, const PeriodicMap &periodicEdges=NullPeriodicMap, const PeriodicMap &periodicFaces=NullPeriodicMap)
 General constructor for expansions of all dimensions without boundary conditions. More...
 
virtual ~AssemblyMapCG ()
 Destructor. More...
 
std::set< ExtraDirDof > & GetCopyLocalDirDofs ()
 
std::set< int > & GetParallelDirBndSign ()
 
- Public Member Functions inherited from Nektar::MultiRegions::AssemblyMap
 AssemblyMap ()
 Default constructor. More...
 
 AssemblyMap (const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &comm, const std::string variable="DefaultVar")
 Constructor with a communicator. More...
 
 AssemblyMap (AssemblyMap *oldLevelMap, const BottomUpSubStructuredGraphSharedPtr &multiLevelGraph)
 Constructor for next level in multi-level static condensation. More...
 
virtual ~AssemblyMap ()
 Destructor. More...
 
LibUtilities::CommSharedPtr GetComm ()
 Retrieves the communicator. More...
 
std::string GetVariable ()
 Retrieves the variable string. More...
 
size_t GetHash () const
 Retrieves the hash of this map. More...
 
int GetLocalToGlobalMap (const int i) const
 
int GetGlobalToUniversalMap (const int i) const
 
int GetGlobalToUniversalMapUnique (const int i) const
 
const Array< OneD, const int > & GetLocalToGlobalMap ()
 
const Array< OneD, const int > & GetGlobalToUniversalMap ()
 
const Array< OneD, const int > & GetGlobalToUniversalMapUnique ()
 
NekDouble GetLocalToGlobalSign (const int i) const
 
const Array< OneD, NekDouble > & GetLocalToGlobalSign () const
 
void LocalToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm=true) const
 
void LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, bool useComm=true) const
 
void GlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
void GlobalToLocal (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
 
void Assemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
void Assemble (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
void UniversalAssemble (Array< OneD, NekDouble > &pGlobal) const
 
void UniversalAssemble (NekVector< NekDouble > &pGlobal) const
 
void UniversalAssemble (Array< OneD, NekDouble > &pGlobal, int offset) const
 
void UniversalAbsMaxBnd (Array< OneD, NekDouble > &bndvals)
 
int GetLocalToGlobalBndMap (const int i) const
 Retrieve the global index of a given local boundary mode. More...
 
const Array< OneD, const int > & GetLocalToGlobalBndMap ()
 Retrieve the global indices of the local boundary modes. More...
 
const Array< OneD, const int > & GetGlobalToUniversalBndMap ()
 
const Array< OneD, const int > & GetGlobalToUniversalBndMapUnique ()
 
bool GetSignChange ()
 Returns true if using a modal expansion requiring a change of sign of some modes. More...
 
NekDouble GetLocalToGlobalBndSign (const int i) const
 Retrieve the sign change of a given local boundary mode. More...
 
Array< OneD, const NekDoubleGetLocalToGlobalBndSign () const
 Retrieve the sign change for all local boundary modes. More...
 
const Array< OneD, const int > & GetBndCondCoeffsToLocalCoeffsMap ()
 Retrieves the local indices corresponding to the boundary expansion modes. More...
 
const Array< OneD, NekDouble > & GetBndCondCoeffsToLocalCoeffsSign ()
 Returns the modal sign associated with a given boundary expansion mode. More...
 
const Array< OneD, const int > & GetBndCondCoeffsToLocalTraceMap ()
 Retrieves the local indices corresponding to the boundary expansion modes to global trace. More...
 
int GetBndCondIDToGlobalTraceID (const int i)
 Returns the global index of the boundary trace giving the index on the boundary expansion. More...
 
const Array< OneD, const int > & GetBndCondIDToGlobalTraceID ()
 
int GetNumGlobalDirBndCoeffs () const
 Returns the number of global Dirichlet boundary coefficients. More...
 
int GetNumLocalDirBndCoeffs () const
 Returns the number of local Dirichlet boundary coefficients. More...
 
int GetNumGlobalBndCoeffs () const
 Returns the total number of global boundary coefficients. More...
 
int GetNumLocalBndCoeffs () const
 Returns the total number of local boundary coefficients. More...
 
int GetNumLocalCoeffs () const
 Returns the total number of local coefficients. More...
 
int GetNumGlobalCoeffs () const
 Returns the total number of global coefficients. More...
 
bool GetSingularSystem () const
 Retrieves if the system is singular (true) or not (false) More...
 
void GlobalToLocalBnd (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc, int offset) const
 
void GlobalToLocalBnd (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
 
void GlobalToLocalBnd (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc, int offset) const
 
void GlobalToLocalBnd (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
void LocalBndToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, int offset, bool UseComm=true) const
 
void LocalBndToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool UseComm=true) const
 
void LocalToLocalBnd (const Array< OneD, const NekDouble > &local, Array< OneD, NekDouble > &locbnd) const
 
void LocalToLocalInt (const Array< OneD, const NekDouble > &local, Array< OneD, NekDouble > &locint) const
 
void LocalBndToLocal (const Array< OneD, const NekDouble > &locbnd, Array< OneD, NekDouble > &local) const
 
void LocalIntToLocal (const Array< OneD, const NekDouble > &locbnd, Array< OneD, NekDouble > &local) const
 
void AssembleBnd (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
 
void AssembleBnd (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
void AssembleBnd (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, int offset) const
 
void AssembleBnd (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
void UniversalAssembleBnd (Array< OneD, NekDouble > &pGlobal) const
 
void UniversalAssembleBnd (NekVector< NekDouble > &pGlobal) const
 
void UniversalAssembleBnd (Array< OneD, NekDouble > &pGlobal, int offset) const
 
int GetFullSystemBandWidth () const
 
int GetNumNonDirVertexModes () const
 
int GetNumNonDirEdgeModes () const
 
int GetNumNonDirFaceModes () const
 
int GetNumDirEdges () const
 
int GetNumDirFaces () const
 
int GetNumNonDirEdges () const
 
int GetNumNonDirFaces () const
 
void PrintStats (std::ostream &out, std::string variable, bool printHeader=true) const
 
const Array< OneD, const int > & GetExtraDirEdges ()
 
std::shared_ptr< AssemblyMapLinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType)
 
int GetBndSystemBandWidth () const
 Returns the bandwidth of the boundary system. More...
 
int GetStaticCondLevel () const
 Returns the level of static condensation for this map. More...
 
int GetNumPatches () const
 Returns the number of patches in this static condensation level. More...
 
const Array< OneD, const unsigned int > & GetNumLocalBndCoeffsPerPatch ()
 Returns the number of local boundary coefficients in each patch. More...
 
const Array< OneD, const unsigned int > & GetNumLocalIntCoeffsPerPatch ()
 Returns the number of local interior coefficients in each patch. More...
 
const AssemblyMapSharedPtr GetNextLevelLocalToGlobalMap () const
 Returns the local to global mapping for the next level in the multi-level static condensation. More...
 
void SetNextLevelLocalToGlobalMap (AssemblyMapSharedPtr pNextLevelLocalToGlobalMap)
 
const PatchMapSharedPtrGetPatchMapFromPrevLevel (void) const
 Returns the patch map from the previous level of the multi-level static condensation. More...
 
bool AtLastLevel () const
 Returns true if this is the last level in the multi-level static condensation. More...
 
GlobalSysSolnType GetGlobalSysSolnType () const
 Returns the method of solving global systems. More...
 
std::string GetPreconType () const
 
NekDouble GetIterativeTolerance () const
 
bool IsAbsoluteTolerance () const
 
int GetMaxIterations () const
 
int GetSuccessiveRHS () const
 
std::string GetLinSysIterSolver () const
 
int GetLowestStaticCondLevel () const
 
void PatchLocalToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
void PatchGlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
void PatchAssemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 

Additional Inherited Members

- Protected Member Functions inherited from Nektar::MultiRegions::AssemblyMapCG
int CreateGraph (const ExpList &locExp, const BndCondExp &bndCondExp, const Array< OneD, const BndCond > &bndConditions, const bool checkIfSystemSingular, const PeriodicMap &periodicVerts, const PeriodicMap &periodicEdges, const PeriodicMap &periodicFaces, DofGraph &graph, BottomUpSubStructuredGraphSharedPtr &bottomUpGraph, std::set< int > &extraDirVerts, std::set< int > &extraDirEdges, int &firstNonDirGraphVertId, int &nExtraDirichlet, int mdswitch=1)
 
void SetUpUniversalC0ContMap (const ExpList &locExp, const PeriodicMap &perVerts=NullPeriodicMap, const PeriodicMap &perEdges=NullPeriodicMap, const PeriodicMap &perFaces=NullPeriodicMap)
 
void CalculateFullSystemBandWidth ()
 Calculate the bandwith of the full matrix system. More...
 
virtual int v_GetLocalToGlobalMap (const int i) const override
 
virtual int v_GetGlobalToUniversalMap (const int i) const override
 
virtual int v_GetGlobalToUniversalMapUnique (const int i) const override
 
virtual const Array< OneD, const int > & v_GetLocalToGlobalMap () override
 
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMap () override
 
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMapUnique () override
 
virtual NekDouble v_GetLocalToGlobalSign (const int i) const override
 
virtual const Array< OneD, NekDouble > & v_GetLocalToGlobalSign () const override
 
virtual void v_LocalToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm) const override
 
virtual void v_LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, bool useComm) const override
 
virtual void v_GlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const override
 
virtual void v_GlobalToLocal (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const override
 
virtual void v_Assemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const override
 
virtual void v_Assemble (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const override
 
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal) const override
 
virtual void v_UniversalAssemble (NekVector< NekDouble > &pGlobal) const override
 
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal, int offset) const override
 
virtual int v_GetFullSystemBandWidth () const override
 
virtual int v_GetNumNonDirVertexModes () const override
 
virtual int v_GetNumNonDirEdgeModes () const override
 
virtual int v_GetNumNonDirFaceModes () const override
 
virtual int v_GetNumDirEdges () const override
 
virtual int v_GetNumDirFaces () const override
 
virtual int v_GetNumNonDirEdges () const override
 
virtual int v_GetNumNonDirFaces () const override
 
virtual const Array< OneD, const int > & v_GetExtraDirEdges () override
 
virtual AssemblyMapSharedPtr v_LinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType) override
 Construct an AssemblyMapCG object which corresponds to the linear space of the current object. More...
 
- Protected Member Functions inherited from Nektar::MultiRegions::AssemblyMap
void CalculateBndSystemBandWidth ()
 Calculates the bandwidth of the boundary system. More...
 
void GlobalToLocalBndWithoutSign (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc)
 
virtual int v_GetLocalToGlobalMap (const int i) const
 
virtual int v_GetGlobalToUniversalMap (const int i) const
 
virtual int v_GetGlobalToUniversalMapUnique (const int i) const
 
virtual const Array< OneD, const int > & v_GetLocalToGlobalMap ()
 
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMap ()
 
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMapUnique ()
 
virtual NekDouble v_GetLocalToGlobalSign (const int i) const
 
virtual const Array< OneD, NekDouble > & v_GetLocalToGlobalSign () const
 
virtual void v_LocalToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm) const
 
virtual void v_LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, bool useComm) const
 
virtual void v_GlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
virtual void v_GlobalToLocal (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
 
virtual void v_Assemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
virtual void v_Assemble (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal) const
 
virtual void v_UniversalAssemble (NekVector< NekDouble > &pGlobal) const
 
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal, int offset) const
 
virtual int v_GetFullSystemBandWidth () const
 
virtual int v_GetNumNonDirVertexModes () const
 
virtual int v_GetNumNonDirEdgeModes () const
 
virtual int v_GetNumNonDirFaceModes () const
 
virtual int v_GetNumDirEdges () const
 
virtual int v_GetNumDirFaces () const
 
virtual int v_GetNumNonDirEdges () const
 
virtual int v_GetNumNonDirFaces () const
 
virtual const Array< OneD, const int > & v_GetExtraDirEdges ()
 
virtual std::shared_ptr< AssemblyMapv_LinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType)
 Generate a linear space mapping from existing mapping. More...
 
- Protected Attributes inherited from Nektar::MultiRegions::AssemblyMapCG
Array< OneD, int > m_localToGlobalMap
 Integer map of local coeffs to global space. More...
 
Array< OneD, NekDoublem_localToGlobalSign
 Integer sign of local coeffs to global space. More...
 
int m_fullSystemBandWidth
 Bandwith of the full matrix system (no static condensation). More...
 
Array< OneD, int > m_globalToUniversalMap
 Integer map of process coeffs to universal space. More...
 
Array< OneD, int > m_globalToUniversalMapUnique
 Integer map of unique process coeffs to universal space (signed) More...
 
int m_numNonDirVertexModes
 Number of non Dirichlet vertex modes. More...
 
int m_numNonDirEdgeModes
 Number of non Dirichlet edge modes. More...
 
int m_numNonDirFaceModes
 Number of non Dirichlet face modes. More...
 
int m_numDirEdges
 Number of Dirichlet edges. More...
 
int m_numDirFaces
 Number of Dirichlet faces. More...
 
int m_numNonDirEdges
 Number of Dirichlet edges. More...
 
int m_numNonDirFaces
 Number of Dirichlet faces. More...
 
int m_numLocalBndCondCoeffs
 Number of local boundary condition coefficients. More...
 
Array< OneD, int > m_extraDirEdges
 Extra dirichlet edges in parallel. More...
 
int m_numLocDirBndCondDofs
 Number of local boundary condition degrees of freedom. More...
 
int m_maxStaticCondLevel
 Maximum static condensation level. More...
 
std::set< ExtraDirDofm_copyLocalDirDofs
 Set indicating degrees of freedom which are Dirichlet but whose value is stored on another processor. More...
 
std::set< int > m_parallelDirBndSign
 Set indicating the local coeffs just touching parallel dirichlet boundary that have a sign change. More...
 
- Protected Attributes inherited from Nektar::MultiRegions::AssemblyMap
LibUtilities::SessionReaderSharedPtr m_session
 Session object. More...
 
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
std::string m_variable
 Variable string identifier. More...
 
size_t m_hash
 Hash for map. More...
 
int m_numLocalBndCoeffs
 Number of local boundary coefficients. More...
 
int m_numGlobalBndCoeffs
 Total number of global boundary coefficients. More...
 
int m_numLocalDirBndCoeffs
 Number of Local Dirichlet Boundary Coefficients. More...
 
int m_numGlobalDirBndCoeffs
 Number of Global Dirichlet Boundary Coefficients. More...
 
bool m_systemSingular
 Flag indicating if the system is singular or not. More...
 
int m_numLocalCoeffs
 Total number of local coefficients. More...
 
int m_numGlobalCoeffs
 Total number of global coefficients. More...
 
bool m_signChange
 Flag indicating if modes require sign reversal. More...
 
Array< OneD, int > m_localToGlobalBndMap
 Integer map of local coeffs to global Boundary Dofs. More...
 
Array< OneD, NekDoublem_localToGlobalBndSign
 Integer sign of local boundary coeffs to global space. More...
 
Array< OneD, int > m_localToLocalBndMap
 Integer map of local boundary coeffs to local boundary system numbering. More...
 
Array< OneD, int > m_localToLocalIntMap
 Integer map of local boundary coeffs to local interior system numbering. More...
 
Array< OneD, int > m_bndCondCoeffsToLocalCoeffsMap
 Integer map of bnd cond coeffs to local coefficients. More...
 
Array< OneD, NekDoublem_bndCondCoeffsToLocalCoeffsSign
 Integer map of sign of bnd cond coeffs to local coefficients. More...
 
Array< OneD, int > m_bndCondCoeffsToLocalTraceMap
 Integer map of bnd cond coeff to local trace coeff. More...
 
Array< OneD, int > m_bndCondIDToGlobalTraceID
 Integer map of bnd cond trace number to global trace number. More...
 
Array< OneD, int > m_globalToUniversalBndMap
 Integer map of process coeffs to universal space. More...
 
Array< OneD, int > m_globalToUniversalBndMapUnique
 Integer map of unique process coeffs to universal space (signed) More...
 
GlobalSysSolnType m_solnType
 The solution type of the global system. More...
 
int m_bndSystemBandWidth
 The bandwith of the global bnd system. More...
 
std::string m_preconType
 Type type of preconditioner to use in iterative solver. More...
 
int m_maxIterations
 Maximum iterations for iterative solver. More...
 
NekDouble m_iterativeTolerance
 Tolerance for iterative solver. More...
 
bool m_isAbsoluteTolerance
 
int m_successiveRHS
 sucessive RHS for iterative solver More...
 
std::string m_linSysIterSolver
 Iterative solver: Conjugate Gradient, GMRES. More...
 
Gs::gs_datam_gsh
 
Gs::gs_datam_bndGsh
 
Gs::gs_datam_dirBndGsh
 gs gather communication to impose Dirhichlet BCs. More...
 
int m_staticCondLevel
 The level of recursion in the case of multi-level static condensation. More...
 
int m_numPatches
 The number of patches (~elements) in the current level. More...
 
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
 The number of bnd dofs per patch. More...
 
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
 The number of int dofs per patch. More...
 
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
 Map from the patches of the previous level to the patches of the current level. More...
 
int m_lowestStaticCondLevel
 Lowest static condensation level. More...
 

Detailed Description

Definition at line 45 of file CoupledLocalToGlobalC0ContMap.h.

Constructor & Destructor Documentation

◆ CoupledLocalToGlobalC0ContMap()

Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr graph,
const SpatialDomains::BoundaryConditionsSharedPtr boundaryConditions,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
const MultiRegions::ExpListSharedPtr pressure,
const size_t  nz_loc,
const bool  CheckforSingularSys = true 
)

This is an vector extension of MultiRegions::AssemblyMapCG::SetUp2DExpansionC0ContMap related to the Linearised Navier Stokes problem

STEP 1: Wrap boundary conditions vector in an array (since routine is set up for multiple fields) and call the graph re-odering subroutine to obtain the reordered values

@TODO: Fix this so that we can extract normals from edges

STEP 2a: Set the mean pressure modes to edges depending on type of direct solver technique;

STEP 2: Count out the number of Dirichlet vertices and edges first

STEP 3: Set up an array which contains the offset information of the different graph vertices.

This basically means to identify how many global degrees of freedom the individual graph vertices correspond. Obviously, the graph vertices corresponding to the mesh-vertices account for a single global DOF. However, the graph vertices corresponding to the element edges correspond to 2*(N-2) global DOF where N is equal to the number of boundary modes on this edge.

STEP 4: Now, all ingredients are ready to set up the actual local to global mapping.

The remainder of the map consists of the element-interior degrees of freedom. This leads to the block-diagonal submatrix as each element-interior mode is globally orthogonal to modes in all other elements.

STEP 5: The boundary condition mapping is generated from the same vertex renumbering and fill in a unique interior map.

Definition at line 53 of file CoupledLocalToGlobalC0ContMap.cpp.

60 : AssemblyMapCG(pSession, fields[0]->GetComm())
61{
62 boost::ignore_unused(graph, boundaryConditions);
63
64 size_t i, j, k, n;
65 size_t cnt = 0, offset = 0;
66 size_t meshVertId;
67 size_t meshEdgeId;
68 size_t bndEdgeCnt;
69 size_t globalId;
70 size_t nEdgeCoeffs;
71 size_t nEdgeInteriorCoeffs;
72 int firstNonDirGraphVertId;
73 size_t nLocBndCondDofs = 0;
74 size_t nLocDirBndCondDofs = 0;
75 int nExtraDirichlet = 0;
79 StdRegions::Orientation edgeOrient;
80 Array<OneD, unsigned int> edgeInteriorMap;
81 Array<OneD, int> edgeInteriorSign;
82 size_t nvel = fields.size();
83
84 const LocalRegions::ExpansionVector &locExpVector = *(fields[0]->GetExp());
85 int id, diff;
86 size_t nel = fields[0]->GetNumElmts();
87
88 MultiRegions::PeriodicMap periodicVerts;
89 MultiRegions::PeriodicMap periodicEdges;
90 MultiRegions::PeriodicMap periodicFaces;
91 vector<map<int, int>> ReorderedGraphVertId(3);
93 int staticCondLevel = 0;
94
95 if (CheckforSingularSys) // all singularity checking by setting flag to true
96 {
97 m_systemSingular = true;
98 }
99 else // Turn off singular checking by setting flag to false
100 {
101 m_systemSingular = false;
102 }
103
104 /**
105 * STEP 1: Wrap boundary conditions vector in an array
106 * (since routine is set up for multiple fields) and call
107 * the graph re-odering subroutine to obtain the reordered
108 * values
109 */
110
111 // Obtain any periodic information and allocate default mapping array
112 fields[0]->GetPeriodicEntities(periodicVerts, periodicEdges, periodicFaces);
113
114 const Array<OneD, const MultiRegions::ExpListSharedPtr> bndCondExp =
115 fields[0]->GetBndCondExpansions();
116
117 Array<OneD, Array<OneD, const SpatialDomains::BoundaryConditionShPtr>>
118 bndConditionsVec(nvel);
119 for (i = 0; i < nvel; ++i)
120 {
121 bndConditionsVec[i] = fields[i]->GetBndConditions();
122 }
123
124 map<int, int> IsDirVertDof;
125 map<int, int> IsDirEdgeDof;
126
128 for (j = 0; j < bndCondExp.size(); ++j)
129 {
130 map<int, int> BndExpVids;
131 // collect unique list of vertex ids for this expansion
132 for (k = 0; k < bndCondExp[j]->GetNumElmts(); ++k)
133 {
134 g = bndCondExp[j]
135 ->GetExp(k)
136 ->as<LocalRegions::Expansion1D>()
137 ->GetGeom1D();
138 BndExpVids[g->GetVid(0)] = g->GetVid(0);
139 BndExpVids[g->GetVid(1)] = g->GetVid(1);
140 }
141
142 for (i = 0; i < nvel; ++i)
143 {
144 if (bndConditionsVec[i][j]->GetBoundaryConditionType() ==
146 {
147 // set number of Dirichlet conditions along edge
148 for (k = 0; k < bndCondExp[j]->GetNumElmts(); ++k)
149 {
150 IsDirEdgeDof[bndCondExp[j]
151 ->GetExp(k)
152 ->as<LocalRegions::Expansion1D>()
153 ->GetGeom1D()
154 ->GetGlobalID()] += 1;
155 }
156
157 // Set number of Dirichlet conditions at vertices
158 // with a clamp on its maximum value being nvel to
159 // handle corners between expansions
160 for (auto &mapIt : BndExpVids)
161 {
162 id = IsDirVertDof[mapIt.second] + 1;
163 IsDirVertDof[mapIt.second] = (id > (int)nvel) ? nvel : id;
164 }
165 }
166 else
167 {
168 // Check to see that edge normals have non-zero
169 // component in this direction since otherwise
170 // also can be singular.
171 /// @TODO: Fix this so that we can extract normals from edges
172 for (k = 0; k < bndCondExp[j]->GetNumElmts(); ++k)
173 {
174 Array<OneD, Array<OneD, NekDouble>> locnorm;
176 bndCondExp[j]
177 ->GetExp(k)
178 ->as<LocalRegions::Expansion1D>();
179 locnorm =
180 loc_exp->GetLeftAdjacentElementExp()->GetTraceNormal(
181 loc_exp->GetLeftAdjacentElementTrace());
182
183 size_t ndir = locnorm.size();
184 if (i < ndir) // account for Fourier version where n can be
185 // larger then ndir
186 {
187 for (size_t l = 0; l < locnorm[0].size(); ++l)
188 {
189 if (fabs(locnorm[i][l]) > NekConstants::kNekZeroTol)
190 {
191 m_systemSingular = false;
192 break;
193 }
194 }
195 }
196 if (m_systemSingular == false)
197 {
198 break;
199 }
200 }
201 }
202 }
203 }
204
205 Array<OneD, map<int, int>> Dofs(2);
206
207 Array<OneD, int> AddMeanPressureToEdgeId(nel, -1);
208 size_t edgeId, vertId;
209
210 // special case of singular problem - need to fix one pressure
211 // dof to a dirichlet edge. Since we attached pressure dof to
212 // last velocity component of edge need to make sure this
213 // component is Dirichlet
215 {
216 id = -1;
217 for (i = 0; i < bndConditionsVec[0].size(); ++i)
218 {
219 if (bndConditionsVec[nvel - 1][i]->GetBoundaryConditionType() ==
221 {
222 id = bndCondExp[i]
223 ->GetExp(0)
224 ->as<LocalRegions::Expansion1D>()
225 ->GetGeom1D()
226 ->GetGlobalID();
227 break;
228 }
229 }
230
231 ASSERTL0(id != -1, " Did not find an edge to attach singular pressure "
232 "degree of freedom");
233
234 // determine element with this edge id. There may be a
235 // more direct way of getting element from spatialDomains
236 for (i = 0; i < nel; ++i)
237 {
238 for (j = 0; j < (size_t)locExpVector[i]->GetNverts(); ++j)
239 {
240 edgeId = (locExpVector[i]
241 ->as<LocalRegions::Expansion2D>()
242 ->GetGeom2D())
243 ->GetEid(j);
244
245 if ((int)edgeId == id)
246 {
247 AddMeanPressureToEdgeId[i] = id;
248 break;
249 }
250 }
251
252 if (AddMeanPressureToEdgeId[i] != -1)
253 {
254 break;
255 }
256 }
257 }
258
259 for (i = 0; i < nel; ++i)
260 {
261 for (j = 0; j < (size_t)locExpVector[i]->GetNverts(); ++j)
262 {
263 vertId =
264 (locExpVector[i]->as<LocalRegions::Expansion2D>()->GetGeom2D())
265 ->GetVid(j);
266 if (Dofs[0].count(vertId) == 0)
267 {
268 Dofs[0][vertId] = nvel * nz_loc;
269
270 // Adjust for a Dirichlet boundary condition to give number to
271 // be solved
272 if (IsDirVertDof.count(vertId) != 0)
273 {
274 Dofs[0][vertId] -= IsDirVertDof[vertId] * nz_loc;
275 }
276 }
277
278 edgeId =
279 (locExpVector[i]->as<LocalRegions::Expansion2D>()->GetGeom2D())
280 ->GetEid(j);
281 if (Dofs[1].count(edgeId) == 0)
282 {
283 Dofs[1][edgeId] =
284 nvel * (locExpVector[i]->GetTraceNcoeffs(j) - 2) * nz_loc;
285 }
286
287 // Adjust for Dirichlet boundary conditions to give number to be
288 // solved
289 if (IsDirEdgeDof.count(edgeId) != 0)
290 {
291 Dofs[1][edgeId] -= IsDirEdgeDof[edgeId] * nz_loc *
292 (locExpVector[i]->GetTraceNcoeffs(j) - 2);
293 }
294 }
295 }
296
297 set<int> extraDirVerts, extraDirEdges;
298
299 CreateGraph(*fields[0], bndCondExp, bndConditionsVec, false, periodicVerts,
300 periodicEdges, periodicFaces, ReorderedGraphVertId,
301 bottomUpGraph, extraDirVerts, extraDirEdges,
302 firstNonDirGraphVertId, nExtraDirichlet, 4);
303 /*
304 SetUp2DGraphC0ContMap(*fields[0],
305 bndCondExp,
306 bndConditionsVec,
307 periodicVerts, periodicEdges,
308 Dofs, ReorderedGraphVertId,
309 firstNonDirGraphVertId, nExtraDirichlet,
310 bottomUpGraph, extraDir, false, 4);
311 */
312
313 /**
314 * STEP 2a: Set the mean pressure modes to edges depending on
315 * type of direct solver technique;
316 */
317
318 // determine which edge to add mean pressure dof based on
319 // ensuring that at least one pressure dof from an internal
320 // patch is associated with its boundary system
321 if (m_session->MatchSolverInfoAsEnum(
323 {
324
325 FindEdgeIdToAddMeanPressure(ReorderedGraphVertId, nel, locExpVector,
326 edgeId, vertId, firstNonDirGraphVertId,
327 IsDirEdgeDof, bottomUpGraph,
328 AddMeanPressureToEdgeId);
329 }
330
331 // Set unset elmts to non-Dirichlet edges.
332 // special case of singular problem - need to fix one
333 // pressure dof to a dirichlet edge
334 for (i = 0; i < nel; ++i)
335 {
336 for (j = 0; j < (size_t)locExpVector[i]->GetNverts(); ++j)
337 {
338 edgeId =
339 (locExpVector[i]->as<LocalRegions::Expansion2D>()->GetGeom2D())
340 ->GetEid(j);
341
342 if (IsDirEdgeDof.count(edgeId) == 0) // interior edge
343 {
344 // setup AddMeanPressureToEdgeId to decide where to
345 // put pressure
346 if (AddMeanPressureToEdgeId[i] == -1)
347 {
348 AddMeanPressureToEdgeId[i] = edgeId;
349 }
350 }
351 }
352 ASSERTL0((AddMeanPressureToEdgeId[i] != -1),
353 "Did not determine "
354 "an edge to attach mean pressure dof");
355 // Add the mean pressure degree of freedom to this edge
356 Dofs[1][AddMeanPressureToEdgeId[i]] += nz_loc;
357 }
358
359 map<int, int> pressureEdgeOffset;
360
361 /**
362 * STEP 2: Count out the number of Dirichlet vertices and edges first
363 */
364 for (i = 0; i < bndCondExp.size(); i++)
365 {
366 for (j = 0; j < bndCondExp[i]->GetNumElmts(); j++)
367 {
368 bndSegExp = bndCondExp[i]->GetExp(j)->as<LocalRegions::SegExp>();
369 for (k = 0; k < nvel; ++k)
370 {
371 if (bndConditionsVec[k][i]->GetBoundaryConditionType() ==
373 {
374 nLocDirBndCondDofs += bndSegExp->GetNcoeffs() * nz_loc;
375 }
376
377 if (bndConditionsVec[k][i]->GetBoundaryConditionType() !=
379 {
380 nLocBndCondDofs += bndSegExp->GetNcoeffs() * nz_loc;
381 }
382 }
383 }
384 }
385
387 {
388 m_numLocalDirBndCoeffs = nLocDirBndCondDofs + nExtraDirichlet + nz_loc;
389 }
390 else
391 {
392 m_numLocalDirBndCoeffs = nLocDirBndCondDofs + nExtraDirichlet;
393 }
394
395 /**
396 * STEP 3: Set up an array which contains the offset information of
397 * the different graph vertices.
398 *
399 * This basically means to identify how many global degrees of
400 * freedom the individual graph vertices correspond. Obviously,
401 * the graph vertices corresponding to the mesh-vertices account
402 * for a single global DOF. However, the graph vertices
403 * corresponding to the element edges correspond to 2*(N-2) global DOF
404 * where N is equal to the number of boundary modes on this edge.
405 */
406 Array<OneD, int> graphVertOffset(
407 nvel * nz_loc *
408 (ReorderedGraphVertId[0].size() + ReorderedGraphVertId[1].size()),
409 0);
410 graphVertOffset[0] = 0;
411
412 m_signChange = false;
413
414 for (i = 0; i < nel; ++i)
415 {
416 locExpansion = locExpVector[i]->as<LocalRegions::Expansion2D>();
417
418 for (j = 0; j < (size_t)locExpansion->GetNtraces(); ++j)
419 {
420 nEdgeCoeffs = locExpansion->GetTraceNcoeffs(j);
421 meshEdgeId = (locExpansion->GetGeom2D())->GetEid(j);
422 meshVertId = (locExpansion->GetGeom2D())->GetVid(j);
423
424 for (k = 0; k < nvel * nz_loc; ++k)
425 {
426 graphVertOffset[ReorderedGraphVertId[0][meshVertId] * nvel *
427 nz_loc +
428 k] = 1;
429 graphVertOffset[ReorderedGraphVertId[1][meshEdgeId] * nvel *
430 nz_loc +
431 k] = (nEdgeCoeffs - 2);
432 }
433
434 bType = locExpansion->GetBasisType(0);
435 // need a sign vector for modal expansions if nEdgeCoeffs >=4
436 if ((nEdgeCoeffs >= 4) && (bType == LibUtilities::eModified_A))
437 {
438 m_signChange = true;
439 }
440 }
441 }
442
443 // Add mean pressure modes;
444 for (i = 0; i < nel; ++i)
445 {
446 graphVertOffset[(ReorderedGraphVertId[1][AddMeanPressureToEdgeId[i]] +
447 1) *
448 nvel * nz_loc -
449 1] += nz_loc;
450 // graphVertOffset[(ReorderedGraphVertId[1][AddMeanPressureToEdgeId[i]])*nvel*nz_loc]
451 // += nz_loc;
452 }
453
454 // Negate the vertices and edges with only a partial
455 // Dirichlet conditon. Essentially we check to see if an edge
456 // has a mixed Dirichlet with Neumann/Robin Condition and if
457 // so negate the offset associated with this vertex.
458
459 map<int, int> DirVertChk;
460
461 for (i = 0; i < bndConditionsVec[0].size(); ++i)
462 {
463 cnt = 0;
464 for (j = 0; j < nvel; ++j)
465 {
466 if (bndConditionsVec[j][i]->GetBoundaryConditionType() ==
468 {
469 cnt++;
470 }
471 }
472
473 // Case where partial Dirichlet boundary condition
474 if ((cnt > 0) && (cnt < nvel))
475 {
476 for (j = 0; j < nvel; ++j)
477 {
478 if (bndConditionsVec[j][i]->GetBoundaryConditionType() ==
480 {
481 // negate graph offsets which should be
482 // Dirichlet conditions
483 for (k = 0; k < bndCondExp[i]->GetNumElmts(); ++k)
484 {
485 // vertices with mix condition;
486 id = bndCondExp[i]
487 ->GetExp(k)
488 ->as<LocalRegions::Expansion1D>()
489 ->GetGeom1D()
490 ->GetVid(0);
491 if (DirVertChk.count(id * nvel + j) == 0)
492 {
493 DirVertChk[id * nvel + j] = 1;
494 for (n = 0; n < nz_loc; ++n)
495 {
496 graphVertOffset[ReorderedGraphVertId[0][id] *
497 nvel * nz_loc +
498 j * nz_loc + n] *= -1;
499 }
500 }
501
502 id = bndCondExp[i]
503 ->GetExp(k)
504 ->as<LocalRegions::Expansion1D>()
505 ->GetGeom1D()
506 ->GetVid(1);
507 if (DirVertChk.count(id * nvel + j) == 0)
508 {
509 DirVertChk[id * nvel + j] = 1;
510 for (n = 0; n < nz_loc; ++n)
511 {
512 graphVertOffset[ReorderedGraphVertId[0][id] *
513 nvel * nz_loc +
514 j * nz_loc + n] *= -1;
515 }
516 }
517
518 // edges with mixed id;
519 id = bndCondExp[i]
520 ->GetExp(k)
521 ->as<LocalRegions::Expansion1D>()
522 ->GetGeom1D()
523 ->GetGlobalID();
524 for (n = 0; n < nz_loc; ++n)
525 {
526 graphVertOffset[ReorderedGraphVertId[1][id] * nvel *
527 nz_loc +
528 j * nz_loc + n] *= -1;
529 }
530 }
531 }
532 }
533 }
534 }
535
536 cnt = 0;
537 // assemble accumulative list of full Dirichlet values.
538 for (i = 0; i < firstNonDirGraphVertId * nvel * nz_loc; ++i)
539 {
540 diff = abs(graphVertOffset[i]);
541 graphVertOffset[i] = cnt;
542 cnt += diff;
543 }
544
545 // set Dirichlet values with negative values to Dirichlet value
546 for (i = firstNonDirGraphVertId * nvel * nz_loc; i < graphVertOffset.size();
547 ++i)
548 {
549 if (graphVertOffset[i] < 0)
550 {
551 diff = -graphVertOffset[i];
552 graphVertOffset[i] = -cnt;
553 cnt += diff;
554 }
555 }
556
557 // Accumulate all interior degrees of freedom with positive values
559
560 // offset values
561 for (i = firstNonDirGraphVertId * nvel * nz_loc; i < graphVertOffset.size();
562 ++i)
563 {
564 if (graphVertOffset[i] >= 0)
565 {
566 diff = graphVertOffset[i];
567 graphVertOffset[i] = cnt;
568 cnt += diff;
569 }
570 }
571
572 // Finally set negative entries (corresponding to Dirichlet
573 // values ) to be positive
574 for (i = firstNonDirGraphVertId * nvel * nz_loc; i < graphVertOffset.size();
575 ++i)
576 {
577 if (graphVertOffset[i] < 0)
578 {
579 graphVertOffset[i] = -graphVertOffset[i];
580 }
581 }
582
583 // Allocate the proper amount of space for the class-data and fill
584 // information that is already known
585 cnt = 0;
588
589 for (i = 0; i < nel; ++i)
590 {
592 nz_loc * (nvel * locExpVector[i]->NumBndryCoeffs() + 1);
593 // add these coeffs up separately since
594 // pressure->GetNcoeffs can include the coefficient in
595 // multiple planes.
596 m_numLocalCoeffs += (pressure->GetExp(i)->GetNcoeffs() - 1) * nz_loc;
597 }
598
600
601 m_localToGlobalMap = Array<OneD, int>(m_numLocalCoeffs, -1);
602 m_localToGlobalBndMap = Array<OneD, int>(m_numLocalBndCoeffs, -1);
603 m_bndCondIDToGlobalTraceID = Array<OneD, int>(nLocBndCondDofs, -1);
604
605 // Set default sign array.
606 m_localToGlobalSign = Array<OneD, NekDouble>(m_numLocalCoeffs, 1.0);
607 m_localToGlobalBndSign = Array<OneD, NekDouble>(m_numLocalBndCoeffs, 1.0);
608
609 m_staticCondLevel = staticCondLevel;
610 m_numPatches = nel;
611
612 m_numLocalBndCoeffsPerPatch = Array<OneD, unsigned int>(nel);
613 m_numLocalIntCoeffsPerPatch = Array<OneD, unsigned int>(nel);
614
615 for (i = 0; i < nel; ++i)
616 {
618 (unsigned int)nz_loc *
619 (nvel * locExpVector[i]->NumBndryCoeffs() + 1);
621 (unsigned int)nz_loc * (pressure->GetExp(i)->GetNcoeffs() - 1);
622 }
623
624 // Set up local to local bnd and local int maps
625 m_localToLocalBndMap = Array<OneD, int>(m_numLocalBndCoeffs, -1);
627 Array<OneD, int>(m_numLocalCoeffs - m_numLocalBndCoeffs, -1);
628
629 size_t bndcnt = 0;
630 size_t intcnt = 0;
631 cnt = 0;
632 for (i = 0; i < nel; ++i)
633 {
634 for (j = 0; j < nz_loc * (nvel * locExpVector[i]->NumBndryCoeffs());
635 ++j)
636 {
637 m_localToLocalBndMap[bndcnt++] = cnt++;
638 }
639
640 for (n = 0; n < nz_loc; ++n)
641 {
642 m_localToLocalBndMap[bndcnt++] = cnt++;
643 for (j = 1; j < (size_t)pressure->GetExp(i)->GetNcoeffs(); ++j)
644 {
645 m_localToLocalIntMap[intcnt++] = cnt++;
646 }
647 }
648 }
649
650 /**
651 * STEP 4: Now, all ingredients are ready to set up the actual
652 * local to global mapping.
653 *
654 * The remainder of the map consists of the element-interior
655 * degrees of freedom. This leads to the block-diagonal submatrix
656 * as each element-interior mode is globally orthogonal to modes
657 * in all other elements.
658 */
659 cnt = 0;
660 size_t nv, velnbndry;
661 Array<OneD, unsigned int> bmap;
662
663 // Loop over all the elements in the domain in shuffled
664 // ordering (element type consistency)
665 for (i = 0; i < nel; ++i)
666 {
667 locExpansion = locExpVector[i]->as<LocalRegions::Expansion2D>();
668
669 velnbndry = locExpansion->NumBndryCoeffs();
670
671 // Require an inverse ordering of the bmap system to store
672 // local numbering system. Therefore get hold of elemental
673 // bmap and set up an inverse map
674 map<int, int> inv_bmap;
675 locExpansion->GetBoundaryMap(bmap);
676 for (j = 0; j < bmap.size(); ++j)
677 {
678 inv_bmap[bmap[j]] = j;
679 }
680
681 // Loop over all edges (and vertices) of element i
682 for (j = 0; j < (size_t)locExpansion->GetNtraces(); ++j)
683 {
684 nEdgeInteriorCoeffs = locExpansion->GetTraceNcoeffs(j) - 2;
685 edgeOrient = (locExpansion->GetGeom())->GetEorient(j);
686 meshEdgeId = (locExpansion->GetGeom())->GetEid(j);
687 meshVertId = (locExpansion->GetGeom())->GetVid(j);
688
689 auto pIt = periodicEdges.find(meshEdgeId);
690
691 // See if this edge is periodic. If it is, then we map all
692 // edges to the one with lowest ID, and align all
693 // coefficients to this edge orientation.
694 if (pIt != periodicEdges.end())
695 {
696 pair<int, StdRegions::Orientation> idOrient =
697 DeterminePeriodicEdgeOrientId(meshEdgeId, edgeOrient,
698 pIt->second);
699 edgeOrient = idOrient.second;
700 }
701
702 locExpansion->GetTraceInteriorToElementMap(
703 j, edgeInteriorMap, edgeInteriorSign, edgeOrient);
704
705 // Set the global DOF for vertex j of element i
706 for (nv = 0; nv < nvel * nz_loc; ++nv)
707 {
708 m_localToGlobalMap[cnt + nv * velnbndry +
709 inv_bmap[locExpansion->GetVertexMap(j)]] =
710 graphVertOffset[ReorderedGraphVertId[0][meshVertId] * nvel *
711 nz_loc +
712 nv];
713
714 // Set the global DOF's for the interior modes of edge j
715 for (k = 0; k < nEdgeInteriorCoeffs; ++k)
716 {
717 m_localToGlobalMap[cnt + nv * velnbndry +
718 inv_bmap[edgeInteriorMap[k]]] =
719 graphVertOffset[ReorderedGraphVertId[1][meshEdgeId] *
720 nvel * nz_loc +
721 nv] +
722 k;
723 }
724 }
725
726 // Fill the sign vector if required
727 if (m_signChange)
728 {
729 for (nv = 0; nv < nvel * nz_loc; ++nv)
730 {
731 for (k = 0; k < nEdgeInteriorCoeffs; ++k)
732 {
733 m_localToGlobalSign[cnt + nv * velnbndry +
734 inv_bmap[edgeInteriorMap[k]]] =
735 (NekDouble)edgeInteriorSign[k];
736 }
737 }
738 }
739 }
740
741 // use difference between two edges of the
742 // AddMeanPressureEdgeId to det nEdgeInteriorCoeffs.
743 nEdgeInteriorCoeffs =
744 graphVertOffset[(ReorderedGraphVertId[1]
745 [AddMeanPressureToEdgeId[i]]) *
746 nvel * nz_loc +
747 1] -
748 graphVertOffset[(ReorderedGraphVertId[1]
749 [AddMeanPressureToEdgeId[i]]) *
750 nvel * nz_loc];
751
752 size_t psize = pressure->GetExp(i)->GetNcoeffs();
753 for (n = 0; n < nz_loc; ++n)
754 {
755 m_localToGlobalMap[cnt + nz_loc * nvel * velnbndry + n * psize] =
756 graphVertOffset
757 [(ReorderedGraphVertId[1][AddMeanPressureToEdgeId[i]] + 1) *
758 nvel * nz_loc -
759 1] +
760 nEdgeInteriorCoeffs +
761 pressureEdgeOffset[AddMeanPressureToEdgeId[i]];
762
763 pressureEdgeOffset[AddMeanPressureToEdgeId[i]] += 1;
764 }
765
766 cnt += (velnbndry * nvel + psize) * nz_loc;
767 }
768
769 // Set up the mapping for the boundary conditions
770 offset = cnt = 0;
771 for (nv = 0; nv < nvel; ++nv)
772 {
773 for (i = 0; i < bndCondExp.size(); i++)
774 {
775 if (bndConditionsVec[nv][i]->GetBoundaryConditionType() ==
777 {
778 continue;
779 }
780
781 for (n = 0; n < nz_loc; ++n)
782 {
783 int ncoeffcnt = 0;
784 for (j = 0; j < bndCondExp[i]->GetNumElmts(); j++)
785 {
786 bndSegExp =
787 bndCondExp[i]->GetExp(j)->as<LocalRegions::SegExp>();
788
789 cnt = offset + bndCondExp[i]->GetCoeff_Offset(j);
790 for (k = 0; k < 2; k++)
791 {
792 meshVertId = (bndSegExp->GetGeom1D())->GetVid(k);
794 bndSegExp->GetVertexMap(k)] =
795 graphVertOffset[ReorderedGraphVertId[0]
796 [meshVertId] *
797 nvel * nz_loc +
798 nv * nz_loc + n];
799 }
800
801 meshEdgeId = (bndSegExp->GetGeom1D())->GetGlobalID();
802 bndEdgeCnt = 0;
803 nEdgeCoeffs = bndSegExp->GetNcoeffs();
804 for (k = 0; k < nEdgeCoeffs; k++)
805 {
806 if (m_bndCondIDToGlobalTraceID[cnt + k] == -1)
807 {
809 graphVertOffset
810 [ReorderedGraphVertId[1][meshEdgeId] *
811 nvel * nz_loc +
812 nv * nz_loc + n] +
813 bndEdgeCnt;
814 bndEdgeCnt++;
815 }
816 }
817 ncoeffcnt += nEdgeCoeffs;
818 }
819 // Note: Can not use bndCondExp[i]->GetNcoeffs()
820 // due to homogeneous extension not returning just
821 // the value per plane
822 offset += ncoeffcnt;
823 }
824 }
825 }
826
827 globalId = Vmath::Vmax(m_numLocalCoeffs, &m_localToGlobalMap[0], 1) + 1;
828 m_numGlobalBndCoeffs = globalId;
829
830 /**
831 * STEP 5: The boundary condition mapping is generated from the
832 * same vertex renumbering and fill in a unique interior map.
833 */
834 cnt = 0;
835 for (i = 0; i < (size_t)m_numLocalCoeffs; ++i)
836 {
837 if (m_localToGlobalMap[i] == -1)
838 {
839 m_localToGlobalMap[i] = globalId++;
840 }
841 else
842 {
843 if (m_signChange)
844 {
846 }
848 }
849 }
850 m_numGlobalCoeffs = globalId;
851
852 // Set up the local to global map for the next level when using
853 // multi-level static condensation
854 if (m_session->MatchSolverInfoAsEnum(
856 {
857 if (m_staticCondLevel < (bottomUpGraph->GetNlevels() - 1))
858 {
859 Array<OneD, int> vwgts_perm(Dofs[0].size() + Dofs[1].size() -
860 firstNonDirGraphVertId);
861 for (i = 0; i < locExpVector.size(); ++i)
862 {
863 locExpansion = locExpVector[i]->as<LocalRegions::Expansion2D>();
864 size_t nVert = locExpansion->GetNverts();
865 for (j = 0; j < nVert; ++j)
866 {
867 meshEdgeId = (locExpansion->GetGeom2D())->GetEid(j);
868 meshVertId = (locExpansion->GetGeom2D())->GetVid(j);
869
870 if (ReorderedGraphVertId[0][meshVertId] >=
871 firstNonDirGraphVertId)
872 {
873 vwgts_perm[ReorderedGraphVertId[0][meshVertId] -
874 firstNonDirGraphVertId] =
875 Dofs[0][meshVertId];
876 }
877
878 if (ReorderedGraphVertId[1][meshEdgeId] >=
879 firstNonDirGraphVertId)
880 {
881 vwgts_perm[ReorderedGraphVertId[1][meshEdgeId] -
882 firstNonDirGraphVertId] =
883 Dofs[1][meshEdgeId];
884 }
885 }
886 }
887
888 bottomUpGraph->ExpandGraphWithVertexWeights(vwgts_perm);
889
892 bottomUpGraph);
893 }
894 }
895}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
void FindEdgeIdToAddMeanPressure(std::vector< std::map< int, int > > &ReorderedGraphVertId, size_t &nel, const LocalRegions::ExpansionVector &locExpVector, size_t &edgeId, size_t &vertId, int &firstNonDirGraphVertId, std::map< int, int > &IsDirEdgeDof, MultiRegions::BottomUpSubStructuredGraphSharedPtr &bottomUpGraph, Array< OneD, int > &AddMeanPressureToEdgeId)
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
int CreateGraph(const ExpList &locExp, const BndCondExp &bndCondExp, const Array< OneD, const BndCond > &bndConditions, const bool checkIfSystemSingular, const PeriodicMap &periodicVerts, const PeriodicMap &periodicEdges, const PeriodicMap &periodicFaces, DofGraph &graph, BottomUpSubStructuredGraphSharedPtr &bottomUpGraph, std::set< int > &extraDirVerts, std::set< int > &extraDirEdges, int &firstNonDirGraphVertId, int &nExtraDirichlet, int mdswitch=1)
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
AssemblyMapCG(const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &comm, const std::string variable="DefaultVar")
Default constructor.
Array< OneD, NekDouble > m_localToGlobalSign
Integer sign of local coeffs to global space.
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:367
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:381
Array< OneD, int > m_localToLocalIntMap
Integer map of local boundary coeffs to local interior system numbering.
Definition: AssemblyMap.h:390
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:378
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:386
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:436
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:336
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:348
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:443
int m_staticCondLevel
The level of recursion in the case of multi-level static condensation.
Definition: AssemblyMap.h:432
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:352
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:354
LibUtilities::CommSharedPtr GetComm()
Retrieves the communicator.
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
The number of int dofs per patch.
Definition: AssemblyMap.h:438
bool m_systemSingular
Flag indicating if the system is singular or not.
Definition: AssemblyMap.h:356
Array< OneD, int > m_localToGlobalBndMap
Integer map of local coeffs to global Boundary Dofs.
Definition: AssemblyMap.h:384
Array< OneD, int > m_localToLocalBndMap
Integer map of local boundary coeffs to local boundary system numbering.
Definition: AssemblyMap.h:388
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:434
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:350
Array< OneD, int > m_bndCondIDToGlobalTraceID
Integer map of bnd cond trace number to global trace number.
Definition: AssemblyMap.h:398
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:50
std::shared_ptr< SegExp > SegExpSharedPtr
Definition: SegExp.h:251
std::shared_ptr< Expansion2D > Expansion2DSharedPtr
Definition: Expansion1D.h:48
std::shared_ptr< Expansion1D > Expansion1DSharedPtr
Definition: Expansion1D.h:52
std::vector< ExpansionSharedPtr > ExpansionVector
Definition: Expansion.h:70
std::shared_ptr< BottomUpSubStructuredGraph > BottomUpSubStructuredGraphSharedPtr
pair< int, StdRegions::Orientation > DeterminePeriodicEdgeOrientId(int meshEdgeId, StdRegions::Orientation edgeOrient, const vector< PeriodicEntity > &periodicEdges)
Determine orientation of an edge to its periodic equivalents, as well as the ID of the representative...
std::map< int, std::vector< PeriodicEntity > > PeriodicMap
static const NekDouble kNekZeroTol
std::shared_ptr< Geometry1D > Geometry1DSharedPtr
Definition: Geometry.h:64
double NekDouble
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.cpp:940
scalarT< T > abs(scalarT< T > in)
Definition: scalar.hpp:298

References tinysimd::abs(), Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::MultiRegions::AssemblyMapCG::CreateGraph(), Nektar::MultiRegions::DeterminePeriodicEdgeOrientId(), Nektar::MultiRegions::eDirectMultiLevelStaticCond, Nektar::SpatialDomains::eDirichlet, Nektar::LibUtilities::eModified_A, Nektar::SpatialDomains::ePeriodic, FindEdgeIdToAddMeanPressure(), Nektar::StdRegions::StdExpansion::GetNcoeffs(), Nektar::StdRegions::StdExpansion::GetNverts(), Nektar::StdRegions::StdExpansion::GetTraceNcoeffs(), Nektar::StdRegions::StdExpansion::GetVertexMap(), Nektar::NekConstants::kNekZeroTol, Nektar::MultiRegions::AssemblyMap::m_bndCondIDToGlobalTraceID, Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndMap, Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndSign, Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap, Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign, Nektar::MultiRegions::AssemblyMap::m_localToLocalBndMap, Nektar::MultiRegions::AssemblyMap::m_localToLocalIntMap, Nektar::MultiRegions::AssemblyMap::m_nextLevelLocalToGlobalMap, Nektar::MultiRegions::AssemblyMap::m_numGlobalBndCoeffs, Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs, Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffsPerPatch, Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalDirBndCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalIntCoeffsPerPatch, Nektar::MultiRegions::AssemblyMap::m_numPatches, Nektar::MultiRegions::AssemblyMap::m_session, Nektar::MultiRegions::AssemblyMap::m_signChange, Nektar::MultiRegions::AssemblyMap::m_staticCondLevel, Nektar::MultiRegions::AssemblyMap::m_systemSingular, Nektar::StdRegions::StdExpansion::NumBndryCoeffs(), CG_Iterations::pressure, and Vmath::Vmax().

Member Function Documentation

◆ FindEdgeIdToAddMeanPressure()

void Nektar::CoupledLocalToGlobalC0ContMap::FindEdgeIdToAddMeanPressure ( std::vector< std::map< int, int > > &  ReorderedGraphVertId,
size_t &  nel,
const LocalRegions::ExpansionVector locExpVector,
size_t &  edgeId,
size_t &  vertId,
int &  firstNonDirGraphVertId,
std::map< int, int > &  IsDirEdgeDof,
MultiRegions::BottomUpSubStructuredGraphSharedPtr bottomUpGraph,
Array< OneD, int > &  AddMeanPressureToEdgeId 
)

Definition at line 897 of file CoupledLocalToGlobalC0ContMap.cpp.

903{
904 boost::ignore_unused(vertId);
905
906 size_t i, j, k;
907
908 // Make list of homogeneous graph edges to elmt mappings
909 Array<TwoD, int> EdgeIdToElmts(ReorderedGraphVertId[1].size(), 2, -1);
910 map<int, int> HomGraphEdgeIdToEdgeId;
911
912 for (i = 0; i < nel; ++i)
913 {
914 size_t nVert = locExpVector[i]->GetNverts();
915 for (j = 0; j < nVert; ++j)
916 {
917 edgeId =
918 (locExpVector[i]->as<LocalRegions::Expansion2D>()->GetGeom2D())
919 ->GetEid(j);
920
921 // note second condition stops us using mixed boundary condition
922 if ((ReorderedGraphVertId[1][edgeId] >= firstNonDirGraphVertId) &&
923 (IsDirEdgeDof.count(edgeId) == 0))
924 {
925 HomGraphEdgeIdToEdgeId[ReorderedGraphVertId[1][edgeId] -
926 firstNonDirGraphVertId] = edgeId;
927
928 if (EdgeIdToElmts[edgeId][0] == -1)
929 {
930 EdgeIdToElmts[edgeId][0] = i;
931 }
932 else
933 {
934 EdgeIdToElmts[edgeId][1] = i;
935 }
936 }
937 }
938 }
939
940 // Start at second to last level and find edge on boundary
941 // to attach element
942 size_t nlevels = bottomUpGraph->GetNlevels();
943
944 // determine a default edge to attach pressure modes to
945 // which is part of the inner solve;
946 int defedge = -1;
947
948 vector<MultiRegions::SubGraphSharedPtr> bndgraphs =
949 bottomUpGraph->GetInteriorBlocks(nlevels);
950 for (i = 0; i < bndgraphs.size(); ++i)
951 {
952 int GlobIdOffset = bndgraphs[i]->GetIdOffset();
953 size_t nVert = bndgraphs[i]->GetNverts();
954
955 for (j = 0; j < nVert; ++j)
956 {
957 // find edge in graph vert list
958 if (HomGraphEdgeIdToEdgeId.count(GlobIdOffset + j) != 0)
959 {
960 edgeId = HomGraphEdgeIdToEdgeId[GlobIdOffset + j];
961
962 if (defedge == -1)
963 {
964 defedge = edgeId;
965 break;
966 }
967 }
968 }
969 if (defedge != -1)
970 {
971 break;
972 }
973 }
974
975 for (size_t n = 1; n < nlevels; ++n)
976 {
977 // produce a map with a key that is the element id
978 // that contains which next level patch it belongs to
979 vector<MultiRegions::SubGraphSharedPtr> bndgraphs =
980 bottomUpGraph->GetInteriorBlocks(n + 1);
981
982 // Fill next level graph of adjacent elements and their level
983 map<int, int> ElmtInBndry;
984
985 for (i = 0; i < bndgraphs.size(); ++i)
986 {
987 int GlobIdOffset = bndgraphs[i]->GetIdOffset();
988 size_t nVert = bndgraphs[i]->GetNverts();
989
990 for (j = 0; j < nVert; ++j)
991 {
992 // find edge in graph vert list
993 if (HomGraphEdgeIdToEdgeId.count(GlobIdOffset + j) != 0)
994 {
995 edgeId = HomGraphEdgeIdToEdgeId[GlobIdOffset + j];
996
997 if (EdgeIdToElmts[edgeId][0] != -1)
998 {
999 ElmtInBndry[EdgeIdToElmts[edgeId][0]] = i;
1000 }
1001 if (EdgeIdToElmts[edgeId][1] != -1)
1002 {
1003 ElmtInBndry[EdgeIdToElmts[edgeId][1]] = i;
1004 }
1005 }
1006 }
1007 }
1008
1009 // Now search interior patches in this level for edges
1010 // that share the same element as a boundary edge and
1011 // assign this elmt that boundary edge
1012 vector<MultiRegions::SubGraphSharedPtr> intgraphs =
1013 bottomUpGraph->GetInteriorBlocks(n);
1014 for (i = 0; i < intgraphs.size(); ++i)
1015 {
1016 int GlobIdOffset = intgraphs[i]->GetIdOffset();
1017 bool SetEdge = false;
1018 int elmtid = 0;
1019 size_t nVert = intgraphs[i]->GetNverts();
1020 for (j = 0; j < nVert; ++j)
1021 {
1022 // Check to see if graph vert is an edge
1023 if (HomGraphEdgeIdToEdgeId.count(GlobIdOffset + j) != 0)
1024 {
1025 edgeId = HomGraphEdgeIdToEdgeId[GlobIdOffset + j];
1026
1027 for (k = 0; k < 2; ++k)
1028 {
1029 // relevant edge id
1030 elmtid = EdgeIdToElmts[edgeId][k];
1031
1032 if (elmtid != -1)
1033 {
1034 auto mapIt = ElmtInBndry.find(elmtid);
1035
1036 if (mapIt != ElmtInBndry.end())
1037 {
1038 // now find a edge in the next level boundary
1039 // graph
1040 int GlobIdOffset1 =
1041 bndgraphs[mapIt->second]->GetIdOffset();
1042 for (int l = 0;
1043 l < bndgraphs[mapIt->second]->GetNverts();
1044 ++l)
1045 {
1046 // find edge in graph vert list
1047 if (HomGraphEdgeIdToEdgeId.count(
1048 GlobIdOffset1 + l) != 0)
1049 {
1050 // June 2012: commenting this condition
1051 // apparently solved the bug caused by
1052 // the edge reordering procedure
1053
1054 // if(AddMeanPressureToEdgeId[elmtid] ==
1055 // -1)
1056 //{
1057
1058 // AddMeanPressureToEdgeId[elmtid] =
1059 // HomGraphEdgeIdToEdgeId[GlobIdOffset1+l];
1060 AddMeanPressureToEdgeId[elmtid] =
1061 defedge;
1062
1063 //}
1064 SetEdge = true;
1065 break;
1066 }
1067 }
1068 }
1069 }
1070 }
1071 }
1072 }
1073
1074 // if we have failed to find matching edge in next
1075 // level patch boundary then set last found elmt
1076 // associated to this interior patch to the
1077 // default edget value
1078 if (SetEdge == false)
1079 {
1080 if (elmtid == -1) // find an elmtid in patch
1081 {
1082 size_t nVert = intgraphs[i]->GetNverts();
1083 for (j = 0; j < nVert; ++j)
1084 {
1085 if (HomGraphEdgeIdToEdgeId.count(GlobIdOffset + j) != 0)
1086 {
1087 edgeId = HomGraphEdgeIdToEdgeId[GlobIdOffset + j];
1088 for (k = 0; k < 2; ++k)
1089 {
1090 // relevant edge id
1091 elmtid = EdgeIdToElmts[edgeId][k];
1092 if (elmtid != -1)
1093 {
1094 break;
1095 }
1096 }
1097 }
1098 if (elmtid != -1)
1099 {
1100 break;
1101 }
1102 }
1103 }
1104 if (AddMeanPressureToEdgeId[elmtid] == -1)
1105 {
1106 AddMeanPressureToEdgeId[elmtid] = defedge;
1107 }
1108 }
1109 }
1110 }
1111}

Referenced by CoupledLocalToGlobalC0ContMap().