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

Base class for constructing local to global mapping of degrees of freedom. More...

#include <AssemblyMap.h>

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

Public Member Functions

 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, bool useComm=true) const
 
void LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, bool useComm=true) const
 
void GlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
void GlobalToLocal (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
 
void Assemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
void Assemble (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
void UniversalAssemble (Array< OneD, NekDouble > &pGlobal) const
 
void UniversalAssemble (NekVector< NekDouble > &pGlobal) const
 
void UniversalAssemble (Array< OneD, NekDouble > &pGlobal, int offset) const
 
int GetLocalToGlobalBndMap (const int i) const
 Retrieve the global index of a given local boundary mode. More...
 
const Array< OneD, const int > & GetLocalToGlobalBndMap ()
 Retrieve the global indices of the local boundary modes. More...
 
const Array< OneD, const int > & GetGlobalToUniversalBndMap ()
 
const Array< OneD, const int > & GetGlobalToUniversalBndMapUnique ()
 
bool GetSignChange ()
 Returns true if using a modal expansion requiring a change of sign of some modes. More...
 
NekDouble GetLocalToGlobalBndSign (const int i) const
 Retrieve the sign change of a given local boundary mode. More...
 
Array< OneD, const NekDoubleGetLocalToGlobalBndSign () const
 Retrieve the sign change for all local boundary modes. More...
 
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, bool printHeader=true) const
 
const Array< OneD, const int > & GetExtraDirEdges ()
 
boost::shared_ptr< AssemblyMapLinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType)
 
int GetBndSystemBandWidth () const
 Returns the bandwidth of the boundary system. More...
 
int GetStaticCondLevel () const
 Returns the level of static condensation for this map. More...
 
int GetNumPatches () const
 Returns the number of patches in this static condensation level. More...
 
const Array< OneD, const
unsigned int > & 
GetNumLocalBndCoeffsPerPatch ()
 Returns the number of local boundary coefficients in each patch. More...
 
const Array< OneD, const
unsigned int > & 
GetNumLocalIntCoeffsPerPatch ()
 Returns the number of local interior coefficients in each patch. More...
 
const AssemblyMapSharedPtr GetNextLevelLocalToGlobalMap () const
 Returns the local to global mapping for the next level in the multi-level static condensation. More...
 
void SetNextLevelLocalToGlobalMap (AssemblyMapSharedPtr pNextLevelLocalToGlobalMap)
 
const PatchMapSharedPtrGetPatchMapFromPrevLevel (void) const
 Returns the patch map from the previous level of the multi-level static condensation. More...
 
bool AtLastLevel () const
 Returns true if this is the last level in the multi-level static condensation. More...
 
GlobalSysSolnType GetGlobalSysSolnType () const
 Returns the method of solving global systems. More...
 
PreconditionerType GetPreconType () const
 
NekDouble GetIterativeTolerance () const
 
int GetMaxIterations () const
 
int GetSuccessiveRHS () const
 
int GetLowestStaticCondLevel () const
 

Protected Member Functions

void CalculateBndSystemBandWidth ()
 Calculates the bandwidth of the boundary system. More...
 
void GlobalToLocalBndWithoutSign (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc)
 

Protected Attributes

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, NekDoublem_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, NekDoublem_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_datam_gsh
 
Gs::gs_datam_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 Member Functions

virtual int v_GetLocalToGlobalMap (const int i) const
 
virtual int v_GetGlobalToUniversalMap (const int i) const
 
virtual int v_GetGlobalToUniversalMapUnique (const int i) const
 
virtual const Array< OneD,
const int > & 
v_GetLocalToGlobalMap ()
 
virtual const Array< OneD,
const int > & 
v_GetGlobalToUniversalMap ()
 
virtual const Array< OneD,
const int > & 
v_GetGlobalToUniversalMapUnique ()
 
virtual NekDouble v_GetLocalToGlobalSign (const int i) const
 
virtual const Array< OneD,
NekDouble > & 
v_GetLocalToGlobalSign () const
 
virtual void v_LocalToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm) const
 
virtual void v_LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, bool useComm) const
 
virtual void v_GlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
virtual void v_GlobalToLocal (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
 
virtual void v_Assemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
virtual void v_Assemble (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal) const
 
virtual void v_UniversalAssemble (NekVector< NekDouble > &pGlobal) const
 
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal, int offset) const
 
virtual int v_GetFullSystemBandWidth () const
 
virtual int v_GetNumNonDirVertexModes () const
 
virtual int v_GetNumNonDirEdgeModes () const
 
virtual int v_GetNumNonDirFaceModes () const
 
virtual int v_GetNumDirEdges () const
 
virtual int v_GetNumDirFaces () const
 
virtual int v_GetNumNonDirEdges () const
 
virtual int v_GetNumNonDirFaces () const
 
virtual const Array< OneD,
const int > & 
v_GetExtraDirEdges ()
 
virtual boost::shared_ptr
< AssemblyMap
v_LinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType)
 Generate a linear space mapping from existing mapping. More...
 

Private Attributes

PatchMapSharedPtr m_patchMapFromPrevLevel
 Mapping information for previous level in MultiLevel Solver. More...
 

Detailed Description

Base class for constructing local to global mapping of degrees of freedom.

This class acts as a base class for constructing mappings between local, global and boundary degrees of freedom. It holds the storage for the maps and provides the accessors needed to retrieve them.

There are two derived classes: AssemblyMapCG and AssemblyMapDG. These perform the actual construction of the maps within their specific contexts.

Definition at line 59 of file AssemblyMap.h.

Constructor & Destructor Documentation

Nektar::MultiRegions::AssemblyMap::AssemblyMap ( )

Default constructor.

Initialises an empty mapping.

Definition at line 79 of file AssemblyMap.cpp.

79  :
80  m_session(),
81  m_comm(),
82  m_hash(0),
89  m_successiveRHS(0),
90  m_gsh(0),
91  m_bndGsh(0)
92  {
93  }
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:308
int m_successiveRHS
sucessive RHS for iterative solver
Definition: AssemblyMap.h:379
size_t m_hash
Hash for map.
Definition: AssemblyMap.h:311
int m_bndSystemBandWidth
The bandwith of the global bnd system.
Definition: AssemblyMap.h:367
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:365
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:320
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:305
No Solution type specified.
Nektar::MultiRegions::AssemblyMap::AssemblyMap ( const LibUtilities::SessionReaderSharedPtr pSession,
const std::string  variable = "DefaultVar" 
)

Constructor with a communicator.

Definition at line 95 of file AssemblyMap.cpp.

References Nektar::NekConstants::kNekIterativeTol, m_iterativeTolerance, m_maxIterations, m_preconType, m_solnType, and m_successiveRHS.

