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

#include <Preconditioner.h>

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

Public Member Functions

 Preconditioner (const std::shared_ptr< GlobalLinSys > &plinsys, const AssemblyMapSharedPtr &pLocToGloMap)
 
virtual ~Preconditioner ()
 
void DoPreconditioner (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
 
void DoPreconditionerWithNonVertOutput (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const Array< OneD, NekDouble > &pNonVertOutput, Array< OneD, NekDouble > &pVertForce=NullNekDouble1DArray)
 
void DoTransformBasisToLowEnergy (Array< OneD, NekDouble > &pInOut)
 
void DoTransformCoeffsFromLowEnergy (Array< OneD, NekDouble > &pInOut)
 
void DoTransformCoeffsToLowEnergy (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
 
void DoTransformBasisFromLowEnergy (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
 
void BuildPreconditioner ()
 
void InitObject ()
 
Array< OneD, NekDoubleAssembleStaticCondGlobalDiagonals ()
 Performs global assembly of diagonal entries to global Schur complement matrix. More...
 
const DNekScalBlkMatSharedPtrGetBlockTransformedSchurCompl () const
 
const DNekScalBlkMatSharedPtrGetBlockCMatrix () const
 
const DNekScalBlkMatSharedPtrGetBlockInvDMatrix () const
 
const DNekScalBlkMatSharedPtrGetBlockSchurCompl () const
 
const DNekScalBlkMatSharedPtrGetBlockTransformationMatrix () const
 
const DNekScalBlkMatSharedPtrGetBlockTransposedTransformationMatrix () const
 
DNekScalMatSharedPtr TransformedSchurCompl (int offset, int bndoffset, const std::shared_ptr< DNekScalMat > &loc_mat)
 

Protected Member Functions

virtual DNekScalMatSharedPtr v_TransformedSchurCompl (int offset, int bndoffset, const std::shared_ptr< DNekScalMat > &loc_mat)
 Get block elemental transposed transformation matrix \(\mathbf{R}^{T}\). More...
 

Protected Attributes

const std::weak_ptr< GlobalLinSysm_linsys
 
PreconditionerType m_preconType
 
DNekMatSharedPtr m_preconditioner
 
std::weak_ptr< AssemblyMapm_locToGloMap
 
LibUtilities::CommSharedPtr m_comm
 

Private Member Functions

void NullPreconditioner (void)
 
virtual void v_InitObject ()
 
virtual void v_DoPreconditioner (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
 Apply a preconditioner to the conjugate gradient method. More...
 
virtual void v_DoPreconditionerWithNonVertOutput (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const Array< OneD, NekDouble > &pNonVertOutput, Array< OneD, NekDouble > &pVertForce)
 Apply a preconditioner to the conjugate gradient method with an output for non-vertex degrees of freedom. More...
 
virtual void v_DoTransformBasisToLowEnergy (Array< OneD, NekDouble > &pInOut)
 Transform from original basis to low energy basis. More...
 
virtual void v_DoTransformCoeffsFromLowEnergy (Array< OneD, NekDouble > &pInOut)
 Transform from low energy coeffs to orignal basis. More...
 
virtual void v_DoTransformCoeffsToLowEnergy (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
 Multiply by the block inverse transformation matrix. More...
 
virtual void v_DoTransformBasisFromLowEnergy (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
 Multiply by the block transposed inverse transformation matrix. More...
 
virtual void v_BuildPreconditioner ()
 

Static Private Attributes

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

Detailed Description

Definition at line 64 of file Preconditioner.h.

Constructor & Destructor Documentation

◆ Preconditioner()

Nektar::MultiRegions::Preconditioner::Preconditioner ( const std::shared_ptr< GlobalLinSys > &  plinsys,
const AssemblyMapSharedPtr pLocToGloMap 
)
@class Preconditioner

This class implements preconditioning for the conjugate

gradient matrix solver.

Definition at line 77 of file Preconditioner.cpp.

80  : m_linsys(plinsys),
81  m_preconType(pLocToGloMap->GetPreconType()),
82  m_locToGloMap(pLocToGloMap)
83  {
84  }
const std::weak_ptr< GlobalLinSys > m_linsys
std::weak_ptr< AssemblyMap > m_locToGloMap

◆ ~Preconditioner()

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

Definition at line 72 of file Preconditioner.h.

72 {}

Member Function Documentation

◆ AssembleStaticCondGlobalDiagonals()

Array< OneD, NekDouble > Nektar::MultiRegions::Preconditioner::AssembleStaticCondGlobalDiagonals ( )

Performs global assembly of diagonal entries to global Schur complement matrix.

Definition at line 185 of file Preconditioner.cpp.

186  {
187  auto asmMap = m_locToGloMap.lock();
188 
189  int nGlobalBnd = asmMap->GetNumGlobalBndCoeffs();
190  int nDirBnd = asmMap->GetNumGlobalDirBndCoeffs();
191  int rows = nGlobalBnd - nDirBnd;
192 
193  DNekScalBlkMatSharedPtr loc_mat;
194  DNekScalMatSharedPtr bnd_mat;
195  int sign1, sign2, gid1, gid2, i, j, n, cnt;
196  Array<OneD, NekDouble> diagonals(rows,0.0);
197 
198  // Extract diagonal contributions of globally assembled
199  // schur complement matrix
200  for (cnt = n = 0; n < m_linsys.lock()->GetNumBlocks(); ++n)
201  {
202  // Get statically condensed local matrix.
203  loc_mat = (m_linsys.lock())->GetStaticCondBlock(n);
204 
205  // Extract boundary block.
206  bnd_mat = loc_mat->GetBlock(0, 0);
207 
208  // Offset by number of rows.
209  int bnd_row = bnd_mat->GetRows();
210 
211  for (i = 0; i < bnd_row; ++i)
212  {
213  gid1 = asmMap->GetLocalToGlobalBndMap (cnt + i)
214  - nDirBnd;
215  sign1 = asmMap->GetLocalToGlobalBndSign(cnt + i);
216 
217  if (gid1 < 0)
218  {
219  continue;
220  }
221 
222  for (j = 0; j < bnd_row; ++j)
223  {
224  gid2 = asmMap->GetLocalToGlobalBndMap (cnt + j)
225  - nDirBnd;
226  sign2 = asmMap->GetLocalToGlobalBndSign(cnt + j);
227 
228  if (gid2 == gid1)
229  {
230  diagonals[gid1] += sign1 * sign2 * (*bnd_mat)(i, j);
231  }
232  }
233  }
234  cnt += bnd_row;
235  }
236 
237  return diagonals;
238  }
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:73

References m_linsys, and m_locToGloMap.

Referenced by Nektar::MultiRegions::PreconditionerDiagonal::StaticCondDiagonalPreconditionerSum().

◆ BuildPreconditioner()

void Nektar::MultiRegions::Preconditioner::BuildPreconditioner ( )
inline

Definition at line 260 of file Preconditioner.h.

261  {
263  }

References v_BuildPreconditioner().

◆ DoPreconditioner()

void Nektar::MultiRegions::Preconditioner::DoPreconditioner ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput 
)
inline

Definition at line 198 of file Preconditioner.h.

201  {
202  v_DoPreconditioner(pInput,pOutput);
203  }
virtual void v_DoPreconditioner(const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
Apply a preconditioner to the conjugate gradient method.

References v_DoPreconditioner().

◆ DoPreconditionerWithNonVertOutput()

void Nektar::MultiRegions::Preconditioner::DoPreconditionerWithNonVertOutput ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput,
const Array< OneD, NekDouble > &  pNonVertOutput,
Array< OneD, NekDouble > &  pVertForce = NullNekDouble1DArray 
)
inline

Definition at line 209 of file Preconditioner.h.

214  {
215  v_DoPreconditionerWithNonVertOutput(pInput,pOutput,pNonVertOutput,
216  pVertForce);
217  }
virtual void v_DoPreconditionerWithNonVertOutput(const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const Array< OneD, NekDouble > &pNonVertOutput, Array< OneD, NekDouble > &pVertForce)
Apply a preconditioner to the conjugate gradient method with an output for non-vertex degrees of free...

References v_DoPreconditionerWithNonVertOutput().

◆ DoTransformBasisFromLowEnergy()

void Nektar::MultiRegions::Preconditioner::DoTransformBasisFromLowEnergy ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput 
)
inline

Definition at line 250 of file Preconditioner.h.

253  {
254  v_DoTransformBasisFromLowEnergy(pInput,pOutput);
255  }
virtual void v_DoTransformBasisFromLowEnergy(const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
Multiply by the block transposed inverse transformation matrix.

References v_DoTransformBasisFromLowEnergy().

◆ DoTransformBasisToLowEnergy()

void Nektar::MultiRegions::Preconditioner::DoTransformBasisToLowEnergy ( Array< OneD, NekDouble > &  pInOut)
inline

Definition at line 222 of file Preconditioner.h.

224  {
226  }
virtual void v_DoTransformBasisToLowEnergy(Array< OneD, NekDouble > &pInOut)
Transform from original basis to low energy basis.

References v_DoTransformBasisToLowEnergy().

◆ DoTransformCoeffsFromLowEnergy()

void Nektar::MultiRegions::Preconditioner::DoTransformCoeffsFromLowEnergy ( Array< OneD, NekDouble > &  pInOut)
inline

Definition at line 231 of file Preconditioner.h.

233  {
235  }
virtual void v_DoTransformCoeffsFromLowEnergy(Array< OneD, NekDouble > &pInOut)
Transform from low energy coeffs to orignal basis.

References v_DoTransformCoeffsFromLowEnergy().

◆ DoTransformCoeffsToLowEnergy()

void Nektar::MultiRegions::Preconditioner::DoTransformCoeffsToLowEnergy ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput 
)
inline

Definition at line 240 of file Preconditioner.h.

243  {
244  v_DoTransformCoeffsToLowEnergy(pInput,pOutput);
245  }
virtual void v_DoTransformCoeffsToLowEnergy(const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
Multiply by the block inverse transformation matrix.

References v_DoTransformCoeffsToLowEnergy().

◆ GetBlockCMatrix()

const DNekScalBlkMatSharedPtr& Nektar::MultiRegions::Preconditioner::GetBlockCMatrix ( ) const
inline

◆ GetBlockInvDMatrix()

const DNekScalBlkMatSharedPtr& Nektar::MultiRegions::Preconditioner::GetBlockInvDMatrix ( ) const
inline

◆ GetBlockSchurCompl()

const DNekScalBlkMatSharedPtr& Nektar::MultiRegions::Preconditioner::GetBlockSchurCompl ( ) const
inline

◆ GetBlockTransformationMatrix()

const DNekScalBlkMatSharedPtr& Nektar::MultiRegions::Preconditioner::GetBlockTransformationMatrix ( ) const
inline

◆ GetBlockTransformedSchurCompl()

const DNekScalBlkMatSharedPtr& Nektar::MultiRegions::Preconditioner::GetBlockTransformedSchurCompl ( ) const
inline

◆ GetBlockTransposedTransformationMatrix()

const DNekScalBlkMatSharedPtr& Nektar::MultiRegions::Preconditioner::GetBlockTransposedTransformationMatrix ( ) const
inline

◆ InitObject()

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

Definition at line 180 of file Preconditioner.h.

181  {
182  v_InitObject();
183  }

References v_InitObject().

◆ NullPreconditioner()

void Nektar::MultiRegions::Preconditioner::NullPreconditioner ( void  )
private

◆ TransformedSchurCompl()

DNekScalMatSharedPtr Nektar::MultiRegions::Preconditioner::TransformedSchurCompl ( int  offset,
int  bndoffset,
const std::shared_ptr< DNekScalMat > &  loc_mat 
)
inline

Definition at line 188 of file Preconditioner.h.

191  {
192  return v_TransformedSchurCompl(offset,bndoffset,loc_mat);
193  }
virtual DNekScalMatSharedPtr v_TransformedSchurCompl(int offset, int bndoffset, const std::shared_ptr< DNekScalMat > &loc_mat)
Get block elemental transposed transformation matrix .

References v_TransformedSchurCompl().

◆ v_BuildPreconditioner()

void Nektar::MultiRegions::Preconditioner::v_BuildPreconditioner ( )
privatevirtual

◆ v_DoPreconditioner()

void Nektar::MultiRegions::Preconditioner::v_DoPreconditioner ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput 
)
privatevirtual

Apply a preconditioner to the conjugate gradient method.

Reimplemented in Nektar::MultiRegions::PreconditionerLowEnergy, Nektar::MultiRegions::PreconditionerLinearWithLowEnergy, Nektar::MultiRegions::PreconditionerLinearWithDiag, Nektar::MultiRegions::PreconditionerLinearWithBlock, Nektar::MultiRegions::PreconditionerLinear, Nektar::MultiRegions::PreconditionerNull, Nektar::MultiRegions::PreconditionerDiagonal, and Nektar::MultiRegions::PreconditionerBlock.

Definition at line 103 of file Preconditioner.cpp.

106  {
107  boost::ignore_unused(pInput, pOutput);
108  NEKERROR(ErrorUtil::efatal,"Method does not exist" );
109  }
#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 DoPreconditioner().

◆ v_DoPreconditionerWithNonVertOutput()

void Nektar::MultiRegions::Preconditioner::v_DoPreconditionerWithNonVertOutput ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput,
const Array< OneD, NekDouble > &  pNonVertOutput,
Array< OneD, NekDouble > &  pVertForce 
)
privatevirtual

Apply a preconditioner to the conjugate gradient method with an output for non-vertex degrees of freedom.

Reimplemented in Nektar::MultiRegions::PreconditionerLinear.

Definition at line 115 of file Preconditioner.cpp.

120  {
121  boost::ignore_unused(pInput, pOutput, pNonVertOutput, pVertForce);
122  NEKERROR(ErrorUtil::efatal, "Method does not exist");
123  }

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

Referenced by DoPreconditionerWithNonVertOutput().

◆ v_DoTransformBasisFromLowEnergy()

void Nektar::MultiRegions::Preconditioner::v_DoTransformBasisFromLowEnergy ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput 
)
privatevirtual

Multiply by the block transposed inverse transformation matrix.

Reimplemented in Nektar::MultiRegions::PreconditionerLowEnergy, and Nektar::MultiRegions::PreconditionerLinearWithLowEnergy.

Definition at line 156 of file Preconditioner.cpp.

159  {
160  boost::ignore_unused(pInput, pOutput);
161  NEKERROR(ErrorUtil::efatal,"Method does not exist" );
162  }

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

Referenced by DoTransformBasisFromLowEnergy().

◆ v_DoTransformBasisToLowEnergy()

void Nektar::MultiRegions::Preconditioner::v_DoTransformBasisToLowEnergy ( Array< OneD, NekDouble > &  pInOut)
privatevirtual

Transform from original basis to low energy basis.

Reimplemented in Nektar::MultiRegions::PreconditionerLowEnergy, and Nektar::MultiRegions::PreconditionerLinearWithLowEnergy.

Definition at line 128 of file Preconditioner.cpp.

130  {
131  boost::ignore_unused(pInOut);
132  }

Referenced by DoTransformBasisToLowEnergy().

◆ v_DoTransformCoeffsFromLowEnergy()

void Nektar::MultiRegions::Preconditioner::v_DoTransformCoeffsFromLowEnergy ( Array< OneD, NekDouble > &  pInOut)
privatevirtual

Transform from low energy coeffs to orignal basis.

Reimplemented in Nektar::MultiRegions::PreconditionerLowEnergy, and Nektar::MultiRegions::PreconditionerLinearWithLowEnergy.

Definition at line 137 of file Preconditioner.cpp.

139  {
140  boost::ignore_unused(pInOut);
141  }

Referenced by DoTransformCoeffsFromLowEnergy().

◆ v_DoTransformCoeffsToLowEnergy()

void Nektar::MultiRegions::Preconditioner::v_DoTransformCoeffsToLowEnergy ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput 
)
privatevirtual

Multiply by the block inverse transformation matrix.

Reimplemented in Nektar::MultiRegions::PreconditionerLowEnergy, and Nektar::MultiRegions::PreconditionerLinearWithLowEnergy.

Definition at line 146 of file Preconditioner.cpp.

149  {
150  boost::ignore_unused(pInput, pOutput);
151  }

Referenced by DoTransformCoeffsToLowEnergy().

◆ v_InitObject()

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

◆ v_TransformedSchurCompl()

DNekScalMatSharedPtr Nektar::MultiRegions::Preconditioner::v_TransformedSchurCompl ( int  offset,
int  bndoffset,
const std::shared_ptr< DNekScalMat > &  loc_mat 
)
protectedvirtual

Get block elemental transposed transformation matrix \(\mathbf{R}^{T}\).

Reimplemented in Nektar::MultiRegions::PreconditionerLowEnergy, and Nektar::MultiRegions::PreconditionerLinearWithLowEnergy.

Definition at line 172 of file Preconditioner.cpp.

175  {
176  boost::ignore_unused(offset, bnd_offset);
177  return loc_mat;
178  }

Referenced by TransformedSchurCompl().

Member Data Documentation

◆ def

std::string Nektar::MultiRegions::Preconditioner::def
staticprivate
Initial value:
=
"Preconditioner", "Diagonal")
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 173 of file Preconditioner.h.

