Go to the documentation of this file.
48 if (psi_.needReference())
50 if (Pstream::master())
52 internalCoeffs_[patchi][facei] +=
53 diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]];
55 boundaryCoeffs_[patchi][facei] +=
56 diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]]
71 if (psi_.mesh().name() != polyMesh::defaultRegion)
80 <<
"fvMatrix<scalar>::solver(const dictionary& solverControls) : "
81 "solver for fvMatrix<scalar>"
86 addBoundaryDiag(
diag(), 0);
89 psi_.boundaryField().scalarInterfaces();
135 fvMat_.addBoundaryDiag(fvMat_.diag(), 0);
138 fvMat_.addBoundarySource(totalSource,
false);
141 solver_->read(solverControls);
145 psi.primitiveFieldRef(),
154 fvMat_.diag() = saveDiag;
156 psi.correctBoundaryConditions();
158 psi.mesh().setSolverPerformance(
psi.name(), solverPerf);
173 <<
"fvMatrix<scalar>::solveSegregated"
174 "(const dictionary& solverControls) : "
175 "solving fvMatrix<scalar>"
191 createOrUpdateLduPrimitiveAssembly();
193 if (psi_.mesh().fluxRequired(psi_.name()))
203 setLduMesh(*lduMeshPtr());
204 transferFvMatrixCoeffs();
205 setBounAndInterCoeffs();
207 manipulateMatrix(cmpt);
211 addBoundaryDiag(
diag(), 0);
214 addBoundarySource(totalSource,
false);
220 interfaces = this->
psi(0).boundaryField().scalarInterfaces();
224 setInterfaces(interfaces, newInterfaces);
243 for (label fieldi = 0; fieldi < nMatrices(); fieldi++)
245 const label cellOffset = lduMeshPtr()->cellOffsets()[fieldi];
246 const auto& psiInternal = this->
psi(fieldi).primitiveField();
248 forAll(psiInternal, localCellI)
250 psi[cellOffset + localCellI] = psiInternal[localCellI];
265 )->solve(
psi, totalSource);
269 for (label fieldi = 0; fieldi < nMatrices(); fieldi++)
277 const label cellOffset = lduMeshPtr()->cellOffsets()[fieldi];
279 forAll(psiInternal, localCellI)
281 psiInternal[localCellI] =
psi[localCellI + cellOffset];
295 if (psi_.mesh().fluxRequired(psi_.name()))
300 lower().setSize(saveLower.size());
304 upper().setSize(saveUpper.size());
308 setLduMesh(psi_.mesh());
311 for (label fieldi = 0; fieldi < nMatrices(); fieldi++)
320 localPsi.mesh().setSolverPerformance(localPsi.name(), solverPerf);
331 addBoundaryDiag(boundaryDiag, 0);
342 source_ - boundaryDiag*psif,
344 psi_.boundaryField().scalarInterfaces(),
350 addBoundarySource(tres_s.
ref());
365 "H("+psi_.name()+
')',
373 extrapolatedCalculatedFvPatchScalarField::typeName
404 dimensions_/(
dimVol*psi_.dimensions()),
405 extrapolatedCalculatedFvPatchScalarField::typeName
410 H1_.primitiveFieldRef() = lduMatrix::H1();
413 H1_.primitiveFieldRef() /= psi_.mesh().V();
414 H1_.correctBoundaryConditions();
int debug
Static debugging option.
tmp< Field< Type > > residual() const
Return the matrix residual.
tmp< GeometricField< Type, fvPatchField, volMesh > > H() const
Return the H operation source.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
OSstream & masterStream(const label communicator)
SolverPerformance< Type > solve()
Solve returning the solution statistics.
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
Ostream & endl(Ostream &os)
Add newline and flush stream.
A dynamically resizable PtrList with allocation management.
conserve primitiveFieldRef()+
#define forAll(list, i)
Loop across all elements in list.
string lower(const std::string &s)
Return string copy transformed with std::tolower on each character.
messageStream Info
Information stream (stdout output on master, null elsewhere)
tmp< volScalarField > H1() const
Return H(1)
#define addProfiling(name, descr)
Define profiling trigger with specified name and description string.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
A scalar instance of fvMatrix.
autoPtr< fvSolver > solver()
Construct and return the solver.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
void correctBoundaryConditions()
Correct boundary field.
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.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
A const Field/List wrapper with possible data conversion.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const dimensionSet dimVol(dimVolume)
Older spelling for dimVolume.
string upper(const std::string &s)
Return string copy transformed with std::toupper on each character.
Container< Type > & ref() const
SolverPerformance< Type > solveSegregated(const dictionary &)
Solve segregated returning the solution statistics.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
const volScalarField & psi
void setComponentReference(const label patchi, const label facei, const direction cmpt, const scalar value)