Nektar++
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Static Private Attributes | List of all members
Nektar::MultiRegions::GlobalLinSys Class Referenceabstract

A global linear system. More...

#include <GlobalLinSys.h>

Inheritance diagram for Nektar::MultiRegions::GlobalLinSys:
[legend]

Public Member Functions

 GlobalLinSys (const GlobalLinSysKey &pKey, const std::weak_ptr< ExpList > &pExpList, const std::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 std::weak_ptr< ExpList > & GetLocMat (void) const
 
void InitObject ()
 
void Initialise (const std::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...
 
std::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...
 

Protected Member Functions

virtual int v_GetNumBlocks ()
 Get the number of blocks in this system. More...
 
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 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 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 Attributes

const GlobalLinSysKey m_linSysKey
 Key associated with this linear system. More...
 
const std::weak_ptr< ExpListm_expList
 Local Matrix System. More...
 
const std::map< int, RobinBCInfoSharedPtrm_robinBCInfo
 Robin boundary info. More...
 
bool m_verbose
 

Private Member Functions

LocalRegions::MatrixKey GetBlockMatrixKey (unsigned int n)
 
virtual void v_Solve (const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out, const AssemblyMapSharedPtr &locToGloMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)=0
 Solve a linear system based on mapping. More...
 
virtual void v_SolveLinearSystem (const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir)=0
 Solve a basic matrix system. More...
 
virtual void v_InitObject ()
 
virtual void v_Initialise (const std::shared_ptr< AssemblyMap > &pLocToGloMap)
 

Static Private Attributes

static std::string lookupIds []
 
static std::string def
 

Detailed Description

A global linear system.

Consider the linear system \(\boldsymbol{M\hat{u}}_g=\boldsymbol{\hat{f}}\). Distinguishing between the boundary and interior components of \(\boldsymbol{\hat{u}}_g\) and \(\boldsymbol{\hat{f}}\) using \(\boldsymbol{\hat{u}}_b\), \(\boldsymbol{\hat{u}}_i\) and \(\boldsymbol{\hat{f}}_b\), \(\boldsymbol{\hat{f}}_i\) respectively, this system can be split into its constituent parts as

\[\left[\begin{array}{cc} \boldsymbol{M}_b&\boldsymbol{M}_{c1}\\ \boldsymbol{M}_{c2}&\boldsymbol{M}_i\\ \end{array}\right] \left[\begin{array}{c} \boldsymbol{\hat{u}_b}\\ \boldsymbol{\hat{u}_i}\\ \end{array}\right]= \left[\begin{array}{c} \boldsymbol{\hat{f}_b}\\ \boldsymbol{\hat{f}_i}\\ \end{array}\right]\]

where \(\boldsymbol{M}_b\) represents the components of \(\boldsymbol{M}\) resulting from boundary-boundary mode interactions, \(\boldsymbol{M}_{c1}\) and \(\boldsymbol{M}_{c2}\) represent the components resulting from coupling between the boundary-interior modes, and \(\boldsymbol{M}_i\) represents the components of \(\boldsymbol{M}\) resulting from interior-interior mode interactions.

The solution of the linear system can now be determined in two steps:

\begin{eqnarray*} \mathrm{step 1:}&\quad&(\boldsymbol{M}_b-\boldsymbol{M}_{c1} \boldsymbol{M}_i^{-1}\boldsymbol{M}_{c2}) \boldsymbol{\hat{u}_b} = \boldsymbol{\hat{f}}_b - \boldsymbol{M}_{c1}\boldsymbol{M}_i^{-1} \boldsymbol{\hat{f}}_i,\nonumber \\ \mathrm{step 2:}&\quad&\boldsymbol{\hat{u}_i}=\boldsymbol{M}_i^{-1} \left( \boldsymbol{\hat{f}}_i - \boldsymbol{M}_{c2}\boldsymbol{\hat{u}_b} \right). \nonumber \\ \end{eqnarray*}

As the inverse of \(\boldsymbol{M}_i^{-1}\) is

\[ \boldsymbol{M}_i^{-1} = \left [\underline{\boldsymbol{M}^e_i} \right ]^{-1} = \underline{[\boldsymbol{M}^e_i]}^{-1} \]

and the following operations can be evaluated as,

\begin{eqnarray*} \boldsymbol{M}_{c1}\boldsymbol{M}_i^{-1}\boldsymbol{\hat{f}}_i & =& \boldsymbol{\mathcal{A}}_b^T \underline{\boldsymbol{M}^e_{c1}} \underline{[\boldsymbol{M}^e_i]}^{-1} \boldsymbol{\hat{f}}_i \\ \boldsymbol{M}_{c2} \boldsymbol{\hat{u}_b} &=& \underline{\boldsymbol{M}^e_{c2}} \boldsymbol{\mathcal{A}}_b \boldsymbol{\hat{u}_b}.\end{eqnarray*}

where \(\boldsymbol{\mathcal{A}}_b \) is the permutation matrix which scatters from global to local degrees of freedom, only the following four matrices should be constructed:

The first three matrices are just a concatenation of the corresponding local matrices and they can be created as such. They also allow for an elemental evaluation of the operations concerned.

The global Schur complement however should be assembled from the concatenation of the local elemental Schur complements, that is,

\[ \boldsymbol{M}_{\mathrm{Schur}}=\boldsymbol{M}_b - \boldsymbol{M}_{c1} \boldsymbol{M}_i^{-1} \boldsymbol{M}_{c2} = \boldsymbol{\mathcal{A}}_b^T \left [\underline{\boldsymbol{M}^e_b - \boldsymbol{M}^e_{c1} [\boldsymbol{M}^e_i]^{-1} (\boldsymbol{M}^e_{c2})} \right ] \boldsymbol{\mathcal{A}}_b \]

and it is the only matrix operation that need to be evaluated on a global level when using static condensation. However, due to the size and sparsity of the matrix \(\boldsymbol{\mathcal{A}}_b\), it is more efficient to assemble the global Schur matrix using the mapping array bmap \([e][i]\) contained in the input argument locToGloMap. The global Schur complement is then constructed as:

\[\boldsymbol{M}_{\mathrm{Schur}}\left[\mathrm{\a bmap}[e][i]\right] \left[\mathrm{\a bmap}[e][j]\right]=\mathrm{\a bsign}[e][i]\cdot \mathrm{\a bsign}[e][j] \cdot\boldsymbol{M}^e_{\mathrm{Schur}}[i][j]\]

All four matrices are stored in the GlobalLinSys returned by this function.

Definition at line 71 of file GlobalLinSys.h.

Constructor & Destructor Documentation

◆ GlobalLinSys()

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

Constructor for full direct matrix solve.

Given a block matrix, construct a global matrix system according to a local to global mapping. #m_linSys is constructed by AssembleFullMatrix().

Parameters
pkeyAssociated linear system key.
locToGloMapLocal to global mapping.

Definition at line 188 of file GlobalLinSys.cpp.

191  : m_linSysKey(pKey), m_expList(pExpList),
192  m_robinBCInfo(m_expList.lock()->GetRobinBCInfo()),
193  m_verbose(
194  m_expList.lock()->GetSession()->DefinesCmdLineArgument("verbose"))
195 {
196  boost::ignore_unused(pLocToGloMap);
197 }
const std::weak_ptr< ExpList > m_expList
Local Matrix System.
Definition: GlobalLinSys.h:123
const std::map< int, RobinBCInfoSharedPtr > m_robinBCInfo
Robin boundary info.
Definition: GlobalLinSys.h:125
const GlobalLinSysKey m_linSysKey
Key associated with this linear system.
Definition: GlobalLinSys.h:121

◆ ~GlobalLinSys()

virtual Nektar::MultiRegions::GlobalLinSys::~GlobalLinSys ( )
inlinevirtual

Definition at line 80 of file GlobalLinSys.h.

81  {
82  }

Member Function Documentation

◆ CreatePrecon()

PreconditionerSharedPtr Nektar::MultiRegions::GlobalLinSys::CreatePrecon ( AssemblyMapSharedPtr  asmMap)
protected

Create a preconditioner object from the parameters defined in the supplied assembly map.

Parameters
asmMapAssembly map used to construct the global system.

Definition at line 214 of file GlobalLinSys.cpp.

215 {
216  PreconditionerType pType = asmMap->GetPreconType();
217  std::string PreconType = MultiRegions::PreconditionerTypeMap[pType];
218  return GetPreconFactory().CreateInstance(PreconType, GetSharedThisPtr(),
219  asmMap);
220 }
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:144
std::shared_ptr< GlobalLinSys > GetSharedThisPtr()
Returns a shared pointer to the current object.
Definition: GlobalLinSys.h:102
const char *const PreconditionerTypeMap[]
PreconFactory & GetPreconFactory()

References Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::MultiRegions::GetPreconFactory(), GetSharedThisPtr(), and Nektar::MultiRegions::PreconditionerTypeMap.

Referenced by Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_AssembleSchurComplement(), Nektar::MultiRegions::GlobalLinSysPETScStaticCond::v_AssembleSchurComplement(), Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_InitObject(), Nektar::MultiRegions::GlobalLinSysPETScStaticCond::v_InitObject(), Nektar::MultiRegions::GlobalLinSysPETScStaticCond::v_PreSolve(), Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_PreSolve(), Nektar::MultiRegions::GlobalLinSysIterative::v_SolveLinearSystem(), and Nektar::MultiRegions::GlobalLinSysPETSc::v_SolveLinearSystem().

◆ DropStaticCondBlock()

void Nektar::MultiRegions::GlobalLinSys::DropStaticCondBlock ( unsigned int  n)
inline

Definition at line 219 of file GlobalLinSys.h.

220 {
221  return v_DropStaticCondBlock(n);
222 }
virtual void v_DropStaticCondBlock(unsigned int n)
Releases the static condensation block matrices from NekManager of n-th expansion using the matrix ke...

References v_DropStaticCondBlock().

◆ GetBlock()

DNekScalMatSharedPtr Nektar::MultiRegions::GlobalLinSys::GetBlock ( unsigned int  n)
inline

Definition at line 209 of file GlobalLinSys.h.

210 {
211  return v_GetBlock(n);
212 }
virtual DNekScalMatSharedPtr v_GetBlock(unsigned int n)
Retrieves the block matrix from n-th expansion using the matrix key provided by the m_linSysKey.

References v_GetBlock().

Referenced by Nektar::MultiRegions::GlobalLinSysDirectFull::AssembleFullMatrix(), Nektar::MultiRegions::GlobalLinSysXxtFull::AssembleMatrixArrays(), and Nektar::MultiRegions::GlobalLinSysPETScFull::GlobalLinSysPETScFull().

◆ GetBlockMatrixKey()

LocalRegions::MatrixKey Nektar::MultiRegions::GlobalLinSys::GetBlockMatrixKey ( unsigned int  n)
private

Assemble the matrix key for each block n

Definition at line 237 of file GlobalLinSys.cpp.

238 {
239 
240  std::shared_ptr<MultiRegions::ExpList> expList = m_expList.lock();
241  int cnt = 0;
242 
243  LocalRegions::ExpansionSharedPtr vExp = expList->GetExp(n);
244 
245  // need to be initialised with zero size for non variable
246  // coefficient case
247  StdRegions::VarCoeffMap vVarCoeffMap;
248 
250 
251  // setup variable factors
252  if (m_linSysKey.GetNVarFactors() > 0)
253  {
254  if (m_linSysKey.GetVarFactors().count(
256  {
257  vConstFactorMap[StdRegions::eFactorSVVDiffCoeff] =
259 
262  "VarCoeffSVVCuroffRatio is set but "
263  " not FactorSVVCutoffRatio");
264 
265  vConstFactorMap[StdRegions::eFactorSVVCutoffRatio] =
267  }
268 
269  if (m_linSysKey.GetVarFactors().count(
271  {
272  vConstFactorMap[StdRegions::eFactorSVVPowerKerDiffCoeff] =
275  }
276 
277  if (m_linSysKey.GetVarFactors().count(
279  {
280  vConstFactorMap[StdRegions::eFactorSVVDGKerDiffCoeff] =
283  }
284  }
285 
286  // retrieve variable coefficients
287  if (m_linSysKey.GetNVarCoeffs() > 0)
288  {
289  cnt = expList->GetPhys_Offset(n);
290 
291  for (auto &x : m_linSysKey.GetVarCoeffs())
292  {
293  vVarCoeffMap[x.first] = x.second + cnt;
294  }
295  }
296 
297  LocalRegions::MatrixKey matkey(m_linSysKey.GetMatrixType(),
298  vExp->DetShapeType(), *vExp, vConstFactorMap,
299  vVarCoeffMap);
300  return matkey;
301 }
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:249
const Array< OneD, const NekDouble > & GetVarFactors(const StdRegions::ConstFactorType &coeff) const
const StdRegions::ConstFactorMap & GetConstFactors() const
Returns all the constants.
const StdRegions::VarCoeffMap & GetVarCoeffs() const
StdRegions::MatrixType GetMatrixType() const
Return the matrix type.
std::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:240
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:282

References ASSERTL1, Nektar::StdRegions::eFactorSVVCutoffRatio, Nektar::StdRegions::eFactorSVVDGKerDiffCoeff, Nektar::StdRegions::eFactorSVVDiffCoeff, Nektar::StdRegions::eFactorSVVPowerKerDiffCoeff, Nektar::MultiRegions::GlobalMatrixKey::GetConstFactors(), Nektar::MultiRegions::GlobalMatrixKey::GetMatrixType(), Nektar::MultiRegions::GlobalMatrixKey::GetNVarCoeffs(), Nektar::MultiRegions::GlobalLinSysKey::GetNVarFactors(), Nektar::MultiRegions::GlobalMatrixKey::GetVarCoeffs(), Nektar::MultiRegions::GlobalLinSysKey::GetVarFactors(), m_expList, and m_linSysKey.

Referenced by v_DropStaticCondBlock(), v_GetBlock(), and v_GetStaticCondBlock().

◆ GetKey()

const GlobalLinSysKey & Nektar::MultiRegions::GlobalLinSys::GetKey ( void  ) const
inline

Returns the key associated with the system.

Definition at line 163 of file GlobalLinSys.h.

164 {
165  return m_linSysKey;
166 }

References m_linSysKey.

◆ GetLocMat()

const std::weak_ptr< ExpList > & Nektar::MultiRegions::GlobalLinSys::GetLocMat ( void  ) const
inline

Definition at line 171 of file GlobalLinSys.h.

172 {
173  return m_expList;
174 }

References m_expList.

◆ GetNumBlocks()

int Nektar::MultiRegions::GlobalLinSys::GetNumBlocks ( )
inline

Definition at line 224 of file GlobalLinSys.h.

225 {
226  return v_GetNumBlocks();
227 }
virtual int v_GetNumBlocks()
Get the number of blocks in this system.

References v_GetNumBlocks().

◆ GetSharedThisPtr()

std::shared_ptr<GlobalLinSys> Nektar::MultiRegions::GlobalLinSys::GetSharedThisPtr ( )
inline

Returns a shared pointer to the current object.

Definition at line 102 of file GlobalLinSys.h.

103  {
104  return shared_from_this();
105  }

Referenced by CreatePrecon().

◆ GetStaticCondBlock()

DNekScalBlkMatSharedPtr Nektar::MultiRegions::GlobalLinSys::GetStaticCondBlock ( unsigned int  n)
inline

Definition at line 214 of file GlobalLinSys.h.

215 {
216  return v_GetStaticCondBlock(n);
217 }
virtual DNekScalBlkMatSharedPtr v_GetStaticCondBlock(unsigned int n)
Retrieves a the static condensation block matrices from n-th expansion using the matrix key provided ...

References v_GetStaticCondBlock().

◆ Initialise()

void Nektar::MultiRegions::GlobalLinSys::Initialise ( const std::shared_ptr< AssemblyMap > &  pLocToGloMap)
inline

◆ InitObject()

void Nektar::MultiRegions::GlobalLinSys::InitObject ( )
inline

Definition at line 198 of file GlobalLinSys.h.

199 {
200  v_InitObject();
201 }

References v_InitObject().

◆ Solve()

void Nektar::MultiRegions::GlobalLinSys::Solve ( const Array< OneD, const NekDouble > &  in,
Array< OneD, NekDouble > &  out,
const AssemblyMapSharedPtr locToGloMap,
const Array< OneD, const NekDouble > &  dirForcing = NullNekDouble1DArray 
)
inline

Solve the linear system for given input and output vectors using a specified local to global map.

Definition at line 179 of file GlobalLinSys.h.

183 {
184  v_Solve(in, out, locToGloMap, dirForcing);
185 }
virtual void v_Solve(const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out, const AssemblyMapSharedPtr &locToGloMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)=0
Solve a linear system based on mapping.

References v_Solve().

◆ SolveLinearSystem()

void Nektar::MultiRegions::GlobalLinSys::SolveLinearSystem ( const int  pNumRows,
const Array< OneD, const NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput,
const AssemblyMapSharedPtr locToGloMap,
const int  pNumDir = 0 
)
inline

Solve the linear system for given input and output vectors.

Definition at line 190 of file GlobalLinSys.h.

194 {
195  v_SolveLinearSystem(pNumRows, pInput, pOutput, locToGloMap, pNumDir);
196 }
virtual void v_SolveLinearSystem(const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir)=0
Solve a basic matrix system.

References v_SolveLinearSystem().

Referenced by Nektar::MultiRegions::GlobalLinSysIterativeFull::v_Solve(), Nektar::MultiRegions::GlobalLinSysPETScFull::v_Solve(), Nektar::MultiRegions::GlobalLinSysStaticCond::v_Solve(), Nektar::MultiRegions::GlobalLinSysXxtFull::v_Solve(), and Nektar::MultiRegions::GlobalLinSysDirectFull::v_Solve().

◆ v_DropStaticCondBlock()

void Nektar::MultiRegions::GlobalLinSys::v_DropStaticCondBlock ( unsigned int  n)
protectedvirtual

Releases the static condensation block matrices from NekManager of n-th expansion using the matrix key provided by the m_linSysKey.

Parameters
nNumber of the expansion

Definition at line 406 of file GlobalLinSys.cpp.

407 {
408  LocalRegions::ExpansionSharedPtr vExp = m_expList.lock()->GetExp(n);
409  vExp->DropLocStaticCondMatrix(GetBlockMatrixKey(n));
410 }
LocalRegions::MatrixKey GetBlockMatrixKey(unsigned int n)

References GetBlockMatrixKey(), and m_expList.

Referenced by DropStaticCondBlock(), and Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::PrepareLocalSchurComplement().

◆ v_GetBlock()

DNekScalMatSharedPtr Nektar::MultiRegions::GlobalLinSys::v_GetBlock ( unsigned int  n)
protectedvirtual

Retrieves the block matrix from n-th expansion using the matrix key provided by the m_linSysKey.

Parameters
nNumber of the expansion.
Returns
Block matrix for the specified expansion.

Definition at line 310 of file GlobalLinSys.cpp.

311 {
312  LocalRegions::ExpansionSharedPtr vExp = m_expList.lock()->GetExp(n);
313  DNekScalMatSharedPtr loc_mat;
314  loc_mat = vExp->GetLocMatrix(GetBlockMatrixKey(n));
315 
316  // apply robin boundary conditions to the matrix.
317  if (m_robinBCInfo.count(n) != 0) // add robin mass matrix
318  {
320 
321  // declare local matrix from scaled matrix.
322  int rows = loc_mat->GetRows();
323  int cols = loc_mat->GetColumns();
324  const NekDouble *dat = loc_mat->GetRawPtr();
325  DNekMatSharedPtr new_mat =
327  Blas::Dscal(rows * cols, loc_mat->Scale(), new_mat->GetRawPtr(), 1);
328 
329  // add local matrix contribution
330  for (rBC = m_robinBCInfo.find(n)->second; rBC; rBC = rBC->next)
331  {
332  vExp->AddRobinMassMatrix(rBC->m_robinID,
333  rBC->m_robinPrimitiveCoeffs, new_mat);
334  }
335 
336  // redeclare loc_mat to point to new_mat plus the scalar.
337  loc_mat = MemoryManager<DNekScalMat>::AllocateSharedPtr(1.0, new_mat);
338  }
339 
340  // finally return the matrix.
341  return loc_mat;
342 }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
static void Dscal(const int &n, const double &alpha, double *x, const int &incx)
BLAS level 1: x = alpha x.
Definition: Blas.hpp:168
std::shared_ptr< RobinBCInfo > RobinBCInfoSharedPtr
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
double NekDouble

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Blas::Dscal(), GetBlockMatrixKey(), m_expList, and m_robinBCInfo.

Referenced by GetBlock(), and Nektar::MultiRegions::GlobalLinSysStaticCond::SetupTopLevel().

◆ v_GetNumBlocks()

int Nektar::MultiRegions::GlobalLinSys::v_GetNumBlocks ( )
protectedvirtual

Get the number of blocks in this system.

At the top level this corresponds to the number of elements in the expansion list.

Reimplemented in Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 228 of file GlobalLinSys.cpp.

229 {
230  return m_expList.lock()->GetExpSize();
231 }

References m_expList.

Referenced by GetNumBlocks().

◆ v_GetStaticCondBlock()

DNekScalBlkMatSharedPtr Nektar::MultiRegions::GlobalLinSys::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 in Nektar::MultiRegions::GlobalLinSysPETScStaticCond, and Nektar::MultiRegions::GlobalLinSysIterativeStaticCond.

Definition at line 352 of file GlobalLinSys.cpp.

353 {
354 
355  LocalRegions::ExpansionSharedPtr vExp = m_expList.lock()->GetExp(n);
356  DNekScalBlkMatSharedPtr loc_mat;
357  loc_mat = vExp->GetLocStaticCondMatrix(GetBlockMatrixKey(n));
358 
359  if (m_robinBCInfo.count(n) != 0) // add robin mass matrix
360  {
361  DNekScalMatSharedPtr tmp_mat;
363 
364  tmp_mat = loc_mat->GetBlock(0, 0);
365 
366  // declare local matrix from scaled matrix.
367  int rows = tmp_mat->GetRows();
368  int cols = tmp_mat->GetColumns();
369  const NekDouble *dat = tmp_mat->GetRawPtr();
370  DNekMatSharedPtr new_mat =
372  Blas::Dscal(rows * cols, tmp_mat->Scale(), new_mat->GetRawPtr(), 1);
373 
374  // add local matrix contribution
375  for (rBC = m_robinBCInfo.find(n)->second; rBC; rBC = rBC->next)
376  {
377  vExp->AddRobinMassMatrix(rBC->m_robinID,
378  rBC->m_robinPrimitiveCoeffs, new_mat);
379  }
380 
381  // redeclare loc_mat to point to new_mat plus the scalar.
382  tmp_mat = MemoryManager<DNekScalMat>::AllocateSharedPtr(1.0, new_mat);
383  DNekScalBlkMatSharedPtr new_loc_mat;
384  unsigned int exp_size[] = {tmp_mat->GetRows(),
385  loc_mat->GetBlock(1, 1)->GetRows()};
386  unsigned int nblks = 2;
388  nblks, nblks, exp_size, exp_size);
389 
390  new_loc_mat->SetBlock(0, 0, tmp_mat);
391  new_loc_mat->SetBlock(0, 1, loc_mat->GetBlock(0, 1));
392  new_loc_mat->SetBlock(1, 0, loc_mat->GetBlock(1, 0));
393  new_loc_mat->SetBlock(1, 1, loc_mat->GetBlock(1, 1));
394  loc_mat = new_loc_mat;
395  }
396 
397  return loc_mat;
398 }
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:79

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Blas::Dscal(), GetBlockMatrixKey(), m_expList, and m_robinBCInfo.

Referenced by GetStaticCondBlock(), and Nektar::MultiRegions::GlobalLinSysStaticCond::SetupTopLevel().

◆ v_Initialise()

void Nektar::MultiRegions::GlobalLinSys::v_Initialise ( const std::shared_ptr< AssemblyMap > &  pLocToGloMap)
privatevirtual

Reimplemented in Nektar::MultiRegions::GlobalLinSysStaticCond.

Definition at line 417 of file GlobalLinSys.cpp.

419 {
420  boost::ignore_unused(pLocToGloMap);
421  NEKERROR(ErrorUtil::efatal, "Method does not exist");
422 }
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:209

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

Referenced by Initialise().

◆ v_InitObject()

void Nektar::MultiRegions::GlobalLinSys::v_InitObject ( )
privatevirtual

◆ v_Solve()

virtual void Nektar::MultiRegions::GlobalLinSys::v_Solve ( const Array< OneD, const NekDouble > &  in,
Array< OneD, NekDouble > &  out,
const AssemblyMapSharedPtr locToGloMap,
const Array< OneD, const NekDouble > &  dirForcing = NullNekDouble1DArray 
)
privatepure virtual

◆ v_SolveLinearSystem()

virtual void Nektar::MultiRegions::GlobalLinSys::v_SolveLinearSystem ( const int  pNumRows,
const Array< OneD, const NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput,
const AssemblyMapSharedPtr locToGloMap,
const int  pNumDir 
)
privatepure virtual

Member Data Documentation

◆ def

std::string Nektar::MultiRegions::GlobalLinSys::def
staticprivate
Initial value:
=
"DirectStaticCond")
static std::string RegisterDefaultSolverInfo(const std::string &pName, const std::string &pValue)
Registers the default string value of a solver info property.

Definition at line 157 of file GlobalLinSys.h.

◆ lookupIds

std::string Nektar::MultiRegions::GlobalLinSys::lookupIds
staticprivate
Initial value:
= {
"GlobalSysSoln", "DirectFull", MultiRegions::eDirectFullMatrix),
"GlobalSysSoln", "DirectStaticCond", MultiRegions::eDirectStaticCond),
"GlobalSysSoln", "DirectMultiLevelStaticCond",
"GlobalSysSoln", "IterativeFull", MultiRegions::eIterativeFull),
"GlobalSysSoln", "IterativeStaticCond",
"GlobalSysSoln", "IterativeMultiLevelStaticCond",
"GlobalSysSoln", "XxtFull", MultiRegions::eXxtFullMatrix),
"GlobalSysSoln", "XxtStaticCond", MultiRegions::eXxtStaticCond),
"GlobalSysSoln", "XxtMultiLevelStaticCond",
"GlobalSysSoln", "PETScFull", MultiRegions::ePETScFullMatrix),
"GlobalSysSoln", "PETScStaticCond", MultiRegions::ePETScStaticCond),
"GlobalSysSoln", "PETScMultiLevelStaticCond",
static std::string RegisterEnumValue(std::string pEnum, std::string pString, int pEnumValue)
Registers an enumeration value.

