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

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

#include <AssemblyMap.h>

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

Public Member Functions

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

Protected Member Functions

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

Protected Attributes

LibUtilities::SessionReaderSharedPtr m_session
 Session object.
LibUtilities::CommSharedPtr m_comm
 Communicator.
size_t m_hash
 Hash for map.
int m_numLocalBndCoeffs
 Number of local boundary coefficients.
int m_numGlobalBndCoeffs
 Total number of global boundary coefficients.
int m_numLocalDirBndCoeffs
 Number of Local Dirichlet Boundary Coefficients.
int m_numGlobalDirBndCoeffs
 Number of Global Dirichlet Boundary Coefficients.
bool m_systemSingular
 Flag indicating if the system is singular or not.
int m_numLocalCoeffs
 Total number of local coefficients.
int m_numGlobalCoeffs
 Total number of global coefficients.
bool m_signChange
 Flag indicating if modes require sign reversal.
Array< OneD, int > m_localToGlobalBndMap
 Integer map of local boundary coeffs to global space.
Array< OneD, NekDoublem_localToGlobalBndSign
 Integer sign of local boundary coeffs to global space.
Array< OneD, int > m_bndCondCoeffsToGlobalCoeffsMap
 Integer map of bnd cond coeffs to global coefficients.
Array< OneD, NekDoublem_bndCondCoeffsToGlobalCoeffsSign
 Integer map of bnd cond coeffs to global coefficients.
Array< OneD, int > m_bndCondTraceToGlobalTraceMap
 Integer map of bnd cond trace number to global trace number.
Array< OneD, int > m_globalToUniversalBndMap
 Integer map of process coeffs to universal space.
Array< OneD, int > m_globalToUniversalBndMapUnique
 Integer map of unique process coeffs to universal space (signed)
GlobalSysSolnType m_solnType
 The solution type of the global system.
int m_bndSystemBandWidth
 The bandwith of the global bnd system.
PreconditionerType m_preconType
 Type type of preconditioner to use in iterative solver.
NekDouble m_iterativeTolerance
 Tolerance for iterative solver.
int m_successiveRHS
 sucessive RHS for iterative solver
Gs::gs_datam_gsh
Gs::gs_datam_bndGsh
int m_staticCondLevel
 The level of recursion in the case of multi-level static condensation.
int m_numPatches
 The number of patches (~elements) in the current level.
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
 The number of bnd dofs per patch.
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
 The number of int dofs per patch.
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
 Map from the patches of the previous level to the patches of the current level.
