44 namespace MultiRegions
57 string GlobalLinSysIterativeFull::className
60 GlobalLinSysIterativeFull::create,
61 "Iterative solver for full matrix system.");
71 GlobalLinSysIterativeFull::GlobalLinSysIterativeFull(
73 const boost::weak_ptr<ExpList> &pExp,
74 const boost::shared_ptr<AssemblyMap> &pLocToGloMap)
79 "This routine should only be used when using an Iterative "
80 "conjugate gradient matrix solve.");
117 boost::shared_ptr<MultiRegions::ExpList> expList =
m_expList.lock();
119 if ((
m_locToGloMap = boost::dynamic_pointer_cast<AssemblyMapCG>(
131 ASSERTL0(
false,
"Unknown map type");
134 bool dirForcCalculated = (bool) pDirForcing.num_elements();
135 int nDirDofs = pLocToGloMap->GetNumGlobalDirBndCoeffs();
136 int nGlobDofs = pLocToGloMap->GetNumGlobalCoeffs();
137 int nDirTotal = nDirDofs;
139 expList->GetComm()->GetRowComm()
147 if(dirForcCalculated)
150 pDirForcing.get(), 1,
157 expList->GeneralMatrixOp(
170 nGlobDofs, tmp, out, pLocToGloMap, nDirDofs);
171 Vmath::Vadd(nGlobDofs-nDirDofs, &out [nDirDofs], 1,
172 &pOutput[nDirDofs], 1, &pOutput[nDirDofs], 1);
176 ASSERTL0(
false,
"Need DG solve if using Dir BCs");
194 boost::shared_ptr<MultiRegions::ExpList> expList =
m_expList.lock();
203 "Robin boundaries not set up in IterativeFull solver.");
207 int nNonDir = nGlobal - nDir;
220 for (
int n = 0; n < expList->GetNumElmts(); ++n)
223 int offset = expList->GetCoeff_Offset(n);
224 int ncoeffs = expList->GetExp(nel)->GetNcoeffs();
233 for(rBC =
m_robinBCInfo.find(nel)->second;rBC; rBC = rBC->next)
235 vExp->AddRobinEdgeContribution(rBC->m_robinID,rBC->m_robinPrimitiveCoeffs, tmp = robin_l + offset);
247 Vmath::Vadd(nGlobal, pOutput, 1, robin_A, 1, pOutput, 1);
GlobalSysSolnType GetGlobalSysSolnType() const
Return the associated solution type.
#define ASSERTL0(condition, msg)
boost::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
virtual void v_UniqueMap()
boost::shared_ptr< RobinBCInfo > RobinBCInfoSharedPtr
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.
virtual ~GlobalLinSysIterativeFull()
Array< OneD, int > m_map
Global to universal unique map.
Describe a linear system.
const GlobalLinSysKey m_linSysKey
Key associated with this linear system.
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.
boost::shared_ptr< AssemblyMap > m_locToGloMap
virtual void v_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...
void Zero(int n, T *x, const int incx)
Zero vector.
GlobalLinSysFactory & GetGlobalLinSysFactory()
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
const std::map< int, RobinBCInfoSharedPtr > m_robinBCInfo
Robin boundary info.
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.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
virtual void v_DoMatrixMultiply(const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)
const boost::weak_ptr< ExpList > m_expList
Local Matrix System.