perfectFluidI.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) 2012-2023 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 "perfectFluid.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Specie>
32 (
33  const Specie& sp,
34  const scalar R,
35  const scalar rho0
36 )
37 :
38  Specie(sp),
39  R_(R),
40  rho0_(rho0)
41 {}
42 
43 
44 template<class Specie>
46 (
47  const word& name,
48  const perfectFluid<Specie>& pf
49 )
50 :
51  Specie(name, pf),
52  R_(pf.R_),
53  rho0_(pf.rho0_)
54 {}
55 
56 
57 template<class Specie>
60 {
62 }
63 
64 
65 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66 
67 template<class Specie>
68 inline Foam::scalar Foam::perfectFluid<Specie>::R() const
69 {
70  return R_;
71 }
72 
73 
74 template<class Specie>
75 inline Foam::scalar Foam::perfectFluid<Specie>::rho(scalar p, scalar T) const
76 {
77  return rho0_ + p/(this->R()*T);
78 }
79 
80 
81 template<class Specie>
82 inline Foam::scalar Foam::perfectFluid<Specie>::h(scalar p, scalar T) const
83 {
84  return p/rho(p, T) - Pstd/rho(Pstd, T);
85 }
86 
87 
88 template<class Specie>
89 inline Foam::scalar Foam::perfectFluid<Specie>::Cp(scalar p, scalar T) const
90 {
91  return
92  sqr(p/(rho(p, T)*T))/this->R()
93  - sqr(Pstd/(rho(Pstd, T)*T))/this->R();
94 }
95 
96 
97 template<class Specie>
98 inline Foam::scalar Foam::perfectFluid<Specie>::e(scalar p, scalar T) const
99 {
100  return 0;
101 }
102 
103 
104 template<class Specie>
105 inline Foam::scalar Foam::perfectFluid<Specie>::Cv(scalar p, scalar T) const
106 {
107  return 0;
108 }
109 
110 
111 template<class Specie>
112 inline Foam::scalar Foam::perfectFluid<Specie>::sp(scalar p, scalar T) const
113 {
114  return -this->R()*log(p/Pstd);
115 }
116 
117 
118 template<class Specie>
119 inline Foam::scalar Foam::perfectFluid<Specie>::sv(scalar p, scalar T) const
120 {
122  return 0;
123 }
124 
125 
126 template<class Specie>
127 inline Foam::scalar Foam::perfectFluid<Specie>::psi(scalar p, scalar T) const
128 {
129  return 1.0/(this->R()*T);
130 }
131 
132 
133 template<class Specie>
134 inline Foam::scalar Foam::perfectFluid<Specie>::Z(scalar p, scalar T) const
135 {
136  return p/(rho(p, T)*this->R()*T);
137 }
138 
139 
140 template<class Specie>
141 inline Foam::scalar Foam::perfectFluid<Specie>::CpMCv(scalar p, scalar T) const
142 {
143  const scalar R = this->R();
144  const scalar rho = this->rho(p, T);
145 
146  return R*sqr(p/(rho*R*T));
147 }
148 
149 
150 template<class Specie>
151 inline Foam::scalar Foam::perfectFluid<Specie>::alphav(scalar p, scalar T) const
152 {
153  return (1 - rho0_/this->rho(p, T))/T;
154 }
155 
156 
157 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
158 
159 template<class Specie>
161 (
162  const perfectFluid<Specie>& pf
163 )
164 {
166 }
167 
168 
169 template<class Specie>
170 inline void Foam::perfectFluid<Specie>::operator*=(const scalar s)
171 {
172  Specie::operator*=(s);
173 }
174 
175 
176 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
177 
178 template<class Specie>
179 inline Foam::perfectFluid<Specie> Foam::operator+
180 (
181  const perfectFluid<Specie>& pf1,
182  const perfectFluid<Specie>& pf2
183 )
184 {
186  return pf1;
187 }
188 
189 
190 template<class Specie>
191 inline Foam::perfectFluid<Specie> Foam::operator*
192 (
193  const scalar s,
194  const perfectFluid<Specie>& pf
195 )
196 {
197  return perfectFluid<Specie>
198  (
199  s*static_cast<const Specie&>(pf),
200  pf.R_,
201  pf.rho0_
202  );
203 }
204 
205 
206 template<class Specie>
207 inline Foam::perfectFluid<Specie> Foam::operator==
208 (
209  const perfectFluid<Specie>& pf1,
210  const perfectFluid<Specie>& pf2
211 )
212 {
213  noCoefficientMixing(perfectFluid);
214  return pf1;
215 }
216 
217 
218 // ************************************************************************* //
#define R(A, B, C, D, E, F, K, M)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Simple extension of the perfect gas equation of state to liquids by the addition of a constant densit...
Definition: perfectFluid.H:122
scalar Cv(scalar p, scalar T) const
Return Cv contribution [J/(kg K].
scalar psi(scalar p, scalar T) const
Return compressibility [s^2/m^2].
perfectFluid(const Specie &sp, const scalar R, const scalar rho0)
Construct from components.
Definition: perfectFluidI.H:32
scalar alphav(const scalar p, const scalar T) const
Return volumetric coefficient of thermal expansion [1/T].
scalar e(const scalar p, const scalar T) const
Return internal energy contribution [J/kg].
Definition: perfectFluidI.H:98
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
Definition: perfectFluidI.H:75
autoPtr< perfectFluid > clone() const
Construct and return a clone.
Definition: perfectFluidI.H:59
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
scalar h(const scalar p, const scalar T) const
Return enthalpy contribution [J/kg].
Definition: perfectFluidI.H:82
scalar Cp(scalar p, scalar T) const
Return Cp contribution [J/(kg K].
Definition: perfectFluidI.H:89
scalar sv(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cv/T [J/kg/K].
scalar sp(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cp/T [J/kg/K].
scalar Z(scalar p, scalar T) const
Return compression factor [].
scalar R() const
Return fluid constant [J/kg/K].
Definition: perfectFluidI.H:68
void operator*=(const scalar)
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const dimensionedScalar Pstd
Standard pressure.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensionedScalar log(const dimensionedScalar &ds)
static scalar R(const scalar a, const scalar x)
Definition: invIncGamma.C:102
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar rho0
volScalarField & p
#define noCoefficientMixing(Type)
Definition: specie.H:159