liquidPropertiesI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 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 "thermophysicalFunction.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 inline Foam::scalar Foam::liquidProperties::limit(const scalar T) const
31 {
32  return T;
33 }
34 
35 
36 inline Foam::scalar Foam::liquidProperties::Y() const
37 {
38  return 1;
39 }
40 
41 
42 inline Foam::scalar Foam::liquidProperties::Tc() const
43 {
44  return Tc_;
45 }
46 
47 
48 inline Foam::scalar Foam::liquidProperties::Pc() const
49 {
50  return Pc_;
51 }
52 
53 
54 inline Foam::scalar Foam::liquidProperties::Vc() const
55 {
56  return Vc_;
57 }
58 
59 
60 inline Foam::scalar Foam::liquidProperties::Zc() const
61 {
62  return Zc_;
63 }
64 
65 
66 inline Foam::scalar Foam::liquidProperties::Tt() const
67 {
68  return Tt_;
69 }
70 
71 
72 inline Foam::scalar Foam::liquidProperties::Pt() const
73 {
74  return Pt_;
75 }
76 
77 
78 inline Foam::scalar Foam::liquidProperties::Tb() const
79 {
80  return Tb_;
81 }
82 
83 
84 inline Foam::scalar Foam::liquidProperties::dipm() const
85 {
86  return dipm_;
87 }
88 
89 
90 inline Foam::scalar Foam::liquidProperties::omega() const
91 {
92  return omega_;
93 }
94 
95 
96 inline Foam::scalar Foam::liquidProperties::delta() const
97 {
98  return delta_;
99 }
100 
101 
102 inline Foam::scalar Foam::liquidProperties::psi(scalar p, scalar T) const
103 {
104  return 0;
105 }
106 
107 
108 inline Foam::scalar Foam::liquidProperties::CpMCv(scalar p, scalar T) const
109 {
110  return 0;
111 }
112 
113 
114 inline Foam::scalar Foam::liquidProperties::Ha(scalar p, scalar T) const
115 {
116  return h(p, T);
117 }
118 
119 
120 inline Foam::scalar Foam::liquidProperties::Hs(scalar p, scalar T) const
121 {
122  return h(p, T);
123 }
124 
125 
126 inline Foam::scalar Foam::liquidProperties::Hf() const
127 {
128  return 0;
129 }
130 
131 
132 inline Foam::scalar Foam::liquidProperties::alphah(scalar p, scalar T) const
133 {
134  return kappa(p, T)/Cp(p, T);
135 }
136 
137 
138 template<class Func>
140 (
141  Func& f,
142  const word& name,
143  const dictionary& dict
144 )
145 {
146  if (dict.found(name))
147  {
148  f = Func(dict.subDict(name));
149  }
150 }
151 
152 
153 template<class Liquid>
155 (
156  Liquid& l,
157  const dictionary& dict
158 )
159 {
160  l.liquidProperties::readIfPresent(dict);
161  readIfPresent(l.rho_, "rho", dict);
162  readIfPresent(l.pv_, "pv", dict);
163  readIfPresent(l.hl_, "hl", dict);
164  readIfPresent(l.Cp_, "Cp", dict);
165  readIfPresent(l.h_, "h", dict);
166  readIfPresent(l.Cpg_, "Cpg", dict);
167  readIfPresent(l.B_, "B", dict);
168  readIfPresent(l.mu_, "mu", dict);
169  readIfPresent(l.mug_, "mug", dict);
170  readIfPresent(l.kappa_, "kappa", dict);
171  readIfPresent(l.kappag_, "kappag", dict);
172  readIfPresent(l.sigma_, "sigma", dict);
173  readIfPresent(l.D_, "D", dict);
174 }
175 
176 
177 template<class Liquid>
179 (
180  const Liquid& l,
181  Ostream& os
182 ) const
183 {
184  l.liquidProperties::write(os);
185  l.rho_.write(os, "rho");
186  l.pv_.write(os, "pv");
187  l.hl_.write(os, "hl");
188  l.Cp_.write(os, "Cp");
189  l.h_.write(os, "h");
190  l.Cpg_.write(os, "Cpg");
191  l.B_.write(os, "B");
192  l.mu_.write(os, "mu");
193  l.mug_.write(os, "mug");
194  l.kappa_.write(os, "kappa");
195  l.kappag_.write(os, "kappag");
196  l.sigma_.write(os, "sigma");
197  l.D_.write(os, "D");
198 }
199 
200 
201 // ************************************************************************* //
scalar delta() const
Solubility parameter [(J/m^3)^(1/2)].
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:667
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
void readIfPresent(const dictionary &dict)
Read and set the properties present it the given dictionary.
scalar Y() const
Mass fraction of this specie in mixture.
scalar Ha(const scalar p, const scalar T) const
Absolute enthalpy [J/kg].
scalar Pt() const
Triple point pressure [Pa].
virtual scalar kappa(scalar p, scalar T) const =0
Liquid thermal conductivity [W/m/K].
scalar dipm() const
Dipole moment [].
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:934
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
scalar Tc() const
Critical temperature [K].
scalar Zc() const
Critical compressibility factor.
scalar Hf() const
Enthalpy of formation [J/kg].
scalar psi(scalar p, scalar T) const
Liquid compressibility [s^2/m^2].
A class for handling words, derived from string.
Definition: word.H:59
virtual scalar h(scalar p, scalar T) const =0
Liquid enthalpy [J/kg] - reference to 298.15 K.
virtual scalar Cp(const scalar p, const scalar T) const =0
Heat capacity at constant pressure [J/kg/K].
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
scalar Vc() const
Critical volume [m^3/kmol].
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volScalarField & T
scalar omega() const
Pitzer&#39;s acentric factor [].
scalar Pc() const
Critical pressure [Pa].
volScalarField & p
scalar alphah(const scalar p, const scalar T) const
Liquid thermal diffusivity of enthalpy [kg/m/s].
virtual void write(Ostream &os) const =0
Write the function coefficients.
scalar Tb() const
Normal boiling temperature [K].
scalar Tt() const
Triple point temperature [K].