97  :
98  m_session(pSession),
99  m_comm(pSession->GetComm()),
100  m_hash(0),
106  m_successiveRHS(0),
107  m_gsh(0),
108  m_bndGsh(0)
109  {
110  // Default value from Solver Info
111  m_solnType = pSession->GetSolverInfoAsEnum<GlobalSysSolnType>(
112  "GlobalSysSoln");
113  m_preconType = pSession->GetSolverInfoAsEnum<PreconditionerType>(
114  "Preconditioner");
115 
116  // Override values with data from GlobalSysSolnInfo section
117  if(pSession->DefinesGlobalSysSolnInfo(variable, "GlobalSysSoln"))
118  {
119  std::string sysSoln = pSession->GetGlobalSysSolnInfo(variable,
120  "GlobalSysSoln");
121  m_solnType = pSession->GetValueAsEnum<GlobalSysSolnType>(
122  "GlobalSysSoln", sysSoln);
123  }
124 
125  if(pSession->DefinesGlobalSysSolnInfo(variable, "Preconditioner"))
126  {
127  std::string precon = pSession->GetGlobalSysSolnInfo(variable,
128  "Preconditioner");
129  m_preconType = pSession->GetValueAsEnum<PreconditionerType>(
130  "Preconditioner", precon);
131  }
132 
133  if(pSession->DefinesGlobalSysSolnInfo(variable,
134  "IterativeSolverTolerance"))
135  {
136  m_iterativeTolerance = boost::lexical_cast<NekDouble>(
137  pSession->GetGlobalSysSolnInfo(variable,
138  "IterativeSolverTolerance").c_str());
139  }
140  else
141  {
142  pSession->LoadParameter("IterativeSolverTolerance",
145  }
146 
147 
148  if(pSession->DefinesGlobalSysSolnInfo(variable,
149  "MaxIterations"))
150  {
151  m_maxIterations = boost::lexical_cast<int>(
152  pSession->GetGlobalSysSolnInfo(variable,
153  "MaxIterations").c_str());
154  }
155  else
156  {
157  pSession->LoadParameter("MaxIterations",
159  5000);
160  }
161 
162 
163  if(pSession->DefinesGlobalSysSolnInfo(variable,"SuccessiveRHS"))
164  {
165  m_successiveRHS = boost::lexical_cast<int>(
166  pSession->GetGlobalSysSolnInfo(variable,
167  "SuccessiveRHS").c_str());
168  }
169  else
170  {
171  pSession->LoadParameter("SuccessiveRHS",
172  m_successiveRHS,0);
173  }
174 
175  }
PreconditionerType m_preconType
Type type of preconditioner to use in iterative solver.
Definition: AssemblyMap.h:370
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:308
static const NekDouble kNekIterativeTol
int m_successiveRHS
sucessive RHS for iterative solver
Definition: AssemblyMap.h:379
size_t m_hash
Hash for map.
Definition: AssemblyMap.h:311
int m_bndSystemBandWidth
The bandwith of the global bnd system.
Definition: AssemblyMap.h:367
NekDouble m_iterativeTolerance
Tolerance for iterative solver.
Definition: AssemblyMap.h:376
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:365
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:320
double NekDouble
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
int m_maxIterations
Maximum iterations for iterative solver.
Definition: AssemblyMap.h:373
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:305
Nektar::MultiRegions::AssemblyMap::AssemblyMap ( AssemblyMap oldLevelMap,
const BottomUpSubStructuredGraphSharedPtr multiLevelGraph 
)

Constructor for next level in multi-level static condensation.

Create a new level of mapping using the information in multiLevelGraph and performing the following steps:

  • STEP 1: setup a mask array to determine to which patch of the new level every patch of the current level belongs. To do so we make four arrays, #gloPatchMask, #globHomPatchMask, #locPatchMask_NekDouble and #locPatchMask. These arrays are then used to check which local dofs of the old level belong to which patch of the new level
  • STEP 2: We calculate how many local bnd dofs of the old level belong to the boundaries of each patch at the new level. We need this information to set up the mapping between different levels.

Definition at line 181 of file AssemblyMap.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, ASSERTL1, CalculateBndSystemBandWidth(), Nektar::MultiRegions::eDirectMultiLevelStaticCond, Nektar::MultiRegions::eIterativeMultiLevelStaticCond, Nektar::MultiRegions::ePETScMultiLevelStaticCond, Nektar::MultiRegions::eXxtMultiLevelStaticCond, GetGlobalSysSolnType(), GetGlobalToUniversalBndMap(), GetGlobalToUniversalBndMapUnique(), GetLocalToGlobalBndMap(), GetLocalToGlobalBndSign(), GetNumGlobalBndCoeffs(), GetNumGlobalDirBndCoeffs(), GetNumLocalBndCoeffs(), GetNumLocalBndCoeffsPerPatch(), GetNumLocalDirBndCoeffs(), GetNumPatches(), GetSignChange(), GetStaticCondLevel(), GlobalToLocalBndWithoutSign(), m_globalToUniversalBndMap, m_globalToUniversalBndMapUnique, m_localToGlobalBndMap, m_localToGlobalBndSign, m_nextLevelLocalToGlobalMap, m_numGlobalBndCoeffs, m_numGlobalCoeffs, m_numGlobalDirBndCoeffs, m_numLocalBndCoeffs, m_numLocalBndCoeffsPerPatch, m_numLocalDirBndCoeffs, m_numLocalIntCoeffsPerPatch, m_numPatches, m_patchMapFromPrevLevel, m_signChange, m_solnType, m_staticCondLevel, Nektar::MultiRegions::RoundNekDoubleToInt(), and sign.

183  :
184  m_session(oldLevelMap->m_session),
185  m_comm(oldLevelMap->GetComm()),
186  m_hash(0),
187  m_solnType(oldLevelMap->m_solnType),
188  m_preconType(oldLevelMap->m_preconType),
189  m_maxIterations(oldLevelMap->m_maxIterations),
190  m_iterativeTolerance(oldLevelMap->m_iterativeTolerance),
191  m_successiveRHS(oldLevelMap->m_successiveRHS),
192  m_gsh(oldLevelMap->m_gsh),
193  m_bndGsh(oldLevelMap->m_bndGsh),
194  m_lowestStaticCondLevel(oldLevelMap->m_lowestStaticCondLevel)
195  {
196  int i;
197  int j;
198  int cnt;
199 
200  //--------------------------------------------------------------
201  // -- Extract information from the input argument
202  int numGlobalBndCoeffsOld = oldLevelMap->GetNumGlobalBndCoeffs();
203  int numGlobalDirBndCoeffsOld = oldLevelMap->GetNumGlobalDirBndCoeffs();
204  int numLocalBndCoeffsOld = oldLevelMap->GetNumLocalBndCoeffs();
205  int numLocalDirBndCoeffsOld = oldLevelMap->GetNumLocalDirBndCoeffs();
206  bool signChangeOld = oldLevelMap->GetSignChange();
207 
208  int staticCondLevelOld = oldLevelMap->GetStaticCondLevel();
209  int numPatchesOld = oldLevelMap->GetNumPatches();
210  GlobalSysSolnType solnTypeOld = oldLevelMap->GetGlobalSysSolnType();
211  const Array<OneD, const unsigned int>& numLocalBndCoeffsPerPatchOld = oldLevelMap->GetNumLocalBndCoeffsPerPatch();
212  //--------------------------------------------------------------
213 
214  //--------------------------------------------------------------
215  int newLevel = staticCondLevelOld+1;
216  /** - STEP 1: setup a mask array to determine to which patch
217  * of the new level every patch of the current
218  * level belongs. To do so we make four arrays,
219  * #gloPatchMask, #globHomPatchMask,
220  * #locPatchMask_NekDouble and #locPatchMask.
221  * These arrays are then used to check which local
222  * dofs of the old level belong to which patch of
223  * the new level
224  */
225  Array<OneD, NekDouble> globPatchMask (numGlobalBndCoeffsOld,-1.0);
226  Array<OneD, NekDouble> globHomPatchMask (globPatchMask+numGlobalDirBndCoeffsOld);
227  Array<OneD, NekDouble> locPatchMask_NekDouble(numLocalBndCoeffsOld,-3.0);
228  Array<OneD, int> locPatchMask (numLocalBndCoeffsOld);
229 
230  // Fill the array globPatchMask as follows:
231  // - The first part (i.e. the glob bnd dofs) is filled with the
232  // value -1
233  // - The second part (i.e. the glob interior dofs) is numbered
234  // according to the patch it belongs to (i.e. dofs in first block
235  // all are numbered 0, the second block numbered are 1, etc...)
236  multiLevelGraph->MaskPatches(newLevel,globHomPatchMask);
237 
238  // Map from Global Dofs to Local Dofs
239  // As a result, we know for each local dof whether
240  // it is mapped to the boundary of the next level, or to which
241  // patch. Based upon this, we can than later associate every patch
242  // of the current level with a patch in the next level.
243  oldLevelMap->GlobalToLocalBndWithoutSign(globPatchMask,locPatchMask_NekDouble);
244 
245  // Convert the result to an array of integers rather than doubles
246  RoundNekDoubleToInt(locPatchMask_NekDouble,locPatchMask);
247 
248  /** - STEP 2: We calculate how many local bnd dofs of the
249  * old level belong to the boundaries of each patch at
250  * the new level. We need this information to set up the
251  * mapping between different levels.
252  */
253 
254  // Retrieve the number of patches at the next level
255  int numPatchesWithIntNew = multiLevelGraph->GetNpatchesWithInterior(newLevel);
256  int numPatchesNew = numPatchesWithIntNew;
257 
258  // Allocate memory to store the number of local dofs associated to
259  // each of elemental boundaries of these patches
260  std::map<int, int> numLocalBndCoeffsPerPatchNew;
261  for(int i = 0; i < numPatchesNew; i++)
262  {
263  numLocalBndCoeffsPerPatchNew[i] = 0;
264  }
265 
266  int minval;
267  int maxval;
268  int curPatch;
269  for(i = cnt = 0; i < numPatchesOld; i++)
270  {
271  // For every patch at the current level, the mask array
272  // locPatchMask should be filled with
273  // - the same (positive) number for each entry (which will
274  // correspond to the patch at the next level it belongs to)
275  // - the same (positive) number for each entry, except some
276  // entries that are -1 (the enties correspond to -1, will be
277  // mapped to the local boundary of the next level patch given
278  // by the positive number)
279  // - -1 for all entries. In this case, we will make an
280  // additional patch only consisting of boundaries at the next
281  // level
282  minval = *min_element(&locPatchMask[cnt],
283  &locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
284  maxval = *max_element(&locPatchMask[cnt],
285  &locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
286  ASSERTL0((minval==maxval)||(minval==-1),"These values should never be the same");
287 
288  if(maxval == -1)
289  {
290  curPatch = numPatchesNew;
291  numLocalBndCoeffsPerPatchNew[curPatch] = 0;
292  numPatchesNew++;
293  }
294  else
295  {
296  curPatch = maxval;
297  }
298 
299  for(j = 0; j < numLocalBndCoeffsPerPatchOld[i]; j++ )
300  {
301  ASSERTL0((locPatchMask[cnt]==maxval)||(locPatchMask[cnt]==minval),
302  "These values should never be the same");
303  if(locPatchMask[cnt] == -1)
304  {
305  ++numLocalBndCoeffsPerPatchNew[curPatch];
306  }
307  cnt++;
308  }
309  }
310 
311  // Count how many local dofs of the old level are mapped
312  // to the local boundary dofs of the new level
314  m_numPatches = numLocalBndCoeffsPerPatchNew.size();
315  m_numLocalBndCoeffsPerPatch = Array<OneD, unsigned int>(m_numPatches);
316  m_numLocalIntCoeffsPerPatch = Array<OneD, unsigned int>(m_numPatches,0u);
317  for(int i = 0; i < m_numPatches; i++)
318  {
319  m_numLocalBndCoeffsPerPatch[i] = (unsigned int) numLocalBndCoeffsPerPatchNew[i];
320  m_numLocalBndCoeffs += numLocalBndCoeffsPerPatchNew[i];
321  }
322  multiLevelGraph->GetNintDofsPerPatch(newLevel,m_numLocalIntCoeffsPerPatch);
323 
324  // Also initialise some more data members
325  m_solnType = solnTypeOld;
330  "This method should only be called for in "
331  "case of multi-level static condensation.");
332  m_staticCondLevel = newLevel;
333  m_signChange = signChangeOld;
334  m_numLocalDirBndCoeffs = numLocalDirBndCoeffsOld;
335  m_numGlobalDirBndCoeffs = numGlobalDirBndCoeffsOld;
336  m_numGlobalBndCoeffs = multiLevelGraph->GetInteriorOffset(newLevel) +
338  m_numGlobalCoeffs = multiLevelGraph->GetNumGlobalDofs(newLevel) +
340  m_localToGlobalBndMap = Array<OneD,int>(m_numLocalBndCoeffs);
341  if(m_signChange)
342  {
343  m_localToGlobalBndSign = Array<OneD,NekDouble>(m_numLocalBndCoeffs);
344  }
345 
347 
348  m_globalToUniversalBndMap = Array<OneD, int>(
349  m_numGlobalBndCoeffs, oldLevelMap->GetGlobalToUniversalBndMap());
350  m_globalToUniversalBndMapUnique = Array<OneD, int>(
351  m_numGlobalBndCoeffs, oldLevelMap->GetGlobalToUniversalBndMapUnique());
352 
353  // Set up an offset array that denotes the offset of the local
354  // boundary degrees of freedom of the next level
355  Array<OneD, int> numLocalBndCoeffsPerPatchOffset(m_numPatches+1,0);
356  for(int i = 1; i < m_numPatches+1; i++)
357  {
358  numLocalBndCoeffsPerPatchOffset[i] += numLocalBndCoeffsPerPatchOffset[i-1] + numLocalBndCoeffsPerPatchNew[i-1];
359  }
360 
361  int additionalPatchCnt = numPatchesWithIntNew;
362  int newid;
363  int blockid;
364  bool isBndDof;
365  NekDouble sign;
366  Array<OneD, int> bndDofPerPatchCnt(m_numPatches,0);
367  for(i = cnt = 0; i < numPatchesOld; i++)
368  {
369  minval = *min_element(&locPatchMask[cnt],&locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
370  maxval = *max_element(&locPatchMask[cnt],&locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
371  ASSERTL0((minval==maxval)||(minval==-1),"These values should never be the same");
372 
373  if(maxval == -1)
374  {
375  curPatch = additionalPatchCnt;
376  additionalPatchCnt++;
377  }
378  else
379  {
380  curPatch = maxval;
381  }
382 
383  for(j = 0; j < numLocalBndCoeffsPerPatchOld[i]; j++ )
384  {
385  ASSERTL0((locPatchMask[cnt]==maxval)||(locPatchMask[cnt]==minval),
386  "These values should never be the same");
387 
388  sign = oldLevelMap->GetLocalToGlobalBndSign(cnt);
389 
390  if(locPatchMask[cnt] == -1)
391  {
392  newid = numLocalBndCoeffsPerPatchOffset[curPatch];
393 
394  m_localToGlobalBndMap[newid] = oldLevelMap->GetLocalToGlobalBndMap(cnt);
395  if(m_signChange)
396  {
397  m_localToGlobalBndSign[ newid ] = sign;
398  }
399 
400  blockid = bndDofPerPatchCnt[curPatch];
401  isBndDof = true;
402 
403 
404  numLocalBndCoeffsPerPatchOffset[curPatch]++;
405  bndDofPerPatchCnt[curPatch]++;
406  }
407  else
408  {
409  newid = oldLevelMap->GetLocalToGlobalBndMap(cnt) -
410  m_numGlobalBndCoeffs+m_numLocalBndCoeffs;
411 
412  blockid = oldLevelMap->GetLocalToGlobalBndMap(cnt)-
413  m_numGlobalDirBndCoeffs - multiLevelGraph->GetInteriorOffset(newLevel,curPatch);
414  isBndDof = false;
415  }
416 
417  sign = isBndDof?1.0:sign;
418 
419  m_patchMapFromPrevLevel->SetPatchMap(cnt,curPatch,blockid,isBndDof,sign);
420 
421  cnt++;
422  }
423  }
425 
426  // Postprocess the computed information - Update the old
427  // level with the mapping to new evel
428  // oldLevelMap->SetLocalBndToNextLevelMap(oldLocalBndToNewLevelMap,oldLocalBndToNewLevelSign);
429  // - Construct the next level mapping object
430  if(m_staticCondLevel < (multiLevelGraph->GetNlevels()-1))
431  {
433  }
434  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
PreconditionerType m_preconType
Type type of preconditioner to use in iterative solver.
Definition: AssemblyMap.h:370
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:308
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
#define sign(a, b)
return the sign(b)*a
Definition: Polylib.cpp:27
PatchMapSharedPtr m_patchMapFromPrevLevel
Mapping information for previous level in MultiLevel Solver.
Definition: AssemblyMap.h:410
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:397
int m_successiveRHS
sucessive RHS for iterative solver
Definition: AssemblyMap.h:379
size_t m_hash
Hash for map.
Definition: AssemblyMap.h:311
NekDouble m_iterativeTolerance
Tolerance for iterative solver.
Definition: AssemblyMap.h:376
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:390
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:365
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:320
double NekDouble
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
The number of int dofs per patch.
Definition: AssemblyMap.h:392
int m_lowestStaticCondLevel
Lowest static condensation level.
Definition: AssemblyMap.h:399
void CalculateBndSystemBandWidth()
Calculates the bandwidth of the boundary system.
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
int m_staticCondLevel
The level of recursion in the case of multi-level static condensation.
Definition: AssemblyMap.h:386
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
int m_maxIterations
Maximum iterations for iterative solver.
Definition: AssemblyMap.h:373
int RoundNekDoubleToInt(NekDouble x)
Rounds a double precision number to an integer.
Definition: AssemblyMap.cpp:57
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:305
Array< OneD, int > m_globalToUniversalBndMap
Integer map of process coeffs to universal space.
Definition: AssemblyMap.h:360
Array< OneD, int > m_globalToUniversalBndMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMap.h:362
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:344
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:388
Nektar::MultiRegions::AssemblyMap::~AssemblyMap ( void  )
virtual

Destructor.

Definition at line 437 of file AssemblyMap.cpp.

438  {
439  }

Member Function Documentation

void Nektar::MultiRegions::AssemblyMap::Assemble ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const

Definition at line 754 of file AssemblyMap.cpp.

References v_Assemble().

Referenced by Nektar::MultiRegions::AssemblyMapCG::v_Assemble().

757  {
758  v_Assemble(loc,global);
759  }
virtual void v_Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
void Nektar::MultiRegions::AssemblyMap::Assemble ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const

Definition at line 761 of file AssemblyMap.cpp.

References v_Assemble().

764  {
765  v_Assemble(loc,global);
766  }
virtual void v_Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
int  offset 
) const

Definition at line 1074 of file AssemblyMap.cpp.

References Nektar::NekVector< DataType >::GetPtr().

Referenced by AssembleBnd(), Nektar::MultiRegions::AssemblyMapDG::v_Assemble(), and Nektar::MultiRegions::AssemblyMapDG::v_LocalToGlobal().

1077  {
1078  AssembleBnd(loc.GetPtr(), global.GetPtr(), offset);
1079  }
void AssembleBnd(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:230
void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const

Definition at line 1082 of file AssemblyMap.cpp.

References AssembleBnd(), and Nektar::NekVector< DataType >::GetPtr().

1085  {
1086  AssembleBnd(loc.GetPtr(), global.GetPtr());
1087  }
void AssembleBnd(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:230
void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
int  offset 
) const

Definition at line 1090 of file AssemblyMap.cpp.

References ASSERTL1, Vmath::Assmb(), m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, m_signChange, UniversalAssembleBnd(), and Vmath::Vcopy().

1093  {
1094  ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local array is not of correct dimension");
1095  ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs-offset,"Global array is not of correct dimension");
1096  Array<OneD,NekDouble> tmp(m_numGlobalBndCoeffs,0.0);
1097 
1098  if(m_signChange)
1099  {
1101  }
1102  else
1103  {
1104  Vmath::Assmb(m_numLocalBndCoeffs,loc.get(), m_localToGlobalBndMap.get(), tmp.get());
1105  }
1106  UniversalAssembleBnd(tmp);
1107  Vmath::Vcopy(m_numGlobalBndCoeffs-offset, tmp.get() + offset, 1, global.get(), 1);
1108  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
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:709
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const

Definition at line 1111 of file AssemblyMap.cpp.

References ASSERTL1, Vmath::Assmb(), m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, m_signChange, UniversalAssembleBnd(), and Vmath::Zero().

1114  {
1115  ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
1116  ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
1117 
1118  Vmath::Zero(m_numGlobalBndCoeffs, global.get(), 1);
1119 
1120  if(m_signChange)
1121  {
1123  loc.get(), m_localToGlobalBndMap.get(), global.get());
1124  }
1125  else
1126  {
1127  Vmath::Assmb(m_numLocalBndCoeffs,loc.get(), m_localToGlobalBndMap.get(), global.get());
1128  }
1129  UniversalAssembleBnd(global);
1130  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
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:709
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:373
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
bool Nektar::MultiRegions::AssemblyMap::AtLastLevel ( ) const

Returns true if this is the last level in the multi-level static condensation.

Definition at line 1197 of file AssemblyMap.cpp.

References m_nextLevelLocalToGlobalMap.

1198  {
1199  return !( m_nextLevelLocalToGlobalMap.get() );
1200  }
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:397
void Nektar::MultiRegions::AssemblyMap::CalculateBndSystemBandWidth ( )
protected

Calculates the bandwidth of the boundary system.

The bandwidth calculated 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 here calculate the bandwith of the global boundary system (as used for static condensation).

Definition at line 456 of file AssemblyMap.cpp.

References m_bndSystemBandWidth, m_localToGlobalBndMap, m_numGlobalDirBndCoeffs, m_numLocalBndCoeffs, m_numLocalBndCoeffsPerPatch, and m_numPatches.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), and Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG().

457  {
458  int i,j;
459  int cnt = 0;
460  int locSize;
461  int maxId;
462  int minId;
463  int bwidth = -1;
464  for(i = 0; i < m_numPatches; ++i)
465  {
466  locSize = m_numLocalBndCoeffsPerPatch[i];
467  maxId = -1;
468  minId = m_numLocalBndCoeffs+1;
469  for(j = 0; j < locSize; j++)
470  {
472  {
473  if(m_localToGlobalBndMap[cnt+j] > maxId)
474  {
475  maxId = m_localToGlobalBndMap[cnt+j];
476  }
477 
478  if(m_localToGlobalBndMap[cnt+j] < minId)
479  {
480  minId = m_localToGlobalBndMap[cnt+j];
481  }
482  }
483  }
484  bwidth = (bwidth>(maxId-minId))?bwidth:(maxId-minId);
485 
486  cnt+=locSize;
487  }
488 
489  m_bndSystemBandWidth = bwidth;
490  }
int m_bndSystemBandWidth
The bandwith of the global bnd system.
Definition: AssemblyMap.h:367
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:390
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:320
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:388
int Nektar::MultiRegions::AssemblyMap::GetBndCondCoeffsToGlobalCoeffsMap ( const int  i)

Retrieves the global index corresponding to a boundary expansion mode.

Definition at line 883 of file AssemblyMap.cpp.

References m_bndCondCoeffsToGlobalCoeffsMap.

885  {
887  }
Array< OneD, int > m_bndCondCoeffsToGlobalCoeffsMap
Integer map of bnd cond coeffs to global coefficients.
Definition: AssemblyMap.h:354
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetBndCondCoeffsToGlobalCoeffsMap ( )

Retrieves the global indices corresponding to the boundary expansion modes.

Definition at line 917 of file AssemblyMap.cpp.

References m_bndCondCoeffsToGlobalCoeffsMap.

Referenced by Nektar::CoupledAssemblyMap::CoupledAssemblyMap().

918  {
920  }
Array< OneD, int > m_bndCondCoeffsToGlobalCoeffsMap
Integer map of bnd cond coeffs to global coefficients.
Definition: AssemblyMap.h:354
NekDouble Nektar::MultiRegions::AssemblyMap::GetBndCondCoeffsToGlobalCoeffsSign ( const int  i)

Returns the modal sign associated with a given boundary expansion mode.

Definition at line 904 of file AssemblyMap.cpp.

References m_bndCondCoeffsToGlobalCoeffsSign, and m_signChange.

905  {
906  if(m_signChange)
907  {
909  }
910  else
911  {
912  return 1.0;
913  }
914  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
Array< OneD, NekDouble > m_bndCondCoeffsToGlobalCoeffsSign
Integer map of bnd cond coeffs to global coefficients.
Definition: AssemblyMap.h:356
int Nektar::MultiRegions::AssemblyMap::GetBndCondTraceToGlobalTraceMap ( const int  i)

Returns the global index of the boundary trace giving the index on the boundary expansion.

Definition at line 890 of file AssemblyMap.cpp.

References ASSERTL1, and m_bndCondTraceToGlobalTraceMap.

892  {
893  ASSERTL1(i < m_bndCondTraceToGlobalTraceMap.num_elements(),
894  "Index out of range.");
896  }
Array< OneD, int > m_bndCondTraceToGlobalTraceMap
Integer map of bnd cond trace number to global trace number.
Definition: AssemblyMap.h:358
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetBndCondTraceToGlobalTraceMap ( )

Definition at line 899 of file AssemblyMap.cpp.

900  {
902  }
Array< OneD, int > m_bndCondTraceToGlobalTraceMap
Integer map of bnd cond trace number to global trace number.
Definition: AssemblyMap.h:358
int Nektar::MultiRegions::AssemblyMap::GetBndSystemBandWidth ( ) const

Returns the bandwidth of the boundary system.

Definition at line 1156 of file AssemblyMap.cpp.

References m_bndSystemBandWidth.

Referenced by Nektar::MultiRegions::AssemblyMapDG::v_GetFullSystemBandWidth().

1157  {
1158  return m_bndSystemBandWidth;
1159  }
int m_bndSystemBandWidth
The bandwith of the global bnd system.
Definition: AssemblyMap.h:367
LibUtilities::CommSharedPtr Nektar::MultiRegions::AssemblyMap::GetComm ( )

Retrieves the communicator.

Definition at line 674 of file AssemblyMap.cpp.

References m_comm.

675  {
676  return m_comm;
677  }
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:308
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetExtraDirEdges ( )

Definition at line 827 of file AssemblyMap.cpp.

References v_GetExtraDirEdges().

828  {
829  return v_GetExtraDirEdges();
830  }
virtual const Array< OneD, const int > & v_GetExtraDirEdges()
int Nektar::MultiRegions::AssemblyMap::GetFullSystemBandWidth ( ) const

Definition at line 787 of file AssemblyMap.cpp.

References v_GetFullSystemBandWidth().

788  {
789  return v_GetFullSystemBandWidth();
790  }
virtual int v_GetFullSystemBandWidth() const
GlobalSysSolnType Nektar::MultiRegions::AssemblyMap::GetGlobalSysSolnType ( ) const

Returns the method of solving global systems.

Definition at line 1203 of file AssemblyMap.cpp.

References m_solnType.

Referenced by AssemblyMap().

1204  {
1205  return m_solnType;
1206  }
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:365
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalBndMap ( )

Definition at line 860 of file AssemblyMap.cpp.

References m_globalToUniversalBndMap.

Referenced by AssemblyMap().

861  {
863  }
Array< OneD, int > m_globalToUniversalBndMap
Integer map of process coeffs to universal space.
Definition: AssemblyMap.h:360
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalBndMapUnique ( )

Definition at line 865 of file AssemblyMap.cpp.

References m_globalToUniversalBndMapUnique.

Referenced by AssemblyMap().

866  {
868  }
Array< OneD, int > m_globalToUniversalBndMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMap.h:362
int Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMap ( const int  i) const

Definition at line 689 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMap().

690  {
691  return v_GetGlobalToUniversalMap(i);
692  }
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMap()
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMap ( )

Definition at line 704 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMap().

705  {
706  return v_GetGlobalToUniversalMap();
707  }
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMap()
int Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMapUnique ( const int  i) const

Definition at line 694 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMapUnique().

695  {
697  }
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMapUnique()
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMapUnique ( )

Definition at line 709 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMapUnique().

710  {
712  }
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMapUnique()
size_t Nektar::MultiRegions::AssemblyMap::GetHash ( ) const

Retrieves the hash of this map.

Definition at line 679 of file AssemblyMap.cpp.

References m_hash.

680  {
681  return m_hash;
682  }
size_t m_hash
Hash for map.
Definition: AssemblyMap.h:311
NekDouble Nektar::MultiRegions::AssemblyMap::GetIterativeTolerance ( ) const

Definition at line 1213 of file AssemblyMap.cpp.

References m_iterativeTolerance.

1214  {
1215  return m_iterativeTolerance;
1216  }
NekDouble m_iterativeTolerance
Tolerance for iterative solver.
Definition: AssemblyMap.h:376
int Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndMap ( const int  i) const

Retrieve the global index of a given local boundary mode.

Definition at line 837 of file AssemblyMap.cpp.

References m_localToGlobalBndMap.

Referenced by AssemblyMap().

838  {
839  return m_localToGlobalBndMap[i];
840  }
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndMap ( void  )

Retrieve the global indices of the local boundary modes.

Definition at line 843 of file AssemblyMap.cpp.

References m_localToGlobalBndMap.

844  {
845  return m_localToGlobalBndMap;
846  }
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
NekDouble Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndSign ( const int  i) const

Retrieve the sign change of a given local boundary mode.

Definition at line 870 of file AssemblyMap.cpp.

References m_localToGlobalBndSign, and m_signChange.

Referenced by AssemblyMap().

871  {
872  if(m_signChange)
873  {
874  return m_localToGlobalBndSign[i];
875  }
876  else
877  {
878  return 1.0;
879  }
880  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
Array< OneD, const NekDouble > Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndSign ( void  ) const

Retrieve the sign change for all local boundary modes.

Definition at line 855 of file AssemblyMap.cpp.

References m_localToGlobalBndSign.

Referenced by Nektar::MultiRegions::AssemblyMapDG::v_GetLocalToGlobalSign().

856  {
857  return m_localToGlobalBndSign;
858  }
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
int Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalMap ( const int  i) const

Definition at line 684 of file AssemblyMap.cpp.

References v_GetLocalToGlobalMap().

685  {
686  return v_GetLocalToGlobalMap(i);
687  }
virtual const Array< OneD, const int > & v_GetLocalToGlobalMap()
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalMap ( )

Definition at line 699 of file AssemblyMap.cpp.

References v_GetLocalToGlobalMap().

700  {
701  return v_GetLocalToGlobalMap();
702  }
virtual const Array< OneD, const int > & v_GetLocalToGlobalMap()
NekDouble Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalSign ( const int  i) const

Definition at line 714 of file AssemblyMap.cpp.

References v_GetLocalToGlobalSign().

715  {
716  return v_GetLocalToGlobalSign(i);
717  }
virtual const Array< OneD, NekDouble > & v_GetLocalToGlobalSign() const
const Array< OneD, NekDouble > & Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalSign ( ) const

Definition at line 719 of file AssemblyMap.cpp.

References v_GetLocalToGlobalSign().

720  {
721  return v_GetLocalToGlobalSign();
722  }
virtual const Array< OneD, NekDouble > & v_GetLocalToGlobalSign() const
int Nektar::MultiRegions::AssemblyMap::GetLowestStaticCondLevel ( ) const
inline

Definition at line 298 of file AssemblyMap.h.

References m_lowestStaticCondLevel.

299  {
301  }
int m_lowestStaticCondLevel
Lowest static condensation level.
Definition: AssemblyMap.h:399
int Nektar::MultiRegions::AssemblyMap::GetMaxIterations ( ) const

Definition at line 1218 of file AssemblyMap.cpp.

References m_maxIterations.

1219  {
1220  return m_maxIterations;
1221  }
int m_maxIterations
Maximum iterations for iterative solver.
Definition: AssemblyMap.h:373
const AssemblyMapSharedPtr Nektar::MultiRegions::AssemblyMap::GetNextLevelLocalToGlobalMap ( ) const

Returns the local to global mapping for the next level in the multi-level static condensation.

Definition at line 1185 of file AssemblyMap.cpp.

References m_nextLevelLocalToGlobalMap.

1186  {
1188  }
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:397
int Nektar::MultiRegions::AssemblyMap::GetNumDirEdges ( ) const

Definition at line 807 of file AssemblyMap.cpp.

References v_GetNumDirEdges().

808  {
809  return v_GetNumDirEdges();
810  }
virtual int v_GetNumDirEdges() const
int Nektar::MultiRegions::AssemblyMap::GetNumDirFaces ( ) const

Definition at line 812 of file AssemblyMap.cpp.

References v_GetNumDirFaces().

813  {
814  return v_GetNumDirFaces();
815  }
virtual int v_GetNumDirFaces() const
int Nektar::MultiRegions::AssemblyMap::GetNumGlobalBndCoeffs ( ) const

Returns the total number of global boundary coefficients.

Definition at line 939 of file AssemblyMap.cpp.

References m_numGlobalBndCoeffs.

Referenced by AssemblyMap().

940  {
941  return m_numGlobalBndCoeffs;
942  }
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
int Nektar::MultiRegions::AssemblyMap::GetNumGlobalCoeffs ( ) const

Returns the total number of global coefficients.

Definition at line 949 of file AssemblyMap.cpp.

References m_numGlobalCoeffs.

950  {
951  return m_numGlobalCoeffs;
952  }
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:344
int Nektar::MultiRegions::AssemblyMap::GetNumGlobalDirBndCoeffs ( ) const

Returns the number of global Dirichlet boundary coefficients.

Definition at line 923 of file AssemblyMap.cpp.

References m_numGlobalDirBndCoeffs.

Referenced by AssemblyMap().

924  {
926  }
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:320
int Nektar::MultiRegions::AssemblyMap::GetNumLocalBndCoeffs ( ) const

Returns the total number of local boundary coefficients.

Definition at line 934 of file AssemblyMap.cpp.

References m_numLocalBndCoeffs.

Referenced by AssemblyMap().

935  {
936  return m_numLocalBndCoeffs;
937  }
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
const Array< OneD, const unsigned int > & Nektar::MultiRegions::AssemblyMap::GetNumLocalBndCoeffsPerPatch ( )

Returns the number of local boundary coefficients in each patch.

Definition at line 1172 of file AssemblyMap.cpp.

References m_numLocalBndCoeffsPerPatch.

Referenced by AssemblyMap().

1173  {
1175  }
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:390
int Nektar::MultiRegions::AssemblyMap::GetNumLocalCoeffs ( ) const

Returns the total number of local coefficients.

Definition at line 944 of file AssemblyMap.cpp.

References m_numLocalCoeffs.

945  {
946  return m_numLocalCoeffs;
947  }
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:333
int Nektar::MultiRegions::AssemblyMap::GetNumLocalDirBndCoeffs ( ) const

Returns the number of local Dirichlet boundary coefficients.

Definition at line 929 of file AssemblyMap.cpp.

References m_numLocalDirBndCoeffs.

Referenced by AssemblyMap().

930  {
931  return m_numLocalDirBndCoeffs;
932  }
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
const Array< OneD, const unsigned int > & Nektar::MultiRegions::AssemblyMap::GetNumLocalIntCoeffsPerPatch ( )

Returns the number of local interior coefficients in each patch.

Definition at line 1179 of file AssemblyMap.cpp.

References m_numLocalIntCoeffsPerPatch.

1180  {
1182  }
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
The number of int dofs per patch.
Definition: AssemblyMap.h:392
int Nektar::MultiRegions::AssemblyMap::GetNumNonDirEdgeModes ( ) const

Definition at line 797 of file AssemblyMap.cpp.

References v_GetNumNonDirEdgeModes().

798  {
799  return v_GetNumNonDirEdgeModes();
800  }
virtual int v_GetNumNonDirEdgeModes() const
int Nektar::MultiRegions::AssemblyMap::GetNumNonDirEdges ( ) const

Definition at line 817 of file AssemblyMap.cpp.

References v_GetNumNonDirEdges().

818  {
819  return v_GetNumNonDirEdges();
820  }
virtual int v_GetNumNonDirEdges() const
int Nektar::MultiRegions::AssemblyMap::GetNumNonDirFaceModes ( ) const

Definition at line 802 of file AssemblyMap.cpp.

References v_GetNumNonDirFaceModes().

803  {
804  return v_GetNumNonDirFaceModes();
805  }
virtual int v_GetNumNonDirFaceModes() const
int Nektar::MultiRegions::AssemblyMap::GetNumNonDirFaces ( ) const

Definition at line 822 of file AssemblyMap.cpp.

References v_GetNumNonDirFaces().

823  {
824  return v_GetNumNonDirFaces();
825  }
virtual int v_GetNumNonDirFaces() const
int Nektar::MultiRegions::AssemblyMap::GetNumNonDirVertexModes ( ) const

Definition at line 792 of file AssemblyMap.cpp.

References v_GetNumNonDirVertexModes().

793  {
794  return v_GetNumNonDirVertexModes();
795  }
virtual int v_GetNumNonDirVertexModes() const
int Nektar::MultiRegions::AssemblyMap::GetNumPatches ( ) const

Returns the number of patches in this static condensation level.

Definition at line 1166 of file AssemblyMap.cpp.

References m_numPatches.

Referenced by AssemblyMap().

1167  {
1168  return m_numPatches;
1169  }
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:388
const PatchMapSharedPtr & Nektar::MultiRegions::AssemblyMap::GetPatchMapFromPrevLevel ( void  ) const

Returns the patch map from the previous level of the multi-level static condensation.

Definition at line 1191 of file AssemblyMap.cpp.

References m_patchMapFromPrevLevel.

1193  {
1194  return m_patchMapFromPrevLevel;
1195  }
PatchMapSharedPtr m_patchMapFromPrevLevel
Mapping information for previous level in MultiLevel Solver.
Definition: AssemblyMap.h:410
PreconditionerType Nektar::MultiRegions::AssemblyMap::GetPreconType ( ) const

Definition at line 1208 of file AssemblyMap.cpp.

References m_preconType.

1209  {
1210  return m_preconType;
1211  }
PreconditionerType m_preconType
Type type of preconditioner to use in iterative solver.
Definition: AssemblyMap.h:370
bool Nektar::MultiRegions::AssemblyMap::GetSignChange ( )

Returns true if using a modal expansion requiring a change of sign of some modes.

Definition at line 848 of file AssemblyMap.cpp.

References m_signChange.

Referenced by AssemblyMap().

849  {
850  return m_signChange;
851  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
bool Nektar::MultiRegions::AssemblyMap::GetSingularSystem ( ) const

Retrieves if the system is singular (true) or not (false)

Definition at line 954 of file AssemblyMap.cpp.

References m_systemSingular.

955  {
956  return m_systemSingular;
957  }
bool m_systemSingular
Flag indicating if the system is singular or not.
Definition: AssemblyMap.h:322
int Nektar::MultiRegions::AssemblyMap::GetStaticCondLevel ( ) const

Returns the level of static condensation for this map.

Definition at line 1161 of file AssemblyMap.cpp.

References m_staticCondLevel.

Referenced by AssemblyMap().

1162  {
1163  return m_staticCondLevel;
1164  }
int m_staticCondLevel
The level of recursion in the case of multi-level static condensation.
Definition: AssemblyMap.h:386
int Nektar::MultiRegions::AssemblyMap::GetSuccessiveRHS ( ) const

Definition at line 1223 of file AssemblyMap.cpp.

References m_successiveRHS.

1224  {
1225  return m_successiveRHS;
1226  }
int m_successiveRHS
sucessive RHS for iterative solver
Definition: AssemblyMap.h:379
void Nektar::MultiRegions::AssemblyMap::GlobalToLocal ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc 
) const

Definition at line 740 of file AssemblyMap.cpp.

References v_GlobalToLocal().

Referenced by Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal().

743  {
744  v_GlobalToLocal(global,loc);
745  }
virtual void v_GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
void Nektar::MultiRegions::AssemblyMap::GlobalToLocal ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const

Definition at line 747 of file AssemblyMap.cpp.

References v_GlobalToLocal().

750  {
751  v_GlobalToLocal(global,loc);
752  }
virtual void v_GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc,
int  offset 
) const

Definition at line 959 of file AssemblyMap.cpp.

References Nektar::NekVector< DataType >::GetPtr().

Referenced by GlobalToLocalBnd(), and Nektar::MultiRegions::AssemblyMapDG::v_GlobalToLocal().

963  {
964  GlobalToLocalBnd(global.GetPtr(), loc.GetPtr(), offset);
965  }
void GlobalToLocalBnd(const NekVector< NekDouble > &global, NekVector< NekDouble > &loc, int offset) const
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:230
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const

Definition at line 968 of file AssemblyMap.cpp.

References Nektar::NekVector< DataType >::GetPtr(), and GlobalToLocalBnd().

971  {
972  GlobalToLocalBnd(global.GetPtr(), loc.GetPtr());
973  }
void GlobalToLocalBnd(const NekVector< NekDouble > &global, NekVector< NekDouble > &loc, int offset) const
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:230
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc,
int  offset 
) const

Definition at line 976 of file AssemblyMap.cpp.

References ASSERTL1, Vmath::Gathr(), m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, m_signChange, and Vmath::Vcopy().

979  {
980  ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
981  ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs-offset,"Global vector is not of correct dimension");
982 
983  // offset input data by length "offset" for Dirichlet boundary conditions.
984  Array<OneD,NekDouble> tmp(m_numGlobalBndCoeffs,0.0);
985  Vmath::Vcopy(m_numGlobalBndCoeffs-offset, global.get(), 1, tmp.get() + offset, 1);
986 
987  if(m_signChange)
988  {
990  }
991  else
992  {
993  Vmath::Gathr(m_numLocalBndCoeffs, tmp.get(), m_localToGlobalBndMap.get(), loc.get());
994  }
995  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
void Gathr(int n, const T *x, const int *y, T *z)
Gather vector z[i] = x[y[i]].
Definition: Vmath.cpp:644
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc 
) const

Definition at line 998 of file AssemblyMap.cpp.

References ASSERTL1, Vmath::Gathr(), m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, and m_signChange.

1001  {
1002  ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
1003  ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
1004 
1005  if(m_signChange)
1006  {
1007  Vmath::Gathr(m_numLocalBndCoeffs, m_localToGlobalBndSign.get(), global.get(), m_localToGlobalBndMap.get(), loc.get());
1008  }
1009  else
1010  {
1011  Vmath::Gathr(m_numLocalBndCoeffs, global.get(), m_localToGlobalBndMap.get(), loc.get());
1012  }
1013  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
void Gathr(int n, const T *x, const int *y, T *z)
Gather vector z[i] = x[y[i]].
Definition: Vmath.cpp:644
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBndWithoutSign ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc 
)
protected

Definition at line 1228 of file AssemblyMap.cpp.

References ASSERTL1, Vmath::Gathr(), m_localToGlobalBndMap, m_numGlobalBndCoeffs, and m_numLocalBndCoeffs.

Referenced by AssemblyMap().

1231  {
1232  ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
1233  ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
1234 
1235  Vmath::Gathr(m_numLocalBndCoeffs, global.get(), m_localToGlobalBndMap.get(), loc.get());
1236  }
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
void Gathr(int n, const T *x, const int *y, T *z)
Gather vector z[i] = x[y[i]].
Definition: Vmath.cpp:644
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
boost::shared_ptr< AssemblyMap > Nektar::MultiRegions::AssemblyMap::LinearSpaceMap ( const ExpList locexp,
GlobalSysSolnType  solnType 
)

Definition at line 832 of file AssemblyMap.cpp.

References v_LinearSpaceMap().

833  {
834  return v_LinearSpaceMap(locexp, solnType);
835  }
virtual boost::shared_ptr< AssemblyMap > v_LinearSpaceMap(const ExpList &locexp, GlobalSysSolnType solnType)
Generate a linear space mapping from existing mapping.
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
int  offset 
) const

Definition at line 1015 of file AssemblyMap.cpp.

References Nektar::NekVector< DataType >::GetPtr().

Referenced by LocalBndToGlobal().

1019  {
1020  LocalBndToGlobal(loc.GetPtr(), global.GetPtr(), offset);
1021  }
void LocalBndToGlobal(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:230
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const

Definition at line 1048 of file AssemblyMap.cpp.

References Nektar::NekVector< DataType >::GetPtr(), and LocalBndToGlobal().

1051  {
1052  LocalBndToGlobal(loc.GetPtr(), global.GetPtr());
1053  }
void LocalBndToGlobal(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:230
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
int  offset 
) const

Definition at line 1024 of file AssemblyMap.cpp.

References ASSERTL1, m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, m_signChange, Vmath::Scatr(), UniversalAssembleBnd(), and Vmath::Vcopy().

1028  {
1029  ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
1030  ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs-offset,"Global vector is not of correct dimension");
1031 
1032  // offset input data by length "offset" for Dirichlet boundary conditions.
1033  Array<OneD,NekDouble> tmp(m_numGlobalBndCoeffs,0.0);
1034 
1035  if(m_signChange)
1036  {
1038  }
1039  else
1040  {
1041  Vmath::Scatr(m_numLocalBndCoeffs, loc.get(), m_localToGlobalBndMap.get(), tmp.get());
1042  }
1043 
1044  UniversalAssembleBnd(tmp);
1045  Vmath::Vcopy(m_numGlobalBndCoeffs-offset, tmp.get()+offset, 1, global.get(), 1);
1046  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
void Scatr(int n, const T *x, const int *y, T *z)
Scatter vector z[y[i]] = x[i].
Definition: Vmath.cpp:673
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const

Definition at line 1056 of file AssemblyMap.cpp.

References ASSERTL1, m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, m_signChange, and Vmath::Scatr().

1059  {
1060  ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
1061  ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
1062 
1063  if(m_signChange)
1064  {
1065  Vmath::Scatr(m_numLocalBndCoeffs, m_localToGlobalBndSign.get(), loc.get(), m_localToGlobalBndMap.get(), global.get());
1066  }
1067  else
1068  {
1069  Vmath::Scatr(m_numLocalBndCoeffs, loc.get(), m_localToGlobalBndMap.get(), global.get());
1070  }
1071  }
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:347
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
void Scatr(int n, const T *x, const int *y, T *z)
Scatter vector z[y[i]] = x[i].
Definition: Vmath.cpp:673
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:352
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Nektar::MultiRegions::AssemblyMap::LocalToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
bool  useComm = true 
) const

Definition at line 724 of file AssemblyMap.cpp.

References v_LocalToGlobal().

Referenced by Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal().

728  {
729  v_LocalToGlobal(loc,global,useComm);
730  }
virtual void v_LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm) const
void Nektar::MultiRegions::AssemblyMap::LocalToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
bool  useComm = true 
) const