int m_lowestStaticCondLevel
 Lowest static condensation level.

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) const
virtual void v_LocalToGlobal (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
virtual void v_GlobalToLocal (const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
virtual void v_GlobalToLocal (const NekVector< NekDouble > &global, NekVector< NekDouble > &loc) const
virtual void v_Assemble (const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
virtual void v_Assemble (const NekVector< NekDouble > &loc, NekVector< NekDouble > &global) const
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal) const
virtual void v_UniversalAssemble (NekVector< NekDouble > &pGlobal) const
virtual void v_UniversalAssemble (Array< OneD, NekDouble > &pGlobal, int offset) const
virtual int v_GetFullSystemBandWidth () const
virtual int v_GetNumNonDirVertexModes () const
virtual int v_GetNumNonDirEdgeModes () const
virtual int v_GetNumNonDirFaceModes () const
virtual int v_GetNumDirEdges () const
virtual int v_GetNumDirFaces () const
virtual int v_GetNumNonDirEdges () const
virtual int v_GetNumNonDirFaces () const
virtual const Array< OneD,
const int > & 
v_GetExtraDirEdges ()
virtual boost::shared_ptr
< AssemblyMap
v_LinearSpaceMap (const ExpList &locexp, GlobalSysSolnType solnType)
 Generate a linear space mapping from existing mapping.

Private Attributes

PatchMapSharedPtr m_patchMapFromPrevLevel
 Mapping information for previous level in MultiLevel Solver.

Detailed Description

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

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

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

Definition at line 59 of file AssemblyMap.h.

Constructor & Destructor Documentation

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

Default constructor.

Initialises an empty mapping.

Definition at line 77 of file AssemblyMap.cpp.

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

Constructor with a communicator.

Definition at line 93 of file AssemblyMap.cpp.

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

:
m_session(pSession),
m_comm(pSession->GetComm()),
m_hash(0),
m_gsh(0),
{
// Default value from Solver Info
m_solnType = pSession->GetSolverInfoAsEnum<GlobalSysSolnType>(
"GlobalSysSoln");
m_preconType = pSession->GetSolverInfoAsEnum<PreconditionerType>(
"Preconditioner");
// Override values with data from GlobalSysSolnInfo section
if(pSession->DefinesGlobalSysSolnInfo(variable, "GlobalSysSoln"))
{
std::string sysSoln = pSession->GetGlobalSysSolnInfo(variable,
"GlobalSysSoln");
m_solnType = pSession->GetValueAsEnum<GlobalSysSolnType>(
"GlobalSysSoln", sysSoln);
}
if(pSession->DefinesGlobalSysSolnInfo(variable, "Preconditioner"))
{
std::string precon = pSession->GetGlobalSysSolnInfo(variable,
"Preconditioner");
m_preconType = pSession->GetValueAsEnum<PreconditionerType>(
"Preconditioner", precon);
}
if(pSession->DefinesGlobalSysSolnInfo(variable,
"IterativeSolverTolerance"))
{
m_iterativeTolerance = boost::lexical_cast<NekDouble>(
pSession->GetGlobalSysSolnInfo(variable,
"IterativeSolverTolerance").c_str());
}
else
{
pSession->LoadParameter("IterativeSolverTolerance",
}
if(pSession->DefinesGlobalSysSolnInfo(variable,"SuccessiveRHS"))
{
m_successiveRHS = boost::lexical_cast<int>(
pSession->GetGlobalSysSolnInfo(variable,
"SuccessiveRHS").c_str());
}
else
{
pSession->LoadParameter("SuccessiveRHS",
}
}
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:

Definition at line 164 of file AssemblyMap.cpp.

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

:
m_session(oldLevelMap->m_session),
m_comm(oldLevelMap->GetComm()),
m_hash(0),
m_globalToUniversalBndMap(oldLevelMap->GetGlobalToUniversalBndMap()),
m_globalToUniversalBndMapUnique(oldLevelMap->GetGlobalToUniversalBndMapUnique()),
m_solnType(oldLevelMap->m_solnType),
m_preconType(oldLevelMap->m_preconType),
m_iterativeTolerance(oldLevelMap->m_iterativeTolerance),
m_successiveRHS(oldLevelMap->m_successiveRHS),
m_gsh(oldLevelMap->m_gsh),
m_bndGsh(oldLevelMap->m_bndGsh),
m_lowestStaticCondLevel(oldLevelMap->m_lowestStaticCondLevel)
{
int i;
int j;
int cnt;
//--------------------------------------------------------------
// -- Extract information from the input argument
int numGlobalBndCoeffsOld = oldLevelMap->GetNumGlobalBndCoeffs();
int numGlobalDirBndCoeffsOld = oldLevelMap->GetNumGlobalDirBndCoeffs();
int numLocalBndCoeffsOld = oldLevelMap->GetNumLocalBndCoeffs();
int numLocalDirBndCoeffsOld = oldLevelMap->GetNumLocalDirBndCoeffs();
bool signChangeOld = oldLevelMap->GetSignChange();
int staticCondLevelOld = oldLevelMap->GetStaticCondLevel();
int numPatchesOld = oldLevelMap->GetNumPatches();
GlobalSysSolnType solnTypeOld = oldLevelMap->GetGlobalSysSolnType();
const Array<OneD, const unsigned int>& numLocalBndCoeffsPerPatchOld = oldLevelMap->GetNumLocalBndCoeffsPerPatch();
//--------------------------------------------------------------
//--------------------------------------------------------------
int newLevel = staticCondLevelOld+1;
/** - 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
*/
Array<OneD, NekDouble> globPatchMask (numGlobalBndCoeffsOld,-1.0);
Array<OneD, NekDouble> globHomPatchMask (globPatchMask+numGlobalDirBndCoeffsOld);
Array<OneD, NekDouble> locPatchMask_NekDouble(numLocalBndCoeffsOld,-3.0);
Array<OneD, int> locPatchMask (numLocalBndCoeffsOld);
// Fill the array globPatchMask as follows:
// - The first part (i.e. the glob bnd dofs) is filled with the
// value -1
// - The second part (i.e. the glob interior dofs) is numbered
// according to the patch it belongs to (i.e. dofs in first block
// all are numbered 0, the second block numbered are 1, etc...)
multiLevelGraph->MaskPatches(newLevel,globHomPatchMask);
// Map from Global Dofs to Local Dofs
// As a result, we know for each local dof whether
// it is mapped to the boundary of the next level, or to which
// patch. Based upon this, we can than later associate every patch
// of the current level with a patch in the next level.
oldLevelMap->GlobalToLocalBndWithoutSign(globPatchMask,locPatchMask_NekDouble);
// Convert the result to an array of integers rather than doubles
RoundNekDoubleToInt(locPatchMask_NekDouble,locPatchMask);
/** - 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.
*/
// Retrieve the number of patches at the next level
int numPatchesWithIntNew = multiLevelGraph->GetNpatchesWithInterior(newLevel);
int numPatchesNew = numPatchesWithIntNew;
// Allocate memory to store the number of local dofs associated to
// each of elemental boundaries of these patches
map<int, int> numLocalBndCoeffsPerPatchNew;
for(int i = 0; i < numPatchesNew; i++)
{
numLocalBndCoeffsPerPatchNew[i] = 0;
}
int minval;
int maxval;
int curPatch;
for(i = cnt = 0; i < numPatchesOld; i++)
{
// For every patch at the current level, the mask array
// locPatchMask should be filled with
// - the same (positive) number for each entry (which will
// correspond to the patch at the next level it belongs to)
// - the same (positive) number for each entry, except some
// entries that are -1 (the enties correspond to -1, will be
// mapped to the local boundary of the next level patch given
// by the positive number)
// - -1 for all entries. In this case, we will make an
// additional patch only consisting of boundaries at the next
// level
minval = *min_element(&locPatchMask[cnt],
&locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
maxval = *max_element(&locPatchMask[cnt],
&locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
ASSERTL0((minval==maxval)||(minval==-1),"These values should never be the same");
if(maxval == -1)
{
curPatch = numPatchesNew;
numLocalBndCoeffsPerPatchNew[curPatch] = 0;
numPatchesNew++;
}
else
{
curPatch = maxval;
}
for(j = 0; j < numLocalBndCoeffsPerPatchOld[i]; j++ )
{
ASSERTL0((locPatchMask[cnt]==maxval)||(locPatchMask[cnt]==minval),
"These values should never be the same");
if(locPatchMask[cnt] == -1)
{
++numLocalBndCoeffsPerPatchNew[curPatch];
}
cnt++;
}
}
// Count how many local dofs of the old level are mapped
// to the local boundary dofs of the new level
m_numPatches = numLocalBndCoeffsPerPatchNew.size();
m_numLocalBndCoeffsPerPatch = Array<OneD, unsigned int>(m_numPatches);
m_numLocalIntCoeffsPerPatch = Array<OneD, unsigned int>(m_numPatches,0u);
for(int i = 0; i < m_numPatches; i++)
{
m_numLocalBndCoeffsPerPatch[i] = (unsigned int) numLocalBndCoeffsPerPatchNew[i];
m_numLocalBndCoeffs += numLocalBndCoeffsPerPatchNew[i];
}
multiLevelGraph->GetNintDofsPerPatch(newLevel,m_numLocalIntCoeffsPerPatch);
// Also initialise some more data members
m_solnType = solnTypeOld;
"This method should only be called for in "
"case of multi-level static condensation.");
m_staticCondLevel = newLevel;
m_signChange = signChangeOld;
m_numLocalDirBndCoeffs = numLocalDirBndCoeffsOld;
m_numGlobalDirBndCoeffs = numGlobalDirBndCoeffsOld;
m_numGlobalBndCoeffs = multiLevelGraph->GetInteriorOffset(newLevel) +
m_numGlobalCoeffs = multiLevelGraph->GetNumGlobalDofs(newLevel) +
{
m_localToGlobalBndSign = Array<OneD,NekDouble>(m_numLocalBndCoeffs);
}
// Set up an offset array that denotes the offset of the local
// boundary degrees of freedom of the next level
Array<OneD, int> numLocalBndCoeffsPerPatchOffset(m_numPatches+1,0);
for(int i = 1; i < m_numPatches+1; i++)
{
numLocalBndCoeffsPerPatchOffset[i] += numLocalBndCoeffsPerPatchOffset[i-1] + numLocalBndCoeffsPerPatchNew[i-1];
}
int additionalPatchCnt = numPatchesWithIntNew;
int newid;
int blockid;
bool isBndDof;
Array<OneD, int> bndDofPerPatchCnt(m_numPatches,0);
for(i = cnt = 0; i < numPatchesOld; i++)
{
minval = *min_element(&locPatchMask[cnt],&locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
maxval = *max_element(&locPatchMask[cnt],&locPatchMask[cnt]+numLocalBndCoeffsPerPatchOld[i]);
ASSERTL0((minval==maxval)||(minval==-1),"These values should never be the same");
if(maxval == -1)
{
curPatch = additionalPatchCnt;
additionalPatchCnt++;
}
else
{
curPatch = maxval;
}
for(j = 0; j < numLocalBndCoeffsPerPatchOld[i]; j++ )
{
ASSERTL0((locPatchMask[cnt]==maxval)||(locPatchMask[cnt]==minval),
"These values should never be the same");
sign = oldLevelMap->GetLocalToGlobalBndSign(cnt);
if(locPatchMask[cnt] == -1)
{
newid = numLocalBndCoeffsPerPatchOffset[curPatch];
m_localToGlobalBndMap[newid] = oldLevelMap->GetLocalToGlobalBndMap(cnt);
{
}
blockid = bndDofPerPatchCnt[curPatch];
isBndDof = true;
numLocalBndCoeffsPerPatchOffset[curPatch]++;
bndDofPerPatchCnt[curPatch]++;
}
else
{
newid = oldLevelMap->GetLocalToGlobalBndMap(cnt) -
blockid = oldLevelMap->GetLocalToGlobalBndMap(cnt)-
m_numGlobalDirBndCoeffs - multiLevelGraph->GetInteriorOffset(newLevel,curPatch);
isBndDof = false;
}
sign = isBndDof?1.0:sign;
m_patchMapFromPrevLevel->SetPatchMap(cnt,curPatch,blockid,isBndDof,sign);
cnt++;
}
}
// Postprocess the computed information - Update the old
// level with the mapping to new evel
// oldLevelMap->SetLocalBndToNextLevelMap(oldLocalBndToNewLevelMap,oldLocalBndToNewLevelSign);
// - Construct the next level mapping object
if(m_staticCondLevel < (multiLevelGraph->GetNlevels()-1))
{
}
}
Nektar::MultiRegions::AssemblyMap::~AssemblyMap ( void  )
virtual

Destructor.

Definition at line 416 of file AssemblyMap.cpp.

{
}

Member Function Documentation

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

Definition at line 729 of file AssemblyMap.cpp.

References v_Assemble().

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

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

Definition at line 736 of file AssemblyMap.cpp.

References v_Assemble().

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

Definition at line 1049 of file AssemblyMap.cpp.

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

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

{
AssembleBnd(loc.GetPtr(), global.GetPtr(), offset);
}
void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const

Definition at line 1057 of file AssemblyMap.cpp.

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

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

Definition at line 1065 of file AssemblyMap.cpp.

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

{
ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local array is not of correct dimension");
ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs-offset,"Global array is not of correct dimension");
Array<OneD,NekDouble> tmp(m_numGlobalBndCoeffs,0.0);
{
}
else
{
}
Vmath::Vcopy(m_numGlobalBndCoeffs-offset, tmp.get() + offset, 1, global.get(), 1);
}
void Nektar::MultiRegions::AssemblyMap::AssembleBnd ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const

Definition at line 1086 of file AssemblyMap.cpp.

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

{
ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
Vmath::Zero(m_numGlobalBndCoeffs, global.get(), 1);
{
loc.get(), m_localToGlobalBndMap.get(), global.get());
}
else
{
Vmath::Assmb(m_numLocalBndCoeffs,loc.get(), m_localToGlobalBndMap.get(), global.get());
}
}
bool Nektar::MultiRegions::AssemblyMap::AtLastLevel ( ) const

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

Definition at line 1172 of file AssemblyMap.cpp.

References m_nextLevelLocalToGlobalMap.

{
return !( m_nextLevelLocalToGlobalMap.get() );
}
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 435 of file AssemblyMap.cpp.

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

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapCG1D::AssemblyMapCG1D(), Nektar::MultiRegions::AssemblyMapCG2D::AssemblyMapCG2D(), Nektar::MultiRegions::AssemblyMapCG3D::AssemblyMapCG3D(), and Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG().

{
int i,j;
int cnt = 0;
int locSize;
int maxId;
int minId;
int bwidth = -1;
for(i = 0; i < m_numPatches; ++i)
{
maxId = -1;
for(j = 0; j < locSize; j++)
{
{
if(m_localToGlobalBndMap[cnt+j] > maxId)
{
maxId = m_localToGlobalBndMap[cnt+j];
}
if(m_localToGlobalBndMap[cnt+j] < minId)
{
minId = m_localToGlobalBndMap[cnt+j];
}
}
}
bwidth = (bwidth>(maxId-minId))?bwidth:(maxId-minId);
cnt+=locSize;
}
}
int Nektar::MultiRegions::AssemblyMap::GetBndCondCoeffsToGlobalCoeffsMap ( const int  i)

Retrieves the global index corresponding to a boundary expansion mode.

Definition at line 858 of file AssemblyMap.cpp.

References m_bndCondCoeffsToGlobalCoeffsMap.

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

Retrieves the global indices corresponding to the boundary expansion modes.

Definition at line 892 of file AssemblyMap.cpp.

References m_bndCondCoeffsToGlobalCoeffsMap.

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

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

Definition at line 879 of file AssemblyMap.cpp.

References m_bndCondCoeffsToGlobalCoeffsSign, and m_signChange.

{
{
}
else
{
return 1.0;
}
}
int Nektar::MultiRegions::AssemblyMap::GetBndCondTraceToGlobalTraceMap ( const int  i)

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

Definition at line 865 of file AssemblyMap.cpp.

References ASSERTL1, and m_bndCondTraceToGlobalTraceMap.

{
"Index out of range.");
}
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetBndCondTraceToGlobalTraceMap ( )

Definition at line 874 of file AssemblyMap.cpp.

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

Returns the bandwidth of the boundary system.

Definition at line 1131 of file AssemblyMap.cpp.

References m_bndSystemBandWidth.

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

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

Retrieves the communicator.

Definition at line 651 of file AssemblyMap.cpp.

References m_comm.

{
return m_comm;
}
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetExtraDirEdges ( )

Definition at line 802 of file AssemblyMap.cpp.

References v_GetExtraDirEdges().

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

Definition at line 762 of file AssemblyMap.cpp.

References v_GetFullSystemBandWidth().

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

Returns the method of solving global systems.

Definition at line 1178 of file AssemblyMap.cpp.

References m_solnType.

Referenced by AssemblyMap().

{
return m_solnType;
}
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::GetGlobalToUniversalBndMap ( )

Definition at line 835 of file AssemblyMap.cpp.

References m_globalToUniversalBndMap.

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

Definition at line 840 of file AssemblyMap.cpp.

References m_globalToUniversalBndMapUnique.

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

Definition at line 666 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMap().

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

Definition at line 681 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMap().

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

Definition at line 671 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMapUnique().

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

Definition at line 686 of file AssemblyMap.cpp.

References v_GetGlobalToUniversalMapUnique().

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

Retrieves the hash of this map.

Definition at line 656 of file AssemblyMap.cpp.

References m_hash.

{
return m_hash;
}
NekDouble Nektar::MultiRegions::AssemblyMap::GetIterativeTolerance ( ) const

Definition at line 1188 of file AssemblyMap.cpp.

References m_iterativeTolerance.

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

Retrieve the global index of a given local boundary mode.

Definition at line 812 of file AssemblyMap.cpp.

References m_localToGlobalBndMap.

Referenced by AssemblyMap().

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

Retrieve the global indices of the local boundary modes.

Definition at line 818 of file AssemblyMap.cpp.

References m_localToGlobalBndMap.

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

Retrieve the sign change of a given local boundary mode.

Definition at line 845 of file AssemblyMap.cpp.

References m_localToGlobalBndSign, and m_signChange.

Referenced by AssemblyMap().

{
{
}
else
{
return 1.0;
}
}
Array< OneD, const NekDouble > Nektar::MultiRegions::AssemblyMap::GetLocalToGlobalBndSign ( void  ) const

Retrieve the sign change for all local boundary modes.

Definition at line 830 of file AssemblyMap.cpp.

References m_localToGlobalBndSign.

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

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

Definition at line 661 of file AssemblyMap.cpp.

References v_GetLocalToGlobalMap().

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

Definition at line 676 of file AssemblyMap.cpp.

References v_GetLocalToGlobalMap().

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

Definition at line 691 of file AssemblyMap.cpp.

References v_GetLocalToGlobalSign().

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

Definition at line 696 of file AssemblyMap.cpp.

References v_GetLocalToGlobalSign().

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

Definition at line 295 of file AssemblyMap.h.

References m_lowestStaticCondLevel.

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

References m_nextLevelLocalToGlobalMap.

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

Definition at line 782 of file AssemblyMap.cpp.

References v_GetNumDirEdges().

{
return v_GetNumDirEdges();
}
int Nektar::MultiRegions::AssemblyMap::GetNumDirFaces ( ) const

Definition at line 787 of file AssemblyMap.cpp.

References v_GetNumDirFaces().

{
return v_GetNumDirFaces();
}
int Nektar::MultiRegions::AssemblyMap::GetNumGlobalBndCoeffs ( ) const

Returns the total number of global boundary coefficients.

Definition at line 914 of file AssemblyMap.cpp.

References m_numGlobalBndCoeffs.

Referenced by AssemblyMap().

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

Returns the total number of global coefficients.

Definition at line 924 of file AssemblyMap.cpp.

References m_numGlobalCoeffs.

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

Returns the number of global Dirichlet boundary coefficients.

Definition at line 898 of file AssemblyMap.cpp.

References m_numGlobalDirBndCoeffs.

Referenced by AssemblyMap().

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

Returns the total number of local boundary coefficients.

Definition at line 909 of file AssemblyMap.cpp.

References m_numLocalBndCoeffs.

Referenced by AssemblyMap().

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

Returns the number of local boundary coefficients in each patch.

Definition at line 1147 of file AssemblyMap.cpp.

References m_numLocalBndCoeffsPerPatch.

Referenced by AssemblyMap().

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

Returns the total number of local coefficients.

Definition at line 919 of file AssemblyMap.cpp.

References m_numLocalCoeffs.

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

Returns the number of local Dirichlet boundary coefficients.

Definition at line 904 of file AssemblyMap.cpp.

References m_numLocalDirBndCoeffs.

Referenced by AssemblyMap().

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

Returns the number of local interior coefficients in each patch.

Definition at line 1154 of file AssemblyMap.cpp.

References m_numLocalIntCoeffsPerPatch.

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

Definition at line 772 of file AssemblyMap.cpp.

References v_GetNumNonDirEdgeModes().

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

Definition at line 792 of file AssemblyMap.cpp.

References v_GetNumNonDirEdges().

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

Definition at line 777 of file AssemblyMap.cpp.

References v_GetNumNonDirFaceModes().

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

Definition at line 797 of file AssemblyMap.cpp.

References v_GetNumNonDirFaces().

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

Definition at line 767 of file AssemblyMap.cpp.

References v_GetNumNonDirVertexModes().

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

Returns the number of patches in this static condensation level.

Definition at line 1141 of file AssemblyMap.cpp.

References m_numPatches.

Referenced by AssemblyMap().

{
return m_numPatches;
}
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 1166 of file AssemblyMap.cpp.

References m_patchMapFromPrevLevel.

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

Definition at line 1183 of file AssemblyMap.cpp.

References m_preconType.

{
return m_preconType;
}
bool Nektar::MultiRegions::AssemblyMap::GetSignChange ( )

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

Definition at line 823 of file AssemblyMap.cpp.

References m_signChange.

Referenced by AssemblyMap().

{
return m_signChange;
}
bool Nektar::MultiRegions::AssemblyMap::GetSingularSystem ( ) const

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

Definition at line 929 of file AssemblyMap.cpp.

References m_systemSingular.

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

Returns the level of static condensation for this map.

Definition at line 1136 of file AssemblyMap.cpp.

References m_staticCondLevel.

Referenced by AssemblyMap().

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

Definition at line 1193 of file AssemblyMap.cpp.

References m_successiveRHS.

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

Definition at line 715 of file AssemblyMap.cpp.

References v_GlobalToLocal().

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

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

Definition at line 722 of file AssemblyMap.cpp.

References v_GlobalToLocal().

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

Definition at line 934 of file AssemblyMap.cpp.

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

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

{
GlobalToLocalBnd(global.GetPtr(), loc.GetPtr(), offset);
}
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const NekVector< NekDouble > &  global,
NekVector< NekDouble > &  loc 
) const

Definition at line 943 of file AssemblyMap.cpp.

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

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

Definition at line 951 of file AssemblyMap.cpp.

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

{
ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs-offset,"Global vector is not of correct dimension");
// offset input data by length "offset" for Dirichlet boundary conditions.
Array<OneD,NekDouble> tmp(m_numGlobalBndCoeffs,0.0);
Vmath::Vcopy(m_numGlobalBndCoeffs-offset, global.get(), 1, tmp.get() + offset, 1);
{
}
else
{
}
}
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBnd ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc 
) const

