Go to the documentation of this file.
45 if (psi_.needReference())
47 if (Pstream::master())
49 internalCoeffs_[patchi][facei] +=
50 diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]];
52 boundaryCoeffs_[patchi][facei] +=
53 diag()[psi_.mesh().boundary()[patchi].faceCells()[facei]]
68 if (psi_.mesh().name() != polyMesh::defaultRegion)
77 <<
"fvMatrix<scalar>::solver(const dictionary& solverControls) : "
78 "solver for fvMatrix<scalar>"
83 addBoundaryDiag(
diag(), 0);
96 psi_.boundaryField().scalarInterfaces(),
119 fvMat_.addBoundaryDiag(fvMat_.diag(), 0);
122 fvMat_.addBoundarySource(totalSource,
false);
125 solver_->read(solverControls);
129 psi.primitiveFieldRef(),
138 fvMat_.diag() = saveDiag;
140 psi.correctBoundaryConditions();
142 psi.mesh().setSolverPerformance(
psi.name(), solverPerf);
157 <<
"fvMatrix<scalar>::solveSegregated"
158 "(const dictionary& solverControls) : "
159 "solving fvMatrix<scalar>"
167 addBoundaryDiag(
diag(), 0);
170 addBoundarySource(totalSource,
false);
179 psi_.boundaryField().scalarInterfaces(),
181 )->
solve(
psi.primitiveFieldRef(), totalSource);
190 psi.correctBoundaryConditions();
192 psi.mesh().setSolverPerformance(
psi.name(), solverPerf);
202 addBoundaryDiag(boundaryDiag, 0);
213 source_ - boundaryDiag*psif,
215 psi_.boundaryField().scalarInterfaces(),
221 addBoundarySource(tres_s.
ref());
236 "H("+psi_.name()+
')',
244 extrapolatedCalculatedFvPatchScalarField::typeName
275 dimensions_/(
dimVol*psi_.dimensions()),
276 extrapolatedCalculatedFvPatchScalarField::typeName
281 H1_.primitiveFieldRef() = lduMatrix::H1();
284 H1_.primitiveFieldRef() /= psi_.mesh().V();
285 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)
Convert to OSstream.
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.
messageStream Info
Information stream (uses stdout - output is on the master only)
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 wrapper with possible data conversion.
Solver class returned by the solver function.
const dimensionSet dimVol(dimVolume)
Older spelling for dimVolume.
SolverPerformance< Type > solveSegregated(const dictionary &)
Solve segregated returning the solution statistics.
const volScalarField & psi
Container< Type > & ref() const
void setComponentReference(const label patchi, const label facei, const direction cmpt, const scalar value)
Set reference level for a component of the solution.