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

Constructs mappings for the C0 scalar continuous Galerkin formulation. More...

#include <AssemblyMapCG.h>

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

Public Member Functions

 AssemblyMapCG (const LibUtilities::SessionReaderSharedPtr &pSession, 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::map< int, std::vector
< ExtraDirDof > > & 		GetExtraDirDofs ()
 
- Public Member Functions inherited from Nektar::MultiRegions::AssemblyMap
 AssemblyMap ()
 Default constructor. More...
 
 AssemblyMap (const LibUtilities::SessionReaderSharedPtr &pSession, 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...
 
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) const
 
void LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) 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
 
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 NekDouble > GetLocalToGlobalBndSign () const
 Retrieve the sign change for all local boundary modes. More...
 
int GetBndCondCoeffsToGlobalCoeffsMap (const int i)
 Retrieves the global index corresponding to a boundary expansion mode. More...
 
const Array< OneD, const int > & GetBndCondCoeffsToGlobalCoeffsMap ()
 Retrieves the global indices corresponding to the boundary expansion modes. More...
 
NekDouble GetBndCondCoeffsToGlobalCoeffsSign (const int i)
 Returns the modal sign associated with a given boundary expansion mode. More...
 
int GetBndCondTraceToGlobalTraceMap (const int i)
 Returns the global index of the boundary trace giving the index on the boundary expansion. More...
 
const Array< OneD, const int > & GetBndCondTraceToGlobalTraceMap ()
 
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 NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
 
void LocalBndToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
void LocalBndToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, int offset) const
 
void LocalBndToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) 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) const
 
const Array< OneD, const int > & GetExtraDirEdges ()
 
boost::shared_ptr< AssemblyMap > LinearSpaceMap (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 PatchMapSharedPtr & GetPatchMapFromPrevLevel (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...
 
PreconditionerType GetPreconType () const
 
NekDouble GetIterativeTolerance () const
 
int GetMaxIterations () const
 
int GetSuccessiveRHS () const
 
int GetLowestStaticCondLevel () const
 

Protected Member Functions

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
 
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) const
 
virtual void v_LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) 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 AssemblyMapSharedPtr v_LinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType)
 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)
 

Protected Attributes

Array< OneD, int > m_localToGlobalMap
 Integer map of local coeffs to global space. More...
 
Array< OneD, NekDouble > m_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::map< int, std::vector
< ExtraDirDof > > 		m_extraDirDofs
 Map indicating degrees of freedom which are Dirichlet but whose value is stored on another processor. More...
 
- Protected Attributes inherited from Nektar::MultiRegions::AssemblyMap
LibUtilities::SessionReaderSharedPtr m_session
 Session object. More...
 
LibUtilities::CommSharedPtr m_comm
 Communicator. 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 boundary coeffs to global space. More...
 
Array< OneD, NekDouble > m_localToGlobalBndSign
 Integer sign of local boundary coeffs to global space. More...
 
Array< OneD, int > m_bndCondCoeffsToGlobalCoeffsMap
 Integer map of bnd cond coeffs to global coefficients. More...
 
Array< OneD, NekDouble > m_bndCondCoeffsToGlobalCoeffsSign
 Integer map of bnd cond coeffs to global coefficients. More...
 
Array< OneD, int > m_bndCondTraceToGlobalTraceMap
 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...
 
PreconditionerType 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...
 
int m_successiveRHS
 sucessive RHS for iterative solver More...
 
Gs::gs_data * m_gsh
 
Gs::gs_data * m_bndGsh
 
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...
 

Private Types

typedef Array< OneD, const
ExpListSharedPtr > 		BndCondExp
 
typedef Array< OneD, const
SpatialDomains::BoundaryConditionShPtr > 		BndCond
 

Detailed Description

Constructs mappings for the C0 scalar continuous Galerkin formulation.

Mappings are created for three possible global solution types:

These mappings are used by GlobalLinSys to generate the global system.

Definition at line 71 of file AssemblyMapCG.h.

Member Typedef Documentation

typedef Array<OneD, const SpatialDomains::BoundaryConditionShPtr> Nektar::MultiRegions::AssemblyMapCG::BndCond
private

Definition at line 75 of file AssemblyMapCG.h.

typedef Array<OneD, const ExpListSharedPtr> Nektar::MultiRegions::AssemblyMapCG::BndCondExp
private

Definition at line 73 of file AssemblyMapCG.h.

Constructor & Destructor Documentation

Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG ( const LibUtilities::SessionReaderSharedPtr &  pSession,
const std::string  variable = "DefaultVar" 
)

Default constructor.

Definition at line 71 of file AssemblyMapCG.cpp.

References m_maxStaticCondLevel.

73  :
74  AssemblyMap(pSession,variable)
75  {
76  pSession->LoadParameter(
77  "MaxStaticCondLevel", m_maxStaticCondLevel, 100);
78  }
Nektar::MultiRegions::AssemblyMapCG::m_maxStaticCondLevel
int m_maxStaticCondLevel
Maximum static condensation level.
Definition: AssemblyMapCG.h:145
Nektar::MultiRegions::AssemblyMap::AssemblyMap
AssemblyMap()
Default constructor.
Definition: AssemblyMap.cpp:79
Nektar::MultiRegions::AssemblyMapCG::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.

STEP 6: 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.

Definition at line 1288 of file AssemblyMapCG.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::MultiRegions::AssemblyMap::CalculateBndSystemBandWidth(), CalculateFullSystemBandWidth(), CreateGraph(), Nektar::MultiRegions::DeterminePeriodicEdgeOrientId(), Nektar::MultiRegions::DeterminePeriodicFaceOrient(), Nektar::MultiRegions::eDirectMultiLevelStaticCond, Nektar::SpatialDomains::eDirichlet, Nektar::MultiRegions::eIterativeMultiLevelStaticCond, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eModified_B, Nektar::MultiRegions::ePETScMultiLevelStaticCond, Nektar::MultiRegions::eXxtMultiLevelStaticCond, Gs::Finalise(), Gs::Gather(), Nektar::MultiRegions::ExpList::GetCoeff_Offset(), Nektar::MultiRegions::ExpList::GetExp(), Nektar::MultiRegions::ExpList::GetOffset_Elmt_Id(), Gs::gs_min, Gs::Init(), Nektar::iterator, Nektar::MultiRegions::AssemblyMap::m_bndCondCoeffsToGlobalCoeffsMap, Nektar::MultiRegions::AssemblyMap::m_bndCondCoeffsToGlobalCoeffsSign, Nektar::MultiRegions::AssemblyMap::m_comm, m_extraDirDofs, Nektar::MultiRegions::AssemblyMap::m_hash, Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndMap, Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndSign, m_localToGlobalMap, m_localToGlobalSign, 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, m_numLocalBndCondCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalIntCoeffsPerPatch, Nektar::MultiRegions::AssemblyMap::m_numPatches, Nektar::MultiRegions::AssemblyMap::m_signChange, Nektar::MultiRegions::AssemblyMap::m_solnType, Nektar::MultiRegions::AssemblyMap::m_staticCondLevel, Nektar::LibUtilities::ReduceSum, SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd(), and Vmath::Vmax().