◆ lookupIds

std::string Nektar::MultiRegions::Preconditioner::lookupIds
staticprivate
Initial value:
= {
"Preconditioner", "Null", eNull),
"Preconditioner", "Diagonal", eDiagonal),
"Preconditioner", "FullLinearSpaceWithDiagonal",
"Preconditioner", "FullLinearSpace",eLinear),
"Preconditioner", "LowEnergyBlock",eLowEnergy),
"Preconditioner", "FullLinearSpaceWithLowEnergyBlock",
"Preconditioner", "Block",eBlock),
"Preconditioner", "FullLinearSpaceWithBlock",eLinearWithBlock),
}
static std::string RegisterEnumValue(std::string pEnum, std::string pString, int pEnumValue)
Registers an enumeration value.
@ eNull
No Solution type specified.

Definition at line 172 of file Preconditioner.h.

◆ m_comm

LibUtilities::CommSharedPtr Nektar::MultiRegions::Preconditioner::m_comm
protected

◆ m_linsys

const std::weak_ptr<GlobalLinSys> Nektar::MultiRegions::Preconditioner::m_linsys
protected

◆ m_locToGloMap

std::weak_ptr<AssemblyMap> Nektar::MultiRegions::Preconditioner::m_locToGloMap
protected

Definition at line 131 of file Preconditioner.h.

