Nektar++
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:
[legend]

Public Member Functions

 AssemblyMap ()
 Default constructor. More...
 
 AssemblyMap (const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &comm, const std::string variable="DefaultVar")
 Constructor with a communicator. More...
 
 AssemblyMap (AssemblyMap *oldLevelMap, const BottomUpSubStructuredGraphSharedPtr &multiLevelGraph)
 Constructor for next level in multi-level static condensation. More...
 
virtual ~AssemblyMap ()
 Destructor. More...
 
LibUtilities::CommSharedPtr GetComm ()
 Retrieves the communicator. More...
 
size_t GetHash () const
 Retrieves the hash of this map. More...
 
int GetLocalToGlobalMap (const int i) const
 
int GetGlobalToUniversalMap (const int i) const
 
int GetGlobalToUniversalMapUnique (const int i) const
 
const Array< OneD, const int > & GetLocalToGlobalMap ()
 
const Array< OneD, const int > & GetGlobalToUniversalMap ()
 
const Array< OneD, const int > & GetGlobalToUniversalMapUnique ()
 
NekDouble GetLocalToGlobalSign (const int i) const
 
const Array< OneD, NekDouble > & GetLocalToGlobalSign () const
 
void LocalToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm=true) const
 
void LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, bool useComm=true) const
 
void GlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
void GlobalToLocal (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
 
void Assemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
void Assemble (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
void UniversalAssemble (Array< OneD, NekDouble > &pGlobal) const
 
void UniversalAssemble (NekVector< NekDouble > &pGlobal) const
 
void UniversalAssemble (Array< OneD, NekDouble > &pGlobal, int offset) const
 
void UniversalAbsMaxBnd (Array< OneD, NekDouble > &bndvals)
 
int GetLocalToGlobalBndMap (const int i) const
 Retrieve the global index of a given local boundary mode. More...
 
const Array< OneD, const int > & GetLocalToGlobalBndMap ()
 Retrieve the global indices of the local boundary modes. More...
 
const Array< OneD, const int > & GetGlobalToUniversalBndMap ()
 
const Array< OneD, const int > & GetGlobalToUniversalBndMapUnique ()
 
bool GetSignChange ()
 Returns true if using a modal expansion requiring a change of sign of some modes. More...
 
NekDouble GetLocalToGlobalBndSign (const int i) const
 Retrieve the sign change of a given local boundary mode. More...
 
Array< OneD, const NekDoubleGetLocalToGlobalBndSign () const
 Retrieve the sign change for all local boundary modes. More...
 
const Array< OneD, const int > & GetBndCondCoeffsToLocalCoeffsMap ()
 Retrieves the local indices corresponding to the boundary expansion modes. More...
 
const Array< OneD, NekDouble > & GetBndCondCoeffsToLocalCoeffsSign ()
 Returns the modal sign associated with a given boundary expansion mode. More...
 
const Array< OneD, const int > & GetBndCondCoeffsToLocalTraceMap ()
 Retrieves the local indices corresponding to the boundary expansion modes to global trace. More...
 
int GetBndCondIDToGlobalTraceID (const int i)
 Returns the global index of the boundary trace giving the index on the boundary expansion. More...
 
const Array< OneD, const int > & GetBndCondIDToGlobalTraceID ()
 
int GetNumGlobalDirBndCoeffs () const
 Returns the number of global Dirichlet boundary coefficients. More...
 
int GetNumLocalDirBndCoeffs () const
 Returns the number of local Dirichlet boundary coefficients. More...
 
int GetNumGlobalBndCoeffs () const
 Returns the total number of global boundary coefficients. More...
 
int GetNumLocalBndCoeffs () const
 Returns the total number of local boundary coefficients. More...
 
int GetNumLocalCoeffs () const
 Returns the total number of local coefficients. More...
 
int GetNumGlobalCoeffs () const
 Returns the total number of global coefficients. More...
 
bool GetSingularSystem () const
 Retrieves if the system is singular (true) or not (false) More...
 
void GlobalToLocalBnd (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc, int offset) const
 
void GlobalToLocalBnd (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
 
void GlobalToLocalBnd (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc, int offset) const
 
void GlobalToLocalBnd (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
void LocalBndToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, int offset, bool UseComm=true) const
 
void LocalBndToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool UseComm=true) const
 
void LocalToLocalBnd (const Array< OneD, const NekDouble > &local, Array< OneD, NekDouble > &locbnd) const
 
void LocalToLocalInt (const Array< OneD, const NekDouble > &local, Array< OneD, NekDouble > &locint) const
 
void LocalBndToLocal (const Array< OneD, const NekDouble > &locbnd, Array< OneD, NekDouble > &local) const
 
void LocalIntToLocal (const Array< OneD, const NekDouble > &locbnd, Array< OneD, NekDouble > &local) const
 
void AssembleBnd (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
 
void AssembleBnd (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
 
void AssembleBnd (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, int offset) const
 
void AssembleBnd (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
void UniversalAssembleBnd (Array< OneD, NekDouble > &pGlobal) const
 
void UniversalAssembleBnd (NekVector< NekDouble > &pGlobal) const
 
void UniversalAssembleBnd (Array< OneD, NekDouble > &pGlobal, int offset) const
 
int GetFullSystemBandWidth () const
 
int GetNumNonDirVertexModes () const
 
int GetNumNonDirEdgeModes () const
 
int GetNumNonDirFaceModes () const
 
int GetNumDirEdges () const
 
int GetNumDirFaces () const
 
int GetNumNonDirEdges () const
 
int GetNumNonDirFaces () const
 
void PrintStats (std::ostream &out, std::string variable, bool printHeader=true) const
 
const Array< OneD, const int > & GetExtraDirEdges ()
 
std::shared_ptr< AssemblyMapLinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType)
 
int GetBndSystemBandWidth () const
 Returns the bandwidth of the boundary system. More...
 
int GetStaticCondLevel () const
 Returns the level of static condensation for this map. More...
 
int GetNumPatches () const
 Returns the number of patches in this static condensation level. More...
 
const Array< OneD, const unsigned int > & GetNumLocalBndCoeffsPerPatch ()
 Returns the number of local boundary coefficients in each patch. More...
 
const Array< OneD, const unsigned int > & GetNumLocalIntCoeffsPerPatch ()
 Returns the number of local interior coefficients in each patch. More...
 
const AssemblyMapSharedPtr GetNextLevelLocalToGlobalMap () const
 Returns the local to global mapping for the next level in the multi-level static condensation. More...
 
void SetNextLevelLocalToGlobalMap (AssemblyMapSharedPtr pNextLevelLocalToGlobalMap)
 
const PatchMapSharedPtrGetPatchMapFromPrevLevel (void) const
 Returns the patch map from the previous level of the multi-level static condensation. More...
 
bool AtLastLevel () const
 Returns true if this is the last level in the multi-level static condensation. More...
 
GlobalSysSolnType GetGlobalSysSolnType () const
 Returns the method of solving global systems. More...
 
PreconditionerType GetPreconType () const
 
NekDouble GetIterativeTolerance () const
 
int GetMaxIterations () const
 
int GetSuccessiveRHS () const
 
std::string GetLinSysIterSolver () const
 
int GetLowestStaticCondLevel () const
 
void PatchLocalToGlobal (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 
void PatchGlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
 
void PatchAssemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
 

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 coeffs to global Boundary Dofs. More...
 
Array< OneD, NekDoublem_localToGlobalBndSign
 Integer sign of local boundary coeffs to global space. More...
 
Array< OneD, int > m_localToLocalBndMap
 Integer map of local boundary coeffs to local boundary system numbering. More...
 
Array< OneD, int > m_localToLocalIntMap
 Integer map of local boundary coeffs to local interior system numbering. More...
 
Array< OneD, int > m_bndCondCoeffsToLocalCoeffsMap
 Integer map of bnd cond coeffs to local coefficients. More...
 
Array< OneD, NekDoublem_bndCondCoeffsToLocalCoeffsSign
 Integer map of sign of bnd cond coeffs to local coefficients. More...
 
Array< OneD, int > m_bndCondCoeffsToLocalTraceMap
 Integer map of bnd cond coeff to local trace coeff. More...
 
Array< OneD, int > m_bndCondIDToGlobalTraceID
 Integer map of bnd cond trace number to global trace number. More...
 
Array< OneD, int > m_globalToUniversalBndMap
 Integer map of process coeffs to universal space. More...
 
Array< OneD, int > m_globalToUniversalBndMapUnique
 Integer map of unique process coeffs to universal space (signed) More...
 
GlobalSysSolnType m_solnType
 The solution type of the global system. More...
 
int m_bndSystemBandWidth
 The bandwith of the global bnd system. More...
 
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...
 
std::string m_linSysIterSolver
 Iterative solver: Conjugate Gradient, GMRES. More...
 
Gs::gs_datam_gsh
 
Gs::gs_datam_bndGsh
 
Gs::gs_datam_dirBndGsh
 gs gather communication to impose Dirhichlet BCs. More...
 
int m_staticCondLevel
 The level of recursion in the case of multi-level static condensation. More...
 
int m_numPatches
 The number of patches (~elements) in the current level. More...
 
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
 The number of bnd dofs per patch. More...
 
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
 The number of int dofs per patch. More...
 
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
 Map from the patches of the previous level to the patches of the current level. More...
 
int m_lowestStaticCondLevel
 Lowest static condensation level. More...
 

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 std::shared_ptr< AssemblyMapv_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 57 of file AssemblyMap.h.

Constructor & Destructor Documentation

◆ AssemblyMap() [1/3]

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

Default constructor.

Initialises an empty mapping.

Definition at line 81 of file AssemblyMap.cpp.

86  m_linSysIterSolver("ConjugateGradient"), m_gsh(0), m_bndGsh(0)
87 {
88 }
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:398
std::string m_linSysIterSolver
Iterative solver: Conjugate Gradient, GMRES.
Definition: AssemblyMap.h:415
int m_successiveRHS
sucessive RHS for iterative solver
Definition: AssemblyMap.h:412
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:332
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:341
int m_bndSystemBandWidth
The bandwith of the global bnd system.
Definition: AssemblyMap.h:400
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:345
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:347
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:335
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:343
@ eNoSolnType
No Solution type specified.

◆ AssemblyMap() [2/3]

Nektar::MultiRegions::AssemblyMap::AssemblyMap ( const LibUtilities::SessionReaderSharedPtr pSession,
const LibUtilities::CommSharedPtr comm,
const std::string  variable = "DefaultVar" 
)

Constructor with a communicator.

Definition at line 90 of file AssemblyMap.cpp.

93  : m_session(pSession), m_comm(comm), m_hash(0), m_numLocalBndCoeffs(0),
96  m_linSysIterSolver("ConjugateGradient"), m_gsh(0), m_bndGsh(0)
97 {
98  // Default value from Solver Info
99  m_solnType =
100  pSession->GetSolverInfoAsEnum<GlobalSysSolnType>("GlobalSysSoln");
101  m_preconType =
102  pSession->GetSolverInfoAsEnum<PreconditionerType>("Preconditioner");
103 
104  // Override values with data from GlobalSysSolnInfo section
105  if (pSession->DefinesGlobalSysSolnInfo(variable, "GlobalSysSoln"))
106  {
107  std::string sysSoln =
108  pSession->GetGlobalSysSolnInfo(variable, "GlobalSysSoln");
109  m_solnType = pSession->GetValueAsEnum<GlobalSysSolnType>(
110  "GlobalSysSoln", sysSoln);
111  }
112 
113  if (pSession->DefinesGlobalSysSolnInfo(variable, "Preconditioner"))
114  {
115  std::string precon =
116  pSession->GetGlobalSysSolnInfo(variable, "Preconditioner");
117  m_preconType = pSession->GetValueAsEnum<PreconditionerType>(
118  "Preconditioner", precon);
119  }
120 
121  if (pSession->DefinesGlobalSysSolnInfo(variable,
122  "IterativeSolverTolerance"))
123  {
124  m_iterativeTolerance = boost::lexical_cast<NekDouble>(
125  pSession->GetGlobalSysSolnInfo(variable, "IterativeSolverTolerance")
126  .c_str());
127  }
128  else
129  {
130  pSession->LoadParameter("IterativeSolverTolerance",
133  }
134 
135  if (pSession->DefinesGlobalSysSolnInfo(variable, "MaxIterations"))
136  {
137  m_maxIterations = boost::lexical_cast<int>(
138  pSession->GetGlobalSysSolnInfo(variable, "MaxIterations").c_str());
139  }
140  else
141  {
142  pSession->LoadParameter("MaxIterations", m_maxIterations, 5000);
143  }
144 
145  if (pSession->DefinesGlobalSysSolnInfo(variable, "SuccessiveRHS"))
146  {
147  m_successiveRHS = boost::lexical_cast<int>(
148  pSession->GetGlobalSysSolnInfo(variable, "SuccessiveRHS").c_str());
149  }
150  else
151  {
152  pSession->LoadParameter("SuccessiveRHS", m_successiveRHS, 0);
153  }
154 
155  if (pSession->DefinesGlobalSysSolnInfo(variable, "LinSysIterSolver"))
156  {
158  pSession->GetGlobalSysSolnInfo(variable, "LinSysIterSolver");
159  }
160  else if (pSession->DefinesSolverInfo("LinSysIterSolver"))
161  {
162  m_linSysIterSolver = pSession->GetSolverInfo("LinSysIterSolver");
163  }
164  else
165  {
166  m_linSysIterSolver = "ConjugateGradient";
167  }
168 }
PreconditionerType m_preconType
Type type of preconditioner to use in iterative solver.
Definition: AssemblyMap.h:403
int m_maxIterations
Maximum iterations for iterative solver.
Definition: AssemblyMap.h:406
NekDouble m_iterativeTolerance
Tolerance for iterative solver.
Definition: AssemblyMap.h:409
static const NekDouble kNekIterativeTol

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

◆ AssemblyMap() [3/3]

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 174 of file AssemblyMap.cpp.

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

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_localToLocalBndMap, m_localToLocalIntMap, m_nextLevelLocalToGlobalMap, m_numGlobalBndCoeffs, m_numGlobalCoeffs, m_numGlobalDirBndCoeffs, m_numLocalBndCoeffs, m_numLocalBndCoeffsPerPatch, m_numLocalCoeffs, m_numLocalDirBndCoeffs, m_numLocalIntCoeffsPerPatch, m_numPatches, m_patchMapFromPrevLevel, m_signChange, m_solnType, m_staticCondLevel, Nektar::MultiRegions::RoundNekDoubleToInt(), and sign.

◆ ~AssemblyMap()

Nektar::MultiRegions::AssemblyMap::~AssemblyMap ( void  )
virtual

Destructor.

Definition at line 478 of file AssemblyMap.cpp.

479 {
480 }

Member Function Documentation

◆ Assemble() [1/2]

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

Definition at line 794 of file AssemblyMap.cpp.

796 {
797  v_Assemble(loc, global);
798 }
virtual void v_Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const

References CG_Iterations::loc, and v_Assemble().

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

◆ Assemble() [2/2]

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

Definition at line 800 of file AssemblyMap.cpp.

802 {
803  v_Assemble(loc, global);
804 }

References CG_Iterations::loc, and v_Assemble().

◆ AssembleBnd() [1/4]

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

Definition at line 1239 of file AssemblyMap.cpp.

1241 {
1242  ASSERTL1(loc.size() >= m_numLocalBndCoeffs,
1243  "Local vector is not of correct dimension");
1244  ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1245  "Global vector is not of correct dimension");
1246 
1247  Vmath::Zero(m_numGlobalBndCoeffs, global.get(), 1);
1248 
1249  if (m_signChange)
1250  {
1252  loc.get(), m_localToGlobalBndMap.get(), global.get());
1253  }
1254  else
1255  {
1257  m_localToGlobalBndMap.get(), global.get());
1258  }
1259  UniversalAssembleBnd(global);
1260 }
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
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:862
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:492

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

◆ AssembleBnd() [2/4]

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

Definition at line 1215 of file AssemblyMap.cpp.

1217 {
1218  ASSERTL1(loc.size() >= m_numLocalBndCoeffs,
1219  "Local array is not of correct dimension");
1220  ASSERTL1(global.size() >= m_numGlobalBndCoeffs - offset,
1221  "Global array is not of correct dimension");
1222  Array<OneD, NekDouble> tmp(m_numGlobalBndCoeffs, 0.0);
1223 
1224  if (m_signChange)
1225  {
1227  loc.get(), m_localToGlobalBndMap.get(), tmp.get());
1228  }
1229  else
1230  {
1232  m_localToGlobalBndMap.get(), tmp.get());
1233  }
1234  UniversalAssembleBnd(tmp);
1235  Vmath::Vcopy(m_numGlobalBndCoeffs - offset, tmp.get() + offset, 1,
1236  global.get(), 1);
1237 }
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1255

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

◆ AssembleBnd() [3/4]

void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const

Definition at line 1209 of file AssemblyMap.cpp.

1211 {
1212  AssembleBnd(loc.GetPtr(), global.GetPtr());
1213 }
void AssembleBnd(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const

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

◆ AssembleBnd() [4/4]

void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
int  offset 
) const

Definition at line 1203 of file AssemblyMap.cpp.

1205 {
1206  AssembleBnd(loc.GetPtr(), global.GetPtr(), offset);
1207 }

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

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

◆ AtLastLevel()

bool Nektar::MultiRegions::AssemblyMap::AtLastLevel ( ) const

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

Definition at line 1326 of file AssemblyMap.cpp.

1327 {
1328  return !(m_nextLevelLocalToGlobalMap.get());
1329 }

References m_nextLevelLocalToGlobalMap.

◆ CalculateBndSystemBandWidth()

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 496 of file AssemblyMap.cpp.

497 {
498  int i, j;
499  int cnt = 0;
500  int locSize;
501  int maxId;
502  int minId;
503  int bwidth = -1;
504  for (i = 0; i < m_numPatches; ++i)
505  {
506  locSize = m_numLocalBndCoeffsPerPatch[i];
507  maxId = -1;
508  minId = m_numLocalBndCoeffs + 1;
509  for (j = 0; j < locSize; j++)
510  {
512  {
513  if (m_localToGlobalBndMap[cnt + j] > maxId)
514  {
515  maxId = m_localToGlobalBndMap[cnt + j];
516  }
517 
518  if (m_localToGlobalBndMap[cnt + j] < minId)
519  {
520  minId = m_localToGlobalBndMap[cnt + j];
521  }
522  }
523  }
524  bwidth = (bwidth > (maxId - minId)) ? bwidth : (maxId - minId);
525 
526  cnt += locSize;
527  }
528 
529  m_bndSystemBandWidth = bwidth;
530 }

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

◆ GetBndCondCoeffsToLocalCoeffsMap()

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetBndCondCoeffsToLocalCoeffsMap ( )

Retrieves the local indices corresponding to the boundary expansion modes.

Definition at line 979 of file AssemblyMap.cpp.

980 {
982 }
Array< OneD, int > m_bndCondCoeffsToLocalCoeffsMap
Integer map of bnd cond coeffs to local coefficients.
Definition: AssemblyMap.h:385

References m_bndCondCoeffsToLocalCoeffsMap.

◆ GetBndCondCoeffsToLocalCoeffsSign()

const Array< OneD, NekDouble > & Nektar::MultiRegions::AssemblyMap::GetBndCondCoeffsToLocalCoeffsSign ( )

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

Definition at line 984 of file AssemblyMap.cpp.

985 {
987 }
Array< OneD, NekDouble > m_bndCondCoeffsToLocalCoeffsSign
Integer map of sign of bnd cond coeffs to local coefficients.
Definition: AssemblyMap.h:387

References m_bndCondCoeffsToLocalCoeffsSign.

◆ GetBndCondCoeffsToLocalTraceMap()

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetBndCondCoeffsToLocalTraceMap ( )

Retrieves the local indices corresponding to the boundary expansion modes to global trace.

Definition at line 989 of file AssemblyMap.cpp.

990 {
992 }
Array< OneD, int > m_bndCondCoeffsToLocalTraceMap
Integer map of bnd cond coeff to local trace coeff.
Definition: AssemblyMap.h:389

References m_bndCondCoeffsToLocalTraceMap.

◆ GetBndCondIDToGlobalTraceID() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetBndCondIDToGlobalTraceID ( )

Definition at line 1000 of file AssemblyMap.cpp.

1001 {
1003 }
Array< OneD, int > m_bndCondIDToGlobalTraceID
Integer map of bnd cond trace number to global trace number.
Definition: AssemblyMap.h:391

References m_bndCondIDToGlobalTraceID.

◆ GetBndCondIDToGlobalTraceID() [2/2]

int Nektar::MultiRegions::AssemblyMap::GetBndCondIDToGlobalTraceID ( const int  i)

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

Definition at line 994 of file AssemblyMap.cpp.

995 {
996  ASSERTL1(i < m_bndCondIDToGlobalTraceID.size(), "Index out of range.");
997  return m_bndCondIDToGlobalTraceID[i];
998 }

References ASSERTL1, and m_bndCondIDToGlobalTraceID.

◆ GetBndSystemBandWidth()

int Nektar::MultiRegions::AssemblyMap::GetBndSystemBandWidth ( ) const

Returns the bandwidth of the boundary system.

Definition at line 1289 of file AssemblyMap.cpp.

1290 {
1291  return m_bndSystemBandWidth;
1292 }

References m_bndSystemBandWidth.

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

◆ GetComm()

LibUtilities::CommSharedPtr Nektar::MultiRegions::AssemblyMap::GetComm ( )

Retrieves the communicator.

Definition at line 718 of file AssemblyMap.cpp.

719 {
720  return m_comm;
721 }

References m_comm.

◆ GetExtraDirEdges()

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetExtraDirEdges ( )

Definition at line 926 of file AssemblyMap.cpp.

927 {
928  return v_GetExtraDirEdges();
929 }
virtual const Array< OneD, const int > & v_GetExtraDirEdges()

References v_GetExtraDirEdges().

◆ GetFullSystemBandWidth()

int Nektar::MultiRegions::AssemblyMap::GetFullSystemBandWidth ( ) const

Definition at line 886 of file AssemblyMap.cpp.

887 {
888  return v_GetFullSystemBandWidth();
889 }
virtual int v_GetFullSystemBandWidth() const

References v_GetFullSystemBandWidth().

◆ GetGlobalSysSolnType()

GlobalSysSolnType Nektar::MultiRegions::AssemblyMap::GetGlobalSysSolnType ( ) const

Returns the method of solving global systems.

Definition at line 1331 of file AssemblyMap.cpp.

1332 {
1333  return m_solnType;
1334 }

References m_solnType.

Referenced by AssemblyMap().

◆ GetGlobalToUniversalBndMap()

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalBndMap ( )

Definition at line 957 of file AssemblyMap.cpp.

958 {
960 }

References m_globalToUniversalBndMap.

Referenced by AssemblyMap().

◆ GetGlobalToUniversalBndMapUnique()

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalBndMapUnique ( )

Definition at line 962 of file AssemblyMap.cpp.

963 {
965 }

References m_globalToUniversalBndMapUnique.

Referenced by AssemblyMap().

◆ GetGlobalToUniversalMap() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMap ( )

Definition at line 748 of file AssemblyMap.cpp.

749 {
750  return v_GetGlobalToUniversalMap();
751 }
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMap()

References v_GetGlobalToUniversalMap().

◆ GetGlobalToUniversalMap() [2/2]

int Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMap ( const int  i) const

Definition at line 733 of file AssemblyMap.cpp.

734 {
735  return v_GetGlobalToUniversalMap(i);
736 }

References v_GetGlobalToUniversalMap().

◆ GetGlobalToUniversalMapUnique() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMapUnique ( )

Definition at line 753 of file AssemblyMap.cpp.

754 {
756 }
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMapUnique()

References v_GetGlobalToUniversalMapUnique().

◆ GetGlobalToUniversalMapUnique() [2/2]

int Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalMapUnique ( const int  i) const

Definition at line 738 of file AssemblyMap.cpp.

739 {
741 }

References v_GetGlobalToUniversalMapUnique().

◆ GetHash()

size_t Nektar::MultiRegions::AssemblyMap::GetHash ( ) const

Retrieves the hash of this map.

Definition at line 723 of file AssemblyMap.cpp.

724 {
725  return m_hash;
726 }

References m_hash.

◆ GetIterativeTolerance()

NekDouble Nektar::MultiRegions::AssemblyMap::GetIterativeTolerance ( ) const

Definition at line 1341 of file AssemblyMap.cpp.

1342 {
1343  return m_iterativeTolerance;
1344 }

References m_iterativeTolerance.

◆ GetLinSysIterSolver()

std::string Nektar::MultiRegions::AssemblyMap::GetLinSysIterSolver ( ) const

Definition at line 1356 of file AssemblyMap.cpp.

1357 {
1358  return m_linSysIterSolver;
1359 }

References m_linSysIterSolver.

◆ GetLocalToGlobalBndMap() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndMap ( void  )

Retrieve the global indices of the local boundary modes.

Definition at line 942 of file AssemblyMap.cpp.

943 {
944  return m_localToGlobalBndMap;
945 }

References m_localToGlobalBndMap.

◆ GetLocalToGlobalBndMap() [2/2]

int Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndMap ( const int  i) const

Retrieve the global index of a given local boundary mode.

Definition at line 937 of file AssemblyMap.cpp.

938 {
939  return m_localToGlobalBndMap[i];
940 }

References m_localToGlobalBndMap.

Referenced by AssemblyMap().

◆ GetLocalToGlobalBndSign() [1/2]

Array< OneD, const NekDouble > Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndSign ( void  ) const

Retrieve the sign change for all local boundary modes.

Definition at line 952 of file AssemblyMap.cpp.

953 {
954  return m_localToGlobalBndSign;
955 }

References m_localToGlobalBndSign.

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

◆ GetLocalToGlobalBndSign() [2/2]

NekDouble Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndSign ( const int  i) const

Retrieve the sign change of a given local boundary mode.

Definition at line 967 of file AssemblyMap.cpp.

968 {
969  if (m_signChange)
970  {
971  return m_localToGlobalBndSign[i];
972  }
973  else
974  {
975  return 1.0;
976  }
977 }

References m_localToGlobalBndSign, and m_signChange.

Referenced by AssemblyMap().

◆ GetLocalToGlobalMap() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalMap ( )

Definition at line 743 of file AssemblyMap.cpp.

744 {
745  return v_GetLocalToGlobalMap();
746 }
virtual const Array< OneD, const int > & v_GetLocalToGlobalMap()

References v_GetLocalToGlobalMap().

◆ GetLocalToGlobalMap() [2/2]

int Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalMap ( const int  i) const

Definition at line 728 of file AssemblyMap.cpp.

729 {
730  return v_GetLocalToGlobalMap(i);
731 }

References v_GetLocalToGlobalMap().

◆ GetLocalToGlobalSign() [1/2]

const Array< OneD, NekDouble > & Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalSign ( ) const

Definition at line 763 of file AssemblyMap.cpp.

764 {
765  return v_GetLocalToGlobalSign();
766 }
virtual const Array< OneD, NekDouble > & v_GetLocalToGlobalSign() const

References v_GetLocalToGlobalSign().

◆ GetLocalToGlobalSign() [2/2]

NekDouble Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalSign ( const int  i) const

Definition at line 758 of file AssemblyMap.cpp.

759 {
760  return v_GetLocalToGlobalSign(i);
761 }

References v_GetLocalToGlobalSign().

◆ GetLowestStaticCondLevel()

int Nektar::MultiRegions::AssemblyMap::GetLowestStaticCondLevel ( ) const
inline

Definition at line 313 of file AssemblyMap.h.

314  {
316  }

References m_lowestStaticCondLevel.

◆ GetMaxIterations()

int Nektar::MultiRegions::AssemblyMap::GetMaxIterations ( ) const

Definition at line 1346 of file AssemblyMap.cpp.

1347 {
1348  return m_maxIterations;
1349 }

References m_maxIterations.

◆ GetNextLevelLocalToGlobalMap()

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 1316 of file AssemblyMap.cpp.

1317 {
1319 }

References m_nextLevelLocalToGlobalMap.

◆ GetNumDirEdges()

int Nektar::MultiRegions::AssemblyMap::GetNumDirEdges ( ) const

Definition at line 906 of file AssemblyMap.cpp.

907 {
908  return v_GetNumDirEdges();
909 }
virtual int v_GetNumDirEdges() const

References v_GetNumDirEdges().

◆ GetNumDirFaces()

int Nektar::MultiRegions::AssemblyMap::GetNumDirFaces ( ) const

Definition at line 911 of file AssemblyMap.cpp.

912 {
913  return v_GetNumDirFaces();
914 }
virtual int v_GetNumDirFaces() const

References v_GetNumDirFaces().

◆ GetNumGlobalBndCoeffs()

int Nektar::MultiRegions::AssemblyMap::GetNumGlobalBndCoeffs ( ) const

Returns the total number of global boundary coefficients.

Definition at line 1020 of file AssemblyMap.cpp.

1021 {
1022  return m_numGlobalBndCoeffs;
1023 }

References m_numGlobalBndCoeffs.

Referenced by AssemblyMap().

◆ GetNumGlobalCoeffs()

int Nektar::MultiRegions::AssemblyMap::GetNumGlobalCoeffs ( ) const

Returns the total number of global coefficients.

Definition at line 1030 of file AssemblyMap.cpp.

1031 {
1032  return m_numGlobalCoeffs;
1033 }

References m_numGlobalCoeffs.

◆ GetNumGlobalDirBndCoeffs()

int Nektar::MultiRegions::AssemblyMap::GetNumGlobalDirBndCoeffs ( ) const

Returns the number of global Dirichlet boundary coefficients.

Definition at line 1005 of file AssemblyMap.cpp.

1006 {
1007  return m_numGlobalDirBndCoeffs;
1008 }

References m_numGlobalDirBndCoeffs.

Referenced by AssemblyMap().

◆ GetNumLocalBndCoeffs()

int Nektar::MultiRegions::AssemblyMap::GetNumLocalBndCoeffs ( ) const

Returns the total number of local boundary coefficients.

Definition at line 1015 of file AssemblyMap.cpp.

1016 {
1017  return m_numLocalBndCoeffs;
1018 }

References m_numLocalBndCoeffs.

Referenced by AssemblyMap().

◆ GetNumLocalBndCoeffsPerPatch()

const Array< OneD, const unsigned int > & Nektar::MultiRegions::AssemblyMap::GetNumLocalBndCoeffsPerPatch ( )

Returns the number of local boundary coefficients in each patch.

Definition at line 1305 of file AssemblyMap.cpp.

1306 {
1308 }

References m_numLocalBndCoeffsPerPatch.

Referenced by AssemblyMap().

◆ GetNumLocalCoeffs()

int Nektar::MultiRegions::AssemblyMap::GetNumLocalCoeffs ( ) const

Returns the total number of local coefficients.

Definition at line 1025 of file AssemblyMap.cpp.

1026 {
1027  return m_numLocalCoeffs;
1028 }

References m_numLocalCoeffs.

◆ GetNumLocalDirBndCoeffs()

int Nektar::MultiRegions::AssemblyMap::GetNumLocalDirBndCoeffs ( ) const

Returns the number of local Dirichlet boundary coefficients.

Definition at line 1010 of file AssemblyMap.cpp.

1011 {
1012  return m_numLocalDirBndCoeffs;
1013 }

References m_numLocalDirBndCoeffs.

Referenced by AssemblyMap().

◆ GetNumLocalIntCoeffsPerPatch()

const Array< OneD, const unsigned int > & Nektar::MultiRegions::AssemblyMap::GetNumLocalIntCoeffsPerPatch ( )

Returns the number of local interior coefficients in each patch.

Definition at line 1311 of file AssemblyMap.cpp.

1312 {
1314 }

References m_numLocalIntCoeffsPerPatch.

◆ GetNumNonDirEdgeModes()

int Nektar::MultiRegions::AssemblyMap::GetNumNonDirEdgeModes ( ) const

Definition at line 896 of file AssemblyMap.cpp.

897 {
898  return v_GetNumNonDirEdgeModes();
899 }
virtual int v_GetNumNonDirEdgeModes() const

References v_GetNumNonDirEdgeModes().

◆ GetNumNonDirEdges()

int Nektar::MultiRegions::AssemblyMap::GetNumNonDirEdges ( ) const

Definition at line 916 of file AssemblyMap.cpp.

917 {
918  return v_GetNumNonDirEdges();
919 }
virtual int v_GetNumNonDirEdges() const

References v_GetNumNonDirEdges().

◆ GetNumNonDirFaceModes()

int Nektar::MultiRegions::AssemblyMap::GetNumNonDirFaceModes ( ) const

Definition at line 901 of file AssemblyMap.cpp.

902 {
903  return v_GetNumNonDirFaceModes();
904 }
virtual int v_GetNumNonDirFaceModes() const

References v_GetNumNonDirFaceModes().

◆ GetNumNonDirFaces()

int Nektar::MultiRegions::AssemblyMap::GetNumNonDirFaces ( ) const

Definition at line 921 of file AssemblyMap.cpp.

922 {
923  return v_GetNumNonDirFaces();
924 }
virtual int v_GetNumNonDirFaces() const

References v_GetNumNonDirFaces().

◆ GetNumNonDirVertexModes()

int Nektar::MultiRegions::AssemblyMap::GetNumNonDirVertexModes ( ) const

Definition at line 891 of file AssemblyMap.cpp.

892 {
893  return v_GetNumNonDirVertexModes();
894 }
virtual int v_GetNumNonDirVertexModes() const

References v_GetNumNonDirVertexModes().

◆ GetNumPatches()

int Nektar::MultiRegions::AssemblyMap::GetNumPatches ( ) const

Returns the number of patches in this static condensation level.

Definition at line 1299 of file AssemblyMap.cpp.

1300 {
1301  return m_numPatches;
1302 }

References m_numPatches.

Referenced by AssemblyMap().

◆ GetPatchMapFromPrevLevel()

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 1321 of file AssemblyMap.cpp.

1322 {
1323  return m_patchMapFromPrevLevel;
1324 }

References m_patchMapFromPrevLevel.

◆ GetPreconType()

PreconditionerType Nektar::MultiRegions::AssemblyMap::GetPreconType ( ) const

Definition at line 1336 of file AssemblyMap.cpp.

1337 {
1338  return m_preconType;
1339 }

References m_preconType.

◆ GetSignChange()

bool Nektar::MultiRegions::AssemblyMap::GetSignChange ( )

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

Definition at line 947 of file AssemblyMap.cpp.

948 {
949  return m_signChange;
950 }

References m_signChange.

Referenced by AssemblyMap().

◆ GetSingularSystem()

bool Nektar::MultiRegions::AssemblyMap::GetSingularSystem ( ) const

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

Definition at line 1035 of file AssemblyMap.cpp.

1036 {
1037  return m_systemSingular;
1038 }
bool m_systemSingular
Flag indicating if the system is singular or not.
Definition: AssemblyMap.h:349

References m_systemSingular.

◆ GetStaticCondLevel()

int Nektar::MultiRegions::AssemblyMap::GetStaticCondLevel ( ) const

Returns the level of static condensation for this map.

Definition at line 1294 of file AssemblyMap.cpp.

1295 {
1296  return m_staticCondLevel;
1297 }

References m_staticCondLevel.

Referenced by AssemblyMap().

◆ GetSuccessiveRHS()

int Nektar::MultiRegions::AssemblyMap::GetSuccessiveRHS ( ) const

Definition at line 1351 of file AssemblyMap.cpp.

1352 {
1353  return m_successiveRHS;
1354 }

References m_successiveRHS.

◆ GlobalToLocal() [1/2]

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

Definition at line 782 of file AssemblyMap.cpp.

784 {
785  v_GlobalToLocal(global, loc);
786 }
virtual void v_GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const

References CG_Iterations::loc, and v_GlobalToLocal().

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

◆ GlobalToLocal() [2/2]

void Nektar::MultiRegions::AssemblyMap::GlobalToLocal ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const

Definition at line 788 of file AssemblyMap.cpp.

790 {
791  v_GlobalToLocal(global, loc);
792 }

References CG_Iterations::loc, and v_GlobalToLocal().

◆ GlobalToLocalBnd() [1/4]

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

Definition at line 1078 of file AssemblyMap.cpp.

1080 {
1081  ASSERTL1(loc.size() >= m_numLocalBndCoeffs,
1082  "Local vector is not of correct dimension");
1083  ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1084  "Global vector is not of correct dimension");
1085 
1086  if (m_signChange)
1087  {
1089  global.get(), m_localToGlobalBndMap.get(), loc.get());
1090  }
1091  else
1092  {
1093  Vmath::Gathr(m_numLocalBndCoeffs, global.get(),
1094  m_localToGlobalBndMap.get(), loc.get());
1095  }
1096 }
void Gathr(int n, const T *sign, const T *x, const int *y, T *z)
Gather vector z[i] = sign[i]*x[y[i]].
Definition: Vmath.cpp:805

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

