twoPhaseMixtureThermo.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "twoPhaseMixtureThermo.H"
29 
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(twoPhaseMixtureThermo, 0);
36 }
37 
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
42 (
43  const fvMesh& mesh
44 )
45 :
46  psiThermo(mesh, word::null),
47  twoPhaseMixture(mesh, *this),
48  thermo1_(NULL),
49  thermo2_(NULL)
50 {
51  {
52  volScalarField T1(IOobject::groupName("T", phase1Name()), T_);
53  T1.write();
54  }
55 
56  {
57  volScalarField T2(IOobject::groupName("T", phase2Name()), T_);
58  T2.write();
59  }
60 
61  thermo1_ = rhoThermo::New(mesh, phase1Name());
62  thermo2_ = rhoThermo::New(mesh, phase2Name());
63 
64  thermo1_->validate(phase1Name(), "e");
65  thermo2_->validate(phase2Name(), "e");
66 
67  correct();
68 }
69 
70 
71 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
72 
74 {}
75 
76 
77 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
78 
80 {
81  thermo1_->he() = thermo1_->he(p_, T_);
82  thermo1_->correct();
83 
84  thermo2_->he() = thermo2_->he(p_, T_);
85  thermo2_->correct();
86 
87  psi_ = alpha1()*thermo1_->psi() + alpha2()*thermo2_->psi();
88  mu_ = alpha1()*thermo1_->mu() + alpha2()*thermo2_->mu();
89  alpha_ = alpha1()*thermo1_->alpha() + alpha2()*thermo2_->alpha();
90 }
91 
92 
94 {
95  return thermo1_->incompressible() && thermo2_->incompressible();
96 }
97 
98 
100 {
101  return thermo1_->isochoric() && thermo2_->isochoric();
102 }
103 
104 
106 (
107  const volScalarField& p,
108  const volScalarField& T
109 ) const
110 {
111  return alpha1()*thermo1_->he(p, T) + alpha2()*thermo2_->he(p, T);
112 }
113 
114 
116 (
117  const scalarField& p,
118  const scalarField& T,
119  const labelList& cells
120 ) const
121 {
122  return
123  scalarField(alpha1(), cells)*thermo1_->he(p, T, cells)
124  + scalarField(alpha2(), cells)*thermo2_->he(p, T, cells);
125 }
126 
127 
129 (
130  const scalarField& p,
131  const scalarField& T,
132  const label patchi
133 ) const
134 {
135  return
136  alpha1().boundaryField()[patchi]*thermo1_->he(p, T, patchi)
137  + alpha2().boundaryField()[patchi]*thermo2_->he(p, T, patchi);
138 }
139 
140 
142 {
143  return alpha1()*thermo1_->hc() + alpha2()*thermo2_->hc();
144 }
145 
146 
148 (
149  const scalarField& h,
150  const scalarField& p,
151  const scalarField& T0,
152  const labelList& cells
153 ) const
154 {
156  return T0;
157 }
158 
159 
161 (
162  const scalarField& h,
163  const scalarField& p,
164  const scalarField& T0,
165  const label patchi
166 ) const
167 {
169  return T0;
170 }
171 
172 
174 {
175  return alpha1()*thermo1_->Cp() + alpha2()*thermo2_->Cp();
176 }
177 
178 
180 (
181  const scalarField& p,
182  const scalarField& T,
183  const label patchi
184 ) const
185 {
186  return
187  alpha1().boundaryField()[patchi]*thermo1_->Cp(p, T, patchi)
188  + alpha2().boundaryField()[patchi]*thermo2_->Cp(p, T, patchi);
189 }
190 
191 
193 {
194  return alpha1()*thermo1_->Cv() + alpha2()*thermo2_->Cv();
195 }
196 
197 
199 (
200  const scalarField& p,
201  const scalarField& T,
202  const label patchi
203 ) const
204 {
205  return
206  alpha1().boundaryField()[patchi]*thermo1_->Cv(p, T, patchi)
207  + alpha2().boundaryField()[patchi]*thermo2_->Cv(p, T, patchi);
208 }
209 
210 
212 {
213  return alpha1()*thermo1_->gamma() + alpha2()*thermo2_->gamma();
214 }
215 
216 
218 (
219  const scalarField& p,
220  const scalarField& T,
221  const label patchi
222 ) const
223 {
224  return
225  alpha1().boundaryField()[patchi]*thermo1_->gamma(p, T, patchi)
226  + alpha2().boundaryField()[patchi]*thermo2_->gamma(p, T, patchi);
227 }
228 
229 
231 {
232  return alpha1()*thermo1_->Cpv() + alpha2()*thermo2_->Cpv();
233 }
234 
235 
237 (
238  const scalarField& p,
239  const scalarField& T,
240  const label patchi
241 ) const
242 {
243  return
244  alpha1().boundaryField()[patchi]*thermo1_->Cpv(p, T, patchi)
245  + alpha2().boundaryField()[patchi]*thermo2_->Cpv(p, T, patchi);
246 }
247 
248 
250 {
251  return
252  alpha1()*thermo1_->CpByCpv()
253  + alpha2()*thermo2_->CpByCpv();
254 }
255 
256 
258 (
259  const scalarField& p,
260  const scalarField& T,
261  const label patchi
262 ) const
263 {
264  return
265  alpha1().boundaryField()[patchi]*thermo1_->CpByCpv(p, T, patchi)
266  + alpha2().boundaryField()[patchi]*thermo2_->CpByCpv(p, T, patchi);
267 }
268 
269 
271 {
272  return mu()/(alpha1()*thermo1_->rho() + alpha2()*thermo2_->rho());
273 }
274 
275 
277 (
278  const label patchi
279 ) const
280 {
281  return
282  mu(patchi)
283  /(
284  alpha1().boundaryField()[patchi]*thermo1_->rho(patchi)
285  + alpha2().boundaryField()[patchi]*thermo2_->rho(patchi)
286  );
287 }
288 
289 
291 {
292  return alpha1()*thermo1_->kappa() + alpha2()*thermo2_->kappa();
293 }
294 
295 
297 (
298  const label patchi
299 ) const
300 {
301  return
302  alpha1().boundaryField()[patchi]*thermo1_->kappa(patchi)
303  + alpha2().boundaryField()[patchi]*thermo2_->kappa(patchi);
304 }
305 
306 
308 (
309  const volScalarField& alphat
310 ) const
311 {
312  return
313  alpha1()*thermo1_->kappaEff(alphat)
314  + alpha2()*thermo2_->kappaEff(alphat);
315 }
316 
317 
319 (
320  const scalarField& alphat,
321  const label patchi
322 ) const
323 {
324  return
325  alpha1().boundaryField()[patchi]*thermo1_->kappaEff(alphat, patchi)
326  + alpha2().boundaryField()[patchi]*thermo2_->kappaEff(alphat, patchi)
327  ;
328 }
329 
330 
332 (
333  const volScalarField& alphat
334 ) const
335 {
336  return
337  alpha1()*thermo1_->alphaEff(alphat)
338  + alpha2()*thermo2_->alphaEff(alphat);
339 }
340 
341 
343 (
344  const scalarField& alphat,
345  const label patchi
346 ) const
347 {
348  return
349  alpha1().boundaryField()[patchi]*thermo1_->alphaEff(alphat, patchi)
350  + alpha2().boundaryField()[patchi]*thermo2_->alphaEff(alphat, patchi)
351  ;
352 }
353 
354 
355 // ************************************************************************* //
virtual bool isochoric() const
Return true if the equation of state is isochoric.
virtual bool incompressible() const
Return true if the equation of state is incompressible.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
virtual tmp< volScalarField > Cpv() const
Heat capacity at constant pressure/volume [J/kg/K].
virtual tmp< volScalarField > Cp() const
Heat capacity at constant pressure [J/kg/K].
virtual tmp< volScalarField > Cv() const
Heat capacity at constant volume [J/kg/K].
volScalarField T_
Temperature [K].
Definition: basicThermo.H:71
virtual tmp< volScalarField > CpByCpv() const
Heat capacity ratio [].
virtual tmp< volScalarField > kappaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
const volScalarField & alpha1() const
Return the phase-fraction of phase 1.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
virtual tmp< volScalarField > mu() const
Dynamic viscosity of mixture [kg/m/s].
Definition: psiThermo.C:111
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual tmp< volScalarField > kappa() const
Thermal diffusivity for temperature of mixture [J/m/s/K].
virtual tmp< volScalarField > gamma() const
Gamma = Cp/Cv [].
const Boundary & boundaryField() const
Return const-reference to the boundary field.
volScalarField & p_
Pressure [Pa].
Definition: basicThermo.H:68
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual ~twoPhaseMixtureThermo()
Destructor.
List< label > labelList
A List of labels.
Definition: labelList.H:56
virtual tmp< volScalarField > hc() const
Chemical enthalpy [J/kg].
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.
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil)-psi)*pSat, rhoMin);#1"/home/ubuntu/OpenFOAM-4.1/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()
Definition: pEqn.H:72
defineTypeNameAndDebug(combustionModel, 0)
virtual volScalarField & he()
Enthalpy/Internal energy [J/kg].
const dimensionedScalar h
Planck constant.
Definition: createFields.H:6
volScalarField alpha_
Laminar thermal diffusuvity [kg/m/s].
Definition: basicThermo.H:74
label patchi
const volScalarField & alpha2() const
Return the phase-fraction of phase 2.
virtual void correct()
Update properties.
volScalarField mu_
Dynamic viscosity [kg/m/s].
Definition: psiThermo.H:62
virtual tmp< volScalarField > alphaEff(const volScalarField &alphat) const
Effective thermal diffusivity of mixture [J/m/s/K].
volScalarField psi_
Compressibility [s^2/m^2].
Definition: psiThermo.H:59
A class for managing temporary objects.
Definition: PtrList.H:54
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:366
twoPhaseMixtureThermo(const fvMesh &mesh)
Construct from mesh.
Namespace for OpenFOAM.
virtual tmp< volScalarField > nu() const
Kinematic viscosity of mixture [m^2/s].