38 lduMatrix::solver::addasymMatrixConstructorToTable<PBiCG>
47 const word& fieldName,
86 const label nCells =
psi.
size();
88 solveScalar* __restrict__ psiPtr =
psi.
begin();
91 solveScalar* __restrict__ pAPtr = pA.begin();
94 solveScalar* __restrict__ wAPtr = wA.begin();
97 matrix_.Amul(wA,
psi, interfaceBouCoeffs_, interfaces_, cmpt);
102 solveScalar* __restrict__ rAPtr = rA.begin();
112 const solveScalar normFactor = this->normFactor(
psi, tsource(), wA, pA);
116 Info<<
" Normalisation factor = " << normFactor <<
endl;
120 solverPerf.initialResidual() =
123 solverPerf.finalResidual() = solverPerf.initialResidual();
129 || !solverPerf.checkConvergence(tolerance_, relTol_, log_)
133 solveScalar* __restrict__ pTPtr = pT.begin();
136 solveScalar* __restrict__ wTPtr = wT.begin();
139 matrix_.Tmul(wT,
psi, interfaceIntCoeffs_, interfaces_, cmpt);
143 solveScalar* __restrict__ rTPtr = rT.begin();
146 solveScalar wArT = 0;
162 const solveScalar wArTold = wArT;
165 preconPtr_->precondition(wA, rA, cmpt);
166 preconPtr_->preconditionT(wT, rT, cmpt);
171 if (solverPerf.nIterations() == 0)
173 for (label cell=0; cell<nCells; cell++)
175 pAPtr[cell] = wAPtr[cell];
176 pTPtr[cell] = wTPtr[cell];
181 const solveScalar
beta = wArT/wArTold;
183 for (label cell=0; cell<nCells; cell++)
185 pAPtr[cell] = wAPtr[cell] +
beta*pAPtr[cell];
186 pTPtr[cell] = wTPtr[cell] +
beta*pTPtr[cell];
192 matrix_.Amul(wA, pA, interfaceBouCoeffs_, interfaces_, cmpt);
193 matrix_.Tmul(wT, pT, interfaceIntCoeffs_, interfaces_, cmpt);
195 const solveScalar wApT =
gSumProd(wA, pT, matrix().
mesh().comm());
198 if (solverPerf.checkSingularity(
mag(wApT)/normFactor))
206 const solveScalar
alpha = wArT/wApT;
208 for (label cell=0; cell<nCells; cell++)
210 psiPtr[cell] +=
alpha*pAPtr[cell];
211 rAPtr[cell] -=
alpha*wAPtr[cell];
212 rTPtr[cell] -=
alpha*wTPtr[cell];
215 solverPerf.finalResidual() =
221 ++solverPerf.nIterations() < maxIter_
222 && !solverPerf.checkConvergence(tolerance_, relTol_, log_)
224 || solverPerf.nIterations() < minIter_
231 if (solverPerf.nIterations() > upperMaxIters)
234 <<
"PBiCG has failed to converge within the maximum iterations: " 235 << upperMaxIters <<
nl 236 <<
" Please try the more robust PBiCGStab solver." 242 preconPtr_->setFinished(solverPerf);
247 ConstPrecisionAdaptor<scalar, solveScalar>(rA)(),
void size(const label n)
Older name for setAddressableSize.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Base class for solution control classes.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Field< solveScalar > solveScalarField
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
void setResidualField(const scalarField &residual, const word &fieldName, const bool initial) const
Set the residual field using an IOField on the object registry if it exists.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
constexpr char nl
The newline '\n' character (0x0a)
scalarProduct< Type, Type >::type gSumProd(const UList< Type > &f1, const UList< Type > &f2, const label comm)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A non-const Field/List wrapper with possible data conversion.
static constexpr const label defaultMaxIter
Default maximum number of iterations for solvers (1000)
A field of fields is a PtrList of fields with reference counting.
virtual solverPerformance solve(scalarField &psi, const scalarField &source, const direction cmpt=0) const
Solve the matrix with this solver.
A class for handling words, derived from Foam::string.
static autoPtr< preconditioner > New(const solver &sol, const dictionary &solverControls)
Return a new preconditioner.
A const Field/List wrapper with possible data conversion.
lduMatrix::solver::addasymMatrixConstructorToTable< PBiCG > addPBiCGAsymMatrixConstructorToTable_
iterator begin() noexcept
Return an iterator to begin traversing the UList.
Preconditioned bi-conjugate gradient solver for asymmetric lduMatrices using a run-time selectable pr...
int debug
Static debugging option.
defineTypeNameAndDebug(combustionModel, 0)
typeOfMag< Type >::type gSumMag(const FieldField< Field, Type > &f)
lduMatrix is a general matrix class in which the coefficients are stored as three arrays...
messageStream Info
Information stream (stdout output on master, null elsewhere)
Internal & ref(const bool updateAccessTime=true)
Same as internalFieldRef()
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
static word getName(const dictionary &)
Find the preconditioner name (directly or from a sub-dictionary)
const volScalarField & psi
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].