◆ GlobalToLocalBnd() [2/4]

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

Definition at line 1052 of file AssemblyMap.cpp.

1055 {
1056  ASSERTL1(loc.size() >= m_numLocalBndCoeffs,
1057  "Local vector is not of correct dimension");
1058  ASSERTL1(global.size() >= m_numGlobalBndCoeffs - offset,
1059  "Global vector is not of correct dimension");
1060 
1061  // offset input data by length "offset" for Dirichlet boundary conditions.
1062  Array<OneD, NekDouble> tmp(m_numGlobalBndCoeffs, 0.0);
1063  Vmath::Vcopy(m_numGlobalBndCoeffs - offset, global.get(), 1,
1064  tmp.get() + offset, 1);
1065 
1066  if (m_signChange)
1067  {
1069  tmp.get(), m_localToGlobalBndMap.get(), loc.get());
1070  }
1071  else
1072  {
1073  Vmath::Gathr(m_numLocalBndCoeffs, tmp.get(),
1074  m_localToGlobalBndMap.get(), loc.get());
1075  }
1076 }

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

◆ GlobalToLocalBnd() [3/4]

void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const

Definition at line 1046 of file AssemblyMap.cpp.

1048 {
1049  GlobalToLocalBnd(global.GetPtr(), loc.GetPtr());
1050 }
void GlobalToLocalBnd(const NekVector< NekDouble > &global, NekVector< NekDouble > &loc, int offset) const

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