1299  : AssemblyMap(pSession, variable)
1300  {
1301  int i, j, k, l;
1302  int cnt = 0;
1303  int intDofCnt;
1304  int meshVertId, meshEdgeId, meshFaceId;
1305  int globalId;
1306  int nEdgeInteriorCoeffs;
1307  int nFaceInteriorCoeffs;
1308  int firstNonDirGraphVertId;
1309  LibUtilities::CommSharedPtr vComm = m_comm->GetRowComm();
1310  LocalRegions::ExpansionSharedPtr exp, bndExp;
1311  LibUtilities::BasisType bType;
1312  StdRegions::Orientation edgeOrient;
1313  StdRegions::Orientation faceOrient;
1314  Array<OneD, unsigned int> edgeInteriorMap;
1315  Array<OneD, int> edgeInteriorSign;
1316  Array<OneD, unsigned int> faceInteriorMap;
1317  Array<OneD, int> faceInteriorSign;
1318  PeriodicMap::const_iterator pIt;
1319 
1320  const LocalRegions::ExpansionVector &locExpVector = *(locExp.GetExp());
1321 
1322  m_signChange = false;
1323 
1324  // Stores vertex, edge and face reordered vertices.
1325  DofGraph graph(3);
1326  DofGraph dofs(3);
1327 
1328  set<int> extraDirVerts, extraDirEdges;
1329  BottomUpSubStructuredGraphSharedPtr bottomUpGraph;
1330 
1331  // Construct list of number of degrees of freedom for each vertex,
1332  // edge and face.
1333  for (i = 0; i < locExpVector.size(); ++i)
1334  {
1335  exp = locExpVector[i];
1336 
1337  for(j = 0; j < locExpVector[i]->GetNverts(); ++j)
1338  {
1339  dofs[0][exp->GetGeom()->GetVid(j)] = 1;
1340  }
1341 
1342  for(j = 0; j < locExpVector[i]->GetNedges(); ++j)
1343  {
1344  if (dofs[1].count(exp->GetGeom()->GetEid(j)) > 0)
1345  {
1346  if (dofs[1][exp->GetGeom()->GetEid(j)] !=
1347  locExpVector[i]->GetEdgeNcoeffs(j)-2)
1348  {
1349  ASSERTL0( (exp->GetEdgeBasisType(j) == LibUtilities::eModified_A) ||
1350  (exp->GetEdgeBasisType(j) == LibUtilities::eModified_B),
1351  "CG with variable order only available with modal expansion");
1352  }
1353  dofs[1][exp->GetGeom()->GetEid(j)] =
1354  min(dofs[1][exp->GetGeom()->GetEid(j)],
1355  locExpVector[i]->GetEdgeNcoeffs(j)-2);
1356  }
1357  else
1358  {
1359  dofs[1][exp->GetGeom()->GetEid(j)] =
1360  exp->GetEdgeNcoeffs(j) - 2;
1361  }
1362  }
1363 
1364  for(j = 0; j < locExpVector[i]->GetNfaces(); ++j)
1365  {
1366  if (dofs[2].count(exp->GetGeom()->GetFid(j)) > 0)
1367  {
1368  if (dofs[2][exp->GetGeom()->GetFid(j)] !=
1369  exp->GetFaceIntNcoeffs(j))
1370  {
1371  ASSERTL0( false,
1372  "CG with variable order not available in 3D");
1373  }
1374  dofs[2][exp->GetGeom()->GetFid(j)] =
1375  min(dofs[2][exp->GetGeom()->GetFid(j)],
1376  exp->GetFaceIntNcoeffs(j));
1377  }
1378  else
1379  {
1380  dofs[2][exp->GetGeom()->GetFid(j)] =
1381  exp->GetFaceIntNcoeffs(j);
1382  }
1383  }
1384  }
1385  // Now use information from all partitions to determine
1386  // the correct size
1387  map<int, int>::iterator dofIt;
1388  // edges
1389  Array<OneD, long> edgeId (dofs[1].size());
1390  Array<OneD, NekDouble> edgeDof (dofs[1].size());
1391  for(dofIt = dofs[1].begin(), i=0; dofIt != dofs[1].end(); dofIt++, i++)
1392  {
1393  edgeId[i] = dofIt->first;
1394  edgeDof[i] = (NekDouble) dofIt->second;
1395  }
1396  Gs::gs_data *tmp = Gs::Init(edgeId, vComm);
1397  Gs::Gather(edgeDof, Gs::gs_min, tmp);
1398  Gs::Finalise(tmp);
1399  for (i=0; i < dofs[1].size(); i++)
1400  {
1401  dofs[1][edgeId[i]] = (int) (edgeDof[i]+0.5);
1402  }
1403  // faces
1404  Array<OneD, long> faceId (dofs[2].size());
1405  Array<OneD, NekDouble> faceDof (dofs[2].size());
1406  for(dofIt = dofs[2].begin(), i=0; dofIt != dofs[2].end(); dofIt++, i++)
1407  {
1408  faceId[i] = dofIt->first;
1409  faceDof[i] = (NekDouble) dofIt->second;
1410  }
1411  Gs::gs_data *tmp2 = Gs::Init(faceId, vComm);
1412  Gs::Gather(faceDof, Gs::gs_min, tmp2);
1413  Gs::Finalise(tmp2);
1414  for (i=0; i < dofs[2].size(); i++)
1415  {
1416  dofs[2][faceId[i]] = (int) (faceDof[i]+0.5);
1417  }
1418 
1419  Array<OneD, const BndCond> bndCondVec(1, bndConditions);
1420 
1421  // Note that nExtraDirichlet is not used in the logic below; it just
1422  // needs to be set so that the coupled solver in
1423  // IncNavierStokesSolver can work.
1424  int nExtraDirichlet;
1425  int nGraphVerts =
1426  CreateGraph(locExp, bndCondExp, bndCondVec,
1427  checkIfSystemSingular, periodicVerts, periodicEdges,
1428  periodicFaces, graph, bottomUpGraph, extraDirVerts,
1429  extraDirEdges, firstNonDirGraphVertId,
1430  nExtraDirichlet);
1431 
1432  /*
1433  * Set up an array which contains the offset information of the
1434  * different graph vertices.
1435  *
1436  * This basically means to identify to how many global degrees of
1437  * freedom the individual graph vertices correspond. Obviously,
1438  * the graph vertices corresponding to the mesh-vertices account
1439  * for a single global DOF. However, the graph vertices
1440  * corresponding to the element edges correspond to N-2 global DOF
1441  * where N is equal to the number of boundary modes on this edge.
1442  */
1443  Array<OneD, int> graphVertOffset(
1444  graph[0].size() + graph[1].size() + graph[2].size() + 1);
1445 
1446  graphVertOffset[0] = 0;
1447 
1448  for(i = 0; i < locExpVector.size(); ++i)
1449  {
1450  exp = locExpVector[locExp.GetOffset_Elmt_Id(i)];
1451 
1452  for(j = 0; j < exp->GetNverts(); ++j)
1453  {
1454  meshVertId = exp->GetGeom()->GetVid(j);
1455  graphVertOffset[graph[0][meshVertId]+1] = 1;
1456  }
1457 
1458  for(j = 0; j < exp->GetNedges(); ++j)
1459  {
1460  nEdgeInteriorCoeffs = exp->GetEdgeNcoeffs(j) - 2;
1461  meshEdgeId = exp->GetGeom()->GetEid(j);
1462  graphVertOffset[graph[1][meshEdgeId]+1]
1463  = dofs[1][meshEdgeId];
1464 
1465  bType = exp->GetEdgeBasisType(j);
1466 
1467  // need a sign vector for modal expansions if nEdgeCoeffs
1468  // >=3 (not 4 because of variable order case)
1469  if(nEdgeInteriorCoeffs &&
1470  (bType == LibUtilities::eModified_A ||
1471  bType == LibUtilities::eModified_B))
1472  {
1473  m_signChange = true;
1474  }
1475  }
1476 
1477  for(j = 0; j < exp->GetNfaces(); ++j)
1478  {
1479  nFaceInteriorCoeffs = exp->GetFaceIntNcoeffs(j);
1480  meshFaceId = exp->GetGeom()->GetFid(j);
1481  graphVertOffset[graph[2][meshFaceId]+1] = dofs[2][meshFaceId];
1482  }
1483  }
1484 
1485  for(i = 1; i < graphVertOffset.num_elements(); i++)
1486  {
1487  graphVertOffset[i] += graphVertOffset[i-1];
1488  }
1489 
1490  // Allocate the proper amount of space for the class-data
1491  m_numLocalCoeffs = numLocalCoeffs;
1492  m_numGlobalDirBndCoeffs = graphVertOffset[firstNonDirGraphVertId];
1493  m_localToGlobalMap = Array<OneD, int>(m_numLocalCoeffs,-1);
1494  m_localToGlobalBndMap = Array<OneD, int>(m_numLocalBndCoeffs,-1);
1495  m_bndCondCoeffsToGlobalCoeffsMap = Array<OneD, int>(m_numLocalBndCondCoeffs,-1);
1496  // If required, set up the sign-vector
1497  if(m_signChange)
1498  {
1499  m_localToGlobalSign = Array<OneD, NekDouble>(m_numLocalCoeffs,1.0);
1500  m_localToGlobalBndSign = Array<OneD, NekDouble>(m_numLocalBndCoeffs,1.0);
1501  m_bndCondCoeffsToGlobalCoeffsSign = Array<OneD,NekDouble>(m_numLocalBndCondCoeffs,1.0);
1502  }
1503 
1504  m_staticCondLevel = 0;
1505  m_numPatches = locExpVector.size();
1506  m_numLocalBndCoeffsPerPatch = Array<OneD, unsigned int>(m_numPatches);
1507  m_numLocalIntCoeffsPerPatch = Array<OneD, unsigned int>(m_numPatches);
1508  for(i = 0; i < m_numPatches; ++i)
1509  {
1510  m_numLocalBndCoeffsPerPatch[i] = (unsigned int)
1511  locExpVector[locExp.GetOffset_Elmt_Id(i)]->NumBndryCoeffs();
1512  m_numLocalIntCoeffsPerPatch[i] = (unsigned int)
1513  locExpVector[locExp.GetOffset_Elmt_Id(i)]->GetNcoeffs() -
1514  locExpVector[locExp.GetOffset_Elmt_Id(i)]->NumBndryCoeffs();
1515  }
1516 
1517  /**
1518  * STEP 6: Now, all ingredients are ready to set up the actual
1519  * local to global mapping.
1520  *
1521  * The remainder of the map consists of the element-interior
1522  * degrees of freedom. This leads to the block-diagonal submatrix
1523  * as each element-interior mode is globally orthogonal to modes
1524  * in all other elements.
1525  */
1526  cnt = 0;
1527 
1528  // Loop over all the elements in the domain
1529  for(i = 0; i < locExpVector.size(); ++i)
1530  {
1531  exp = locExpVector[i];
1532  cnt = locExp.GetCoeff_Offset(i);
1533  for(j = 0; j < exp->GetNverts(); ++j)
1534  {
1535  meshVertId = exp->GetGeom()->GetVid(j);
1536 
1537  // Set the global DOF for vertex j of element i
1538  m_localToGlobalMap[cnt+exp->GetVertexMap(j)] =
1539  graphVertOffset[graph[0][meshVertId]];
1540  }
1541 
1542  for(j = 0; j < exp->GetNedges(); ++j)
1543  {
1544  nEdgeInteriorCoeffs = exp->GetEdgeNcoeffs(j)-2;
1545  edgeOrient = exp->GetGeom()->GetEorient(j);
1546  meshEdgeId = exp->GetGeom()->GetEid(j);
1547 
1548  pIt = periodicEdges.find(meshEdgeId);
1549 
1550  // See if this edge is periodic. If it is, then we map all
1551  // edges to the one with lowest ID, and align all
1552  // coefficients to this edge orientation.
1553  if (pIt != periodicEdges.end())
1554  {
1555  pair<int, StdRegions::Orientation> idOrient =
1556  DeterminePeriodicEdgeOrientId(
1557  meshEdgeId, edgeOrient, pIt->second);
1558  edgeOrient = idOrient.second;
1559  }
1560 
1561  exp->GetEdgeInteriorMap(j,edgeOrient,edgeInteriorMap,edgeInteriorSign);
1562 
1563  // Set the global DOF's for the interior modes of edge j
1564  for(k = 0; k < dofs[1][exp->GetGeom()->GetEid(j)]; ++k)
1565  {
1566  m_localToGlobalMap[cnt+edgeInteriorMap[k]] =
1567  graphVertOffset[graph[1][meshEdgeId]]+k;
1568  }
1569  for(k = dofs[1][exp->GetGeom()->GetEid(j)]; k < nEdgeInteriorCoeffs; ++k)
1570  {
1571  m_localToGlobalMap[cnt+edgeInteriorMap[k]] =
1572  graphVertOffset[graph[1][meshEdgeId]];
1573  }
1574 
1575  // Fill the sign vector if required
1576  if(m_signChange)
1577  {
1578  for(k = 0; k < dofs[1][exp->GetGeom()->GetEid(j)]; ++k)
1579  {
1580  m_localToGlobalSign[cnt+edgeInteriorMap[k]] = (NekDouble) edgeInteriorSign[k];
1581  }
1582  for(k = dofs[1][exp->GetGeom()->GetEid(j)]; k < nEdgeInteriorCoeffs; ++k)
1583  {
1584  m_localToGlobalSign[cnt+edgeInteriorMap[k]] = 0.0;
1585  }
1586  }
1587  }
1588 
1589  for(j = 0; j < exp->GetNfaces(); ++j)
1590  {
1591  nFaceInteriorCoeffs = exp->GetFaceIntNcoeffs(j);
1592  faceOrient = exp->GetGeom()->GetForient(j);
1593  meshFaceId = exp->GetGeom()->GetFid(j);
1594 
1595  pIt = periodicFaces.find(meshFaceId);
1596 
1597  if (pIt != periodicFaces.end() &&
1598  meshFaceId == min(meshFaceId, pIt->second[0].id))
1599  {
1600  faceOrient = DeterminePeriodicFaceOrient(faceOrient,pIt->second[0].orient);
1601  }
1602 
1603  exp->GetFaceInteriorMap(j,faceOrient,faceInteriorMap,faceInteriorSign);
1604 
1605  // Set the global DOF's for the interior modes of face j
1606  for(k = 0; k < dofs[2][exp->GetGeom()->GetFid(j)]; ++k)
1607  {
1608  m_localToGlobalMap[cnt+faceInteriorMap[k]] =
1609  graphVertOffset[graph[2][meshFaceId]]+k;
1610  }
1611  for(k = dofs[2][exp->GetGeom()->GetFid(j)]; k < nFaceInteriorCoeffs; ++k)
1612  {
1613  m_localToGlobalMap[cnt+faceInteriorMap[k]] =
1614  graphVertOffset[graph[2][meshFaceId]];
1615  }
1616 
1617  if(m_signChange)
1618  {
1619  for(k = 0; k < dofs[2][exp->GetGeom()->GetFid(j)]; ++k)
1620  {
1621  m_localToGlobalSign[cnt+faceInteriorMap[k]] = (NekDouble) faceInteriorSign[k];
1622  }
1623  for(k = dofs[2][exp->GetGeom()->GetFid(j)]; k < nFaceInteriorCoeffs; ++k)
1624  {
1625  m_localToGlobalSign[cnt+faceInteriorMap[k]] = 0;
1626  }
1627  }
1628 
1629  }
1630  }
1631 
1632  // Set up the mapping for the boundary conditions
1633  cnt = 0;
1634  int offset = 0;
1635  for(i = 0; i < bndCondExp.num_elements(); i++)
1636  {
1637  set<int> foundExtraVerts, foundExtraEdges;
1638  for(j = 0; j < bndCondExp[i]->GetNumElmts(); j++)
1639  {
1640  bndExp = bndCondExp[i]->GetExp(j);
1641  cnt = offset + bndCondExp[i]->GetCoeff_Offset(j);
1642  for(k = 0; k < bndExp->GetNverts(); k++)
1643  {
1644  meshVertId = bndExp->GetGeom()->GetVid(k);
1645  m_bndCondCoeffsToGlobalCoeffsMap[cnt+bndExp->GetVertexMap(k)] = graphVertOffset[graph[0][meshVertId]];
1646 
1647  if (bndConditions[i]->GetBoundaryConditionType() !=
1648  SpatialDomains::eDirichlet)
1649  {
1650  continue;
1651  }
1652 
1653  set<int>::iterator iter = extraDirVerts.find(meshVertId);
1654  if (iter != extraDirVerts.end() &&
1655  foundExtraVerts.count(meshVertId) == 0)
1656  {
1657  int loc = bndCondExp[i]->GetCoeff_Offset(j) +
1658  bndExp->GetVertexMap(k);
1659  int gid = graphVertOffset[
1660  graph[0][meshVertId]];
1661  ExtraDirDof t(loc, gid, 1.0);
1662  m_extraDirDofs[i].push_back(t);
1663  foundExtraVerts.insert(meshVertId);
1664  }
1665  }
1666 
1667  for(k = 0; k < bndExp->GetNedges(); k++)
1668  {
1669  nEdgeInteriorCoeffs = bndExp->GetEdgeNcoeffs(k)-2;
1670  edgeOrient = bndExp->GetGeom()->GetEorient(k);
1671  meshEdgeId = bndExp->GetGeom()->GetEid(k);
1672 
1673  pIt = periodicEdges.find(meshEdgeId);
1674 
1675  // See if this edge is periodic. If it is, then we map
1676  // all edges to the one with lowest ID, and align all
1677  // coefficients to this edge orientation.
1678  if (pIt != periodicEdges.end())
1679  {
1680  pair<int, StdRegions::Orientation> idOrient =
1681  DeterminePeriodicEdgeOrientId(
1682  meshEdgeId, edgeOrient, pIt->second);
1683  edgeOrient = idOrient.second;
1684  }
1685 
1686  bndExp->GetEdgeInteriorMap(
1687  k,edgeOrient,edgeInteriorMap,edgeInteriorSign);
1688 
1689  for(l = 0; l < nEdgeInteriorCoeffs; ++l)
1690  {
1691  m_bndCondCoeffsToGlobalCoeffsMap[cnt+edgeInteriorMap[l]] =
1692  graphVertOffset[graph[1][meshEdgeId]]+l;
1693  }
1694 
1695  // Fill the sign vector if required
1696  if(m_signChange)
1697  {
1698  for(l = 0; l < nEdgeInteriorCoeffs; ++l)
1699  {
1700  m_bndCondCoeffsToGlobalCoeffsSign[cnt+edgeInteriorMap[l]] = (NekDouble) edgeInteriorSign[l];
1701  }
1702  }
1703 
1704  if (bndConditions[i]->GetBoundaryConditionType() !=
1705  SpatialDomains::eDirichlet)
1706  {
1707  continue;
1708  }
1709 
1710  set<int>::iterator iter = extraDirEdges.find(meshEdgeId);
1711  if (iter != extraDirEdges.end() &&
1712  foundExtraEdges.count(meshEdgeId) == 0 &&
1713  nEdgeInteriorCoeffs > 0)
1714  {
1715  for(l = 0; l < nEdgeInteriorCoeffs; ++l)
1716  {
1717  int loc = bndCondExp[i]->GetCoeff_Offset(j) +
1718  edgeInteriorMap[l];
1719  int gid = graphVertOffset[
1720  graph[1][meshEdgeId]]+l;
1721  ExtraDirDof t(loc, gid, edgeInteriorSign[l]);
1722  m_extraDirDofs[i].push_back(t);
1723  }
1724  foundExtraEdges.insert(meshEdgeId);
1725  }
1726  }
1727 
1728  meshFaceId = bndExp->GetGeom()->GetGlobalID();
1729  intDofCnt = 0;
1730  for(k = 0; k < bndExp->GetNcoeffs(); k++)
1731  {
1732  if(m_bndCondCoeffsToGlobalCoeffsMap[cnt+k] == -1)
1733  {
1734  m_bndCondCoeffsToGlobalCoeffsMap[cnt+k] =
1735  graphVertOffset[graph[bndExp->GetNumBases()][meshFaceId]]+intDofCnt;
1736  intDofCnt++;
1737  }
1738  }
1739  }
1740  offset += bndCondExp[i]->GetNcoeffs();
1741  }
1742 
1743  globalId = Vmath::Vmax(m_numLocalCoeffs,&m_localToGlobalMap[0],1)+1;
1744  m_numGlobalBndCoeffs = globalId;
1745 
1746 
1747  /*
1748  * The boundary condition mapping is generated from the same vertex
1749  * renumbering.
1750  */
1751  cnt=0;
1752  for(i = 0; i < m_numLocalCoeffs; ++i)
1753  {
1754  if(m_localToGlobalMap[i] == -1)
1755  {
1756  m_localToGlobalMap[i] = globalId++;
1757  }
1758  else
1759  {
1760  if(m_signChange)
1761  {
1762  m_localToGlobalBndSign[cnt]=m_localToGlobalSign[i];
1763  }
1764  m_localToGlobalBndMap[cnt++]=m_localToGlobalMap[i];
1765  }
1766  }
1767 
1768  m_numGlobalCoeffs = globalId;
1769 
1770  SetUpUniversalC0ContMap(locExp, periodicVerts, periodicEdges, periodicFaces);
1771 
1772  // Since we can now have multiple entries to m_extraDirDofs due to
1773  // periodic boundary conditions we make a call to work out the
1774  // multiplicity of all entries and invert (Need to be after
1775  // SetupUniversalC0ContMap)
1776  Array<OneD, NekDouble> valence(m_numGlobalBndCoeffs,0.0);
1777 
1778  // Fill in Dirichlet coefficients that are to be sent to other
1779  // processors with a value of 1
1780  map<int, vector<ExtraDirDof> >::iterator Tit;
1781 
1782  // Generate valence for extraDirDofs
1783  for (Tit = m_extraDirDofs.begin(); Tit != m_extraDirDofs.end(); ++Tit)
1784  {
1785  for (i = 0; i < Tit->second.size(); ++i)
1786  {
1787  valence[Tit->second[i].get<1>()] = 1.0;
1788  }
1789  }
1790 
1791  UniversalAssembleBnd(valence);
1792 
1793  // Set third argument of tuple to inverse of valence.
1794  for (Tit = m_extraDirDofs.begin(); Tit != m_extraDirDofs.end(); ++Tit)
1795  {
1796  for (i = 0; i < Tit->second.size(); ++i)
1797  {
1798  boost::get<2>(Tit->second.at(i)) /= valence[Tit->second.at(i).get<1>()];
1799  }
1800  }
1801 
1802  // Set up the local to global map for the next level when using
1803  // multi-level static condensation
1804  if ((m_solnType == eDirectMultiLevelStaticCond ||
1805  m_solnType == eIterativeMultiLevelStaticCond ||
1806  m_solnType == eXxtMultiLevelStaticCond ||
1807  m_solnType == ePETScMultiLevelStaticCond) && nGraphVerts)
1808  {
1809  if (m_staticCondLevel < (bottomUpGraph->GetNlevels()-1))
1810  {
1811  Array<OneD, int> vwgts_perm(
1812  dofs[0].size() + dofs[1].size() + dofs[2].size()
1813  - firstNonDirGraphVertId);
1814 
1815  for (i = 0; i < locExpVector.size(); ++i)
1816  {
1817  exp = locExpVector[locExp.GetOffset_Elmt_Id(i)];
1818 
1819  for (j = 0; j < exp->GetNverts(); ++j)
1820  {
1821  meshVertId = exp->GetGeom()->GetVid(j);
1822 
1823  if (graph[0][meshVertId] >= firstNonDirGraphVertId)
1824  {
1825  vwgts_perm[graph[0][meshVertId] -
1826  firstNonDirGraphVertId] =
1827  dofs[0][meshVertId];
1828  }
1829  }
1830 
1831  for (j = 0; j < exp->GetNedges(); ++j)
1832  {
1833  meshEdgeId = exp->GetGeom()->GetEid(j);
1834 
1835  if (graph[1][meshEdgeId] >= firstNonDirGraphVertId)
1836  {
1837  vwgts_perm[graph[1][meshEdgeId] -
1838  firstNonDirGraphVertId] =
1839  dofs[1][meshEdgeId];
1840  }
1841  }
1842 
1843  for (j = 0; j < exp->GetNfaces(); ++j)
1844  {
1845  meshFaceId = exp->GetGeom()->GetFid(j);
1846 
1847  if (graph[2][meshFaceId] >= firstNonDirGraphVertId)
1848  {
1849  vwgts_perm[graph[2][meshFaceId] -
1850  firstNonDirGraphVertId] =
1851  dofs[2][meshFaceId];
1852  }
1853  }
1854  }
1855 
1856  bottomUpGraph->ExpandGraphWithVertexWeights(vwgts_perm);
1857  m_nextLevelLocalToGlobalMap = MemoryManager<AssemblyMap>::
1858  AllocateSharedPtr(this, bottomUpGraph);
1859  }
1860  }
1861 
1862  m_hash = boost::hash_range(m_localToGlobalMap.begin(),
1863  m_localToGlobalMap.end());
1864 
1865  // Add up hash values if parallel
1866  int hash = m_hash;
1867  vComm->AllReduce(hash, LibUtilities::ReduceSum);
1868  m_hash = hash;
1869 
1870  CalculateBndSystemBandWidth();
1871  CalculateFullSystemBandWidth();
1872  }
ASSERTL0
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
Nektar::MultiRegions::AssemblyMap::m_signChange
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:345
Nektar::MultiRegions::AssemblyMap::m_numGlobalBndCoeffs
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:314
Nektar::MultiRegions::AssemblyMap::m_comm
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:306
Gs::Gather
static void Gather(Nektar::Array< OneD, NekDouble > pU, gs_op pOp, gs_data *pGsh, Nektar::Array< OneD, NekDouble > pBuffer=NullNekDouble1DArray)
Performs a gather-scatter operation of the provided values.
Definition: GsLib.hpp:224
Nektar::MemoryManager::AllocateSharedPtr
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
Definition: NekMemoryManager.hpp:235
Gs::Finalise
static void Finalise(gs_data *pGsh)
Deallocates the GSLib mapping data.
Definition: GsLib.hpp:206
Vmath::Vmax
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:765
Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap
void SetUpUniversalC0ContMap(const ExpList &locExp, const PeriodicMap &perVerts=NullPeriodicMap, const PeriodicMap &perEdges=NullPeriodicMap, const PeriodicMap &perFaces=NullPeriodicMap)
Definition: AssemblyMapCG.cpp:1993
Nektar::MultiRegions::AssemblyMapCG::CreateGraph
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)
Definition: AssemblyMapCG.cpp:80
Gs::Init
static gs_data * Init(const Nektar::Array< OneD, long > pId, const LibUtilities::CommSharedPtr &pComm)
Initialise Gather-Scatter map.
Definition: GsLib.hpp:151
Nektar::LibUtilities::eModified_A
Principle Modified Functions .
Definition: BasisType.h:49
Nektar::MultiRegions::BottomUpSubStructuredGraphSharedPtr
boost::shared_ptr< BottomUpSubStructuredGraph > BottomUpSubStructuredGraphSharedPtr
Definition: SubStructuredGraph.h:60
Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:331
Nektar::MultiRegions::AssemblyMapCG::m_numLocalBndCondCoeffs
int m_numLocalBndCondCoeffs
Number of local boundary condition coefficients.
Definition: AssemblyMapCG.h:139
Nektar::MultiRegions::AssemblyMapCG::m_extraDirDofs
std::map< int, std::vector< ExtraDirDof > > m_extraDirDofs
Map indicating degrees of freedom which are Dirichlet but whose value is stored on another processor...
Definition: AssemblyMapCG.h:148
Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth
void CalculateFullSystemBandWidth()
Calculate the bandwith of the full matrix system.
Definition: AssemblyMapCG.cpp:2338
Nektar::LocalRegions::ExpansionVector
std::vector< ExpansionSharedPtr > ExpansionVector
Definition: Expansion.h:70
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
Nektar::MultiRegions::AssemblyMap::m_nextLevelLocalToGlobalMap
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:395
Nektar::MultiRegions::AssemblyMap::m_hash
size_t m_hash
Hash for map.
Definition: AssemblyMap.h:309
Nektar::LibUtilities::CommSharedPtr
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffsPerPatch
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:388
Nektar::MultiRegions::AssemblyMap::m_solnType
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:363
Nektar::MultiRegions::AssemblyMap::m_bndCondCoeffsToGlobalCoeffsSign
Array< OneD, NekDouble > m_bndCondCoeffsToGlobalCoeffsSign
Integer map of bnd cond coeffs to global coefficients.
Definition: AssemblyMap.h:354
Nektar::MultiRegions::ExtraDirDof
boost::tuple< int, int, NekDouble > ExtraDirDof
Definition: AssemblyMapCG.h:54
Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
Nektar::LibUtilities::eModified_B
Principle Modified Functions .
Definition: BasisType.h:50
Nektar::MultiRegions::DeterminePeriodicFaceOrient
StdRegions::Orientation DeterminePeriodicFaceOrient(StdRegions::Orientation faceOrient, StdRegions::Orientation perFaceOrient)
Determine relative orientation between two faces.
Definition: AssemblyMapCG.cpp:1950
Nektar::NekDouble
double NekDouble
Definition: NektarUnivTypeDefs.hpp:44
Nektar::LocalRegions::ExpansionSharedPtr
boost::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
Nektar::MultiRegions::AssemblyMap::m_numLocalIntCoeffsPerPatch
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
The number of int dofs per patch.
Definition: AssemblyMap.h:390
Nektar::MultiRegions::AssemblyMap::CalculateBndSystemBandWidth
void CalculateBndSystemBandWidth()
Calculates the bandwidth of the boundary system.
Definition: AssemblyMap.cpp:456
Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndMap
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:348
Nektar::MultiRegions::AssemblyMap::m_bndCondCoeffsToGlobalCoeffsMap
Array< OneD, int > m_bndCondCoeffsToGlobalCoeffsMap
Integer map of bnd cond coeffs to global coefficients.
Definition: AssemblyMap.h:352
Nektar::iterator
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
Definition: MatrixOperations.cpp:96
Nektar::MultiRegions::DofGraph
std::vector< std::map< int, int > > DofGraph
Definition: AssemblyMapCG.h:56
Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffs
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:312
Nektar::MultiRegions::AssemblyMap::m_staticCondLevel
int m_staticCondLevel
The level of recursion in the case of multi-level static condensation.
Definition: AssemblyMap.h:384
Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndSign
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign
Array< OneD, NekDouble > m_localToGlobalSign
Integer sign of local coeffs to global space.
Definition: AssemblyMapCG.h:117
Nektar::MultiRegions::AssemblyMap::AssemblyMap
AssemblyMap()
Default constructor.
Definition: AssemblyMap.cpp:79
Gs::gs_data
Definition: GsLib.hpp:117
Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
Definition: AssemblyMap.cpp:1128
Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:342
Nektar::MultiRegions::DeterminePeriodicEdgeOrientId
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...
Definition: AssemblyMapCG.cpp:1902
Nektar::MultiRegions::AssemblyMap::m_numPatches
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:386
Nektar::MultiRegions::AssemblyMapCG::~AssemblyMapCG ( )
virtual

