Go to the documentation of this file.
38 lduMatrix::solver::addasymMatrixConstructorToTable<PBiCG>
47 const word& fieldName,
82 lduMatrix::preconditioner::getName(controlDict_) + typeName,
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;
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;
160 const solveScalar wArTold = wArT;
163 preconPtr->precondition(wA, rA, cmpt);
164 preconPtr->preconditionT(wT, rT, cmpt);
179 const solveScalar
beta = wArT/wArTold;
190 matrix_.Amul(wA, pA, interfaceBouCoeffs_, interfaces_, cmpt);
191 matrix_.Tmul(wT, pT, interfaceIntCoeffs_, interfaces_, cmpt);
193 const solveScalar wApT =
gSumProd(wA, pT, matrix().
mesh().comm());
204 const solveScalar
alpha = wArT/wApT;
230 <<
"PBiCG has failed to converge within the maximum number"
231 " of iterations " <<
max(defaultMaxIter_, maxIter_) <<
nl
232 <<
" Please try the more robust PBiCGStab solver."
int debug
Static debugging option.
typeOfMag< Type >::type gSumMag(const FieldField< Field, Type > &f)
A class for handling words, derived from Foam::string.
A field of fields is a PtrList of fields with reference counting.
lduMatrix::solver::addasymMatrixConstructorToTable< PBiCG > addPBiCGAsymMatrixConstructorToTable_
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
lduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
Abstract base-class for lduMatrix solvers.
Ostream & endl(Ostream &os)
Add newline and flush stream.
void setResidualField(const scalarField &residual, const word &fieldName, const bool initial) const
messageStream Info
Information stream (uses stdout - output is on the master only)
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
scalarProduct< Type, Type >::type gSumProd(const UList< Type > &f1, const UList< Type > &f2, const label comm)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A Field wrapper with possible data conversion.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
FieldType & ref()
Allow modification without const-ref check.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A const Field wrapper with possible data conversion.
const volScalarField & psi
virtual solverPerformance solve(scalarField &psi, const scalarField &source, const direction cmpt=0) const
Solve the matrix with this solver.
defineTypeNameAndDebug(combustionModel, 0)
A cell is defined as a list of faces with extra functionality.