◆ GlobalToLocalBnd() [4/4]

void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc,
int  offset 
) const

Definition at line 1040 of file AssemblyMap.cpp.

1042 {
1043  GlobalToLocalBnd(global.GetPtr(), loc.GetPtr(), offset);
1044 }

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

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

◆ GlobalToLocalBndWithoutSign()

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

Definition at line 1361 of file AssemblyMap.cpp.

1363 {
1364  ASSERTL1(loc.size() >= m_numLocalBndCoeffs,
1365  "Local vector is not of correct dimension");
1366  ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1367  "Global vector is not of correct dimension");
1368 
1370  loc.get());
1371 }

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

Referenced by AssemblyMap().

◆ LinearSpaceMap()

std::shared_ptr< AssemblyMap > Nektar::MultiRegions::AssemblyMap::LinearSpaceMap ( const ExpList locexp,
GlobalSysSolnType  solnType 
)

Definition at line 931 of file AssemblyMap.cpp.

933 {
934  return v_LinearSpaceMap(locexp, solnType);
935 }
virtual std::shared_ptr< AssemblyMap > v_LinearSpaceMap(const ExpList &locexp, GlobalSysSolnType solnType)
Generate a linear space mapping from existing mapping.

References v_LinearSpaceMap().

◆ LocalBndToGlobal() [1/2]