Destructor.

Definition at line 1877 of file AssemblyMapCG.cpp.

References Gs::Finalise(), Nektar::MultiRegions::AssemblyMap::m_bndGsh, and Nektar::MultiRegions::AssemblyMap::m_gsh.

1878  {
1879  Gs::Finalise(m_gsh);
1880  Gs::Finalise(m_bndGsh);
1881  }
Gs::Finalise
static void Finalise(gs_data *pGsh)
Deallocates the GSLib mapping data.
Definition: GsLib.hpp:206

Member Function Documentation

void Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth ( )
protected

Calculate the bandwith of the full matrix system.

The bandwidth calculated here corresponds to what is referred to as half-bandwidth. If the elements of the matrix are designated as a_ij, it corresponds to the maximum value of |i-j| for non-zero a_ij. As a result, the value also corresponds to the number of sub- or super-diagonals.

The bandwith can be calculated elementally as it corresponds to the maximal elemental bandwith (i.e. the maximal difference in global DOF index for every element).

We caluclate here the bandwith of the full global system.

Definition at line 2338 of file AssemblyMapCG.cpp.

References m_fullSystemBandWidth, m_localToGlobalMap, Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffsPerPatch, Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalIntCoeffsPerPatch, and Nektar::MultiRegions::AssemblyMap::m_numPatches.

Referenced by AssemblyMapCG().

2339  {
2340  int i,j;
2341  int cnt = 0;
2342  int locSize;
2343  int maxId;
2344  int minId;
2345  int bwidth = -1;
2346  for(i = 0; i < m_numPatches; ++i)
2347  {
2348  locSize = m_numLocalBndCoeffsPerPatch[i]+m_numLocalIntCoeffsPerPatch[i];
2349  maxId = -1;
2350  minId = m_numLocalCoeffs+1;
2351  for(j = 0; j < locSize; j++)
2352  {
2353  if(m_localToGlobalMap[cnt+j] >= m_numGlobalDirBndCoeffs)
2354  {
2355  if(m_localToGlobalMap[cnt+j] > maxId)
2356  {
2357  maxId = m_localToGlobalMap[cnt+j];
2358  }
2359 
2360  if(m_localToGlobalMap[cnt+j] < minId)
2361  {
2362  minId = m_localToGlobalMap[cnt+j];
2363  }
2364  }
2365  }
2366  bwidth = (bwidth>(maxId-minId))?bwidth:(maxId-minId);
2367 
2368  cnt+=locSize;
2369  }
2370 
2371  m_fullSystemBandWidth = bwidth;
2372  }
Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:331
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffsPerPatch
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:388
Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
Nektar::MultiRegions::AssemblyMap::m_numLocalIntCoeffsPerPatch
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
The number of int dofs per patch.
Definition: AssemblyMap.h:390
Nektar::MultiRegions::AssemblyMapCG::m_fullSystemBandWidth
int m_fullSystemBandWidth
Bandwith of the full matrix system (no static condensation).
Definition: AssemblyMapCG.h:119
Nektar::MultiRegions::AssemblyMap::m_numPatches
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:386
int Nektar::MultiRegions::AssemblyMapCG::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 
)
protected
  • Count verts, edges, face and add up edges and face sizes
  • Periodic vertices
  • Periodic edges
  • Periodic faces

STEP 4: Fill the #graph[0] and #graph[1] with the optimal ordering from boost.

Definition at line 80 of file AssemblyMapCG.cpp.

References ASSERTL0, ASSERTL1, Nektar::MultiRegions::CuthillMckeeReordering(), Nektar::MultiRegions::eDirectFullMatrix, Nektar::MultiRegions::eDirectMultiLevelStaticCond, Nektar::MultiRegions::eDirectStaticCond, Nektar::SpatialDomains::eDirichlet, Nektar::MultiRegions::eIterativeFull, Nektar::MultiRegions::eIterativeMultiLevelStaticCond, Nektar::MultiRegions::eIterativeStaticCond, Nektar::SpatialDomains::eNeumann, Nektar::MultiRegions::ePETScFullMatrix, Nektar::MultiRegions::ePETScMultiLevelStaticCond, Nektar::MultiRegions::ePETScStaticCond, Nektar::MultiRegions::eXxtFullMatrix, Nektar::MultiRegions::eXxtMultiLevelStaticCond, Nektar::MultiRegions::eXxtStaticCond, Gs::Finalise(), Gs::Gather(), Nektar::MultiRegions::ExpList::GetExp(), Nektar::LocalRegions::Expansion::GetGeom(), Nektar::MultiRegions::ExpList::GetOffset_Elmt_Id(), Nektar::MultiRegions::GlobalSysSolnTypeMap, Gs::gs_add, Vmath::Imax(), Gs::Init(), Nektar::iterator, Nektar::MultiRegions::AssemblyMap::m_comm, m_extraDirEdges, Nektar::MultiRegions::AssemblyMap::m_lowestStaticCondLevel, m_numDirEdges, m_numDirFaces, Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffs, m_numLocalBndCondCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalDirBndCoeffs, m_numNonDirEdgeModes, m_numNonDirEdges, m_numNonDirFaceModes, m_numNonDirFaces, m_numNonDirVertexModes, Nektar::MultiRegions::AssemblyMap::m_session, Nektar::MultiRegions::AssemblyMap::m_solnType, Nektar::MultiRegions::AssemblyMap::m_systemSingular, Nektar::MultiRegions::MultiLevelBisectionReordering(), Nektar::MultiRegions::NoReordering(), Nektar::LibUtilities::ReduceMax, Nektar::LibUtilities::ReduceMin, Nektar::LibUtilities::ReduceSum, and Vmath::Vsum().

Referenced by AssemblyMapCG(), and Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap().

