Go to the documentation of this file.
37 namespace incompressible
51 Info<<
"Storing initial values of turbulence variables" <<
endl;
92 Info<<
"Allocating mean values of turbulence variables" <<
endl;
160 <<
"Cloning " << sf.name() <<
endl;
181 const word name1 = f1.name();
182 const word name2 = f2.name();
200 solverControl_(SolverControl),
212 TMVar1InitPtr_(
nullptr),
213 TMVar2InitPtr_(
nullptr),
214 nutInitPtr_(
nullptr),
215 TMVar1MeanPtr_(
nullptr),
216 TMVar2MeanPtr_(
nullptr),
234 nutPtr_(cloneAutoTmp(rmv.
nutPtr_)),
235 dPtr_(cloneAutoTmp(rmv.
dPtr_)),
239 TMVar1InitPtr_(
nullptr),
240 TMVar2InitPtr_(
nullptr),
241 nutInitPtr_(
nullptr),
242 TMVar1MeanPtr_(
nullptr),
243 TMVar2MeanPtr_(
nullptr),
267 mesh.time().constant(),
277 const word modelType(
dict.getOrDefault<
word>(
"RASModel",
"laminar"));
279 Info<<
"Creating references for RASModel variables : " << modelType <<
endl;
281 auto cstrIter = dictionaryConstructorTablePtr_->cfind(modelType);
283 if (!cstrIter.found())
290 *dictionaryConstructorTablePtr_
426 return dPtr_().constCast();
472 <<
"jutJacobianVar1 not implemented for the current turbulence model."
473 <<
"Returning zero field" <<
endl;
482 mesh_.time().timeName(),
501 <<
"nutJacobianVar2 not implemented for the current turbulence model."
502 <<
"Returning zero field" <<
endl;
511 mesh_.time().timeName(),
547 Info<<
"Reseting mean turbulent fields to zero" <<
endl;
573 scalar avIter(iAverageIter);
574 scalar oneOverItP1 = 1./(avIter + 1);
575 scalar mult = avIter*oneOverItP1;
607 mesh_.time().timeName(),
625 TMVar1Ptr_().constCast().correctBoundaryConditions();
626 if (solverControl_.average())
628 TMVar1MeanPtr_().correctBoundaryConditions();
634 TMVar2Ptr_().constCast().correctBoundaryConditions();
635 if (solverControl_.average())
637 TMVar2MeanPtr_().correctBoundaryConditions();
643 nutPtr_().constCast().correctBoundaryConditions();
644 if (solverControl_.average())
646 nutMeanPtr_().correctBoundaryConditions();
static autoPtr< T > New(Args &&... args)
Construct autoPtr of T with forwarding arguments.
autoPtr< volScalarField > nutInitPtr_
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
void reset(T *p=nullptr) noexcept
Delete managed object and set to new given pointer.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & nutBaseName() const
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
bool doAverageIter() const
A class for handling words, derived from Foam::string.
const volScalarField & d() const
const volScalarField & TMVar1Inst() const
return references to instantaneous turbulence fields
const solverControl & solverControl_
defineTypeNameAndDebug(adjointEikonalSolver, 0)
Info<< "Reading field U\n"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\n"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Abstract base class for objective functions. No point in making this runTime selectable since its chi...
autoPtr< volScalarField > TMVar1InitPtr_
static constexpr const zero Zero
Global zero.
autoPtr< volScalarField > TMVar1MeanPtr_
static word timeName(const scalar t, const int precision=precision_)
static const word propertiesName
Default name of the turbulence properties dictionary.
void restoreInitValues()
Restore turbulent fields to their initial values.
bool valid() const noexcept
True if the managed pointer is non-null.
Base class for solver control classes.
void allocateInitValues()
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool storeInitValues() const
Re-initialize.
A simple single-phase transport model based on viscosityModel.
virtual void correctBoundaryConditions(const incompressible::turbulenceModel &turbulence)
correct bounday conditions of turbulent fields
defineRunTimeSelectionTable(adjointSensitivity, dictionary)
RASModelVariables(const fvMesh &mesh, const solverControl &SolverControl)
Construct from components.
void computeMeanFields()
Compute mean fields on the fly.
autoPtr< volScalarField > nutMeanPtr_
virtual void transfer(RASModelVariables &rmv)
Transfer turbulence fields from an another object.
const volScalarField & nutRefInst() const
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void allocateMeanFields()
void resetMeanFields()
Reset mean fields to zero.
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
virtual tmp< volScalarField > nutJacobianVar2(const singlePhaseTransportModel &laminarTransport) const
bool useAveragedFields() const
static autoPtr< RASModelVariables > New(const fvMesh &mesh, const solverControl &SolverControl)
Return a reference to the selected turbulence model.
autoPtr< volScalarField > TMVar2MeanPtr_
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
label & averageIter()
Return average iteration index reference.
const volScalarField & TMVar2Inst() const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
const volScalarField & nutRef() const
bool hasTMVar1() const
Bools to idenify which turbulent fields are present.
Mesh data needed to do the Finite Volume discretisation.
dictionary subOrEmptyDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX, const bool mandatory=false) const
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< volSymmTensorField > devReff(const singlePhaseTransportModel &laminarTransport, const volVectorField &U) const
Return stress tensor based on the mean flow variables.
#define DebugInfo
Report an information message using Foam::Info.
autoTmp cloneAutoTmp(const autoTmp &source) const
const word & TMVar1BaseName() const
Turbulence field names.
void copyAndRename(volScalarField &f1, volScalarField &f2)
const word & TMVar2BaseName() const
static const word null
An empty word.
const Time & time() const
Return the top-level database.
const volScalarField & TMVar1() const
Return references to turbulence fields.
Templated abstract base class for single-phase incompressible turbulence models.
virtual tmp< volScalarField > nutJacobianVar1(const singlePhaseTransportModel &laminarTransport) const
Return nut Jacobian wrt the TM vars.
autoPtr< volScalarField > TMVar2InitPtr_
const volScalarField & TMVar2() const
#define WarningInFunction
Report a warning using Foam::Warning.
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
bool average() const
Whether averaging is enabled or not.
singlePhaseTransportModel laminarTransport(U, phi)
autoPtr< RASModelVariables > clone() const
Clone.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)