void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
bool  UseComm = true 
) const

Definition at line 1130 of file AssemblyMap.cpp.

1133 {
1134  ASSERTL1(loc.size() >= m_numLocalBndCoeffs,
1135  "Local vector is not of correct dimension");
1136  ASSERTL1(global.size() >= m_numGlobalBndCoeffs,
1137  "Global vector is not of correct dimension");
1138 
1139  if (m_signChange)
1140  {
1142  loc.get(), m_localToGlobalBndMap.get(), global.get());
1143  }
1144  else
1145  {
1147  m_localToGlobalBndMap.get(), global.get());
1148  }
1149  if (UseComm)
1150  {
1151  Gs::Gather(global, Gs::gs_max, m_bndGsh);
1152  }
1153 }
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:270
@ gs_max
Definition: GsLib.hpp:65
void Scatr(int n, const T *x, const int *y, T *z)
Scatter vector z[y[i]] = x[i].
Definition: Vmath.cpp:822

References ASSERTL1, Gs::Gather(), Gs::gs_max, CG_Iterations::loc, m_bndGsh, m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, m_signChange, and Vmath::Scatr().

◆ LocalBndToGlobal() [2/2]

void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
int  offset,
bool  UseComm = true 
) const

Definition at line 1098 of file AssemblyMap.cpp.

1101 {
1102  ASSERTL1(loc.size() >= m_numLocalBndCoeffs,
1103  "Local vector is not of correct dimension");
1104  ASSERTL1(global.size() >= m_numGlobalBndCoeffs - offset,
1105  "Global vector is not of correct dimension");
1106 
1107  // offset input data by length "offset" for Dirichlet boundary conditions.
1108  Array<OneD, NekDouble> tmp(m_numGlobalBndCoeffs, 0.0);
1109 
1110  if (m_signChange)
1111  {
1113  loc.get(), m_localToGlobalBndMap.get(), tmp.get());
1114  }
1115  else
1116  {
1118  m_localToGlobalBndMap.get(), tmp.get());
1119  }
1120 
1121  // Ensure each processor has unique value with a max gather.
1122  if (UseComm)
1123  {
1125  }
1126  Vmath::Vcopy(m_numGlobalBndCoeffs - offset, tmp.get() + offset, 1,
1127  global.get(), 1);
1128 }

References ASSERTL1, Gs::Gather(), Gs::gs_max, CG_Iterations::loc, m_bndGsh, m_localToGlobalBndMap, m_localToGlobalBndSign, m_numGlobalBndCoeffs, m_numLocalBndCoeffs, m_signChange, Vmath::Scatr(), and Vmath::Vcopy().

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

◆ LocalBndToLocal()

void Nektar::MultiRegions::AssemblyMap::LocalBndToLocal ( const Array< OneD, const NekDouble > &  locbnd,
Array< OneD, NekDouble > &  local 
) const

Definition at line 1179 of file AssemblyMap.cpp.

1181 {
1182  ASSERTL1(locbnd.size() >= m_numLocalBndCoeffs,
1183  "LocBnd vector is not of correct dimension");
1184  ASSERTL1(local.size() >= m_numLocalCoeffs,
1185  "Local vector is not of correct dimension");
1186 
1188  local.get());
1189 }

References ASSERTL1, m_localToLocalBndMap, m_numLocalBndCoeffs, m_numLocalCoeffs, and Vmath::Scatr().

◆ LocalIntToLocal()

void Nektar::MultiRegions::AssemblyMap::LocalIntToLocal ( const Array< OneD, const NekDouble > &  locbnd,
Array< OneD, NekDouble > &  local 
) const

Definition at line 1191 of file AssemblyMap.cpp.

1193 {
1194  ASSERTL1(locint.size() >= m_numLocalCoeffs - m_numLocalBndCoeffs,
1195  "LocBnd vector is not of correct dimension");
1196  ASSERTL1(local.size() >= m_numLocalCoeffs,
1197  "Local vector is not of correct dimension");
1198 
1200  m_localToLocalIntMap.get(), local.get());
1201 }

References ASSERTL1, m_localToLocalIntMap, m_numLocalBndCoeffs, m_numLocalCoeffs, and Vmath::Scatr().

◆ LocalToGlobal() [1/2]

void Nektar::MultiRegions::AssemblyMap::LocalToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
bool  useComm = true 
) const

Definition at line 768 of file AssemblyMap.cpp.

771 {
772  v_LocalToGlobal(loc, global, useComm);
773 }
virtual void v_LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global, bool useComm) const

References CG_Iterations::loc, and v_LocalToGlobal().

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

◆ LocalToGlobal() [2/2]

void Nektar::MultiRegions::AssemblyMap::LocalToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
bool  useComm = true 
) const

Definition at line 775 of file AssemblyMap.cpp.

778 {
779  v_LocalToGlobal(loc, global, useComm);
780 }

References CG_Iterations::loc, and v_LocalToGlobal().

◆ LocalToLocalBnd()

void Nektar::MultiRegions::AssemblyMap::LocalToLocalBnd ( const Array< OneD, const NekDouble > &  local,
Array< OneD, NekDouble > &  locbnd 
) const

Definition at line 1155 of file AssemblyMap.cpp.

1157 {
1158  ASSERTL1(locbnd.size() >= m_numLocalBndCoeffs,
1159  "LocBnd vector is not of correct dimension");
1160  ASSERTL1(local.size() >= m_numLocalCoeffs,
1161  "Local vector is not of correct dimension");
1162 
1164  locbnd.get());
1165 }

References ASSERTL1, Vmath::Gathr(), m_localToLocalBndMap, m_numLocalBndCoeffs, and m_numLocalCoeffs.

◆ LocalToLocalInt()

void Nektar::MultiRegions::AssemblyMap::LocalToLocalInt ( const Array< OneD, const NekDouble > &  local,
Array< OneD, NekDouble > &  locint 
) const

Definition at line 1167 of file AssemblyMap.cpp.

