liquidPropertiesI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28inline Foam::scalar Foam::liquidProperties::limit(const scalar T) const
29{
30 return T;
31}
32
33
34inline Foam::scalar Foam::liquidProperties::Y() const
35{
36 return 1;
37}
38
39
40inline Foam::scalar Foam::liquidProperties::Tc() const
41{
42 return Tc_;
43}
44
45
46inline Foam::scalar Foam::liquidProperties::Pc() const
47{
48 return Pc_;
49}
50
51
52inline Foam::scalar Foam::liquidProperties::Vc() const
53{
54 return Vc_;
55}
56
57
58inline Foam::scalar Foam::liquidProperties::Zc() const
59{
60 return Zc_;
61}
62
63
64inline Foam::scalar Foam::liquidProperties::Tt() const
65{
66 return Tt_;
67}
68
69
70inline Foam::scalar Foam::liquidProperties::Pt() const
71{
72 return Pt_;
73}
74
75
76inline Foam::scalar Foam::liquidProperties::Tb() const
77{
78 return Tb_;
79}
80
81
82inline Foam::scalar Foam::liquidProperties::dipm() const
83{
84 return dipm_;
85}
86
87
88inline Foam::scalar Foam::liquidProperties::omega() const
89{
90 return omega_;
91}
92
93
94inline Foam::scalar Foam::liquidProperties::delta() const
95{
96 return delta_;
97}
98
99
100inline Foam::scalar Foam::liquidProperties::psi(scalar p, scalar T) const
101{
102 return 0;
103}
104
105
106inline Foam::scalar Foam::liquidProperties::CpMCv(scalar p, scalar T) const
107{
108 return 0;
109}
110
111
112inline Foam::scalar Foam::liquidProperties::Ha(scalar p, scalar T) const
113{
114 return h(p, T);
115}
116
117
118inline Foam::scalar Foam::liquidProperties::Hs(scalar p, scalar T) const
119{
120 return h(p, T);
121}
122
123
124inline Foam::scalar Foam::liquidProperties::Hc() const
125{
126 return 0;
127}
128
129
130inline Foam::scalar Foam::liquidProperties::alphah(scalar p, scalar T) const
131{
132 return kappa(p, T)/Cp(p, T);
133}
134
135
136template<class Func>
138(
139 Func& f,
140 const word& name,
141 const dictionary& dict
142)
143{
144 if (dict.found(name))
145 {
146 f = Func(dict.subDict(name));
147 }
148}
149
150
151template<class Liquid>
153(
154 Liquid& l,
155 const dictionary& dict
156)
157{
158 l.liquidProperties::readIfPresent(dict);
159 readIfPresent(l.rho_, "rho", dict);
160 readIfPresent(l.pv_, "pv", dict);
161 readIfPresent(l.hl_, "hl", dict);
162 readIfPresent(l.Cp_, "Cp", dict);
163 readIfPresent(l.h_, "h", dict);
164 readIfPresent(l.Cpg_, "Cpg", dict);
165 readIfPresent(l.B_, "B", dict);
166 readIfPresent(l.mu_, "mu", dict);
167 readIfPresent(l.mug_, "mug", dict);
168 readIfPresent(l.kappa_, "K", dict);
169 readIfPresent(l.kappag_, "kappag", dict);
170 readIfPresent(l.sigma_, "sigma", dict);
171 readIfPresent(l.D_, "D", dict);
172}
173
174
175template<class Liquid>
177(
178 const Liquid& l,
179 Ostream& os
180) const
181{
182 l.liquidProperties::writeData(os); os << nl;
183 l.rho_.writeData(os); os << nl;
184 l.pv_.writeData(os); os << nl;
185 l.hl_.writeData(os); os << nl;
186 l.Cp_.writeData(os); os << nl;
187 l.h_.writeData(os); os << nl;
188 l.Cpg_.writeData(os); os << nl;
189 l.B_.writeData(os); os << nl;
190 l.mu_.writeData(os); os << nl;
191 l.mug_.writeData(os); os << nl;
192 l.kappa_.writeData(os); os << nl;
193 l.kappag_.writeData(os); os << nl;
194 l.sigma_.writeData(os); os << nl;
195 l.D_.writeData(os); os << endl;
196}
197
198
199// ************************************************************************* //
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
friend complex limit(const complex &c1, const complex &c2)
Definition: complexI.H:263
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Definition: dictionary.C:460
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Definition: dictionaryI.H:87
scalar Vc() const
Critical volume [m^3/kmol].
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
scalar Zc() const
Critical compressibility factor.
scalar delta() const
Solubility parameter [(J/m^3)^(1/2)].
scalar Tb() const
Normal boiling temperature [K].
scalar Pc() const
Critical pressure [Pa].
scalar Hc() const
Chemical enthalpy [J/kg].
scalar Pt() const
Triple point pressure [Pa].
scalar omega() const
Pitzer's acentric factor [].
scalar Ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kg].
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
scalar Tt() const
Triple point temperature [K].
scalar alphah(const scalar p, const scalar T) const
Liquid thermal diffusivity of enthalpy [kg/ms].
scalar dipm() const
Dipole moment [].
scalar Y() const
No of moles of this species in mixture.
virtual void writeData(Ostream &os) const =0
Write the function coefficients.
scalar Tc() const
Critical temperature [K].
void readIfPresent(const dictionary &dict)
Read and set the properties present it the given dictionary.
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
const volScalarField & psi
const volScalarField & T
const volScalarField & Cp
Definition: EEqn.H:7
OBJstream os(runTime.globalPath()/outputName)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
labelList f(nPoints)
dictionary dict
volScalarField & h
propsDict readIfPresent("fields", acceptFields)