Definition at line 156 of file GlobalLinSys.h.

◆ m_expList

const std::weak_ptr<ExpList> Nektar::MultiRegions::GlobalLinSys::m_expList
protected

Local Matrix System.

Definition at line 123 of file GlobalLinSys.h.

Referenced by Nektar::MultiRegions::GlobalLinSysDirectFull::AssembleFullMatrix(), Nektar::MultiRegions::GlobalLinSysXxtFull::AssembleMatrixArrays(), Nektar::MultiRegions::GlobalLinSysPETSc::CalculateReordering(), Nektar::MultiRegions::GlobalLinSysStaticCond::ConstructNextLevelCondensedSystem(), Nektar::MultiRegions::GlobalLinSysIterative::DoProjection(), GetBlockMatrixKey(), GetLocMat(), Nektar::MultiRegions::GlobalLinSysIterative::GlobalLinSysIterative(), Nektar::MultiRegions::GlobalLinSysPETSc::GlobalLinSysPETSc(), Nektar::MultiRegions::GlobalLinSysPETScFull::GlobalLinSysPETScFull(), Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::PrepareLocalSchurComplement(), Nektar::MultiRegions::GlobalLinSysIterative::Set_Rhs_Magnitude(), Nektar::MultiRegions::GlobalLinSysStaticCond::SetupTopLevel(), Nektar::MultiRegions::GlobalLinSysIterative::UpdateKnownSolutions(), Nektar::MultiRegions::GlobalLinSysXxtStaticCond::v_AssembleSchurComplement(), Nektar::MultiRegions::GlobalLinSysPETScFull::v_DoMatrixMultiply(), Nektar::MultiRegions::GlobalLinSysIterativeFull::v_DoMatrixMultiply(), v_DropStaticCondBlock(), v_GetBlock(), v_GetNumBlocks(), v_GetStaticCondBlock(), Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_InitObject(), Nektar::MultiRegions::GlobalLinSysPETScStaticCond::v_InitObject(), Nektar::MultiRegions::GlobalLinSysIterativeStaticCond::v_PreSolve(), Nektar::MultiRegions::GlobalLinSysIterativeFull::v_Solve(), Nektar::MultiRegions::GlobalLinSysPETScFull::v_Solve(), Nektar::MultiRegions::GlobalLinSysXxtFull::v_Solve(), Nektar::MultiRegions::GlobalLinSysDirectFull::v_Solve(), and Nektar::MultiRegions::GlobalLinSysIterative::v_SolveLinearSystem().

◆ m_linSysKey

const GlobalLinSysKey Nektar::MultiRegions::GlobalLinSys::m_linSysKey
protected

◆ m_robinBCInfo

const std::map<int, RobinBCInfoSharedPtr> Nektar::MultiRegions::GlobalLinSys::m_robinBCInfo
protected

◆ m_verbose

bool Nektar::MultiRegions::GlobalLinSys::m_verbose
protected