Definition at line 732 of file AssemblyMap.cpp.

References v_LocalToGlobal().

736  {
737  v_LocalToGlobal(loc,global,useComm);
738  }
virtual void v_LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm) const
void Nektar::MultiRegions::AssemblyMap::PrintStats ( std::ostream &  out,
std::string  variable,
bool  printHeader = true 
) const

Definition at line 1238 of file AssemblyMap.cpp.

References Vmath::Assmb(), m_globalToUniversalBndMapUnique, m_localToGlobalBndMap, m_nextLevelLocalToGlobalMap, m_numGlobalBndCoeffs, m_numGlobalCoeffs, m_numGlobalDirBndCoeffs, m_numLocalBndCoeffs, m_numLocalCoeffs, m_numLocalDirBndCoeffs, m_session, Nektar::LibUtilities::ReduceMax, Nektar::LibUtilities::ReduceMin, Nektar::LibUtilities::ReduceSum, and UniversalAssembleBnd().

1240  {
1242  = m_session->GetComm()->GetRowComm();
1243  bool isRoot = vRowComm->GetRank() == 0;
1244  int n = vRowComm->GetSize();
1245  int i;
1246 
1247  // Determine number of global degrees of freedom.
1248  int globBndCnt = 0, globDirCnt = 0;
1249 
1250  for (i = 0; i < m_numGlobalBndCoeffs; ++i)
1251  {
1253  {
1254  globBndCnt++;
1255 
1256  if (i < m_numGlobalDirBndCoeffs)
1257  {
1258  globDirCnt++;
1259  }
1260  }
1261  }
1262 
1263  int globCnt = m_numGlobalCoeffs - m_numGlobalBndCoeffs + globBndCnt;
1264 
1265  // Calculate maximum valency
1266  Array<OneD, NekDouble> tmpLoc (m_numLocalBndCoeffs, 1.0);
1267  Array<OneD, NekDouble> tmpGlob(m_numGlobalBndCoeffs, 0.0);
1268 
1269  Vmath::Assmb(
1270  m_numLocalBndCoeffs, tmpLoc.get(),
1271  m_localToGlobalBndMap.get(), tmpGlob.get());
1272  UniversalAssembleBnd(tmpGlob);
1273 
1274  int totGlobDof = globCnt;
1275  int totGlobBndDof = globBndCnt;
1276  int totGlobDirDof = globDirCnt;
1277  int totLocalDof = m_numLocalCoeffs;
1278  int totLocalBndDof = m_numLocalBndCoeffs;
1279  int totLocalDirDof = m_numLocalDirBndCoeffs;
1280 
1281  int meanValence = 0;
1282  int maxValence = 0;
1283  int minValence = 10000000;
1284  for (int i = 0; i < m_numGlobalBndCoeffs; ++i)
1285  {
1287  {
1288  continue;
1289  }
1290 
1291  if (tmpGlob[i] > maxValence)
1292  {
1293  maxValence = tmpGlob[i];
1294  }
1295  if (tmpGlob[i] < minValence)
1296  {
1297  minValence = tmpGlob[i];
1298  }
1299  meanValence += tmpGlob[i];
1300  }
1301 
1302  vRowComm->AllReduce(maxValence, LibUtilities::ReduceMax);
1303  vRowComm->AllReduce(minValence, LibUtilities::ReduceMin);
1304  vRowComm->AllReduce(meanValence, LibUtilities::ReduceSum);
1305  vRowComm->AllReduce(totGlobDof, LibUtilities::ReduceSum);
1306  vRowComm->AllReduce(totGlobBndDof, LibUtilities::ReduceSum);
1307  vRowComm->AllReduce(totGlobDirDof, LibUtilities::ReduceSum);
1308  vRowComm->AllReduce(totLocalDof, LibUtilities::ReduceSum);
1309  vRowComm->AllReduce(totLocalBndDof, LibUtilities::ReduceSum);
1310  vRowComm->AllReduce(totLocalDirDof, LibUtilities::ReduceSum);
1311 
1312  meanValence /= totGlobBndDof;
1313 
1314  if (isRoot)
1315  {
1316  if (printHeader)
1317  {
1318  out << "Assembly map statistics for field " << variable
1319  << ":" << endl;
1320  }
1321 
1322  out << " - Number of local/global dof : "
1323  << totLocalDof << " " << totGlobDof << endl;
1324  out << " - Number of local/global boundary dof : "
1325  << totLocalBndDof << " " << totGlobBndDof << endl;
1326  out << " - Number of local/global Dirichlet dof : "
1327  << totLocalDirDof << " " << totGlobDirDof << endl;
1328  out << " - dof valency (min/max/mean) : "
1329  << minValence << " " << maxValence << " " << meanValence
1330  << endl;
1331 
1332  if (n > 1)
1333  {
1334  NekDouble mean = m_numLocalCoeffs, mean2 = mean * mean;
1335  NekDouble minval = mean, maxval = mean;
1336  Array<OneD, NekDouble> tmp(1);
1337 
1338  for (i = 1; i < n; ++i)
1339  {
1340  vRowComm->Recv(i, tmp);
1341  mean += tmp[0];
1342  mean2 += tmp[0]*tmp[0];
1343 
1344  if (tmp[0] > maxval)
1345  {
1346  maxval = tmp[0];
1347  }
1348  if (tmp[0] < minval)
1349  {
1350  minval = tmp[0];
1351  }
1352  }
1353 
1354  if (maxval > 0.1)
1355  {
1356  out << " - Local dof dist. (min/max/mean/dev) : "
1357  << minval << " " << maxval << " " << (mean / n)
1358  << " " << sqrt(mean2/n - mean*mean/n/n) << endl;
1359  }
1360 
1361  vRowComm->Block();
1362 
1363  mean = minval = maxval = m_numLocalBndCoeffs;
1364  mean2 = mean * mean;
1365 
1366  for (i = 1; i < n; ++i)
1367  {
1368  vRowComm->Recv(i, tmp);
1369  mean += tmp[0];
1370  mean2 += tmp[0]*tmp[0];
1371 
1372  if (tmp[0] > maxval)
1373  {
1374  maxval = tmp[0];
1375  }
1376  if (tmp[0] < minval)
1377  {
1378  minval = tmp[0];
1379  }
1380  }
1381 
1382  out << " - Local bnd dof dist. (min/max/mean/dev) : "
1383  << minval << " " << maxval << " " << (mean / n) << " "
1384  << sqrt(mean2/n - mean*mean/n/n) << endl;
1385  }
1386  }
1387  else
1388  {
1389  Array<OneD, NekDouble> tmp(1);
1390  tmp[0] = m_numLocalCoeffs;
1391  vRowComm->Send(0, tmp);
1392  vRowComm->Block();
1393  tmp[0] = m_numLocalBndCoeffs;
1394  vRowComm->Send(0, tmp);
1395  }
1396 
1397  // Either we have no more levels in the static condensation, or we
1398  // are not multi-level.
1400  {
1401  return;
1402  }
1403 
1404  int level = 2;
1406  while (tmp->m_nextLevelLocalToGlobalMap)
1407  {
1408  tmp = tmp->m_nextLevelLocalToGlobalMap;
1409  ++level;
1410  }
1411 
1412  // Print out multi-level static condensation information.
1413  if (n > 1)
1414  {
1415  if (isRoot)
1416  {
1417  NekDouble mean = level, mean2 = mean * mean;
1418  int minval = level, maxval = level;
1419 
1420  Array<OneD, NekDouble> tmpRecv(1);
1421  for (i = 1; i < n; ++i)
1422  {
1423  vRowComm->Recv(i, tmpRecv);
1424  mean += tmpRecv[0];
1425  mean2 += tmpRecv[0]*tmpRecv[0];
1426 
1427  if (tmpRecv[0] > maxval)
1428  {
1429  maxval = (int)(tmpRecv[0] + 0.5);
1430  }
1431  if (tmpRecv[0] < minval)
1432  {
1433  minval = (int)(tmpRecv[0] + 0.5);
1434  }
1435  }
1436 
1437  out << " - M-level sc. dist. (min/max/mean/dev) : "
1438  << minval << " " << maxval << " " << (mean / n) << " "
1439  << sqrt(mean2/n - mean*mean/n/n) << endl;
1440  }
1441  else
1442  {
1443  Array<OneD, NekDouble> tmpSend(1);
1444  tmpSend[0] = level;
1445  vRowComm->Send(0, tmpSend);
1446  }
1447  }
1448  else
1449  {
1450  out << " - Number of static cond. levels : "
1451  << level << endl;
1452  }
1453 
1454  if (isRoot)
1455  {
1456  out << "Stats at lowest static cond. level:" << endl;
1457  }
1458  tmp->PrintStats(out, variable, false);
1459  }
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
boost::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
Definition: AssemblyMap.h:53
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:333
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:397
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:55
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:320
double NekDouble
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:709
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:314
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:305
Array< OneD, int > m_globalToUniversalBndMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMap.h:362
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:344
void Nektar::MultiRegions::AssemblyMap::SetNextLevelLocalToGlobalMap ( AssemblyMapSharedPtr  pNextLevelLocalToGlobalMap)
void Nektar::MultiRegions::AssemblyMap::UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal) const
void Nektar::MultiRegions::AssemblyMap::UniversalAssemble ( NekVector< NekDouble > &  pGlobal) const

