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, const bool &IsLocal=false)
 
void DoAssembleLoc (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const bool &ZeroDir)
 Apply an assembly and scatter back to lcoal array. More...
 
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...
 
virtual void v_InitObject ()
 
virtual void v_DoPreconditioner (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const bool &isLocal=false)
 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 ()
 

Protected Attributes

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

Private Member Functions

void NullPreconditioner (void)
 

Static Private Attributes

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

Detailed Description

This class implements preconditioning for the conjugate gradient matrix solver.

This class implements diagonal preconditioning for the conjugate gradient matrix solver.

Definition at line 63 of file Preconditioner.h.

Constructor & Destructor Documentation

◆ Preconditioner()

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

Definition at line 54 of file Preconditioner.cpp.

56 : m_linsys(plinsys), m_preconType(pLocToGloMap->GetPreconType()),
57 m_locToGloMap(pLocToGloMap)
58{
59}
const std::weak_ptr< GlobalLinSys > m_linsys
std::weak_ptr< AssemblyMap > m_locToGloMap

◆ ~Preconditioner()

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

Definition at line 71 of file Preconditioner.h.

72 {
73 }

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
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) - nDirBnd;
214 sign1 = asmMap->GetLocalToGlobalBndSign(cnt + i);
215
216 if (gid1 < 0)
217 {
218 continue;
219 }
220
221 for (j = 0; j < bnd_row; ++j)
222 {
223 gid2 = asmMap->GetLocalToGlobalBndMap(cnt + j) - nDirBnd;
224 sign2 = asmMap->GetLocalToGlobalBndSign(cnt + j);
225
226 if (gid2 == gid1)
227 {
228 diagonals[gid1] += sign1 * sign2 * (*bnd_mat)(i, j);
229 }
230 }
231 }
232 cnt += bnd_row;
233 }
234
235 return diagonals;
236}
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:79

References m_linsys, and m_locToGloMap.

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

◆ BuildPreconditioner()

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

Definition at line 239 of file Preconditioner.h.

240{
242}

References v_BuildPreconditioner().

◆ DoAssembleLoc()

void Nektar::MultiRegions::Preconditioner::DoAssembleLoc ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput,
const bool &  ZeroDir 
)

Apply an assembly and scatter back to lcoal array.

Definition at line 89 of file Preconditioner.cpp.

92{
93 auto assMap = m_locToGloMap.lock();
94 GlobalSysSolnType solvertype = assMap->GetGlobalSysSolnType();
95
96 if (solvertype == eIterativeFull)
97 {
98 assMap->Assemble(pInput, pOutput);
99 if (ZeroDir)
100 {
101 int nDir = assMap->GetNumGlobalDirBndCoeffs();
102 Vmath::Zero(nDir, pOutput, 1);
103 }
104 assMap->GlobalToLocal(pOutput, pOutput);
105 }
106 else // bnd version.
107 {
108 assMap->AssembleBnd(pInput, pOutput);
109 if (ZeroDir)
110 {
111 int nDir = assMap->GetNumGlobalDirBndCoeffs();
112 Vmath::Zero(nDir, pOutput, 1);
113 }
114 assMap->GlobalToLocalBnd(pOutput, pOutput);
115 }
116}
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273

References Nektar::MultiRegions::eIterativeFull, m_locToGloMap, and Vmath::Zero().

◆ DoPreconditioner()

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

Definition at line 181 of file Preconditioner.h.

184{
185 v_DoPreconditioner(pInput, pOutput, ZeroDir);
186}
virtual void v_DoPreconditioner(const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const bool &isLocal=false)
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 191 of file Preconditioner.h.

195{
196 v_DoPreconditionerWithNonVertOutput(pInput, pOutput, pNonVertOutput,
197 pVertForce);
198}
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 230 of file Preconditioner.h.

232{
233 v_DoTransformBasisFromLowEnergy(pInput, pOutput);
234}
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 203 of file Preconditioner.h.

205{
207}
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 212 of file Preconditioner.h.

