45 namespace functionObjects
74 bool Foam::functionObjects::age::converged
77 const scalar initialResidual
80 if (initialResidual < tolerance_)
82 Info<<
"Field " << typeName
83 <<
" converged in " <<
nCorr <<
" correctors" 93 template<
class GeoField>
95 Foam::functionObjects::age::newField
105 scopedName(baseName),
139 phiName_ =
dict.getOrDefault<
word>(
"phi",
"phi");
140 rhoName_ =
dict.getOrDefault<
word>(
"rho",
"rho");
141 schemesField_ =
dict.getOrDefault<
word>(
"schemesField", typeName);
142 tolerance_ =
dict.getOrDefault<scalar>(
"tolerance", 1
e-5);
143 nCorr_ =
dict.getOrDefault<
int>(
"nCorr", 5);
144 diffusion_ =
dict.getOrDefault<
bool>(
"diffusion",
false);
160 mesh_.time().timeName(),
172 const word divScheme(
"div(phi," + schemesField_ +
")");
175 scalar relaxCoeff = 0;
176 if (mesh_.relaxEquation(schemesField_))
178 relaxCoeff = mesh_.equationRelaxationFactor(schemesField_);
196 tmp<volScalarField> tmuEff;
197 word laplacianScheme;
208 "laplacian(" + tmuEff().name() +
',' + schemesField_ +
")";
211 for (
int i = 0; i <= nCorr_; ++i)
223 ageEqn.relax(relaxCoeff);
227 if (converged(i, ageEqn.solve().initialResidual()))
237 tmp<volScalarField> tnuEff;
238 word laplacianScheme;
249 "laplacian(" + tnuEff().name() +
',' + schemesField_ +
")";
252 for (
int i = 0; i <= nCorr_; ++i)
265 ageEqn.relax(relaxCoeff);
269 if (converged(i, ageEqn.solve().initialResidual()))
278 Info<<
"Min/max age:" 279 <<
min(age).value() <<
' ' 284 word fieldName = typeName;
286 return store(fieldName, tage);
292 return writeObject(typeName);
virtual bool execute()
Execute.
virtual bool read(const dictionary &)
Read the data.
defineTypeNameAndDebug(ObukhovLength, 0)
fvMatrix< scalar > fvScalarMatrix
A list of keyword definitions, which are a keyword followed by a number of values (eg...
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
age(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
constexpr char nl
The newline '\n' character (0x0a)
wordList patchTypes(nPatches)
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
Ostream & endl(Ostream &os)
Add newline and flush stream.
Calculate the matrix for the laplacian of the field.
Generic dimensioned Type class.
IncompressibleTurbulenceModel< transportModel > turbulenceModel
Ignore writing from objectRegistry::writeObject()
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Macros for easy insertion into run-time selection tables.
bool read(const char *buf, int32_t &val)
Same as readInt32.
#define forAll(list, i)
Loop across all elements in list.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
static const char *const typeName
Typename for Field.
const dimensionedScalar e
Elementary charge.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from Foam::string.
label size() const noexcept
The number of elements in the list.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
Calculate the matrix for the divergence of the given field and flux.
List< word > wordList
A List of words.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
const polyBoundaryMesh & patches
Automatically write from objectRegistry::writeObject()
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Internal & ref(const bool updateAccessTime=true)
Same as internalFieldRef()
virtual bool read(const dictionary &dict)
Read optional controls.
virtual bool write()
Write.
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
Defines the attributes of an object for which implicit objectRegistry management is supported...
const fvMesh & mesh_
Reference to the fvMesh.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
static constexpr const zero Zero
Global zero (0)