1169 {
1170  ASSERTL1(locint.size() >= m_numLocalCoeffs - m_numLocalBndCoeffs,
1171  "Locint vector is not of correct dimension");
1172  ASSERTL1(local.size() >= m_numLocalCoeffs,
1173  "Local vector is not of correct dimension");
1174 
1176  m_localToLocalIntMap.get(), locint.get());
1177 }

References ASSERTL1, Vmath::Gathr(), m_localToLocalIntMap, m_numLocalBndCoeffs, and m_numLocalCoeffs.

◆ PatchAssemble()

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

Definition at line 862 of file AssemblyMap.cpp.

864 {
865  Array<OneD, const NekDouble> local;
866  Array<OneD, const int> map = m_patchMapFromPrevLevel->GetNewLevelMap();
867  Array<OneD, const NekDouble> sign = m_patchMapFromPrevLevel->GetSign();
868 
869  if (global.data() == loc.data())
870  {
871  local = Array<OneD, NekDouble>(map.size(), loc.data());
872  }
873  else
874  {
875  local = loc; // create reference
876  }
877 
878  // since we are calling mapping from level down from array
879  // the m_numLocaBndCoeffs represents the size of the
880  // boundary elements we need to assemble into
881  Vmath::Zero(m_numLocalCoeffs, global.get(), 1);
882 
883  Vmath::Assmb(map.size(), sign.get(), local.get(), map.get(), global.get());
884 }

References Vmath::Assmb(), CG_Iterations::loc, m_numLocalCoeffs, m_patchMapFromPrevLevel, sign, and Vmath::Zero().

◆ PatchGlobalToLocal()

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

Definition at line 842 of file AssemblyMap.cpp.

844 {
845  Array<OneD, const NekDouble> glo;
846 
847  Array<OneD, const int> map = m_patchMapFromPrevLevel->GetNewLevelMap();
848  Array<OneD, const NekDouble> sign = m_patchMapFromPrevLevel->GetSign();
849 
850  if (global.data() == loc.data())
851  {
852  glo = Array<OneD, NekDouble>(global.size(), global.data());
853  }
854  else
855  {
856  glo = global; // create reference
857  }
858 
859  Vmath::Gathr(map.size(), sign.get(), glo.get(), map.get(), loc.get());
860 }

References Vmath::Gathr(), CG_Iterations::loc, m_patchMapFromPrevLevel, and sign.

◆ PatchLocalToGlobal()

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

Definition at line 822 of file AssemblyMap.cpp.

824 {
825  Array<OneD, const NekDouble> local;
826 
827  Array<OneD, const int> map = m_patchMapFromPrevLevel->GetNewLevelMap();
828  Array<OneD, const NekDouble> sign = m_patchMapFromPrevLevel->GetSign();
829 
830  if (global.data() == loc.data())
831  {
832  local = Array<OneD, NekDouble>(map.size(), loc.data());
833  }
834  else
835  {
836  local = loc; // create reference
837  }
838 
839  Vmath::Scatr(map.size(), sign.get(), local.get(), map.get(), global.get());
840 }

References CG_Iterations::loc, m_patchMapFromPrevLevel, Vmath::Scatr(), and sign.

◆ PrintStats()

void Nektar::MultiRegions::AssemblyMap::PrintStats ( std::ostream &  out,
std::string  variable,
bool  printHeader = true 
) const

Definition at line 1373 of file AssemblyMap.cpp.

1375 {
1376  LibUtilities::CommSharedPtr vRowComm = m_session->GetComm()->GetRowComm();
1377  bool isRoot = vRowComm->GetRank() == 0;
1378  int n = vRowComm->GetSize();
1379  int i;
1380 
1381  // Determine number of global degrees of freedom.
1382  int globBndCnt = 0, globDirCnt = 0;
1383 
1384  for (i = 0; i < m_numGlobalBndCoeffs; ++i)
1385  {
1387  {
1388  globBndCnt++;
1389 
1390  if (i < m_numGlobalDirBndCoeffs)
1391  {
1392  globDirCnt++;
1393  }
1394  }
1395  }
1396 
1397  int globCnt = m_numGlobalCoeffs - m_numGlobalBndCoeffs + globBndCnt;
1398 
1399  // Calculate maximum valency
1400  Array<OneD, NekDouble> tmpLoc(m_numLocalBndCoeffs, 1.0);
1401  Array<OneD, NekDouble> tmpGlob(m_numGlobalBndCoeffs, 0.0);
1402 
1404  tmpGlob.get());
1405  UniversalAssembleBnd(tmpGlob);
1406 
1407  int totGlobDof = globCnt;
1408  int totGlobBndDof = globBndCnt;
1409  int totGlobDirDof = globDirCnt;
1410  int totLocalDof = m_numLocalCoeffs;
1411  int totLocalBndDof = m_numLocalBndCoeffs;
1412  int totLocalDirDof = m_numLocalDirBndCoeffs;
1413 
1414  int meanValence = 0;
1415  int maxValence = 0;
1416  int minValence = 10000000;
1417  for (int i = 0; i < m_numGlobalBndCoeffs; ++i)
1418  {
1420  {
1421  continue;
1422  }
1423 
1424  if (tmpGlob[i] > maxValence)
1425  {
1426  maxValence = tmpGlob[i];
1427  }
1428  if (tmpGlob[i] < minValence)
1429  {
1430  minValence = tmpGlob[i];
1431  }
1432  meanValence += tmpGlob[i];
1433  }
1434 
1435  vRowComm->AllReduce(maxValence, LibUtilities::ReduceMax);
1436  vRowComm->AllReduce(minValence, LibUtilities::ReduceMin);
1437  vRowComm->AllReduce(meanValence, LibUtilities::ReduceSum);
1438  vRowComm->AllReduce(totGlobDof, LibUtilities::ReduceSum);
1439  vRowComm->AllReduce(totGlobBndDof, LibUtilities::ReduceSum);
1440  vRowComm->AllReduce(totGlobDirDof, LibUtilities::ReduceSum);
1441  vRowComm->AllReduce(totLocalDof, LibUtilities::ReduceSum);
1442  vRowComm->AllReduce(totLocalBndDof, LibUtilities::ReduceSum);
1443  vRowComm->AllReduce(totLocalDirDof, LibUtilities::ReduceSum);
1444 
1445  meanValence /= totGlobBndDof;
1446 
1447  if (isRoot)
1448  {
1449  if (printHeader)
1450  {
1451  out << "Assembly map statistics for field " << variable << ":"
1452  << endl;
1453  }
1454 
1455  out << " - Number of local/global dof : " << totLocalDof
1456  << " " << totGlobDof << endl;
1457  out << " - Number of local/global boundary dof : " << totLocalBndDof
1458  << " " << totGlobBndDof << endl;
1459  out << " - Number of local/global Dirichlet dof : " << totLocalDirDof
1460  << " " << totGlobDirDof << endl;
1461  out << " - dof valency (min/max/mean) : " << minValence
1462  << " " << maxValence << " " << meanValence << endl;
1463 
1464  if (n > 1)
1465  {
1466  NekDouble mean = m_numLocalCoeffs, mean2 = mean * mean;
1467  NekDouble minval = mean, maxval = mean;
1468  Array<OneD, NekDouble> tmp(1);
1469 
1470  for (i = 1; i < n; ++i)
1471  {
1472  vRowComm->Recv(i, tmp);
1473  mean += tmp[0];
1474  mean2 += tmp[0] * tmp[0];
1475 
1476  if (tmp[0] > maxval)
1477  {
1478  maxval = tmp[0];
1479  }
1480  if (tmp[0] < minval)
1481  {
1482  minval = tmp[0];
1483  }
1484  }
1485 
1486  if (maxval > 0.1)
1487  {
1488  out << " - Local dof dist. (min/max/mean/dev) : " << minval
1489  << " " << maxval << " " << (mean / n) << " "
1490  << sqrt(mean2 / n - mean * mean / n / n) << endl;
1491  }
1492 
1493  vRowComm->Block();
1494 
1495  mean = minval = maxval = m_numLocalBndCoeffs;
1496  mean2 = mean * mean;
1497 
1498  for (i = 1; i < n; ++i)
1499  {
1500  vRowComm->Recv(i, tmp);
1501  mean += tmp[0];
1502  mean2 += tmp[0] * tmp[0];
1503 
1504  if (tmp[0] > maxval)
1505  {
1506  maxval = tmp[0];
1507  }
1508  if (tmp[0] < minval)
1509  {
1510  minval = tmp[0];
1511  }
1512  }
1513 
1514  out << " - Local bnd dof dist. (min/max/mean/dev) : " << minval
1515  << " " << maxval << " " << (mean / n) << " "
1516  << sqrt(mean2 / n - mean * mean / n / n) << endl;
1517  }
1518  }
1519  else
1520  {
1521  Array<OneD, NekDouble> tmp(1);
1522  tmp[0] = m_numLocalCoeffs;
1523  vRowComm->Send(0, tmp);
1524  vRowComm->Block();
1525  tmp[0] = m_numLocalBndCoeffs;
1526  vRowComm->Send(0, tmp);
1527  }
1528 
1529  // Either we have no more levels in the static condensation, or we
1530  // are not multi-level.
1532  {
1533  return;
1534  }
1535 
1536  int level = 2;
1538  while (tmp->m_nextLevelLocalToGlobalMap)
1539  {
1540  tmp = tmp->m_nextLevelLocalToGlobalMap;
1541  ++level;
1542  }
1543 
1544  // Print out multi-level static condensation information.
1545  if (n > 1)
1546  {
1547  if (isRoot)
1548  {
1549  NekDouble mean = level, mean2 = mean * mean;
1550  int minval = level, maxval = level;
1551 
1552  Array<OneD, NekDouble> tmpRecv(1);
1553  for (i = 1; i < n; ++i)
1554  {
1555  vRowComm->Recv(i, tmpRecv);
1556  mean += tmpRecv[0];
1557  mean2 += tmpRecv[0] * tmpRecv[0];
1558 
1559  if (tmpRecv[0] > maxval)
1560  {
1561  maxval = (int)(tmpRecv[0] + 0.5);
1562  }
1563  if (tmpRecv[0] < minval)
1564  {
1565  minval = (int)(tmpRecv[0] + 0.5);
1566  }
1567  }
1568 
1569  out << " - M-level sc. dist. (min/max/mean/dev) : " << minval
1570  << " " << maxval << " " << (mean / n) << " "
1571  << sqrt(mean2 / n - mean * mean / n / n) << endl;
1572  }
1573  else
1574  {
1575  Array<OneD, NekDouble> tmpSend(1);
1576  tmpSend[0] = level;
1577  vRowComm->Send(0, tmpSend);
1578  }
1579  }
1580  else
1581  {
1582  out << " - Number of static cond. levels : " << level << endl;
1583  }
1584 
1585  if (isRoot)
1586  {
1587  out << "Stats at lowest static cond. level:" << endl;
1588  }
1589  tmp->PrintStats(out, variable, false);
1590 }
std::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:54
std::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
Definition: AssemblyMap.h:51
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:291

