31#include "phaseCompressibleTurbulenceModel.H"
42#include "surfaceInterpolate.H"
47template<
class BasePhaseModel>
51 word phiName(IOobject::groupName(
"phi", this->
name()));
56 U.mesh().time().timeName(),
61 if (phiHeader.typeHeaderOk<surfaceScalarField>(
true))
63 Info<<
"Reading face flux field " << phiName <<
endl;
65 return tmp<surfaceScalarField>
67 new surfaceScalarField
72 U.mesh().time().timeName(),
83 Info<<
"Calculating face flux field " << phiName <<
endl;
87 U.boundaryField().size(),
88 calculatedFvPatchScalarField::typeName
95 isA<fixedValueFvPatchVectorField>(
U.boundaryField()[i])
96 || isA<slipFvPatchVectorField>(
U.boundaryField()[i])
97 || isA<partialSlipFvPatchVectorField>(
U.boundaryField()[i])
104 return tmp<surfaceScalarField>
106 new surfaceScalarField
111 U.mesh().time().timeName(),
126template<
class BasePhaseModel>
130 const word& phaseName,
134 BasePhaseModel(
fluid, phaseName, index),
189 IOobject::groupName(
"continuityErrorFlow", this->
name()),
196 continuityErrorSources_
200 IOobject::groupName(
"continuityErrorSources", this->
name()),
217template<
class BasePhaseModel>
220 BasePhaseModel::correct();
232template<
class BasePhaseModel>
235 BasePhaseModel::correctKinematics();
251 K_.ref() = 0.5*
magSqr(this->
U());
256template<
class BasePhaseModel>
259 BasePhaseModel::correctTurbulence();
261 turbulence_->correct();
265template<
class BasePhaseModel>
268 BasePhaseModel::correctEnergyTransport();
270 turbulence_->correctEnergyTransport();
274template<
class BasePhaseModel>
281template<
class BasePhaseModel>
292 +
fvm::SuSp(- this->continuityError(), U_)
294 + turbulence_->divDevRhoReff(U_)
299template<
class BasePhaseModel>
312 +
fvm::SuSp(- this->continuityErrorSources(), U_)
314 + turbulence_->divDevRhoReff(U_)
319template<
class BasePhaseModel>
327template<
class BasePhaseModel>
335template<
class BasePhaseModel>
343template<
class BasePhaseModel>
351template<
class BasePhaseModel>
359template<
class BasePhaseModel>
367template<
class BasePhaseModel>
375template<
class BasePhaseModel>
383template<
class BasePhaseModel>
396template<
class BasePhaseModel>
402 DUDtf_ =
byDt(phi_ - phi_.oldTime());
409template<
class BasePhaseModel>
413 return continuityErrorFlow_ + continuityErrorSources_;
417template<
class BasePhaseModel>
421 return continuityErrorFlow_;
425template<
class BasePhaseModel>
429 return continuityErrorSources_;
433template<
class BasePhaseModel>
450template<
class BasePhaseModel>
458template<
class BasePhaseModel>
465template<
class BasePhaseModel>
469 return turbulence_->mut();
473template<
class BasePhaseModel>
477 return turbulence_->muEff();
481template<
class BasePhaseModel>
485 return turbulence_->nut();
489template<
class BasePhaseModel>
493 return turbulence_->nuEff();
497template<
class BasePhaseModel>
501 return turbulence_->kappaEff();
505template<
class BasePhaseModel>
509 return turbulence_->kappaEff(patchi);
513template<
class BasePhaseModel>
517 return turbulence_->alphaEff();
521template<
class BasePhaseModel>
525 return turbulence_->alphaEff(patchi);
529template<
class BasePhaseModel>
533 return turbulence_->k();
537template<
class BasePhaseModel>
541 return turbulence_->pPrime();
Defines the attributes of an object for which implicit objectRegistry management is supported,...
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
writeOption writeOpt() const noexcept
The write option.
void correctBoundaryVelocity(volVectorField &U) const
Correct the boundary velocity for the rotation of the MRF region.
Class which represents a moving fluid phase. Holds the velocity, fluxes and turbulence model....
virtual void correctKinematics()
Correct the kinematics.
virtual tmp< volScalarField > continuityErrorSources() const
Return the continuity error due to any sources.
virtual void correctEnergyTransport()
Correct the energy transport e.g. alphat.
virtual tmp< fvVectorMatrix > UEqn()
Return the momentum equation.
virtual tmp< surfaceScalarField > DUDtf() const
Return the substantive acceleration on the faces.
virtual tmp< volScalarField > muEff() const
Return the effective dynamic viscosity.
virtual void correct()
Correct the phase properties other than the thermo and turbulence.
virtual tmp< volScalarField > k() const
Return the turbulent kinetic energy.
virtual tmp< volVectorField > U() const
Access const reference to U.
virtual tmp< volScalarField > mut() const
Return the turbulent dynamic viscosity.
virtual surfaceScalarField & alphaRhoPhiRef()
Access the mass flux of the phase.
virtual void correctTurbulence()
Correct the turbulence.
virtual surfaceScalarField & phiRef()
Access the volumetric flux.
virtual tmp< volScalarField > alphaEff() const
Return the effective thermal diffusivity.
virtual tmp< surfaceScalarField > alphaPhi() const
Constant access the volumetric flux of the phase.
virtual tmp< volScalarField > K() const
Return the phase kinetic energy.
virtual tmp< volScalarField > pPrime() const
Return the phase-pressure'.
virtual tmp< volScalarField > continuityErrorFlow() const
Return the continuity error due to the flow field.
virtual tmp< volScalarField > nuEff() const
Return the effective kinematic viscosity.
virtual tmp< volScalarField > divU() const
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
virtual surfaceScalarField & alphaPhiRef()
Access the volumetric flux of the phase.
virtual tmp< surfaceScalarField > alphaRhoPhi() const
Return the mass flux of the phase.
virtual tmp< surfaceScalarField > phi() const
Constant access the volumetric flux.
virtual tmp< fvVectorMatrix > UfEqn()
Return the momentum equation for the face-based algorithm.
virtual tmp< volScalarField > nut() const
Return the turbulent kinematic viscosity.
virtual volVectorField & URef()
Access the velocity.
virtual tmp< volVectorField > DUDt() const
Return the substantive acceleration.
virtual tmp< volScalarField > kappaEff() const
Return the effective thermal conductivity.
virtual bool stationary() const
Return whether the phase is stationary.
virtual tmp< volScalarField > continuityError() const
Return the continuity error.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Class to represent a system of phases and model interfacial transfers between them.
const IOMRFZoneList & MRF() const
Return MRF zones.
fv::options & fvOptions() const
Access the fvOptions.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Calculate the face-flux of the given field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the finiteVolume matrix for implicit and explicit sources.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
zeroField SuSp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
List< word > wordList
A List of words.
tmp< volScalarField > byDt(const volScalarField &vf)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const dimensionSet dimDensity
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
#define forAll(list, i)
Loop across all elements in list.
static const char *const typeName
The type name used in ensight case files.