Go to the documentation of this file.
44 namespace functionObjects
62 if (!foundObject<volScalarField>(fieldName_))
76 store(fieldName_, tfldPtr);
78 if (phaseName_ !=
"none")
84 return lookupObjectRef<volScalarField>(fieldName_);
94 const word Dname(
"D" +
s.name());
103 mesh_.time().timeName(),
113 if (nutName_ !=
"none")
124 findObject<incompressible::turbulenceModel>
134 alphaD_ *
turb->nu() + alphaDt_ *
turb->nut()
142 findObject<compressible::turbulenceModel>
152 alphaD_ *
turb->mu() + alphaDt_ *
turb->mut()
163 mesh_.time().timeName(),
176 Foam::functionObjects::scalarTransport::scalarTransport
184 fieldName_(
dict.getOrDefault<
word>(
"field",
"s")),
185 phiName_(
dict.getOrDefault<
word>(
"phi",
"phi")),
186 rhoName_(
dict.getOrDefault<
word>(
"rho",
"rho")),
187 nutName_(
dict.getOrDefault<
word>(
"nut",
"none")),
188 phaseName_(
dict.getOrDefault<
word>(
"phase",
"none")),
189 phasePhiCompressedName_
191 dict.getOrDefault<
word>(
"phasePhiCompressed",
"alphaPhiUn")
196 resetOnStartUp_(
false),
197 schemesField_(
"unknown-schemesField"),
199 bounded01_(
dict.getOrDefault(
"bounded01",
true))
226 dict.readIfPresent(
"phi", phiName_);
227 dict.readIfPresent(
"rho", rhoName_);
228 dict.readIfPresent(
"nut", nutName_);
229 dict.readIfPresent(
"phase", phaseName_);
230 dict.readIfPresent(
"bounded01", bounded01_);
232 schemesField_ =
dict.getOrDefault(
"schemesField", fieldName_);
233 constantD_ =
dict.readIfPresent(
"D", D_);
234 alphaD_ =
dict.getOrDefault<scalar>(
"alphaD", 1);
235 alphaDt_ =
dict.getOrDefault<scalar>(
"alphaDt", 1);
237 dict.readIfPresent(
"nCorr", nCorr_);
238 dict.readIfPresent(
"resetOnStartUp", resetOnStartUp_);
240 if (
dict.found(
"fvOptions"))
242 fvOptions_.reset(
dict.subDict(
"fvOptions"));
261 word divScheme(
"div(phi," + schemesField_ +
")");
262 word laplacianScheme(
"laplacian(" +
D.name() +
"," + schemesField_ +
")");
265 scalar relaxCoeff = 0.0;
266 if (mesh_.relaxEquation(schemesField_))
268 relaxCoeff = mesh_.equationRelaxationFactor(schemesField_);
272 if (phaseName_ !=
"none")
283 D.dimensions().reset(limitedPhiAlpha.dimensions()/
dimLength);
287 for (label i = 0; i <= nCorr_; i++)
292 +
fvm::div(limitedPhiAlpha,
s, divScheme)
298 sEqn.
relax(relaxCoeff);
299 fvOptions_.constrain(sEqn);
300 sEqn.
solve(mesh_.solverDict(schemesField_));
302 tTPhiUD = sEqn.
flux();
322 for (label i = 0; i <= nCorr_; i++)
334 sEqn.
relax(relaxCoeff);
336 fvOptions_.constrain(sEqn);
338 sEqn.
solve(mesh_.solverDict(schemesField_));
343 for (label i = 0; i <= nCorr_; i++)
354 sEqn.
relax(relaxCoeff);
356 fvOptions_.constrain(sEqn);
358 sEqn.
solve(mesh_.solverDict(schemesField_));
364 <<
"Incompatible dimensions for phi: " <<
phi.dimensions() <<
nl
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
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 dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
bool read(const char *buf, int32_t &val)
Same as readInt32.
static word timeName(const scalar t, const int precision=precision_)
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual bool execute()
Calculate the scalarTransport.
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool store(word &fieldName, const tmp< ObjectType > &tfield, bool cacheable=false)
Store the field in the (sub) objectRegistry under the given name.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual bool read(const dictionary &)
Read the scalarTransport data.
virtual ~scalarTransport()
Destructor.
Calculate the matrix for the divergence of the given field and flux.
void setFluxRequired(const word &name) const
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
word name(const complex &c)
Return string representation of complex.
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Calculate the matrix for the laplacian of the field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
virtual bool write()
Do nothing.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool read(const dictionary &dict)
Read optional controls.
Macros for easy insertion into run-time selection tables.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Calculate the matrix for implicit and explicit sources.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const fvMesh & mesh_
Reference to the fvMesh.
A class representing the concept of a field of 1 used to avoid unnecessary manipulations for objects ...
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
const dimensionedScalar & D
defineTypeNameAndDebug(ObukhovLength, 0)
A special matrix type and solver, designed for finite volume solutions of scalar equations....
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
const Time & time() const
Return the top-level database.
const dimensionSet dimVolume(pow3(dimLength))
compressible::turbulenceModel & turb
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Calculate the matrix for the first temporal derivative.
dimensionedScalar pos(const dimensionedScalar &ds)