References Vmath::Assmb(), CellMLToNektar.pycml::level, 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, tinysimd::sqrt(), and UniversalAssembleBnd().

◆ SetNextLevelLocalToGlobalMap()

void Nektar::MultiRegions::AssemblyMap::SetNextLevelLocalToGlobalMap ( AssemblyMapSharedPtr  pNextLevelLocalToGlobalMap)

◆ UniversalAbsMaxBnd()

void Nektar::MultiRegions::AssemblyMap::UniversalAbsMaxBnd ( Array< OneD, NekDouble > &  bndvals)

Definition at line 1284 of file AssemblyMap.cpp.

1285 {
1286  Gs::Gather(bndvals, Gs::gs_amax, m_dirBndGsh);
1287 }
Gs::gs_data * m_dirBndGsh
gs gather communication to impose Dirhichlet BCs.
Definition: AssemblyMap.h:420
@ gs_amax
Definition: GsLib.hpp:66

References Gs::Gather(), Gs::gs_amax, and m_dirBndGsh.

◆ UniversalAssemble() [1/3]

void Nektar::MultiRegions::AssemblyMap::UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal) const

◆ UniversalAssemble() [2/3]

void Nektar::MultiRegions::AssemblyMap::UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal,
int  offset 
) const

Definition at line 816 of file AssemblyMap.cpp.

818 {
819  v_UniversalAssemble(pGlobal, offset);
820 }

References v_UniversalAssemble().

◆ UniversalAssemble() [3/3]

void Nektar::MultiRegions::AssemblyMap::UniversalAssemble ( NekVector< NekDouble > &  pGlobal) const

Definition at line 811 of file AssemblyMap.cpp.

812 {
813  v_UniversalAssemble(pGlobal);
814 }

References v_UniversalAssemble().

◆ UniversalAssembleBnd() [1/3]

void Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd ( Array< OneD, NekDouble > &  pGlobal) const

Definition at line 1262 of file AssemblyMap.cpp.

1263 {
1264  ASSERTL1(pGlobal.size() >= m_numGlobalBndCoeffs, "Wrong size.");
1265  Gs::Gather(pGlobal, Gs::gs_add, m_bndGsh);
1266 }
@ gs_add
Definition: GsLib.hpp:62

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

Referenced by AssembleBnd(), PrintStats(), and UniversalAssembleBnd().

◆ UniversalAssembleBnd() [2/3]

void Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd ( Array< OneD, NekDouble > &  pGlobal,
int  offset 
) const

Definition at line 1273 of file AssemblyMap.cpp.

1275 {
1276  Array<OneD, NekDouble> tmp(offset);
1277  if (offset > 0)
1278  Vmath::Vcopy(offset, pGlobal, 1, tmp, 1);
1279  UniversalAssembleBnd(pGlobal);
1280  if (offset > 0)
1281  Vmath::Vcopy(offset, tmp, 1, pGlobal, 1);
1282 }

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

◆ UniversalAssembleBnd() [3/3]

void Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd ( NekVector< NekDouble > &  pGlobal) const

Definition at line 1268 of file AssemblyMap.cpp.

1269 {
1270  UniversalAssembleBnd(pGlobal.GetPtr());
1271 }

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

◆ v_Assemble() [1/2]

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

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

Definition at line 618 of file AssemblyMap.cpp.

620 {
621  boost::ignore_unused(loc, global);
622  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
623 }
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:209

References Nektar::ErrorUtil::efatal, CG_Iterations::loc, and NEKERROR.

Referenced by Assemble().

◆ v_Assemble() [2/2]

void Nektar::MultiRegions::AssemblyMap::v_Assemble ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const
privatevirtual

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

Definition at line 625 of file AssemblyMap.cpp.

627 {
628  boost::ignore_unused(loc, global);
629  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
630 }

References Nektar::ErrorUtil::efatal, CG_Iterations::loc, and NEKERROR.

◆ v_GetExtraDirEdges()

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetExtraDirEdges ( )
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 702 of file AssemblyMap.cpp.