214{
216}
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 221 of file Preconditioner.h.

223{
224 v_DoTransformCoeffsToLowEnergy(pInput, pOutput);
225}
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 164 of file Preconditioner.h.

165{
166 v_InitObject();
167}

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 172 of file Preconditioner.h.

174{
175 return v_TransformedSchurCompl(offset, bndoffset, loc_mat);
176}
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 ( )
protectedvirtual

◆ v_DoPreconditioner()

void Nektar::MultiRegions::Preconditioner::v_DoPreconditioner ( const Array< OneD, NekDouble > &  pInput,
Array< OneD, NekDouble > &  pOutput,
const bool &  isLocal = false 
)
protectedvirtual

◆ 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 
)
protectedvirtual

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 122 of file Preconditioner.cpp.

127{
128 NEKERROR(ErrorUtil::efatal, "Method does not exist");
129}

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 
)
protectedvirtual

Multiply by the block transposed inverse transformation matrix.

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

Definition at line 159 of file Preconditioner.cpp.

162{
163 NEKERROR(ErrorUtil::efatal, "Method does not exist");
164}

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

Referenced by DoTransformBasisFromLowEnergy().

◆ v_DoTransformBasisToLowEnergy()

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

Transform from original basis to low energy basis.

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

Definition at line 134 of file Preconditioner.cpp.

136{
137}

Referenced by DoTransformBasisToLowEnergy().

◆ v_DoTransformCoeffsFromLowEnergy()

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

Transform from low energy coeffs to orignal basis.

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

Definition at line 142 of file Preconditioner.cpp.

144{
145}

Referenced by DoTransformCoeffsFromLowEnergy().

◆ v_DoTransformCoeffsToLowEnergy()

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

Multiply by the block inverse transformation matrix.

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

Definition at line 150 of file Preconditioner.cpp.

153{
154}

Referenced by DoTransformCoeffsToLowEnergy().

◆ v_InitObject()

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

◆ 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::PreconditionerLinearWithLowEnergy, and Nektar::MultiRegions::PreconditionerLowEnergy.

Definition at line 174 of file Preconditioner.cpp.

177{
178 return loc_mat;
179}

Referenced by TransformedSchurCompl().

Member Data Documentation

◆ def

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

◆ lookupIds

std::string Nektar::MultiRegions::Preconditioner::lookupIds[]
staticprivate

Definition at line 156 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 123 of file Preconditioner.h.

Referenced by AssembleStaticCondGlobalDiagonals(), Nektar::MultiRegions::PreconditionerBlock::BlockPreconditionerCG(), Nektar::MultiRegions::PreconditionerBlock::BlockPreconditionerHDG(), Nektar::MultiRegions::PreconditionerLowEnergy::CreateVariablePMask(), Nektar::MultiRegions::PreconditionerDiagonal::DiagonalPreconditionerSum(), DoAssembleLoc(), Nektar::MultiRegions::PreconditionerLowEnergy::SetupBlockTransformationMatrix(), Nektar::MultiRegions::PreconditionerDiagonal::StaticCondDiagonalPreconditionerSum(), Nektar::MultiRegions::PreconditionerDiagonal::v_BuildPreconditioner(), Nektar::MultiRegions::PreconditionerJacobi::v_BuildPreconditioner(), Nektar::MultiRegions::PreconditionerLinear::v_BuildPreconditioner(), Nektar::MultiRegions::PreconditionerLowEnergy::v_BuildPreconditioner(), Nektar::MultiRegions::PreconditionerBlock::v_DoPreconditioner(), Nektar::MultiRegions::PreconditionerDiagonal::v_DoPreconditioner(), Nektar::MultiRegions::PreconditionerNull::v_DoPreconditioner(), Nektar::MultiRegions::PreconditionerJacobi::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 122 of file Preconditioner.h.

◆ m_preconType

std::string Nektar::MultiRegions::Preconditioner::m_preconType
protected

Definition at line 121 of file Preconditioner.h.