Definition at line 774 of file AssemblyMap.cpp.

References v_UniversalAssemble().

776  {
777  v_UniversalAssemble(pGlobal);
778  }
virtual void v_UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
void Nektar::MultiRegions::AssemblyMap::UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal,
int  offset 
) const

Definition at line 780 of file AssemblyMap.cpp.

References v_UniversalAssemble().

783  {
784  v_UniversalAssemble(pGlobal, offset);
785  }
virtual void v_UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
void Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd ( Array< OneD, NekDouble > &  pGlobal) const

Definition at line 1132 of file AssemblyMap.cpp.

References ASSERTL1, Gs::Gather(), Gs::gs_add, m_bndGsh, and m_numGlobalBndCoeffs.

Referenced by AssembleBnd(), Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), LocalBndToGlobal(), PrintStats(), and UniversalAssembleBnd().

1134  {
1135  ASSERTL1(pGlobal.num_elements() == m_numGlobalBndCoeffs,
1136  "Wrong size.");
1137  Gs::Gather(pGlobal, Gs::gs_add, m_bndGsh);
1138  }
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:316
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:239
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd ( NekVector< NekDouble > &  pGlobal) const

Definition at line 1140 of file AssemblyMap.cpp.

References Nektar::NekVector< DataType >::GetPtr(), and UniversalAssembleBnd().