Definition at line 973 of file AssemblyMap.cpp.

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

{
ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
{
}
else
{
Vmath::Gathr(m_numLocalBndCoeffs, global.get(), m_localToGlobalBndMap.get(), loc.get());
}
}
void Nektar::MultiRegions::AssemblyMap::GlobalToLocalBndWithoutSign ( const Array< OneD, const NekDouble > &  global,
Array< OneD, NekDouble > &  loc 
)
protected

Definition at line 1198 of file AssemblyMap.cpp.

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

Referenced by AssemblyMap().

{
ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
Vmath::Gathr(m_numLocalBndCoeffs, global.get(), m_localToGlobalBndMap.get(), loc.get());
}
boost::shared_ptr< AssemblyMap > Nektar::MultiRegions::AssemblyMap::LinearSpaceMap ( const ExpList locexp,
GlobalSysSolnType  solnType 
)

Definition at line 807 of file AssemblyMap.cpp.

References v_LinearSpaceMap().

{
return v_LinearSpaceMap(locexp, solnType);
}
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global,
int  offset 
) const

Definition at line 990 of file AssemblyMap.cpp.

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

Referenced by LocalBndToGlobal().

{
LocalBndToGlobal(loc.GetPtr(), global.GetPtr(), offset);
}
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const

Definition at line 1023 of file AssemblyMap.cpp.

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