95  {
96  int graphVertId = 0;
97  int vMaxVertId = -1;
98  int i, j, k, l, cnt;
99  int meshVertId, meshEdgeId, meshFaceId;
100  int meshVertId2, meshEdgeId2;
101 
102  LocalRegions::ExpansionSharedPtr exp, bndExp;
103  const LocalRegions::ExpansionVector &locExpVector =
104  *(locExp.GetExp());
105  LibUtilities::CommSharedPtr vComm = m_comm->GetRowComm();
106  PeriodicMap::const_iterator pIt;
107 
108  m_numLocalBndCondCoeffs = 0;
109  m_systemSingular = checkIfSystemSingular;
110 
111  for(i = 0; i < bndCondExp.num_elements(); i++)
112  {
113  // Check to see if any value on boundary has Dirichlet value.
114  cnt = 0;
115  for(k = 0; k < bndConditions.num_elements(); ++k)
116  {
117  if (bndConditions[k][i]->GetBoundaryConditionType() ==
118  SpatialDomains::eDirichlet)
119  {
120  cnt++;
121  }
122  if (bndConditions[k][i]->GetBoundaryConditionType() !=
123  SpatialDomains::eNeumann)
124  {
125  m_systemSingular = false;
126  }
127  }
128 
129  // Find the maximum boundary vertex ID on this process. This is
130  // used later to pin a vertex if the system is singular.
131  for (j = 0; j < bndCondExp[i]->GetNumElmts(); ++j)
132  {
133  bndExp = bndCondExp[i]->GetExp(j)->as<LocalRegions::Expansion>();
134  for (k = 0; k < bndExp->GetNverts(); ++k)
135  {
136  if (vMaxVertId < bndExp->GetGeom()->GetVid(k))
137  {
138  vMaxVertId = bndExp->GetGeom()->GetVid(k);
139  }
140  }
141 
142  }
143 
144  // If all boundaries are Dirichlet take out of mask
145  if(cnt == bndConditions.num_elements())
146  {
147  for(j = 0; j < bndCondExp[i]->GetNumElmts(); j++)
148  {
149  bndExp = bndCondExp[i]->GetExp(j);
150 
151  for (k = 0; k < bndExp->GetNverts(); k++)
152  {
153  meshVertId = bndExp->GetGeom()->GetVid(k);
154  if (graph[0].count(meshVertId) == 0)
155  {
156  graph[0][meshVertId] = graphVertId++;
157  }
158  }
159 
160  for (k = 0; k < bndExp->GetNedges(); k++)
161  {
162  meshEdgeId = bndExp->GetGeom()->GetEid(k);
163  if (graph[1].count(meshEdgeId) == 0)
164  {
165  graph[1][meshEdgeId] = graphVertId++;
166  }
167  }
168 
169  // Possibly not a face.
170  meshFaceId = bndExp->GetGeom()->GetGlobalID();
171  const int bndDim = bndExp->GetNumBases();
172  if (graph[bndDim].count(meshFaceId) == 0)
173  {
174  graph[bndDim][meshFaceId] = graphVertId++;
175  }
176  m_numLocalDirBndCoeffs += bndExp->GetNcoeffs();
177  }
178  }
179 
180  m_numLocalBndCondCoeffs += bndCondExp[i]->GetNcoeffs();
181  }
182 
183  // Number of dirichlet edges and faces (not considering periodic
184  // BCs)
185  m_numDirEdges = graph[1].size();
186  m_numDirFaces = graph[2].size();
187 
188  /*
189  * The purpose of this routine is to deal with those degrees of
190  * freedom that are Dirichlet, but do not have a local Dirichlet
191  * boundary condition expansion set.
192  *
193  * For example, in 2D, consider a triangulation of a square into two
194  * triangles. Now imagine one edge of the square is Dirichlet and
195  * the problem is run on two processors. On one processor, one
196  * triangle vertex is Dirichlet, but doesn't know this since the
197  * Dirichlet composite lives on the other processor.
198  *
199  * When the global linear system is solved therefore, there is an
200  * inconsistency that at best leads to an inaccurate answer or a
201  * divergence of the system.
202  *
203  * This routine identifies such cases for 2D, and also for 3D where
204  * e.g. edges may have the same problem (consider an extrusion of
205  * the case above, for example).
206  */
207 
208  // Collate information on Dirichlet vertices from all processes
209  int n = vComm->GetSize();
210  int p = vComm->GetRank();
211 
212  // At this point, graph only contains information from Dirichlet
213  // boundaries. Therefore make a global list of the vert and edge
214  // information on all processors.
215  Array<OneD, int> vertcounts (n, 0);
216  Array<OneD, int> vertoffsets(n, 0);
217  Array<OneD, int> edgecounts (n, 0);
218  Array<OneD, int> edgeoffsets(n, 0);
219  vertcounts[p] = graph[0].size();
220  edgecounts[p] = graph[1].size();
221  vComm->AllReduce(vertcounts, LibUtilities::ReduceSum);
222  vComm->AllReduce(edgecounts, LibUtilities::ReduceSum);
223 
224  for (i = 1; i < n; ++i)
225  {
226  vertoffsets[i] = vertoffsets[i-1] + vertcounts[i-1];
227  edgeoffsets[i] = edgeoffsets[i-1] + edgecounts[i-1];
228  }
229 
230  int nTotVerts = Vmath::Vsum(n,vertcounts,1);
231  int nTotEdges = Vmath::Vsum(n,edgecounts,1);
232 
233  Array<OneD, int> vertlist(nTotVerts, 0);
234  Array<OneD, int> edgelist(nTotEdges, 0);
235  std::map<int, int>::iterator it;
236 
237  // construct list of global ids of global vertices
238  for (it = graph[0].begin(), i = 0;
239  it != graph[0].end();
240  ++it, ++i)
241  {
242  vertlist[vertoffsets[p] + i] = it->first;
243  }
244 
245  // construct list of global ids of global edges
246  for (it = graph[1].begin(), i = 0;
247  it != graph[1].end();
248  ++it, ++i)
249  {
250  edgelist[edgeoffsets[p] + i] = it->first;
251  }
252  vComm->AllReduce(vertlist, LibUtilities::ReduceSum);
253  vComm->AllReduce(edgelist, LibUtilities::ReduceSum);
254 
255  // Now we have a list of all Dirichlet vertices and edges on all
256  // processors.
257  nExtraDirichlet = 0;
258  map<int, int> extraDirVertIds, extraDirEdgeIds;
259 
260  // Ensure Dirchlet vertices are consistently recorded between
261  // processes (e.g. Dirichlet region meets Neumann region across a
262  // partition boundary requires vertex on partition to be Dirichlet).
263  //
264  // To do this we look over all elements and vertices in local
265  // partition and see if they match the values stored in the vertlist
266  // from othe processors and if so record the meshVertId/meshEdgeId
267  // and the processor it comes from.
268  for (i = 0; i < n; ++i)
269  {
270  if (i == p)
271  {
272  continue;
273  }
274 
275  for(j = 0; j < locExpVector.size(); j++)
276  {
277  exp = locExpVector[locExp.GetOffset_Elmt_Id(j)];
278 
279  for(k = 0; k < exp->GetNverts(); k++)
280  {
281  meshVertId = exp->GetGeom()->GetVid(k);
282  if(graph[0].count(meshVertId) == 0)
283  {
284  for (l = 0; l < vertcounts[i]; ++l)
285  {
286  if (vertlist[vertoffsets[i]+l] == meshVertId)
287  {
288  extraDirVertIds[meshVertId] = i;
289  graph[0][meshVertId] = graphVertId++;
290  nExtraDirichlet++;
291  }
292  }
293  }
294  }
295 
296  for(k = 0; k < exp->GetNedges(); k++)
297  {
298  meshEdgeId = exp->GetGeom()->GetEid(k);
299  if(graph[1].count(meshEdgeId) == 0)
300  {
301  for (l = 0; l < edgecounts[i]; ++l)
302  {
303  if (edgelist[edgeoffsets[i]+l] == meshEdgeId)
304  {
305  extraDirEdgeIds[meshEdgeId] = i;
306  graph[1][meshEdgeId] = graphVertId++;
307  nExtraDirichlet += exp->GetEdgeNcoeffs(k)-2;
308  }
309  }
310  }
311  }
312  }
313  }
314 
315  // Low Energy preconditioner needs to know how many extra Dirichlet
316  // edges are on this process so store map in array.
317  map<int, int>::const_iterator mapConstIt;
318  m_extraDirEdges = Array<OneD, int>(extraDirEdgeIds.size(), -1);
319  for (mapConstIt = extraDirEdgeIds.begin(), i = 0;
320  mapConstIt != extraDirEdgeIds.end(); mapConstIt++)
321  {
322  meshEdgeId = mapConstIt->first;
323  m_extraDirEdges[i++] = meshEdgeId;
324  }
325 
326  // Now we have a list of all vertices and edges that are Dirichlet
327  // and not defined on the local partition as well as which processor
328  // they are stored on.
329  //
330  // Make a full list of all such entities on all processors and which
331  // processor they belong to.
332  for (i = 0; i < n; ++i)
333  {
334  vertcounts [i] = 0;
335  vertoffsets[i] = 0;
336  edgecounts [i] = 0;
337  edgeoffsets[i] = 0;
338  }
339 
340  vertcounts[p] = extraDirVertIds.size();
341  edgecounts[p] = extraDirEdgeIds.size();
342  vComm->AllReduce(vertcounts, LibUtilities::ReduceSum);
343  vComm->AllReduce(edgecounts, LibUtilities::ReduceSum);
344  nTotVerts = Vmath::Vsum(n, vertcounts, 1);
345  nTotEdges = Vmath::Vsum(n, edgecounts, 1);
346 
347  vertoffsets[0] = edgeoffsets[0] = 0;
348 
349  for (i = 1; i < n; ++i)
350  {
351  vertoffsets[i] = vertoffsets[i-1] + vertcounts[i-1];
352  edgeoffsets[i] = edgeoffsets[i-1] + edgecounts[i-1];
353  }
354 
355  Array<OneD, int> vertids (nTotVerts, 0);
356  Array<OneD, int> edgeids (nTotEdges, 0);
357  Array<OneD, int> vertprocs(nTotVerts, 0);
358  Array<OneD, int> edgeprocs(nTotEdges, 0);
359 
360  for (it = extraDirVertIds.begin(), i = 0;
361  it != extraDirVertIds.end(); ++it, ++i)
362  {
363  vertids [vertoffsets[p]+i] = it->first;
364  vertprocs[vertoffsets[p]+i] = it->second;
365  }
366 
367  for (it = extraDirEdgeIds.begin(), i = 0;
368  it != extraDirEdgeIds.end(); ++it, ++i)
369  {
370  edgeids [edgeoffsets[p]+i] = it->first;
371  edgeprocs[edgeoffsets[p]+i] = it->second;
372  }
373 
374  vComm->AllReduce(vertids, LibUtilities::ReduceSum);
375  vComm->AllReduce(vertprocs, LibUtilities::ReduceSum);
376  vComm->AllReduce(edgeids, LibUtilities::ReduceSum);
377  vComm->AllReduce(edgeprocs, LibUtilities::ReduceSum);
378 
379  // Set up list of vertices that need to be shared to other
380  // partitions
381  for (i = 0; i < nTotVerts; ++i)
382  {
383  if (vComm->GetRank() == vertprocs[i])
384  {
385  extraDirVerts.insert(vertids[i]);
386  }
387  }
388 
389  // Set up list of edges that need to be shared to other partitions
390  for (i = 0; i < nTotEdges; ++i)
391  {
392  if (vComm->GetRank() == edgeprocs[i])
393  {
394  extraDirEdges.insert(edgeids[i]);
395  }
396  }
397 
398  // Check between processes if the whole system is singular
399  int s = m_systemSingular ? 1 : 0;
400  vComm->AllReduce(s, LibUtilities::ReduceMin);
401  m_systemSingular = s == 1 ? true : false;
402 
403  // Find the minimum boundary vertex ID on each process
404  Array<OneD, int> bcminvertid(n, 0);
405  bcminvertid[p] = vMaxVertId;
406  vComm->AllReduce(bcminvertid, LibUtilities::ReduceMax);
407 
408  // Find the process rank with the minimum boundary vertex ID
409  int maxIdx = Vmath::Imax(n, bcminvertid, 1);
410 
411  // If the system is singular, the process with the maximum number of
412  // BCs will set a Dirichlet vertex to make system non-singular.
413  // Note: we find the process with maximum boundary regions to ensure
414  // we do not try to set a Dirichlet vertex on a partition with no
415  // intersection with the boundary.
416  meshVertId = 0;
417 
418  if (m_systemSingular && checkIfSystemSingular && maxIdx == p)
419  {
420  if (m_session->DefinesParameter("SingularVertex"))
421  {
422  m_session->LoadParameter("SingularVertex", meshVertId);
423  }
424  else if (vMaxVertId == -1)
425  {
426  // All boundaries are periodic.
427  meshVertId = locExpVector[0]->GetGeom()->GetVid(0);
428  }
429  else
430  {
431  // Set pinned vertex to that with minimum vertex ID to
432  // ensure consistency in parallel.
433  meshVertId = bcminvertid[p];
434  }
435 
436  if (graph[0].count(meshVertId) == 0)
437  {
438  graph[0][meshVertId] = graphVertId++;
439  }
440  }
441 
442  vComm->AllReduce(meshVertId, LibUtilities::ReduceSum);
443 
444  // When running in parallel, we need to ensure that the singular
445  // mesh vertex is communicated to any periodic vertices, otherwise
446  // the system may diverge.
447  if(m_systemSingular && checkIfSystemSingular)
448  {
449  // Firstly, we check that no other processors have this
450  // vertex. If they do, then we mark the vertex as also being
451  // Dirichlet.
452  if (maxIdx != p)
453  {
454  for (i = 0; i < locExpVector.size(); ++i)
455  {
456  for (j = 0; j < locExpVector[i]->GetNverts(); ++j)
457  {
458  if (locExpVector[i]->GetGeom()->GetVid(j) !=
459  meshVertId)
460  {
461  continue;
462  }
463 
464  if (graph[0].count(meshVertId) == 0)
465  {
466  graph[0][meshVertId] =
467  graphVertId++;
468  }
469  }
470  }
471  }
472 
473  // In the case that meshVertId is periodic with other vertices,
474  // this process and all other processes need to make sure that
475  // the periodic vertices are also marked as Dirichlet.
476  int gId;
477 
478  // At least one process (maxBCidx) will have already associated
479  // a graphVertId with meshVertId. Others won't even have any of
480  // the vertices. The logic below is designed to handle both
481  // cases.
482  if (graph[0].count(meshVertId) == 0)
483  {
484  gId = -1;
485  }
486  else
487  {
488  gId = graph[0][meshVertId];
489  }
490 
491  for (pIt = periodicVerts.begin();
492  pIt != periodicVerts.end(); ++pIt)
493  {
494  // Either the vertex is local to this processor (in which
495  // case it will be in the pIt->first position) or else
496  // meshVertId might be contained within another processor's
497  // vertex list. The if statement below covers both cases. If
498  // we find it, set as Dirichlet with the vertex id gId.
499  if (pIt->first == meshVertId)
500  {
501  graph[0][meshVertId] = gId < 0 ? graphVertId++ : gId;
502 
503  for (i = 0; i < pIt->second.size(); ++i)
504  {
505  if (pIt->second[i].isLocal)
506  {
507  graph[0][pIt->second[i].id] = graph[0][meshVertId];
508  }
509  }
510  }
511  else
512  {
513  bool found = false;
514  for (i = 0; i < pIt->second.size(); ++i)
515  {
516  if (pIt->second[i].id == meshVertId)
517  {
518  found = true;
519  break;
520  }
521  }
522 
523  if (found)
524  {
525  graph[0][pIt->first] = gId < 0 ? graphVertId++ : gId;
526 
527  for (i = 0; i < pIt->second.size(); ++i)
528  {
529  if (pIt->second[i].isLocal)
530  {
531  graph[0][pIt->second[i].id] = graph[0][pIt->first];
532  }
533  }
534  }
535  }
536  }
537  }
538 
539  // Add extra dirichlet boundary conditions to count.
540  m_numLocalDirBndCoeffs += nExtraDirichlet;
541  firstNonDirGraphVertId = graphVertId;
542 
543  typedef boost::adjacency_list<
544  boost::setS, boost::vecS, boost::undirectedS> BoostGraph;
545  BoostGraph boostGraphObj;
546 
547  vector<map<int,int> > tempGraph(3);
548  map<int, int> vwgts_map;
549  Array<OneD, int> localVerts;
550  Array<OneD, int> localEdges;
551  Array<OneD, int> localFaces;
552 
553  int tempGraphVertId = 0;
554  int localVertOffset = 0;
555  int localEdgeOffset = 0;
556  int localFaceOffset = 0;
557  int nTotalVerts = 0;
558  int nTotalEdges = 0;
559  int nTotalFaces = 0;
560  int nVerts;
561  int nEdges;
562  int nFaces;
563  int vertCnt;
564  int edgeCnt;
565  int faceCnt;
566 
567  m_numNonDirVertexModes = 0;
568  m_numNonDirEdges = 0;
569  m_numNonDirFaces = 0;
570  m_numNonDirFaceModes = 0;
571  m_numNonDirFaceModes = 0;
572  m_numLocalBndCoeffs = 0;
573 
574  map<int,int> EdgeSize;
575  map<int,int> FaceSize;
576 
577  /// - Count verts, edges, face and add up edges and face sizes
578  for(i = 0; i < locExpVector.size(); ++i)
579  {
580  exp = locExpVector[locExp.GetOffset_Elmt_Id(i)];
581  nTotalVerts += exp->GetNverts();
582  nTotalEdges += exp->GetNedges();
583  nTotalFaces += exp->GetNfaces();
584 
585  nEdges = exp->GetNedges();
586  for(j = 0; j < nEdges; ++j)
587  {
588  meshEdgeId = exp->GetGeom()->GetEid(j);
589  if (EdgeSize.count(meshEdgeId) > 0)
590  {
591  EdgeSize[meshEdgeId] =
592  min(EdgeSize[meshEdgeId],
593  exp->GetEdgeNcoeffs(j) - 2);
594  }
595  else
596  {
597  EdgeSize[meshEdgeId] = exp->GetEdgeNcoeffs(j) - 2;
598  }
599  }
600 
601  nFaces = exp->GetNfaces();
602  faceCnt = 0;
603  for(j = 0; j < nFaces; ++j)
604  {
605  meshFaceId = exp->GetGeom()->GetFid(j);
606  if (FaceSize.count(meshFaceId) > 0)
607  {
608  FaceSize[meshFaceId] =
609  min(FaceSize[meshFaceId],
610  exp->GetFaceIntNcoeffs(j));
611  }
612  else
613  {
614  FaceSize[meshFaceId] = exp->GetFaceIntNcoeffs(j);
615  }
616  FaceSize[meshFaceId] = exp->GetFaceIntNcoeffs(j);
617  }
618  }
619 
620  /// - Periodic vertices
621  for (pIt = periodicVerts.begin(); pIt != periodicVerts.end(); ++pIt)
622  {
623  meshVertId = pIt->first;
624 
625  // This periodic vertex is joined to a Dirichlet condition.
626  if (graph[0].count(pIt->first) != 0)
627  {
628  for (i = 0; i < pIt->second.size(); ++i)
629  {
630  meshVertId2 = pIt->second[i].id;
631  if (graph[0].count(meshVertId2) == 0 &&
632  pIt->second[i].isLocal)
633  {
634  graph[0][meshVertId2] =
635  graph[0][meshVertId];
636  }
637  }
638  continue;
639  }
640 
641  // One of the attached vertices is Dirichlet.
642  bool isDirichlet = false;
643  for (i = 0; i < pIt->second.size(); ++i)
644  {
645  if (!pIt->second[i].isLocal)
646  {
647  continue;
648  }
649 
650  meshVertId2 = pIt->second[i].id;
651  if (graph[0].count(meshVertId2) > 0)
652  {
653  isDirichlet = true;
654  break;
655  }
656  }
657 
658  if (isDirichlet)
659  {
660  graph[0][meshVertId] =
661  graph[0][pIt->second[i].id];
662 
663  for (j = 0; j < pIt->second.size(); ++j)
664  {
665  meshVertId2 = pIt->second[i].id;
666  if (j == i || !pIt->second[j].isLocal ||
667  graph[0].count(meshVertId2) > 0)
668  {
669  continue;
670  }
671 
672  graph[0][meshVertId2] =
673  graph[0][pIt->second[i].id];
674  }
675 
676  continue;
677  }
678 
679  // Otherwise, see if a vertex ID has already been set.
680  for (i = 0; i < pIt->second.size(); ++i)
681  {
682  if (!pIt->second[i].isLocal)
683  {
684  continue;
685  }
686 
687  if (tempGraph[0].count(pIt->second[i].id) > 0)
688  {
689  break;
690  }
691  }
692 
693  if (i == pIt->second.size())
694  {
695  boost::add_vertex(boostGraphObj);
696  tempGraph[0][meshVertId] = tempGraphVertId++;
697  m_numNonDirVertexModes++;
698  }
699  else
700  {
701  tempGraph[0][meshVertId] = tempGraph[0][pIt->second[i].id];
702  }
703  }
704 
705  // Store the temporary graph vertex id's of all element edges and
706  // vertices in these 3 arrays below
707  localVerts = Array<OneD, int>(nTotalVerts,-1);
708  localEdges = Array<OneD, int>(nTotalEdges,-1);
709  localFaces = Array<OneD, int>(nTotalFaces,-1);
710 
711  // Set up vertex numbering
712  for(i = 0; i < locExpVector.size(); ++i)
713  {
714  exp = locExpVector[i];
715  vertCnt = 0;
716  nVerts = exp->GetNverts();
717  for(j = 0; j < nVerts; ++j)
718  {
719  meshVertId = exp->GetGeom()->GetVid(j);
720  if(graph[0].count(meshVertId) == 0)
721  {
722  if(tempGraph[0].count(meshVertId) == 0)
723  {
724  boost::add_vertex(boostGraphObj);
725  tempGraph[0][meshVertId] = tempGraphVertId++;
726  m_numNonDirVertexModes+=1;
727  }
728  localVerts[localVertOffset+vertCnt++] = tempGraph[0][meshVertId];
729  vwgts_map[ tempGraph[0][meshVertId] ] = 1;
730  }
731  }
732 
733  localVertOffset+=nVerts;
734  }
735 
736  /// - Periodic edges
737  for (pIt = periodicEdges.begin(); pIt != periodicEdges.end(); ++pIt)
738  {
739  meshEdgeId = pIt->first;
740 
741  // This periodic edge is joined to a Dirichlet condition.
742  if (graph[1].count(pIt->first) != 0)
743  {
744  for (i = 0; i < pIt->second.size(); ++i)
745  {
746  meshEdgeId2 = pIt->second[i].id;
747  if (graph[1].count(meshEdgeId2) == 0 &&
748  pIt->second[i].isLocal)
749  {
750  graph[1][meshEdgeId2] =
751  graph[1][meshEdgeId];
752  }
753  }
754  continue;
755  }
756 
757  // One of the attached edges is Dirichlet.
758  bool isDirichlet = false;
759  for (i = 0; i < pIt->second.size(); ++i)
760  {
761  if (!pIt->second[i].isLocal)
762  {
763  continue;
764  }
765 
766  meshEdgeId2 = pIt->second[i].id;
767  if (graph[1].count(meshEdgeId2) > 0)
768  {
769  isDirichlet = true;
770  break;
771  }
772  }
773 
774  if (isDirichlet)
775  {
776  graph[1][meshEdgeId] =
777  graph[1][pIt->second[i].id];
778 
779  for (j = 0; j < pIt->second.size(); ++j)
780  {
781  meshEdgeId2 = pIt->second[i].id;
782  if (j == i || !pIt->second[j].isLocal ||
783  graph[1].count(meshEdgeId2) > 0)
784  {
785  continue;
786  }
787 
788  graph[1][meshEdgeId2] =
789  graph[1][pIt->second[i].id];
790  }
791 
792  continue;
793  }
794 
795  // Otherwise, see if a edge ID has already been set.
796  for (i = 0; i < pIt->second.size(); ++i)
797  {
798  if (!pIt->second[i].isLocal)
799  {
800  continue;
801  }
802 
803  if (tempGraph[1].count(pIt->second[i].id) > 0)
804  {
805  break;
806  }
807  }
808 
809  if (i == pIt->second.size())
810  {
811  boost::add_vertex(boostGraphObj);
812  tempGraph[1][meshEdgeId] = tempGraphVertId++;
813  m_numNonDirEdgeModes += EdgeSize[meshEdgeId];
814  m_numNonDirEdges++;
815  }
816  else
817  {
818  tempGraph[1][meshEdgeId] = tempGraph[1][pIt->second[i].id];
819  }
820  }
821 
822  int nEdgeIntCoeffs, nFaceIntCoeffs;
823 
824  // Set up edge numbering
825  for(i = 0; i < locExpVector.size(); ++i)
826  {
827  exp = locExpVector[i];
828  edgeCnt = 0;
829  nEdges = exp->GetNedges();
830 
831  for(j = 0; j < nEdges; ++j)
832  {
833  nEdgeIntCoeffs = exp->GetEdgeNcoeffs(j) - 2;
834  meshEdgeId = exp->GetGeom()->GetEid(j);
835  if(graph[1].count(meshEdgeId) == 0)
836  {
837  if(tempGraph[1].count(meshEdgeId) == 0)
838  {
839  boost::add_vertex(boostGraphObj);
840  tempGraph[1][meshEdgeId] = tempGraphVertId++;
841  m_numNonDirEdgeModes+=nEdgeIntCoeffs;
842 
843  m_numNonDirEdges++;
844  }
845  localEdges[localEdgeOffset+edgeCnt++] = tempGraph[1][meshEdgeId];
846  vwgts_map[ tempGraph[1][meshEdgeId] ] = nEdgeIntCoeffs;
847  }
848  }
849 
850  localEdgeOffset+=nEdges;
851  }
852 
853  /// - Periodic faces
854  for (pIt = periodicFaces.begin(); pIt != periodicFaces.end(); ++pIt)
855  {
856  if (!pIt->second[0].isLocal)
857  {
858  // The face mapped to is on another process.
859  meshFaceId = pIt->first;
860  ASSERTL0(graph[2].count(meshFaceId) == 0,
861  "This periodic boundary edge has been specified before");
862  boost::add_vertex(boostGraphObj);
863  tempGraph[2][meshFaceId] = tempGraphVertId++;
864  nFaceIntCoeffs = FaceSize[meshFaceId];
865  m_numNonDirFaceModes+=nFaceIntCoeffs;
866  m_numNonDirFaces++;
867  }
868  else if (pIt->first < pIt->second[0].id)
869  {
870  ASSERTL0(graph[2].count(pIt->first) == 0,
871  "This periodic boundary face has been specified before");
872  ASSERTL0(graph[2].count(pIt->second[0].id) == 0,
873  "This periodic boundary face has been specified before");
874 
875  boost::add_vertex(boostGraphObj);
876  tempGraph[2][pIt->first] = tempGraphVertId;
877  tempGraph[2][pIt->second[0].id] = tempGraphVertId++;
878  nFaceIntCoeffs = FaceSize[pIt->first];
879  m_numNonDirFaceModes+=nFaceIntCoeffs;
880  m_numNonDirFaces++;
881  }
882  }
883 
884  // setup face numbering
885  for(i = 0; i < locExpVector.size(); ++i)
886  {
887  exp = locExpVector[i];
888  nFaces = exp->GetNfaces();
889  faceCnt = 0;
890  for(j = 0; j < nFaces; ++j)
891  {
892  nFaceIntCoeffs = exp->GetFaceIntNcoeffs(j);
893  meshFaceId = exp->GetGeom()->GetFid(j);
894  if(graph[2].count(meshFaceId) == 0)
895  {
896  if(tempGraph[2].count(meshFaceId) == 0)
897  {
898  boost::add_vertex(boostGraphObj);
899  tempGraph[2][meshFaceId] = tempGraphVertId++;
900  m_numNonDirFaceModes+=nFaceIntCoeffs;
901 
902  m_numNonDirFaces++;
903  }
904  localFaces[localFaceOffset+faceCnt++] = tempGraph[2][meshFaceId];
905  vwgts_map[ tempGraph[2][meshFaceId] ] = nFaceIntCoeffs;
906  }
907  }
908  m_numLocalBndCoeffs += exp->NumBndryCoeffs();
909 
910  localFaceOffset+=nFaces;
911  }
912 
913  localVertOffset=0;
914  localEdgeOffset=0;
915  localFaceOffset=0;
916  for(i = 0; i < locExpVector.size(); ++i)
917  {
918  exp = locExpVector[i];
919  nVerts = exp->GetNverts();
920  nEdges = exp->GetNedges();
921  nFaces = exp->GetNfaces();
922 
923  // Now loop over all local faces, edges and vertices of this
924  // element and define that all other faces, edges and verices of
925  // this element are adjacent to them.
926 
927  // Vertices
928  for(j = 0; j < nVerts; j++)
929  {
930  if(localVerts[j+localVertOffset]==-1)
931  {
932  break;
933  }
934  // associate to other vertices
935  for(k = 0; k < nVerts; k++)
936  {
937  if(localVerts[k+localVertOffset]==-1)
938  {
939  break;
940  }
941  if(k!=j)
942  {
943  boost::add_edge( (size_t) localVerts[j+localVertOffset],
944  (size_t) localVerts[k+localVertOffset],boostGraphObj);
945  }
946  }
947  // associate to other edges
948  for(k = 0; k < nEdges; k++)
949  {
950  if(localEdges[k+localEdgeOffset]==-1)
951  {
952  break;
953  }
954  boost::add_edge( (size_t) localVerts[j+localVertOffset],
955  (size_t) localEdges[k+localEdgeOffset],boostGraphObj);
956  }
957  // associate to other faces
958  for(k = 0; k < nFaces; k++)
959  {
960  if(localFaces[k+localFaceOffset]==-1)
961  {
962  break;
963  }
964  boost::add_edge( (size_t) localVerts[j+localVertOffset],
965  (size_t) localFaces[k+localFaceOffset],boostGraphObj);
966  }
967  }
968 
969  // Edges
970  for(j = 0; j < nEdges; j++)
971  {
972  if(localEdges[j+localEdgeOffset]==-1)
973  {
974  break;
975  }
976  // Associate to other edges
977  for(k = 0; k < nEdges; k++)
978  {
979  if(localEdges[k+localEdgeOffset]==-1)
980  {
981  break;
982  }
983  if(k!=j)
984  {
985  boost::add_edge( (size_t) localEdges[j+localEdgeOffset],
986  (size_t) localEdges[k+localEdgeOffset],boostGraphObj);
987  }
988  }
989  // Associate to vertices
990  for(k = 0; k < nVerts; k++)
991  {
992  if(localVerts[k+localVertOffset]==-1)
993  {
994  break;
995  }
996  boost::add_edge( (size_t) localEdges[j+localEdgeOffset],
997  (size_t) localVerts[k+localVertOffset],boostGraphObj);
998  }
999  // Associate to faces
1000  for(k = 0; k < nFaces; k++)
1001  {
1002  if(localFaces[k+localFaceOffset]==-1)
1003  {
1004  break;
1005  }
1006  boost::add_edge( (size_t) localEdges[j+localEdgeOffset],
1007  (size_t) localFaces[k+localFaceOffset],boostGraphObj);
1008  }
1009  }
1010 
1011  // Faces
1012  for(j = 0; j < nFaces; j++)
1013  {
1014  if(localFaces[j+localFaceOffset]==-1)
1015  {
1016  break;
1017  }
1018  // Associate to other faces
1019  for(k = 0; k < nFaces; k++)
1020  {
1021  if(localFaces[k+localFaceOffset]==-1)
1022  {
1023  break;
1024  }
1025  if(k!=j)
1026  {
1027  boost::add_edge( (size_t) localFaces[j+localFaceOffset],
1028  (size_t) localFaces[k+localFaceOffset],boostGraphObj);
1029  }
1030  }
1031  // Associate to vertices
1032  for(k = 0; k < nVerts; k++)
1033  {
1034  if(localVerts[k+localVertOffset]==-1)
1035  {
1036  break;
1037  }
1038  boost::add_edge( (size_t) localFaces[j+localFaceOffset],
1039  (size_t) localVerts[k+localVertOffset],boostGraphObj);
1040  }
1041  // Associate to edges
1042  for(k = 0; k < nEdges; k++)
1043  {
1044  if(localEdges[k+localEdgeOffset]==-1)
1045  {
1046  break;
1047  }
1048  boost::add_edge( (size_t) localFaces[j+localFaceOffset],
1049  (size_t) localEdges[k+localEdgeOffset],boostGraphObj);
1050  }
1051  }
1052 
1053  localVertOffset+=nVerts;
1054  localEdgeOffset+=nEdges;
1055  localFaceOffset+=nFaces;
1056  }
1057 
1058  // Container to store vertices of the graph which correspond to
1059  // degrees of freedom along the boundary and periodic BCs.
1060  set<int> partVerts;
1061 
1062  if (m_solnType == eIterativeMultiLevelStaticCond ||
1063  m_solnType == eXxtMultiLevelStaticCond)
1064  {
1065  vector<long> procVerts, procEdges, procFaces;
1066  set <int> foundVerts, foundEdges, foundFaces;
1067 
1068  // Loop over element and construct the procVerts and procEdges
1069  // vectors, which store the geometry IDs of mesh vertices and
1070  // edges respectively which are local to this process.
1071  for(i = cnt = 0; i < locExpVector.size(); ++i)
1072  {
1073  int elmtid = locExp.GetOffset_Elmt_Id(i);
1074  exp = locExpVector[elmtid];
1075  for (j = 0; j < exp->GetNverts(); ++j)
1076  {
1077  int vid = exp->GetGeom()->GetVid(j)+1;
1078  if (foundVerts.count(vid) == 0)
1079  {
1080  procVerts.push_back(vid);
1081  foundVerts.insert(vid);
1082  }
1083  }
1084 
1085  for (j = 0; j < exp->GetNedges(); ++j)
1086  {
1087  int eid = exp->GetGeom()->GetEid(j)+1;
1088 
1089  if (foundEdges.count(eid) == 0)
1090  {
1091  procEdges.push_back(eid);
1092  foundEdges.insert(eid);
1093  }
1094  }
1095 
1096  for (j = 0; j < exp->GetNfaces(); ++j)
1097  {
1098  int fid = exp->GetGeom()->GetFid(j)+1;
1099 
1100  if (foundFaces.count(fid) == 0)
1101  {
1102  procFaces.push_back(fid);
1103  foundFaces.insert(fid);
1104  }
1105  }
1106  }
1107 
1108  int unique_verts = foundVerts.size();
1109  int unique_edges = foundEdges.size();
1110  int unique_faces = foundFaces.size();
1111 
1112  // Now construct temporary GS objects. These will be used to
1113  // populate the arrays tmp3 and tmp4 with the multiplicity of
1114  // the vertices and edges respectively to identify those
1115  // vertices and edges which are located on partition boundary.
1116  Array<OneD, long> vertArray(unique_verts, &procVerts[0]);
1117  Gs::gs_data *tmp1 = Gs::Init(vertArray, vComm);
1118  Array<OneD, NekDouble> tmp4(unique_verts, 1.0);
1119  Array<OneD, NekDouble> tmp5(unique_edges, 1.0);
1120  Array<OneD, NekDouble> tmp6(unique_faces, 1.0);
1121  Gs::Gather(tmp4, Gs::gs_add, tmp1);
1122  Gs::Finalise(tmp1);
1123 
1124  if (unique_edges > 0)
1125  {
1126  Array<OneD, long> edgeArray(unique_edges, &procEdges[0]);
1127  Gs::gs_data *tmp2 = Gs::Init(edgeArray, vComm);
1128  Gs::Gather(tmp5, Gs::gs_add, tmp2);
1129  Gs::Finalise(tmp2);
1130  }
1131 
1132  if (unique_faces > 0)
1133  {
1134  Array<OneD, long> faceArray(unique_faces, &procFaces[0]);
1135  Gs::gs_data *tmp3 = Gs::Init(faceArray, vComm);
1136  Gs::Gather(tmp6, Gs::gs_add, tmp3);
1137  Gs::Finalise(tmp3);
1138  }
1139 
1140  // Finally, fill the partVerts set with all non-Dirichlet
1141  // vertices which lie on a partition boundary.
1142  for (i = 0; i < unique_verts; ++i)
1143  {
1144  if (tmp4[i] > 1.0)
1145  {
1146  if (graph[0].count(procVerts[i]-1) == 0)
1147  {
1148  partVerts.insert(tempGraph[0][procVerts[i]-1]);
1149  }
1150  }
1151  }
1152 
1153  for (i = 0; i < unique_edges; ++i)
1154  {
1155  if (tmp5[i] > 1.0)
1156  {
1157  if (graph[1].count(procEdges[i]-1) == 0)
1158  {
1159  partVerts.insert(tempGraph[1][procEdges[i]-1]);
1160  }
1161  }
1162  }
1163 
1164  for (i = 0; i < unique_faces; ++i)
1165  {
1166  if (tmp6[i] > 1.0)
1167  {
1168  if (graph[2].count(procFaces[i]-1) == 0)
1169  {
1170  partVerts.insert(tempGraph[2][procFaces[i]-1]);
1171  }
1172  }
1173  }
1174 
1175  // Now fill with all vertices on periodic BCs
1176  for (pIt = periodicVerts.begin(); pIt != periodicVerts.end(); ++pIt)
1177  {
1178  if (graph[0].count(pIt->first) == 0)
1179  {
1180  partVerts.insert(tempGraph[0][pIt->first]);
1181  }
1182  }
1183  for (pIt = periodicEdges.begin(); pIt != periodicEdges.end(); ++pIt)
1184  {
1185  if (graph[1].count(pIt->first) == 0)
1186  {
1187  partVerts.insert(tempGraph[1][pIt->first]);
1188  }
1189  }
1190  for (pIt = periodicFaces.begin(); pIt != periodicFaces.end(); ++pIt)
1191  {
1192  if (graph[2].count(pIt->first) == 0)
1193  {
1194  partVerts.insert(tempGraph[2][pIt->first]);
1195  }
1196  }
1197  }
1198 
1199  int nGraphVerts = tempGraphVertId;
1200  Array<OneD, int> perm(nGraphVerts);
1201  Array<OneD, int> iperm(nGraphVerts);
1202  Array<OneD, int> vwgts(nGraphVerts);
1203  ASSERTL1(vwgts_map.size()==nGraphVerts,"Non matching dimensions");
1204  for(i = 0; i < nGraphVerts; ++i)
1205  {
1206  vwgts[i] = vwgts_map[i];
1207  }
1208 
1209  if(nGraphVerts)
1210  {
1211  switch(m_solnType)
1212  {
1213  case eDirectFullMatrix:
1214  case eIterativeFull:
1215  case eIterativeStaticCond:
1216  case ePETScStaticCond:
1217  case ePETScFullMatrix:
1218  case eXxtFullMatrix:
1219  case eXxtStaticCond:
1220  {
1221  NoReordering(boostGraphObj,perm,iperm);
1222  break;
1223  }
1224 
1225  case eDirectStaticCond:
1226  {
1227  CuthillMckeeReordering(boostGraphObj,perm,iperm);
1228  break;
1229  }
1230 
1231  case ePETScMultiLevelStaticCond:
1232  case eDirectMultiLevelStaticCond:
1233  case eIterativeMultiLevelStaticCond:
1234  case eXxtMultiLevelStaticCond:
1235  {
1236  MultiLevelBisectionReordering(
1237  boostGraphObj, perm, iperm, bottomUpGraph,
1238  partVerts, mdswitch);
1239  break;
1240  }
1241  default:
1242  {
1243  ASSERTL0(false,
1244  "Unrecognised solution type: " + std::string(
1245  GlobalSysSolnTypeMap[m_solnType]));
1246  }
1247  }
1248  }
1249 
1250  // For parallel multi-level static condensation determine the lowest
1251  // static condensation level amongst processors.
1252  if (m_solnType == eIterativeMultiLevelStaticCond ||
1253  m_solnType == eXxtMultiLevelStaticCond)
1254  {
1255  m_lowestStaticCondLevel = bottomUpGraph->GetNlevels()-1;
1256  vComm->AllReduce(m_lowestStaticCondLevel,
1257  LibUtilities::ReduceMax);
1258  }
1259  else
1260  {
1261  m_lowestStaticCondLevel = 0;
1262  }
1263 
1264  /**
1265  * STEP 4: Fill the #graph[0] and
1266  * #graph[1] with the optimal ordering from boost.
1267  */
1268  map<int, int>::iterator mapIt;
1269  for(mapIt = tempGraph[0].begin(); mapIt != tempGraph[0].end(); mapIt++)
1270  {
1271  graph[0][mapIt->first] = iperm[mapIt->second] + graphVertId;
1272  }
1273  for(mapIt = tempGraph[1].begin(); mapIt != tempGraph[1].end(); mapIt++)
1274  {
1275  graph[1][mapIt->first] = iperm[mapIt->second] + graphVertId;
1276  }
1277  for(mapIt = tempGraph[2].begin(); mapIt != tempGraph[2].end(); mapIt++)
1278  {
1279  graph[2][mapIt->first] = iperm[mapIt->second] + graphVertId;
1280  }
1281 
1282  return nGraphVerts;
1283  }
ASSERTL0
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
Nektar::MultiRegions::AssemblyMap::m_systemSingular
bool m_systemSingular
Flag indicating if the system is singular or not.
Definition: AssemblyMap.h:320
Nektar::MultiRegions::MultiLevelBisectionReordering
void MultiLevelBisectionReordering(const BoostGraph &graph, Array< OneD, int > &perm, Array< OneD, int > &iperm, BottomUpSubStructuredGraphSharedPtr &substructgraph, std::set< int > partVerts, int mdswitch)
Definition: SubStructuredGraph.cpp:862
Nektar::MultiRegions::AssemblyMap::m_comm
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:306
Gs::Gather
static void Gather(Nektar::Array< OneD, NekDouble > pU, gs_op pOp, gs_data *pGsh, Nektar::Array< OneD, NekDouble > pBuffer=NullNekDouble1DArray)
Performs a gather-scatter operation of the provided values.
Definition: GsLib.hpp:224
Gs::Finalise
static void Finalise(gs_data *pGsh)
Deallocates the GSLib mapping data.
Definition: GsLib.hpp:206
Nektar::MultiRegions::AssemblyMapCG::m_extraDirEdges
Array< OneD, int > m_extraDirEdges
Extra dirichlet edges in parallel.
Definition: AssemblyMapCG.h:141
Gs::Init
static gs_data * Init(const Nektar::Array< OneD, long > pId, const LibUtilities::CommSharedPtr &pComm)
Initialise Gather-Scatter map.
Definition: GsLib.hpp:151
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirVertexModes
int m_numNonDirVertexModes
Number of non Dirichlet vertex modes.
Definition: AssemblyMapCG.h:125
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirEdges
int m_numNonDirEdges
Number of Dirichlet edges.
Definition: AssemblyMapCG.h:135
Nektar::MultiRegions::CuthillMckeeReordering
void CuthillMckeeReordering(const BoostGraph &graph, Array< OneD, int > &perm, Array< OneD, int > &iperm)
Definition: SubStructuredGraph.cpp:839
Nektar::MultiRegions::AssemblyMapCG::m_numLocalBndCondCoeffs
int m_numLocalBndCondCoeffs
Number of local boundary condition coefficients.
Definition: AssemblyMapCG.h:139
Vmath::Imax
int Imax(int n, const T *x, const int incx)
Return the index of the maximum element in x.
Definition: Vmath.cpp:741
Nektar::LocalRegions::ExpansionVector
std::vector< ExpansionSharedPtr > ExpansionVector
Definition: Expansion.h:70
Nektar::MultiRegions::GlobalSysSolnTypeMap
const char *const GlobalSysSolnTypeMap[]
Definition: MultiRegions.hpp:96
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirFaceModes
int m_numNonDirFaceModes
Number of non Dirichlet face modes.
Definition: AssemblyMapCG.h:129
Nektar::LibUtilities::CommSharedPtr
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
Nektar::MultiRegions::AssemblyMap::m_solnType
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:363
Nektar::MultiRegions::AssemblyMapCG::m_numDirFaces
int m_numDirFaces
Number of Dirichlet faces.
Definition: AssemblyMapCG.h:133
Nektar::LocalRegions::ExpansionSharedPtr
boost::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
Nektar::MultiRegions::AssemblyMap::m_lowestStaticCondLevel
int m_lowestStaticCondLevel
Lowest static condensation level.
Definition: AssemblyMap.h:397
Nektar::MultiRegions::AssemblyMap::m_numLocalDirBndCoeffs
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:316
Nektar::iterator
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
Definition: MatrixOperations.cpp:96
Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffs
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:312
Nektar::MultiRegions::AssemblyMapCG::m_numDirEdges
int m_numDirEdges
Number of Dirichlet edges.
Definition: AssemblyMapCG.h:131
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirEdgeModes
int m_numNonDirEdgeModes
Number of non Dirichlet edge modes.
Definition: AssemblyMapCG.h:127
Gs::gs_data
Definition: GsLib.hpp:117
Nektar::MultiRegions::AssemblyMap::m_session
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:303
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirFaces
int m_numNonDirFaces
Number of Dirichlet faces.
Definition: AssemblyMapCG.h:137
Vmath::Vsum
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.cpp:723
ASSERTL1
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
Nektar::MultiRegions::NoReordering
void NoReordering(const BoostGraph &graph, Array< OneD, int > &perm, Array< OneD, int > &iperm)
Definition: SubStructuredGraph.cpp:1074
std::map<int, std::vector<ExtraDirDof> >& Nektar::MultiRegions::AssemblyMapCG::GetExtraDirDofs ( )
inline

