constantFilmThermo.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-2017 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 "constantFilmThermo.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace regionModels
35 {
36 namespace surfaceFilmModels
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 defineTypeNameAndDebug(constantFilmThermo, 0);
42 
44 (
45  filmThermoModel,
46  constantFilmThermo,
47  dictionary
48 );
49 
50 
51 void constantFilmThermo::init(thermoData& td)
52 {
53  if (coeffDict_.readIfPresent(td.name_, td.value_))
54  {
55  td.set_ = true;
56  }
57 }
58 
59 
60 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
61 
62 constantFilmThermo::constantFilmThermo
63 (
64  surfaceFilmModel& film,
65  const dictionary& dict
66 )
67 :
68  filmThermoModel(typeName, film, dict),
69  name_(coeffDict_.lookup("specie")),
70  rho0_("rho0"),
71  mu0_("mu0"),
72  sigma0_("sigma0"),
73  Cp0_("Cp0"),
74  kappa0_("kappa0"),
75  hl0_("hl0"),
76  pv0_("pv0"),
77  W0_("W0"),
78  Tb0_("Tb0")
79 {
80  init(rho0_);
81  init(mu0_);
82  init(sigma0_);
83  init(Cp0_);
84  init(kappa0_);
85  init(hl0_);
86  init(pv0_);
87  init(W0_);
88  init(Tb0_);
89 }
90 
91 
92 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
93 
95 {}
96 
97 
98 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
99 
101 {
102  return name_;
103 }
104 
105 
107 (
108  const scalar p,
109  const scalar T
110 ) const
111 {
112  if (!rho0_.set_)
113  {
114  coeffDict_.lookup(rho0_.name_) >> rho0_.value_;
115  rho0_.set_ = true;
116  }
117 
118  return rho0_.value_;
119 }
120 
121 
123 (
124  const scalar p,
125  const scalar T
126 ) const
127 {
128  if (!mu0_.set_)
129  {
130  coeffDict_.lookup(mu0_.name_) >> mu0_.value_;
131  mu0_.set_ = true;
132  }
133 
134  return mu0_.value_;
135 }
136 
137 
139 (
140  const scalar p,
141  const scalar T
142 ) const
143 {
144  if (!sigma0_.set_)
145  {
146  coeffDict_.lookup(sigma0_.name_) >> sigma0_.value_;
147  sigma0_.set_ = true;
148  }
149 
150  return sigma0_.value_;
151 }
152 
153 
155 (
156  const scalar p,
157  const scalar T
158 ) const
159 {
160  if (!Cp0_.set_)
161  {
162  coeffDict_.lookup(Cp0_.name_) >> Cp0_.value_;
163  Cp0_.set_ = true;
164  }
165 
166  return Cp0_.value_;
167 }
168 
169 
171 (
172  const scalar p,
173  const scalar T
174 ) const
175 {
176  if (!kappa0_.set_)
177  {
178  coeffDict_.lookup(kappa0_.name_) >> kappa0_.value_;
179  kappa0_.set_ = true;
180  }
181 
182  return kappa0_.value_;
183 }
184 
185 
187 (
188  const scalar p,
189  const scalar T
190 ) const
191 {
192  if (!D0_.set_)
193  {
194  coeffDict_.lookup(D0_.name_) >> D0_.value_;
195  D0_.set_ = true;
196  }
197 
198  return D0_.value_;
199 }
200 
201 
203 (
204  const scalar p,
205  const scalar T
206 ) const
207 {
208  if (!hl0_.set_)
209  {
210  coeffDict_.lookup(hl0_.name_) >> hl0_.value_;
211  hl0_.set_ = true;
212  }
213 
214  return hl0_.value_;
215 }
216 
217 
219 (
220  const scalar p,
221  const scalar T
222 ) const
223 {
224  if (!pv0_.set_)
225  {
226  coeffDict_.lookup(pv0_.name_) >> pv0_.value_;
227  pv0_.set_ = true;
228  }
229 
230  return pv0_.value_;
231 }
232 
233 
234 scalar constantFilmThermo::W() const
235 {
236  if (!W0_.set_)
237  {
238  coeffDict_.lookup(W0_.name_) >> W0_.value_;
239  W0_.set_ = true;
240  }
241 
242  return W0_.value_;
243 }
244 
245 
246 scalar constantFilmThermo::Tb(const scalar p) const
247 {
248  if (!Tb0_.set_)
249  {
250  coeffDict_.lookup(Tb0_.name_) >> Tb0_.value_;
251  Tb0_.set_ = true;
252  }
253 
254  return Tb0_.value_;
255 }
256 
257 
259 {
261  (
262  new volScalarField
263  (
264  IOobject
265  (
266  type() + ':' + rho0_.name_,
267  film().time().timeName(),
268  film().regionMesh(),
271  ),
272  film().regionMesh(),
273  dimensionedScalar("0", dimDensity, 0.0),
274  extrapolatedCalculatedFvPatchScalarField::typeName
275  )
276  );
277 
278  trho.ref().primitiveFieldRef() = this->rho(0, 0);
279  trho.ref().correctBoundaryConditions();
280 
281  return trho;
282 }
283 
284 
286 {
288  (
289  new volScalarField
290  (
291  IOobject
292  (
293  type() + ':' + mu0_.name_,
294  film().time().timeName(),
295  film().regionMesh(),
298  ),
299  film().regionMesh(),
301  extrapolatedCalculatedFvPatchScalarField::typeName
302  )
303  );
304 
305  tmu.ref().primitiveFieldRef() = this->mu(0, 0);
306  tmu.ref().correctBoundaryConditions();
307 
308  return tmu;
309 }
310 
311 
313 {
314  tmp<volScalarField> tsigma
315  (
316  new volScalarField
317  (
318  IOobject
319  (
320  type() + ':' + sigma0_.name_,
321  film().time().timeName(),
322  film().regionMesh(),
325  ),
326  film().regionMesh(),
327  dimensionedScalar("0", dimMass/sqr(dimTime), 0.0),
328  extrapolatedCalculatedFvPatchScalarField::typeName
329  )
330  );
331 
332  tsigma.ref().primitiveFieldRef() = this->sigma(0, 0);
333  tsigma.ref().correctBoundaryConditions();
334 
335  return tsigma;
336 }
337 
338 
340 {
342  (
343  new volScalarField
344  (
345  IOobject
346  (
347  type() + ':' + Cp0_.name_,
348  film().time().timeName(),
349  film().regionMesh(),
352  ),
353  film().regionMesh(),
355  extrapolatedCalculatedFvPatchScalarField::typeName
356  )
357  );
358 
359  tCp.ref().primitiveFieldRef() = this->Cp(0, 0);
360  tCp.ref().correctBoundaryConditions();
361 
362  return tCp;
363 }
364 
365 
367 {
368  tmp<volScalarField> tkappa
369  (
370  new volScalarField
371  (
372  IOobject
373  (
374  type() + ':' + kappa0_.name_,
375  film().time().timeName(),
376  film().regionMesh(),
379  ),
380  film().regionMesh(),
382  extrapolatedCalculatedFvPatchScalarField::typeName
383  )
384  );
385 
386  tkappa.ref().primitiveFieldRef() = this->kappa(0, 0);
387  tkappa.ref().correctBoundaryConditions();
388 
389  return tkappa;
390 }
391 
392 
393 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
394 
395 } // End namespace surfaceFilmModels
396 } // End namespace regionModels
397 } // End namespace Foam
398 
399 // ************************************************************************* //
const surfaceFilmModel & film() const
Return const access to the film surface film model.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< volScalarField > trho
virtual const word & name() const
Return the specie name.
virtual scalar Tb(const scalar p) const
Return boiling temperature [K].
const dictionary coeffDict_
Coefficients dictionary.
Definition: subModelBase.H:79
Macros for easy insertion into run-time selection tables.
virtual scalar D(const scalar p, const scalar T) const
Return diffusivity [m2/s].
addToRunTimeSelectionTable(surfaceFilmModel, kinematicSingleLayer, mesh)
virtual tmp< volScalarField > Cp() const
Return specific heat capacity [J/kg/K].
A class for handling words, derived from string.
Definition: word.H:59
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
virtual scalar W() const
Return molecular weight [kg/kmol].
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
word timeName
Definition: getTimeIndex.H:3
const dimensionSet dimPressure
const dimensionSet dimPower
virtual tmp< volScalarField > sigma() const
Return surface tension [kg/s2].
virtual tmp< volScalarField > mu() const
Return dynamic viscosity [Pa.s].
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
const dimensionSet dimEnergy
virtual tmp< volScalarField > rho() const
Return density [kg/m3].
const dimensionSet dimDensity
virtual scalar hl(const scalar p, const scalar T) const
Return latent heat [J/kg].
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual tmp< volScalarField > kappa() const
Return thermal conductivity [W/m/K].
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
defineTypeNameAndDebug(kinematicSingleLayer, 0)
Namespace for OpenFOAM.
virtual scalar pv(const scalar p, const scalar T) const
Return vapour pressure [Pa].
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576