{
LocalBndToGlobal(loc.GetPtr(), global.GetPtr());
}
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global,
int  offset 
) const

Definition at line 999 of file AssemblyMap.cpp.

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

{
ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs-offset,"Global vector is not of correct dimension");
// offset input data by length "offset" for Dirichlet boundary conditions.
Array<OneD,NekDouble> tmp(m_numGlobalBndCoeffs,0.0);
{
}
else
{
}
Vmath::Vcopy(m_numGlobalBndCoeffs-offset, tmp.get()+offset, 1, global.get(), 1);
}
void Nektar::MultiRegions::AssemblyMap::LocalBndToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const

Definition at line 1031 of file AssemblyMap.cpp.

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

{
ASSERTL1(loc.num_elements() >= m_numLocalBndCoeffs,"Local vector is not of correct dimension");
ASSERTL1(global.num_elements() >= m_numGlobalBndCoeffs,"Global vector is not of correct dimension");
{
}
else
{
Vmath::Scatr(m_numLocalBndCoeffs, loc.get(), m_localToGlobalBndMap.get(), global.get());
}
}
void Nektar::MultiRegions::AssemblyMap::LocalToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const

Definition at line 701 of file AssemblyMap.cpp.

References v_LocalToGlobal().

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

{
v_LocalToGlobal(loc,global);
}
void Nektar::MultiRegions::AssemblyMap::LocalToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const

