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

A global linear system. More...

#include <GlobalLinSysIterativeStaticCond.h>

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

Public Types

typedef NekSparseDiagBlkMatrix
< StorageSmvBsr< NekDouble > > 
DNekSmvBsrDiagBlkMat
 
typedef boost::shared_ptr
< DNekSmvBsrDiagBlkMat
DNekSmvBsrDiagBlkMatSharedPtr
 

Public Member Functions

 GlobalLinSysIterativeStaticCond (const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &locToGloMap)
 Constructor for full direct matrix solve. More...
 
 GlobalLinSysIterativeStaticCond (const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const DNekScalBlkMatSharedPtr pSchurCompl, const DNekScalBlkMatSharedPtr pBinvD, const DNekScalBlkMatSharedPtr pC, const DNekScalBlkMatSharedPtr pInvD, const boost::shared_ptr< AssemblyMap > &locToGloMap, const PreconditionerSharedPtr pPrecon)
 Constructor for full direct matrix solve. More...
 
virtual ~GlobalLinSysIterativeStaticCond ()
 
- Public Member Functions inherited from Nektar::MultiRegions::GlobalLinSysIterative
 GlobalLinSysIterative (const GlobalLinSysKey &pKey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
 Constructor for full direct matrix solve. More...
 
virtual ~GlobalLinSysIterative ()
 
- Public Member Functions inherited from Nektar::MultiRegions::GlobalLinSys
 GlobalLinSys (const GlobalLinSysKey &pKey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
 Constructor for full direct matrix solve. More...
 
virtual ~GlobalLinSys ()
 
const GlobalLinSysKeyGetKey (void) const
 Returns the key associated with the system. More...
 
const boost::weak_ptr< ExpList > & GetLocMat (void) const
 
void InitObject ()
 
void Initialise (const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
 
void Solve (const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out, const AssemblyMapSharedPtr &locToGloMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
 Solve the linear system for given input and output vectors using a specified local to global map. More...
 
boost::shared_ptr< GlobalLinSysGetSharedThisPtr ()
 Returns a shared pointer to the current object. More...
 
int GetNumBlocks ()
 
DNekScalMatSharedPtr GetBlock (unsigned int n)
 
DNekScalBlkMatSharedPtr GetStaticCondBlock (unsigned int n)
 
void DropStaticCondBlock (unsigned int n)
 
void SolveLinearSystem (const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir=0)
 Solve the linear system for given input and output vectors. More...
 
- Public Member Functions inherited from Nektar::MultiRegions::GlobalLinSysStaticCond
 GlobalLinSysStaticCond (const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &locToGloMap)
 Constructor for full direct matrix solve. More...
 
virtual ~GlobalLinSysStaticCond ()
 

Static Public Member Functions

static GlobalLinSysSharedPtr create (const GlobalLinSysKey &pLinSysKey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of class. More...
 
static std::string className2
 

Protected Member Functions

virtual DNekScalBlkMatSharedPtr v_GetStaticCondBlock (unsigned int n)
 Retrieves a the static condensation block matrices from n-th expansion using the matrix key provided by the m_linSysKey. More...
 
virtual
GlobalLinSysStaticCondSharedPtr 
v_Recurse (const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const DNekScalBlkMatSharedPtr pSchurCompl, const DNekScalBlkMatSharedPtr pBinvD, const DNekScalBlkMatSharedPtr pC, const DNekScalBlkMatSharedPtr pInvD, const boost::shared_ptr< AssemblyMap > &locToGloMap)
 
virtual DNekScalBlkMatSharedPtr v_PreSolve (int scLevel, NekVector< NekDouble > &F_GlobBnd)
 
virtual void v_BasisTransform (Array< OneD, NekDouble > &pInOut, int offset)
 
virtual void v_BasisInvTransform (Array< OneD, NekDouble > &pInOut)
 
- Protected Member Functions inherited from Nektar::MultiRegions::GlobalLinSysIterative
void DoAconjugateProjection (const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir)
 A-conjugate projection technique. More...
 
void DoConjugateGradient (const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir)
 Actual iterative solve. More...
 
void Set_Rhs_Magnitude (const NekVector< NekDouble > &pIn)
 
- Protected Member Functions inherited from Nektar::MultiRegions::GlobalLinSys
virtual DNekScalMatSharedPtr v_GetBlock (unsigned int n)
 Retrieves the block matrix from n-th expansion using the matrix key provided by the m_linSysKey. More...
 
virtual void v_DropStaticCondBlock (unsigned int n)
 Releases the static condensation block matrices from NekManager of n-th expansion using the matrix key provided by the m_linSysKey. More...
 
PreconditionerSharedPtr CreatePrecon (AssemblyMapSharedPtr asmMap)
 Create a preconditioner object from the parameters defined in the supplied assembly map. More...
 
- Protected Member Functions inherited from Nektar::MultiRegions::GlobalLinSysStaticCond
virtual int v_GetNumBlocks ()
 Get the number of blocks in this system. More...
 
virtual void v_Solve (const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out, const AssemblyMapSharedPtr &locToGloMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
 Solve the linear system for given input and output vectors using a specified local to global map. More...
 
virtual void v_Initialise (const boost::shared_ptr< AssemblyMap > &locToGloMap)
 Initialise this object. More...
 
void SetupTopLevel (const boost::shared_ptr< AssemblyMap > &locToGloMap)
 Set up the storage for the Schur complement or the top level of the multi-level Schur complement. More...
 
void ConstructNextLevelCondensedSystem (const boost::shared_ptr< AssemblyMap > &locToGloMap)
 

Private Member Functions

virtual void v_InitObject ()
 
void v_AssembleSchurComplement (const boost::shared_ptr< AssemblyMap > locToGloMap)
 Assemble the Schur complement matrix. More...
 
void PrepareLocalSchurComplement ()
 Prepares local representation of Schur complement stored as a sparse block-diagonal matrix. More...
 
virtual void v_DoMatrixMultiply (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
 Perform a Shur-complement matrix multiply operation. More...
 
virtual void v_UniqueMap ()
 

Private Attributes

DNekScalBlkMatSharedPtr m_S1Blk
 
std::vector< double > m_storage
 Dense storage for block Schur complement matrix. More...
 
std::vector< const double * > m_denseBlocks
 Vector of pointers to local matrix data. More...
 
Array< OneD, unsigned int > m_rows
 Ranks of local matrices. More...
 
Array< OneD, NekDoublem_scale
 Scaling factors for local matrices. More...
 
DNekSmvBsrDiagBlkMatSharedPtr m_sparseSchurCompl
 Sparse representation of Schur complement matrix at this level. More...
 

Static Private Attributes

static std::string storagedef
 Utility strings. More...
 
static std::string storagelookupIds []
 

Additional Inherited Members

- Protected Attributes inherited from Nektar::MultiRegions::GlobalLinSysIterative
Array< OneD, int > m_map
 Global to universal unique map. More...
 
int m_maxiter
 maximum iterations More...
 
NekDouble m_tolerance
 Tolerance of iterative solver. More...
 
NekDouble m_rhs_magnitude
 dot product of rhs to normalise stopping criterion More...
 
NekDouble m_rhs_mag_sm
 cnt to how many times rhs_magnitude is called More...
 
PreconditionerSharedPtr m_precon
 
MultiRegions::PreconditionerType m_precontype
 
int m_totalIterations
 
bool m_useProjection
 Whether to apply projection technique. More...
 
bool m_root
 Provide verbose output and root if parallel. More...
 
bool m_verbose
 
boost::circular_buffer< Array
< OneD, NekDouble > > 
m_prevLinSol
 Storage for solutions to previous linear problems. More...
 
int m_numPrevSols
 Total counter of previous solutions. More...
 
- Protected Attributes inherited from Nektar::MultiRegions::GlobalLinSys
const GlobalLinSysKey m_linSysKey
 Key associated with this linear system. More...
 
const boost::weak_ptr< ExpListm_expList
 Local Matrix System. More...
 
const map< int,
RobinBCInfoSharedPtr
m_robinBCInfo
 Robin boundary info. More...
 
- Protected Attributes inherited from Nektar::MultiRegions::GlobalLinSysStaticCond
GlobalLinSysStaticCondSharedPtr m_recursiveSchurCompl
 Schur complement for Direct Static Condensation. More...
 
DNekScalBlkMatSharedPtr m_schurCompl
 Block Schur complement matrix. More...
 
DNekScalBlkMatSharedPtr m_BinvD
 Block $ BD^{-1} $ matrix. More...
 
DNekScalBlkMatSharedPtr m_C
 Block $ C $ matrix. More...
 
DNekScalBlkMatSharedPtr m_invD
 Block $ D^{-1} $ matrix. More...
 
boost::shared_ptr< AssemblyMapm_locToGloMap
 Local to global map. More...
 
Array< OneD, NekDoublem_wsp
 Workspace array for matrix multiplication. More...
 

Detailed Description

A global linear system.

Solves a linear system iteratively using single- or multi-level static condensation.

Definition at line 72 of file GlobalLinSysIterativeStaticCond.h.

Member Typedef Documentation

Definition at line 77 of file GlobalLinSysIterativeStaticCond.h.

Definition at line 79 of file GlobalLinSysIterativeStaticCond.h.

Constructor & Destructor Documentation

Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::GlobalLinSysIterativeStaticCond ( const GlobalLinSysKey pKey,
const boost::weak_ptr< ExpList > &  pExpList,
const boost::shared_ptr< AssemblyMap > &  pLocToGloMap 
)

Constructor for full direct matrix solve.

For a matrix system of the form

\[ \left[ \begin{array}{cc} \boldsymbol{A} & \boldsymbol{B}\\ \boldsymbol{C} & \boldsymbol{D} \end{array} \right] \left[ \begin{array}{c} \boldsymbol{x_1}\\ \boldsymbol{x_2} \end{array}\right] = \left[ \begin{array}{c} \boldsymbol{y_1}\\ \boldsymbol{y_2} \end{array}\right], \]

where $\boldsymbol{D}$ and $(\boldsymbol{A-BD^{-1}C})$ are invertible, store and assemble a static condensation system, according to a given local to global mapping. #m_linSys is constructed by AssembleSchurComplement().

Parameters
mKeyAssociated matrix key.
pLocMatSysLocalMatrixSystem
locToGloMapLocal to global mapping.

Definition at line 109 of file GlobalLinSysIterativeStaticCond.cpp.

References ASSERTL1, Nektar::MultiRegions::eIterativeMultiLevelStaticCond, Nektar::MultiRegions::eIterativeStaticCond, and Nektar::MultiRegions::GlobalLinSysKey::GetGlobalSysSolnType().

113  : GlobalLinSys (pKey, pExpList, pLocToGloMap),
114  GlobalLinSysIterative (pKey, pExpList, pLocToGloMap),
115  GlobalLinSysStaticCond(pKey, pExpList, pLocToGloMap)
116  {
117  ASSERTL1((pKey.GetGlobalSysSolnType()==eIterativeStaticCond)||
118  (pKey.GetGlobalSysSolnType()==eIterativeMultiLevelStaticCond),
119  "This constructor is only valid when using static "
120  "condensation");
121  ASSERTL1(pKey.GetGlobalSysSolnType()
122  == pLocToGloMap->GetGlobalSysSolnType(),
123  "The local to global map is not set up for the requested "
124  "solution type");
125  }
GlobalLinSysIterative(const GlobalLinSysKey &pKey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
Constructor for full direct matrix solve.
GlobalLinSys(const GlobalLinSysKey &pKey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
Constructor for full direct matrix solve.
GlobalLinSysStaticCond(const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &locToGloMap)
Constructor for full direct matrix solve.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::GlobalLinSysIterativeStaticCond ( const GlobalLinSysKey mkey,
const boost::weak_ptr< ExpList > &  pExpList,
const DNekScalBlkMatSharedPtr  pSchurCompl,
const DNekScalBlkMatSharedPtr  pBinvD,
const DNekScalBlkMatSharedPtr  pC,
const DNekScalBlkMatSharedPtr  pInvD,
const boost::shared_ptr< AssemblyMap > &  locToGloMap,
const PreconditionerSharedPtr  pPrecon 
)

Constructor for full direct matrix solve.

Definition at line 131 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MultiRegions::GlobalLinSysStaticCond::m_BinvD, Nektar::MultiRegions::GlobalLinSysStaticCond::m_C, Nektar::MultiRegions::GlobalLinSysStaticCond::m_invD, Nektar::MultiRegions::GlobalLinSysIterative::m_precon, m_S1Blk, and Nektar::MultiRegions::GlobalLinSysStaticCond::m_schurCompl.

140  : GlobalLinSys (pKey, pExpList, pLocToGloMap),
141  GlobalLinSysIterative (pKey, pExpList, pLocToGloMap),
142  GlobalLinSysStaticCond(pKey, pExpList, pLocToGloMap)
143  {
144  m_schurCompl = pSchurCompl;
145  m_S1Blk = pSchurCompl;
146  m_BinvD = pBinvD;
147  m_C = pC;
148  m_invD = pInvD;
149  m_precon = pPrecon;
150  }
DNekScalBlkMatSharedPtr m_schurCompl
Block Schur complement matrix.
DNekScalBlkMatSharedPtr m_C
Block matrix.
DNekScalBlkMatSharedPtr m_invD
Block matrix.
GlobalLinSysIterative(const GlobalLinSysKey &pKey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
Constructor for full direct matrix solve.
GlobalLinSys(const GlobalLinSysKey &pKey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
Constructor for full direct matrix solve.
GlobalLinSysStaticCond(const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &locToGloMap)
Constructor for full direct matrix solve.
DNekScalBlkMatSharedPtr m_BinvD
Block matrix.
Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::~GlobalLinSysIterativeStaticCond ( )
virtual

Definition at line 200 of file GlobalLinSysIterativeStaticCond.cpp.

201  {
202 
203  }

Member Function Documentation

static GlobalLinSysSharedPtr Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::create ( const GlobalLinSysKey pLinSysKey,
const boost::weak_ptr< ExpList > &  pExpList,
const boost::shared_ptr< AssemblyMap > &  pLocToGloMap 
)
inlinestatic

Creates an instance of this class.

Definition at line 82 of file GlobalLinSysIterativeStaticCond.h.

86  {
89  AllocateSharedPtr(pLinSysKey, pExpList, pLocToGloMap);
90  p->InitObject();
91  return p;
92  }
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
GlobalLinSysIterativeStaticCond(const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &locToGloMap)
Constructor for full direct matrix solve.
boost::shared_ptr< GlobalLinSys > GlobalLinSysSharedPtr
Pointer to a GlobalLinSys object.
Definition: GlobalLinSys.h:52
void Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::PrepareLocalSchurComplement ( )
private

Prepares local representation of Schur complement stored as a sparse block-diagonal matrix.

Populates sparse block-diagonal schur complement matrix from the block matrices stored in #m_blkMatrices.

Definition at line 321 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL1, Nektar::MultiRegions::eContiguous, Nektar::eFULL, Nektar::MultiRegions::eNonContiguous, Nektar::MultiRegions::eSparse, Nektar::NekConstants::kNekZeroTol, m_denseBlocks, Nektar::MultiRegions::GlobalLinSys::m_expList, m_rows, m_scale, Nektar::MultiRegions::GlobalLinSysStaticCond::m_schurCompl, m_sparseSchurCompl, m_storage, Vmath::Smul(), and Nektar::MultiRegions::GlobalLinSys::v_DropStaticCondBlock().

Referenced by v_AssembleSchurComplement().

322  {
323  LocalMatrixStorageStrategy storageStrategy =
324  m_expList.lock()->GetSession()->
325  GetSolverInfoAsEnum<LocalMatrixStorageStrategy>(
326  "LocalMatrixStorageStrategy");
327 
328  switch(storageStrategy)
329  {
332  {
333  size_t storageSize = 0;
334  int nBlk = m_schurCompl->GetNumberOfBlockRows();
335 
336  m_scale = Array<OneD, NekDouble> (nBlk, 1.0);
337  m_rows = Array<OneD, unsigned int> (nBlk, 0U);
338 
339  // Determine storage requirements for dense blocks.
340  for (int i = 0; i < nBlk; ++i)
341  {
342  m_rows[i] = m_schurCompl->GetBlock(i,i)->GetRows();
343  m_scale[i] = m_schurCompl->GetBlock(i,i)->Scale();
344  storageSize += m_rows[i] * m_rows[i];
345  }
346 
347  // Assemble dense storage blocks.
348  DNekScalMatSharedPtr loc_mat;
349  m_denseBlocks.resize(nBlk);
350  double *ptr = 0;
351 
352  if (MultiRegions::eContiguous == storageStrategy)
353  {
354  m_storage.resize (storageSize);
355  ptr = &m_storage[0];
356  }
357 
358  for (unsigned int n = 0; n < nBlk; ++n)
359  {
360  loc_mat = m_schurCompl->GetBlock(n,n);
361 
362  if (MultiRegions::eContiguous == storageStrategy)
363  {
364  int loc_lda = loc_mat->GetRows();
365  int blockSize = loc_lda * loc_lda;
366  m_denseBlocks[n] = ptr;
367  for(int i = 0; i < loc_lda; ++i)
368  {
369  for(int j = 0; j < loc_lda; ++j)
370  {
371  ptr[j*loc_lda+i] = (*loc_mat)(i,j);
372  }
373  }
374  ptr += blockSize;
376  m_expList.lock()->GetOffset_Elmt_Id(n));
377  }
378  else
379  {
380  m_denseBlocks[n] = loc_mat->GetRawPtr();
381  }
382  }
383  break;
384  }
386  {
387  DNekScalMatSharedPtr loc_mat;
388  int loc_lda;
389  int blockSize = 0;
390 
391  // First run through to split the set of local matrices into
392  // partitions of fixed block size, and count number of local
393  // matrices that belong to each partition.
394  std::vector<std::pair<int,int> > partitions;
395  for(int n = 0; n < m_schurCompl->GetNumberOfBlockRows(); ++n)
396  {
397  loc_mat = m_schurCompl->GetBlock(n,n);
398  loc_lda = loc_mat->GetRows();
399 
400  ASSERTL1(loc_lda >= 0,
401  boost::lexical_cast<std::string>(n) + "-th "
402  "matrix block in Schur complement has "
403  "rank 0!");
404 
405  if (blockSize == loc_lda)
406  {
407  partitions[partitions.size()-1].first++;
408  }
409  else
410  {
411  blockSize = loc_lda;
412  partitions.push_back(make_pair(1,loc_lda));
413  }
414  }
415 
416  MatrixStorage matStorage = eFULL;
417 
418  // Create a vector of sparse storage holders
420  sparseStorage (partitions.size());
421 
422  for (int part = 0, n = 0; part < partitions.size(); ++part)
423  {
424  BCOMatType partMat;
425 
426  for(int k = 0; k < partitions[part].first; ++k, ++n)
427  {
428  loc_mat = m_schurCompl->GetBlock(n,n);
429  loc_lda = loc_mat->GetRows();
430 
431  ASSERTL1(loc_lda == partitions[part].second,
432  boost::lexical_cast<std::string>(n) + "-th"
433  " matrix block in Schur complement has "
434  "unexpected rank");
435 
436  NekDouble scale = loc_mat->Scale();
437  if(fabs(scale-1.0) > NekConstants::kNekZeroTol)
438  {
439  Array<OneD, NekDouble> matarray(loc_lda*loc_lda);
440  Vmath::Smul(loc_lda*loc_lda,scale,
441  loc_mat->GetRawPtr(),1,&matarray[0],1);
442  partMat[make_pair(k,k)] = BCOEntryType(matarray);
443  }
444  else // scale factor is 1.0
445  {
446  partMat[make_pair(k,k)] = BCOEntryType(
447  loc_lda*loc_lda, loc_mat->GetRawPtr());
448  }
449 
451  m_expList.lock()->GetOffset_Elmt_Id(n));
452  }
453 
454  sparseStorage[part] =
457  partitions[part].first, partitions[part].first,
458  partitions[part].second, partMat, matStorage );
459  }
460 
461  // Create block diagonal matrix
463  AllocateSharedPtr(sparseStorage);
464 
465  break;
466  }
467  default:
468  ErrorUtil::NekError("Solver info property \
469  LocalMatrixStorageStrategy takes values \
470  Contiguous, Non-contiguous and Sparse");
471  }
472  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
DNekScalBlkMatSharedPtr m_schurCompl
Block Schur complement matrix.
std::vector< double > m_storage
Dense storage for block Schur complement matrix.
boost::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
Array< OneD, NekDouble > BCOEntryType
static const NekDouble kNekZeroTol
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:199
Array< OneD, NekDouble > m_scale
Scaling factors for local matrices.
std::map< CoordType, BCOEntryType > BCOMatType
Array< OneD, unsigned int > m_rows
Ranks of local matrices.
double NekDouble
std::vector< const double * > m_denseBlocks
Vector of pointers to local matrix data.
DNekSmvBsrDiagBlkMatSharedPtr m_sparseSchurCompl
Sparse representation of Schur complement matrix at this level.
virtual void v_DropStaticCondBlock(unsigned int n)
Releases the static condensation block matrices from NekManager of n-th expansion using the matrix ke...
Array< OneD, SparseStorageSharedPtr > SparseStorageSharedPtrVector
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
const boost::weak_ptr< ExpList > m_expList
Local Matrix System.
Definition: GlobalLinSys.h:129
void Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_AssembleSchurComplement ( const boost::shared_ptr< AssemblyMap locToGloMap)
privatevirtual

Assemble the Schur complement matrix.

Assemble the schur complement matrix from the block matrices stored in #m_blkMatrices and the given local to global mapping information.

Parameters
locToGloMapLocal to global mapping information.

Reimplemented from Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 228 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::convertCooToBco(), Nektar::MultiRegions::GlobalLinSys::CreatePrecon(), Nektar::eFULL, Nektar::MultiRegions::eIterativeMultiLevelStaticCond, Nektar::MultiRegions::GlobalLinSysKey::GetGlobalSysSolnType(), Nektar::MultiRegions::GlobalMatrixKey::GetMatrixType(), Nektar::MultiRegions::GlobalLinSys::m_expList, Nektar::MultiRegions::GlobalLinSys::m_linSysKey, Nektar::MultiRegions::GlobalLinSysStaticCond::m_locToGloMap, Nektar::MultiRegions::GlobalLinSysIterative::m_precon, Nektar::MultiRegions::GlobalLinSysStaticCond::m_schurCompl, m_sparseSchurCompl, PrepareLocalSchurComplement(), and v_UniqueMap().

230  {
231  int i,j,n,cnt,gid1,gid2;
232  NekDouble sign1,sign2;
233 
234  bool doGlobalOp = m_expList.lock()->GetGlobalOptParam()->
235  DoGlobalMatOp(m_linSysKey.GetMatrixType());
236 
237  // Set up unique map
238  v_UniqueMap();
239 
240  // Build precon again if we in multi-level static condensation (a
241  // bit of a hack)
244  {
246  m_precon->BuildPreconditioner();
247  }
248 
249  if (!doGlobalOp)
250  {
252  return;
253  }
254 
255  int nBndDofs = pLocToGloMap->GetNumGlobalBndCoeffs();
256  int NumDirBCs = pLocToGloMap->GetNumGlobalDirBndCoeffs();
257  unsigned int rows = nBndDofs - NumDirBCs;
258  unsigned int cols = nBndDofs - NumDirBCs;
259 
260  // COO sparse storage to assist in assembly
261  COOMatType gmat_coo;
262 
263  // Get the matrix storage structure
264  // (whether to store only one triangular part, if symmetric)
265  MatrixStorage matStorage = eFULL;
266 
267  // assemble globally
268  DNekScalMatSharedPtr loc_mat;
269  int loc_lda;
270  for(n = cnt = 0; n < m_schurCompl->GetNumberOfBlockRows(); ++n)
271  {
272  loc_mat = m_schurCompl->GetBlock(n,n);
273  loc_lda = loc_mat->GetRows();
274 
275  // Set up Matrix;
276  for(i = 0; i < loc_lda; ++i)
277  {
278  gid1 = pLocToGloMap->GetLocalToGlobalBndMap (cnt + i)
279  - NumDirBCs;
280  sign1 = pLocToGloMap->GetLocalToGlobalBndSign(cnt + i);
281 
282  if(gid1 >= 0)
283  {
284  for(j = 0; j < loc_lda; ++j)
285  {
286  gid2 = pLocToGloMap->GetLocalToGlobalBndMap(cnt+j)
287  - NumDirBCs;
288  sign2 = pLocToGloMap->GetLocalToGlobalBndSign(cnt+j);
289 
290  if (gid2 >= 0)
291  {
292  gmat_coo[std::make_pair(gid1,gid2)] +=
293  sign1*sign2*(*loc_mat)(i,j);
294  }
295  }
296  }
297  }
298  cnt += loc_lda;
299  }
300 
302  sparseStorage (1);
303 
304  BCOMatType partMat;
305  convertCooToBco(rows, cols, 1, gmat_coo, partMat);
306 
307  sparseStorage[0] =
309  AllocateSharedPtr(rows, cols, 1, partMat, matStorage );
310 
311  // Create block diagonal matrix
313  AllocateSharedPtr(sparseStorage);
314  }
GlobalSysSolnType GetGlobalSysSolnType() const
Return the associated solution type.
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
std::map< CoordType, NekDouble > COOMatType
DNekScalBlkMatSharedPtr m_schurCompl
Block Schur complement matrix.
boost::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::map< CoordType, BCOEntryType > BCOMatType
double NekDouble
void PrepareLocalSchurComplement()
Prepares local representation of Schur complement stored as a sparse block-diagonal matrix...
StdRegions::MatrixType GetMatrixType() const
Return the matrix type.
const GlobalLinSysKey m_linSysKey
Key associated with this linear system.
Definition: GlobalLinSys.h:127
DNekSmvBsrDiagBlkMatSharedPtr m_sparseSchurCompl
Sparse representation of Schur complement matrix at this level.
Array< OneD, SparseStorageSharedPtr > SparseStorageSharedPtrVector
void convertCooToBco(const unsigned int blkRows, const unsigned int blkColumns, const unsigned int blkDim, const COOMatType &cooMat, BCOMatType &bcoMat)
Definition: SparseUtils.cpp:48
boost::shared_ptr< AssemblyMap > m_locToGloMap
Local to global map.
PreconditionerSharedPtr CreatePrecon(AssemblyMapSharedPtr asmMap)
Create a preconditioner object from the parameters defined in the supplied assembly map...
const boost::weak_ptr< ExpList > m_expList
Local Matrix System.
Definition: GlobalLinSys.h:129
void Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_BasisInvTransform ( Array< OneD, NekDouble > &  pInOut)
protectedvirtual

Reimplemented from Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 562 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MultiRegions::GlobalLinSysIterative::m_precon.

564  {
565  m_precon->DoTransformFromLowEnergy(pInOut);
566  }
void Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_BasisTransform ( Array< OneD, NekDouble > &  pInOut,
int  offset 
)
protectedvirtual

Reimplemented from Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 555 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MultiRegions::GlobalLinSysIterative::m_precon.

558  {
559  m_precon->DoTransformToLowEnergy(pInOut, offset);
560  }
void Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_DoMatrixMultiply ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput 
)
privatevirtual

Perform a Shur-complement matrix multiply operation.

Implements Nektar::MultiRegions::GlobalLinSysIterative.

Definition at line 477 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MultiRegions::GlobalMatrixKey::GetMatrixType(), m_denseBlocks, Nektar::MultiRegions::GlobalLinSys::m_expList, Nektar::MultiRegions::GlobalLinSys::m_linSysKey, Nektar::MultiRegions::GlobalLinSysStaticCond::m_locToGloMap, m_rows, m_scale, m_sparseSchurCompl, and Nektar::MultiRegions::GlobalLinSysStaticCond::m_wsp.

480  {
481  int nLocal = m_locToGloMap->GetNumLocalBndCoeffs();
482  int nDir = m_locToGloMap->GetNumGlobalDirBndCoeffs();
483  bool doGlobalOp = m_expList.lock()->GetGlobalOptParam()->
484  DoGlobalMatOp(m_linSysKey.GetMatrixType());
485 
486  if(doGlobalOp)
487  {
488  // Do matrix multiply globally
489  Array<OneD, NekDouble> in = pInput + nDir;
490  Array<OneD, NekDouble> out = pOutput + nDir;
491 
492  m_sparseSchurCompl->Multiply(in,out);
493  m_locToGloMap->UniversalAssembleBnd(pOutput, nDir);
494  }
495  else if (m_sparseSchurCompl)
496  {
497  // Do matrix multiply locally using block-diagonal sparse matrix
498  Array<OneD, NekDouble> tmp = m_wsp + nLocal;
499 
500  m_locToGloMap->GlobalToLocalBnd(pInput, m_wsp);
501  m_sparseSchurCompl->Multiply(m_wsp,tmp);
502  m_locToGloMap->AssembleBnd(tmp, pOutput);
503  }
504  else
505  {
506  // Do matrix multiply locally, using direct BLAS calls
507  m_locToGloMap->GlobalToLocalBnd(pInput, m_wsp);
508  int i, cnt;
509  Array<OneD, NekDouble> tmpout = m_wsp + nLocal;
510  for (i = cnt = 0; i < m_denseBlocks.size(); cnt += m_rows[i], ++i)
511  {
512  const int rows = m_rows[i];
513  Blas::Dgemv('N', rows, rows,
514  m_scale[i], m_denseBlocks[i], rows,
515  m_wsp.get()+cnt, 1,
516  0.0, tmpout.get()+cnt, 1);
517  }
518  m_locToGloMap->AssembleBnd(tmpout, pOutput);
519  }
520  }
Array< OneD, NekDouble > m_wsp
Workspace array for matrix multiplication.
Array< OneD, NekDouble > m_scale
Scaling factors for local matrices.
Array< OneD, unsigned int > m_rows
Ranks of local matrices.
std::vector< const double * > m_denseBlocks
Vector of pointers to local matrix data.
StdRegions::MatrixType GetMatrixType() const
Return the matrix type.
const GlobalLinSysKey m_linSysKey
Key associated with this linear system.
Definition: GlobalLinSys.h:127
DNekSmvBsrDiagBlkMatSharedPtr m_sparseSchurCompl
Sparse representation of Schur complement matrix at this level.
boost::shared_ptr< AssemblyMap > m_locToGloMap
Local to global map.
const boost::weak_ptr< ExpList > m_expList
Local Matrix System.
Definition: GlobalLinSys.h:129
DNekScalBlkMatSharedPtr Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_GetStaticCondBlock ( unsigned int  n)
protectedvirtual

Retrieves a the static condensation block matrices from n-th expansion using the matrix key provided by the m_linSysKey.

Parameters
nNumber of the expansion
Returns
2x2 Block matrix holding the static condensation matrices for the n-th expansion.

Reimplemented from Nektar::MultiRegions::GlobalLinSys.

Definition at line 206 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::MultiRegions::GlobalLinSysStaticCond::m_locToGloMap, m_S1Blk, and Nektar::MultiRegions::GlobalLinSysStaticCond::m_schurCompl.

207  {
208  DNekScalBlkMatSharedPtr schurComplBlock;
209  int scLevel = m_locToGloMap->GetStaticCondLevel();
210  DNekScalBlkMatSharedPtr sc = scLevel == 0 ? m_S1Blk : m_schurCompl;
211  DNekScalMatSharedPtr localMat = sc->GetBlock(n,n);
212  unsigned int nbdry = localMat->GetRows();
213  unsigned int nblks = 1;
214  unsigned int esize[1] = {nbdry};
215 
216  schurComplBlock = MemoryManager<DNekScalBlkMat>
217  ::AllocateSharedPtr(nblks, nblks, esize, esize);
218  schurComplBlock->SetBlock(0, 0, localMat);
219 
220  return schurComplBlock;
221  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
DNekScalBlkMatSharedPtr m_schurCompl
Block Schur complement matrix.
boost::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
boost::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:74
boost::shared_ptr< AssemblyMap > m_locToGloMap
Local to global map.
void Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_InitObject ( )
privatevirtual

Reimplemented from Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 153 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::MultiRegions::GlobalLinSys::CreatePrecon(), Nektar::eDIAGONAL, Nektar::StdRegions::eHybridDGHelmBndLam, Nektar::MultiRegions::GlobalMatrixKey::GetMatrixType(), Nektar::MultiRegions::GlobalLinSys::Initialise(), Nektar::MultiRegions::GlobalLinSys::m_expList, Nektar::MultiRegions::GlobalLinSys::m_linSysKey, Nektar::MultiRegions::GlobalLinSysStaticCond::m_locToGloMap, Nektar::MultiRegions::GlobalLinSysIterative::m_precon, m_S1Blk, Nektar::MultiRegions::GlobalLinSysStaticCond::m_schurCompl, and Nektar::MultiRegions::GlobalLinSysStaticCond::SetupTopLevel().

154  {
156 
157  // Allocate memory for top-level structure
159 
160  // Setup Block Matrix systems
161  int n, n_exp = m_expList.lock()->GetNumElmts();
162 
163  MatrixStorage blkmatStorage = eDIAGONAL;
164  const Array<OneD,const unsigned int>& nbdry_size
165  = m_locToGloMap->GetNumLocalBndCoeffsPerPatch();
166 
168  ::AllocateSharedPtr(nbdry_size, nbdry_size, blkmatStorage);
169 
170  // Preserve original matrix in m_S1Blk
171  for (n = 0; n < n_exp; ++n)
172  {
173  DNekScalMatSharedPtr mat = m_schurCompl->GetBlock(n, n);
174  m_S1Blk->SetBlock(n, n, mat);
175  }
176 
177  // Build preconditioner
178  m_precon->BuildPreconditioner();
179 
180  // Do transform of Schur complement matrix
181  for (n = 0; n < n_exp; ++n)
182  {
183  if (m_linSysKey.GetMatrixType() !=
185  {
186  DNekScalMatSharedPtr mat = m_S1Blk->GetBlock(n, n);
187  DNekScalMatSharedPtr t = m_precon->TransformedSchurCompl(
188  m_expList.lock()->GetOffset_Elmt_Id(n), mat);
189  m_schurCompl->SetBlock(n, n, t);
190  }
191  }
192 
193  // Construct this level
195  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
DNekScalBlkMatSharedPtr m_schurCompl
Block Schur complement matrix.
void Initialise(const boost::shared_ptr< AssemblyMap > &pLocToGloMap)
Definition: GlobalLinSys.h:214
boost::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
void SetupTopLevel(const boost::shared_ptr< AssemblyMap > &locToGloMap)
Set up the storage for the Schur complement or the top level of the multi-level Schur complement...
StdRegions::MatrixType GetMatrixType() const
Return the matrix type.
const GlobalLinSysKey m_linSysKey
Key associated with this linear system.
Definition: GlobalLinSys.h:127
boost::shared_ptr< AssemblyMap > m_locToGloMap
Local to global map.
PreconditionerSharedPtr CreatePrecon(AssemblyMapSharedPtr asmMap)
Create a preconditioner object from the parameters defined in the supplied assembly map...
const boost::weak_ptr< ExpList > m_expList
Local Matrix System.
Definition: GlobalLinSys.h:129
DNekScalBlkMatSharedPtr Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_PreSolve ( int  scLevel,
NekVector< NekDouble > &  F_GlobBnd 
)
protectedvirtual

Reimplemented from Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 527 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MultiRegions::GlobalLinSys::CreatePrecon(), Nektar::NekConstants::kNekUnsetDouble, Nektar::MultiRegions::GlobalLinSysStaticCond::m_locToGloMap, Nektar::MultiRegions::GlobalLinSysIterative::m_precon, Nektar::MultiRegions::GlobalLinSysIterative::m_rhs_magnitude, m_S1Blk, Nektar::MultiRegions::GlobalLinSysStaticCond::m_schurCompl, and Nektar::MultiRegions::GlobalLinSysIterative::Set_Rhs_Magnitude().

530  {
531  if (scLevel == 0)
532  {
533  // When matrices are supplied to the constructor at the top
534  // level, the preconditioner is never set up.
535  if (!m_precon)
536  {
538  m_precon->BuildPreconditioner();
539  }
540 
541  Set_Rhs_Magnitude(F_GlobBnd);
542 
543  return m_S1Blk;
544  }
545  else
546  {
547  // for multilevel iterative solver always use rhs
548  // vector value with no weighting
550 
551  return m_schurCompl;
552  }
553  }
DNekScalBlkMatSharedPtr m_schurCompl
Block Schur complement matrix.
void Set_Rhs_Magnitude(const NekVector< NekDouble > &pIn)
static const NekDouble kNekUnsetDouble
NekDouble m_rhs_magnitude
dot product of rhs to normalise stopping criterion
boost::shared_ptr< AssemblyMap > m_locToGloMap
Local to global map.
PreconditionerSharedPtr CreatePrecon(AssemblyMapSharedPtr asmMap)
Create a preconditioner object from the parameters defined in the supplied assembly map...
GlobalLinSysStaticCondSharedPtr Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_Recurse ( const GlobalLinSysKey mkey,
const boost::weak_ptr< ExpList > &  pExpList,
const DNekScalBlkMatSharedPtr  pSchurCompl,
const DNekScalBlkMatSharedPtr  pBinvD,
const DNekScalBlkMatSharedPtr  pC,
const DNekScalBlkMatSharedPtr  pInvD,
const boost::shared_ptr< AssemblyMap > &  locToGloMap 
)
protectedvirtual

Implements Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 568 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MultiRegions::GlobalLinSysIterative::m_precon.

576  {
578  GlobalLinSysIterativeStaticCond>::AllocateSharedPtr(
579  mkey, pExpList, pSchurCompl, pBinvD, pC, pInvD, l2gMap,
580  m_precon);
581  sys->Initialise(l2gMap);
582  return sys;
583  }
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
GlobalLinSysIterativeStaticCond(const GlobalLinSysKey &mkey, const boost::weak_ptr< ExpList > &pExpList, const boost::shared_ptr< AssemblyMap > &locToGloMap)
Constructor for full direct matrix solve.
boost::shared_ptr< GlobalLinSysIterativeStaticCond > GlobalLinSysIterativeStaticCondSharedPtr
void Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_UniqueMap ( )
privatevirtual

Implements Nektar::MultiRegions::GlobalLinSysIterative.

Definition at line 522 of file GlobalLinSysIterativeStaticCond.cpp.

References Nektar::MultiRegions::GlobalLinSysStaticCond::m_locToGloMap, and Nektar::MultiRegions::GlobalLinSysIterative::m_map.

Referenced by v_AssembleSchurComplement().

523  {
524  m_map = m_locToGloMap->GetGlobalToUniversalBndMapUnique();
525  }
Array< OneD, int > m_map
Global to universal unique map.
boost::shared_ptr< AssemblyMap > m_locToGloMap
Local to global map.

Member Data Documentation

string Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::className
static
Initial value:
"IterativeStaticCond",
"Iterative static condensation.")

Name of class.

Registers the class with the Factory.

Definition at line 95 of file GlobalLinSysIterativeStaticCond.h.

string Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::className2
static
Initial value:
"IterativeMultiLevelStaticCond",
"Iterative multi-level static condensation.")

Definition at line 96 of file GlobalLinSysIterativeStaticCond.h.

std::vector<const double*> Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::m_denseBlocks
private

Vector of pointers to local matrix data.

Definition at line 142 of file GlobalLinSysIterativeStaticCond.h.

Referenced by PrepareLocalSchurComplement(), and v_DoMatrixMultiply().

Array<OneD, unsigned int> Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::m_rows
private

Ranks of local matrices.

Definition at line 144 of file GlobalLinSysIterativeStaticCond.h.

Referenced by PrepareLocalSchurComplement(), and v_DoMatrixMultiply().

DNekScalBlkMatSharedPtr Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::m_S1Blk
private
Array<OneD, NekDouble> Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::m_scale
private

Scaling factors for local matrices.

Definition at line 146 of file GlobalLinSysIterativeStaticCond.h.

Referenced by PrepareLocalSchurComplement(), and v_DoMatrixMultiply().

DNekSmvBsrDiagBlkMatSharedPtr Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::m_sparseSchurCompl
private

Sparse representation of Schur complement matrix at this level.

Definition at line 148 of file GlobalLinSysIterativeStaticCond.h.

Referenced by PrepareLocalSchurComplement(), v_AssembleSchurComplement(), and v_DoMatrixMultiply().

std::vector<double> Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::m_storage
private

Dense storage for block Schur complement matrix.

Definition at line 140 of file GlobalLinSysIterativeStaticCond.h.

Referenced by PrepareLocalSchurComplement().

std::string Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::storagedef
staticprivate
Initial value:
=
"LocalMatrixStorageStrategy",
"Sparse")

Utility strings.

Definition at line 150 of file GlobalLinSysIterativeStaticCond.h.

std::string Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::storagelookupIds
staticprivate
Initial value:
= {
"LocalMatrixStorageStrategy",
"Contiguous",
"LocalMatrixStorageStrategy",
"Non-contiguous",
"LocalMatrixStorageStrategy",
"Sparse",
}

Definition at line 151 of file GlobalLinSysIterativeStaticCond.h.