Definition at line 108 of file AssemblyMapCG.h.

References m_extraDirDofs.

109  {
110  return m_extraDirDofs;
111  }
Nektar::MultiRegions::AssemblyMapCG::m_extraDirDofs
std::map< int, std::vector< ExtraDirDof > > m_extraDirDofs
Map indicating degrees of freedom which are Dirichlet but whose value is stored on another processor...
Definition: AssemblyMapCG.h:148
void Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap ( const ExpList &  locExp,
const PeriodicMap &  perVerts = NullPeriodicMap,
const PeriodicMap &  perEdges = NullPeriodicMap,
const PeriodicMap &  perFaces = NullPeriodicMap 
)
protected

Sets up the global to universal mapping of degrees of freedom across processors.

Definition at line 1993 of file AssemblyMapCG.cpp.

References Nektar::MultiRegions::DeterminePeriodicEdgeOrientId(), Nektar::MultiRegions::DeterminePeriodicFaceOrient(), Nektar::MultiRegions::ExpList::GetCoeff_Offset(), Nektar::MultiRegions::ExpList::GetExp(), Gs::Init(), Nektar::MultiRegions::AssemblyMap::m_bndGsh, Nektar::MultiRegions::AssemblyMap::m_comm, Nektar::MultiRegions::AssemblyMap::m_globalToUniversalBndMap, Nektar::MultiRegions::AssemblyMap::m_globalToUniversalBndMapUnique, m_globalToUniversalMap, m_globalToUniversalMapUnique, Nektar::MultiRegions::AssemblyMap::m_gsh, m_localToGlobalMap, Nektar::MultiRegions::AssemblyMap::m_numGlobalBndCoeffs, Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs, Nektar::LibUtilities::ReduceMax, Nektar::LibUtilities::ReduceSum, Gs::Unique(), and Vmath::Zero().