Definition at line 708 of file AssemblyMap.cpp.

References v_LocalToGlobal().

{
v_LocalToGlobal(loc,global);
}
void Nektar::MultiRegions::AssemblyMap::PrintStats ( std::ostream &  out,
std::string  variable 
) const

Definition at line 1208 of file AssemblyMap.cpp.

References AssembleBnd(), m_globalToUniversalBndMapUnique, m_numGlobalBndCoeffs, m_numGlobalCoeffs, m_numGlobalDirBndCoeffs, m_numLocalBndCoeffs, m_numLocalCoeffs, m_numLocalDirBndCoeffs, m_session, Nektar::LibUtilities::ReduceMax, Nektar::LibUtilities::ReduceMin, and Nektar::LibUtilities::ReduceSum.

{
= m_session->GetComm()->GetRowComm();
bool isRoot = vRowComm->GetRank() == 0;
int n = vRowComm->GetSize();
int i;
// Determine number of global degrees of freedom.
int globBndCnt = 0, globDirCnt = 0;
for (i = 0; i < m_numGlobalBndCoeffs; ++i)
{
{
globBndCnt++;
{
globDirCnt++;
}
}
}
int globCnt = m_numGlobalCoeffs - m_numGlobalBndCoeffs + globBndCnt;
// Calculate maximum valency
Array<OneD, NekDouble> tmpLoc (m_numLocalBndCoeffs, 1.0);
Array<OneD, NekDouble> tmpGlob(m_numGlobalBndCoeffs, 0.0);
AssembleBnd(tmpLoc, tmpGlob);
int totGlobDof = globCnt;
int totGlobBndDof = globBndCnt;
int totGlobDirDof = globDirCnt;
int totLocalDof = m_numLocalCoeffs;
int totLocalBndDof = m_numLocalBndCoeffs;
int totLocalDirDof = m_numLocalDirBndCoeffs;
int meanValence = 0;
int maxValence = 0;
int minValence = 10000000;
for (int i = 0; i < m_numGlobalBndCoeffs; ++i)
{
{
continue;
}
if (tmpGlob[i] > maxValence)
{
maxValence = tmpGlob[i];
}
if (tmpGlob[i] < minValence)
{
minValence = tmpGlob[i];
}
meanValence += tmpGlob[i];
}
vRowComm->AllReduce(maxValence, LibUtilities::ReduceMax);
vRowComm->AllReduce(minValence, LibUtilities::ReduceMin);
vRowComm->AllReduce(meanValence, LibUtilities::ReduceSum);
vRowComm->AllReduce(totGlobDof, LibUtilities::ReduceSum);
vRowComm->AllReduce(totGlobBndDof, LibUtilities::ReduceSum);
vRowComm->AllReduce(totGlobDirDof, LibUtilities::ReduceSum);
vRowComm->AllReduce(totLocalDof, LibUtilities::ReduceSum);
vRowComm->AllReduce(totLocalBndDof, LibUtilities::ReduceSum);
vRowComm->AllReduce(totLocalDirDof, LibUtilities::ReduceSum);
meanValence /= totGlobBndDof;
if (isRoot)
{
out << "Assembly map statistics for field " << variable << ":"
<< endl;
out << " - Number of local/global dof : "
<< totLocalDof << " " << totGlobDof << endl;
out << " - Number of local/global boundary dof : "
<< totLocalBndDof << " " << totGlobBndDof << endl;
out << " - Number of local/global Dirichlet dof : "
<< totLocalDirDof << " " << totGlobDirDof << endl;
out << " - dof valency (min/max/mean) : "
<< minValence << " " << maxValence << " " << meanValence
<< endl;
if (n > 1)
{
NekDouble mean = m_numLocalCoeffs, mean2 = mean * mean;
NekDouble minval = mean, maxval = mean;
Array<OneD, NekDouble> tmp(1);
for (i = 1; i < n; ++i)
{
vRowComm->Recv(i, tmp);
mean += tmp[0];
mean2 += tmp[0]*tmp[0];
if (tmp[0] > maxval)
{
maxval = tmp[0];
}
if (tmp[0] < minval)
{
minval = tmp[0];
}
}
out << " - Local dof dist. (min/max/mean/dev) : "
<< minval << " " << maxval << " " << (mean / n) << " "
<< sqrt(mean2/n - mean*mean/n/n) << endl;
vRowComm->Block();
mean = minval = maxval = m_numLocalBndCoeffs;
mean2 = mean * mean;
for (i = 1; i < n; ++i)
{
vRowComm->Recv(i, tmp);
mean += tmp[0];
mean2 += tmp[0]*tmp[0];
if (tmp[0] > maxval)
{
maxval = tmp[0];
}
if (tmp[0] < minval)
{
minval = tmp[0];
}
}
out << " - Local bnd dof dist. (min/max/mean/dev) : "
<< minval << " " << maxval << " " << (mean / n) << " "
<< sqrt(mean2/n - mean*mean/n/n) << endl;
}
}
else
{
Array<OneD, NekDouble> tmp(1);
tmp[0] = m_numLocalCoeffs;
vRowComm->Send(0, tmp);
vRowComm->Block();
vRowComm->Send(0, tmp);
}
}
void Nektar::MultiRegions::AssemblyMap::SetNextLevelLocalToGlobalMap ( AssemblyMapSharedPtr  pNextLevelLocalToGlobalMap)
void Nektar::MultiRegions::AssemblyMap::UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal) const

