59 const std::shared_ptr<AssemblyMap> &pLocToGloMap)
63 const int nDirDofs = pLocToGloMap->GetNumGlobalDirBndCoeffs();
65 int i, j, n, cnt, gid1, gid2, loc_lda;
71 pLocToGloMap->GetGlobalToUniversalMapUnique(),
75 SetUpMatVec(pLocToGloMap->GetNumGlobalCoeffs(), nDirDofs);
84 for (n = cnt = 0; n <
m_expList.lock()->GetNumElmts(); ++n)
87 loc_lda = loc_mat->GetRows();
89 for (i = 0; i < loc_lda; ++i)
91 gid1 = pLocToGloMap->GetLocalToGlobalMap(cnt + i) - nDirDofs;
92 sign1 = pLocToGloMap->GetLocalToGlobalSign(cnt + i);
96 for (j = 0; j < loc_lda; ++j)
99 pLocToGloMap->GetLocalToGlobalMap(cnt + j) - nDirDofs;
100 sign2 = pLocToGloMap->GetLocalToGlobalSign(cnt + j);
104 value = sign1 * sign2 * (*loc_mat)(i, j);
105 MatSetValue(
m_matrix, gid1ro, gid2ro, value,
115 MatAssemblyBegin(
m_matrix, MAT_FINAL_ASSEMBLY);
116 MatAssemblyEnd(
m_matrix, MAT_FINAL_ASSEMBLY);
130 bool dirForcCalculated = (bool)pDirForcing.size();
131 int nDirDofs = pLocToGloMap->GetNumGlobalDirBndCoeffs();
132 int nGlobDofs = pLocToGloMap->GetNumGlobalCoeffs();
133 int nLocDofs = pLocToGloMap->GetNumLocalCoeffs();
135 int nDirTotal = nDirDofs;
136 std::shared_ptr<MultiRegions::ExpList> expList =
m_expList.lock();
137 expList->GetComm()->GetRowComm()->AllReduce(nDirTotal,
145 if (dirForcCalculated)
149 pDirForcing.size() >= nLocDofs,
150 "DirForcing is not of sufficient size. Is it in local space?");
151 Vmath::Vsub(nLocDofs, pLocInput, 1, pDirForcing, 1, rhs, 1);
157 expList->GeneralMatrixOp(
m_linSysKey, pLocOutput, rhs);
167 int offset = expList->GetCoeff_Offset(n);
171 for (rBC = r.second; rBC; rBC = rBC->next)
173 vExp->AddRobinTraceContribution(
174 rBC->m_robinID, rBC->m_robinPrimitiveCoeffs,
175 pLocOutput + offset, rhsloc = rhs + offset);
178 Vmath::Vsub(nLocDofs, pLocInput, 1, rhs, 1, rhs, 1);
187 Vmath::Vadd(nLocDofs, diff, 1, pLocOutput, 1, pLocOutput, 1);
206 std::shared_ptr<MultiRegions::ExpList> expList =
m_expList.lock();
#define ASSERTL0(condition, msg)
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
const std::weak_ptr< ExpList > m_expList
Local Matrix System.
const std::map< int, RobinBCInfoSharedPtr > m_robinBCInfo
Robin boundary info.
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.
const GlobalLinSysKey m_linSysKey
Key associated with this linear system.
DNekScalMatSharedPtr GetBlock(unsigned int n)
Describe a linear system.
GlobalLinSysPETScFull(const GlobalLinSysKey &pLinSysKey, const std::weak_ptr< ExpList > &pExpList, const std::shared_ptr< AssemblyMap > &pLocToGloMap)
Constructor for full direct matrix solve.
static std::string className
Name of class.
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_Solve(const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out, const AssemblyMapSharedPtr &locToGloMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray) override
Solve the linear system for given input and output vectors using a specified local to global map.
std::shared_ptr< AssemblyMap > m_locToGloMap
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.
A PETSc global linear system.
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.
void SetUpSolver(NekDouble tolerance)
Set up KSP solver object.
void SetUpMatVec(int nGlobal, int nDir)
Construct PETSc matrix and vector handles.
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.
std::shared_ptr< Expansion > ExpansionSharedPtr
std::shared_ptr< RobinBCInfo > RobinBCInfoSharedPtr
GlobalLinSysFactory & GetGlobalLinSysFactory()
std::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.