1142  {
1143  UniversalAssembleBnd(pGlobal.GetPtr());
1144  }
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
Array< OneD, DataType > & GetPtr()
Definition: NekVector.cpp:230
void Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd ( Array< OneD, NekDouble > &  pGlobal,
int  offset 
) const

Definition at line 1146 of file AssemblyMap.cpp.

References UniversalAssembleBnd(), and Vmath::Vcopy().

1149  {
1150  Array<OneD, NekDouble> tmp(offset);
1151  if (offset > 0) Vmath::Vcopy(offset, pGlobal, 1, tmp, 1);
1152  UniversalAssembleBnd(pGlobal);
1153  if (offset > 0) Vmath::Vcopy(offset, tmp, 1, pGlobal, 1);
1154  }
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::MultiRegions::AssemblyMap::v_Assemble ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 575 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by Assemble().

578  {
579  ASSERTL0(false, "Not defined for this type of mapping.");
580  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::MultiRegions::AssemblyMap::v_Assemble ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 582 of file AssemblyMap.cpp.

References ASSERTL0.

585  {
586  ASSERTL0(false, "Not defined for this type of mapping.");
587  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetExtraDirEdges ( )
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 659 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetExtraDirEdges().

660  {
661  ASSERTL0(false, "Not defined for this type of mapping.");
662  static Array<OneD, const int> result;
663  return result;
664  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetFullSystemBandWidth ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 611 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetFullSystemBandWidth().

612  {
613  ASSERTL0(false, "Not defined for this type of mapping.");
614  return 0;
615  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMap ( const int  i) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 499 of file AssemblyMap.cpp.

References ASSERTL0.

500  {
501  ASSERTL0(false, "Not defined for this type of mapping.");
502  return 0;
503  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMap ( void  )
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 518 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetGlobalToUniversalMap().

519  {
520  ASSERTL0(false, "Not defined for this type of mapping.");
521  static Array<OneD, const int> result;
522  return result;
523  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMapUnique ( const int  i) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 505 of file AssemblyMap.cpp.

References ASSERTL0.

506  {
507  ASSERTL0(false, "Not defined for this type of mapping.");
508  return 0;
509  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMapUnique ( void  )
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 525 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetGlobalToUniversalMapUnique().

526  {
527  ASSERTL0(false, "Not defined for this type of mapping.");
528  static Array<OneD, const int> result;
529  return result;
530  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalMap ( const int  i) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 493 of file AssemblyMap.cpp.

References ASSERTL0.

494  {
495  ASSERTL0(false, "Not defined for this type of mapping.");
496  return 0;
497  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalMap ( void  )
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 511 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetLocalToGlobalMap().

512  {
513  ASSERTL0(false, "Not defined for this type of mapping.");
514  static Array<OneD,const int> result;
515  return result;
516  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
NekDouble Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalSign ( const int  i) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 532 of file AssemblyMap.cpp.

References ASSERTL0.

533  {
534  ASSERTL0(false, "Not defined for this type of mapping.");
535  return 0.0;
536  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
const Array< OneD, NekDouble > & Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalSign ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 538 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetLocalToGlobalSign().

539  {
540  ASSERTL0(false, "Not defined for this type of mapping.");
541  static Array<OneD, NekDouble> result;
542  return result;
543  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetNumDirEdges ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 635 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumDirEdges().

636  {
637  ASSERTL0(false, "Not defined for this type of mapping.");
638  return 0;
639  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetNumDirFaces ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 641 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumDirFaces().

642  {
643  ASSERTL0(false, "Not defined for this type of mapping.");
644  return 0;
645  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirEdgeModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 623 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirEdgeModes().

624  {
625  ASSERTL0(false, "Not defined for this type of mapping.");
626  return 0;
627  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirEdges ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 647 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirEdges().

648  {
649  ASSERTL0(false, "Not defined for this type of mapping.");
650  return 0;
651  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirFaceModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 629 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirFaceModes().

630  {
631  ASSERTL0(false, "Not defined for this type of mapping.");
632  return 0;
633  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirFaces ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 653 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirFaces().

654  {
655  ASSERTL0(false, "Not defined for this type of mapping.");
656  return 0;
657  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirVertexModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 617 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirVertexModes().

618  {
619  ASSERTL0(false, "Not defined for this type of mapping.");
620  return 0;
621  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::MultiRegions::AssemblyMap::v_GlobalToLocal ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 561 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GlobalToLocal().

564  {
565  ASSERTL0(false, "Not defined for this type of mapping.");
566  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::MultiRegions::AssemblyMap::v_GlobalToLocal ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 568 of file AssemblyMap.cpp.

References ASSERTL0.

571  {
572  ASSERTL0(false, "Not defined for this type of mapping.");
573  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
boost::shared_ptr< AssemblyMap > Nektar::MultiRegions::AssemblyMap::v_LinearSpaceMap ( const ExpList locexp,
GlobalSysSolnType  solnType 
)
privatevirtual

Generate a linear space mapping from existing mapping.

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 666 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by LinearSpaceMap().

668  {
669  ASSERTL0(false, "Not defined for this sub class");
670  static boost::shared_ptr<AssemblyMap> result;
671  return result;
672  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::MultiRegions::AssemblyMap::v_LocalToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
bool  useComm 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 545 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by LocalToGlobal().

549  {
550  ASSERTL0(false, "Not defined for this type of mapping.");
551  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::MultiRegions::AssemblyMap::v_LocalToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
bool  useComm 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 553 of file AssemblyMap.cpp.

References ASSERTL0.

557  {
558  ASSERTL0(false, "Not defined for this type of mapping.");
559  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::MultiRegions::AssemblyMap::v_UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 589 of file AssemblyMap.cpp.

Referenced by UniversalAssemble().

591  {
592  // Do nothing here since multi-level static condensation uses a
593  // AssemblyMap and thus will call this routine in serial.
594  }
void Nektar::MultiRegions::AssemblyMap::v_UniversalAssemble ( NekVector< NekDouble > &  pGlobal) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG, and Nektar::MultiRegions::AssemblyMapDG.

Definition at line 596 of file AssemblyMap.cpp.

598  {
599  // Do nothing here since multi-level static condensation uses a
600  // AssemblyMap and thus will call this routine in serial.
601  }
void Nektar::MultiRegions::AssemblyMap::v_UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal,
int  offset 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 603 of file AssemblyMap.cpp.

606  {
607  // Do nothing here since multi-level static condensation uses a
608  // AssemblyMap and thus will call this routine in serial.
609  }

Member Data Documentation

Array<OneD,int> Nektar::MultiRegions::AssemblyMap::m_bndCondCoeffsToGlobalCoeffsMap
protected
Array<OneD,NekDouble> Nektar::MultiRegions::AssemblyMap::m_bndCondCoeffsToGlobalCoeffsSign
protected

Integer map of bnd cond coeffs to global coefficients.

Definition at line 356 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), and GetBndCondCoeffsToGlobalCoeffsSign().

Array<OneD,int> Nektar::MultiRegions::AssemblyMap::m_bndCondTraceToGlobalTraceMap
protected

Integer map of bnd cond trace number to global trace number.

Definition at line 358 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), and GetBndCondTraceToGlobalTraceMap().

Gs::gs_data* Nektar::MultiRegions::AssemblyMap::m_bndGsh
protected
int Nektar::MultiRegions::AssemblyMap::m_bndSystemBandWidth
protected

The bandwith of the global bnd system.

Definition at line 367 of file AssemblyMap.h.

Referenced by CalculateBndSystemBandWidth(), and GetBndSystemBandWidth().

LibUtilities::CommSharedPtr Nektar::MultiRegions::AssemblyMap::m_comm
protected
Array<OneD,int> Nektar::MultiRegions::AssemblyMap::m_globalToUniversalBndMap
protected
Array<OneD,int> Nektar::MultiRegions::AssemblyMap::m_globalToUniversalBndMapUnique
protected
Gs::gs_data* Nektar::MultiRegions::AssemblyMap::m_gsh
protected
size_t Nektar::MultiRegions::AssemblyMap::m_hash
protected
NekDouble Nektar::MultiRegions::AssemblyMap::m_iterativeTolerance
protected

Tolerance for iterative solver.

Definition at line 376 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetIterativeTolerance().

Array<OneD,int> Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndMap
protected
Array<OneD,NekDouble> Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndSign
protected
int Nektar::MultiRegions::AssemblyMap::m_lowestStaticCondLevel
protected
int Nektar::MultiRegions::AssemblyMap::m_maxIterations
protected

Maximum iterations for iterative solver.

Definition at line 373 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetMaxIterations().

AssemblyMapSharedPtr Nektar::MultiRegions::AssemblyMap::m_nextLevelLocalToGlobalMap
protected

Map from the patches of the previous level to the patches of the current level.

The local to global mapping of the next level of recursion

Definition at line 397 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), AtLastLevel(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNextLevelLocalToGlobalMap(), and PrintStats().

int Nektar::MultiRegions::AssemblyMap::m_numGlobalBndCoeffs
protected
int Nektar::MultiRegions::AssemblyMap::m_numGlobalCoeffs
protected

Total number of global coefficients.

This corresponds to the number of total number of coefficients

  • For CG this corresponds to the total of bnd + int DOFs.
  • For DG this corresponds to the number of bnd DOFs. This means that m_numGlobalCoeffs = m_numGlobalBndCoeffs This way, we can consider the trace-system solve as a statically condensed solve without interior DOFs. This allows us to use the same global system solver for both cases.

Definition at line 344 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumGlobalCoeffs(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapCG::v_Assemble(), and Nektar::MultiRegions::AssemblyMapCG::v_LinearSpaceMap().

int Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs
protected
int Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffs
protected
Array<OneD, unsigned int> Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffsPerPatch
protected
int Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs
protected

Total number of local coefficients.

This corresponds to the number of total number of coefficients

  • For CG this corresponds to the total of bnd + int DOFs
  • For DG this corresponds to the number of bnd DOFs. This means that m_numLocalCoeffs = m_numLocalBndCoeffs This way, we can consider the trace-system solve as a statically condensed solve without interior DOFs. This allows us to use the same global system solver for both cases.

Definition at line 333 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth(), Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumLocalCoeffs(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG::v_Assemble(), Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal(), and Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal().

int Nektar::MultiRegions::AssemblyMap::m_numLocalDirBndCoeffs
protected
Array<OneD, unsigned int> Nektar::MultiRegions::AssemblyMap::m_numLocalIntCoeffsPerPatch
protected
int Nektar::MultiRegions::AssemblyMap::m_numPatches
protected
PatchMapSharedPtr Nektar::MultiRegions::AssemblyMap::m_patchMapFromPrevLevel
private

Mapping information for previous level in MultiLevel Solver.

Definition at line 410 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetPatchMapFromPrevLevel().

PreconditionerType Nektar::MultiRegions::AssemblyMap::m_preconType
protected

Type type of preconditioner to use in iterative solver.

Definition at line 370 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetPreconType().

LibUtilities::SessionReaderSharedPtr Nektar::MultiRegions::AssemblyMap::m_session
protected
bool Nektar::MultiRegions::AssemblyMap::m_signChange
protected
GlobalSysSolnType Nektar::MultiRegions::AssemblyMap::m_solnType
protected
int Nektar::MultiRegions::AssemblyMap::m_staticCondLevel
protected
int Nektar::MultiRegions::AssemblyMap::m_successiveRHS
protected

sucessive RHS for iterative solver

Definition at line 379 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetSuccessiveRHS().

bool Nektar::MultiRegions::AssemblyMap::m_systemSingular
protected