58 "PETSc static condensation.");
63 "PETSc multi-level static condensation.");
86 const std::shared_ptr<AssemblyMap> &pLocToGloMap)
93 "This constructor is only valid when using static "
96 pLocToGloMap->GetGlobalSysSolnType(),
97 "The local to global map is not set up for the requested "
109 const std::shared_ptr<AssemblyMap> &pLocToGloMap,
139 int n, n_exp =
m_expList.lock()->GetNumElmts();
146 for (n = 0; n < n_exp; ++n)
152 m_precon->TransformedSchurCompl(n, cnt, mat);
154 cnt += mat->GetRows();
170 int i, j, n, cnt, gid1, gid2, loc_lda;
173 const int nDirDofs = pLocToGloMap->GetNumGlobalDirBndCoeffs();
191 pLocToGloMap->GetGlobalToUniversalBndMapUnique(),
195 SetUpMatVec(pLocToGloMap->GetNumGlobalBndCoeffs(), nDirDofs);
203 string variable = pLocToGloMap->GetVariable();
205 if (pSession->DefinesGlobalSysSolnInfo(variable,
206 "IterativeSolverTolerance"))
208 IterativeSolverTolerance = boost::lexical_cast<double>(
209 pSession->GetGlobalSysSolnInfo(variable,
"IterativeSolverTolerance")
212 else if (pSession->DefinesParameter(
"IterativeSolverTolerance"))
214 IterativeSolverTolerance =
215 pSession->GetParameter(
"IterativeSolverTolerance");
227 for (n = cnt = 0; n <
m_schurCompl->GetNumberOfBlockRows(); ++n)
230 loc_lda = loc_mat->GetRows();
232 for (i = 0; i < loc_lda; ++i)
234 gid1 = pLocToGloMap->GetLocalToGlobalBndMap(cnt + i) - nDirDofs;
235 sign1 = pLocToGloMap->GetLocalToGlobalBndSign(cnt + i);
239 for (j = 0; j < loc_lda; ++j)
241 gid2 = pLocToGloMap->GetLocalToGlobalBndMap(cnt + j) -
243 sign2 = pLocToGloMap->GetLocalToGlobalBndSign(cnt + j);
247 value = sign1 * sign2 * (*loc_mat)(i, j);
248 MatSetValue(
m_matrix, gid1ro, gid2ro, value,
258 MatAssemblyBegin(
m_matrix, MAT_FINAL_ASSEMBLY);
259 MatAssemblyEnd(
m_matrix, MAT_FINAL_ASSEMBLY);
267 unsigned int nbdry = localMat->GetRows();
268 unsigned int nblks = 1;
269 unsigned int esize[1] = {nbdry};
272 nblks, nblks, esize, esize);
273 schurComplBlock->SetBlock(0, 0, localMat);
275 return schurComplBlock;
308 const int nHomDofs = pNumRows - pNumDir;
317 locToGloMap->AssembleBnd(pInput, Glo);
324 VecAssemblyBegin(
m_b);
330 KSPConvergedReason reason;
331 KSPGetConvergedReason(
m_ksp, &reason);
332 ASSERTL0(reason > 0,
"PETSc solver diverged, reason is: " +
333 std::string(KSPConvergedReasons[reason]));
344 locToGloMap->GlobalToLocalBnd(Glo, pOutput);
351 m_precon->DoTransformBasisToLowEnergy(pInOut);
357 m_precon->DoTransformCoeffsFromLowEnergy(pInOut);
363 m_precon->DoTransformCoeffsToLowEnergy(pInput, pOutput);
382 int nLocBndDofs = asmMap->GetNumLocalBndCoeffs();
383 int nDirDofs = asmMap->GetNumGlobalDirBndCoeffs();
386 asmMap->GlobalToLocalBnd(input, in.GetPtr(), nDirDofs);
387 out = (*m_schurCompl) * in;
388 asmMap->AssembleBnd(out.
GetPtr(), output, nDirDofs);
396 const std::shared_ptr<AssemblyMap> &l2gMap)
400 mkey, pExpList, pSchurCompl, pBinvD, pC, pInvD, l2gMap,
m_precon);
401 sys->Initialise(l2gMap);
#define ASSERTL0(condition, msg)
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
const std::weak_ptr< ExpList > m_expList
Local Matrix System.
const GlobalLinSysKey m_linSysKey
Key associated with this linear system.
PreconditionerSharedPtr CreatePrecon(AssemblyMapSharedPtr asmMap)
Create a preconditioner object from the parameters defined in the supplied assembly map.
void Initialise(const std::shared_ptr< AssemblyMap > &pLocToGloMap)
Describe a linear system.
GlobalSysSolnType GetGlobalSysSolnType() const
Return the associated solution type.
A PETSc global linear system.
Vec m_x
PETSc vector objects used for local storage.
PETScMatMult m_matMult
Enumerator to select matrix multiplication type.
std::vector< int > m_reorderedMap
Reordering that takes universal IDs to a unique row in the PETSc matrix.
void SetUpScatter()
Set up PETSc local (equivalent to Nektar++ global) and global (equivalent to universal) scatter maps.
Mat m_matrix
PETSc matrix object.
PreconditionerSharedPtr m_precon
void SetUpSolver(NekDouble tolerance)
Set up KSP solver object.
void SetUpMatVec(int nGlobal, int nDir)
Construct PETSc matrix and vector handles.
VecScatter m_ctx
PETSc scatter context that takes us between Nektar++ global ordering and PETSc vector ordering.
void CalculateReordering(const Array< OneD, const int > &glo2uniMap, const Array< OneD, const int > &glo2unique, const AssemblyMapSharedPtr &pLocToGloMap)
Calculate a reordering of universal IDs for PETSc.
KSP m_ksp
KSP object that represents solver system.
void v_CoeffsFwdTransform(const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput) override
static std::string className2
void v_AssembleSchurComplement(std::shared_ptr< AssemblyMap > locToGloMap) override
Assemble the Schur complement matrix.
static GlobalLinSysSharedPtr create(const GlobalLinSysKey &pLinSysKey, const std::weak_ptr< ExpList > &pExpList, const std::shared_ptr< AssemblyMap > &pLocToGloMap)
Creates an instance of this class.
void v_SolveLinearSystem(const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir) override
Solve linear system using PETSc.
GlobalLinSysStaticCondSharedPtr v_Recurse(const GlobalLinSysKey &mkey, const std::weak_ptr< ExpList > &pExpList, const DNekScalBlkMatSharedPtr pSchurCompl, const DNekScalBlkMatSharedPtr pBinvD, const DNekScalBlkMatSharedPtr pC, const DNekScalBlkMatSharedPtr pInvD, const std::shared_ptr< AssemblyMap > &locToGloMap) override
~GlobalLinSysPETScStaticCond() override
void v_CoeffsBwdTransform(Array< OneD, NekDouble > &pInOut) override
void v_InitObject() override
void v_DoMatrixMultiply(const Array< OneD, const NekDouble > &input, Array< OneD, NekDouble > &output) override
Apply matrix-vector multiplication using local approach and the assembly map.
DNekScalBlkMatSharedPtr v_GetStaticCondBlock(unsigned int n) override
Retrieves a the static condensation block matrices from n-th expansion using the matrix key provided ...
static std::string className
Name of class.
void v_PreSolve(int scLevel, Array< OneD, NekDouble > &F_bBnd) override
GlobalLinSysPETScStaticCond(const GlobalLinSysKey &mkey, const std::weak_ptr< ExpList > &pExpList, const std::shared_ptr< AssemblyMap > &locToGloMap)
Constructor for full direct matrix solve.
void v_BasisFwdTransform(Array< OneD, NekDouble > &pInOut) override
DNekScalBlkMatSharedPtr m_schurCompl
Block Schur complement matrix.
std::weak_ptr< AssemblyMap > m_locToGloMap
Local to global map.
void SetupTopLevel(const std::shared_ptr< AssemblyMap > &locToGloMap)
Set up the storage for the Schur complement or the top level of the multi-level Schur complement.
DNekScalBlkMatSharedPtr m_BinvD
Block matrix.
DNekScalBlkMatSharedPtr m_C
Block matrix.
DNekScalBlkMatSharedPtr m_invD
Block matrix.
StdRegions::MatrixType GetMatrixType() const
Return the matrix type.
Array< OneD, DataType > & GetPtr()
std::shared_ptr< SessionReader > SessionReaderSharedPtr
@ ePETScMultiLevelStaticCond
std::shared_ptr< GlobalLinSysStaticCond > GlobalLinSysStaticCondSharedPtr
GlobalLinSysFactory & GetGlobalLinSysFactory()
std::shared_ptr< Preconditioner > PreconditionerSharedPtr
std::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
std::shared_ptr< GlobalLinSysPETScStaticCond > GlobalLinSysPETScStaticCondSharedPtr
static const NekDouble kNekIterativeTol
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
void Zero(int n, T *x, const int incx)
Zero vector.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)