Referenced by AssemblyMapCG().

1998  {
1999  LocalRegions::ExpansionSharedPtr exp;
2000  int nVert = 0;
2001  int nEdge = 0;
2002  int nFace = 0;
2003  int maxEdgeDof = 0;
2004  int maxFaceDof = 0;
2005  int maxIntDof = 0;
2006  int dof = 0;
2007  int cnt;
2008  int i,j,k;
2009  int meshVertId;
2010  int meshEdgeId;
2011  int meshFaceId;
2012  int elementId;
2013  int vGlobalId;
2014  int maxBndGlobalId = 0;
2015  StdRegions::Orientation edgeOrient;
2016  StdRegions::Orientation faceOrient;
2017  Array<OneD, unsigned int> edgeInteriorMap;
2018  Array<OneD, int> edgeInteriorSign;
2019  Array<OneD, unsigned int> faceInteriorMap;
2020  Array<OneD, int> faceInteriorSign;
2021  Array<OneD, unsigned int> interiorMap;
2022  PeriodicMap::const_iterator pIt;
2023 
2024  const LocalRegions::ExpansionVector &locExpVector = *(locExp.GetExp());
2025  LibUtilities::CommSharedPtr vCommRow = m_comm->GetRowComm();
2026 
2027  m_globalToUniversalMap = Nektar::Array<OneD, int>(m_numGlobalCoeffs, -1);
2028  m_globalToUniversalMapUnique = Nektar::Array<OneD, int>(m_numGlobalCoeffs, -1);
2029  m_globalToUniversalBndMap = Nektar::Array<OneD, int>(m_numGlobalBndCoeffs, -1);
2030  m_globalToUniversalBndMapUnique = Nektar::Array<OneD, int>(m_numGlobalBndCoeffs, -1);
2031 
2032  // Loop over all the elements in the domain to gather mesh data
2033  for(i = 0; i < locExpVector.size(); ++i)
2034  {
2035  exp = locExpVector[i];
2036  nVert += exp->GetNverts();
2037  nEdge += exp->GetNedges();
2038  nFace += exp->GetNfaces();
2039  // Loop over all edges (and vertices) of element i
2040  for(j = 0; j < exp->GetNedges(); ++j)
2041  {
2042  dof = exp->GetEdgeNcoeffs(j)-2;
2043  maxEdgeDof = (dof > maxEdgeDof ? dof : maxEdgeDof);
2044  }
2045  for(j = 0; j < exp->GetNfaces(); ++j)
2046  {
2047  dof = exp->GetFaceIntNcoeffs(j);
2048  maxFaceDof = (dof > maxFaceDof ? dof : maxFaceDof);
2049  }
2050  exp->GetInteriorMap(interiorMap);
2051  dof = interiorMap.num_elements();
2052  maxIntDof = (dof > maxIntDof ? dof : maxIntDof);
2053  }
2054 
2055  // Tell other processes about how many dof we have
2056  vCommRow->AllReduce(nVert, LibUtilities::ReduceSum);
2057  vCommRow->AllReduce(nEdge, LibUtilities::ReduceSum);
2058  vCommRow->AllReduce(nFace, LibUtilities::ReduceSum);
2059  vCommRow->AllReduce(maxEdgeDof, LibUtilities::ReduceMax);
2060  vCommRow->AllReduce(maxFaceDof, LibUtilities::ReduceMax);
2061  vCommRow->AllReduce(maxIntDof, LibUtilities::ReduceMax);
2062 
2063  // Assemble global to universal mapping for this process
2064  for(i = 0; i < locExpVector.size(); ++i)
2065  {
2066  exp = locExpVector[i];
2067  cnt = locExp.GetCoeff_Offset(i);
2068 
2069  // Loop over all vertices of element i
2070  for(j = 0; j < exp->GetNverts(); ++j)
2071  {
2072  meshVertId = exp->GetGeom()->GetVid(j);
2073  vGlobalId = m_localToGlobalMap[cnt+exp->GetVertexMap(j)];
2074 
2075  pIt = perVerts.find(meshVertId);
2076  if (pIt != perVerts.end())
2077  {
2078  for (k = 0; k < pIt->second.size(); ++k)
2079  {
2080  meshVertId = min(meshVertId, pIt->second[k].id);
2081  }
2082  }
2083 
2084  m_globalToUniversalMap[vGlobalId] = meshVertId + 1;
2085  m_globalToUniversalBndMap[vGlobalId]=m_globalToUniversalMap[vGlobalId];
2086  maxBndGlobalId = (vGlobalId > maxBndGlobalId ? vGlobalId : maxBndGlobalId);
2087  }
2088 
2089  // Loop over all edges of element i
2090  for(j = 0; j < exp->GetNedges(); ++j)
2091  {
2092  meshEdgeId = exp->GetGeom()->GetEid(j);
2093  pIt = perEdges.find(meshEdgeId);
2094  edgeOrient = exp->GetGeom()->GetEorient(j);
2095 
2096  if (pIt != perEdges.end())
2097  {
2098  pair<int, StdRegions::Orientation> idOrient =
2099  DeterminePeriodicEdgeOrientId(
2100  meshEdgeId, edgeOrient, pIt->second);
2101  meshEdgeId = idOrient.first;
2102  edgeOrient = idOrient.second;
2103  }
2104 
2105  exp->GetEdgeInteriorMap(j,edgeOrient,edgeInteriorMap,edgeInteriorSign);
2106  dof = exp->GetEdgeNcoeffs(j)-2;
2107 
2108  // Set the global DOF's for the interior modes of edge j
2109  // run backwards because of variable P case "ghost" modes
2110  for(k = dof-1; k >= 0; --k)
2111  {
2112  vGlobalId = m_localToGlobalMap[cnt+edgeInteriorMap[k]];
2113  m_globalToUniversalMap[vGlobalId]
2114  = nVert + meshEdgeId * maxEdgeDof + k + 1;
2115  m_globalToUniversalBndMap[vGlobalId]=m_globalToUniversalMap[vGlobalId];
2116  maxBndGlobalId = (vGlobalId > maxBndGlobalId ? vGlobalId : maxBndGlobalId);
2117  }
2118  }
2119 
2120  // Loop over all faces of element i
2121  for(j = 0; j < exp->GetNfaces(); ++j)
2122  {
2123  faceOrient = exp->GetGeom()->GetForient(j);
2124 
2125  meshFaceId = exp->GetGeom()->GetFid(j);
2126 
2127  pIt = perFaces.find(meshFaceId);
2128  if (pIt != perFaces.end())
2129  {
2130  if(meshFaceId == min(meshFaceId, pIt->second[0].id))
2131  {
2132  faceOrient = DeterminePeriodicFaceOrient(faceOrient,pIt->second[0].orient);
2133  }
2134  meshFaceId = min(meshFaceId, pIt->second[0].id);
2135  }
2136 
2137 
2138  exp->GetFaceInteriorMap(j,faceOrient,faceInteriorMap,faceInteriorSign);
2139  dof = exp->GetFaceIntNcoeffs(j);
2140 
2141 
2142  for(k = dof-1; k >= 0; --k)
2143  {
2144  vGlobalId = m_localToGlobalMap[cnt+faceInteriorMap[k]];
2145  m_globalToUniversalMap[vGlobalId]
2146  = nVert + nEdge*maxEdgeDof + meshFaceId * maxFaceDof
2147  + k + 1;
2148  m_globalToUniversalBndMap[vGlobalId]=m_globalToUniversalMap[vGlobalId];
2149 
2150  maxBndGlobalId = (vGlobalId > maxBndGlobalId ? vGlobalId : maxBndGlobalId);
2151  }
2152  }
2153 
2154  // Add interior DOFs to complete universal numbering
2155  exp->GetInteriorMap(interiorMap);
2156  dof = interiorMap.num_elements();
2157  elementId = (exp->GetGeom())->GetGlobalID();
2158  for (k = 0; k < dof; ++k)
2159  {
2160  vGlobalId = m_localToGlobalMap[cnt+interiorMap[k]];
2161  m_globalToUniversalMap[vGlobalId]
2162  = nVert + nEdge*maxEdgeDof + nFace*maxFaceDof + elementId*maxIntDof + k + 1;
2163  }
2164  }
2165 
2166  // Set up the GSLib universal assemble mapping
2167  // Internal DOF do not participate in any data
2168  // exchange, so we keep these set to the special GSLib id=0 so
2169  // they are ignored.
2170  Nektar::Array<OneD, long> tmp(m_numGlobalCoeffs);
2171  Vmath::Zero(m_numGlobalCoeffs, tmp, 1);
2172  Nektar::Array<OneD, long> tmp2(m_numGlobalBndCoeffs, tmp);
2173  for (unsigned int i = 0; i < m_numGlobalBndCoeffs; ++i)
2174  {
2175  tmp[i] = m_globalToUniversalMap[i];
2176  }
2177 
2178  m_gsh = Gs::Init(tmp, vCommRow);
2179  m_bndGsh = Gs::Init(tmp2, vCommRow);
2180  Gs::Unique(tmp, vCommRow);
2181  for (unsigned int i = 0; i < m_numGlobalCoeffs; ++i)
2182  {
2183  m_globalToUniversalMapUnique[i] = (tmp[i] >= 0 ? 1 : 0);
2184  }
2185  for (unsigned int i = 0; i < m_numGlobalBndCoeffs; ++i)
2186  {
2187  m_globalToUniversalBndMapUnique[i] = (tmp2[i] >= 0 ? 1 : 0);
2188  }
2189  }
Nektar::MultiRegions::AssemblyMap::m_numGlobalBndCoeffs
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:314
Nektar::MultiRegions::AssemblyMap::m_comm
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:306
Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMapUnique
Array< OneD, int > m_globalToUniversalMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMapCG.h:123
Gs::Init
static gs_data * Init(const Nektar::Array< OneD, long > pId, const LibUtilities::CommSharedPtr &pComm)
Initialise Gather-Scatter map.
Definition: GsLib.hpp:151
Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMap
Array< OneD, int > m_globalToUniversalMap
Integer map of process coeffs to universal space.
Definition: AssemblyMapCG.h:121
Nektar::LocalRegions::ExpansionVector
std::vector< ExpansionSharedPtr > ExpansionVector
Definition: Expansion.h:70
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
Nektar::LibUtilities::CommSharedPtr
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
Nektar::MultiRegions::DeterminePeriodicFaceOrient
StdRegions::Orientation DeterminePeriodicFaceOrient(StdRegions::Orientation faceOrient, StdRegions::Orientation perFaceOrient)
Determine relative orientation between two faces.
Definition: AssemblyMapCG.cpp:1950
Gs::Unique
static void Unique(const Nektar::Array< OneD, long > pId, const LibUtilities::CommSharedPtr &pComm)
Updates pId to negate all-but-one references to each universal ID.
Definition: GsLib.hpp:184
Nektar::LocalRegions::ExpansionSharedPtr
boost::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
Nektar::MultiRegions::AssemblyMap::m_globalToUniversalBndMap
Array< OneD, int > m_globalToUniversalBndMap
Integer map of process coeffs to universal space.
Definition: AssemblyMap.h:358
Nektar::MultiRegions::AssemblyMap::m_globalToUniversalBndMapUnique
Array< OneD, int > m_globalToUniversalBndMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMap.h:360
Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:342
Nektar::MultiRegions::DeterminePeriodicEdgeOrientId
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...
Definition: AssemblyMapCG.cpp:1902
Vmath::Zero
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:359
void Nektar::MultiRegions::AssemblyMapCG::v_Assemble ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2493 of file AssemblyMapCG.cpp.