Definition at line 743 of file AssemblyMap.cpp.

References v_UniversalAssemble().

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

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

Definition at line 749 of file AssemblyMap.cpp.

References v_UniversalAssemble().

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

Definition at line 755 of file AssemblyMap.cpp.

References v_UniversalAssemble().

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

Definition at line 1107 of file AssemblyMap.cpp.

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

Referenced by AssembleBnd(), LocalBndToGlobal(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), and UniversalAssembleBnd().

{
ASSERTL1(pGlobal.num_elements() == m_numGlobalBndCoeffs,
"Wrong size.");
}
void Nektar::MultiRegions::AssemblyMap::UniversalAssembleBnd ( NekVector< NekDouble > &  pGlobal) const

Definition at line 1115 of file AssemblyMap.cpp.

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

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

Definition at line 1121 of file AssemblyMap.cpp.

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

{
Array<OneD, NekDouble> tmp(offset);
if (offset > 0) Vmath::Vcopy(offset, pGlobal, 1, tmp, 1);
if (offset > 0) Vmath::Vcopy(offset, tmp, 1, pGlobal, 1);
}
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 552 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by Assemble().

{
ASSERTL0(false, "Not defined for this type of mapping.");
}
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 559 of file AssemblyMap.cpp.

References ASSERTL0.

