Go to the documentation of this file.
50 psiThermo(
U.
mesh(), word::null),
51 twoPhaseMixture(
U.
mesh(), *this),
52 interfaceProperties(
alpha1(),
U, *this),
90 thermo1_->he() = thermo1_->he(p_, T_);
93 thermo2_->he() = thermo2_->he(p_, T_);
100 psi_ =
alpha1()*thermo1_->psi() +
alpha2()*thermo2_->psi();
101 mu_ =
alpha1()*thermo1_->mu() +
alpha2()*thermo2_->mu();
102 alpha_ =
alpha1()*thermo1_->alpha() +
alpha2()*thermo2_->alpha();
110 return thermo1_->thermoName() +
',' + thermo2_->thermoName();
116 return thermo1_->incompressible() && thermo2_->incompressible();
122 return thermo1_->isochoric() && thermo2_->isochoric();
157 alpha1().boundaryField()[patchi]*thermo1_->he(
p,
T, patchi)
158 +
alpha2().boundaryField()[patchi]*thermo2_->he(
p,
T, patchi);
164 return alpha1()*thermo1_->hc() +
alpha2()*thermo2_->hc();
196 return alpha1()*thermo1_->Cp() +
alpha2()*thermo2_->Cp();
208 alpha1().boundaryField()[patchi]*thermo1_->Cp(
p,
T, patchi)
209 +
alpha2().boundaryField()[patchi]*thermo2_->Cp(
p,
T, patchi);
215 return alpha1()*thermo1_->Cv() +
alpha2()*thermo2_->Cv();
227 alpha1().boundaryField()[patchi]*thermo1_->Cv(
p,
T, patchi)
228 +
alpha2().boundaryField()[patchi]*thermo2_->Cv(
p,
T, patchi);
234 return alpha1()*thermo1_->gamma() +
alpha2()*thermo2_->gamma();
246 alpha1().boundaryField()[patchi]*thermo1_->gamma(
p,
T, patchi)
247 +
alpha2().boundaryField()[patchi]*thermo2_->gamma(
p,
T, patchi);
253 return alpha1()*thermo1_->Cpv() +
alpha2()*thermo2_->Cpv();
265 alpha1().boundaryField()[patchi]*thermo1_->Cpv(
p,
T, patchi)
266 +
alpha2().boundaryField()[patchi]*thermo2_->Cpv(
p,
T, patchi);
273 alpha1()*thermo1_->CpByCpv()
274 +
alpha2()*thermo2_->CpByCpv();
286 alpha1().boundaryField()[patchi]*thermo1_->CpByCpv(
p,
T, patchi)
287 +
alpha2().boundaryField()[patchi]*thermo2_->CpByCpv(
p,
T, patchi);
299 return mu()/(
alpha1()*thermo1_->rho() +
alpha2()*thermo2_->rho());
311 alpha1().boundaryField()[patchi]*thermo1_->rho(patchi)
312 +
alpha2().boundaryField()[patchi]*thermo2_->rho(patchi)
319 return alpha1()*thermo1_->kappa() +
alpha2()*thermo2_->kappa();
329 alpha1().boundaryField()[patchi]*thermo1_->kappa(patchi)
330 +
alpha2().boundaryField()[patchi]*thermo2_->kappa(patchi);
337 alpha1()*thermo1_->alphahe()
338 +
alpha2()*thermo2_->alphahe();
348 alpha1().boundaryField()[patchi]*thermo1_->alphahe(patchi)
349 +
alpha2().boundaryField()[patchi]*thermo2_->alphahe(patchi);
359 alpha1()*thermo1_->kappaEff(alphat)
360 +
alpha2()*thermo2_->kappaEff(alphat);
371 alpha1().boundaryField()[patchi]*thermo1_->kappaEff(alphat, patchi)
372 +
alpha2().boundaryField()[patchi]*thermo2_->kappaEff(alphat, patchi);
382 alpha1()*thermo1_->alphaEff(alphat)
383 +
alpha2()*thermo2_->alphaEff(alphat);
394 alpha1().boundaryField()[patchi]*thermo1_->alphaEff(alphat, patchi)
395 +
alpha2().boundaryField()[patchi]*thermo2_->alphaEff(alphat, patchi);
List< label > labelList
A List of labels.
virtual tmp< volScalarField > alphahe() const
Thermal diffusivity for energy of mixture [kg/m/s].
virtual bool read()
Read base transportProperties dictionary.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual tmp< volScalarField > kappaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
A class for handling words, derived from Foam::string.
const dimensionedScalar mu
Atomic mass unit.
A class for managing temporary objects.
virtual tmp< volScalarField > Cv() const
Heat capacity at constant volume [J/kg/K].
virtual tmp< volScalarField > gamma() const
Gamma = Cp/Cv [].
const volScalarField & alpha2
virtual bool read()
Read thermophysical properties dictionary.
virtual void flush() const
Forcibly wait until all output done. Flush any cached data.
virtual bool isochoric() const
Return true if the equation of state is isochoric.
const fileOperation & fileHandler()
Get current file handler.
const volScalarField & alpha1
virtual tmp< volScalarField > Cp() const
Heat capacity at constant pressure [J/kg/K].
virtual tmp< volScalarField > CpByCpv() const
Heat capacity ratio [].
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
virtual tmp< volScalarField > hc() const
Chemical enthalpy [J/kg].
virtual void correctThermo()
Correct the thermodynamics of each phase.
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
const dimensionedScalar h
Planck constant.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
static autoPtr< rhoThermo > New(const fvMesh &, const word &phaseName=word::null)
Selector.
virtual bool incompressible() const
Return true if the equation of state is incompressible.
virtual tmp< volScalarField > Cpv() const
Heat capacity at constant pressure/volume [J/kg/K].
GeometricField< scalar, fvPatchField, volMesh > volScalarField
virtual tmp< volScalarField > kappa() const
Thermal diffusivity for temperature of mixture [J/m/s/K].
virtual word thermoName() const
Return the name of the thermo physics.
virtual void correct()
Update mixture properties.
GeometricField< vector, fvPatchField, volMesh > volVectorField
virtual tmp< scalarField > THE(const scalarField &h, const scalarField &p, const scalarField &T0, const labelList &cells) const
Temperature from enthalpy/internal energy for cell-set.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual tmp< volScalarField > W() const
Molecular weight [kg/kmol].
virtual ~twoPhaseMixtureThermo()
Destructor.
twoPhaseMixtureThermo(const volVectorField &U, const surfaceScalarField &phi)
Construct from components.
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
defineTypeNameAndDebug(combustionModel, 0)
bool read()
Read transportProperties dictionary.
virtual tmp< volScalarField > nu() const
Kinematic viscosity of mixture [m^2/s].