liquidFilmThermo.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 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 "liquidFilmThermo.H"
27 #include "demandDrivenData.H"
28 #include "thermoSingleLayer.H"
29 #include "SLGThermo.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace regionModels
37 {
38 namespace surfaceFilmModels
39 {
40 
41 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42 
43 defineTypeNameAndDebug(liquidFilmThermo, 0);
44 
46 (
47  filmThermoModel,
48  liquidFilmThermo,
49  dictionary
50 );
51 
52 
53 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
54 
56 {
57  if (!isA<thermoSingleLayer>(owner_))
58  {
60  (
61  "const thermoSingleLayer& liquidFilmThermo::thermoFilm() const"
62  )
63  << "Thermo model requires a " << thermoSingleLayer::typeName
64  << " film to supply the pressure and temperature, but "
65  << owner_.type() << " film model selected. "
66  << "Use the 'useReferenceValues' flag to employ reference "
67  << "pressure and temperature" << exit(FatalError);
68  }
69 
70  return refCast<const thermoSingleLayer>(owner_);
71 }
72 
73 
75 {
76  if (liquidPtr_ != NULL)
77  {
78  return;
79  }
80 
81  dict.lookup("liquid") >> name_;
82 
83  if (owner_.primaryMesh().foundObject<SLGThermo>("SLGThermo"))
84  {
85  // retrieve from film thermo
86  ownLiquid_ = false;
87 
88  const SLGThermo& thermo =
89  owner_.primaryMesh().lookupObject<SLGThermo>("SLGThermo");
90  label id = thermo.liquidId(name_);
91  liquidPtr_ = &thermo.liquids().properties()[id];
92  }
93  else
94  {
95  // new liquid create
96  ownLiquid_ = true;
97 
98  liquidPtr_ = new liquidProperties(dict.subDict(name_ + "Coeffs"));
99  }
100 }
101 
102 
103 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
104 
106 (
107  surfaceFilmModel& owner,
108  const dictionary& dict
109 )
110 :
111  filmThermoModel(typeName, owner, dict),
112  name_("unknown_liquid"),
113  liquidPtr_(NULL),
114  ownLiquid_(false),
115  useReferenceValues_(readBool(coeffDict_.lookup("useReferenceValues"))),
116  pRef_(0.0),
117  TRef_(0.0)
118 {
120 
122  {
123  coeffDict_.lookup("pRef") >> pRef_;
124  coeffDict_.lookup("TRef") >> TRef_;
125  }
126 }
127 
128 
129 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
130 
132 {
133  if (ownLiquid_)
134  {
136  }
137 }
138 
139 
140 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
141 
143 {
144  return name_;
145 }
146 
147 
149 (
150  const scalar p,
151  const scalar T
152 ) const
153 {
154  return liquidPtr_->rho(p, T);
155 }
156 
157 
159 (
160  const scalar p,
161  const scalar T
162 ) const
163 {
164  return liquidPtr_->mu(p, T);
165 }
166 
167 
169 (
170  const scalar p,
171  const scalar T
172 ) const
173 {
174  return liquidPtr_->sigma(p, T);
175 }
176 
177 
179 (
180  const scalar p,
181  const scalar T
182 ) const
183 {
184  return liquidPtr_->Cp(p, T);
185 }
186 
187 
189 (
190  const scalar p,
191  const scalar T
192 ) const
193 {
194  return liquidPtr_->K(p, T);
195 }
196 
197 
198 scalar liquidFilmThermo::D
199 (
200  const scalar p,
201  const scalar T
202 ) const
203 {
204  return liquidPtr_->D(p, T);
205 }
206 
207 
209 (
210  const scalar p,
211  const scalar T
212 ) const
213 {
214  return liquidPtr_->hl(p, T);
215 }
216 
217 
219 (
220  const scalar p,
221  const scalar T
222 ) const
223 {
224  return liquidPtr_->pv(p, T);
225 }
226 
227 
228 scalar liquidFilmThermo::W() const
229 {
230  return liquidPtr_->W();
231 }
232 
233 
234 scalar liquidFilmThermo::Tb(const scalar p) const
235 {
236  return liquidPtr_->pvInvert(p);
237 }
238 
239 
241 {
243  (
244  new volScalarField
245  (
246  IOobject
247  (
248  type() + ":rho",
249  owner().time().timeName(),
250  owner().regionMesh(),
253  ),
254  owner().regionMesh(),
255  dimensionedScalar("0", dimDensity, 0.0),
256  zeroGradientFvPatchScalarField::typeName
257  )
258  );
259 
260  scalarField& rho = trho().internalField();
261 
263  {
264  forAll(rho, cellI)
265  {
266  rho[cellI] = this->rho(pRef_, TRef_);
267  }
268  }
269  else
270  {
271  const thermoSingleLayer& film = thermoFilm();
272 
273  const volScalarField& T = film.T();
274  const volScalarField& p = film.pPrimary();
275 
276  forAll(rho, cellI)
277  {
278  rho[cellI] = this->rho(p[cellI], T[cellI]);
279  }
280  }
281 
282  trho().correctBoundaryConditions();
283 
284  return trho;
285 }
286 
287 
289 {
291  (
292  new volScalarField
293  (
294  IOobject
295  (
296  type() + ":mu",
297  owner().time().timeName(),
298  owner().regionMesh(),
301  ),
302  owner().regionMesh(),
304  zeroGradientFvPatchScalarField::typeName
305  )
306  );
307 
308  scalarField& mu = tmu().internalField();
309 
311  {
312  forAll(mu, cellI)
313  {
314  mu[cellI] = this->mu(pRef_, TRef_);
315  }
316  }
317  else
318  {
319  const thermoSingleLayer& film = thermoFilm();
320 
321  const volScalarField& T = film.T();
322  const volScalarField& p = film.pPrimary();
323 
324  forAll(mu, cellI)
325  {
326  mu[cellI] = this->mu(p[cellI], T[cellI]);
327  }
328  }
329 
330  tmu().correctBoundaryConditions();
331 
332  return tmu;
333 }
334 
335 
337 {
338  tmp<volScalarField> tsigma
339  (
340  new volScalarField
341  (
342  IOobject
343  (
344  type() + ":sigma",
345  owner().time().timeName(),
346  owner().regionMesh(),
349  ),
350  owner().regionMesh(),
351  dimensionedScalar("0", dimMass/sqr(dimTime), 0.0),
352  zeroGradientFvPatchScalarField::typeName
353  )
354  );
355 
356  scalarField& sigma = tsigma().internalField();
357 
359  {
360  forAll(sigma, cellI)
361  {
362  sigma[cellI] = this->sigma(pRef_, TRef_);
363  }
364  }
365  else
366  {
367  const thermoSingleLayer& film = thermoFilm();
368 
369  const volScalarField& T = film.T();
370  const volScalarField& p = film.pPrimary();
371 
372  forAll(sigma, cellI)
373  {
374  sigma[cellI] = this->sigma(p[cellI], T[cellI]);
375  }
376  }
377 
378  tsigma().correctBoundaryConditions();
379 
380  return tsigma;
381 }
382 
383 
385 {
387  (
388  new volScalarField
389  (
390  IOobject
391  (
392  type() + ":Cp",
393  owner().time().timeName(),
394  owner().regionMesh(),
397  ),
398  owner().regionMesh(),
400  zeroGradientFvPatchScalarField::typeName
401  )
402  );
403 
404  scalarField& Cp = tCp().internalField();
405 
407  {
408  forAll(Cp, cellI)
409  {
410  Cp[cellI] = this->Cp(pRef_, TRef_);
411  }
412  }
413  else
414  {
415  const thermoSingleLayer& film = thermoFilm();
416 
417  const volScalarField& T = film.T();
418  const volScalarField& p = film.pPrimary();
419 
420  forAll(Cp, cellI)
421  {
422  Cp[cellI] = this->Cp(p[cellI], T[cellI]);
423  }
424  }
425 
426  tCp().correctBoundaryConditions();
427 
428  return tCp;
429 }
430 
431 
433 {
434  tmp<volScalarField> tkappa
435  (
436  new volScalarField
437  (
438  IOobject
439  (
440  type() + ":kappa",
441  owner().time().timeName(),
442  owner().regionMesh(),
445  ),
446  owner().regionMesh(),
448  zeroGradientFvPatchScalarField::typeName
449  )
450  );
451 
452  scalarField& kappa = tkappa().internalField();
453 
455  {
456  forAll(kappa, cellI)
457  {
458  kappa[cellI] = this->kappa(pRef_, TRef_);
459  }
460  }
461  else
462  {
463  const thermoSingleLayer& film = thermoFilm();
464 
465  const volScalarField& T = film.T();
466  const volScalarField& p = film.pPrimary();
467 
468  forAll(kappa, cellI)
469  {
470  kappa[cellI] = this->kappa(p[cellI], T[cellI]);
471  }
472  }
473 
474  tkappa().correctBoundaryConditions();
475 
476  return tkappa;
477 }
478 
479 
480 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
481 
482 } // End namespace surfaceFilmModels
483 } // End namespace regionModels
484 } // End namespace Foam
485 
486 // ************************************************************************* //
const dimensionSet dimPressure
surfaceFilmModel & owner_
Reference to the owner surface film model.
The thermophysical properties of a liquidProperties.
virtual scalar pv(const scalar p, const scalar T) const
Return vapour pressure [Pa].
bool readBool(Istream &)
Definition: boolIO.C:60
virtual scalar W() const
Return molecular weight [kg/kmol].
virtual tmp< volScalarField > Cp() const
Return specific heat capacity [J/kg/K].
virtual tmp< volScalarField > sigma() const
Return surface tension [kg/s2].
virtual scalar Cp(scalar p, scalar T) const
Liquid heat capacity [J/(kg K)].
liquidFilmThermo(const liquidFilmThermo &)
Disallow default bitwise copy construct.
defineTypeNameAndDebug(kinematicSingleLayer, 0)
virtual const volScalarField & T() const
Return the film mean temperature [K].
bool foundObject(const word &name) const
Is the named Type found?
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
void deleteDemandDrivenData(DataPtr &dataPtr)
virtual tmp< volScalarField > mu() const
Return dynamic viscosity [Pa.s].
A class for handling words, derived from string.
Definition: word.H:59
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
const dictionary coeffDict_
Coefficients dictionary.
Definition: subModelBase.H:79
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
addToRunTimeSelectionTable(surfaceFilmModel, kinematicSingleLayer, mesh)
virtual tmp< volScalarField > kappa() const
Return thermal conductivity [W/m/K].
const dimensionSet dimEnergy
scalar W() const
Molecular weight [kg/kmol].
const dimensionSet dimDensity
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const dimensionSet dimPower
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
virtual scalar D(scalar p, scalar T) const
Vapour diffussivity [m2/s].
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
Namespace for OpenFOAM.
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:62
const thermoSingleLayer & thermoFilm() const
Return a reference to a thermo film.
bool useReferenceValues_
Flag to indicate that reference values of p and T should be used.
bool ownLiquid_
Flag to indicate that model owns the liquid object.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual scalar K(scalar p, scalar T) const
Liquid thermal conductivity [W/(m K)].
const surfaceFilmModel & owner() const
Return const access to the owner surface film model.
#define forAll(list, i)
Definition: UList.H:421
void initLiquid(const dictionary &dict)
Initialise the liquid pointer.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
Macros for easy insertion into run-time selection tables.
Template functions to aid in the implementation of demand driven data.
const liquidProperties * liquidPtr_
Pointer to the liquid properties.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
virtual scalar Tb(const scalar p) const
Return boiling temperature [K].
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual scalar pv(scalar p, scalar T) const
Vapour pressure [Pa].
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
Definition: regionModelI.H:31
const PtrList< liquidProperties > & properties() const
Return the liquid properties.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
virtual scalar D(const scalar p, const scalar T) const
Return diffusivity [m2/s].
error FatalError
virtual scalar pvInvert(scalar p) const
Invert the vapour pressure relationship to retrieve the.
virtual scalar rho(scalar p, scalar T) const
Liquid rho [kg/m^3].
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
virtual tmp< volScalarField > rho() const
Return density [kg/m3].
virtual scalar mu(scalar p, scalar T) const
Liquid viscosity [Pa s].
const volScalarField & pPrimary() const
Pressure / [Pa].
virtual scalar hl(const scalar p, const scalar T) const
Return latent heat [J/kg].
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
virtual const word & name() const
Return the specie name.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
virtual scalar hl(scalar p, scalar T) const
Heat of vapourisation [J/kg].
virtual scalar sigma(scalar p, scalar T) const
Surface tension [N/m].
const liquidMixtureProperties & liquids() const
Return reference to the global (additional) liquids.
Definition: SLGThermo.C:124
word timeName
Definition: getTimeIndex.H:3
tmp< volScalarField > trho
label liquidId(const word &cmptName, bool allowNotFound=false) const
Index of liquid component.
Definition: SLGThermo.C:185