Referenced by AssembleStaticCondGlobalDiagonals(), Nektar::MultiRegions::PreconditionerBlock::BlockPreconditionerCG(), Nektar::MultiRegions::PreconditionerBlock::BlockPreconditionerHDG(), Nektar::MultiRegions::PreconditionerLowEnergy::CreateVariablePMask(), Nektar::MultiRegions::PreconditionerDiagonal::DiagonalPreconditionerSum(), Nektar::MultiRegions::PreconditionerLowEnergy::SetupBlockTransformationMatrix(), Nektar::MultiRegions::PreconditionerDiagonal::StaticCondDiagonalPreconditionerSum(), Nektar::MultiRegions::PreconditionerDiagonal::v_BuildPreconditioner(), Nektar::MultiRegions::PreconditionerLinear::v_BuildPreconditioner(), Nektar::MultiRegions::PreconditionerLowEnergy::v_BuildPreconditioner(), Nektar::MultiRegions::PreconditionerBlock::v_DoPreconditioner(), Nektar::MultiRegions::PreconditionerDiagonal::v_DoPreconditioner(), Nektar::MultiRegions::PreconditionerLinearWithLowEnergy::v_DoPreconditioner(), Nektar::MultiRegions::PreconditionerLowEnergy::v_DoPreconditioner(), Nektar::MultiRegions::PreconditionerLinear::v_DoPreconditionerWithNonVertOutput(), Nektar::MultiRegions::PreconditionerLowEnergy::v_DoTransformBasisFromLowEnergy(), Nektar::MultiRegions::PreconditionerLowEnergy::v_DoTransformBasisToLowEnergy(), Nektar::MultiRegions::PreconditionerLowEnergy::v_DoTransformCoeffsFromLowEnergy(), Nektar::MultiRegions::PreconditionerLowEnergy::v_DoTransformCoeffsToLowEnergy(), Nektar::MultiRegions::PreconditionerBlock::v_InitObject(), Nektar::MultiRegions::PreconditionerLinearWithBlock::v_InitObject(), Nektar::MultiRegions::PreconditionerLinearWithDiag::v_InitObject(), Nektar::MultiRegions::PreconditionerLinearWithLowEnergy::v_InitObject(), and Nektar::MultiRegions::PreconditionerLowEnergy::v_InitObject().

◆ m_preconditioner

DNekMatSharedPtr Nektar::MultiRegions::Preconditioner::m_preconditioner
protected

Definition at line 130 of file Preconditioner.h.

◆ m_preconType

PreconditionerType Nektar::MultiRegions::Preconditioner::m_preconType
protected

Definition at line 129 of file Preconditioner.h.