{
ASSERTL0(false, "Not defined for this type of mapping.");
}
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetExtraDirEdges ( )
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 636 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetExtraDirEdges().

{
ASSERTL0(false, "Not defined for this type of mapping.");
static Array<OneD, const int> result;
return result;
}
int Nektar::MultiRegions::AssemblyMap::v_GetFullSystemBandWidth ( ) const
privatevirtual

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

Definition at line 588 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetFullSystemBandWidth().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
int Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMap ( const int  i) const
privatevirtual

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

Definition at line 478 of file AssemblyMap.cpp.

References ASSERTL0.

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMap ( void  )
privatevirtual

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

Definition at line 497 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetGlobalToUniversalMap().

{
ASSERTL0(false, "Not defined for this type of mapping.");
static Array<OneD, const int> result;
return result;
}
int Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMapUnique ( const int  i) const
privatevirtual

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

Definition at line 484 of file AssemblyMap.cpp.

References ASSERTL0.

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetGlobalToUniversalMapUnique ( void  )
privatevirtual

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

Definition at line 504 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetGlobalToUniversalMapUnique().

{
ASSERTL0(false, "Not defined for this type of mapping.");
static Array<OneD, const int> result;
return result;
}
int Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalMap ( const int  i) const
privatevirtual

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

Definition at line 472 of file AssemblyMap.cpp.

References ASSERTL0.

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
const Array< OneD, const int > & Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalMap ( void  )
privatevirtual

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

Definition at line 490 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetLocalToGlobalMap().

{
ASSERTL0(false, "Not defined for this type of mapping.");
static Array<OneD,const int> result;
return result;
}
NekDouble Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalSign ( const int  i) const
privatevirtual

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

Definition at line 511 of file AssemblyMap.cpp.

References ASSERTL0.

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0.0;
}
const Array< OneD, NekDouble > & Nektar::MultiRegions::AssemblyMap::v_GetLocalToGlobalSign ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 517 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetLocalToGlobalSign().

{
ASSERTL0(false, "Not defined for this type of mapping.");
static Array<OneD, NekDouble> result;
return result;
}
int Nektar::MultiRegions::AssemblyMap::v_GetNumDirEdges ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 612 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumDirEdges().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
int Nektar::MultiRegions::AssemblyMap::v_GetNumDirFaces ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 618 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumDirFaces().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirEdgeModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 600 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirEdgeModes().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirEdges ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 624 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirEdges().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirFaceModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 606 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirFaceModes().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirFaces ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 630 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirFaces().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
int Nektar::MultiRegions::AssemblyMap::v_GetNumNonDirVertexModes ( ) const
privatevirtual

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 594 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GetNumNonDirVertexModes().

{
ASSERTL0(false, "Not defined for this type of mapping.");
return 0;
}
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 538 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by GlobalToLocal().

{
ASSERTL0(false, "Not defined for this type of mapping.");
}
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 545 of file AssemblyMap.cpp.

References ASSERTL0.

{
ASSERTL0(false, "Not defined for this type of mapping.");
}
boost::shared_ptr< AssemblyMap > Nektar::MultiRegions::AssemblyMap::v_LinearSpaceMap ( const ExpList locexp,
GlobalSysSolnType  solnType 
)
privatevirtual

Generate a linear space mapping from existing mapping.

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 643 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by LinearSpaceMap().

{
ASSERTL0(false, "Not defined for this sub class");
static boost::shared_ptr<AssemblyMap> result;
return result;
}
void Nektar::MultiRegions::AssemblyMap::v_LocalToGlobal ( const Array< OneD, const NekDouble > &  loc,
Array< OneD, NekDouble > &  global 
) const
privatevirtual

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

Definition at line 524 of file AssemblyMap.cpp.

References ASSERTL0.

Referenced by LocalToGlobal().

{
ASSERTL0(false, "Not defined for this type of mapping.");
}
void Nektar::MultiRegions::AssemblyMap::v_LocalToGlobal ( const NekVector< NekDouble > &  loc,
NekVector< NekDouble > &  global 
) const
privatevirtual

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

Definition at line 531 of file AssemblyMap.cpp.

References ASSERTL0.

{
ASSERTL0(false, "Not defined for this type of mapping.");
}
void Nektar::MultiRegions::AssemblyMap::v_UniversalAssemble ( Array< OneD, NekDouble > &  pGlobal) const
privatevirtual

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

Definition at line 566 of file AssemblyMap.cpp.

Referenced by UniversalAssemble().

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

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

Definition at line 573 of file AssemblyMap.cpp.

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

Reimplemented in Nektar::MultiRegions::AssemblyMapCG.

Definition at line 580 of file AssemblyMap.cpp.

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

Member Data Documentation

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

Integer map of bnd cond coeffs to global coefficients.

Definition at line 351 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetBndCondCoeffsToGlobalCoeffsMap(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

Integer map of bnd cond coeffs to global coefficients.

Definition at line 353 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), GetBndCondCoeffsToGlobalCoeffsSign(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

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

Definition at line 355 of file AssemblyMap.h.

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

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

Definition at line 376 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), and UniversalAssembleBnd().

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

The bandwith of the global bnd system.

Definition at line 364 of file AssemblyMap.h.

Referenced by CalculateBndSystemBandWidth(), and GetBndSystemBandWidth().

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

Communicator.

Definition at line 305 of file AssemblyMap.h.