703 {
704  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
705  static Array<OneD, const int> result;
706  return result;
707 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetExtraDirEdges().

◆ v_GetFullSystemBandWidth()

int Nektar::MultiRegions::AssemblyMap::v_GetFullSystemBandWidth ( ) const
privatevirtual

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

Definition at line 654 of file AssemblyMap.cpp.

655 {
656  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
657  return 0;
658 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetFullSystemBandWidth().

◆ v_GetGlobalToUniversalMap() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMap ( void  )
privatevirtual

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

Definition at line 560 of file AssemblyMap.cpp.

561 {
562  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
563  static Array<OneD, const int> result;
564  return result;
565 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetGlobalToUniversalMap().

◆ v_GetGlobalToUniversalMap() [2/2]

int Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMap ( const int  i) const
privatevirtual

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

Definition at line 539 of file AssemblyMap.cpp.

540 {
541  boost::ignore_unused(i);
542  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
543  return 0;
544 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

◆ v_GetGlobalToUniversalMapUnique() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMapUnique ( void  )
privatevirtual

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

Definition at line 567 of file AssemblyMap.cpp.

568 {
569  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
570  static Array<OneD, const int> result;
571  return result;
572 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetGlobalToUniversalMapUnique().

◆ v_GetGlobalToUniversalMapUnique() [2/2]

int Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMapUnique ( const int  i) const
privatevirtual

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

Definition at line 546 of file AssemblyMap.cpp.

547 {
548  boost::ignore_unused(i);
549  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
550  return 0;
551 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

◆ v_GetLocalToGlobalMap() [1/2]

const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalMap ( void  )
privatevirtual

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

Definition at line 553 of file AssemblyMap.cpp.

554 {
555  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
556  static Array<OneD, const int> result;
557  return result;
558 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetLocalToGlobalMap().

◆ v_GetLocalToGlobalMap() [2/2]

int Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalMap ( const int  i) const
privatevirtual

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

Definition at line 532 of file AssemblyMap.cpp.

533 {
534  boost::ignore_unused(i);
535  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
536  return 0;
537 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

◆ v_GetLocalToGlobalSign() [1/2]

const Array< OneD, NekDouble > & Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalSign ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 581 of file AssemblyMap.cpp.

582 {
583  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
584  static Array<OneD, NekDouble> result;
585  return result;
586 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetLocalToGlobalSign().

◆ v_GetLocalToGlobalSign() [2/2]

NekDouble Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalSign ( const int  i) const
privatevirtual

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

Definition at line 574 of file AssemblyMap.cpp.

575 {
576  boost::ignore_unused(i);
577  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
578  return 0.0;
579 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

◆ v_GetNumDirEdges()

int Nektar::MultiRegions::AssemblyMap::v_GetNumDirEdges ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 678 of file AssemblyMap.cpp.

679 {
680  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
681  return 0;
682 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumDirEdges().

◆ v_GetNumDirFaces()

int Nektar::MultiRegions::AssemblyMap::v_GetNumDirFaces ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 684 of file AssemblyMap.cpp.

685 {
686  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
687  return 0;
688 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumDirFaces().

◆ v_GetNumNonDirEdgeModes()

int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirEdgeModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 666 of file AssemblyMap.cpp.

667 {
668  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
669  return 0;
670 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumNonDirEdgeModes().

◆ v_GetNumNonDirEdges()

int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirEdges ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 690 of file AssemblyMap.cpp.

691 {
692  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
693  return 0;
694 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumNonDirEdges().

◆ v_GetNumNonDirFaceModes()

int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirFaceModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 672 of file AssemblyMap.cpp.

673 {
674  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
675  return 0;
676 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumNonDirFaceModes().

◆ v_GetNumNonDirFaces()

int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirFaces ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 696 of file AssemblyMap.cpp.

697 {
698  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
699  return 0;
700 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumNonDirFaces().

◆ v_GetNumNonDirVertexModes()

int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirVertexModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 660 of file AssemblyMap.cpp.

661 {
662  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
663  return 0;
664 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetNumNonDirVertexModes().

◆ v_GlobalToLocal() [1/2]

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

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

Definition at line 604 of file AssemblyMap.cpp.

606 {
607  boost::ignore_unused(loc, global);
608  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
609 }

References Nektar::ErrorUtil::efatal, CG_Iterations::loc, and NEKERROR.

Referenced by GlobalToLocal().

◆ v_GlobalToLocal() [2/2]

void Nektar::MultiRegions::AssemblyMap::v_GlobalToLocal ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const
privatevirtual

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

Definition at line 611 of file AssemblyMap.cpp.

613 {
614  boost::ignore_unused(loc, global);
615  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
616 }

References Nektar::ErrorUtil::efatal, CG_Iterations::loc, and NEKERROR.

◆ v_LinearSpaceMap()

std::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 709 of file AssemblyMap.cpp.

711 {
712  boost::ignore_unused(locexp, solnType);
713  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
714  static std::shared_ptr<AssemblyMap> result;
715  return result;
716 }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by LinearSpaceMap().

◆ v_LocalToGlobal() [1/2]

void Nektar::MultiRegions::AssemblyMap::v_LocalToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
bool  useComm 
) const
privatevirtual

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

Definition at line 588 of file AssemblyMap.cpp.

591 {
592  boost::ignore_unused(loc, global, useComm);
593  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
594 }

References Nektar::ErrorUtil::efatal, CG_Iterations::loc, and NEKERROR.

Referenced by LocalToGlobal().

◆ v_LocalToGlobal() [2/2]

void Nektar::MultiRegions::AssemblyMap::v_LocalToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
bool  useComm 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 596 of file AssemblyMap.cpp.

599 {
600  boost::ignore_unused(loc, global, useComm);
601  NEKERROR(ErrorUtil::efatal, "Not defined for this type of mapping.");
602 }

References Nektar::ErrorUtil::efatal, CG_Iterations::loc, and NEKERROR.

◆ v_UniversalAssemble() [1/3]

void Nektar::MultiRegions::AssemblyMap::v_UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal) const
privatevirtual

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

Definition at line 632 of file AssemblyMap.cpp.

633 {
634  boost::ignore_unused(pGlobal);
635  // Do nothing here since multi-level static condensation uses a
636  // AssemblyMap and thus will call this routine in serial.
637 }

Referenced by UniversalAssemble().

◆ v_UniversalAssemble() [2/3]

void Nektar::MultiRegions::AssemblyMap::v_UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal,
int  offset 
) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 646 of file AssemblyMap.cpp.

648 {
649  boost::ignore_unused(pGlobal, offset);
650  // Do nothing here since multi-level static condensation uses a
651  // AssemblyMap and thus will call this routine in serial.
652 }

◆ v_UniversalAssemble() [3/3]

void Nektar::MultiRegions::AssemblyMap::v_UniversalAssemble ( NekVector< NekDouble > &  pGlobal) const
privatevirtual

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

Definition at line 639 of file AssemblyMap.cpp.

640 {
641  boost::ignore_unused(pGlobal);
642  // Do nothing here since multi-level static condensation uses a
643  // AssemblyMap and thus will call this routine in serial.
644 }

Member Data Documentation

◆ m_bndCondCoeffsToLocalCoeffsMap

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

Integer map of bnd cond coeffs to local coefficients.

Definition at line 385 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), and GetBndCondCoeffsToLocalCoeffsMap().

◆ m_bndCondCoeffsToLocalCoeffsSign

Array<OneD, NekDouble> Nektar::MultiRegions::AssemblyMap::m_bndCondCoeffsToLocalCoeffsSign
protected

Integer map of sign of bnd cond coeffs to local coefficients.

Definition at line 387 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), and GetBndCondCoeffsToLocalCoeffsSign().

◆ m_bndCondCoeffsToLocalTraceMap

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

Integer map of bnd cond coeff to local trace coeff.

Definition at line 389 of file AssemblyMap.h.

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

◆ m_bndCondIDToGlobalTraceID

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

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

Definition at line 391 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), and GetBndCondIDToGlobalTraceID().

◆ m_bndGsh

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

◆ m_bndSystemBandWidth

int Nektar::MultiRegions::AssemblyMap::m_bndSystemBandWidth
protected

The bandwith of the global bnd system.

Definition at line 400 of file AssemblyMap.h.

Referenced by CalculateBndSystemBandWidth(), and GetBndSystemBandWidth().

◆ m_comm

LibUtilities::CommSharedPtr Nektar::MultiRegions::AssemblyMap::m_comm
protected

◆ m_dirBndGsh

Gs::gs_data* Nektar::MultiRegions::AssemblyMap::m_dirBndGsh
protected

gs gather communication to impose Dirhichlet BCs.

Definition at line 420 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), and UniversalAbsMaxBnd().

◆ m_globalToUniversalBndMap

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

◆ m_globalToUniversalBndMapUnique

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

◆ m_gsh

Gs::gs_data* Nektar::MultiRegions::AssemblyMap::m_gsh
protected

◆ m_hash

size_t Nektar::MultiRegions::AssemblyMap::m_hash
protected

◆ m_iterativeTolerance

NekDouble Nektar::MultiRegions::AssemblyMap::m_iterativeTolerance
protected

Tolerance for iterative solver.

Definition at line 409 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetIterativeTolerance().

◆ m_linSysIterSolver

std::string Nektar::MultiRegions::AssemblyMap::m_linSysIterSolver
protected

Iterative solver: Conjugate Gradient, GMRES.

Definition at line 415 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetLinSysIterSolver().

◆ m_localToGlobalBndMap

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

◆ m_localToGlobalBndSign

Array<OneD, NekDouble> Nektar::MultiRegions::AssemblyMap::m_localToGlobalBndSign
protected

◆ m_localToLocalBndMap

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

◆ m_localToLocalIntMap

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

◆ m_lowestStaticCondLevel

int Nektar::MultiRegions::AssemblyMap::m_lowestStaticCondLevel
protected

◆ m_maxIterations

int Nektar::MultiRegions::AssemblyMap::m_maxIterations
protected

Maximum iterations for iterative solver.

Definition at line 406 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetMaxIterations().

◆ m_nextLevelLocalToGlobalMap

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 435 of file AssemblyMap.h.

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

◆ m_numGlobalBndCoeffs

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

◆ m_numGlobalCoeffs

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 371 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(), Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal(), and Nektar::MultiRegions::AssemblyMapCG::v_LinearSpaceMap().

◆ m_numGlobalDirBndCoeffs

int Nektar::MultiRegions::AssemblyMap::m_numGlobalDirBndCoeffs
protected

◆ m_numLocalBndCoeffs

int Nektar::MultiRegions::AssemblyMap::m_numLocalBndCoeffs
protected

◆ m_numLocalBndCoeffsPerPatch

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

◆ m_numLocalCoeffs

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 360 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapCG::AssemblyMapCG(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth(), Nektar::CoupledAssemblyMap::CoupledAssemblyMap(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumLocalCoeffs(), LocalBndToLocal(), LocalIntToLocal(), LocalToLocalBnd(), LocalToLocalInt(), PatchAssemble(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG::v_Assemble(), Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal(), and Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal().

◆ m_numLocalDirBndCoeffs

int Nektar::MultiRegions::AssemblyMap::m_numLocalDirBndCoeffs
protected

◆ m_numLocalIntCoeffsPerPatch

Array<OneD, unsigned int> Nektar::MultiRegions::AssemblyMap::m_numLocalIntCoeffsPerPatch
protected

◆ m_numPatches

int Nektar::MultiRegions::AssemblyMap::m_numPatches
protected

◆ m_patchMapFromPrevLevel

PatchMapSharedPtr Nektar::MultiRegions::AssemblyMap::m_patchMapFromPrevLevel
private

Mapping information for previous level in MultiLevel Solver.

Definition at line 447 of file AssemblyMap.h.

Referenced by AssemblyMap(), GetPatchMapFromPrevLevel(), PatchAssemble(), PatchGlobalToLocal(), and PatchLocalToGlobal().

◆ m_preconType

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

Type type of preconditioner to use in iterative solver.

Definition at line 403 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetPreconType().

◆ m_session

LibUtilities::SessionReaderSharedPtr Nektar::MultiRegions::AssemblyMap::m_session
protected

◆ m_signChange

bool Nektar::MultiRegions::AssemblyMap::m_signChange
protected

◆ m_solnType

GlobalSysSolnType Nektar::MultiRegions::AssemblyMap::m_solnType
protected

◆ m_staticCondLevel

int Nektar::MultiRegions::AssemblyMap::m_staticCondLevel
protected

◆ m_successiveRHS

int Nektar::MultiRegions::AssemblyMap::m_successiveRHS
protected

sucessive RHS for iterative solver

Definition at line 412 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetSuccessiveRHS().

◆ m_systemSingular

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