References Vmath::Assmb(), m_localToGlobalMap, m_localToGlobalSign, Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs, Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs, Nektar::MultiRegions::AssemblyMap::m_signChange, Nektar::MultiRegions::AssemblyMap::UniversalAssemble(), and Vmath::Zero().

2496  {
2497  Array<OneD, const NekDouble> local;
2498  if(global.data() == loc.data())
2499  {
2500  local = Array<OneD, NekDouble>(loc.num_elements(),loc.data());
2501  }
2502  else
2503  {
2504  local = loc; // create reference
2505  }
2506  //ASSERTL1(loc.get() != global.get(),"Local and Global Arrays cannot be the same");
2507 
2508  Vmath::Zero(m_numGlobalCoeffs, global.get(), 1);
2509 
2510  if(m_signChange)
2511  {
2512  Vmath::Assmb(m_numLocalCoeffs, m_localToGlobalSign.get(), local.get(), m_localToGlobalMap.get(), global.get());
2513  }
2514  else
2515  {
2516  Vmath::Assmb(m_numLocalCoeffs, local.get(), m_localToGlobalMap.get(), global.get());
2517  }
2518  UniversalAssemble(global);
2519  }
Nektar::MultiRegions::AssemblyMap::m_signChange
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:345
Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:331
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
Nektar::MultiRegions::AssemblyMap::UniversalAssemble
void UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
Definition: AssemblyMap.cpp:764
Vmath::Assmb
void Assmb(int n, const T *x, const int *y, T *z)
Assemble z[y[i]] += x[i]; z should be zero'd first.
Definition: Vmath.cpp:695
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign
Array< OneD, NekDouble > m_localToGlobalSign
Integer sign of local coeffs to global space.
Definition: AssemblyMapCG.h:117
Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:342
Vmath::Zero
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:359
void Nektar::MultiRegions::AssemblyMapCG::v_Assemble ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2521 of file AssemblyMapCG.cpp.

References Nektar::MultiRegions::AssemblyMap::Assemble(), and Nektar::NekVector< DataType >::GetPtr().

2524  {
2525  Assemble(loc.GetPtr(),global.GetPtr());
2526  }
Nektar::MultiRegions::AssemblyMap::Assemble
void Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
Definition: AssemblyMap.cpp:750
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMapCG::v_GetExtraDirEdges ( )
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2590 of file AssemblyMapCG.cpp.

References m_extraDirEdges.

2591  {
2592  return m_extraDirEdges;
2593  }
Nektar::MultiRegions::AssemblyMapCG::m_extraDirEdges
Array< OneD, int > m_extraDirEdges
Extra dirichlet edges in parallel.
Definition: AssemblyMapCG.h:141
int Nektar::MultiRegions::AssemblyMapCG::v_GetFullSystemBandWidth ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2550 of file AssemblyMapCG.cpp.

References m_fullSystemBandWidth.

2551  {
2552  return m_fullSystemBandWidth;
2553  }
Nektar::MultiRegions::AssemblyMapCG::m_fullSystemBandWidth
int m_fullSystemBandWidth
Bandwith of the full matrix system (no static condensation).
Definition: AssemblyMapCG.h:119
int Nektar::MultiRegions::AssemblyMapCG::v_GetGlobalToUniversalMap ( const int  i) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2380 of file AssemblyMapCG.cpp.

References m_globalToUniversalMap.

2381  {
2382  return m_globalToUniversalMap[i];
2383  }
Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMap
Array< OneD, int > m_globalToUniversalMap
Integer map of process coeffs to universal space.
Definition: AssemblyMapCG.h:121
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMapCG::v_GetGlobalToUniversalMap ( void  )
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2397 of file AssemblyMapCG.cpp.

References m_globalToUniversalMap.

2398  {
2399  return m_globalToUniversalMap;
2400  }
Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMap
Array< OneD, int > m_globalToUniversalMap
Integer map of process coeffs to universal space.
Definition: AssemblyMapCG.h:121
int Nektar::MultiRegions::AssemblyMapCG::v_GetGlobalToUniversalMapUnique ( const int  i) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2385 of file AssemblyMapCG.cpp.

References m_globalToUniversalMapUnique.

2386  {
2387  return m_globalToUniversalMapUnique[i];
2388  }
Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMapUnique
Array< OneD, int > m_globalToUniversalMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMapCG.h:123
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMapCG::v_GetGlobalToUniversalMapUnique ( void  )
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2403 of file AssemblyMapCG.cpp.

References m_globalToUniversalMapUnique.

2404  {
2405  return m_globalToUniversalMapUnique;
2406  }
Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMapUnique
Array< OneD, int > m_globalToUniversalMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMapCG.h:123
int Nektar::MultiRegions::AssemblyMapCG::v_GetLocalToGlobalMap ( const int  i) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2375 of file AssemblyMapCG.cpp.

References m_localToGlobalMap.

2376  {
2377  return m_localToGlobalMap[i];
2378  }
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMapCG::v_GetLocalToGlobalMap ( void  )
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2391 of file AssemblyMapCG.cpp.

References m_localToGlobalMap.

2392  {
2393  return m_localToGlobalMap;
2394  }
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
NekDouble Nektar::MultiRegions::AssemblyMapCG::v_GetLocalToGlobalSign ( const int  i) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2408 of file AssemblyMapCG.cpp.

References m_localToGlobalSign, and Nektar::MultiRegions::AssemblyMap::m_signChange.

2410  {
2411  if(m_signChange)
2412  {
2413  return m_localToGlobalSign[i];
2414  }
2415  else
2416  {
2417  return 1.0;
2418  }
2419  }
Nektar::MultiRegions::AssemblyMap::m_signChange
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:345
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign
Array< OneD, NekDouble > m_localToGlobalSign
Integer sign of local coeffs to global space.
Definition: AssemblyMapCG.h:117
const Array< OneD, NekDouble > & Nektar::MultiRegions::AssemblyMapCG::v_GetLocalToGlobalSign ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2421 of file AssemblyMapCG.cpp.

References m_localToGlobalSign.

2422  {
2423  return m_localToGlobalSign;
2424  }
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign
Array< OneD, NekDouble > m_localToGlobalSign
Integer sign of local coeffs to global space.
Definition: AssemblyMapCG.h:117
int Nektar::MultiRegions::AssemblyMapCG::v_GetNumDirEdges ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2570 of file AssemblyMapCG.cpp.

References m_numDirEdges.

2571  {
2572  return m_numDirEdges;
2573  }
Nektar::MultiRegions::AssemblyMapCG::m_numDirEdges
int m_numDirEdges
Number of Dirichlet edges.
Definition: AssemblyMapCG.h:131
int Nektar::MultiRegions::AssemblyMapCG::v_GetNumDirFaces ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2575 of file AssemblyMapCG.cpp.

References m_numDirFaces.

2576  {
2577  return m_numDirFaces;
2578  }
Nektar::MultiRegions::AssemblyMapCG::m_numDirFaces
int m_numDirFaces
Number of Dirichlet faces.
Definition: AssemblyMapCG.h:133
int Nektar::MultiRegions::AssemblyMapCG::v_GetNumNonDirEdgeModes ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2560 of file AssemblyMapCG.cpp.

References m_numNonDirEdgeModes.

2561  {
2562  return m_numNonDirEdgeModes;
2563  }
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirEdgeModes
int m_numNonDirEdgeModes
Number of non Dirichlet edge modes.
Definition: AssemblyMapCG.h:127
int Nektar::MultiRegions::AssemblyMapCG::v_GetNumNonDirEdges ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2580 of file AssemblyMapCG.cpp.

References m_numNonDirEdges.

2581  {
2582  return m_numNonDirEdges;
2583  }
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirEdges
int m_numNonDirEdges
Number of Dirichlet edges.
Definition: AssemblyMapCG.h:135
int Nektar::MultiRegions::AssemblyMapCG::v_GetNumNonDirFaceModes ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2565 of file AssemblyMapCG.cpp.

References m_numNonDirFaceModes.

2566  {
2567  return m_numNonDirFaceModes;
2568  }
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirFaceModes
int m_numNonDirFaceModes
Number of non Dirichlet face modes.
Definition: AssemblyMapCG.h:129
int Nektar::MultiRegions::AssemblyMapCG::v_GetNumNonDirFaces ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2585 of file AssemblyMapCG.cpp.

References m_numNonDirFaces.

2586  {
2587  return m_numNonDirFaces;
2588  }
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirFaces
int m_numNonDirFaces
Number of Dirichlet faces.
Definition: AssemblyMapCG.h:137
int Nektar::MultiRegions::AssemblyMapCG::v_GetNumNonDirVertexModes ( ) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2555 of file AssemblyMapCG.cpp.

References m_numNonDirVertexModes.

2556  {
2557  return m_numNonDirVertexModes;
2558  }
Nektar::MultiRegions::AssemblyMapCG::m_numNonDirVertexModes
int m_numNonDirVertexModes
Number of non Dirichlet vertex modes.
Definition: AssemblyMapCG.h:125
void Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc 
) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2461 of file AssemblyMapCG.cpp.

References Vmath::Gathr(), m_localToGlobalMap, m_localToGlobalSign, Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs, and Nektar::MultiRegions::AssemblyMap::m_signChange.

2464  {
2465  Array<OneD, const NekDouble> glo;
2466  if(global.data() == loc.data())
2467  {
2468  glo = Array<OneD, NekDouble>(global.num_elements(),global.data());
2469  }
2470  else
2471  {
2472  glo = global; // create reference
2473  }
2474 
2475 
2476  if(m_signChange)
2477  {
2478  Vmath::Gathr(m_numLocalCoeffs, m_localToGlobalSign.get(), glo.get(), m_localToGlobalMap.get(), loc.get());
2479  }
2480  else
2481  {
2482  Vmath::Gathr(m_numLocalCoeffs, glo.get(), m_localToGlobalMap.get(), loc.get());
2483  }
2484  }
Nektar::MultiRegions::AssemblyMap::m_signChange
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:345
Vmath::Gathr
void Gathr(int n, const T *x, const int *y, T *z)
Gather vector z[i] = x[y[i]].
Definition: Vmath.cpp:630
Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:331
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign
Array< OneD, NekDouble > m_localToGlobalSign
Integer sign of local coeffs to global space.
Definition: AssemblyMapCG.h:117
void Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2486 of file AssemblyMapCG.cpp.

