Go to the documentation of this file.
38 lduMatrix::solver::addsymMatrixConstructorToTable<PPCG>
45 void Foam::PPCG::gSumMagProd
47 FixedList<solveScalar, 3>& globalSum,
52 label& outstandingRequest,
56 const label nCells = a.size();
59 for (label cell=0; cell<nCells; ++cell)
61 globalSum[0] += a[cell]*
b[cell];
62 globalSum[1] += a[cell]*
c[cell];
96 const label comm = matrix().
mesh().
comm();
97 const label nCells =
psi.size();
101 matrix_.Amul(w,
psi, interfaceBouCoeffs_, interfaces_, cmpt);
108 const solveScalar normFactor = this->normFactor(
psi, source, w,
p);
112 Info<<
" Normalisation factor = " << normFactor <<
endl;
128 matrix_.Amul(w, u, interfaceBouCoeffs_, interfaces_, cmpt);
139 label outstandingRequest = -1;
143 gSumMagProd(globalSum, u, r, w, r, outstandingRequest, comm);
154 gSumMagProd(globalSum, w, u, m, r, outstandingRequest, comm);
159 matrix_.Amul(
n, m, interfaceBouCoeffs_, interfaces_, cmpt);
161 solveScalar
alpha = 0.0;
162 solveScalar
gamma = 0.0;
176 outstandingRequest = -1;
179 const solveScalar gammaOld =
gamma;
180 gamma = globalSum[0];
181 const solveScalar
delta = globalSum[1];
192 (minIter_ <= 0 || solverPerf.
nIterations() >= minIter_)
233 gSumMagProd(globalSum, u, r, w, r, outstandingRequest, comm);
244 gSumMagProd(globalSum, w, u, m, r, outstandingRequest, comm);
248 matrix_.Amul(
n, m, interfaceBouCoeffs_, interfaces_, cmpt);
260 const word& fieldName,
int debug
Static debugging option.
static word getName(const dictionary &)
Find the preconditioner name (directly or from a sub-dictionary)
A class for handling words, derived from Foam::string.
A field of fields is a PtrList of fields with reference counting.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
lduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
static bool & parRun()
Is this a parallel run?
static autoPtr< preconditioner > New(const solver &sol, const dictionary &solverControls)
Return a new preconditioner.
dimensioned< typename typeOfMag< Type >::type > sumMag(const DimensionedField< Type, GeoMesh > &df)
Abstract base-class for lduMatrix solvers.
Ostream & endl(Ostream &os)
Add newline and flush stream.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
messageStream Info
Information stream (uses stdout - output is on the master only)
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
static void waitRequest(const label i)
Wait until request i has finished.
lduMatrix::solver::addsymMatrixConstructorToTable< PPCG > addPPCGSymMatrixConstructorToTable_
Field< solveScalar > solveScalarField
virtual void precondition(solveScalarField &wA, const solveScalarField &rA, const direction cmpt=0) const =0
Return wA the preconditioned form of residual rA.
virtual solverPerformance solve(scalarField &psi, const scalarField &source, const direction cmpt=0) const
Solve the matrix with this solver.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A Field wrapper with possible data conversion.
static int & msgType()
Message tag of standard messages.
FieldType & ref()
Allow modification without const-ref check.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
A 1D vector of objects of type <T> with a fixed length <N>.
A const Field wrapper with possible data conversion.
const lduMesh & mesh() const
Return the LDU mesh from which the addressing is obtained.
virtual label comm() const =0
Return communicator used for parallel communication.
const dimensionedScalar c
Speed of light in a vacuum.
solverPerformance scalarSolveCG(solveScalarField &psi, const solveScalarField &source, const direction cmpt, const bool cgMode) const
const volScalarField & psi
defineTypeNameAndDebug(combustionModel, 0)
A cell is defined as a list of faces with extra functionality.