Referenced by GetComm(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DGraphC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalTraceMap(), and Nektar::MultiRegions::AssemblyMapCG::v_LinearSpaceMap().

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

Integer map of process coeffs to universal space.

Definition at line 357 of file AssemblyMap.h.

Referenced by GetGlobalToUniversalBndMap(), Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), and Nektar::MultiRegions::AssemblyMapDG::v_GetGlobalToUniversalMap().

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

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

Definition at line 359 of file AssemblyMap.h.

Referenced by GetGlobalToUniversalBndMapUnique(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), and Nektar::MultiRegions::AssemblyMapDG::v_GetGlobalToUniversalMapUnique().

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

Definition at line 375 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal(), Nektar::MultiRegions::AssemblyMapCG::v_UniversalAssemble(), and Nektar::MultiRegions::AssemblyMapDG::v_UniversalAssemble().

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

Hash for map.

Definition at line 308 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), GetHash(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

Tolerance for iterative solver.

Definition at line 370 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetIterativeTolerance().

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

Integer map of local boundary coeffs to global space.

Definition at line 347 of file AssemblyMap.h.

Referenced by AssembleBnd(), AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), CalculateBndSystemBandWidth(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetLocalToGlobalBndMap(), GlobalToLocalBnd(), GlobalToLocalBndWithoutSign(), LocalBndToGlobal(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), and Nektar::MultiRegions::AssemblyMapDG::v_GetLocalToGlobalMap().

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

Integer sign of local boundary coeffs to global space.

Definition at line 349 of file AssemblyMap.h.

Referenced by AssembleBnd(), AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetLocalToGlobalBndSign(), GlobalToLocalBnd(), LocalBndToGlobal(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

Lowest static condensation level.

Definition at line 393 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), GetLowestStaticCondLevel(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DGraphC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), AtLastLevel(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNextLevelLocalToGlobalMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

Total number of global boundary coefficients.

Definition at line 313 of file AssemblyMap.h.

Referenced by AssembleBnd(), AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumGlobalBndCoeffs(), GlobalToLocalBnd(), GlobalToLocalBndWithoutSign(), LocalBndToGlobal(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapDG::SetUpUniversalDGMap(), and UniversalAssembleBnd().

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

Total number of global coefficients.

This corresponds to the number of total number of coefficients

Definition at line 341 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumGlobalCoeffs(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG::SetUpUniversalC0ContMap(), Nektar::MultiRegions::AssemblyMapCG::v_Assemble(), and Nektar::MultiRegions::AssemblyMapCG::v_LinearSpaceMap().

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

Number of Global Dirichlet Boundary Coefficients.

Definition at line 317 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), CalculateBndSystemBandWidth(), Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumGlobalDirBndCoeffs(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG::v_LinearSpaceMap().

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

Number of local boundary coefficients.

Definition at line 311 of file AssemblyMap.h.

Referenced by AssembleBnd(), AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), CalculateBndSystemBandWidth(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumLocalBndCoeffs(), GlobalToLocalBnd(), GlobalToLocalBndWithoutSign(), LocalBndToGlobal(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DGraphC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

The number of bnd dofs per patch.

Definition at line 384 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), CalculateBndSystemBandWidth(), Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumLocalBndCoeffsPerPatch(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

int Nektar::MultiRegions::AssemblyMap::m_numLocalCoeffs
protected

Total number of local coefficients.

This corresponds to the number of total number of coefficients

Definition at line 330 of file AssemblyMap.h.

Referenced by Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumLocalCoeffs(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG::v_Assemble(), Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal(), and Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal().

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

Number of Local Dirichlet Boundary Coefficients.

Definition at line 315 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumLocalDirBndCoeffs(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

The number of int dofs per patch.

Definition at line 386 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumLocalIntCoeffsPerPatch(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

The number of patches (~elements) in the current level.

Definition at line 382 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), CalculateBndSystemBandWidth(), Nektar::MultiRegions::AssemblyMapCG::CalculateFullSystemBandWidth(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetNumPatches(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

Mapping information for previous level in MultiLevel Solver.

Definition at line 404 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetPatchMapFromPrevLevel().

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

Type type of preconditioner to use in iterative solver.

Definition at line 367 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetPreconType().

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

Session object.

Definition at line 302 of file AssemblyMap.h.

Referenced by Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), PrintStats(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DGraphC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG::v_LinearSpaceMap().

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

Flag indicating if modes require sign reversal.

Definition at line 344 of file AssemblyMap.h.

Referenced by AssembleBnd(), AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetBndCondCoeffsToGlobalCoeffsSign(), GetLocalToGlobalBndSign(), GetSignChange(), GlobalToLocalBnd(), LocalBndToGlobal(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG::v_Assemble(), Nektar::MultiRegions::AssemblyMapCG::v_GetLocalToGlobalSign(), Nektar::MultiRegions::AssemblyMapCG::v_GlobalToLocal(), and Nektar::MultiRegions::AssemblyMapCG::v_LocalToGlobal().

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

The solution type of the global system.

Definition at line 362 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), GetGlobalSysSolnType(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DGraphC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

The level of recursion in the case of multi-level static condensation.

Definition at line 380 of file AssemblyMap.h.

Referenced by AssemblyMap(), Nektar::MultiRegions::AssemblyMapDG::AssemblyMapDG(), Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetStaticCondLevel(), Nektar::MultiRegions::AssemblyMapCG1D::SetUp1DExpansionC0ContMap(), Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap(), and Nektar::MultiRegions::AssemblyMapCG3D::SetUp3DExpansionC0ContMap().

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

sucessive RHS for iterative solver

Definition at line 373 of file AssemblyMap.h.

Referenced by AssemblyMap(), and GetSuccessiveRHS().

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

Flag indicating if the system is singular or not.

Definition at line 319 of file AssemblyMap.h.

Referenced by Nektar::CoupledLocalToGlobalC0ContMap::CoupledLocalToGlobalC0ContMap(), GetSingularSystem(), and Nektar::MultiRegions::AssemblyMapCG2D::SetUp2DExpansionC0ContMap().