28 #include "phaseCompressibleTurbulenceModel.H" 39 #include "surfaceInterpolate.H" 44 template<
class BasePhaseModel>
48 word phiName(IOobject::groupName(
"phi", this->
name()));
53 U.mesh().time().timeName(),
60 Info<<
"Reading face flux field " << phiName <<
endl;
62 return tmp<surfaceScalarField>
69 U.mesh().time().timeName(),
80 Info<<
"Calculating face flux field " << phiName <<
endl;
84 U.boundaryField().size(),
85 calculatedFvPatchScalarField::typeName
88 forAll(U.boundaryField(), i)
92 isA<fixedValueFvPatchVectorField>(U.boundaryField()[i])
93 || isA<slipFvPatchVectorField>(U.boundaryField()[i])
94 || isA<partialSlipFvPatchVectorField>(U.boundaryField()[i])
97 phiTypes[i] = fixedValueFvPatchScalarField::typeName;
101 return tmp<surfaceScalarField>
108 U.mesh().time().timeName(),
123 template<
class BasePhaseModel>
126 const phaseSystem& fluid,
127 const word& phaseName,
131 BasePhaseModel(fluid, phaseName, index),
136 IOobject::groupName(
"U", this->
name()),
149 IOobject::groupName(
"alphaPhi", this->
name()),
160 IOobject::groupName(
"alphaRhoPhi", this->
name()),
186 IOobject::groupName(
"continuityErrorFlow", this->
name()),
193 continuityErrorSources_
197 IOobject::groupName(
"continuityErrorSources", this->
name()),
214 template<
class BasePhaseModel>
221 template<
class BasePhaseModel>
228 continuityErrorSources_ = - (this->
fluid().fvOptions()(*
this,
rho)&rho);
232 template<
class BasePhaseModel>
236 this->
fluid().MRF().correctBoundaryVelocity(U_);
237 correctContinuityErrors();
241 template<
class BasePhaseModel>
244 BasePhaseModel::correctKinematics();
266 template<
class BasePhaseModel>
269 correctContinuityErrors();
273 template<
class BasePhaseModel>
276 BasePhaseModel::correctTurbulence();
278 turbulence_->correct();
282 template<
class BasePhaseModel>
285 BasePhaseModel::correctEnergyTransport();
287 turbulence_->correctEnergyTransport();
291 template<
class BasePhaseModel>
298 template<
class BasePhaseModel>
309 +
fvm::SuSp(- this->continuityError(), U_)
311 + turbulence_->divDevRhoReff(U_)
316 template<
class BasePhaseModel>
329 +
fvm::SuSp(- this->continuityErrorSources(), U_)
331 + turbulence_->divDevRhoReff(U_)
336 template<
class BasePhaseModel>
344 template<
class BasePhaseModel>
352 template<
class BasePhaseModel>
360 template<
class BasePhaseModel>
368 template<
class BasePhaseModel>
376 template<
class BasePhaseModel>
384 template<
class BasePhaseModel>
392 template<
class BasePhaseModel>
400 template<
class BasePhaseModel>
409 return tmp<volVectorField>(DUDt_());
413 template<
class BasePhaseModel>
419 DUDtf_ =
byDt(phi_ - phi_.oldTime());
422 return tmp<surfaceScalarField>(DUDtf_());
426 template<
class BasePhaseModel>
430 return continuityErrorFlow_ + continuityErrorSources_;
434 template<
class BasePhaseModel>
438 return continuityErrorFlow_;
442 template<
class BasePhaseModel>
446 return continuityErrorSources_;
450 template<
class BasePhaseModel>
463 return tmp<volScalarField>(K_());
467 template<
class BasePhaseModel>
471 return divU_.valid() ? tmp<volScalarField>(divU_()) : tmp<volScalarField>();
475 template<
class BasePhaseModel>
482 template<
class BasePhaseModel>
486 return turbulence_->mut();
490 template<
class BasePhaseModel>
494 return turbulence_->muEff();
498 template<
class BasePhaseModel>
502 return turbulence_->nut();
506 template<
class BasePhaseModel>
510 return turbulence_->nuEff();
514 template<
class BasePhaseModel>
518 return turbulence_->kappaEff();
522 template<
class BasePhaseModel>
526 return turbulence_->kappaEff(patchi);
530 template<
class BasePhaseModel>
534 return turbulence_->alphaEff();
538 template<
class BasePhaseModel>
542 return turbulence_->alphaEff(patchi);
546 template<
class BasePhaseModel>
550 return turbulence_->k();
554 template<
class BasePhaseModel>
558 return turbulence_->pPrime();
#define forAll(list, i)
Loop across all elements in list.
virtual tmp< volScalarField > divU() const
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
tmp< GeometricField< Type, fvPatchField, volMesh > > DDt(const surfaceScalarField &phi, const GeometricField< Type, fvPatchField, volMesh > &psi)
ThermalDiffusivity< PhaseCompressibleTurbulenceModel< phaseModel > > phaseCompressibleTurbulenceModel
Typedef for phaseCompressibleTurbulenceModel.
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
virtual surfaceScalarField & alphaRhoPhiRef()
Access the mass flux of the phase.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual tmp< volScalarField > alphaEff() const
Return the effective thermal diffusivity.
virtual bool stationary() const
Return whether the phase is stationary.
MovingPhaseModel(const phaseSystem &fluid, const word &phaseName, const label index)
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual tmp< volScalarField > mut() const
Return the turbulent dynamic viscosity.
virtual tmp< volVectorField > DUDt() const
Return the substantive acceleration.
virtual tmp< surfaceScalarField > alphaPhi() const
Return the volumetric flux of the phase.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
rhoReactionThermo & thermo
virtual void correctTurbulence()
Correct the turbulence.
GeometricField< vector, fvPatchField, volMesh > volVectorField
CGAL::Exact_predicates_exact_constructions_kernel K
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvPatchField< scalar >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
virtual tmp< volScalarField > k() const
Return the turbulent kinetic energy.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Calculate the first temporal derivative.
Class which represents a moving fluid phase. Holds the velocity, fluxes and turbulence model and can ...
virtual void correctThermo()
Correct the thermodynamics.
virtual tmp< surfaceScalarField > alphaRhoPhi() const
Return the mass flux of the phase.
virtual tmp< fvVectorMatrix > UfEqn()
Return the momentum equation for the face-based algorithm.
virtual ~MovingPhaseModel()
Destructor.
Calculate the face-flux of the given field.
static word groupName(Name name, const word &group)
virtual surfaceScalarField & phiRef()
Access the volumetric flux.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Calculate the matrix for the first temporal derivative.
virtual tmp< volScalarField > muEff() const
Return the effective dynamic viscosity.
virtual tmp< volScalarField > continuityErrorSources() const
Return the continuity error due to any sources.
virtual tmp< volScalarField > continuityErrorFlow() const
Return the continuity error due to the flow field.
virtual volVectorField & URef()
Access the velocity.
virtual tmp< volScalarField > pPrime() const
Return the phase-pressure'.
Calculate the divergence of the given field.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual tmp< volScalarField > nut() const
Return the turbulent kinematic viscosity.
virtual tmp< volScalarField > kappaEff() const
Return the effective thermal conductivity.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-7/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
virtual tmp< volScalarField > nuEff() const
Return the effective kinematic viscosity.
virtual void correct()
Correct the phase properties other than the thermo and turbulence.
virtual tmp< volVectorField > U() const
Return the velocity.
word name(const complex &)
Return a string representation of a complex.
virtual void correctEnergyTransport()
Correct the energy transport e.g. alphat.
tmp< volScalarField > byDt(const volScalarField &vf)
const dimensionSet dimDensity
virtual tmp< volScalarField > continuityError() const
Return the continuity error.
Calculate the matrix for the divergence of the given field and flux.
List< word > wordList
A List of words.
virtual tmp< volScalarField > K() const
Return the phase kinetic energy.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual surfaceScalarField & alphaPhiRef()
Access the volumetric flux of the phase.
virtual tmp< fvVectorMatrix > UEqn()
Return the momentum equation.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
virtual void correctKinematics()
Correct the kinematics.
virtual tmp< surfaceScalarField > phi() const
Return the volumetric flux.
A class for managing temporary objects.
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual tmp< surfaceScalarField > DUDtf() const
Return the substantive acceleration on the faces.
Calculate the matrix for implicit and explicit sources.