26 #include "twoPhaseMixtureThermo.H" 27 #include "gradientEnergyFvPatchScalarField.H" 28 #include "mixedEnergyFvPatchScalarField.H" 35 defineTypeNameAndDebug(twoPhaseMixtureThermo, 0);
41 Foam::twoPhaseMixtureThermo::twoPhaseMixtureThermo
46 psiThermo(mesh, word::null),
47 twoPhaseMixture(mesh, *this),
52 volScalarField T1(IOobject::groupName(
"T", phase1Name()), T_);
57 volScalarField T2(IOobject::groupName(
"T", phase2Name()), T_);
61 thermo1_ = rhoThermo::New(mesh, phase1Name());
62 thermo2_ = rhoThermo::New(mesh, phase2Name());
64 thermo1_->validate(phase1Name(),
"e");
65 thermo2_->validate(phase2Name(),
"e");
73 Foam::twoPhaseMixtureThermo::~twoPhaseMixtureThermo()
79 void Foam::twoPhaseMixtureThermo::correct()
81 thermo1_->he() = thermo1_->he(p_, T_);
84 thermo2_->he() = thermo2_->he(p_, T_);
87 psi_ = alpha1()*thermo1_->psi() + alpha2()*thermo2_->psi();
88 mu_ = alpha1()*thermo1_->mu() + alpha2()*thermo2_->mu();
89 alpha_ = alpha1()*thermo1_->alpha() + alpha2()*thermo2_->alpha();
93 bool Foam::twoPhaseMixtureThermo::incompressible()
const 95 return thermo1_->incompressible() && thermo2_->incompressible();
99 bool Foam::twoPhaseMixtureThermo::isochoric()
const 101 return thermo1_->isochoric() && thermo2_->isochoric();
105 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::he
107 const volScalarField& p,
108 const volScalarField& T
111 return alpha1()*thermo1_->he(p, T) + alpha2()*thermo2_->he(p, T);
115 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::he
117 const scalarField& p,
118 const scalarField& T,
119 const labelList& cells
123 scalarField(alpha1(), cells)*thermo1_->he(p, T, cells)
124 + scalarField(alpha2(), cells)*thermo2_->he(p, T, cells);
128 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::he
130 const scalarField& p,
131 const scalarField& T,
136 alpha1().boundaryField()[patchi]*thermo1_->he(p, T, patchi)
137 + alpha2().boundaryField()[patchi]*thermo2_->he(p, T, patchi);
141 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::hc()
const 143 return alpha1()*thermo1_->hc() + alpha2()*thermo2_->hc();
147 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::THE
149 const scalarField& h,
150 const scalarField& p,
151 const scalarField& T0,
152 const labelList& cells
155 notImplemented(
"twoPhaseMixtureThermo::THE(...)");
160 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::THE
162 const scalarField& h,
163 const scalarField& p,
164 const scalarField& T0,
168 notImplemented(
"twoPhaseMixtureThermo::THE(...)");
173 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::Cp()
const 175 return alpha1()*thermo1_->Cp() + alpha2()*thermo2_->Cp();
179 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::Cp
181 const scalarField& p,
182 const scalarField& T,
187 alpha1().boundaryField()[patchi]*thermo1_->Cp(p, T, patchi)
188 + alpha2().boundaryField()[patchi]*thermo2_->Cp(p, T, patchi);
192 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::Cv()
const 194 return alpha1()*thermo1_->Cv() + alpha2()*thermo2_->Cv();
198 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::Cv
200 const scalarField& p,
201 const scalarField& T,
206 alpha1().boundaryField()[patchi]*thermo1_->Cv(p, T, patchi)
207 + alpha2().boundaryField()[patchi]*thermo2_->Cv(p, T, patchi);
211 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::gamma()
const 213 return alpha1()*thermo1_->gamma() + alpha2()*thermo2_->gamma();
217 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::gamma
219 const scalarField& p,
220 const scalarField& T,
225 alpha1().boundaryField()[patchi]*thermo1_->gamma(p, T, patchi)
226 + alpha2().boundaryField()[patchi]*thermo2_->gamma(p, T, patchi);
230 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::Cpv()
const 232 return alpha1()*thermo1_->Cpv() + alpha2()*thermo2_->Cpv();
236 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::Cpv
238 const scalarField& p,
239 const scalarField& T,
244 alpha1().boundaryField()[patchi]*thermo1_->Cpv(p, T, patchi)
245 + alpha2().boundaryField()[patchi]*thermo2_->Cpv(p, T, patchi);
249 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::CpByCpv()
const 252 alpha1()*thermo1_->CpByCpv()
253 + alpha2()*thermo2_->CpByCpv();
257 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::CpByCpv
259 const scalarField& p,
260 const scalarField& T,
265 alpha1().boundaryField()[patchi]*thermo1_->CpByCpv(p, T, patchi)
266 + alpha2().boundaryField()[patchi]*thermo2_->CpByCpv(p, T, patchi);
270 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::kappa()
const 272 return alpha1()*thermo1_->kappa() + alpha2()*thermo2_->kappa();
276 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::kappa
282 alpha1().boundaryField()[patchi]*thermo1_->kappa(patchi)
283 + alpha2().boundaryField()[patchi]*thermo2_->kappa(patchi);
287 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::kappaEff
289 const volScalarField& alphat
293 alpha1()*thermo1_->kappaEff(alphat)
294 + alpha2()*thermo2_->kappaEff(alphat);
298 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::kappaEff
300 const scalarField& alphat,
305 alpha1().boundaryField()[patchi]*thermo1_->kappaEff(alphat, patchi)
306 + alpha2().boundaryField()[patchi]*thermo2_->kappaEff(alphat, patchi)
311 Foam::tmp<Foam::volScalarField> Foam::twoPhaseMixtureThermo::alphaEff
313 const volScalarField& alphat
317 alpha1()*thermo1_->alphaEff(alphat)
318 + alpha2()*thermo2_->alphaEff(alphat);
322 Foam::tmp<Foam::scalarField> Foam::twoPhaseMixtureThermo::alphaEff
324 const scalarField& alphat,
329 alpha1().boundaryField()[patchi]*thermo1_->alphaEff(alphat, patchi)
330 + alpha2().boundaryField()[patchi]*thermo2_->alphaEff(alphat, patchi)