References Nektar::NekVector< DataType >::GetPtr(), and Nektar::MultiRegions::AssemblyMap::GlobalToLocal().

2489  {
2490  GlobalToLocal(global.GetPtr(),loc.GetPtr());
2491  }
Nektar::MultiRegions::AssemblyMap::GlobalToLocal
void GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
Definition: AssemblyMap.cpp:736
AssemblyMapSharedPtr Nektar::MultiRegions::AssemblyMapCG::v_LinearSpaceMap ( const ExpList &  locexp,
GlobalSysSolnType  solnType 
)
protectedvirtual

Construct an AssemblyMapCG object which corresponds to the linear space of the current object.

This function is used to create a linear-space assembly map, which is then used in the linear space preconditioner in the conjugate gradient solve.

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2199 of file AssemblyMapCG.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::MultiRegions::eNull, Nektar::MultiRegions::ExpList::GetExp(), Gs::Init(), Nektar::MultiRegions::AssemblyMap::m_comm, m_globalToUniversalMap, m_localToGlobalMap, Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs, Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs, Nektar::MultiRegions::AssemblyMap::m_session, Gs::Unique(), and Vmath::Zero().

2201  {
2202  AssemblyMapCGSharedPtr returnval;
2203 
2204  int i, j;
2205  int nverts = 0;
2206  const boost::shared_ptr<LocalRegions::ExpansionVector> exp
2207  = locexp.GetExp();
2208  int nelmts = exp->size();
2209 
2210  // Get Default Map and turn off any searched values.
2211  returnval = MemoryManager<AssemblyMapCG>
2212  ::AllocateSharedPtr(m_session);
2213  returnval->m_solnType = solnType;
2214  returnval->m_preconType = eNull;
2215  returnval->m_maxStaticCondLevel = 0;
2216  returnval->m_signChange = false;
2217  returnval->m_comm = m_comm;
2218 
2219  // Count the number of vertices
2220  for (i = 0; i < nelmts; ++i)
2221  {
2222  nverts += (*exp)[i]->GetNverts();
2223  }
2224 
2225  returnval->m_numLocalCoeffs = nverts;
2226  returnval->m_localToGlobalMap = Array<OneD, int>(nverts, -1);
2227 
2228  // Store original global ids in this map
2229  returnval->m_localToGlobalBndMap = Array<OneD, int>(nverts, -1);
2230 
2231  int cnt = 0;
2232  int cnt1 = 0;
2233  Array<OneD, int> GlobCoeffs(m_numGlobalCoeffs, -1);
2234 
2235  // Set up local to global map;
2236  for (i = 0; i < nelmts; ++i)
2237  {
2238  for (j = 0; j < (*exp)[i]->GetNverts(); ++j)
2239  {
2240  returnval->m_localToGlobalMap[cnt] =
2241  returnval->m_localToGlobalBndMap[cnt] =
2242  m_localToGlobalMap[cnt1 + (*exp)[i]->GetVertexMap(j,true)];
2243  GlobCoeffs[returnval->m_localToGlobalMap[cnt]] = 1;
2244 
2245  // Set up numLocalDirBndCoeffs
2246  if ((returnval->m_localToGlobalMap[cnt]) <
2247  m_numGlobalDirBndCoeffs)
2248  {
2249  returnval->m_numLocalDirBndCoeffs++;
2250  }
2251  cnt++;
2252  }
2253  cnt1 += (*exp)[i]->GetNcoeffs();
2254  }
2255 
2256  cnt = 0;
2257  // Reset global numbering and count number of dofs
2258  for (i = 0; i < m_numGlobalCoeffs; ++i)
2259  {
2260  if (GlobCoeffs[i] != -1)
2261  {
2262  GlobCoeffs[i] = cnt++;
2263  }
2264  }
2265 
2266  // Set up number of globalCoeffs;
2267  returnval->m_numGlobalCoeffs = cnt;
2268 
2269  // Set up number of global Dirichlet boundary coefficients
2270  for (i = 0; i < m_numGlobalDirBndCoeffs; ++i)
2271  {
2272  if (GlobCoeffs[i] != -1)
2273  {
2274  returnval->m_numGlobalDirBndCoeffs++;
2275  }
2276  }
2277 
2278  // Set up global to universal map
2279  if (m_globalToUniversalMap.num_elements())
2280  {
2281  LibUtilities::CommSharedPtr vCommRow
2282  = m_session->GetComm()->GetRowComm();
2283  int nglocoeffs = returnval->m_numGlobalCoeffs;
2284  returnval->m_globalToUniversalMap
2285  = Array<OneD, int> (nglocoeffs);
2286  returnval->m_globalToUniversalMapUnique
2287  = Array<OneD, int> (nglocoeffs);
2288 
2289  // Reset local to global map and setup universal map
2290  for (i = 0; i < nverts; ++i)
2291  {
2292  cnt = returnval->m_localToGlobalMap[i];
2293  returnval->m_localToGlobalMap[i] = GlobCoeffs[cnt];
2294 
2295  returnval->m_globalToUniversalMap[GlobCoeffs[cnt]] =
2296  m_globalToUniversalMap[cnt];
2297  }
2298 
2299  Nektar::Array<OneD, long> tmp(nglocoeffs);
2300  Vmath::Zero(nglocoeffs, tmp, 1);
2301  for (unsigned int i = 0; i < nglocoeffs; ++i)
2302  {
2303  tmp[i] = returnval->m_globalToUniversalMap[i];
2304  }
2305  returnval->m_gsh = Gs::Init(tmp, vCommRow);
2306  Gs::Unique(tmp, vCommRow);
2307  for (unsigned int i = 0; i < nglocoeffs; ++i)
2308  {
2309  returnval->m_globalToUniversalMapUnique[i]
2310  = (tmp[i] >= 0 ? 1 : 0);
2311  }
2312  }
2313  else // not sure this option is ever needed.
2314  {
2315  for (i = 0; i < nverts; ++i)
2316  {
2317  cnt = returnval->m_localToGlobalMap[i];
2318  returnval->m_localToGlobalMap[i] = GlobCoeffs[cnt];
2319  }
2320  }
2321 
2322  return returnval;
2323  }
Nektar::MultiRegions::AssemblyMap::m_comm
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:306
Nektar::MemoryManager::AllocateSharedPtr
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
Definition: NekMemoryManager.hpp:235
Gs::Init
static gs_data * Init(const Nektar::Array< OneD, long > pId, const LibUtilities::CommSharedPtr &pComm)
Initialise Gather-Scatter map.
Definition: GsLib.hpp:151
Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMap
Array< OneD, int > m_globalToUniversalMap
Integer map of process coeffs to universal space.
Definition: AssemblyMapCG.h:121
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
Nektar::LibUtilities::CommSharedPtr
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
Gs::Unique
static void Unique(const Nektar::Array< OneD, long > pId, const LibUtilities::CommSharedPtr &pComm)
Updates pId to negate all-but-one references to each universal ID.
Definition: GsLib.hpp:184
Nektar::MultiRegions::eNull
No Solution type specified.
Definition: MultiRegions.hpp:123
Nektar::MultiRegions::AssemblyMap::m_session
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:303
Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:342
Vmath::Zero
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:359
Nektar::MultiRegions::AssemblyMapCGSharedPtr
boost::shared_ptr< AssemblyMapCG > AssemblyMapCGSharedPtr
Definition: AssemblyMapCG.h:52
void Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2426 of file AssemblyMapCG.cpp.

References Gs::Gather(), Gs::gs_max, Nektar::MultiRegions::AssemblyMap::m_gsh, m_localToGlobalMap, m_localToGlobalSign, Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs, Nektar::MultiRegions::AssemblyMap::m_signChange, and Vmath::Scatr().

2429  {
2430  Array<OneD, const NekDouble> local;
2431  if(global.data() == loc.data())
2432  {
2433  local = Array<OneD, NekDouble>(loc.num_elements(),loc.data());
2434  }
2435  else
2436  {
2437  local = loc; // create reference
2438  }
2439 
2440 
2441  if(m_signChange)
2442  {
2443  Vmath::Scatr(m_numLocalCoeffs, m_localToGlobalSign.get(), local.get(), m_localToGlobalMap.get(), global.get());
2444  }
2445  else
2446  {
2447  Vmath::Scatr(m_numLocalCoeffs, local.get(), m_localToGlobalMap.get(), global.get());
2448  }
2449 
2450  // ensure all values are unique by calling a max
2451  Gs::Gather(global, Gs::gs_max, m_gsh);
2452  }
Nektar::MultiRegions::AssemblyMap::m_signChange
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:345
Gs::Gather
static void Gather(Nektar::Array< OneD, NekDouble > pU, gs_op pOp, gs_data *pGsh, Nektar::Array< OneD, NekDouble > pBuffer=NullNekDouble1DArray)
Performs a gather-scatter operation of the provided values.
Definition: GsLib.hpp:224
Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:331
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
Array< OneD, int > m_localToGlobalMap
Integer map of local coeffs to global space.
Definition: AssemblyMapCG.h:115
Vmath::Scatr
void Scatr(int n, const T *x, const int *y, T *z)
Scatter vector z[y[i]] = x[i].
Definition: Vmath.cpp:659
Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign
Array< OneD, NekDouble > m_localToGlobalSign
Integer sign of local coeffs to global space.
Definition: AssemblyMapCG.h:117
void Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2454 of file AssemblyMapCG.cpp.

References Nektar::NekVector< DataType >::GetPtr(), and Nektar::MultiRegions::AssemblyMap::LocalToGlobal().

2457  {
2458  LocalToGlobal(loc.GetPtr(),global.GetPtr());
2459  }
Nektar::MultiRegions::AssemblyMap::LocalToGlobal
void LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
Definition: AssemblyMap.cpp:722
void Nektar::MultiRegions::AssemblyMapCG::v_UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2528 of file AssemblyMapCG.cpp.

References Gs::Gather(), Gs::gs_add, and Nektar::MultiRegions::AssemblyMap::m_gsh.

2530  {
2531  Gs::Gather(pGlobal, Gs::gs_add, m_gsh);
2532  }
Gs::Gather
static void Gather(Nektar::Array< OneD, NekDouble > pU, gs_op pOp, gs_data *pGsh, Nektar::Array< OneD, NekDouble > pBuffer=NullNekDouble1DArray)
Performs a gather-scatter operation of the provided values.
Definition: GsLib.hpp:224
void Nektar::MultiRegions::AssemblyMapCG::v_UniversalAssemble ( NekVector< NekDouble > &  pGlobal) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2534 of file AssemblyMapCG.cpp.

References Nektar::NekVector< DataType >::GetPtr(), and Nektar::MultiRegions::AssemblyMap::UniversalAssemble().

2536  {
2537  UniversalAssemble(pGlobal.GetPtr());
2538  }
Nektar::MultiRegions::AssemblyMap::UniversalAssemble
void UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
Definition: AssemblyMap.cpp:764
void Nektar::MultiRegions::AssemblyMapCG::v_UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal,
int  offset 
) const
protectedvirtual

Reimplemented from Nektar::MultiRegions::AssemblyMap.

Definition at line 2540 of file AssemblyMapCG.cpp.

References Nektar::MultiRegions::AssemblyMap::UniversalAssemble(), and Vmath::Vcopy().

2543  {
2544  Array<OneD, NekDouble> tmp(offset);
2545  Vmath::Vcopy(offset, pGlobal, 1, tmp, 1);
2546  UniversalAssemble(pGlobal);
2547  Vmath::Vcopy(offset, tmp, 1, pGlobal, 1);
2548  }
Nektar::MultiRegions::AssemblyMap::UniversalAssemble
void UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
Definition: AssemblyMap.cpp:764
Vmath::Vcopy
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047

Member Data Documentation

std::map<int, std::vector<ExtraDirDof> > Nektar::MultiRegions::AssemblyMapCG::m_extraDirDofs
protected

Map indicating degrees of freedom which are Dirichlet but whose value is stored on another processor.

Definition at line 148 of file AssemblyMapCG.h.

Referenced by AssemblyMapCG(), Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), and GetExtraDirDofs().

Array<OneD, int> Nektar::MultiRegions::AssemblyMapCG::m_extraDirEdges
protected

Extra dirichlet edges in parallel.

Definition at line 141 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetExtraDirEdges().

int Nektar::MultiRegions::AssemblyMapCG::m_fullSystemBandWidth
protected

Bandwith of the full matrix system (no static condensation).

Definition at line 119 of file AssemblyMapCG.h.

Referenced by CalculateFullSystemBandWidth(), and v_GetFullSystemBandWidth().

Array<OneD,int> Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMap
protected

Integer map of process coeffs to universal space.

Definition at line 121 of file AssemblyMapCG.h.

Referenced by Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), SetUpUniversalC0ContMap(), v_GetGlobalToUniversalMap(), and v_LinearSpaceMap().

Array<OneD,int> Nektar::MultiRegions::AssemblyMapCG::m_globalToUniversalMapUnique
protected

Integer map of unique process coeffs to universal space (signed)

Definition at line 123 of file AssemblyMapCG.h.

Referenced by Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), SetUpUniversalC0ContMap(), and v_GetGlobalToUniversalMapUnique().

Array<OneD,int> Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalMap
protected

Integer map of local coeffs to global space.

Definition at line 115 of file AssemblyMapCG.h.

Referenced by AssemblyMapCG(), CalculateFullSystemBandWidth(), Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), SetUpUniversalC0ContMap(), v_Assemble(), v_GetLocalToGlobalMap(), v_GlobalToLocal(), v_LinearSpaceMap(), and v_LocalToGlobal().

Array<OneD,NekDouble> Nektar::MultiRegions::AssemblyMapCG::m_localToGlobalSign
protected

Integer sign of local coeffs to global space.

Definition at line 117 of file AssemblyMapCG.h.

Referenced by AssemblyMapCG(), Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), v_Assemble(), v_GetLocalToGlobalSign(), v_GlobalToLocal(), and v_LocalToGlobal().

int Nektar::MultiRegions::AssemblyMapCG::m_maxStaticCondLevel
protected

Maximum static condensation level.

Definition at line 145 of file AssemblyMapCG.h.

Referenced by AssemblyMapCG().

int Nektar::MultiRegions::AssemblyMapCG::m_numDirEdges
protected

Number of Dirichlet edges.

Definition at line 131 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetNumDirEdges().

int Nektar::MultiRegions::AssemblyMapCG::m_numDirFaces
protected

Number of Dirichlet faces.

Definition at line 133 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetNumDirFaces().

int Nektar::MultiRegions::AssemblyMapCG::m_numLocalBndCondCoeffs
protected

Number of local boundary condition coefficients.

Definition at line 139 of file AssemblyMapCG.h.

Referenced by AssemblyMapCG(), and CreateGraph().

int Nektar::MultiRegions::AssemblyMapCG::m_numLocDirBndCondDofs
protected

Number of local boundary condition degrees of freedom.

Definition at line 143 of file AssemblyMapCG.h.

int Nektar::MultiRegions::AssemblyMapCG::m_numNonDirEdgeModes
protected

Number of non Dirichlet edge modes.

Definition at line 127 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetNumNonDirEdgeModes().

int Nektar::MultiRegions::AssemblyMapCG::m_numNonDirEdges
protected

Number of Dirichlet edges.

Definition at line 135 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetNumNonDirEdges().

int Nektar::MultiRegions::AssemblyMapCG::m_numNonDirFaceModes
protected

Number of non Dirichlet face modes.

Definition at line 129 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetNumNonDirFaceModes().

int Nektar::MultiRegions::AssemblyMapCG::m_numNonDirFaces
protected

Number of Dirichlet faces.

Definition at line 137 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetNumNonDirFaces().

int Nektar::MultiRegions::AssemblyMapCG::m_numNonDirVertexModes
protected

Number of non Dirichlet vertex modes.

Definition at line 125 of file AssemblyMapCG.h.

Referenced by CreateGraph(), and v_GetNumNonDirVertexModes().

Generated on Wed May 4 2016 17:54:31 for Nektar++ by   doxygen 1.8.8